一维有限差分法之稳定性分析范文.docVIP

WORD格式 可编辑 PAGE 专业知识整理分享 20 13 ～20 14 学年第 一 学期实验报告 班级　学号：姓名：　实验时间： 11月29日 　成绩： 实验 项目 一维有限差分法之稳定性分析 所属 课程 微分方程数值解 实 验 目 的 了解一维有限差分法之算法构造，如Euler法，改进Euler等。 实 验 环 境 Matlab 7.0 实 验 内 容 掌握一维有限差分法之算法构造，用Euler法、Taylor法求解课本上45页实习题，并分析稳定性 实 验 过 程 第二题： （1）利用Euler法计算出U（2），用Milne迭代法+Jacobi法： 定义函数f2 function y=f2(t,u) y=t^3*u^3-t*u; % Milne迭代法+Jacobi方法，例2 %Milne迭代法+Jacobi迭代法，例2 n=input(n=);%n等分 T=input(T=);%时间上界 H=input(H=);%小区间H等分 IT=input(IT=);%Jacobi迭代次数，内迭代次数 h=T/n; hh=h/H; u=zeros(n+1,1); u(1)=0.5; x(1)=0.5; for mm=1:H x(mm+1)=x(mm)+hh*f2(mm*hh,x(mm)); end u(2)=x(H+1); for m=1:n-1%Milne开始 u0=u(m)+h*f2(m*h,u(m)); for k=1:IT%Jacobi%Jacobi迭代开始 w=-u(m)-4/3*h*f2((m+1)*h,u(m+1))-1/3*h*f2(m*h,u(m)); u1=1/3*h*f2((m+2)*h,u0)-w; u0=u1; end%Jacobi%Jacobi迭代结束 u(m+2)=u0; end%Milne结束 %u t=0:h:T; uu=1./sqrt(t.^2+1+3*exp(t.^2)); subplot(2,1,1) u-uu plot(t,u,bo-,t,uu,r*-) hold on u2=zeros(n+1,1); u2(1)=0.5; x(1)=0.5; for mm=1:H x(mm+1)=x(mm)+hh*f2(mm*hh,x(mm)); end u2(2)=x(H+1); for m=1:n-1%Milne开始 u0=u2(m)+h*f2(m*h,u2(m)); for k=1:IT%Jacobi迭代开始 w=-u2(m)-4/3*h*f2((m+1)*h,u2(m+1))-1/3*h*f2(m*h,u2(m)); u1=u0-(u0-1/3*h*f2((m+2)*h,u0)+w)/(1-1/3*h*(3*(m+2)^3*h^3*u0^2-(m+2)*h)); u0=u1; end%Jacobi%Jacobi迭代结束 u2(m+2)=u0; end%Milne结束 %u t=0:h:T; uu=1./sqrt(t.^2+1+3*exp(t.^2)); subplot(2,1,2) u2-uu plot(t,u2,go-,t,uu,r*-) 结果：（画图如下） （2）用Adams二步内插法+Newton法： n=input(n=); T=input(T=); H=input(H=); IT=input(IT=); h=T/n; hh=h/H; u=zeros(n+1,1); u(1)=0.5; x(1)=0.5; for mm=1:H x(mm+1)=x(mm)+hh*f2(mm*hh,x(mm)); end u(2)=x(H+1); for m=1:n-1 u0=u(m)+h*f2(m*h,u(m)); w=-u(m+1)-2/3*h*f2((m+1)*h,u(m+1))-1/12*h*f2(m*h,u(m)); for k=1:IT u1=u0-(u0-5/12*h*f2((m+2)*h,u0)+w)/(1-5/12*h*(3*(m+2)^3*h^3*u0^2-(m+2)*h)); u0=u1; end u(m+2)=u0; end t=0:h:T; uu=1./sqrt(t.^2+1+3*exp(t.^2)); plot(t,u,r*-,