M文档下载】ATLAB)课后实验答案.pdf
61页
实验一 MATLAB运算基础2 l + 2z-0.45 51 . 先 求下列表达式的值，然后显示MATLAB工作空间的使用情况并保存全部变量/八 2 s i n 8 5°⑴(2) . = ；l n (x + Jl + x ?) ,其中 x =0 . 3 a - - 0 . 3 a(3) Z 3=——Y—s i n (a + 0.3) + l n；, a = — 3.0, — 2.9 ,…， 2.9 ,3.0? 0
2 )将矩阵C的右下角3 X 2子矩阵赋给Do(3 )查看MATLAB工作空间的使用情况解：.运算结果：E=(reshape(1:1:25,5,5)),;F=[3 0 16;17 -6 9;0 23 -4;9 7 0;413 11];C= E*FH=C(3:5,2:3)C =93 15077258 335237423 520397588 705557753 890717 H =520 397705 557890 7174 .完成下列操作：( 1 )求［ 100, 999］之间能被21整除的数的个数 2 )建立一个字符串向量，删除其中的大写字母解：(1 )结果: m=100:999;n=find(mod(m,21 )==0);length(n)ans =43(2).建立一个字符串向量 例如：ch='ABC123d4e56Fg9';则要求结果是:ch='ABC123d4e56Fg9';k=find(ch>='A'&ch<='Z');ch(k)=[]ch =123d4e56g9 实验二MATLAB矩阵分析与处理1 .设有分块矩阵4^F3x302x3 §2x2，其 中E、R、0、S分别为单位矩阵、随机矩阵、零矩阵和对角阵，试通过数值计算验证R + RSS2解：M文件如下;E0输出结果:S =1 00 2A1.0000000.53830.442701.000000.99610.1067001.00000.07820.96190001.0000000002.0000 a =1.0000001.07671.328001.000001.99230.3200001.00000.15642.88570001.0000000004.0000ans =0000000000000000000000000S a n s ,所以 A?E R + RSO S22 .产 生5阶希尔伯特矩阵H和5阶帕斯卡矩阵P ,且求其行列式的值Hh和Hp以及它们的条件数Th和T p ,判断哪个矩阵性能更好。
为什么？解：M文件如下：输出结果:H =1.00000.50000.33330.25000.20000.50000.33330.25000.20000.16670.33330.25000.20000.16670.1429 因为它们的条件数Th»Tp,所以pascal矩阵性能更好0.25000.20000.20000.16670.16670.14290.14290.12500.12500.1111P =1 111 11 234 51 3610 151 41020 351 51535 70Hh =3.7493e-012Hp =1Th =4.7661 e+005Tp =8.5175e+0033 .建立一个5 X 5 矩阵，求它的行列式值、迹、秩和范数解： M 文件如下: ■ Editor - Untitled2»输出结果为:17241815235714164613202210121921311182529d =5070000t =65c1 =6.8500c2 = 5.4618cinf =6.85004 . 已知--29 6 18一A = 20 5 12_ -8 8 5求 A 的特征值及特征向量，并分析其数学意义。
解：M 文件如图:输出结果为:v =0.71300.28030.2733-0.6084-0.78670.87250.34870.55010.4050 D =-25.3169 0 00 -10.5182 00 0 16.8351数学意义：V 的 3 个列向量是A 的特征向量，D 的主对角线上3 个是A 的特征值，特别的,V 的 3 个列向量分别是D 的 3 个特征值的特征向量5 . 下面是一个线性方程组：1-411-311-21-51-6-41-5-31-4(1 )求方程的解2 )将方程右边向量元素b3改为0.53再求解，并比较b3的变化和解的相对变化3 )计算系数矩阵A 的条件数并分析结论解： M 文件如下：目回区I■ E d ito r - U n titled 2 »输出结果: X =1.20000.60000.6000X2 =1.20000.60000.6000c =1.3533e+003由结果，X 和 X 2的 值 •样 ，这表示b 的微小变化对方程解也影响较小，而 A 的条件数算得较小，所以数值稳定性较好，A 是较好的矩阵6 .建立A 矩阵，试比较sqrtm(A)和 sqrt(A),分析它们的区别。
解：M 文件如下：运行结果有:A =1620965818125 b1 =3.8891-0.11023.21033.29172.14360.36980.38552.07601.7305b2 =4.00002.44954.24264.47212.23613.46413.00002.82842.2361b =16.00006.000018.000020.00005.000012.00009.00008.00005.0000分析结果知：sqrtm(A)是类似A 的数值平方根(这可山b1*b1=A的结果看出) ，而 sqrt(A)则是对A 中的每个元素开根号，两则区别就在于此实验三选择结构程序设计1 .求分段函数的值x2 + x-6 x < OJLx -3y = ' x2 -5x + 6 0< x< 5且x 丰 2及x 丰 3X2 - X - 1 其他用 if语句实现，分别输出x=-5.0,-3.0,1,0,2.0,2.5,3.0,5.0时的y 值解：M 文件如下： 运算结果有:f(-5)y =14» f(-3 )y =ii» f(Dy =2 » f(2)y =i» f(2.5)y =-0.2500» f(3 )y =5» f(5 )y =192 .输入一个百分制成绩，要求输出成绩等级A、B、C、D、E o其 中9 0分~100分 为A,8 0分~8 9分为B, 7 9分~7 9分为C, 6 0分〜6 9分为D, 6 0分以下为E。
要求：(1 )分别用if语句和switch语句实现2 )输入百分制成绩后要判断该成绩的合理性，对不合理的成绩应输出出错信息解：M文件如下 试算结果:score=88grade =Bscore= 123错误：输入的成绩不是百分制成绩3 .硅谷公司员工的工资计算方法如下:(1 )工作时数超过120小时者，超过部分加发15%2 )工作时数低于60小时者，扣发700元3 )其余按每小时84元计发试编程按输入的工号和该号员工的工时数，计算应发工资 解：M文件下4 .设计程序，完成两位数的加、减、乘、除四则运算，即产生两个两位随机整数，再输入一个运算符号，做相应的运算，并显示相应的结果解：M文件如下； 运算结果例:a =38b =33输入一个运算符: 八false a =92b =40输入一个运算符: +c =1325 . 建 立 5 X 6 矩阵，要求输出矩阵第n 行元素当 n 值超过矩阵的行数时，自动转为输出矩阵最后一行元素，并给出出错信息解：M 文件如下：运算结果如卜: 输入一个 5 彳 亍 6 歹 I 」 矩阵 A=[l 2 3 4 5 5;2 3 4 5 7 6;2 2 2 2 2 3;11 2 3 9 7 3;2 3 4 5 6 7]输入一正整数n=411 2 3 9 7 3输 入 ♦ 个 5 行 6 歹 U 矩阵 A=[l 2 3 45 5;2 345 7 6;2 222 2 3;11 2 39 7 3;2 3 45 67]输入一正整数n=62 3 4 5 6 7ans =Error using ==> dispToo many input arguments.实验四循环结构程序设计1 . 根据［■ =《 + / +J + … + 4, 求 n 的近似值。
当 n 分别取100、1000、10000时，结果是多少？要求：分别用循环结构和向量运算（ 使用sum 函数）来实现解：M 文件如下: 运行结果如下：K » %循环结构计算pi值y=o；n=input('n=*);for i=1 :ny=y+1/i/i；endpi=sqrt(6*y)n=100Pi = 3.1321n=1000Pi =3.1406n=10000pi =3.1415%向量方法计算P i值n=inputfn=');i=1./(1:n).A2;s=sum(i);pi=sqrt(6*s)n=100Pi =3.1321n=1000Pi =3.1406n=10000Pi =3.14152 .根据y = l + ' + L + ・・・ +」一 ，求:3 5 2/2-1 ⑴ y v 3时的最大n值⑵ 与 ⑴ 的n值对应的y值解：M一文件如下：运行结果如下：K » y=0;n=0;while y<3n=n+1;y=y+1/(2*n-1);endynif y>3n=n-1;endny =3.0033 n =57n =563 . 考虑以下迭代公式:a其中a 、b为正的学数 1 ) 编写程序求迭代的结果，迭代的终止条件为| X n + 「 X n | W 1 0 - 5 , 迭代初值X °= 1 . 0 , 迭代次数不超过5 0 0 次。
2 ) 如果迭代过程收敛于r , 那 么 r 的准确_b 值+ J b2是 + 4 .,,当( a , b ) 的值取( 1 , 1 ) 、( 8 , 3 ) 、( 1 0 , 0 . 1 ) 时，分别对迭代结果和准确值进行比较解：M文件如下： 运算结果如下;请输入正数a=1请输入正数b=1x =0.6180r =0.6180-4.70160.6180-1.6180s = -0.0000 -2.2361请输入正数a=8请输入正数b=3x =1.70161.7016-1.61801.7016-4.7016s =0.0-6.4031请输入正数a=10请输入正数b=0.1x =3.11273.1127-4.70163.1127-3.2127s =-0.0000-6.3254 4 . 已知/ i = 1 〃 = 1力 =0 " = 2力 =1 〃 = 3= 工1一 2 / 一2+ 力. 3 〃> 3求 fM o o 中：(1 )最大值、最小值、各数之和2 )正数、零、负数的个数解：M一文件以下是运算结果：max(f)=437763282635min(f)=-899412113528sum(f)=-742745601951c1=49c2=2c3=49 5 .若两个连续自然数的乘积减1 是素数，则称这两个边疆自然数是亲密数对，该素数是亲密素数。
例如， 2 X 3 7 = 5 ,由于5 是素数， 所以2 和 3 是亲密数， 5 是亲密素数 求[2,50]区间内：(1 )亲密数对的对数2 )与上述亲密数对对应的所有亲密素数之和解：M 文件：运算结果为:j =29S =23615 实验五函数文件1 . 定义一个函数文件，求给定复数的指数、对数、正弦和余弦，并在命令文件中调用该函数文件解：M 文件如下：函数fushu.M文件:function [e,l,s,c] = fushu(z)% fu sh u 复数的指数，对数，正弦，余弦的计算%e 复数的指数函数值%1 复数的对数函数值%s 复数的正弦函数值%c 复数的余弦函数值e=exp(z);l=log(z);s=sin(z);c=cos(z);命令文件M：z=input( ' 请输入一个复数z=');[a,b,c,d]=fushu(z)运算结果如下：z=input 请输入一个复数z=');[a,b,c,d]=fushu(z)请输入一个复数z=1+ia =1.4687 + 2.2874ib =c =0.3466 + 0.7854i 1.2985 + 0.6350id =0.8337 - 0.9889i2. 一物理系统可用下列方程组来表示:叫mx scions。
0-sin < 9 0%00cos0a2町g0m2一sin0000一cos1.m2S.从键盘输入3 、m2和 的值，求 a 、a2、M 和 N2的值其 中 g 取 9 .8 ,输入8 时以角度为单位要求：定义一个求解线性方程组AX =B的函数文件，然后在命令文件中调用该函数文件解： M 文件函数fc.M文件：function X= fc(A,B)%fc fc 是求解线性方程的函数%A A 是未知矩阵的系数矩阵X=A\B；命令M 文件：cic;m1=input( ' 输入 m1 = );m2=input( 输入 m2=);theta=input 输入 theta=');x=theta*pi/180;g=9.8;A=[m1*cos(x) -ml -sin(x) 0m1*sin(x) 0 cos(x) 00 m2 -sin(x) 00 0 -cos(x) 1];B=[0;m1*g;0;m2*g];X=fc(A,B)运算结果：输入m1=1输入m2=1输入 theta=30X =7.8400 3.39486.789615.68003. 一个自然数是素数，且它的数字位置经过任意对换后仍为素数例 如 13是绝对素数。
试求所有两位绝对素数要求：定义一个判断素数的函数文件解：M 文件：函数prime.m文件function [p] = prime(p)% 输入p 的范围，找出其中的素数m=p(length(p));for i=2:sqrt(m)n=find(rem(p,i)==O&p-=i);P(n)=[]；% 将 p 中能被i 整除，而却不等于i 的元素，即下标为n 的元素剔除，其余的即为素数endP；命令文件：cic;p=10:99;p=prime(p); %找 出 10到 9 9 内的所有素数p=10*rem(p,10)+(p-rem(p,10))/10;% 将 p 素数矩阵每个元素个位十位调换顺序p=prime(p)% 再对对换后的素数矩阵找出所有的素数运算结果：p =1131711373173797794. 设 /( %) =- - - - - - - - - - -+- - - - - -- - - - - - -( 尤 -2)2+0.1 *-3 )4 + 0 .0 1, 编写一个MATLAB函数文件仅.m ,使得调用f(x)时，x 可用矩阵代入，得出的f(x)为同阶矩阵解：函数仅.m 文件：function f= fx(x)%fx 仅求算x 矩阵下的f(x)的函数值A=0.1+(x-2).A2; B=0.01+(x-3).A4;f=1./A+1./B;命令文件：cic;x=input。
输入矩阵 x=');f=fx(x)运算结果：» x=input( ' 输入矩阵 x=')；f=fx(x)输入矩阵x=[7 2;12 5]f =0.0437 10.99010.0101 0.1724/(30) + /(20)(1) Sf(n)=n+10ln(n2+ 5)B t,求 y 的值2 )当 f(n)=1X2+2X3+3X4+...+nX(n+1)时，求 y 的值解：⑴函数f.m 文件：function f=f(x)f=x+10*log(xA2+5);命令文件：clc;nl=input(* nl=');n2=input(* n2=');n3=input('n3=');yl=f(nl);y2=f(n2);y3=f(n3);y=yl/(y2+y3)运算结果如下:n1=40n2=30n3=20 y =0.6390(2).函数g.m 文件function s= g(n)for i=l:ng(i)=i* (i + 1);ends=sum(g);命令文件：clc;nl=input('nl=');n2=input(* n2=');n3=input(* n3=');yi=g (ni)；y2=g(n2);y3=g (n3);y=yi/(y2+y3)运算结果如下:n1=40n2=30n3=20y =1.7662实验六高层绘图操作1. 设 y = 0.5 +3 sin xl + x2c o s x ,在 x=o~2n区间取101点，绘制函数的曲线。
解：M 文件如下：clc;x=linspace(0,2*pi r 101); y=(0.5 + 3*sin(x) ./ (1+x.A2));plot(xz y)运行结果有：2 .已知＞ 1=¥, y2=COS(2t), y3=)T X y 2 ,完成下列操作：(1 )在同一坐标系下用不同的颜色和线型绘制三条曲线2 )以子图形式绘制三条曲线3 )分别用条形图、阶梯图、杆图和填充图绘制三条曲线解 ：(1) M文件：cic;x=-pi:pi/100:pi;y1=x.A2;y2=cos(2*x);y3=y1 .*y2;plot(x,y1,,b-',x,y2,'r:',x,y3,,k-')运行结果：(2) M文件：clc;x=-pi:pi/100:pi;yl=x. 人 2;y2=cos(2*x);y3=yl.*y2;subplot(1,3,1);plot (xz ylz fb-');title('yl=xA2 ');subplot(1,3,2);plot(x,y2,'r:');title(T y2 = cos(2x)');subplot(1, 3,3);plot(x,y3, fk--');title(1y3=yl*y21);. 运行结果:(3) M文件: clc;x=-pi:pi/100:pi;yl=x.A2;y2=cos(2 *x);y3=yl.*y2;subplot (2,2,1);plot (x,yl, 'b-',x,y2,'r :',x,y3, 'k--1);subplot (2,2,2);bar(xf yl, * b1);title(1yl=xA2 ');subplot(2,2,3);bar(xzy2,'r');title(1y2=cos (2x)');subplot (2,2,4);bar(x,y3, * k ');title(1y3=yl*y2');由上面的M 文件，只要依次将“ bar” 改 为 “ stairs"、“ stem”、“ fill” ， 再适当更改区间取的点数，运行程序即可，即有下面的结果：3 . 已知y = I _ ____—ln(x + J l + f)x > 012在-5WxW5区间绘制函数曲线。
解：M 文件：clc;x=-5:0.01:5;y=(x+sqrt(pi))/ (exp(2)).*(x<=0)+0.5*1og(x+sqrt(1+x.A2)).* (x>0);plot(x,y)运行结果：由图可看出，函数在零点不连续4 .绘制极坐标曲线 「= 2$的8+ 110),并分析参数a、b、n 对曲线形状的影响解：M 文件如下：cic;theta=0:pi/100:2*pi; a=input( ' 输入 a=');b=input( ' 输入 b=');n=input( ' 输入 n=');rho=a*sin(b+n*theta);polar(theta,rho,'m')采用控制变量法的办法，固定两个参数，变动第三个参数观察输出图象的变化分析结果：由这8 个图知道，当 a,n固定时，图形的形状也就固定了，b 只影响图形的旋转的角度；当 a,b固定时，n 只影响图形的扇形数，特别地，当 n 是奇数时，扇叶数就是n,当是偶数时，扇叶数则是2n个；当 b,n固定时，a 影响的是图形大小，特别地，当 a 是整数时，图形半径大小就是a5 . 绘制函数的曲线图和等高线z = cos x cos ye 4其 中 x 的 2 1 个值均匀分布卜5, 5］ 范围，y 的 3 1 个值均匀分布在［ 0, 10］ , 要求使用subplot(2,1,1)和 subplot(2,1,2)将产生的曲面图和等高线图画在同一个窗口上。
解：M 文件：cic;x=linspace(-5,5,21);y=linspace(0,10,31);[x,y]=meshgrid(x,y)；z=cos(x).*cos(y).*exp(-sqrt(x.A2+y.A2)/4);subplot(2,1,1);surf(x,y,z);title( ' 曲面图' ) ；subplot(2,1,2);surfc(x,y,z);title( ' 等高线图力运行结果：6 . 绘制曲面图形，并进行插值着色处理X = cos s cos t.兀y = cos 5 sin/ 0 < 5 < —,0 < f< ——2 2. sins解：M 文件：clc;s=0:pi/100:pi/2; t=0:pi/100:3*pi/2;[s,t]=meshgrid(s,t);x=cos(s).*cos(t);y=cos(s).*sin(t);z=sin(s);subplot(2,2,l);mesh(x,y,z);title( 沫着色的图形) ;subplot(2,2,2);surf(x,y,z);title('shading faceted ( 缺 省 ) , ) ;subplot(2,2,3)；surf(x,y,z);shading flat;titleC^hading flat1);subplot(2,2,4)；surf(x,y,z);shading interp;titleC'shading interp*);运行结果有：实验七低层绘图操作二、实验内容1 . 建立一个图形窗口，使之背景颜色为红色，并在窗口上保留原有的菜单项，而且在按下鼠标器的左键之后显示出Left Button Pressed字样。
解：M 文件如下：cic;hf=figure('color',[1 0 0],...,WindowButtonDownFcn','disp("Left Button Pressed.")');运行结果:左击鼠标后：2 . 先利用默认属性绘制曲线y=x?e2x,然后通过图形句柄操作来改变曲线的颜色、线型和线宽，并利用文件对象给曲线添加文字标注解：M 文件：cic;x=-2:0.01:2; y=x.A2.*exp(2*x);h=plot(x,y);set(h,*color',[0.4,0.2.0.5],linestyle',.'linewidth',2);text(1.5,1.5A2*exp(2*1.5),'\leftarrow xA2exp(2x),,'fontsize',9);运行结果:3 . 利用曲面对象绘制曲面 v(x,t)=10e 0°1xsin(2000 n t-0.2x+ n )解：M 文件：cic;x=0:0.1:2*pi;[x,t]=meshgrid(x);v=10*exp(-0.01*x).*sin(2000*pi*t-0.2*x+pi);axes(,view',[-37,30]);hs=surface(x,t,v,'facecolor',...[0.2Q.3Q3],'edgecolor?flat');grid on;xlabel('x-axis'); ylabel('y-axis');zlabelfz-axis*);titleCmesh-surf*);pause %按任意键继续set(hs,'FaceColor7flat');text(O,O,O；曲面力运行结果：按任意键继续：4 . 以任意位置子图形式绘制出正弦、余弦、正切和余切函数曲线。
5 . 生成一个圆柱体，并进行光照和材质处理解：M 文件：[x,y,z]=cylinder(3,500); %cylinder 是生成柱体的函数surf(x,y,z);title，圆柱体的光照和材料处理' ) ;XlabelCX-axis*);YlabelCY-axis*);ZlabeKZ-axis*);axis([-5,5,-5,5,0,l])grid off;lightCColor'/r',Position',[-4,0,0],'style',infinite*);shading interp; material shiny;view(0,10);lighting phong;axis off;运行结果：实验八 数据处理与多项式计算1 .利用MATLAB提供的rand函数生成30000个符合均匀分布的随机数，然后检验随机数的性质：(1 )均值和标准方差2 )最大元素和最小元素3 )大于0.5的随机数个数占总数的百分比解：M 文件：cic;x=rand(1,30000);mu=mean(x) %求这30000个均匀分布随机数的平均值sig=std(x) %求其标准差1y=length(find(x>0.5)); % 找 出大于 0.5 数的个数p=y/30000 %大于0.5的所占百分比运行结果:mu =0.499488553231043sig =0.288599933559786P =0.4994000000000002 .将 100个学生5 门功课的成绩存入矩阵P 中，进行如下处理:(1 )分别求每门课的最高分、最低分及相应学生序号。
2 )分别求每门课的平均分和标准方差 ( 3) 5 门课总分的最高分、最低分及相应学生序号4 )将 5 门课总分按从大到小顺序存入zcj中，相应学生序号存入xsxh提示：上机调试时，为避免输入学生成绩的麻烦，可用取值范围在[45,95]之间的随机矩阵来表示学生成绩解：M 文件：cic;t=45+50*rand(100,5);P=fix(t); %生 成 100个学生5 门功课成绩[x,l]=max(P)% x为每门课最高分行向量, I 为相应学生序号[y,k]=min(P)% y为每门课最低分行向列, k 为相应学生序号mu=mean(P) % 每门课的平均值行向量sig=std(P) % 每门课的标准差行向量s=sum(P,2) % 5门课总分的列向量[X,m]=max(s)%5门课总分的最高分X 与相应学生序号m[Y,n]=min(s)%5门课总分的最低分Y 与相应学生序号n[zcj,xsxh]=sort(s)%zcj为 5 门课总分从大到小排序，相应学生序号xsxh运行结果：实验表1 室内外温度观测结果(°C)3 . 某气象观测得某|16:00~18:00之间每隔2h的室内外温度(0 C )如实验表1 所示。
时间h681012141618室内温度t118.020.022.025.030.028.024.0室外温度t215.019.024.028.034.032.030.0试用三次样条插值分别求出该日室内外6:30~18:30之间每隔2h各点的近似温度(°C)解：M 文件：cic;h=6:2:18;t1=[18.0 20.0 22.0 25.0 30.0 28.0 24.0]；t2=[15.0 19.0 24.0 28.0 34.0 32.0 30.0];T1=interp1(h,t1,'spline')%室内的3 次样条插值温度T2=interp1(h,t2,'spline')%室外的3 次样条插值温度运行结果：4 . 已知Igx在[1,101]区间10个整数采样点的函数值如实验表2 所示实验表2 Igx在 10个采样点的函数值F~ ~ 1 ~~11 21 31 41 51 61 71 81 97101_____________________Igx~~0~ ~ 1.0414~ ~ 1.32221.49141.6128 1.7076 1.7853 1.8513 1.9085 1.9510 2.0043试求Igx的 5 次拟合多项式p (x ),并绘制出Igx和 p(x)在[1,101]区间的函数曲线。
解：M 文件：x=1:10:101;y=igio(x)；P=polyfit(x,y,5)y1=polyval(P,x);plot(x,y,':o',x,y1运行结果：Warning: Polynomial is badly conditioned. Add points with distinct Xvalues, reduce the degree of the polynomial, or try centeringand scaling as described in HELP POLYFIT.> In polyfit at 80P =0.0000 -0.0000 0.0001 -0.0058 0.1537 -0.1326(这里出现警告是提示不必用5 价函数就已经可以完美拟合了，是可以降价拟合 )在[1,101]的区间函数图像5 .有 3 个多项式 P1(x)=x4+2x3+4x2+5, P2(x)=x+2, P3(x)=x2+2x+3,试进行下列操作：(1) P(x)=P1(x)+P2(x)P3(x)«(2 )求 P(x)的根3 )当 x 取矩阵A 的每一元素时，求 P(x)的值。
其 中 ：--I1.2-1.4-A =0.7523.5052.5_(4 )当以矩阵A 为自变量时，求 P(x)的值其中A 的值与第(3)题相同解：M 文件：clc;clear;p1 =[1,2,4,0,5]；p2=[1,2]；p3=[1,2,3]；p2=[0,0,0,p2];p3=[0,0,p3]；p4=conv(p2,p3);% p4是 p2与 p3的乘积后的多项式np4=length(p4);np1=length(p1);p=[zeros(1,np4-np1) p1]+p4% 求 p(x)=p1(x)+p2(x)x=roots(p)% 求 p(x)的根 A=[-1 1.2 -1.4;0.75 2 3.5;0 5 2.5];y=polyval(p,A) % x取矩阵A 的每一元素时的p(x)值运行结果:P =00001387 11x =-1.3840+ 1.8317i-1.3840- 1.8317i-0.1160+ 1.4400i-0.1160- 1.4400iy =1.0e+003 *0.0100 0.0382 0.01250.0223 0.0970 0.41220.0110 1.2460 0.1644实验九数值微积分与方程数值求解1 . 求函数在指定点的数值导数。
XX2x3/(x) =I2x3x2,x = l,2,3026x解：M 文件：clc;clear;x=1；i=1 ；f=inline('det([x xA2 xA3;1 2*x 3*xA2;0 2 6*x])');while x<=3.01g(i)=f(x)；』i+1;x=x+0.01; % 以 0.01的步长增加，可再缩小步长提高精度 运行结果:endg;t=1:0.01:3.01;dx=diff(g)/0.01;% 差分法近似求导f1=dx(1)%x=1的数值倒数f2=dx(101)%x=2的数值倒数f3=dx(length(g)-1)%x=3的数值倒数f1 =6.0602f2 =24.1202f3 =54.18022 . 用数值方法求定积分1) /, = £ Jcos产 +4sin⑵y + ldf 的近似值⑵A = f*力解：M 文件：clc;clear;f=inline('sqrt(cos(tA2)+4*sin(2*t).A2+1)');l1=quad(f,0,2*pi)g=inline(1og(1 +x)./(1 +x.A2),);l2=quad(gQ2*pi)运行结果：3 . 分别用3 种不同的数值方法解线性方程组。
6 x + 5 y - 2 z + 5 〃 = - 49 x - y + 4 z - 〃 = 1 3«：3 x + 4y + 2 z -2 〃 = 13 x -9 y + 2 〃 = 1 1解：M 文件：clc;clear; A=[6 5 -2 5;9 -1 4-1;342 -2;3 -9 0 2];b=[-4 13 1 11]1;x=A\by=inv(A)*b[L,U]=lu(A);z=U\(L\b)运行结果：4 .求非齐次线性方程组的通解2xl + lx2 + 3X3 + % = 6< 3X j + 5X2 + 2X3 + 2X4 = 49x] + 4X2 +X3 +7X4 = 2解：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)==ndispC 有唯一解 x*)；x=A\b;elsedisp('有无穷个解，特 解 x , 基础解系y');x=A\b;y=null(A,'r');endelsedisp。
无解) ;x=[]；endelse %齐次方程组dispC有零解x');x=zeros(n,1);if rank(A)
1) f( x ) = -------- - -----匚 在 (0,1)内的最小值e(2) /( xpx2) = 2x： +4x闻 — lox/? + x2在[0, ] 附近的最小值点和最小值解：M 文件：function f=g(u)x=u(1);y=u(2);f=2*x.A3+4*x.*yA3-10*x.*y+y.A2;clc;clear;format longf=inline('(xA3+cos(x)+x*log(x))/exp(x)');[x,fmin1 ]=fminbnd(f,0,1)[U,fmin2]=fminsearch('g,,[0,0])运行结果7 . 求微分方程的数值解’吗-5 生+ y = 0dx~ dx< y ( 0 ) = 0y ' ( o ) = o 解：M 文件：function xdot= sys( x,y)xdot=[y(2);(5*y(2)-y(1))/x];clc;clear;x0=1.0e-9;xf=20;[x,y]=ode45('sys',[x0,xf],[0 0]);[x,y]运行结果：8 .求微分方程组的数值解，并绘制解的曲线。
消 = 为 为y ‘2=f %y '3 =-051% %" 0 ) = 0,%(0) = 1,%( ) = 1解： 令 y1=x,y2=y,y3=z;这样方程变为：x' = yz< V1 = -xz， 自变量是z' = -0.51 盯x(0) = 0,y(0) = l,z(0) = lM 文件：function xdot=sys(x,y)xdot=[y(2)*y(3);-y(1)*y(3);-0.51 *y(1)*y(2)];clc;clear;t0=0;tf=8;[x,y]=ode23Csys',[tO,tf],[0,1,1])Plot(x,y)运行结果:图形:实验十 符号计算基础与符号微积分1 . 已知x=6,y=5,利用符号表达式求 _ x + 1z y/3+x-y[y提示: 定义符号常数x=sym(6), y=sym('5')解：M 文件：clear all;clc;x=sym('6');y=sym('5');z=(1 +x)/(sqrt(3+x)-sqrt(y))运行结果：z =-7/(5A( 1/2) - 3)2 . 分解因式⑴x R ⑵5135解：M 文件：clear all;clc;syms x y;t=sym('5135');a=xA4-yA4;factor(a)facto r(t)运行结果:ans =(x - y)*(x H卜 y)*(xA2 + yA2)ans =5*13*793.化简表达式。
1) sin d cos 夕2 - cos 尸]sin /72o 4x2 + 8x + 3I2-/ _ 42x + l解：M 文件：clear all;clc;syms betal beta2 x;f 1 =sin(beta1 )*cos(beta2)-cos(beta1 )*sin(beta2);simplify(fl) % ( 1 ) 问f2=(4*xA2+8*x+3)/(2*x+1);simplify(f2) % ( 2 ) 问 运行结果:ans =sin(beta1 - beta2)ans =2*x + 3解：M 文件:4 . 已知0 o-1 00 1 _,A =ad_Sb ce fh kP\ =010100()-01，P2 =101完成下列运算：(1) B=Pi • P2 • A 0(3 )包 括 B 矩阵主对角线元素的下三角阵2) B 的逆矩阵并验证结果4) B 的行列式值运行结果:clear all;clc;syms a b c d e f g h k;p1=[0 1 0;1 OO;OO 1];p2=[1 OO;O 1 O;1 0 1];A=[a b c;d e f;g h k];B=p1*p2*AB1=inv(B)% B的逆矩阵B1*B% 验证逆矩阵结果B2=tril(B)d=det(B)B =[ d, e, f][ a, b, c][ a + g, b + h, c + k]B1 =[ -(c*h - b*k)/(a*f*h - b*f*g - c*d*h + c*e*g - a*e*k + b*d*k), (b*f - c*e + f*h - e*k)/(a*f*h -b*f*g - c*d*h + c*e*g - a*e*k + b*d*k), -(b*f - c*e)/(a*f*h - b*f*g ・ c*d*h + c*e*g - a*e*k +b*d*k)][(c*g - a*k)/(a*f*h - b*f*g - c*d*h + c*e*g - a*e*k + b*d*k), -(a*f - c*d + f*g - d*k)/(a*f*h -b*f*g - c*d*h + c*e*g - a*e*k + b*d*k), (a*f - c*d)/(a*f*h - b*f*g - c*d*h + c*e*g - a*e*k + b*d*k)][(a*h - b*g)/(a*f*h - b*f*g - c*d*h + c*e*g - a*e*k + b*d*k), (a*e - b*d - d*h +e*g)/(a*f*h - b*f*g - c*d*h + c*e*g - a*e*k + b*d*k), -(a*e - b*d)/(a*f*h - b*f*g - c*d*h +c*e*g - a*e*k + b*d*k)]ans =[1,0,0][0,1,0][0,0, 11B2 =o叫ab ,“a ,d =a*f*h - b*f*g - c*d*h + c*e*g - a*e*k + b*d*k5 .用符号方法求卜列极限或导数。
x(esinr + l)-2(elant-l) G Ja rc c o sx(1) hm--------------z------------ (2) lim------------------I sin x iT-r VX + 1g 、 l-cos(2x) _p, , „ / 八「 优 / ~ ] dA d2A d2A(3)y = -------求 y',y" (4)已知 A= ， 分别求丁, 丁 , 丁7x tcosx Inx dx at dxdt(5)已知八x,y) = ( 一 一 2%成 —, 求 粤 总ox dxdyx=0,y=l解：M 文件:clear all;clc;syms x t a y z;f1 =(x*(exp(sin(x))+1 )-2*(exp(tan(x))-1 ))/sin(x)A3;limit(fl)%(1)f2=(sqrt(pi)-sqrt(acos(x)))/sqrt(x+1);%(2)Iimit(f2,x,-1 /right1)y=(1-cos(2*x))/x;y1=diff(y)y2=diff(y,2)%0)A=[aAx tA3;t*cos(x) log(x)];%G) Ax1=diff(A,x,1)At2=diff(A,t,2)Axt=diff(Ax1 ,t)f=(xA2-2*x)*exp(-xA2-zA2-x*z);Zx=-diff(f,x)/diff(f,z)dfxz=diff(diff(f,x),z);x=sym(,0');z=sym(,1');% 出)eval(dfxz)%符号运算返回数值运行结果：ans =-1/2ans =-Infyi =(2*sin(2*x))/x + (cos(2*x) - 1 )/xA2y2 =(4*cos(2*x))/x - (4*sin(2*x))/xA2 - (2*(cos(2*x) -1 ))/xA3Ax1 =[ a"*log ⑻，0][-t*sin(x), 1/x]At2 =[0, 6*t][0, 0] Axt =[ 0,0][ -sin(x), 0]Zx =-(exp(xA2 + x*z + zA2)*((2*x ・ 2)/exp(xA2 + x*z + zA2) + ((2*x - xA2)*(2*x + z))/exp(xA2 +x*z + zA2)))/((2*x - xA2)*(x + 2*z))ans =4/exp(1)6 . 用符号方法求下列积分。
p dx p dxJ T 7 7 7 7 % rcsinx)2 折 ——⑶(4) (n2e \\ + eyd xN X +1解：M 文件:clear;clc;x=sym(*x');f1=1/(1+xA4+xA8);%(1)f 2=1 /(asi n (x)) A2/sq rt( 1 -x A2);%⑵f3=(xA2+1)/(xA4+1);%⑶f4=exp(x)*(1 +exp(x))A2;F1=int(f1)F2=int(f2)F3=int(f3,0,inf)F4=int(f4,0,log(2))%(4)运行结果:Fl =(3A(l/2)*log(xA2 + 3A(l/2)*x + 1))/12 - (3A(l/2)*(2*atan(- (3A(l/2)*xA3)/3 - (2*3A(l/2)*x)/3)-2*atan((3A(l/2)*x)/3)))/12 - (3A(l/2)*log(xA2 - 3A(l/2)*x + 1 ))/12F2 = -l/asin(x)F3 =(pi*2A(l/2))/2F4 =exp(l 8729944304496077/9007199254740992)/3 + exp(6243314768165359/9007199254740992)+ exp(6243314768165359/4503599627370496) - 7/3实验十一级数与方程符号求解1 .级数符号求和。
10, 1( 1 )计算 s = £ —一2 〃- 1( 2 )求 级 数 的 和 函 数 ，并求X j之和解：M文件：clear all;clc;n=sym(* n *);x=sym(' x ');Sl=symsum(1/(2*n-l)z n, lz10)S2=symsum(nA2 *xA(n-1) , n, 1rinf)S3=symsum (n 人 2/5 人 n,n, 1 z inf) %vpa (S3)可以转化成小数运行结果：S1 =31037876/14549535S2 =piecewise([abs(x) < 1, -(xA2 + x)/(x*(x - 1)A3)]) S3 =15/322 . 将 Inx在 x=1处按5 次多项式展开为泰勒级数解：M 文件：clear all;clc;x=sym(1x1;f=log(x);taylor(f,x, 6, 1)运行结果：ans =x -(x -1 )A2/2 + (x - 1)A3/3 -(x -1 )A4/4 + (x - 1)A5/5 - 13 . 求下列方程的符号解5、 i----------(1) ln(l + x )--------------= 2 (2) x2 + + l -1 = 01 + sin xr \ o x c , ro u 八 / A\ + I'— —100 = 0(3) 3xe + 5 sin x - 78.5 = 0 (4)〈 、 ,3x + 5 y -8 = 0解：M 文件：clear all;clc;xl=solve(* log(x+1)-5/(1+sin(x))=2， )x2=solve('xA2+9*sqrt(x+1)-1')x3=solve(* 3*x*exp(x)+5*sin(x)-78.5*)[x4 y4]=solve(1sqrt(xA2+yA2)-100 \ f 3*x+5*y-8f)运行结果：x1 =521.67926389905839979437366649258 x2 =-1(3A(1/2)*i*(4/(9*(6465A(1/2)/2 + 2171/54)A(1/3)) - (1/2*6465A(1/2) + 2171/54)A(1/3)))/2 -(6465A(1/2)/2 + 2171/54)A(1/3)/2 - 2/(9*(6465A(1/2)/2 + 2171/54)A(1/3)) + 1/31/3 - (6465A(1/2)/2 + 2171/54)A(1/3)/2 - (3A(1/2)*i*(4/(9*(6465A(1/2)/2 + 2171/54)A(1/3))-(1/2*6465A(1/2) + 2171/54)A(1/3)))/2 - 2/(9*(6465A(1/2)/2 + 2171/54)八 (1/3))x3 =2.3599419584772910151699327715486x4 =12/17- (10*21246A(1/2))/17(10*21246A(1/2))/17 + 12/17y4 =(6*21246A(1/2))/17 + 20/1720/17-(6*21246A(1/2))/174 . 求微分方程初值问题的符号解，并与数值解进行比较。
2+ 4 立 + 29) ，= 0- dx2 dx[y(0) = 0,y1(0) = 15解：M 文件：clear all;clc;dsolve('D2y+4*Dy+29*y','y(0)=0'，'Dy(0)=15','x')运行结果：ans =(3*sin(5*x))/exp(2*x)5 . 求微分方程组的通解 dx A 八 八— =2x-3y + 3zdt— = 4x-5y + 3z+ = 4x-4y+ 2z解：M文件：clear all;clc;[x y z]=dsolve(,Dx=2*x-3*y+3*z,r ...» Dy=4 *x-5*y+3* z *z ,Dz=4*x-4*y+2*z'z 't')运行结果:X =C1/exp(t) + C2*exp(2*t)y =C1/exp(t) + C2*exp(2*t) + C3/exp(2*t)z =C2*exp(2*t) + C3/exp(2*t) 。
