MATLAB-实验报告

1、.MATLAB 语言与应用实验报告参考教材高等应用数学问题的MATLAB求解。第一部分 MATLAB语言编程、科学绘图与基本数学问题求解（4学时）主要内容：掌握MATLAB语言编程基础、科学绘图方法、微积分问题、线性代数问题等基本数学问题的求解与应用。1-1 安装MATLAB软件，应用demo命令了解主要功能，熟悉基本功能，会用help命令。>>demo>> 运行结果Figure 1 demo命令弹出窗口>> help demo>> 运行结果 demo Access product demos via Help browser. demo ope

2、ns the Help browser and selects the MATLAB Demos entry in the table of contents. demo SUBTOPIC CATEGORY opens the Demos entry to the specified CATEGORY. CATEGORY is a product or group within SUBTOPIC. SUBTOPIC is 'matlab', 'toolbox', 'simulink', or 'blockset'. When SU

3、BTOPIC is 'matlab' or 'simulink', do not specify CATEGORY to show all demos for the product. Examples: demo 'matlab' demo 'toolbox' 'signal' demo 'matlab' 'getting started' See also echodemo, grabcode, help, helpbrowser. Reference page in Help

4、browser doc demo1-2 用MATLAB语句输入矩阵和 ， 前面给出的是矩阵，如果给出命令将得出什么结果？>>A=1 2 3 4;4 3 2 1;2 3 4 1;3 2 4 1;>>B=1+4j 2+3j 3+2j 4+j;4+j 3+2j 2+3j 1+4j;2+3j 3+2j 4+j 1+4j;3+2j 2+3j 4+j 1+4j;>>A,B>> 运行结果A = 1 2 3 4 4 3 2 1 2 3 4 1 3 2 4 1B = 1.0000 + 4.0000i 2.0000 + 3.0000i 3.0000 + 2.0000

5、i 4.0000 + 1.0000i 4.0000 + 1.0000i 3.0000 + 2.0000i 2.0000 + 3.0000i 1.0000 + 4.0000i 2.0000 + 3.0000i 3.0000 + 2.0000i 4.0000 + 1.0000i 1.0000 + 4.0000i 3.0000 + 2.0000i 2.0000 + 3.0000i 4.0000 + 1.0000i 1.0000 + 4.0000i>> A(5,6)=5 %增补A矩阵的第5行，第6列为5；矩阵空余处均按0增补A = 1 2 3 4 0 0 4 3 2 1 0 0 2 3 4

6、 1 0 0 3 2 4 1 0 0 0 0 0 0 0 51-3 假设已知矩阵，试给出相应的MATLAB命令，将其全部偶数行提取出来，赋给矩阵，用命令生成矩阵，用上述命令检验一下结果是不是正确。>>A=magic(8);>>B=A(2:2:end,:);>>A(B(1:end,1)%提取B中每行首列元素，查看其在矩阵A中的对应位置>> 运行结果A = 64 2 3 61 60 6 7 57 9 55 54 12 13 51 50 16 17 47 46 20 21 43 42 24 40 26 27 37 36 30 31 33 32 34 3

7、5 29 28 38 39 25 41 23 22 44 45 19 18 48 49 15 14 52 53 11 10 56 8 58 59 5 4 62 63 1B = 9 55 54 12 13 51 50 16 40 26 27 37 36 30 31 33 41 23 22 44 45 19 18 48 8 58 59 5 4 62 63 1ans =%提取B中每行首列元素，查看其在矩阵A中的对应位置，其显示结果如下 2 4 6 81-4 用数值方法可以求出，试不采用循环的形式求出和式的数值解。由于数值方法是采用double形式进行计算的，难以保证有效位数字，所以结果不一定精确。试

8、采用运算的方法求该和式的精确值。>>i=1:63;s=sum(2.i);s%对2i进行级数求和的数值解>>j=1:63;s1=sum(sym(2).j);s1%对2i进行级数求和的解析解>>e=s-s1%解析解与数值解的差异>> 运行结果s =%数值解 1.8447e+019 s1 =%解析解 18446744073709551614 e =%误差 21-5 选择合适的步距绘制出下面的图形。（1），其中； （2），其中。>>t1=-1:0.01:-0.2,-0.2:0.001:0.2,0.2:0.01:1; %变步长绘制>>

9、;y1=sin(1./t1);>>t2=-pi:0.01:-2,-2:0.001:-1,-1:0.01:1,1:0.001:2,2:0.01:pi; %变步长绘制>>y2=sin(tan(t2)-tan(sin(t2);>>subplot(2,1,1),plot(t1,y1) %绘制y1函数>>subplot(2,1,2),plot(t2,y2) %绘制y2函数>> 运行结果Figure 2 二维曲线绘制结果1-6 试绘制出二元函数的三维图和三视图。>>clear all clear,clc;>>xx=-2:0

10、.1:-1.2,-1.1:0.02:-0.9,-0.8:0.1:0.8,0.9:0.02:1.1,1.2:0.1:2;>>yy=-1:0.1:-0.2,-0.1:0.02:0.1,0.2:0.1:1;x,y=meshgrid(xx,yy);%变步长绘制二元函数曲线，使曲线显示更为合理，优化显示奇异点附近曲面>>z=1./(sqrt(1-x).2+y.2)+1./(sqrt(1+x).2+y.2);>>subplot(2,2,1),surf(x,y,z),shading flat,zlim(0,15),title('图1 三维视图')>&

11、gt;subplot(2,2,2),surf(x,y,z),shading flat,zlim(0,15),view(0,90),title('图2 俯视图')>>subplot(2,2,3),surf(x,y,z),shading flat,zlim(0,15),view(-90,0),title('图3 左视图')>>subplot(2,2,4),surf(x,y,z),shading flat,zlim(0,15),view(0,0),title('图4 正视图')>> 运行结果Figure 3 曲面绘制

12、结果三维显示图1-7 试求出如下极限。（1）； （2）； （3）。>>syms x;syms y>>f1=(3x+9x)(1/x);>>L1=limit(f1,x,inf)%（1）的极限求取算法>>f2=x*y/(sqrt(x*y+1)-1);>>L2=limit(limit(f2,x,0),y,0) %（2）的极限求取算法>>f3=(1-cos(x2+y2)/(x2+y2)*exp(x2+y2); >>L3=limit(limit(f3,x,0),y,0) %（3）的极限求取算法>> 运行结果（1

13、）；L1 = 9 （2）；L2 = 2 （3）。L3 = 01-8 已知参数方程，试求出和。>>clear all,clear,clc;>>syms t;>>x=log(cos(t);>>y=cos(t)-t*sin(t);>>f1=diff(y,t)/diff(x,t)%dy/dx的运算结果>>f2=diff(y,t,2)/diff(x,t,2) %d2y/dx2的运算结果>>f2=subs(f2,t,sym(pi)/3)%d2y/dx2|t=pi/3的运算结果>>latex(f1);latex(

14、f2); %latex形式显示结果>> 运行结果f1 = (cos(t)*(2*sin(t) + t*cos(t)/sin(t)%运算结果为：= f2 = (3*cos(t) - t*sin(t)/(sin(t)2/cos(t)2 + 1) %运算结果为： f2 = 3/8 - (pi*3(1/2)/24%运算结果为：ans =1-9 假设，试求。>>syms x;syms y;syms t;>>f=int(exp(-t2),t,0,x*y);>>F=(x/y)*diff(f,x,2)-2*diff(diff(f,x),y)+diff(f,y,2

15、);>>latex(F);>> 运行结果F = (2*x2*y2)/exp(x2*y2) - (2*x3*y)/exp(x2*y2) - 2/exp(x2*y2)%的运算结果为：1-10 试求出下面的极限。 （1）；（2）。>>clear all,clear,clc;>>syms n;syms k;>>s1=limit(symsum(1/(2*k)2-1),k,1,n),n,inf)%(1)的计算算法>>s2=limit(n*symsum(1/(n2+k*sym(pi),k,1,n),n,inf)%(2)的计算算法>

16、> 运行结果s1 = 1/2%(1)运行结果 s2 = 1%(2)运行结果1-11 试求出以下的曲线积分。 （1），为曲线， 。 （2），其中为正向上半椭圆。>>syms t1;>>syms a positive;>>x1=a*(cos(t1)+t1*sin(t1);>>y1=a*(sin(t1)-t1*cos(t1);>>I1=int(x12+y12)*sqrt(diff(x1,t1)2+diff(y1,t1)2),t1,0,2*pi)%(1)的计算算法>>syms x y t2;>>syms a b

17、c positive;>>x=c*cos(t2)/a;>>y=c*sin(t2)/b;>>F1=y*x3+exp(y);>>F2=x*y3+x*exp(y)-2*y;>>ds=diff(x,t2);diff(y,t2);>>I2=int(F1 F2*ds,t2,0,pi) %(1)的计算算法>>latex(I1),latex(I2)>> 运行结果I1 = 2*pi2*a3*(2*pi2 + 1)%的曲线积分 I2 = -(2*c*(15*b4 - 2*c4)/(15*a*b4)%的曲线积分 ans

18、= %(1)的计算结果的LATEX显示ans = %(2)的计算结果的LATEX显示1-12 试求出Vandermonde矩阵的行列式，并以最简的形式显示结果。>>clear all,clear,clc>>syms a b c d e;>>A=a4 a3 a2 a 1;b4 b3 b2 b 1;c4 c3 c2 c 1;d4 d3 d2 d 1;e4 e3 e2 e 1;>>D=simple(det(A)>> 运行结果D = (a - b)*(a - c)*(a - d)*(b - c)*(a - e)*(b - d)*(b - e)

19、*(c - d)*(c - e)*(d - e) 1-13 试对矩阵进行Jordan变换，并得出变换矩阵。>>A=-2 0.5 -0.5 0.5;0 -1.5 0.5 -0.5;2 0.5 -4.5 0.5;2 1 -2 -2;>>V,J=jordan(A)>> 运行结果V = 0 0.5000 0.5000 -0.2500 0 0 0.5000 1.0000 0.2500 0.5000 0.5000 -0.2500 0.2500 0.5000 1.0000 -0.2500J = -4 0 0 0 0 -2 1 0 0 0 -2 1 0 0 0 -21-14

20、 试用数值方法和解析方法求取下面的Sylvester方程，并验证得出的结果。>>A=3 -6 -4 0 5;1 4 2 -2 4;-6 3 -6 7 3;-13 10 0 -11 0;0 4 0 3 4;>>B=3 -2 1; -2 -9 2;-2 -1 9;>>C=-2 1 -1;4 1 2;5 -6 1;6 -4 4;-6 6 -3;>>X1=lyap(A,B,C),e1=norm(A*X1+X1*B+C) %X的数值解法及误差>>X2=lyap(sym(A),B,C),e2=norm(double(A*X2+X2*B+C) %X

21、的解析解法及误差>> 运行结果%(1)Sylvester方程的数值解法X1 = 4.0569 14.5128 -1.5653 -0.0356 -25.0743 2.7408 -9.4886 -25.9323 4.4177 -2.6969 -21.6450 2.8851 -7.7229 -31.9100 3.7634e1 =%X的数值解误差 4.3904e-013%(2)Sylvester方程的解析解法X2 = 434641749950/107136516451, 4664546747350/321409549353, -503105815912/321409549353 -3809

22、507498/107136516451, -8059112319373/321409549353, 880921527508/321409549353 -1016580400173/107136516451, -8334897743767/321409549353, 1419901706449/321409549353 -288938859984/107136516451, -6956912657222/321409549353, 927293592476/321409549353 -827401644798/107136516451, -10256166034813/321409549353

23、, 1209595497577/321409549353e2 =%X的解析解误差 01-15 假设已知矩阵如下，试求出，。>>syms t;>>A=-4.5 0 0.5 -1.5;-0.5 -4 0.5 -0.5;1.5 1 -2.5 1.5;0 -1 -1 -3;>>F1=exp(A.*t)>>F2=sin(A.*t)>>F3=exp(A.*t)*sin(A.2.*exp(A.*t).*t)>> 运行结果F1 =%的计算结果 1/exp(9*t)/2), 1, exp(t/2), 1/exp(3*t)/2) 1/exp(

24、t/2), 1/exp(4*t), exp(t/2), 1/exp(t/2) exp(3*t)/2), exp(t), 1/exp(5*t)/2), exp(3*t)/2) 1, 1/exp(t), 1/exp(t), 1/exp(3*t) F2 =%的计算结果 -sin(9*t)/2), 0, sin(t/2), -sin(3*t)/2)sin(t/2), -sin(4*t),sin(t/2), sin(t/2)sin(3*t)/2), sin(t), -sin(5*t)/2), sin(3*t)/2) 0, -sin(t), -sin(t), -sin(3*t) F3 =%的计算结果 si

25、n(t/(4*exp(t/2) + exp(t/2)*sin(9*t*exp(3*t)/2)/4) + sin(81*t)/(4*exp(9*t)/2)/exp(9*t)/2), sin(16*t)/exp(4*t) + sin(t/exp(t)/exp(3*t)/2) + sin(t*exp(t)*exp(t/2),sin(t*exp(t/2)/4) + sin(t/exp(t)/exp(3*t)/2) + sin(t*exp(t/2)/4)/exp(9*t)/2) + exp(t/2)*sin(25*t)/(4*exp(5*t)/2),sin(t/(4*exp(t/2) + sin(9*t

26、)/exp(3*t)/exp(3*t)/2) + exp(t/2)*sin(9*t*exp(3*t)/2)/4) + sin(9*t)/(4*exp(3*t)/2)/exp(9*t)/2) sin(t/(4*exp(t/2)/exp(4*t) + exp(t/2)*sin(9*t*exp(3*t)/2)/4) + sin(81*t)/(4*exp(9*t)/2)/exp(t/2), sin(t/exp(t)/exp(t/2) + sin(16*t)/exp(4*t)/exp(4*t) + sin(t*exp(t)*exp(t/2), sin(t/exp(t)/exp(t/2) + sin(t*

27、exp(t/2)/4)/exp(t/2) + sin(t*exp(t/2)/4)/exp(4*t) + exp(t/2)*sin(25*t)/(4*exp(5*t)/2), sin(t/(4*exp(t/2)/exp(4*t) + sin(9*t)/exp(3*t)/exp(t/2) + sin(9*t)/(4*exp(3*t)/2)/exp(t/2) + exp(t/2)*sin(9*t*exp(3*t)/2)/4) sin(9*t*exp(3*t)/2)/4)/exp(5*t)/2) + exp(3*t)/2)*sin(81*t)/(4*exp(9*t)/2) + exp(t)*sin(t

28、/(4*exp(t/2), exp(3*t)/2)*sin(t/exp(t) + sin(t*exp(t)/exp(5*t)/2) + exp(t)*sin(16*t)/exp(4*t), exp(3*t)/2)*sin(t/exp(t) + exp(3*t)/2)*sin(t*exp(t/2)/4) + sin(25*t)/(4*exp(5*t)/2)/exp(5*t)/2) + exp(t)*sin(t*exp(t/2)/4), exp(3*t)/2)*sin(9*t)/exp(3*t) + exp(3*t)/2)*sin(9*t)/(4*exp(3*t)/2) + sin(9*t*exp

29、(3*t)/2)/4)/exp(5*t)/2) + exp(t)*sin(t/(4*exp(t/2) sin(81*t)/(4*exp(9*t)/2) + sin(t/(4*exp(t/2)/exp(t) + sin(9*t*exp(3*t)/2)/4)/exp(t), sin(t/exp(t)/exp(3*t) + sin(16*t)/exp(4*t)/exp(t) + sin(t*exp(t)/exp(t), sin(t*exp(t/2)/4) + sin(t/exp(t)/exp(3*t) + sin(t*exp(t/2)/4)/exp(t) + sin(25*t)/(4*exp(5*t

30、)/2)/exp(t), sin(9*t)/(4*exp(3*t)/2) + sin(t/(4*exp(t/2)/exp(t) + sin(9*t)/exp(3*t)/exp(3*t) + sin(9*t*exp(3*t)/2)/4)/exp(t)%的计算结果LATEX显示F3 =第二部分数学问题求解与数据处理（4学时）主要内容：掌握代数方程与最优化问题、微分方程问题、数据处理问题的MATLAB求解方法。2-1 对下列的函数进行Laplace变换。（1）；（2）；（3）。>>syms t;>>syms a;>>fa=sin(a*t)/t;>>fb

31、=t5*sin(a*t);>>fc=t8*cos(a*t);>>Fa=laplace(fa)>>Fb=laplace(fb)>>Fc=laplace(fc)>> 运行结果Fa =atan(a/s) %的LATEX显示结果为 Fb =(720*a*s)/(a2 + s2)4 - (3840*a*s3)/(a2 + s2)5 + (3840*a*s5)/(a2 + s2)6 %的LATEX显示结果为 Fc = (362880*s)/(a2 + s2)5 - (4838400*s3)/(a2 + s2)6 + (17418240*s5)/(

32、a2 + s2)7 - (23224320*s7)/(a2 + s2)8 + (10321920*s9)/(a2 + s2)9%的LATEX显示结果为2-2 对下面的式进行Laplace反变换。（1）；（2）；（3）。>>syms s t;>>syms a b;>>Fa=1/(sqrt(s2)*(s2-a2)*(s+b);>>Fb=sqrt(s-a)-sqrt(s-b);>>Fc=log(s-a)/(s-b);>>fa=ilaplace(Fa,s,t)>>fb=ilaplace(Fb,s,t)>>f

33、c=ilaplace(Fc,s,t)>>latex(fa),latex(fb),latex(fc)>> 运行结果fa = -ilaplace(1/(b + s)*(a2 - s2)*(s2)(1/2), s, t) %LATEX结果 fb = ilaplace(s - a)(1/2), s, t) - ilaplace(s - b)(1/2), s, t) %LATEX结果 fc = exp(b*t)/t - exp(a*t)/t%的LATEX显示结果为2-3 试求出下面函数的Fourier变换，对得出的结果再进行Fourier反变换，观察是否能得出原来函数。（1）；（

34、2）。>>syms w x t;>>fx=x2*(3*pi-2*abs(x); >>ft=t2*(t-2*pi)2;>>Fx=fourier(fx),Ft=fourier(ft) %Fourier正变换>>fx1=ifourier(Fx,w,x),ft1=ifourier(Ft,w,t)%Fourier反变换>>latex(fx1),latex(ft1)>> 运行结果Fx = - 24/w4 - 6*pi2*dirac(w, 2) %Fourier正变换结果 Ft = 2*pi*dirac(w, 4) + pi

35、2*dirac(w, 3)*8*i - 8*pi3*dirac(w, 2) %Fourier正变换结果 fx1 = (6*pi2*x2 - 4*pi*x3*(2*heaviside(x) - 1)/(2*pi)%(1)执行正变换后进行Fourier反变换结果 ft1 = (2*pi*t4 - 8*pi2*t3 + 8*pi3*t2)/(2*pi) %(2) 执行正变换后进行Fourier反变换结果2-4 请将下述时域序列函数进行Z变换，并对结果进行反变换检验。（1）；（2）；（3）。>>syms a T k;>>fa=cos(a*T*k);>>fb=(k*T

36、)2*exp(-a*k*T);>>fc=(a*k*T-1+exp(-a*k*T)/a;>>Za=ztrans(fa);Zb=ztrans(fb);Zc=ztrans(fc);>>za=simple(iztrans(Za),z,k),latex(za)>>zb=simple(iztrans(Zb),z,k),latex(zb)>>zc=simple(iztrans(Zc),z,k),latex(zc)>> 运行结果Za = (z*(z - cos(T*a)/(z2 - 2*cos(T*a)*z + 1) %(1)Z变换结果Z

37、b = (T2*z*exp(T*a)*(z*exp(T*a) + 1)/(z*exp(T*a) - 1)3%(2)Z变换Zc = (T*z)/(z - 1)2 + z/(a*(z - exp(-T*a) - z/(a*(z - 1) %(3)Z变换 za = cos(T*a*k) %(1)Z变换后反变换结果zb =T2*k2*exp(-T*a)k%(2)Z变换后反变换结果zc = T*k + (exp(-T*a)k - 1)/a%(3)Z变换后反变换结果2-5 用数值求解函数求解下述一元和二元方程的根，并对得出的结果进行检验。（1）；（2）。>>X=solve('exp(-

38、(x+1)2+pi/2)*sin(5*x+2)=0') %(1)求解>>eval('exp(-(X+1).2+pi/2)*sin(5*X+2)')>>ezplot('exp(-(x+1)2+pi/2)*sin(5*x+2)')>>syms x; %(2)求解>> Y=solve('(x2+y2+x*y)*exp(-x2-y2-x*y)=0','y')>>subs('(x2+y2+x*y)*exp(-x2-y2-x*y)=0','y'

39、,Y)>>ezplot('(x2+y2+x*y)*exp(-x2-y2-x*y)')>> 运行结果X = -2/5%(1)利用solve方法解求的答案ans = 0%(1)计算结果的检验 Y = %(2)利用solve方法求得的y关于x的表达式 - x/2 + (3(1/2)*x*i)/2 - x/2 - (3(1/2)*x*i)/2ans = %(2)计算结果带入原方程后进行检校的结果为0 exp(x*(x/2 - (3(1/2)*x*i)/2) - (x/2 - (3(1/2)*x*i)/2)2 - x2)*(x2 + (x/2 - (3(1/2)*

40、x*i)/2)2 - x*(x/2 - (3(1/2)*x*i)/2) = 0exp(x*(x/2 + (3(1/2)*x*i)/2) - (x/2 + (3(1/2)*x*i)/2)2 - x2)*(x2 + (x/2 + (3(1/2)*x*i)/2)2 - x*(x/2 + (3(1/2)*x*i)/2) = 0%(2)检校结果计算过程的latex表示在(1)的计算中，应用数值解法solve只能求取出方程的一个根，而由图解法可知方程还应有其他解如-2.9134，-2.285，-1.6567，-1.0283，0.2284及0.8567等。这是算法不允许选择初始搜索点造成的。故对(1)的算法

41、进行数值精确化查找，利用fsolve函数通过以下过程得到所需结果：>> y=(x)exp(-(x+1)2+pi/2).*sin(5*x+2);>> x1,f=fsolve(y,-2.9134)>> x2,f=fsolve(y,-2.285)>> x3,f=fsolve(y,-1.6567)>> x4,f=fsolve(y,-1.0283)>> x5,f=fsolve(y,0.2284)>> x6,f=fsolve(y,0.8567)>> x,f=fsolve(y,-0.4)方程解x1x2x3x4x5

42、x6x图解法-2.9134-2.285-1.6567-1.02830.22840.8567-0.4精确解-2.9133-2.2850-1.6566-1.02830.22830.8566-0.4000误差3.7547e-08-2.3410e-087.6555e-16-4.7028e-108.6807e-08-1.1271e-0802-6 试求出使得取得极小值的值。>> syms x c; %利用解析解法和图解法联合求解函数的取得极小值时的c值 >> y=int(exp(x)-c*x)2,x,0,1);>> y1=diff(y,c);ezplot(y1)%求取函

43、数的一阶导数>> t0=solve(y1);ezplot(y)%求出一阶导数为零的点 >> y2=diff(y1);b=subs(y2,c,t0) %验证二阶导数为正>> 运行结果t0 = 3%一阶导数为零的c值为3，即函数取得极小值时c为3b = 2/3%该点处二阶导数的验证结果为2/3>0，恰为极小值点图1 原函数图2 原函数一阶导数2-7 试求解下面的非线性规划问题。 首先编辑函数opt_con1记录Ceq和C矩阵：function c,ceq = opt_con1(x)ceq=;c=x(1)+x(2);x(1)*x(2)-x(1)-x(2)+1

44、.5;-10-x(1)*x(2);end以下为在Matlab命令行中执行的命令：>>y=(x)exp(x(1)*(4*x(1)2+2*x(2)2+4*x(1)*x(2)+2*x(2)+1); %录入目标函数>>ff=optimset;ff.LargeScale='off' %设置控制选项 >>ff.TolFun=1e-30;ff.TolX=1e-15;ff.TolCon=1e-20;>>x0=0;0;xm=-10;-10;xM=10;10;A=;B=;Aeq=;Beq=;%设置可行解条件矩阵>>x,f_opt,c,d=

45、fmincon(y,x0,A,B,Aeq,Beq,xm,xM,opt_con1,ff) %用fmincon求解 >> 运行结果x =%x(1)和x(2)的最优结果 1.1825 -1.7398f_opt =%目标函数的最优值fopt 3.0608c =%表示目标函数非线性规划求解成功 2*;第38页，共38页d = %求解过程中相关参量 iterations: 23funcCount: 113constrviolation: 0stepsize: 1.5267e-15algorithm: 'interior-point'firstorderopt: 8.3962e-

46、08cgiterations: 27message: 'Local minimum possible. Constraints satisfied.2-8 求解下面的整数线性规划问题。 >>f=-592;381;273;55;48;37;23; %前面的“-”号代表求取最大值>>A=3534,2356,1767,589,528,451,304;>>b=119567;>>lb=zeros(7,1);>>intcon=1,2,3,4,5,6,7;>>x,fval=intlinprog(f,intcon,A,b,lb)&

47、gt;> 运行结果LP: Optimal objective value is -20012.250000. Cut Generation: Applied 2 gomory cuts. Lower bound is -19979.000000. Relative gap is 0.00%. Optimal solution found.Intlinprog stopped at the root node because the objective value is within a gap tolerance of the optimal value; options.TolGapA

48、bs = 0 (the default value). The intcon variables are integer within tolerance,options.TolInteger = 1e-05 (the default value).x =%取得最大值时x的相关取值 32.0000 2.0000 1.0000 0 0 0 0fval =%当x为上述向量时取得的函数的最大值为19979 -199792-9 试求出微分方程的解析解通解，并求出满足边界条件的解析解。>>syms x y1 y2;>>y1=dsolve('D2y+(2-1/x)*Dy+(

49、1-1/x)*y=x2*exp(-5*x)');>>y2=dsolve('D2y+(2-1/x)*Dy+(1-1/x)*y=x2*exp(-5*x)','y(1)=pi','y(pi)=1');>>y1,y2>> 运行结果y1 = C2/exp(t*(x - 1)/x) + C3/exp(t) - x2/(exp(5*x)*(1/x - 1)%解析通解latex表示 y2 = C13/exp(t) - x2/(exp(5*x)*(1/x - 1) - (exp(x - 1)/x)*(C13*x - C1

50、3 + exp(1)*piy(pi) + (x3*exp(1)/exp(5*x) - x*exp(1)*piy(pi)/(exp(1)*exp(t*(x - 1)/x)*(x - 1)%满足边界条件特解latex表示2-10 试求出下面微分方程的通解。（1）；（2）。>>syms x y1 y2 t;>>y1=dsolve('D2x+2*t*Dx+t2*x=t+1');>>y2=dsolve('Dy+2*x*Dy=x*exp(-x2)');>>y1,y2>> 运行结果y1 = (exp(t2/2 - t

51、)/2 - (2(1/2)*pi(1/2)*erf(2(1/2)*t*i)/2 - (2(1/2)*i)/2)*i)/(2*exp(1/2)/exp(t*(t - 2)/2) + C15/exp(t*(t - 2)/2) + C16/exp(t*(t + 2)/2) - exp(t2/2 + t)/(2*exp(t*(t + 2)/2) %(1)通解的latex表示y2 = C18 + (t*x)/(exp(x2)*(2*x + 1) %(2)通解的latex表示2-11 考虑著名的化学反应方程组，选定，且，绘制仿真结果的三维相轨迹，并得出其在x-y平面上的投影。在实际求解中建议将作为附加参数

52、，同样的方程若设，时，绘制出状态变量的二维图和三维图。>>f=(t,x,a,b,c)-x(2)-x(3);x(1)+a*x(2);b+(x(1)-c)*x(3);>>options=odeset;options.RelTol=1e-7;%修改精度>>t_final=100;x0=0;0;0;>>a=0.2;b=0.2;c=5.7;>>t,x=ode45(f,0,t_final,x0,options,a,b,c);>>subplot(2,2,1),plot3(x(:,1),x(:,2),x(:,3),title('图1.1 三维相轨迹(a=b=0.2,c=5.7)');>>subplot(2,2,2),plot3(x(:,1),x(:,2),x(:,3),view(0,90),title('图1.2 x-y平面投影(a=b=0.2,c=5.7)');>>a=0.2;b=0.5;c=10;>>t,x=ode45(f,0,t_final,x0,options,a,b,c);>>subplot(2,2,3),plot3(x(:,1),x(:,2),