第3章 MATLAB應用

1、School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用第三章第三章 MatlabMatlab應用應用 School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用 3.1 解析解與數(shù)值解 解析解符號運算 數(shù)值解近似解 有許多問題沒有解析解，有些問題有解析解但不實用。 許多問題的解決依靠的是數(shù)值解法，工程上更為關心問題的數(shù)值解-非線性方程的根、非線性系統(tǒng)的響應。 仿真本身就是尋找系統(tǒng)的數(shù)值解。Sc

2、hool of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用3.2 線性代數(shù) AX=B X=A-1B 矩陣的特征參數(shù) 秩rank(A) 行列式det(A)-Determinant 特征根eig(A)-Eigenvalue 范數(shù)norm(A) 跡trace(A) 對角元素和School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用3.3 微積分問題q數(shù)值積分數(shù)值積分y, n=quad(F, a, b, t

3、ol)q 參數(shù)參數(shù) F 可通過以下方式給出：可通過以下方式給出：l quad(sin(x),0,pi) 2l F=inline(sin(x); quad(F,0,2*pi) 0l myfun.mfunction y=myfun(x)y=sin(x);l quad(myfun, 0,2*pi)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用3.5 非線性方程q解非線性方程解非線性方程x = fsolve(F,X0);q 例：例：fun = inline(sin(x);x = fsolve

4、(fun,1 4 10, optimset(fsolve)x =0.0000 3.1416 9.4248School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用3.6 數(shù)據(jù)插值與統(tǒng)計分析數(shù)據(jù)插值與統(tǒng)計分析q插值插值yi=interp1(x, y, xi, 方法方法) 方法方法=linear, cubic, splineq 例：例：x=-3:3;y=-1;-1;-1;0;1;1;1;xi=-3:0.1:3;yi1=interp1(x,y,xi,linear)yi2=interp1(x,y,xi,

5、cubic)yi3=interp1(x,y,xi,spline)plot(x,y,o,xi,yi1,-*, xi,yi2,-.x, xi,yi3,-b)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用q擬合-多項式p=polyfit(x,y,3) yy=polyval(p,t);q 例：例：x=-3:3;y=-1 -1 -1 0

6、1 1 1;p=polyfit(x,y,3)t=-3:.1:3;yy=polyval(p,t);plot(x,y,o,t,yy)p=-0.0556 0.0000 0.8175 -0.0000Polyval:多項式求值Y = P(1)*XN + P(2)*X(N-1) + . + P(N)*X + P(N+1)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用擬合-多項式School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機

7、仿真第三章第三章 Matlab應用應用3.6 數(shù)據(jù)插值與統(tǒng)計分析數(shù)據(jù)插值與統(tǒng)計分析q統(tǒng)計統(tǒng)計x=rand(n,m)x=randn(n,m) N(,2) yi= +ximean()std()q 例例x=randn(10000, 1)；x=2+5*x；x_mean=mean(x)：2.0053x_std=std(x)：5.0057School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用 y=randn(30000,1) x=linspace(-3,3,61)%步長0.1 yn=hist(y,x);

8、%y落在以x為中心的桶內的個數(shù)。 yp=yn/(30000*0.1);%歸一化 bar(x,yp)Linspace(a,b,N) a，b內等間距N個點。Linspace(1,2,3) 1 1.5 2School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用bar(1 2 3,4 5 6)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用3.7 符號數(shù)學工具箱q MatlabMatlab 符

9、號運算是通過符號運算是通過符號數(shù)學工具箱符號數(shù)學工具箱（Symbolic Symbolic Math ToolboxMath Toolbox）來實現(xiàn)的。來實現(xiàn)的。MatlabMatlab 符號數(shù)學工具箱是符號數(shù)學工具箱是建立在功能強大的建立在功能強大的 Maple Maple 軟件的基礎上的，當軟件的基礎上的，當 MatlabMatlab 進行符號運算時，它就請求進行符號運算時，它就請求 Maple Maple 軟件去計算并將結軟件去計算并將結果返回給果返回給 MatlabMatlab。 q符號表達式的運算，符號表達式的復合、化簡，符號矩陣的運算，符號微積分、符號作圖，符號代數(shù)方程求解，符號微

10、分方程求解等。School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用q 符號表達式的建立：符號表達式的建立：syms x, y, ty1=sin(x)+cos(x)y2=exp(t)符號表達式的建立School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用微積分微積分求極限Limit()syms x ; limit(sin(x)/x)-1求積分int()syms x; int(1/(1+x2)

11、 -atan(x)求泰勒級數(shù)taylor()syms x; taylor(exp(-x),6)-1-x+1/2*x2-1/6*x3+1/24*x4-1/120*x5求導數(shù)diff()syms x; diff(sin(x)-cos(x)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用線性代數(shù)線性代數(shù)求逆矩陣Inv()syms t; a= t 2*t; 2*t, t b=inv(a) b=-1/3/t, 2/3/t; 2/3/t, -1/3/t求行列式Det()syms t; a= t 2

12、*t; 2*t, t det(a)：-3*t2School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用解方程解方程解方程Solve(eqn1,eqn2,.,eqnN,var1,var2,.,varN)solve(p*sin(x) = r) ：ans =asin(r/p) x,y = solve(x2 + x*y + y = 3,x2 - 4*x + 3 = 0) x = 1， 3 y = 1， -3/2 解微分方程Dsolve(eq1,eq2, . ,ini1,ini2, . ,v)dsolv

13、e(Dx = -a*x) ：ans = exp(-a*t)*C1 解微分方程Dsolve()x = dsolve(Dx = -a*x,x(0) = 1)：x = exp(-a*t)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用變換變換反拉普拉斯變換Ilaplace()syms s; ilaplace(1/(s-1) ：exp(t)拉普拉斯變換Laplace()syms t; laplace(t+sin(t)：1/s2+1/(s2+1)School of Arming Engineer

14、ing沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用3.8控制系統(tǒng)工具箱控制系統(tǒng)工具箱-線性系統(tǒng)線性系統(tǒng) 傳遞函數(shù)G=tf(num, den)零極點輸入G=zpk(z, p, k) 串聯(lián)G=G1*G2 并聯(lián)G=G1+G2 反饋G=feedback(G1, G2)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用 階躍step(G) 沖擊impulse(G) 波特圖bode 幅值和相位裕度margin 根軌跡rlocus School

15、 of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用例題：已知單位反饋系統(tǒng)的開環(huán)傳遞函數(shù)為 （K=10, T=2），求 系統(tǒng)的閉環(huán)傳遞函數(shù)。 系統(tǒng)的單位階躍響應。 繪制系統(tǒng)的波特圖。 繪制系統(tǒng)關于K的根軌跡。) 1)(1()(sTssKsGSchool of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用rlocus(g) g=zpk(,0 -.5 -1,5) gb=feedback(g,1) step(gb,

16、 15) figure margin(g) figure rlocus(g)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用step(gb, 15)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用margin(g)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用rloc

17、us(g)School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用3.4 微分方程的數(shù)值解法T,Y = solver(odefun,tspan,y0)其中其中 y0 為初值條件，為初值條件，tspan為求解區(qū)間；為求解區(qū)間；Matlab在數(shù)值求解在數(shù)值求解時時自動對求解區(qū)間進行分割自動對求解區(qū)間進行分割，T (向量向量) 中返回的是分割點的中返回的是分割點的值值(自變量自變量)，Y (向量向量) 中返回的是解函數(shù)在這些分割點上的中返回的是解函數(shù)在這些分割點上的函數(shù)值。函數(shù)值。solver 為為

18、Matlab的的ODE求解器求解器（可以是（可以是 ode45、ode23、ode113、ode15s、ode23s、ode23t、ode23tb）沒有一種算法可以有效地解決所有的沒有一種算法可以有效地解決所有的 ODE 問題，因此問題，因此MATLAB 提供了多種提供了多種ODE求解器，對于不同的求解器，對于不同的ODE，可以，可以調用不同的調用不同的求解器求解器。School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用Matlab提供的ODE求解器求解器ODE類型特點說明ode45非剛性單

19、步法；4階5級 R-K 方法；累計截斷誤差為 (x)3大部分場合的首選方法ode23非剛性單步法； 2階3級 R-K 方法；累計截斷誤差為 (x)3使用于精度較低的情形ode113非剛性多步法；Adams算法；高低精度均可到 10-310-6計算時間比 ode45 短ode23t適度剛性 采用梯形算法適度剛性情形ode15s剛性多步法；Gears 反向數(shù)值微分；精度中等若 ode45 失效時，可嘗試使用ode23s剛性單步法；2 階Rosebrock 算法；低精度當精度較低時，計算時間比 ode15s 短ode23tb剛性梯形算法；低精度當精度較低時，計算時間比ode15s短School of

20、 Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用微分方程的數(shù)值求解微分方程的數(shù)值求解微分方程t, x=ode45(函數(shù)名, tspan, x0) 函數(shù)文件lorenzeq.mfunction xdot = lorenzeq(t, x)xdot=-8/3*x(1)+x(2)*x(3); -10*x(2)+10*x(3); -x(1)*x(2)+28*x(2)-x(3);腳本文件mymain1.mt,x=ode45(lorenzeq, 0, 10, 0, 0, 1e-10)plot(t,x) 例：洛侖茲方程

21、School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用彈道仿真實例考慮魚雷在水平面內追蹤目標的情況?？紤]魚雷在水平面內追蹤目標的情況。我們需要建立：我們需要建立：目標的運動方程目標的運動方程魚雷的運動方程魚雷的運動方程根據(jù)目標和魚雷的運動狀態(tài)，按照設定的導引根據(jù)目標和魚雷的運動狀態(tài)，按照設定的導引方法，控制舵機使魚雷追蹤目標。方法，控制

22、舵機使魚雷追蹤目標。School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用坐標系目標模型 魚雷模型 導引方法School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用地面坐標系地面坐標系OX0Y0Z0 坐標原點選在地面適當?shù)奈恢谩?X0軸處于地平面內指向某一適當方向；Y0軸垂直于地面并指向上方，此軸稱為鉛垂軸； Z0軸處于地平面內，垂直于OX0Y0平面。雷體坐標系雷體坐標系 ox1y1z1 坐

23、標原點O取在魚雷重心處。 三個坐標軸與雷體固連。 OX1軸沿雷體幾何對稱軸指向雷頭。OY1軸處于魚雷縱對稱面內。 速度坐標系速度坐標系OXYZ 原點O取在魚雷重心處，OX軸為重心運動瞬時速度矢量v。一般情況下，v不在魚雷對稱平面內。OY軸在雷體縱對稱面內，且垂直于OX軸，指向上方，稱為升力軸。School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matla

24、b應用應用School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用航向角：雷體縱軸ox1在水平面的投影與地面坐標系ox0之間的夾角。側滑角：雷體縱軸ox1與魚雷的速度方向ox在水平面內的夾角。彈道偏角：速度方向ox與地面坐標系ox0之間的夾角 。返回School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用目標模型在地面坐標系內：目標運動方程組如下：TTTTTTTVZVXsincos式中，T、T

25、、XT、ZT分別為目標彈道偏角、回旋角速度、縱向距離和側向距離。返回School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用魚雷側向運動的方程組-魚雷側滑角魚雷側滑角y-回旋回旋(偏航偏航)角速度角速度r-偏航偏航(垂直垂直)舵角舵角-航向角航向角-彈道偏角彈道偏角Vm-速度速度Xe、Ze-地面系中坐標。地面系中坐標。sincos84.11915425.9261. 319.8 . 5memeyryyryVZVXSchool of Arming Engineering沈陽理工大學裝備工程學院沈陽

26、理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用魚雷與目標的相對運動關系魚雷與目標的相對運動關系 School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用 x0為基準線；T表示魚雷的位置；t表示目標的位置，Tt是瞄準線，或稱視線；r表示魚雷與目標的相對距離；q是瞄準線和基準線之間的夾角，稱為舷角；v和vt分別是魚雷和目標的速度；和t分別稱為提前角或前置角；和t分別是魚雷和目標的彈道偏角彈道偏角，將魚雷和目標視為質點時，也是魚雷和目標的航向角。School of Arm

27、ing Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用導引方法導引方法 School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用導引方法 魚雷與目標相對距離為 X = XT-Xe，Z=ZT-Ze，視線角定義為 q=atan(-Z/X) q為地球視線角，或=q-，為雷體系中的提前角。魚雷尾追法的操舵規(guī)律為r=-k，k-比例系數(shù)，k=0.5, |r |=70s時，T=-0.1rad/s，目標開始做回旋機動。School

28、of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用定義函數(shù)torpedo2ship.m function dx=torpedo2ship(t,x) delr=0.5*(x(3)-atan(x(5)-x(8)/(x(7)-x(4); if abs(delr)10*0.0175 delr=sign(delr)*0.175; end if t=70 omigat=0; %target move in direct lineelse else omigat=-.1; %target maneuver endS

29、chool of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用續(xù) dx=-5.8*x(1)+0.19*x(2)-3.61*delr; 92.425*x(1)-15*x(2)-119.84*delr; x(2); 25*cos(x(3)-x(1); -25*sin(x(3)-x(1); omigat; 5*cos(x(6); -5*sin(x(6); School of Arming Engineering沈陽理工大學裝備工程學院沈陽理工大學裝備工程學院計算機仿真計算機仿真第三章第三章 Matlab應用應用trajectory.mx0=0 0 0 0 0 0.4 1500 1500; t,x=ode45(torpedo2ship, 0,97, x0); plo