数学实验第三版（主编：李继成赵小艳）课后练习答案（十二）（1）（2）

C.L.L

已于 2024-06-23 09:38:00 修改

阅读量1k

点赞数 25

分类专栏：数学实验答案（matlab）文章标签：算法 matlab

于 2024-02-17 10:17:09 首次发布

本文链接：https://blog.csdn.net/2301_80199493/article/details/136135665

版权

数学实验答案（matlab）专栏收录该内容

31 篇文章 3 订阅

订阅专栏

实验十二：微分方程模型

练习一

1.用MATLAB软件求下列微分方程的通解或者初值解:

(1)

clc;clear;
syms y(x) 
y=x*diff(y)+2*y==1-1/x;
dsolve(y)

ans =

(x - 2)/(2*x) + C1/x^2

(2)

clc;clear;
syms f(x)
h=diff(f)+f*tan(x)==cos(x)^2;
y=dsolve(h)

y =sin(2*x)/2 + C1*cos(x)

(3)

clc;clear;
syms y(x)
h=diff(y,2)==sin(2*x)-y;
dy=diff(y);
y=dsolve(h,y(0)==0,dy(0)==1)

y =(5*sin(x))/3 - sin(2*x)/3

(4)

clc;clear;
syms y(x)
h=x*diff(y)-2*y==x^3*tan(x)*sec(x);
y=dsolve(h,y(pi/3)==2)

y =x^2/cos(x) - (x^2*(2*pi^2 - 18))/pi^2

2.高为16m,半径为5m的正圆柱形容器装满水,在该容器的侧面底部装有一个半径为0.02m的圆形排水孔,已知物理学中的托里拆利(Torriceli)定理:从小孔流出液体的流速正比于水面高度的平方根.完成下面的实验任务:

(1)容器内水面高度与时间t的函数关系是什么?多长时间容器内的水可以完全排空?

clc;clear;
syms h(t) k
m=-25*pi*diff(h)==k*sqrt(h)*pi*4*10^(-4);
h=dsolve(m,h(0)==16)
t=solve(h(2)==0,t)

h =

((k*t)/125000 + 4)^2

((k*t)/125000 - 4)^2

ans =

500000/k

(2)画出水面高度与时间的关系曲线.

经过查找资料我们得知，常数k=4.427188724235731

clc;clear;
syms h(t) k
m=-25*pi*diff(h)==k*sqrt(h)*pi*4*10^(-4);
h=dsolve(m,h(0)==16)
k=4.427188724235731;
f=eval(h(2));
ezplot(f,[0,99999]);

练习二

1.利用示例2中算得的数据,编程分别计算初值问题(12-11)的两种数值解与解析解在结点x=pi:0.1:2*pi处的距离和,看看与理论分析是否一致.

clc;clear;
format default
syms y(x)
m=diff(y)==3/x*y+x^3*(exp(x)+cos(x))-2*x;
j=dsolve(m,y(pi)==(exp(pi)+2/pi)*pi^3);
jj=matlabFunction(j);
x=pi:0.1:2*pi;
sum(jj(x))%解析解
y1=(exp(pi)+2/pi)*pi^3;
f=@(x,y)3/x*y+x^3*(exp(x)+cos(x))-2*x;
jj2=[y1];
for x=pi:0.1:(2*pi-0.1)
    y1=y1+0.1*f(x,y1);
    jj2=[jj2,y1];
end
sum(jj2)%数值解1
y1=(exp(pi)+2/pi)*pi^3;jj3=[y1];
for x=pi:0.1:2*pi-0.1
    jzj=y1+0.1*f(x,y1);
    yy=y1+0.1/2*(f(x,y1)+f(x+0.1,jzj));
    jj3=[jj3,yy];
    y1=yy;
end
sum(jj3) %数值解2

ans =

8.7052e+05

ans =

7.7816e+05

ans =

8.6997e+05

与理论方法相比，数值方法得到的解还是和理论有一定的差距，但是欧拉两步法比欧拉法误差小得多。

2.分别用欧拉方法、欧拉两步法和改进的欧拉方法求下面初值问题的数值解.

clc;clear;f=@(x,y)(x-sin(y))/cos(y);
%欧拉方法（显式）
y0=pi/4;y1=[y0];
for x=0:0.01:3/2-0.01
    y0=y0+0.01*f(x,y0);
    y1=[y1,y0];
end
plot(0:0.01:3/2,y1,'.');
%欧拉两步法
y0=pi/4;y2=[y0,y1(2)];
x=0:0.01:3/2-0.01;
for i=1:length(x)-1
    y=y2(i)+2*0.01*f(x(i+1),y2(i+1));
    y2=[y2,y];
end
hold on
plot(0:0.01:3/2,y2,'.');
%改进的欧拉方法
y1=pi/4;y3=[y1];x=0:0.01:3/2;
for i=1:length(x)-1
    ycz=y3(i)+0.01*f(x(i),y3(i));
    y=y3(i)+0.01/2*(f(x(i),y3(i))+f(x(i+1),ycz));
    y3=[y3,y];
end
hold on
plot(0:0.01:3/2,y3,'.');
legend('欧拉方法（显示）','欧拉两步法','改进的欧拉方法');

3.求下面方程的数值解,并画出数值解曲线

首先，我们要把原式进行换元使其变为一阶微分方程：

令，则原式化为，然后在一步步带回去即可。

clear;clc;%欧拉两步法
f=@(x,t)2*x*t/(1+x^2);
y1=3;y2=y1+0.01*0;
j1=[y1,y2];x=0:0.01:2;
for i=1:length(x)-2
    y=j1(i)+2*0.01*f(x(i+1),j1(i+1));
    j1=[j1,y];
end
y11=1;y22=y11+0.01*j1(1);
j11=[y11,y22];x=0:0.01:2;
for i=1:length(x)-2
    y=j11(i)+2*0.01*j1(i+1);
    j11=[j11,y];
end
j11
plot(0:0.01:2,j11);

j11 =

列 1 至 10

1.0000 1.0300 1.0600 1.0900 1.1200 1.1501 1.1802 1.2103 1.2405 1.2707

列 11 至 20

1.3010 1.3313 1.3617 1.3922 1.4227 1.4534 1.4840 1.5149 1.5458 1.5768

列 21 至 30

1.6079 1.6392 1.6706 1.7021 1.7337 1.7656 1.7975 1.8297 1.8618 1.8944

列 31 至 40

1.9269 1.9598 1.9926 2.0259 2.0592 2.0928 2.1265 2.1606 2.1947 2.2293

列 41 至 50

2.2638 2.2989 2.3339 2.3695 2.4050 2.4411 2.4772 2.5138 2.5504 2.5876

列 51 至 60