積分的數(shù)值求解函數(shù)

1、積分的數(shù)值求解函數(shù),積分: (求面積) trapz(x , y)給定數(shù)據(jù)點(diǎn)x和y，計(jì)算y=f(x)下的梯形面積積分。精度更高的求積分函數(shù) quad,QUADL;近似求原函數(shù)并做圖mm2302二維積分:Dblquad三維積分: Triplequad,微分: 具有其內(nèi)在的敏感性,故對于被積分的數(shù)據(jù)使通過實(shí)驗(yàn)獲得時(shí),避免使用.應(yīng)該相對數(shù)據(jù)做( ),然后對( )進(jìn)行微分. Polyfit Polyder diff

2、(x)數(shù)組元素間的差值(diff(y)./diff(x)近似求微分)，結(jié)果分析mm2305,mm2306 結(jié)論二維微分:梯度gradient 主要用作圖形數(shù)據(jù)可視化(mm2307)quiver曲率del2(mm2308),微分的數(shù)值求解函數(shù),微分方程(ODE)（1）,1.IVP格式2.IVP微分方程解法程序[t , y]=ode45(‘fname’ , [tspan] , yo)用4階/5階龍格-庫塔算法解微分方程組 ot

3、her 6 ODE solvers: ODE23, ODE113, ODE15S, ODE23S, ODE23T, ODE23TB 設(shè)置選項(xiàng):OPTIONS=ODESET (‘name1’,value1,…), ODEGET；,微分方程 (ODE)（2）,output functions: ODEPLOT,ODEPHAS2,ODEPHAS3, ODEPRINT；evaluating solution:

4、 DEVAL； ODEexamples:RIGIDODE,BALLODE, ORBITODE 邊界值問題(BVP)和偏微分方程問題,基本用法（連續(xù)時(shí)間系統(tǒng)仿真1）,1、微分方程變形：ydot=ldefile(t,y)m文件的編寫function ydot=vdpol(t,y,mu)%VDPOL van der Pol equation.% Ydot=VDPOL(t,Y)% Ydot(1) = Y(

5、2)% Ydot(2) = mu*(1-Y(1)^2)*Y(2)-Y(1)% mu = 2%ydot = [y(2,:); mu*(1-y(1,:)^2)*y(2,:)-y(1,:)];ydot = [y(2); mu*(1-y(1)^2)*y(2)-y(1)];2、求解命令的編寫2401.doc,基本用法（連續(xù)時(shí)間系統(tǒng)仿真2）,例 :微分方程的數(shù)值解法是在科學(xué)與工程計(jì)算中經(jīng)常遇到的問題。著名的 Lorenz 模型的狀態(tài)方程表

6、示為： 若令 ，且假設(shè)初值為 ，e 為一個(gè)小常數(shù)，洛倫茲.doc,離散系統(tǒng)仿真,FilterCon,SIMULINK動態(tài)仿真,SIMULINK動態(tài)仿真是MATLAB軟件包之一，用于可視化的動態(tài)系統(tǒng)仿真，適用于連續(xù)、離散系統(tǒng)，也適用于線性和非線性系統(tǒng)。 1）啟動SIMULINK，進(jìn)入

7、SIMULINK窗口； 2）在SIMULINK窗口下，借助模塊庫，創(chuàng)建系統(tǒng) 框圖模型并調(diào)整模塊參數(shù)。 3）兩種方式（ SIMULINK模型窗口、matlab command 窗口sim）下，設(shè)置仿真參數(shù)后， 啟動仿真。 4）輸出仿真結(jié)果：利用scope模塊，利用變量返 回值、利用to workspace模塊。,SIMULINK動態(tài)仿真(2),模塊平衡點(diǎn)（trim)和線性化(linmod、lin

8、mod2、dlinmod等）子系統(tǒng)創(chuàng)建及封裝 S-函數(shù)非線性系統(tǒng)優(yōu)化設(shè)計(jì),作業(yè)2:系統(tǒng)時(shí)間響應(yīng)和動態(tài)仿真,應(yīng)用4階龍格-庫塔法和使用simulink軟件包仿真程序，求如下系統(tǒng)的響應(yīng)。,作業(yè)2 分析求解(1),運(yùn)用四階龍格-塔庫公式編程實(shí)現(xiàn) 公式：x(n+1)=x(n)+h/6(k1+2*k2+2*k3+k4) k1=f[x(n),t(n)] k2

9、=f[x(n)+h/2*k1,t(n)+h/2] k3=f[x(n)+h/2*k2,t(n)+h/2] k4=f[x(n)+h*k3,t(n)+h],編程實(shí)現(xiàn),1.定義函數(shù)dydtfunction ydot=dydt(t,y)%VDPOL van der Pol equation.% Ydot=VDPOL(t,Y)% Ydot(1) = Y(2)% Y

10、dot(2) = Y(3)% Ydot(3) = Y(4)% Ydot(4) = Y(5) % Ydot(5) = Y(6)% Ydot(6) = -a0*Y(1)-a1*Y(2)-a2*Y(3)-a3*Y(4)-a4*Y(5) –a5*Y(6)+1% ydot =[Y(2);Y(3);Y(4);Y(5);Y(6); … -a0*Y(1)-a1*Y(2)-a2*Y(3)-a3*Y(4)-a4*Y

11、(5) –a5*Y(6)+1],編程實(shí)現(xiàn),ydot =[y(2);y(3);y(4);y(5);y(6); … -186.5*y(1)-840*y(2)-904.4*y(3)-237.3*y(4)-76.1*y(5)-8.8*y(6)+1]2.在workspace鍵入 »tspan=[0:0.1:5]; » y0=[0;0;0;0;0;0]; »[t,y

12、]=ode45('dydt',tspan,y0); » z= -0.2*y(:,1)+1187.6*y(:,2)+3699.6*y(:,3)+327*y(:,4)+65*y(:,5); » plot(z) 或者遍寫M文件 參見 work2.m,作業(yè)2 分析求解(2),用simulink系統(tǒng)仿真求解 （1）常微分方程模型

13、（時(shí)域分析） （2）傳遞函數(shù)模型 （3）狀態(tài)方程模型,方法一：常微分方程模型,分析：系統(tǒng)的輸入變量為u(t),輸出變量為y(t)，系統(tǒng)微分方程是： D6y+8.8*D5y+76.1*D4y+237.3*D3y+904.4*D2y+840*Dy+186.5*y=65*D4u+327*D3u+3699.6*D2u+1187.6*Du- 0.2*u 實(shí)現(xiàn)：運(yùn)用微分模塊、加法器、比例器構(gòu)造系統(tǒng)求解

14、 參見work21.mdl,方法二：傳遞函數(shù)模型,分析：將方程兩邊分別取拉普拉斯變換（Laplace),求傳遞函數(shù)G(s)。 y(t)?Y(s), u(t)?U(s), (a6s6+ a5s5+ a4s4+ a3s3+ a2s2+ a1s + a0) Y(s) =（b4s4+ b3s3+ b2s2+ b1s+ b0) U(

15、s) G(s)=Y(s)/U(s) = （b4s4+ b3s3+ b2s2+ b1s+ b0) / (a6s6+ a5s5+ a4s4+ a3s3+ a2s2+ a1s + a0) 實(shí)現(xiàn)：用傳遞函數(shù)模塊建立系統(tǒng) 參見work22.mdl,方法三：狀態(tài)方程模型,思路：復(fù)雜問題簡單化。用n個(gè)一階的微分方程去替換一個(gè)n階的微分方程。形式：X’(t)=AX

16、(t)+BU(t) Y(t)=CX(t)+DU(t) 其中X(t)= [x(1),x(2),…,x(n)] X’(t)=[x’(1),x’(2),…,x’(n)] =[x(2),x(3),…x(n),x’(n)] x’(n)用包含 x(1),x(2),…,x(n),u(t)的式子表示,狀

17、態(tài)方程模型(2),方法: 重在求A、B、C、D矩陣 由常微分方程寫出狀態(tài)方程組 x’(1)=x(2) x’(2)=x(3) x’(3)=x(4) x’(4)=x(5) x’(5)=x(6) x’(6)=-a0x(1)-a1x(2)-a2x(3)-a3x(4)-a4x(5) – a5x(6)+u y(t)=b0x(1)+ b1x(2)+b3x(3)+b4x(4) 則可以得到各矩陣,狀態(tài)方程模型(3)

18、,A=[0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 -a0 -a1 -a2 -a3 -a4 -a5]B=[0 0 0 0