ans =
0.4498 + 0.7623i 0.5526 + 0.2068i 0.6555 - 0.3487i 1.0185 + 0.0842i 1.2515 + 0.0228i 1.4844 - 0.0385i 1.5873 - 0.5940i 1.9503 - 0.1611i 2.3134 + 0.2717i
>> sqrt(A)
ans =
1.0000 1.4142 1.7321 2.0000 2.2361 2.4495 2.6458 2.8284 3.0000
实验四 求余弦的积分并绘出图像
一、实验目的
1、了解绘图工具的使用。 二、实验要求
1．(不定积分)用int计算下列不定积分，并用diff验证
dx?xsinxdx，?1?cosx，，，
2
输入以下指令： >> syms x;
>> f=x*sin(x^2); >> int(f,'x') ans =
-cos(x^2)/2 验证：
>> diff(-cos(x^2)/2) ans =
x*sin(x^2)
dx?1?cosx
输入以下指令： >> syms x;
>> f=1/(1+cos(x)); >> int(f) ans =
tan(x/2) 验证：
2xsinxdx?>> diff(tan(x/2)) ans =
tan(x/2)^2/2 + 1/2
2．(定积分)用trapz,int计算下列定积分(2个) 1sinx2?dx?0x,?0exsin(2x)dx 1sinx?0xdx 输入指令： >> syms x;
>> f=sin(x)/x; >> int(f,'x',0,1) ans =
sinint(1)
>> x=0:0.1:2*pi; f=exp(x).*sin(2*x); s=trapz(x,f) s =
-209.5581
0?2?exsin(2x)dxx2y2??1943．(椭圆的周长) 用定积分的方法计算椭圆的周长
t=0:0.001:2*pi; a=2; b=3;
x=a*sin(t); y=b*cos(t);
>> X=[0 x(1:end-1)]; Y=[0 y(1:end-1)]; x=x-X; y=y-Y;
d=sqrt(x.^2+y.^2); d=sum(d) d =
18.8651
??(1?x?y)dxdy 4．（二重积分）计算积分x?y?2y 指令为：
>> fun=inline('(1+x+y).*(x.^2+y.^2-2*y<=0)','x','y'); >> i=dblquad(fun,-1,1,0,2) i =
6.283
22exp(?x2)dx???1?x45． (广义积分)计算广义积分 指令为： >>syms x;
>> f=exp(-x.^2)./(1+x.^4); >> int(f,'x',-inf,inf) ans =
(4*pi^(1/2)*hypergeom([1], [5/4, 7/4], -1/4))/3 + (2^(1/2)*pi*(cos(1) - sin(1)))/2
?
实验五、六 用matlab求解多项式并用plot绘制函数图象
（常微分方程）
一、实验目的
1、 了解MATLAB中主要用dsolve求符号解析解，ode45,ode23,ode15s求数值解。 2、 s=dsolve(‘方程1’, ‘方程2’,…,’初始条件1’，’初始条件2’ …,’自变量’) 用字符串方程表示，自变量缺省值为t。导数用D表示，2阶导数用D2表示，以此类推。S返回解析解。在方程组情形，s为一个符号结构。 [tout,yout]=ode45(‘yprime’,[t0,tf],y0) 采用变步长四阶Runge-Kutta法和五阶Runge-Kutta-Felhberg法求数值解，yprime是用以表示f(t,y)的M文件名，t0表示自变量的初始值，tf表示自变量的终值，y0表示初始向量值。输出向量tout表示节点(t0,t1, …,tn)T,输出矩阵yout表示数值解，每一列对应y的一个分量。若无输出参数，则自动作出图形。 ode45是最常用的求解微分方程数值解的命令，对于刚性方程组不宜采用。ode23与ode45类似，只是精度低一些。ode12s用来求解刚性方程组，是用格式同ode45。可以用help dsolve, help ode45查阅有关这些命令的详细信息. 3、熟悉plot绘图 二、实验内容
1．求下列微分方程的解析解(2个)
指令为：
y=dsolve('D2y+2*Dy-3*y=exp(-3*x)','x') y =
C2*exp(x) - (x*exp(-3*x))/4 - exp(-3*x)/16 + C3*exp(-3*x)
指令为：
y=dsolve('D2y+Dy+y=cos(x)','y(0)=0','Dy(0)=1.5','x') y =
sin((3^(1/2)*x)/2)*(cos(x - (3^(1/2)*x)/2)/2 - cos(x + (3^(1/2)*x)/2)/2 + (3^(1/2)*cos(x - (3^(1/2)*x)/2))/3 + (3^(1/2)*cos(x +
(3^(1/2)*x)/2))/3 + (3^(1/2)*sin(x - (3^(1/2)*x)/2))/6 + (3^(1/2)*sin(x + (3^(1/2)*x)/2))/6) + (3^(1/2)*exp(-x/2)*sin((3^(1/2)*x)/2))/3 - (3^(1/2)*cos((3^(1/2)*x)/2)*((sin(x*(3^(1/2)/2 - 1))/2 -
cos(x*(3^(1/2)/2 - 1))*(3^(1/2)/2 - 1))/((3^(1/2)/2 - 1)^2 + 1/4) + (sin(x*(3^(1/2)/2 + 1))/2 - cos(x*(3^(1/2)/2 + 1))*(3^(1/2)/2 + 1))/((3^(1/2)/2 + 1)^2 + 1/4)))/3
2．求方程
(1?x2)y\?2xy',y(0)?1,y'(0)?3
的解析解和数值解,并进行比较（用plot绘图） 解析解：
>> s=dsolve('(1+x^2)*D2y-2*x*Dy','y(0)=1','Dy(0)=3','x') s =
x*(x^2 + 3) + 1 数值解：
先建立m文件
function dy=myfun_1(x,y) dy=zeros(2,1); dy(1)=y(2);
dy(2)=2*x*y(1)/(1+x^2); end
后再命令窗口输入：
>> [x,y]=ode23(@myfun_1,[0,1000],[0,1])； >>plot(x,y(:,1),'r+',x,y(:,2),'g*')
得到：
3．分别用ode45和ode15s求解Van-del-Pol方程
?d2x2dx?2?1000(1?x)?x?0dt?dt?x(0)?0,x'?0)?1??的数值解,并进行比较.（用plot绘图）
function dy=vdp1000(t,y)
dy=zeros(2,1); dy(1)=y(2);
dy(2)=1000*(1-y(1)^2)*y(2)-y(1); end
[T,Y]=ode15s('vdp1000',[0 3000],[0 1]); plot(T,Y(:,1),'-')