MATLAB程序设计与应用(第二版)实验答案

MATLAB程序设计与应用（第二版）实验参考答案

%实验一 MATLAB运算基础 %第1题 %（1）

z1=2*sin(85*pi/180)/(1+exp(2))

%（2）

x=[2,1+2i;-0.45,5];

z2=0.5*log(x+sqrt(1+x.^2)) %(3)

a=-3.0:0.1:3.0;

z3=(exp(0.3*a)-exp(-0.3*a))/2.*sin(a+0.3)+log((0.3+a)/2) %(4)

t=0:0.5:2.5;

z4=t.^2.*(t>=0&t<1)+(t.^2-1).*(t>=1&t<2)+(t.^2-2*t+1).*(t>=2&t<3)

%第2题

A=[12 34 -4;34 7 87;3 65 7]; B=[1 3 -1;2 0 3;3 -2 7]; A+6*B

A-B+eye(size(A)) A*B A.*B A^3 A.^3 A/B B\A [A,B]

[A([1,3],:);B^2]

%第3题

A=[1 2 3 4 5;6 7 8 9 10;11 12 13 14 15;16 17 18 19 20;21 22 23 24 25] B=[3 0 16;17 -6 9;0 23 -4;9 7 0;4 13 11] C=A*B

F=size(C)

D=C(F(1)-2:F(1),F(2)-1:F(2)) whos

%第4题 %(1):

A=100:999; B=rem(A,21);

C=length(find(B==0)) %(2):

A='lsdhKSDLKklsdkl';

k=find(A>='A'&A<='Z'); A(k)=[]

%实验二 MATLAB矩阵分析与处理 %第1题 E=eye(3); R=rand(3,2); O=zeros(2,3); S=diag([2,3]); A=[E,R;O,S]; A^2

B=[E,(R+R*S);O,S^2] %第2题 H=hilb(5) P=pascal(5) Hh=det(H) Hp=det(P) Th=cond(H) Tp=cond(P)

%第3题

A=fix(10*rand(5)) H=det(A)

Trace=trace(A) Rank=rank(A) Norm=norm(A)

%第4题

A=[-29,6,18;20,5,12;-8,8,5] [V,D]=eig(A)

%数学意义略 %第5题方法一 %(1):

A=[1/2,1/3,1/4;1/3,1/4,1/5;1/4,1/5,1/6]; b=[0.95,0.67,0.52]'; x=inv(A)*b %(2):

B=[0.95,0.67,0.53]'; x=inv(A)*B %(3): cond(A)

%第5题方法二 A=hilb(4) A(:,1)=[] A(4,:)=[]

B=[0.95,0.67,0.52]';

X=inv(A)*B

B1=[0.95,0.67,0.53]'; X1=inv(A)*B1 N=cond(B) N1=cond(B1)

Na=cond(A) %矩阵A为病态矩阵 %第6题

A=[1,4,9;16,25,36;49,64,81] B=sqrtm(A)

C=sqrt(A) %sqrtm函数是以矩阵为单位进行计算，sqrt函数是以矩阵中的元素进行计算

%实验三 选择程序结构设计 %第1题程序一

x=[-5.0,-3.0,1.0,2.0,2.5,3.0,5.0];

y=[]; %建立存放所有y值的矩阵 for x0=x

if x0<0&x0~=-3

y=[y,x0*x0+x0-6];

elseif x0>=0&x0<5&x0~=2&x0~=3 y=[y,x0*x0-5*x0+6]; else

y=[y,x0*x0-x0-1]; end end

x %输出所有x y %输出所有y %第1题程序二

x=[-5,-3,1,2,2.5,3,5]; y=[];

for a=1:7

if x(a)<0&x(a)~=-3

y=[y,(x(a))^2+x(a)-6];

elseif x(a)>=0&x(a)<5&x(a)~=2&x(a)~=3 y=[y,(x(a))^2-5*x(a)+6]; else

y=[y,x(a)*x(a)-x(a)-1]; end end

%第2题程序一

x=input('请输入一个百分制成绩：'); if x>100|x<0

disp('您输入的成绩不是百分制成绩，请重新输入。'); else

if x<=100&x>=90 disp('A');

elseif x<=89&x>=80 disp('B');

elseif x<=79&x>=70 disp('C');

elseif x<=69&x>60 disp('D'); else

disp('E'); end end

%第2题程序二

s=input('请输入一个成绩（0分到100分之间）：'); %s用于存放成绩 while 1 %判断输入成绩的合理性 if s<0|s>100

disp('输入的成绩需在0到100之间，请重新输入：') s=input('请输入一个成绩（0分到100分之间）：'); else

break; end end

switch fix(s/10) %对成绩做出等级判断 case {9,10} disp('A') case 8

disp('B') case 7

disp('C') case 6

disp('D') otherwise disp('E') end

%第3题

n=input('请输入员工工号：'); h=input('该员工工作时数是：'); if h>120

x=(h-120)*84*(1+0.15)+120*84; elseif h<60 x=h*84-700; else x=h*84; end

disp([num2str(n),'号员工','的应发工资为',num2str(x)]); %第4题（还可以用switch语句实现）

a=fix(10+(99-10)*rand(1,2)) %产生两个随机整数

x=a(1); y=a(2);

t=input('请输入运算符号：','s'); if t=='+' z=x+y;

elseif t=='-' z=x-y;

elseif t=='*' z=x*y;

elseif t=='/' z=x/y; end

disp([num2str(x),t,num2str(y),'=',num2str(z)]) %结果

%第5题

a=rand(5,6) %产生5x6的随机矩阵 n=input('请输入您要输出矩阵的第几行：'); if n>5

disp('超出了矩阵的行数，矩阵的最后一行为：') a(5,:) else

disp(['矩阵的第',num2str(n),'行为：']) a(n,:) end

%实验四 循环结构程序设计 %第1题程序一 s=0;

n=input('n=?'); for i=1:n

s=s+1/i/i; end

PI=sqrt(6*s) pi

%第1题程序二

n=input('n=?'); a=1:n;

输出运算

b=1./a.^2;

PI=sqrt(6*sum(b)) pi

%第2题 y=0; n=1;

while(y<3)

y=y+1/(2*n-1); n=n+1; end

y=y-1/(2*(n-1)-1) n=n-2

%第3题

a=input('a=?'); b=input('b=?'); Xn=1;

Xn1=a/(b+Xn); n=0;

while abs(Xn1-Xn)>1e-5 Xn=Xn1;

Xn1=a/(b+Xn); n=n+1; if n==500 break; end end n Xn1

r1=(-b+sqrt(b*b+4*a))/2 r2=(-b-sqrt(b*b+4*a))/2 %第4题

for i=1:100 if i==1 f(i)=1;

elseif i==2 f(i)=0;

elseif i==3 f(i)=1; else

f(i)=f(i-1)-2*f(i-2)+f(i-3); end end max(f) min(f)

sum(f)

length(find(f>0)) length(find(f==0)) length(find(f<0)) %第5题 s=0;n=0; for i=2:49 b=i*(i+1)-1; m=fix(sqrt(b)); for j=2:m

if rem(b,j)==0 break end end

if j==m n=n+1; s=s+b; end end n s

%实验五 函数文件 %第1题

function y=mat1(x) %建立函数文件mat1.m y=[exp(x),log(x),sin(x),cos(x)]; %在命令窗口调用上述函数文件： y=mat1(1+i)

%第2题程序一

function [a,b,N,M]=shiyanwu2(m,n,t)

A=[m*cos(t*pi/180),-m,-sin(t*pi/180),0;m*sin(t*pi/180),0,cos(t*pi/180),0;0,n,-sin(t*pi/180),0;0,0,-cos(t*pi/180),1]; B=[0,9.8*m,0,9.8*n]; C=inv(A)*B'; a=C(1); b=C(2); N=C(3); M=C(4);

%在命令窗口调用该函数文件: m1=input('m1='); m2=input('m2=');

theta=input('theta=');

[a1,a2,N1,N2]=shiyanwu2(m1,m2,theta)

%第2题程序二

function X=mat2(m1,m2,t) g=9.8;

A=[m1*cos(t*pi/180),-m1,-sin(t*pi/180),0;m1*sin(t*pi/180),0,cos(t*pi/180),0;0,m2,-sin(t*pi/180),0;0,0,-cos(t*pi/180),1]; B=[0;m1*g;0;m2*g]; X=inv(A)*B;

%在命令窗口调用该函数文件： X=mat2(1,1,60)

%第3题

function flag=mat3(x) flag=1;

for i=2:sqrt(x) if rem(x,i)==0 flag=0; break; end end

%在命令窗口调用该函数文件： for i=10:99

j=10*rem(i,10)+fix(i/10); if mat3(i)&mat3(j) disp(i) end end

%第4题

function y=fx(x)

y=1./((x-2).^2+0.1)+1./((x-3).^4+0.01); %在命令窗口调用该函数文件： y=fx(2)

a=[1,2;3,4]; y=fx(a)

%第5题 %(1)

function f1=mat5(n) f1=n+10*log(n*n+5);

%在命令窗口中调用该函数文件：

y=mat5(40)/(mat5(30)+mat5(20)) %(2)方法一

function f2=mat6(n) f2=0;

for i=1:n

f2=f2+i*(i+1); end

%在命令窗口中调用该函数文件如：

y=mat6(40)/(mat6(30)+mat6(20)) %(2)方法二

function f2=mat7(n) i=1:n;

m=i.*(i+1); f2=sum(m); end

%在命令窗口中调用该函数文件如：

y=mat7(40)/(mat7(30)+mat7(20))

%实验六 高层绘图操作 %第1题

x=linspace(0,2*pi,101);

y=(0.5+3*sin(x)./(1+x.^2)).*cos(x); plot(x,y)

%第2题 %（1）

x=linspace(-2*pi,2*pi,100); y1=x.^2;

y2=cos(2*x); y3=y1.*y2;

plot(x,y1,'b-',x,y2,'r:',x,y3,'y--'); text(4,16,'\leftarrow y1=x^2');

text(6*pi/4,-1,'\downarrow y2=cos(2*x)');

text(-1.5*pi,-2.25*pi*pi,'\uparrow y3=y1*y2'); %（2）

x=linspace(-2*pi,2*pi,100); y1=x.^2;

y2=cos(2*x); y3=y1.*y2;

subplot(1,3,1);%分区 plot(x,y1);

title('y1=x^2');%设置标题 subplot(1,3,2); plot(x,y2);

title('y2=cos(2*x)'); subplot(1,3,3); plot(x,y3);

title('y3=x^2*cos(2*x)'); %（3）

x=linspace(-2*pi,2*pi,20); y1=x.^2;

subplot(2,2,1);%分区

bar(x,y1);

title('y1=x^2的条形图');%设置标题 subplot(2,2,2); stairs(x,y1);

title('y1=x^2的阶梯图'); subplot(2,2,3); stem(x,y1);

title('y1=x^2的杆图'); subplot(2,2,4);

fill(x,y1,'r');%如果少了'r'则会出错 title('y1=x^2的填充图'); %其他的函数照样做。 %第3题

x=-5:0.01:5;

y=[];%起始设y为空向量 for x0=x

if x0<=0 %不能写成x0=<0

y=[y,(x0+sqrt(pi))/exp(2)]; %将x对应的函数值放到y中 else

y=[y,0.5*log(x0+sqrt(1+x0^2))]; end end

plot(x,y)

%第4题

a=input('a='); b=input('b='); n=input('n=');

t=-2*pi:0.01:2*pi; r=a*sin(b+n*t); polar(t,r)

%第5题

x=linspace(-5,5,21); y=linspace(0,10,31);

[x,y]=meshgrid(x,y);%在[-5,5]*[0,10]的范围内生成网格坐标 z=cos(x).*cos(y).*exp(-sqrt(x.^2+y.^2)/4); subplot(2,1,1); surf(x,y,z); subplot(2,1,2);

contour3(x,y,z,50);%其中50为高度的等级数，越大越密 %第6题

ezsurf('cos(s)*cos(t)','cos(s)*sin(t)','sin(s)',[0,0.5*pi,0,1.5*pi]); %利用ezsurf隐函数

shading interp %进行插值着色处理

%实验七 低层绘图操作 %第1题

h=figure('MenuBar','figure','color','r','WindowButtonDownFcn','disp(''Left Button Pressed'')')

%第2题

x=-2:0.01:2;

y=x.^2.*exp(2*x); h=line(x,y);

set(h,'color','r','linestyle',':','linewidth',2) text(1,exp(2),'y=x^2*exp(2*x)')

%第3题

t=0:0.00001:0.001; [t,x]=meshgrid(t);

v=10*exp(-0.01*x).*sin(2000*pi*t-0.2*x+pi); axes('view',[-37.5,30]); h=surface(t,x,v);

title('v=10*exp(-0.01*x).*sin(2000*pi*t-0.2*x+pi)'); xlabel(Ct'),ylabel('x'),zlabel('v') %第4题

x=0:0.01:2*pi; y1=sin(x); y2=cos(x); y3=tan(x); y4=cot(x);

subplot(2,2,1); plot(x,y1);

subplot(2,2,2); plot(x,y2);

subplot(2,2,3); plot(x,y3);

subplot(2,2,4); plot(x,y4);

%第5题

cylinder(5);

light('Position',[0,1,1]); material shiny

%实验八 数据处理与多项式运算 %第1题 %(1)

A=rand(1,30000); b=mean(A) std(A,0,2) %(2)

max(A) min(A) %(3) n=0;

for i=1:30000 if A(i)>0.5 n=n+1; end end

p=n/30000

%第2题 %(1)

A=45+51*rand(100,5); [Y,U]=max(A) [a,b]=min(A) %(2)

m=mean(A) s=std(A) %(3)

sum(A,2)

[Y,U]=max(ans) [a,b]=min(ans) %(4)

[zcj,xsxh]=sort(ans) %第3题 h=6:2:18;

x=6.5:2:17.5;

t1=[18,20,22,25,30,28,24]; t2=[15,19,24,28,34,32,30]; T1=spline(h,t1,x) T2=spline(h,t2,x) %第4题

x=1:0.1:101; y1=log10(x);

p=polyfit(x,y1,5) y2=polyval(p,x);

plot(x,y1,':',x,y2,'-') %第5题 %(1)

p1=[1,2,4,0,5]; p2=[1,2]; p3=[1,2,3];

p=p1+[0,conv(p2,p3)] %为使两向量大小相同，所以补0 %(2)

A=roots(p) %(3)

A=[-1,1.2,-1.4;0.75,2,3.5;0,5,2.5]; polyval(p,A) %(4)

polyvalm(p,A)

%实验九 数值微积分与方程数值求解 %第1题

for x=1:3

fx=[x,x^2,x^3;1,2*x,3*x;0,2,6*x]; diff(fx) end

%第2题 %（1）

f1=inline('sqrt(cos(x.*x)+4*sin((2*x).^2)+1)'); I1=quad(f1,0,2*pi)

%（2）

f2=inline('log(1+x)./(1+x.^2)'); I2=quadl(f2,0,1)

%第3题

A=[6,5,-2,5;9,-1,4,-1;3,4,2,-2;3,-9,0,2]; b=[-4;13;1;11]; x1=A\b

x2=inv(A)*b [L,U]=lu(A); x3=U\(L\b)

%第4题

%首先建立求解线性方程组的函数文件line_solution.m。

function [x,y]=line_solution(A,b) [m,n]=size(A); y=[];

if norm(b)>0 %非齐次方程组

if rank(A)==rank([A,b])

if rank(A)==n %有惟一解

disp('原方程组有惟一解x'); x=A\b;

else %方程组有无穷多个解，基础解系

disp('原方程组有无穷个解，特解为x，其齐次方程组的基础解系为y'); x=A\b;

y=null(A,'r'); end else

disp('方程组无解'); %方程组无解 x=[];

end

else %齐次方程组

disp('原方程组有零解x');

x=zeros(n,1); %0解 if rank(A)<n

disp('方程组有无穷个解，基础解系为y'); %非0解 y=null(A,'r'); end end

A=[2,7,3,1;3,5,2,2;9,4,1,7]; b=[6,4,2]';

[x,y]=line_solution(A,b); x,y %第5题 %（1）

z=fzero('3*x+sin(x)-exp(x)',1.5)

%（2）

%首先建立函数文件myfun.m。

function F=myfun(X) x=X(1); y=X(2); z=X(3);

F(1)=sin(x)+y^2+log(z)-7; F(2)=3*x+2^y-z^3+1; F(3)=x+y+z-5;

X=fsolve('myfun',[1,1,1],optimset('Display', 'off'))

%第6题 %（1）

f=inline('(x^3+cos(x)+x*log(x))/exp(x)'); fminbnd(f,0,1)

%（2）

%建立函数文件fx12.m：

function f=fx12(u) x1=u(1);x2=u(2);

f=2*x1.^3+4*x1*x2.^3-10*x1.*x2+x2.^2;

[U,fmin]=fminsearch('fx12',[0,0])

%第7题

%建立函数文件ztf.m

function yy=ztf(x,y)

yy=[(5*y(1)-y(2))/x;y(1)];

[x,y]=ode23(@ztf,[0,5],[0;0])

%第8题

%建立函数文件ztfun.m

function yy=ztfun(t,y)

yy=[y(2)*y(3);-y(1)*y(3);-0.51*y(1)*y(2)];

[t,y]=ode23(@ztfun,[0,5],[0;1;1]) plot3(y(:,1),y(:,2),y(:,3))

%实验十 符号计算基础与符号微积分 %第1题

x=sym('6'); y=sym('5');

z=(x+1)/(sqrt(3+x)-sqrt(y))

%第2题

z=sym('x*x*x*x-y*y*y*y'); factor(z)

num=sym('5135'); factor(num)

%第3题

syms beta1 beta2

f1=sym('sin(beta1)*sin(beta2)-cos(beta1)*cos(beta2)'); simplify(f1) simple(f1) sym x

f2=sym('(4*x*x+8*x+3)/(2*x+1)'); simplify(f2) simple(f2)

%第4题

syms a b c d e f g h i P1=[0,1,0;1,0,0;0,0,1]; P2=[1,0,0;0,1,0;1,0,1]; A=[a,b,c;d,e,f;g,h,i]; B=P1*P2*A inv(B) tril(B) det(B)

%第5题

syms x a t

f1=(x*(exp(sin(x))+1)-2*(exp(tan(x))-1))/sin(x)^3 limit(f1,x,0) %求(1)

f2=(sqrt(pi)-sqrt(acos(x)))/sqrt(x+1); limit(f2,x,-1,'right') %求(2) f3=(1-cos(2*x))/x; diff(f3,x)

diff(f3,x,2) %求(3) A=[power(a,x),t^3;t*cos(x),log(x)]; diff(A,x) diff(A,t,2)

diff(diff(A,x),t) %求(4) f=(x*x-2*x)*exp(-x*x-y*y-x*y);%求(5) df51=-diff(f,y)/diff(f,x) df52=diff(diff(f,x),y)

%第6题 sym x

format short e

int('1/(1+x^4+x^8)',x)

int('1/((asin(x))^2*sqrt(1-x*x))',x) int('(x*x+1)/(x^4+1)',0,inf)

a=int('exp(x)*(1+exp(x))^2',0,log(2)) eval(a)

%实验十一 级数与方程符号求解 %第1题 syms n x;

symsum(1/(2*n-1),n,1,10)

s=symsum(n^2*x^(n-1),n,1,inf) x=sym('5'); x=1/x;

s=symsum(n^2*x^(n-1),n,1,inf)/5

%第2题

x=sym('x'); f=log(x);

taylor(f,6,1) %第3题

x=solve(sym('log(1+x)-5/(1+sin(x))=2')) x=solve(sym('x^2+9*sqrt(x+1)-1')) x=solve('3*x*exp(x)+5*sin(x)-78.5')

[x y]=solve(sym('sqrt(x^2+y^2)-100=0,3*x+5*y-8=0'))

%第4题

y=dsolve('D2y+4*Dy+29*y','y(0)=0,Dy(0)=15','x')

%第5题

[x,y,z]=dsolve('Dx-2*x+3*y-3*z','Dy-4*x+5*y-3*z','Dz-4*x+4*y-2*z','t')

%实验十二 菜单与对话框设计 %第1~4题 （略）

%第5题

%实验十三 Simulink的应用 %第1题

（1）启动Simulink并打开模型编辑窗口。

（2）将所需模块添加到模型中。单击模块库浏览器中的Sources，在右边的窗口中找到Sine Wave模块，然后用鼠标将其拖到模型编辑窗口，再重复4次，得到5个正弦源。同样，在Math Operations中把Add模块拖到模型编辑窗口，在Sinks中把Scope模块拖到模型编辑窗口。

（3）设置模块参数并连接各个模块组成仿真模型。先双击各个正弦源，打开其Block Parameters对话框，分别设置Frequency（频率）为2*pi、6*pi、10*pi、14*pi和18*pi，设置Amplitude（幅值）为1、1/3、1/5、1/7和1/9，其余参数不改变。对于求和模块，将符号列表List of signs设置为+++++。

设置模块参数后，用连线将各个模块连接起来组成仿真模型，如图所示。

图 x( t)的仿真模型

（4）设置系统仿真参数。单击模型编辑窗口Simulink菜单中的Configuration Parameters命令，打开仿真参数设置对话框，选择Solver选项卡。在Star time和Stop time两个编辑框内分别设置起始为0，停止时间为1秒。把算法选择中的Type 设为Fixed-step（固定步长算法），并在其右栏的具体算法框选择ode5（Dormand-Prince），即5阶Runge-Kutta算法，再把Fixed step size设置为0.001秒。

（5）开始系统仿真。单击模型编辑窗口中的Start simulation按钮或选择模型编辑窗口Simulink菜单中的Start命令开始系统仿真。

（6）观察仿真结果。系统仿真结束后，双击仿真模型中的示波器模块，得到仿真结果。单击示波器窗口工具栏上的Autoscale按钮，可以自动调整坐标来使波形刚好完整显示，这是由5次谐波合成的方波。

%第2~6题 （略）

%实验十四 外部程序接口 %第1题

clear all

fid=fopen('SY12_5.m','r'); textall=fscanf(fid,'%c',inf); ch=zeros(26,1); for n=1:size(textall,2) if isletter(textall(n))

if textall(n)>='a' && textall(n)<='z' textall(n)=char(textall(n)-('a'-'A')); end

m=textall(n)-'A'+1; ch(m)=ch(m)+1; end end ch'

%第2~5题 （略）

%实验十五 综合实验 %第1题

%定义函数文件

function zp=fs(z,n)

z=0;

for k=1:n;

z=1+imag(z)-1.4*real(z)^2+i*real(z)*0.3; zp=[zp,z]; end

%在主窗口调用该函数： n=30000; z=0;

zp=fs(z,n);

plot(zp,'p','markersize',2,'color','r') %第2题

quad('1/sqrt(2*pi)*exp(-x.*x/2)',0,1) syms x

int(1/sqrt(2*pi)*exp(-x.*x/2),0,1)

%第3题

%定义函数文件

function dx=apollo(t,x) mu=1/82.45; mustar=1-mu;

r1=sqrt((x(1)+mu)^2+x(3)^2);

r2=sqrt((x(1)-mustar)^2+x(3)^2);

dx=[x(2);2*x(4)+x(1)-mustar*(x(1)+mu)/r1^3-mu*(x(1)-mustar)/r2^3; x(4);-2*x(2)+x(3)-mustar*x(3)/r1^3-mu*x(3)/r2^3]; %在主窗口调用该函数：

x0=[1.2;0;0;-1.04935751];

options=odeset('reltol',1e-8);

%该命令的另一种写法是options=odeset;options.reltol=1e-8; [t,y]=ode45(@apollo,[0,20],x0,options);

plot(y(:,1),y(:,3))%绘制x和y图像，也就是卫星的轨迹 title('Appollo卫星运动轨迹') xlabel('X') ylabel('Y')

%第4题

%这是一个力矩平衡问题，可列方程：500×1.5×cosα=400×2.0×cos(π/3-α) %解得α＝26.8021°。 %程序一：迭代法解方程 x1=input('x1=?');

x=pi/3.0-acos(15.0*cos(x1)/16); while abs(x-x1)>=10e-8 x1=x;

x=pi/3.0-acos(15.0*cos(x1)/16); if x>pi/3.0 disp('error');

end

y=x1*180/pi

%程序二：利用绘图功能

%分别绘制两个小孩所产生力矩随α变化的曲线，两曲线的交点即是跷跷板平衡时的α alpha=linspace(0,pi/2,80000);

m1=500*1.5*cos(alpha); %第一个小孩产生的力矩 m2=400*2*cos(pi/3-alpha); %第二个小孩产生的力矩 k=find(abs(m1-m2)<1e-2); %找平衡点 alpha0=alpha(k);

m0=500*1.5*cos(alpha0);

plot(alpha,m1,alpha,m2,alpha0,m0,'kp'); alpha=alpha0*180/pi %化为角度 %第5题

clear; close;

fplot('18-t^(2/3)',[0,20]);grid on;hold on; fplot('5+t+2*t^(2/3)',[0,20],'r');hold off; %发现t约为4

[t,f,h]=fsolve('18-x^(2/3)-5-x-2*x^(2/3)',4) %求得t=4.6465

t=linspace(0,t,100); y=18-t.^(2/3)-5-t-2*t.^(2/3); trapz(t,y)-20

%最大利润6.3232(百万元)

本文来源：https://www.bwwdw.com/article/9q61.html