Manovellismo in Matlab
Salve a tutti,
ho un problema con Matlab: vorrei riuscire a rappresentare il meccanismo di funzionamento di un pistone sia in termini numerici che
con parte grafica. Quello che per ora sono riuscita a fare riguarda un semplice meccanismo di manovellismo ordinario centrato con sistema biella- manovella. Qualcuno riesce a darmi una mano per creare lo stesso modello per un pistone?
Grazie!!
Posto qui il manovellismo prodotto:
ho un problema con Matlab: vorrei riuscire a rappresentare il meccanismo di funzionamento di un pistone sia in termini numerici che
con parte grafica. Quello che per ora sono riuscita a fare riguarda un semplice meccanismo di manovellismo ordinario centrato con sistema biella- manovella. Qualcuno riesce a darmi una mano per creare lo stesso modello per un pistone?
Grazie!!
Posto qui il manovellismo prodotto:
Sfunction mocen % **************************************************** % programma per la simulazione della cinematica del % manovellismo ordinario centrato % **************************************************** % Inserimento dati stato = 1; while stato == 1 a = input('Assegna raggio manovella in [m]: '); b = input('Assegna lunghezza biella in [m]: '); if a < b stato = 0; else disp('La manovella deve essere più corta della biella') end end alfa = [1:0.1:360].*(pi/180); alfap = input('Assegna velocita'' angolare manovella [rad/s]: '); alfapp = 0; % Analisi di posizione disp('Risoluzione analisi di posizione') c = a*cos(alfa) + b*sqrt(1-((a*sin(alfa))./b).^2); beta = asin(-(a*sin(alfa))./b); %% PRIMA E SECONDA APPROSSIMAZIONE c1_approx=a*cos(alfa) + b; %Buona approx solo se lambda è molto minore di 1 c2_approx=a*cos(alfa) + b + a*a/b/4*(cos(2*alfa)-1) ; %% PIEDE DI BIELLA xb=a+b-c; xb1_approx=a+b-c1_approx; xb2_approx=a+b-c2_approx; % Analisi di velocità disp('Risoluzione analisi di velocità') for i = 1:length(alfa) A = [1 b*sin(beta(i)); 0 -b*cos(beta(i))]; B = [-alfap*a*sin(alfa(i)); alfap*a*cos(alfa(i))]; xp = inv(A)*B; cp(i) = xp(1); betap(i) = xp(2); end % cp = -alfap*a*(sin(alfa) - cos(alfa).*tan(beta)); % betap = -alfap*((a*cos(alfa))./(b*cos(beta))); %% PRIMA E SECONDA APPROSSIMAZIONE cp1_approx=a*alfap * ( -sin(alfa) ); cp2_approx=a*alfap * ( -sin(alfa) -a/b/2*sin(2*alfa) ); % Analisi di accelerazione disp('Risoluzione analisi di accelerazione') for i = 1:length(alfa) A = [1 b*sin(beta(i)); 0 -b*cos(beta(i))]; B = -[ alfapp*a*sin(alfa(i)) + alfap^2*a*cos(alfa(i)) + betap(i)^2*b*cos(beta(i)); -alfapp*a*cos(alfa(i)) + alfap^2*a*sin(alfa(i)) + betap(i)^2*b*sin(beta(i))]; xpp = inv(A)*B; cpp(i) = xpp(1); betapp(i) = xpp(2); end %% PRIMA E SECONDA APPROSSIMAZIONE cpp1_approx=a* alfap^2 * ( -cos(alfa) ); cpp2_approx=a* alfap^2 * ( -cos(alfa) -a/b*cos(2*alfa) ); % Visualizzazione risultati stato = 1; while stato == 1 k = menu('Visualizzazione risultati', ... 'Grafici posizione', ... 'Grafici velocita''', ... 'Grafici accelerazione', ... 'Animazione del moto', ... 'C.I.R.',... 'Fine'); if k == 1 figure(1) subplot(211) %plot(alfa,c) plot(alfa,c,alfa,c1_approx,alfa,c2_approx) xlabel('Rotazione manovella [rad]') ylabel(sprintf('Posizione piede di biella [m]')) legend('Esatta','Prima approx','Seconda approx') ax = axis; axis([0 2*pi ax(3) ax(4)]) grid on subplot(212) plot(alfa,beta) xlabel('Rotazione manovella [rad]') ylabel(sprintf('Rotazione biella [rad]\n')) ax = axis; axis([0 2*pi ax(3) ax(4)]) grid on elseif k == 2 figure(2) subplot(211) %plot(alfa,cp) plot(alfa,cp,alfa,cp1_approx,alfa,cp2_approx) xlabel('Rotazione manovella [rad]') ylabel(sprintf('Velocita'' piede di biella [m/s]')) legend('Esatta','Prima approx','Seconda approx') ax = axis; axis([0 2*pi ax(3) ax(4)]) grid on subplot(212) plot(alfa,betap) xlabel('Rotazione manovella [rad]') ylabel(sprintf('Velocita'' angolare biella [rad/s]\n')) ax = axis; axis([0 2*pi ax(3) ax(4)]) grid on figure(20) plot(a-xb,-cp./alfap,(a-xb1_approx),-cp1_approx./alfap,(a-xb2_approx),-cp2_approx./alfap) xlabel('r-x_B') ylabel('v_B/\omega') legend('Esatta','Prima approx','Seconda approx') grid axis equal text(a*1.1,0,'PMS') text(-a*1.25,0,'PMI') elseif k == 3 figure(3) subplot(211) %plot(alfa,cpp) plot(alfa,cpp,alfa,cpp1_approx,alfa,cpp2_approx) xlabel('Rotazione manovella [rad]') ylabel(sprintf('Accelerazione piede di biella [m/s^2]')) legend('Esatta','Prima approx','Seconda approx') ax = axis; axis([0 2*pi ax(3) ax(4)]) grid on subplot(212) plot(alfa,betapp) xlabel('Rotazione manovella [rad]') ylabel(sprintf('Accelerazione angolare biella [rad/s^2]\n')) ax = axis; axis([0 2*pi ax(3) ax(4)]) grid on figure(20) plot((a-xb),-cpp./alfap^2,(a-xb1_approx),-cpp1_approx./alfap^2,(a-xb2_approx),-cpp2_approx./alfap^2) xlabel('r-x') ylabel('a_B/\omega^2') legend('Esatta','Prima approx','Seconda approx') grid axis equal text(a*1.1,0,'PMS') text(-a*1.25,0,'PMI') elseif k == 4 figure(4) set(gcf,'DoubleBuffer','on') ngiri = 3; nxframe = 5; nframe = ngiri*360/nxframe; for i = 1:nframe igrado = nxframe*i; n = floor(igrado/360); igrado = igrado - 360*n; ic = find(alfa*180/pi >= igrado); x_manov = a*cos(alfa(ic(1))); y_manov = a*sin(alfa(ic(1))); x_biell = c(ic(1)); y_biell = 0; plot([0 x_manov],[0 y_manov],'r', ... [x_manov x_biell],[y_manov y_biell],'b',... [0 x_manov x_biell],[0 y_manov y_biell],'ok') axis([-a-b a+b -a-b a+b]) title(sprintf('Ciclo: %d',n+1)) xlabel('[m]') ylabel('[m]') grid on drawnow end elseif k==5 angolo=10; ic = find(angolo/180*pi==alfa); x_manov = a*cos(alfa(ic)); y_manov = a*sin(alfa(ic)); x_biell = c(ic(1)); y_biell = 0; figure(10) plot([0 x_manov],[0 y_manov],'r', ... [x_manov x_biell],[y_manov y_biell],'b',... [0 x_manov x_biell],[0 y_manov y_biell],'ok') axis([-a a+b -a-b/2 a+b/2]),grid on hold on %plot CIR plot([x_manov x_biell],[y_manov y_manov/x_manov*x_biell],'k:') plot([x_biell x_biell], [0 y_manov/x_manov*x_biell],'k:') plot(x_biell,y_manov/x_manov*x_biell,'ko') text(x_biell*0.9,y_manov/x_manov*x_biell*1.1,'CIR') plot([0 0],[0 y_manov/x_manov*x_biell],'k:') plot([x_manov 0],[y_manov -tan(beta(ic))*x_biell],'k:') plot(0,-tan(beta(ic))*x_biell,'ro') text(-0.3,-tan(beta(ic))*x_biell*1.1,'P') text(-0.3,0.,'O') text(x_manov,y_manov*1.8,'A') text(x_biell*0.95,-0.2,'B') if abs(a+b/2)<y_manov/x_manov*x_biell axis([-a a+b 0 y_manov/x_manov*x_biell]) end figure(11) subplot(2,1,1) plot(alfa,cp), hold on plot(alfa(ic),cp(ic),'ko'), grid xlabel('Rotazione manovella [rad]') ylabel(sprintf('v_B [m/s]')) subplot(2,1,2) plot(alfa,-alfap*a), hold on plot(alfa(ic),-alfap*a,'ko'), grid xlabel('Rotazione manovella [rad]') ylabel(sprintf('v_A [m/s]')) ax = axis; axis([0 2*pi ax(3) ax(4)]) grid on % keyboard else stato = 0; close all end end $