当前位置:首页 > 微分方程作业
P10习题
1.用Euler法和改进的Euler法求u’=-5u (0≤t≤1),u(0)=1的数值解,步长h=0.1,0.05;并比较两个算法的精度。
解:function du=Euler_fun1(t,u) du=-5*u;clear;
h=0.1;tend=1;N=1/h;t(1)=0;u(1)=1; t=h.*(0:N); for n=1:N
u(n+1)=u(n)+h*Euler_fun1(t(n),u(n)); end
plot(t,u,'*');hold on for n=1:N
v(1)=u(n)+h*Euler_fun1(t(n),u(n)); for k=1:6
v(k+1)=u(n)+h/2*(Euler_fun1(t(n),u(n))+Euler_fun1(t(n+1),v(k))); end
u(n+1)=v(k+1); end
plot(t,u,'o');
sol=dsolve('Du=-5*u','u(0)=1'); u_real=eval(sol); plot(t,u_real,'r');
将上述 h 换为0.05得:
由图像知道:
显然改进的Euler法要比Euler法精确度要高;
‘’’
3.将u=-u(0≤t≤1),u(0)=0,u(0)=1化为一阶方程组,并用Euler法和改进的的Euler
法求解,步长h=0.1,0.05;并比较两个算法的精度。
解:
function du=fun31(y) du=y;
function dy=fun32(u) dy=-u; clear;
h=0.1;tend=1;N=1/h;t(1)=0;u(1)=0;y(1)=; t=h.*(0:N); for n=1:N
u(n+1)=u(n)+h*y(n); y(n+1)=y(n)+h*(-u(n)); end
plot(t,u,'*');hold on for n=1:N
v(1)=u(n)+h*fun31(y(n)); w(1)=y(n)+h*fun32(u(n));
for k=1:6
v(k+1)=u(n)+h/2*(fun31(y(n))+fun31(...w(k))); w(k+1)=y(n)+h/2*(fun32(u(n))+fun32(...v(k))); end
u(n+1)=v(k+1); y(n+1)=w(k+1); end
plot(t,u,'o');
sol=dsolve('D2u=-u','u(0)=0','Du(0)=1'; u_real=eval(sol); plot(t,u_real,'r');
将上述 h 换为0.05得:
由图像可以知道:
显然改进的Euler法要比Euler法精确度要高;
实习题(二)
1.取步长 h?0.1 ,分别用Euler 法和改进的Euler 法求下列初值问题的解,并与真解
相比较.
2x??u'?u?,0?x?1,(1)? u??u(0)?1,真解 u(x)?1?2x; 解:
function du=fun1(x,u) du=u-2*x/u; clear;
h=0.1;xend=1;N=1/h;x(1)=0;u(1)=1; x=h.*(0:N);
%——Eluer法——%
for n=1:N
u(n+1)=u(n)+h*fun1(x(n),u(n)); end
plot(x,u,'*'); hold on
%——改进的Eluer法——%
for n=1:N
v(1)=u(n)+h*fun1(x(n),u(n)); for k=1:6
v(k+1)=u(n)+h/2*(fun1(x(n),u(n))+fun1(x(n+1),v(k))); end
u(n+1)=v(k+1); end
plot(x,u,'o'); hold on
%——真解——%
u_real=sqrt(1+2*x); plot(x,u_real,'r');
由图像可以知道:
显然改进的Euler法要比Euler法精确度要高;
共分享92篇相关文档