【正文】 Z=sin(r)./r。 meshz(X,Y,Z) 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 作圖舉例 (surf ) [X,Y]=meshgrid(8::8)。 r=sqrt(X.^2+Y.^2)+eps。 Z=sin(r)./r。 mesh(X,Y,Z) 例： 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 例： 繪制等高線 meshc 空間曲面作圖舉例 (meshc) [X,Y]=meshgrid(8::8)。 r=sqrt(X.^2+Y.^2)+eps。 y=[8::8]。 mesh, meshc, meshz,surf meshc 調用方式與 mesh 相同，在 mesh 基礎上增加等高線 meshz 調用方式與 mesh 相同，在 mesh 基礎上屏蔽邊界面 surf 調用方式與 mesh 相同，繪制三維 表面圖 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 空間三維作圖 (mesh) ? 網格生成函數(shù)： meshgrid x, y 為給定的向量， X, Y 是網格劃分后得到的網格矩陣 ? 繪制由函數(shù) z=z(x,y) 確定的曲面時，首先需產生一個網格矩陣，然后計算函數(shù)在各網格點上的值。 plot3(x,y,z)。 y=sin(t)。 例 ： 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 空間三維作圖 (plot3) ? 三維曲線 ： plot3 設三維曲線的參數(shù)方程為： x=x(t),y=y(t),z=z(t), 則其圖形可由下面的命令繪出 ： 例： 三維螺旋線 plot3的用法與 plot 類似 t=[0::10*pi]。 subplot(2,2,4)。 subplot(2,2,3)。 subplot(2,2,2)。 subplot(2,2,1)。 p 表示第 p 個繪圖子區(qū)域。)。 text(pi/2,cos(pi/2),39。cos(x)39。 plot(x,cos(x),’k*’)。 plot(x,cos(x),’r+:’)。 例 ： y=cos(x) 在 [0, 4*pi] 上的圖像 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty ? 點和線的基本屬性 圖形的屬性 (plot) ? plot(x,y,string) 其中 string 是用 單引號 括起來的字符串，用來指定圖形的屬性 （點、線的形狀和顏色） 屬性可以全部指定，也 可以只指定其中某幾個 排列順序任意 紅色、虛線、離散點用加號 藍色、點劃線 離散點為菱形 黑色、實線 離散點用星號 x=[0::2*pi]。 y=cos(x)。)。 plot(t,x(:,1),39。verderpol39。0]。mu=7。 clear。 mu*(1x(1)^2)*x(2) x(1)]。 ( ) ,d y d yyyd t d tyy???? ? ? ???? ? ? ??例 5： 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 數(shù)值求解舉例 ? 先編寫函數(shù)文件 function xprime=verderpol(t,x) global mu。 %此時 x=[0::] [T,Y] = solver(odefun,tspan,y0) 求初值問題 的數(shù)值解，求解范圍為 [0,] 22 2 201()dyy x xdxy? ? ? ? ???? ??例 4： 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 數(shù)值求解舉例 如果需求解的問題是 高階 常微分方程，則需將其化為 一階常微分方程組 ，此時需用 函數(shù)文件 來定義該常微分方程組。 [x,y]=ode23(fun,[0,],1)。y39。x39。2*y+2*x^2+2*x39。 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty Matlab提供的 ODE求解器 求解器 ODE類型 特點 說明 ode45 非剛性 單步法； 4， 5 階 RK 方法；累計截斷誤差為 (△ x)3 大部分場合的 首選方法 ode23 非剛性 單步法； 2， 3 階 RK 方法；累計截斷誤差為 (△ x)3 使用于精度較低的情形 ode113 非剛性 多步法； Adams算法；高低精度均可到 103～ 106 計算時間比 ode45 短 ode23t 適度剛性 采用梯形算法 適度剛性情形 ode15s 剛性 多步法； Gear’s 反向數(shù)值微分；精度中等 若 ode45 失效時，可嘗試使用 ode23s 剛性 單步法； 2 階 Rosebrock 算法；低精度 當精度較低時，計算時間比 ode15s 短 ode23tb 剛性 梯形算法；低精度 當精度較低時，計算時間比 ode15s短 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 參數(shù)說明 odefun 為 顯式常微分方程 ，可以用命令 inline 定義，或在 函數(shù)文件 中定義，然后通過函數(shù)句柄調用。) 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty RungeKutta 方法 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty Euler 法與 RK法誤差比較 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty Matlab函數(shù)數(shù)值求解 [T,Y] = solver(odefun,tspan,y0) 其中 y0 為初值條件， tspan為求解區(qū)間； Matlab在數(shù)值求解時 自動對求解區(qū)間進行分割 ， T (向量 ) 中返回的是分割點的值 (自變量 )， Y (向量 ) 中返回的是解函數(shù)在這些分割點上的函數(shù)值。 end plot(szj(:,1),szj(:,2), 39。 szj=[szj。 y=y+h*(l1+2*l2+2*l3+l4)/6。y39。x39。},{x+h/2,y+l2*h/2})。,39。 l3=subs(f,{39。y39。x39。},{x,y})。,39。 for i=1:n1 % i=1:n l1=subs(f,{39。 y=1。 % n=(ba)/h。 h=。 a=0。y+2*x/y^239。) 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty Euler折線法舉例（續(xù)） 解析解： 13352 233/xy e x??? ? ?????解析解 近似解 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty RungeKutta 方法 ? 為了減小誤差，可采用以下方法： ? 讓步長 h 取得更小一些； ? 改用具有較高精度的數(shù)值方法： ? 龍格 庫塔方法 RungeKutta (龍格 庫塔 ) 方法 ? 是 一類 求解常微分方程的數(shù)值方法 ? 有多種不同的迭代格式 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty RungeKutta 方法 ? 用得較多的是 四階 RK方法 0 0 11 1 2 3 4 ( 2 2 )/ 6( ), kkkky y x x x hy y h L L L L??? ? ??? ? ? ? ? ??12132432222( , )( / , / )( / , / )( , )kkkkkkkkL f x yL f x h y h LL f x h y h LL f x h y h L???? ? ???? ? ??? ? ? ??其中 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 四階 RK 方法 源程序 clear。 end szj plot(szj(:,1),szj(:,2),39。 szj=[szj。},{x,y})。,39。 for i=1:n1 % i=1:n y=y+h*subs(f,{39。 y=1。 % n=(ba)/h。 h=。 a=0。y+2*x/y^239。 0010||ttxy????? ??數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty ? Euler折線法 數(shù)值求解： Euler 折線法 ? RungeKutta（龍格 庫塔）方法 ? 利用 Matlab函數(shù)求數(shù)值解 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 初值問題的 Euler折線法 ? 具體步驟： 等距剖分： 0 1 2 1nna x x x x x b?? ? ? ? ? ? ?步長： 1 0 1 2 1( ) / , , , , ,kkh x x kn nba? ? ?? ? ? ?? 分割求解區(qū)間 ? 差商代替微商 1( ) ( )kkky x y xdyd x hx? ?? 1 ( ) ( ) 39。t39。y(0)=039。x(0)=139。Dyx3*y=039。Dx+5*x+y=exp(t)39。 12()ye?數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty dsolve 舉例 例 3： 求微分方程組 在初值條件 下的特解，并畫出解函數(shù)的圖形。x39。y(1)=2*exp(1)39。x*Dy+yexp(x)=039。 039。39。； D3y y39。； D2y y39。)。 ％ dy/dx = 2x dsolve(39。x39。Dy=2*x39。) syms x。,39。 22 xdy x y x edx??? y=dsolve(39。) 其中 y 為輸出， eq eq ...為微分方程， condcond ...為初值條件， v 為自變量。, ... ,39。,39。, ... ,39。,39。 微分方程求解 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 符號解法－－ dsolve 求解析解 ? dsolve 的使用 y=dsolve(39。 ? 由于實際應用的需要，人們必須求解微分方程。 S=symsum(f,n,1,inf) 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty 代數(shù)方程和微分方程 求解 ? 代數(shù)方程求解 solve(f,v)： 求方程關于指定自變量的解， f 可以是用字符串表示的方程 、 符號表達式 或 符號方程 ； ? 微分方程求解 數(shù)學建模專題一 － Matlab 基礎 Lxy, China Jiliang Universty ? 自牛頓發(fā)明微積分以來，微分方程在描述事物運動規(guī)律上已發(fā)揮了重要的作用。 S=symsum(f,n,1,inf) S100=symsum(f,n,1,100) symsum(f,v,a,b): 求和 symsum(f,a,b): 關于 默認變量 求和 ()bvafv??例：計算級數(shù) 及其前 100項的部分和 21