Z Skrypty dla studentów Ekonofizyki UPGOW
(Utworzył nową stronę „== Modelowanie: Armata == Mamy: :<math>\vec F = m \vec a</math> więc: :<math>\vec F = m \vec \ddot x</math> w kartezjańskim układzie współrzędnych mamy: :<math>\b…”) |
(→Modelowanie: Armata) |
||
| (Nie pokazano 1 wersji pomiędzy niniejszymi.) | |||
| Linia 24: | Linia 24: | ||
:<math>\begin{cases} \dot x=v_x \\ \dot y=v_y\\ \dot v_x = -\gamma/m \; (v_x-u(t)) \\ \dot v_y=-\gamma/m\; v_y-g \end{cases}</math> | :<math>\begin{cases} \dot x=v_x \\ \dot y=v_y\\ \dot v_x = -\gamma/m \; (v_x-u(t)) \\ \dot v_y=-\gamma/m\; v_y-g \end{cases}</math> | ||
| + | |||
| + | == Kod źródłowy == | ||
| + | <source lang="matlab"> | ||
| + | clear all | ||
| + | close all | ||
| + | disp('This is a gun!') | ||
| + | |||
| + | function dx = gun(X,T) | ||
| + | g=10; m=1; gama=1.5; | ||
| + | u=(T>0.1)*10; | ||
| + | dx = zeros(4,1); | ||
| + | dx(1) = X(3); % x | ||
| + | dx(2) = X(4); % y | ||
| + | dx(3) = -gama/m*( X(3)-u ); % vx | ||
| + | dx(4) = -gama/m*X(4) - g; % vy | ||
| + | return | ||
| + | end | ||
| + | |||
| + | n=125; | ||
| + | t = linspace(0,1.22,n); | ||
| + | i=1; | ||
| + | subplot(2,1,1) | ||
| + | hold on | ||
| + | alfas=linspace(0.01,pi/2,10); | ||
| + | v0=1.; | ||
| + | for alfa=alfas | ||
| + | X0 = [0,0,v0*cos(alfa),v0*sin(alfa)]; | ||
| + | [X,T,MSG]=lsode(@gun,X0,t); | ||
| + | plot(X(:,1),X(:,2),'r-') | ||
| + | i0=find(X(:,2)<0)(1)-1; | ||
| + | i1=find(X(:,2)<0)(1); | ||
| + | x_interp(i)=X(i0,1)+( 0 -X(i0,2))/(X(i1,2)-X(i0,2))*(X(i1,1)-X(i0,1)); | ||
| + | plot(x_interp(i), 0 ,'g*') | ||
| + | i+=1; | ||
| + | endfor | ||
| + | |||
| + | hold off | ||
| + | grid on | ||
| + | ylim([0,.05]) | ||
| + | xlim([0,.1]) | ||
| + | subplot(2,1,2) | ||
| + | plot(alfas*180/pi,x_interp,'bo') | ||
| + | alfa_fine=linspace(0,90,1000); | ||
| + | x_fine=interp1(alfas*180/pi,x_interp,alfa_fine,'cubic'); | ||
| + | [zasieg,alfa_i]=max(x_fine); | ||
| + | hold on | ||
| + | plot(alfa_fine,x_fine) | ||
| + | plot(alfa_fine(alfa_i),zasieg,'ro') | ||
| + | printf("Najwiekszy zasieg dla strzalu pod katem %f stopni\n",alfa_fine(alfa_i)) | ||
| + | </source> | ||
Aktualna wersja na dzień 21:34, 4 sty 2010
Modelowanie: Armata
Mamy: \[\vec F = m \vec a\] więc: \[\vec F = m \vec \ddot x\] w kartezjańskim układzie współrzędnych mamy: \[\begin{cases} F_x = m \ddot x \\ F_y = m \ddot y \end{cases}\] Siły w ogólności zależą od prędkości, czasu i położenia. W przypadku lotu pocisku możemy założyć, że działa na niego siła tarcia oraz siła ciężkości.
\[\vec F = \vec T + \vec Q \] Gdzie Q to ciężar \[\vec Q = \begin{cases} 0 \\ -mg \end{cases}\]
Siła tarcia zależy od prędkości ruchu posicku względem powietrza. Jeżeli będzie wiał wiatr \[\vec u(t) = \begin{cases} u(t) \\0 \end{cases}\] to trzeba jego prękość odjąć od prędkości ruchu pocisku względem ziemi (wiatr pozorny). \[\vec T = - \gamma (\vec v-\vec{u(t)}) = \begin{cases} -\gamma (v_x-u(t) ) \\ -\gamma v_y \end{cases} \]
Czyli mamy układ równań:
\[\begin{cases} \dot x=v_x \\ \dot y=v_y\\ \dot v_x = -\gamma/m \; (v_x-u(t)) \\ \dot v_y=-\gamma/m\; v_y-g \end{cases}\]
Kod źródłowy
clear all close all disp('This is a gun!') function dx = gun(X,T) g=10; m=1; gama=1.5; u=(T>0.1)*10; dx = zeros(4,1); dx(1) = X(3); % x dx(2) = X(4); % y dx(3) = -gama/m*( X(3)-u ); % vx dx(4) = -gama/m*X(4) - g; % vy return end n=125; t = linspace(0,1.22,n); i=1; subplot(2,1,1) hold on alfas=linspace(0.01,pi/2,10); v0=1.; for alfa=alfas X0 = [0,0,v0*cos(alfa),v0*sin(alfa)]; [X,T,MSG]=lsode(@gun,X0,t); plot(X(:,1),X(:,2),'r-') i0=find(X(:,2)<0)(1)-1; i1=find(X(:,2)<0)(1); x_interp(i)=X(i0,1)+( 0 -X(i0,2))/(X(i1,2)-X(i0,2))*(X(i1,1)-X(i0,1)); plot(x_interp(i), 0 ,'g*') i+=1; endfor hold off grid on ylim([0,.05]) xlim([0,.1]) subplot(2,1,2) plot(alfas*180/pi,x_interp,'bo') alfa_fine=linspace(0,90,1000); x_fine=interp1(alfas*180/pi,x_interp,alfa_fine,'cubic'); [zasieg,alfa_i]=max(x_fine); hold on plot(alfa_fine,x_fine) plot(alfa_fine(alfa_i),zasieg,'ro') printf("Najwiekszy zasieg dla strzalu pod katem %f stopni\n",alfa_fine(alfa_i))