常州大学MATLAB实习报告.docx

常州大学MATLAB实习报告.docx

《常州大学MATLAB实习报告.docx》由会员分享，可在线阅读，更多相关《常州大学MATLAB实习报告.docx（17页珍藏版）》请在冰豆网上搜索。

常州大学MATLAB实习报告

详细求解

3．已知t=an2+bn，测得对应数据如下：

（多项式插值interp1）

t=[0,20,40,60,80,100,120,140,160,183.5];

n=[0,1153,2045,2800,3466,4068,4621,5135,5619,6152];

试求a和b的值。

t=[0,20,40,60,80,100,120,140,160,183.5];

n=[0,1153,2045,2800,3466,4068,4621,5135,5619,6152];

t0=[0,20,40,60,80,100,120,140,160,183.5];

n0=[0,1153,2045,2800,3466,4068,4621,5135,5619,6152];

n=0:

0.001:

6152;

t=interp1（n0,t0,n,'spline'）;%interp1函数进行多项式插直

p=polyfit（n,t,2）%polyfit函数进行多项拟合

p=

0.00000.01440.0631

4．请用梯形法、辛普森法分别计算积分值

（trapz、quad）

f=inline（'sqrt（x.^2+x+1）','x'）;%inline定以内置函数

>>quad（f,0,1）%辛普森法

ans=

1.3369

>>x=0:

0.01:

1;y=sqrt（x.^2+x+1）;

trapz（x,y）%trapz为梯形法

ans=

1.3369

5计算二重积分

（使用函数dblquad）

dblquad（'x.^2+y.^2+x*y+2*x+y+1',0,1,0,2）%函数求解二重积分

ans=

10.3333

8用两种方法求解Ax=b的解。

（A为四阶随机矩阵，b为四阶向量，自己创建）。

方法1.A=rand（4）

A=

0.81470.63240.95750.9572

0.90580.09750.96490.4854

0.12700.27850.15760.8003

0.91340.54690.97060.1419

>>B=[1;2;3;4];x=inv（A）*B%求矩阵的逆

x=

73.6009

6.6966

-69.4126

3.4100

方法2.A（:

5）=B;rref（A）%将A化为阶梯状

ans=

1.000000073.6009

01.0000006.6966

001.00000-69.4126

0001.00003.4100

9.

，用两种方法求函数的根，并求其极值与零点。

求根.

solve（'x^3+（x-4/5）^2/（x+5/4）^3-5*x-5/x'）;

x=double（ans）;

>>fori=1:

length（x）

xx（i）=isreal（x（i））;

end

>>x=x（xx）

x=

2.4156

零点fzero（'x.^3+（x-0.8）.^2/（x+1.25）.^3-5*（x+1/x）',3）

ans=

2.4156

极值：

>>symsxy

>>y=x.^3+（x-0.8）.^2/（x+1.25）.^3-5*（x+1/x）;

>>dydx=diff（y）;

>>solve（dydx）;x=double（ans）;

>>fori=1:

length（x）

xx（i）=isreal（x（i））;

end

>>t=x（xx）

t=

-0.4694

-2.4039

>>x=t

（1）;y1=subs（y）

y1=

16.2832

>>x=t

（2）;y2=subs（y）

y2=

-6.4732

>>z=diff（diff（y））;

>>x=t

（1）;z1=subs（z）

z1=

205.8164

>>x=t

（2）;z2=subs（z）

z2=

-53.5382

函数f（x）有一个极大值点x1=-2.4039,极大值为-6.4732；一个极小值点x2=-2.4039,极小值为16.2832

10.f（x）的定义如下：

，写一个matlab函数func1实现该函数。

Functions11461101_10（x）

ifx<0&&x==-4

y=x^2+x-6;

elseifx>=0&x<10&x~=2&x~=3

y=x^2-5*x+6;

else

y=x.^2-x-1;

end

fprintf（'%d\n',y）;

13.写一个MATLAB函式pifun.m来计算下列级数：

f（n）=4*（1-1/3+1/5-1/7+1/9-1/11+...）

其中n为函式的输入，级数和f（n）则是函式的输出。

此外，你必须进行下列事项：

（1）使用tic和toc指令来测量pifun（100000）的计算时间。

如果你不知道如何使用这两个指令，请使用helptic及helptoc来查出它们的用法。

我的计算机是Pentium-450，所得的计算时间约为2秒。

请说明你的计算机规格以及其计算时间。

functionf=11461101_13（n）

l=1;s=0;

fori=1:

n

s=s+l/（2*i-1）;

l=-l;

end

s=4*s;

disp（'f（n）的值是'）;disp（num2str（s））

>>tic,s11461101_13（100000）,toc

f（n）的值是

3.1416

Elapsedtimeis0.028928seconds.

（2）使用flops指令来测量pifun（100000）所用到floatingpointoperations的计算次数。

如果你不知道如何使用这个指令，请使用helpflops来查出它的用法。

14.写一个MATLAB的递归函式fibo.m来计算Fibonacci数列，其定义如下：

fibo（n+2）=fibo（n+1）+fibo（n）

此数列的启始条件如下：

fibo

（1）=0,fibo

（2）=1.

functionf=fibo（n）

ifn==1

f=0;

elseifn==2

f=1;

else

f=fibo（n-1）+fibo（n-2）;

end

15．求下列函数的极小点：

1）

；

f='x

（1）^2+4*x

（2）^2+9*x（3）^2-2*x

（1）+18*x

（2）';

[x,fm]=fminsearch（f,[0,0,0]）

x=

1.0000-2.25000.0000

fm=

-21.2500

2）

；

f='x

（1）^2+3/2*x

（2）^2-2*x

（1）*x

（2）+x

（1）-2*x

（2）';

[x,fm]=fminsearch（f,[0,0]）

x=

0.50001.0000

fm=

-0.7500

3）

.

f='（x

（1）-1）^4+x

（2）^2';

>>[x,fm]=fminsearch（f,[0;1]）

x=

1.0000

-0.0000

fm=

4.2415e-17

第1），2）题的初始点可任意选取，

第3）题的初始点取为

.

16.解线性方程组

并求系数矩阵的行列式。

a=[51-101;103-12;-1-1053;0024-1];

>>rref（a）

ans=

1.00000001.4000

01.000000-5.9000

001.000000.1000

0001.0000-0.3000

17．设f（x,y）=4sin（x3y），求

。

symsfxy;

f=4*sin（x.^3*y）;

z=diff（diff（f,x）,y）;%diff函数进行求导

x=2;y=3;

subs（z）%subs函数进行元素替换

ans=

1.0636e+03

 18．求方程3x4+4x3-20x+5=0的所有解。

c=[340-205];roots（c）

ans=

-1.5003+1.5470i

-1.5003-1.5470i

1.4134

0.2539

19．对于迭代模型

取初值x0=0,y0=0,进行3000次迭代，对于k>1000,在（xk,yk）处亮一点（注意不要连线）可得所谓Henon引力线图。

 x

（1）=0;

>>y

（1）=0;

>>fori=1:

3000;

x（i+1）=1+y（i）-1.4*x（i）^2;

y（i+1）=0.3*x（i）;

holdon

plot（x（i）,y（i）,'*b'）

end

20．：

请设计一个程序，程序完成下列功能：

（1）让用户输入一个矩阵A；

（2）在A中找出小于0的矩阵元素的位置；

（3）在A中统计出等于0的元素的个数；

（4）显示A的行数和列数；

（5）找出矩阵A各元素中的最大值（显示值，不显示元素位置）。

functionf=s11461101_20（）

a=input（'输入一个矩阵A:

'）;

[m,n]=size（a）;

t=0;p=0;

fori=1:

m

forj=1:

n

ifa（i,j）<0p=p+1;

fprintf（'（矩阵中小于0的矩阵元素的位置%d,%d）\n',p,i,j）;

end

ifa（i,j）==0t=t+1;

end

end

end

B=a（:

）;

max=B

（1）;

fori=2:

m*n

ifB（i）>maxmax=B（i）;

end

end

ifp==0disp（'.'）;

end

fprintf（'矩阵中等于0的元素的个数%d\n',t）;

fprintf（'矩阵的行数和列数%d\n',m,n）;

fprintf（'矩阵中各元素的最大值%d\n',max）;

21．请设计一个程序，程序完成下列功能：

（1）让用户依次输入两个字符串s1和s2；

（2）比较两个字符串的长度并显示比较结果；

（3）判断s1与s2有没有长度在3个字符以上的相同子串，显示判断结果。

s1=input（'输入字符s1','s'）;

s2=input（'输入字符s2','s'）;

l1=length（s1）;

l2=length（s2）;

ifl1==l2

disp（'长度相等'）

elseifl1>l2

disp（'s1>s2'）

else

disp（'s1

end

e=strncmp（s1,s2,3）;

ife==0

disp（'没有长度在3个字符以上的相同子串'）

else

disp（'有长度在3个字符以上的相同子串'）

end

22．编写程序模拟杨氏双缝干涉

两相干光源到接收屏上P点距离r1=（D2+（y-a/2）2）1/2,r2=（D2+（y+a/2）2）1/2,相位差

Φ=2π（r2-r1）/λ,光强I=4I0cos2（Φ/2）编写程序

>>lam=500e-9;

a=2e-3;D=1;

ym=5*lam*D/a;xs=ym;

n=101;ys=linspace（-ym,ym,n）;%linspace用于产生-ymym之间的n点矢量

fori=1:

n

r1=sqrt（（ys（i）-a/2）.^2+D^2）;

r2=sqrt（（ys（i）+a/2）.^2+D^2）;

phi=2*pi*（r2-r1）./lam;

B（i,:

）=sum（4*cos（phi/2）.^2）;

end

N=255;

Br=（B/4.0）*N;

subplot（1,2,1）;%产生两幅图中的第一幅

image（xs,ys,Br）;

colormap（gray（N））;

subplot（1,2,2）;%产生第二幅

plot（B,ys）;

24．绘制三元函数

的可视化图形

x=-5:

0.05:

5;y=-5:

0.05:

5;z=-5:

0.05:

5;

[x,y,z]=meshgrid（x,y,z）;

v=x.^2+y.^2+z.^2;

isosurface（x,y,z,v,10）;

axisequal

25．绘制

的图象

[x,y,z,v]=flow;

q=1./sqrt（1-x）.*log（x-y+eps）-z;

p=patch（isosurface（x,y,z,v,0））;

isonormals（x,y,z,q,p）

set（p,'FaceColor','blue','EdgeColor','none'）;

view（3）

camlight

第二部分（8题选择4题）：

28.某公司投资2000万元建成一条生产线。

投产后，在时刻t的追加成本和追加收益分别为G（t）=

（百万元/年）,H（t）=

（百万元/年）。

试确定该生产线在合适何时停产可获最大利润？

最大利润是多少？

要求：

写出数学模型、M函数（如果需要的话）、运算命令和计算结果。

 解：

利润函数

（百万元）

f（t）=H（t）-G（t）=175-t-3t2/3,f（t）单调递减，则当f（t）=0时利润最大,即H（T）=G（T）时，R（t）取得最大利润。

命令行：

clear;close;fplot（'18-t^（2/3）',[0,20]）;gridon;holdon;fplot（'5+t+2*t^（2/3）',[0,20],’r’）;holdoff;[t,f,h]=fsolve（'18-x^（2/3）-5-x-2*x^（2/3）',4）求得t=4.6465t=linspace（0,t,100）;y=18-t.^（2/3）-5-t-2*t.^（2/3）;trapz（t,y）-20最大利润6.3232（百万元）

从图中，我们可以发现t约为4，通过命令行[t,f,h]=fsolve（'18-x^（2/3）-5-x-2*x^（2/3）',4），求得t=4.6465，

>>fplot（'18-t^（2/3）',[0,20]）;

>>gridon;

>>holdon;

>>fplot（'5+t+2*t^（2/3）',[0,20],'r'）;

>>legend（'H（t）','G（t）'）

>>holdoff;

>>[t,f,h]=fsolve（'18-x^（2/3）-5-x-2*x^（2/3）',4）

t=

4.6465

f=

1.1458e-13

h=

1

>>t=linspace（0,t,100）;

>>y=18-t.^（2/3）-5-t-2*t.^（2/3）;

>>trapz（t,y）-20

ans=

6.3232

29.一幢楼房的后面是一个很大的花园。

在花园中紧靠着楼房建有一个温室，温室高10英尺，延伸进花园7尺。

清洁工要打扫温室上方的楼房的窗户。

他只有借助于梯子，一头放在花园中，一头靠在楼房的墙上，攀援上去进行工作。

他只有一架20米长的梯子，你认为他能否成功？

能满足要求的梯子的最小长度是多少？

步骤：

1．先进行问题分析，明确问题；

2．建立模型，并运用Matlab函数求解；

3．对结果进行分析说明；

4．设计程序画出图形，对问题进行直观的分析和了解（主要用画线函数plot，line）

若梯子能恰好放置，则梯子的长度f与倾角x关系如下：

（其中a=

m,b=3.048m，0

）

接下来用Matlab求该函数的最小值：

>L=inline（'3.048/sin（x）+2.333/cos（x）','x'）;%inline为内置函数

x=fminbnd（L,0.01,pi/2-0.01）%fminbnd求最小值

l=subs（L）%subs求函数最小值

ezplot（'3.048/sin（x）+2.333/cos（x）',[0,pi/2]）%ezplot符号绘图函数

x=

0.8299

l=

7.5874

34.一半径为5m的球形水罐充满了水，底部有一半径为b=0.1m的小孔漏水，问多少时间以后，水面下降至离底部0.5m?

>>f=inline（'（（x-5）.^2-25）./（0.62*0.01*sqrt（2*9.8*x））'）;

>>quadl（f,10,0.5）

ans=

2.9889e+003

35.经调查发现，电饭锅销售速度与当时的销量成正比。

现在我们来建立一个数学模型以预测销量。

解：

设x（t）表示t时刻的销量,x0为初始时刻t0的销量,模型1（指数增长模型）

其中k为常数。

解得x（t）=x0exp（k（t-t0））。

当k>0,t时，x（t），这对于销售初期可认为是合适的，长期显然不合适。

设x为全部需要量，那么销售速度与当时的潜在需要量（1-x/x）成正比

模型2（阻滞增长模型）

设t0=0（年）,x0=1（万台）,x=100（万台）

k=0.9（年-1万台-1），

>>fplot（'exp（0.9*x）',[0,10]）;%模型1解析解

>>holdon;

>>[t,x]=ode45（inline（'0.9*x*（1-x/1000）','t','x'）,[010],1）;%模型2数值解

>>plot（t,x）;

>>axis（[01001500]）;

>>holdoff;

模型比较

短期预报二个模型相近,但作为长期预报,后者较前者合理。