中科院matlab课件第6章微分方程问题的解法幻灯片.ppt

编写函数： function [t,y]=nlbound(funcn,funcv,tspan,x0f,tol,varargin) t0=tspan(1);tfinal=tspan(2); ga=x0f(1); gb=x0f(2); m=1; m0=0; while (norm(m-m0)tol), m0=m; [t,v]=ode45(funcv,tspan,[ga;m;0;1],varargin); m=m0-(v(end,1)-gb)/(v(end,3)); end [t,y]=ode45(funcn,tspan,[ga;m],varargin); 例： 编写两个函数： function xdot=c7fun3(t,x) xdot=[x(2); 2*x(1)*x(2); x(4); 2*x(2)*x(3)+2*x(1)*x(4)]; function xdot=c7fun4(t,x) xdot=[x(2); 2*x(1)*x(2)]; [t,y]=nlbound(c7fun4,c7fun3,[0,pi/2],[-1,1],1e-8); plot(t,y); set(gca,xlim,[0,pi/2]); 精确解： 检验： y0=tan(t-pi/4); norm(y(:,1)-y0) ans = 1.6629e-005 norm(y(end,1)-1) ans = 5.2815e-006 6.4.2 线性微分方程的有限差分算法 把等式左边用差商表示。 编写函数： function [x,y]=fdiff(funcs,tspan,x0f,n) t0=tspan(1);tfinal=tspan(2); ga=x0f(1); gb=x0f(2); h=(tfinal-t0)/n; for i=1:n, x(i)=t0+h*(i-1); end, x0=x(1:n-1); t=-2+h^2*feval(funcs,x0,2); tmp=feval(funcs,x0,1); v=1+h*tmp/2; w=1-h*tmp/2; b=h^2*feval(funcs,x0,3); b(1)=b(1)-w(1)*ga; b(n-1)=b(n-1)-v(n-1)*gb; b=b; A=diag(t); for i=1:n-2, A(i,i+1)=v(i); A(i+1,i)=w(i+1); end y=inv(A)*b; x=[x tfinal]; y=[ga; y; gb]; 例： 编写函数： function y=c7fun5(x,key) switch key case 1, y=1+x; case 2, y=1-x; otherwise, y=1+x.^2; end [t,y]=fdiff(c7fun5,[0,1],[1,4],50); plot(t,y) 6.5 偏微分方程求解入门 6.5.1 偏微分方程组求解 函数描述： 边界条件的函数描述： 初值条件的函数描述： u0=pdeic(x) 例： 函数描述： function [c,f,s]=c7mpde(x,t,u,du) c=[1;1]; y=u(1)-u(2); F=exp(5.73*y)-exp(-11.46*y); s=F*[-1; 1]; f=[0.024*du(1); 0.17*du(2)]; 描述边界条件的函数 function [pa,qa,pb,qb]=c7mpbc(xa,ua,xb,ub,t) pa=[0; ua(2)]; qa=[1;0]; pb=[ub(1)-1; 0]; qb=[0;1]; 描述初值： function u0=c7mpic(x) u0=[1; 0]; 求解： x=0:.05:1; t=0:0.05:2; m=0; sol=pdepe(m,@c7mpde,@c7mpic,@c7mpbc,x,t); surf(x,t,sol(:,:,1))， figure; surf(x,t,sol(:,:,2)) 6.5.2 二阶偏微分方程的数学描述 椭圆型偏微分方程： 抛物线型偏微分方程： 双曲型偏微分方程： 特征值型偏微分方程： 6.5.3 偏微分方程的求解界面应用简介 6.5.3.1 偏微分方程求解程序概述 启动偏微