MATLAB上机作业(数值分析方法(上理))
- 格式:pdf
- 大小:592.70 KB
- 文档页数:13
引论习题题目:用二分法求方程x^3-x-1=0在[1,2]内的近似根,要求误差不超过10-3程序:function [root,n]=dichotomy(a,b,eps)% dichotomy:二分法求根函数% f:带求解方程;[a,b]求解区间;n为二分次数if nargin<2 disp('输入变量不够');return;else if nargin>3 disp('输入变量过多');return;else nargin==2eps=1.0e-3;endendif a>b disp('区间输入错误')returnendc=(a+b)/2;if f(c)==0y=cn=1returnendn=0;while abs(b-a)>epsif f(a)*f(c)>0a=c;else b=c;endc=(a+b)/2;n=n+1;endroot=c;end运算结果:定义函数f如下:function [yy] = f(x)%定义了测试函数% x为输入变量yy=x^3-x-1;end在commond window中运行结果如下:>> [root,n]=dichotomy(1,2,1.0e-3)root =1.32471n =10习题一题目一给出概率积分2xxy e dx-=的数据表i 0 1 2 3xi 0.46 0.47 0.48 0.49yi 0.4846555 0.4937452 0.5027498 0.5116683 用二次插值计算,试问:(1)当x=0.472时该积分值等于多少?(2)当x为何值时积分值等于0.5?程序:function [ y,t ] = interpolation( X,Y,x )%拉格朗日插值函数%X为自变量数组,Y为函数值数组,x为插值点,t为插值基函数n=length(X);y=0;t=zeros(1,n);for i=1:nt(i)=1;for j=1:nif i==j continueendt(i)=t(i)*(x-X(j))/(X(i)-X(j));endy=y+t(i)*Y(i);endend运算结果:在commond window中运行结果如下:(1)当x=0.472时该积分值等于多少?>> clear>> X=[0.46 0.47 0.48 0.49];>> Y=[0.4846555 0.4937452 0.5027498 0.5116683];>> x=0.472;>> [ y,t ] = interpolation( X,Y,x )y =0.4956t =-0.0480 0.8640 0.2160 -0.0320(2)当x为何值时积分值等于0.5?>> clear>> X=[0.46 0.47 0.48 0.49];>> Y=[0.4846555 0.4937452 0.5027498 0.5116683];>> y=0.5;>> [ x,t ] = interpolation( Y,X,y )x =0.4769t =-0.0452 0.3356 0.7707 -0.0611题目二构造适合下列数据表的三次插值样条函数S(x):x -1 0 1 3 y -1 1 3 5 y' 6 1程序:function [ y,m ] =spline( X,Y,x,m1,mn )%样条插值函数%X,Y为输入的自变量、函数值数组;%x为插值点(数组),y为插值结果(数组),m为各插值节点的一阶导数值%m1为X(1)的一阶导数;mn为X(n)的一阶导数n=length(X);alpha=zeros(n-2,1);beta=zeros(n-2,1);for i=1:(n-2)alpha(i)=(X(i+1)-X(i))/(X(i+1)-X(i)+X(i+2)-X(i+1));beta(i)=3*((1-alpha(i))*(Y(i+1)-Y(i))/(X(i+1)-X(i))+alpha(i)*(Y(i+2)-Y(i+1))/(X(i+2)-X(i+1) ));endm=zeros(n,1);m(1)=m1;m(n)=mn;A=zeros((n-2),(n-2));B=zeros((n-2),1);for j=1:(n-2)A(j,j)=2;endfor k=2:(n-2)A(k,(k-1))=1-alpha(k);endfor l=1:(n-3)A(l,(l+1))=alpha(l);endB(1)=beta(1)-(1-alpha(1))*m1;B(n-2)=beta(n-2)-alpha(n-2)*mn;for p=2:(n-3)B(p)=beta(p);endm(2:(n-1))=A\B;s=length(x);y=zeros(1,s);for t=1:sr=1;for q=1:nif x(t)>=X(q)&&x(t)<=X(q+1)break;endr=r+1;endxx=(x(t)-X(r))/(X(r+1)-X(r));y(t)=(xx-1)^2*(2*xx+1)*Y(r)+xx^2*(-2*xx+3)*Y(r+1)+(X(r+1)-X(r))*xx*(xx-1)^2*m(r)+( X(r+1)-X(r))*xx^2*(xx-1)*m(r+1);endend运算结果:在commond window中运行结果如下:>> clear>> X=[-1 0 1 3];>> Y=[-1 1 3 5];>> m1=6;mn=1;>> x=-1:0.1:3;>> [ y,m ] =spline( X,Y,x,m1,mn );plot(x,y)求得结果如下(整理之后):其中红色数据表示插值节点将以上结果绘制成图,如下所示:习题二题目编写一通用型复化辛普生公式,能够对任意长度的等间距离散数据进行积分运算。
1.用二分法解方程 x-lnx=2 在区间【2 ,4】内的根方法: 二分法算法:f=inline('x-2-log(x)');a=2;b=4;er=b-a; ya=f(a);er0=.00001;while er>er0x0=.5*(a+b);y0=f(x0);if ya*y0<0b=x0;elsea=x0;ya=y0;enddisp([a,b]);er=b-a;k=k+1;end求解结果:>> answer13 43.0000 3.50003.0000 3.25003.1250 3.25003.1250 3.18753.1250 3.15633.1406 3.15633.1406 3.14843.1445 3.1484 3.1445 3.1465 3.1455 3.1465 3.1460 3.1465 3.1460 3.1462 3.1461 3.1462 3.1462 3.14623.1462 3.1462 3.1462 3.1462 3.1462 3.1462 最终结果为: 3.14622.试编写MATLAB 函数实现Newton 插值,要求能输出插值多项式。
对函数141)(2+=x x f 在区间[-5,5]上实现10次多项式插值。
Matlab 程序代码如下:%此函数实现y=1/(1+4*x^2)的n 次Newton 插值,n 由调用函数时指定 %函数输出为插值结果的系数向量(行向量)和插值多项式 算法:function [t y]=func5(n) x0=linspace(-5,5,n+1)'; y0=1./(1.+4.*x0.^2); b=zeros(1,n+1); for i=1:n+1 s=0; for j=1:i t=1; for k=1:iif k~=jt=(x0(j)-x0(k))*t;end;end;s=s+y0(j)/t;end;b(i)=s;end;t=linspace(0,0,n+1);for i=1:ns=linspace(0,0,n+1);s(n+1-i:n+1)=b(i+1).*poly(x0(1:i));t=t+s;end;t(n+1)=t(n+1)+b(1);y=poly2sym(t);10次插值运行结果:[b Y]=func5(10)b =Columns 1 through 4-0.0000 0.0000 0.0027 -0.0000Columns 5 through 8-0.0514 -0.0000 0.3920 -0.0000Columns 9 through 11-1.1433 0.0000 1.0000Y =- (7319042784910035*x^10)/147573952589676412928 + x^9/18446744073709551616 + (256*x^8)/93425 -x^7/1152921504606846976 -(28947735013693*x^6)/562949953421312 -(3*x^5)/72057594037927936 + (36624*x^4)/93425 -(5*x^3)/36028797018963968 -(5148893614132311*x^2)/4503599627370496 +(7*x)/36028797018963968 + 1b为插值多项式系数向量,Y为插值多项式。
数值分析———Matlab上机作业学院:班级:老师:姓名:学号:第二章解线性方程组的直接解法第14题【解】1、编写一个追赶法的函数输入a,b,c,d输出结果x,均为数组形式function x=Zhuiganfa(a,b,c,d)%首先说明:追赶法是适用于三对角矩阵的线性方程组求解的方法,并不适用于其他类型矩阵。
%定义三对角矩阵A的各组成单元。
方程为Ax=d%b为A的对角线元素(1~n),a为-1对角线元素(2~n),c为+1对角线元素(1~n-1)。
% A=[2 -1 0 0% -1 3 -2 0% 0 -2 4 -3% 0 0 -3 5]% a=[-1 -2 -3];c=[-1 -2 -3];b=[2 3 4 5];d=[6 1 -2 1];n=length(b);u(1)=b(1);y(1)=d(1);for i=2:nl(i)=a(i-1)/u(i-1);%先求l(i)u(i)=b(i)-c(i-1)*l(i);%再求u(i)%A=LU,Ax=LUx=d,y=Ux,%Ly=d,由于L是下三角矩阵,对角线均为1,所以可求y(i)y(i)=d(i)-l(i)*y(i-1);endx(n)=y(n)/u(n);for i=(n-1):-1:1%Ux=y,由于U是上三角矩阵,所以可求x(i)x(i)=(y(i)-c(i)*x(i+1))/u(i);end2、输入已知参数>>a=[2 2 2 2 2 2 2];>>b=[2 5 5 5 5 5 5 5];>>c=[2 2 2 2 2 2 2];>>d=[220/27 0 0 0 0 0 0 0];3、按定义格式调用函数>>x=zhuiganfa(a,b,c,d)4、输出结果x=[8.147775166909105 -4.073701092835030 2.036477565178471 -1.017492820111148 0.507254485099400 -0.250643392637350 0.119353996493976 -0.047741598597591]第15题【解】1、编写一个程序生成题目条件生成线性方程组A x=b 的系数矩阵A 和右端项量b ,分别定义矩阵A 、B 、a 、b 分别表示系数矩阵,其中1(10.1;,1,2,...,)j ij i i a x x i i j n -==+=或1(,1,2,...,)1ij a i j n i j ==+-分别构成A 、B 对应右端项量分别a 、b 。
数值计算上机实习题目(matlab编程)非线性方程求根一、实验目的本次实验通过上机实习,了解迭代法求解非线性方程数值解的过程和步骤。
二、实验要求1、用迭代法求方程230x x e -=的根。
要求:确定迭代函数?(x),使得x=?(x),并求一根。
提示:构造迭代函数2ln(3)x ?=。
2、对上面的方程用牛顿迭代计算。
3、用割线法求方程3()310f x x x =--=在02x =附近的根。
误差限为410-,取012, 1.9x x ==。
三、实验内容1、(1)首先编写迭代函数,记为iterate.mfunction y=iterate(x)x1=g(x); % x 为初始值。
n=1;while(abs(x1-x)>=1.0e-6)&(n<=1000) % 迭代终止的原则。
x=x1;x1=g(x);n=n+1;endx1 %近似根n %迭代步数(2)后编制函数文件?(x),记为g.mfunction y=g(x)y=log(3*x.^2);(3)设初始值为0、3、-3、1000,观察初始值对求解的影响。
将结果记录在文档中。
>>iterate(0)>>iterate(3) 等等2、(1)首先编制牛顿迭代函数如下,记为newton.mfunction y=newton(x0)x1=x0-fc(x0)/df(x0); % 牛顿迭代格式n=1;while(abs(x1-x0)>=1.0e-6)&(n<=1000000) % 迭代终止的原则。
x0=x1;x1=x0-fc(x0)/df(x0);n=n+1;endx1 %近似根n %迭代步数(2)对题目中的方程编制函数文件,记为fc.mfunction y=fc(x)y=3*x.^2-exp(x)编制函数的导数文件,记为df.mfunction y=df(x)y=6*x-exp(x)(3)在MATLAB 命令窗计算,当设初始值为0时,newton(0);给定不同的初始值,观察用牛顿法求解时所需要的迭代步数,并与上面第一题的迭代步数比较。
MATLAB上机内容及作业无约束优化求解函数fminsearch和fminunc求解无约束非线性优化问题的函数包括fminsearch函数和fminunc函数。
fminsearch和fminunc的函数相同,但fminunc函数可以在最优解处获得目标函数的梯度和Hessian矩阵值。
无约束优化数学模型为:minf(x)x∈rn求解无约束非线性优化问题的步骤为:第一步:先写目标函数的m文件;第二步是在命令窗口中调用相应的优化函数。
1.Fminsearch函数调用格式为[x,fval]=fminsearch(@myfun,x0)输出参数的含义:x:返回最优解的设计变量值;fval:在最优设计变量值时,目标函数的最小值;Exitflag:返回算法的终止标志。
在下列情况下,>0表示算法收敛到最优解;=0表示算法已达到最大迭代次数;<0表示算法不收敛。
output:返回优化算法信息的一个数据结构,有以下信息:输出Iteration表示迭代次数output Algorithm表示算法output Funccount表示函数求值的次数输入参数的含义:@Myfun:目标函数的m文件,前面加“@”或表示为“Myfun”,可以任意命名;X0:调用优化函数时,需要先给设计变量赋值;2、fminunc函数的调用格式[x,fval]=fminunc(@myfun,x0)grad:返回目标函数在最优解处的梯度信息;hessian:返回目标函数在最优解处的hessian矩阵信息。
其余含义同上。
3、实例已知某一优化问题的数学模型为:Minf(x)=3x12+2x1x2+x22x∈ r n由MATLAB程序编写的代码为:第一步:首先编写目标函数的.m文件,并保存为examplefsearch.m文件(先单击file菜单,后点击new命令中的m―file,即可打开m文件编辑窗口进行代码的编辑,在英文状态下输入程序代码),代码为:函数f=examplefsearch(x)f=3*x(1)^2+2*x(1)*x(2)+x(2)^2;1.0e-0.08*-0.79140.2260%分别为x1和x2的最优点的值(近似为0)Fval=1.5722e-016%,对应于最佳优势的最佳目标函数值(约为0)4。
数值分析实习报告姓名:gestepoA学号:201*******班级:***班序言随着计算机技术的迅速发展,数值分析在工程技术领域中的应用越来越广泛,并且成为数学与计算机之间的桥梁。
要解决工程问题,往往需要处理很多数学模型,不仅要研究各种数学问题的数值解法,同时也要分析所用的数值解法在理论上的合理性,如解法所产生的误差能否满足精度要求:解法是否稳定、是否收敛及熟练的速度等。
而且还能减少大量的人工计算。
由于工程实际中所遇到的数学模型求解过程迭代次数很多,计算量很大,所以需要借助如MATLAB,C++,VB,JAVA的辅助软件来解决,得到一个满足误差限的解。
本文所计算题目,均采用MATLAB进行编程,MATLAB被称为第四代计算机语言,利用其丰富的函数资源,使编程人员从繁琐的程序代码中解放出来MATLAB最突出的特点就是简洁,它用更直观的、符合人们思维习惯的代码。
它具有以下优点:1友好的工作平台和编程环境。
MATLAB界面精致,人机交互性强,操作简单。
2简单易用的程序语言。
MATLAB是一个高级的矩阵/阵列语言,包含控制语言、函数、数据结构,具有输入、输出和面向对象编程特点。
用户可以在命令窗口中将输入语句与执行命令同步,也可以先编好一个较大的复杂的应用程序(M文件)后再一起运行。
3强大的科学计算机数据处理能力。
包含大量计算算法的集合,拥有600多个工程中要用到的数学运算函数。
4出色的图像处理功能,可以方便地输出二维图像,便于我们绘制函数图像。
目录1 第一题.................................................................. 4 1.1 实验目的ﻩ41.2 实验原理和方法 (4)1.3 实验结果........................................................... 51.3.1 最佳平方逼近法ﻩ51.3.2 拉格朗日插值法ﻩ71.3.3 对比ﻩ82 第二题.................................................................. 92.1实验目的9ﻩ2.2实验原理和方法ﻩ10102.3实验结果ﻩ2.3.1 第一问10ﻩ2.3.2 第二问 (11)2.3.3 第三问1ﻩ13第三题 (12)3.1实验目的 (12)123.2 实验原理和方法ﻩ3.3实验结果1ﻩ24 MATLAB程序 (14)1第一题某过程涉及两变量x和y,拟分别用插值多项式和多项式拟合给出其对应规律的4.65880.3719.64484.272170345.2795 43 7.48478.2392 ⑴请用次数分别为3,4,5,6的多项式拟合并给出最好近似结果f(x)。
电子科技大学电子工程学院标准实验报告(实验)课程名称MATLAB与数值分析学生姓名:李培睿学号:2013020904026指导教师:程建一、实验名称《MATLAB与数值分析》第一次上机实验二、实验目的1. 熟练掌握矩阵的生成、加、减、乘、除、转置、行列式、逆、范数等运算操作。
(用.m文件和Matlab函数编写一个对给定矩阵进行运算操作的程序)2. 熟练掌握算术符号操作和基本运算操作,包括矩阵合并、向量合并、符号转换、展开符号表达式、符号因式分解、符号表达式的化简、代数方程的符号解析解、特征多项式、函数的反函数、函数计算器、微积分、常微分方程的符号解、符号函数的画图等。
(用.m文件编写进行符号因式分解和函数求反的程序)3. 掌握Matlab函数的编写规范。
4、掌握Matlab常用的绘图处理操作,包括:基本平面图、图形注释命令、三维曲线和面的填充、三维等高线等。
(用.m文件编写在一个图形窗口上绘制正弦和余弦函数的图形,并给出充分的图形注释)5. 熟练操作MATLAB软件平台,能利用M文件完成MATLAB的程序设计。
三、实验内容1. 编程实现以下数列的图像,用户能输入不同的初始值以及系数。
并以x,y为坐标显示图像x(n+1) = a*x(n)-b*(y(n)-x(n)^2);y(n+1) = b*x(n)+a*(y(n)-x(n)^2)2. 编程实现奥运5环图,允许用户输入环的直径。
3. 实现对输入任意长度向量元素的冒泡排序的升序排列。
不允许使用sort 函数。
四、实验数据及结果分析题目一:①在Editor窗口编写函数代码如下:并将编写的函数文件用“draw.m”储存在指定地址;②在Command窗口输入如下命令:③得到图形结果如下:题目二:①在Editor窗口编写函数代码如下:并将编写的函数文件用“circle.m”储存在指定地址;②再次在Editor窗口编写代码:并将编写的函数文件用“Olympic.m”储存在指定地址;③在Command窗口输入如下指令(半径可任意输入):④按回车执行,将在图形窗口获得五环旗:题目三:①在Editor窗口编写函数代码如下:并用.将编写的函数文件用“qipaofa.m”储存在指定地址;②在Command窗口输入一组乱序数值,则可以得到升序排序结果如下:五、总结及心得体会1.要熟悉MATLAB编译软件的使用方法,明白有关语法,语句的基本用法,才可以在编写程序的时候游刃有余,不至于寸步难行。
不动点迭代:运行结果:代码截图:牛顿迭代:运行结果:代码截图:1)运行结果:代码截图:a二分法运行结果:代码截图:牛顿法:运行结果:代码截图:弦位法:代码截图:第二次上机作业运行结果:代码截图:A欧拉方法运行结果:代码截图:B改进的欧拉方法运行结果:代码截图:C四阶龙格库塔方法:代码截图:第三次上机作业第6章:(7版)P368: 5(e), 8, 9Using1) Gaussian elimination2) Gaussian elimination with partial pivoting.3) Gaussian elimination with scaled partial pivoting and three-digit chopping arithmetic to solvethe following linear systems, and compare the approximations to the actual solution.1.19x1 +2.11x2 − 100x3 + x4 = 1.12,14.2x1 − 0.122x2 + 12.2x3 − x4 = 3.44,100x2 − 99.9x3 + x4 = 2.15,15.3x1 + 0.110x2 − 13.1x3 − x4 = 4.16.Actual solution [0.176, 0.0126,−0.0206,−1.18]解:1)高斯消元法运行结果代码截图:2)部分选主元高斯消元法运行结果:代码截图:3)具有比例因子部分选主元高斯消元法运行结果:代码截图:2.Solve Ax = b using the Crout factorization for tridiagonal systems. Let A be the 10×10 tridiagonal matrix given byaii = 2, ai,i+1 = ai,i−1 = −1, for each i = 2, . . . , 9and a11 = a10,10 = 2, a12 = a10,9 = −1.Let b be the ten-dimensional column vector given byb1 = b10 = 1 and bi = 0, for each i = 2, 3, . . . , 9.运行结果:代码截图:3.1)雅克比迭代运行结果:程序截图:2)GS迭代运行结果代码截图:3)SOR方法运行结果:代码截图:。
第二次上机作业一. 任务:用MATLAB 语言编写连续函数最佳平方逼近的算法程序(函数式M 文件)。
并用此程序进行数值试验,写出计算实习报告。
二. 程序功能要求:在后面的附一leastp.m 的基础上进行修改,使其更加完善。
要求算法程序可以适应不同的具体函数,具有一定的通用性。
所编程序具有以下功能:1. 用Lengendre 多项式做基,并适合于构造任意次数的最佳平方逼近多项式。
可利用递推关系 0112()1,()()(21)()(1)()/2,3,.....n n n P x P x xP x n xP x n P x n n --===---⎡⎤⎣⎦=2. 被逼近函数f(x)不用内联函数构造,而改用M 文件建立数学函数。
这样,此程序可通过修改建立数学函数的M 文件以适用不同的被逼近函数(要学会用函数句柄)。
3. 要考虑一般的情况]1,1[],[)(+-≠∈b a x f 。
因此,程序中要有变量代换的功能。
4. 计算组合系数时,计算函数的积分采用变步长复化梯形求积法(见附三)。
5. 程序中应包括帮助文本和必要的注释语句。
另外,程序中也要有必要的反馈信息。
6. 程序输入:(1)待求的被逼近函数值的数据点0x (可以是一个数值或向量)(2)区间端点:a,b 。
7. 程序输出:(1)拟合系数:012,,,...,n c c c c(2)待求的被逼近函数值00001102200()()()()()n n s x c P x c P x c P x c P x =++++三:数值试验要求:1. 试验函数:()cos ,[0,4]f x x x x =∈+;也可自选其它的试验函数。
2. 用所编程序直接进行计算,检测程序的正确性,并理解算法。
3. 分别求二次、三次、。
最佳平方逼近函数)x s (。
4. 分别作出逼近函数)x s (和被逼近函数)(x f 的曲线图进行比较。
(分别用绘图函数plot(0x ,s(0x ))和fplot(‘x cos x ’,[x 1 x 2,y 1,y 2]))四:计算实习报告要求:1.简述方法的基本原理,程序功能,使用说明。
数值分析作业MATLAB数值分析是研究用数值方法解决数学问题的一门学科。
它的主要目标是通过计算机编程解决数学问题,尤其是那些无法通过解析方法解决的问题。
MATLAB是一种常用的数值分析软件,它提供了丰富的数值计算函数和工具箱,能够方便地进行各种数值分析方法的实现和计算。
数值分析的研究内容很广泛,包括数值计算方法、数值逼近、数值微分和数值积分等。
在数值计算方法中,最常用的有数值解线性方程组、数值解非线性方程、数值积分、数值微分等。
例如,通过使用MATLAB的线性方程组求解函数或者工具箱中的线性代数函数,可以解决各种形式的线性方程组。
通过MATLAB的非线性方程求解函数,可以解决各种非线性方程的数值解。
而数值积分和数值微分则可以通过MATLAB的积分函数和微分函数来实现,实现对函数的积分和微分操作。
数值逼近是数值分析的重要内容之一,它研究的是如何用简单的函数逼近给定的复杂函数。
在MATLAB中,可以通过多项式逼近、三次样条、拉格朗日插值、最小二乘逼近等方法来实现数值逼近的计算。
例如,使用MATLAB的插值函数interp1可以实现一维函数的插值计算,使用MATLAB 的polyfit函数可以拟合一维数据集合的多项式曲线。
而对于二维函数和三维函数的逼近,可以使用MATLAB的interp2和interp3函数来实现。
数值微分和数值积分是数值分析中的基本操作之一、它们可以根据给定的函数计算函数的导数和积分。
在MATLAB中,使用diff函数可以计算一维函数的导数,使用trapz和quad函数可以计算一维函数的定积分和数值积分。
而对于二维函数和三维函数的微分和积分,可以使用MATLAB 的grad函数和integral2函数来实现。
此外,MATLAB还提供了很多其他的数学函数和工具,包括解微分方程、优化问题、曲线拟合和最小二乘等。
对于一些复杂的数学问题,可以通过使用MATLAB的符号计算工具箱来实现符号计算。
1.用两种方法解sin 5cos 20.31x x +=-,如何一次求出此方程的四个根。
解:方法一:在Matlab 命令窗口输入命令x=solve('sin(5*x)+cos(2*x)=-0.31')运行后输出结果:x =- 0.36685340479904406504913603551901 - 0.089925518994485866153169602644131*i 方法二:在Matlab 命令窗口输入命令E1=sym('sin(5*x)+cos(2*x)=-0.31');[x1]=solve(E1)运行后输出结果:x1 =- 0.36685340479904406504913603551901 - 0.089925518994485866153169602644131*i 2.用三种方法解方程11852912312x x x x -+-=。
解:方法一:在Matlab 命令窗口输入命令x=solve('9*x^11-12*x^8+x^5-3*x^2=12')运行后输出结果:x =1.1945355092112803561808169303487 - 0.25878939596021341127295614138703 + 0.98996423045277565907337492165509*i - 0.91081125361412082546165956220494 + 0.34557668006489020731472236114125*i - 0.6567748898900820562684890790666 - 0.94093805304063227438930031737768*i - 0.6567748898900820562684890790666 + 0.94093805304063227438930031737768*i 0.34695114229720670305614226506505 + 0.85428006847946699100354057810685*i - 0.25878939596021341127295614138703 - 0.98996423045277565907337492165509*i 0.88215664256156941185655405241917 + 0.47468310614563172153151897912015*i 0.34695114229720670305614226506505 - 0.85428006847946699100354057810685*i 0.88215664256156941185655405241917 - 0.47468310614563172153151897912015*i - 0.91081125361412082546165956220494 - 0.34557668006489020731472236114125*i 方法二:在Matlab 命令窗口输入命令fa=[9,0,0,-12,0,0,1,0,0,-3,0,-12];xk=roots(fa)运行后输出结果:xk =-0.9108 + 0.3456i-0.9108 - 0.3456i-0.6568 + 0.9409i-0.6568 - 0.9409i-0.2588 + 0.9900i-0.2588 - 0.9900i1.1945 + 0.0000i0.8822 + 0.4747i0.8822 - 0.4747i0.3470 + 0.8543i0.3470 - 0.8543i方法三:在Matlab命令窗口输入命令E1=sym('9*x^11-12*x^8+x^5-3*x^2-12=0'); [x1]=solve(E1)运行后输出结果:x1 =1.1945355092112803561808169303487 - 0.25878939596021341127295614138703 + 0.98996423045277565907337492165509*i - 0.91081125361412082546165956220494 + 0.34557668006489020731472236114125*i - 0.6567748898900820562684890790666 - 0.94093805304063227438930031737768*i - 0.6567748898900820562684890790666 + 0.94093805304063227438930031737768*i 0.34695114229720670305614226506505 + 0.85428006847946699100354057810685*i - 0.25878939596021341127295614138703 - 0.98996423045277565907337492165509*i 0.88215664256156941185655405241917 + 0.47468310614563172153151897912015*i 0.34695114229720670305614226506505 - 0.85428006847946699100354057810685*i 0.88215664256156941185655405241917 - 0.47468310614563172153151897912015*i - 0.91081125361412082546165956220494 - 0.34557668006489020731472236114125*3 用MATLAB 方法求函数11852()912312f x x x x x =-+--的导数'()f x 。
Newton法及改进的Newton法源程序:clear%%%% 输入函数f=input('请输入需要求解函数>>','s')%%%求解f(x)的导数df=diff(f);%%%改进常数或重根数miu=2;%%%初始值x0x0=input('input initial value x0>>');k=0;%迭代次数max=100;%最大迭代次数R=eval(subs(f,'x0','x'));%求解f(x0),以确定初值x0时否就是解while (abs(R)>1e-8)x1=x0-miu*eval(subs(f,'x0','x'))/eval(subs(df,'x0','x'));R=x1-x0;x0=x1;k=k+1;if (eval(subs(f,'x0','x'))<1e-10);breakendif k>max;%如果迭代次数大于给定值,认为迭代不收敛,重新输入初值ss=input('maybe result is error,choose a new x0,y/n?>>','s');ifstrcmp(ss,'y')x0=input('input initial value x0>>');k=0; elsebreakendendendk;%给出迭代次数x=x0;%给出解Gauss消元法源程序:cleara=input('输入系数阵:>>\n')b=input('输入列阵b:>>\n')n=length(b);A=[a b]x=zeros(n,1);%%%函数主体Yipusilong-1;for k=1:n-1;%%%是否进行主元选取if abs(A(k,k))<yipusilong;%事先给定的认为有必要选主元的小数yzhuyuan=1;elseyzhuyuan=0;endifyzhuyuan;%%%%选主元t=A(k,k);for r=k+1:n;if abs(A(r,k))>abs(t)p=r;else p=k;endend%%%交换元素if p~=k; for q=k:n+1;s=A(k,q);A(k,q)=A(p,q);A(p,q)=s;endendend%%%判断系数矩阵是否奇异或病态非常严重if abs(A(k,k))<yipusilongdisp(‘矩阵奇异,解可能不正确’)end%%%%计算消元,得三角阵for r=k+1:n;m=A(r,k)/A(k,k);for q=k:n+1;A(r,q)=A(r,q)-A(k,q)*m;endendend%%%%求解xx(n)=A(n,n+1)/A(n,n);for k=n-1:-1:1;s=0;for r=k+1:n;s=s+A(k,r)*x(r);endt=(A(k,n+1)-s) x(k)=(A(k,n+1)-s)/A(k,k) end Lagrange插值源程序:n=input('将区间分为的等份数输入:\n');s=[-1+2/n*[0:n]];%%%给定的定点,Rf为给定的函数x=-1:0.01:1;f=0;for q=1:n+1;l=1;%求插值基函数for k=1:n+1;if k~=q;l=l.*(x-s(k))./(s(q)-s(k));elsel=l;endendf=f+Rf(s(q))*l;%求插值函数endplot(x,f,'r')%作出插值函数曲线grid onhold on分段线性插值源程序clearn=input('将区间分为的等份数输入:\n');s=[-1+2/n*[0:n]];%%%给定的定点,Rf为给定的函数m=0;hh=0.001;for x=-1:hh:1;ff=0;for k=1:n+1;%%%求插值基函数switch kcase 1if x<=s(2);l=(x-s(2))./(s(1)-s(2));elsel=0;endcase n+1if x>s(n);l=(x-s(n))./(s(n+1)-s(n)); elsel=0;endotherwiseif x>=s(k-1)&x<=s(k);l=(x-s(k-1))./(s(k)-s(k-1)); else if x>=s(k)&x<=s(k+1);l=(x-s(k+1))./(s(k)-s(k+1)); elsel=0;endendendff=ff+Rf(s(k))*l;%%求插值函数值endm=m+1;f(m)=ff;end%%%作出曲线x=-1:hh:1;plot(x,f,'r');grid onhold on三次样条插值源程序:(采用第一边界条件)clearn=input('将区间分为的等份数输入:\n');%%%插值区间a=-1;b=1;hh=0.001;%画图的步长s=[a+(b-a)/n*[0:n]];%%%给定的定点,Rf为给定的函数%%%%第一边界条件Rf"(-1),Rf"(1)v=5000*1/(1+25*a*a)^3-50/(1+25*a*a)^4;for k=1:n;%取出节点间距h(k)=s(k+1)-s(k);endfor k=1:n-1;%求出系数向量lamuda,miula(k)=h(k+1)/(h(k+1)+h(k));miu(k)=1-la(k);end%%%%赋值系数矩阵Afor k=1:n-1;for p=1:n-1;switch pcase kA(k,p)=2;case k-1A(k,p)=miu(p+1);case k+1A(k,p)=la(p-1);otherwiseA(k,p)=0;endendend%%%%求出d阵for k=1:n-1;switch kcase 1d(k)=6*f2c([s(k) s(k+1) s(k+2)])-miu(k)*v; case n-1d(k)=6*f2c([s(k) s(k+1) s(k+2)])-la(k)*v; otherwised(k)=6*f2c([s(k) s(k+1) s(k+2)]);endend %%%%求解M阵M=A\d';M=[v;M;v];%%%%m=0;f=0;for x=a:hh:b;if x==a;p=1;elsep=ceil((x-s(1))/((b-a)/n));endff1=0;ff2=0;ff3=0;ff4=0;m=m+1;ff1=1/h(p)*(s(p+1)-x)^3*M(p)/6;ff2=1/h(p)*(x-s(p))^3*M(p+1)/6;ff3=((Rf(s(p+1))-Rf(s(p)))/h(p)-h(p)*(M(p+1)-M(p))/6)*(x-s(p));ff4=Rf(s(p))-M(p)*h(p)*h(p)/6;f(m)=ff1+ff2+ff3+ff4 ;end%%%作出插值图形x=a:hh:b;plot(x,f,'k')hold on grid on 多项式最小二乘法源程序clear%%%给定测量数据点(s,f)s=[3 4 5 6 7 8 9];f=[2.01 2.98 3.50 5.02 5.47 6.02 7.05];%%%计算给定的数据点的数目n=length(f);%%%给定需要拟合的数据的最高次多项式的次数m=10;%%%程序主体for k=0:m;g=zeros(1,m+1);for j=0:m;t=0;for i=1:n;%计算内积(fai(si),fai(si))t=t+fai(s(i),j)*fai(s(i),k);endg(j+1)=t;endA(k+1,:)=g;%法方程的系数矩阵t=0;for i=1:n;%计算内积(f(si),fai(si))t=t+f(i)*fai(s(i),k);endb(k+1,1)=t;enda=A\b%求出多项式系数x=[s(1):0.01:s(n)]';y=0;fori=0:m;y=y+a(i+1)*fai(x,i);endplot(x,y)%作出拟合成的多项式的曲线grid onhold onplot(s,f,'rx') %在上图中标记给定的点表中,L10(x)为Lagrang e插值的10次多项式,S10(x),S40(x)分别代表n=10,40的三次样条插值函数,X10(x),X40(x)分别代表n=10,40的线性分段插值函数。
实验报告课程名称:数值分析实验项目:曲线拟合/数值积分专业班级:姓名:学号:实验室号:实验组号:实验时间:20.10.24 批阅时间:指导教师:成绩:工业大学实验报告(适用计算机程序设计类)专业班级:学号:姓名:实验名称:曲线拟合与函数插值附件A 工业大学实验报告(适用计算机程序设计类)专业班级:学号:姓名:实验步骤或程序:附录一:1.利用二次,三次,四次多项式进行拟合:1.1 MATLAB代码如下:clear;clc;close allt=[0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24];y=[14 13 13 13 13 14 15 17 19 21 22 24 27 30 31 30 28 26 24 23 21 19 17 16 15];%输入数据hold on[p2 s2]=polyfit(t,y,2);%对于上面的数据进行2次多项式拟合,其中s2包括R(系数矩阵的QR分解的上三角阵),%df(自由度),normr(拟合误差平方和的算术平方根)。
y2=polyval(p2,t);%返回多项式拟合曲线在t处的值[p3 s3]=polyfit(t,y,3);y3=polyval(p3,t);[p4 s4]=polyfit(t,y,4);y4=polyval(p4,t);plot(t,y,'ro')%画图plot(t,y2,'g-')plot(t,y3,'m^-')plot(t,y4,'bs-')xlabel('t')ylabel('y')legend('原始数据','2次多项式拟合','3次多项式拟合','4次多项式拟合')1.2 二次,三次,四次多项式拟合的结果分别如下:(1)总的拟合结果在工作区的显示如下:(2)其次二次多项式拟合的结果为:(3)其中三次多项式拟合的结果:(4)其中四次多项式拟合的结果为:1.3 拟合的图像为:1.4 拟合的多项式为:根据工作区得出的数据列出最后的拟合多项式为:(1)y=7.416+2.594t-0.094t^2(2)y=12.251-0.102t+0.193t^2-0.008t^3(3)y=15.604-3.526t+0.866t^2-0.052t^3+0.0009t^42.形如2()()b t c y t ae--=的函数,其中,,a b c 为待定常数。
一、实验目的通过本次上机实验,掌握数值分析中常用的算法,如二分法、牛顿法、不动点迭代法、弦截法等,并能够运用这些算法解决实际问题。
同时,提高编程能力,加深对数值分析理论知识的理解。
二、实验环境1. 操作系统:Windows 102. 编程语言:MATLAB3. 实验工具:MATLAB数值分析工具箱三、实验内容1. 二分法求方程根二分法是一种常用的求方程根的方法,适用于连续函数。
其基本思想是:从区间[a, b]中选取中点c,判断f(c)的符号,若f(c)与f(a)同号,则新的区间为[a, c],否则为[c, b]。
重复此过程,直至满足精度要求。
2. 牛顿法求方程根牛顿法是一种迭代法,适用于可导函数。
其基本思想是:利用函数在某点的导数值,求出函数在该点的切线方程,切线与x轴的交点即为方程的近似根。
3. 不动点迭代法求方程根不动点迭代法是一种迭代法,适用于具有不动点的函数。
其基本思想是:从初始值x0开始,不断迭代函数g(x)的值,直至满足精度要求。
4. 弦截法求方程根弦截法是一种线性近似方法,适用于可导函数。
其基本思想是:利用两点间的直线近似代替曲线,求出直线与x轴的交点作为方程的近似根。
四、实验步骤1. 二分法求方程根(1)编写二分法函数:function [root, error] = bisection(a, b, tol)(2)输入初始区间[a, b]和精度要求tol(3)调用函数计算根:[root, error] = bisection(a, b, tol)2. 牛顿法求方程根(1)编写牛顿法函数:function [root, error] = newton(f, df, x0, tol)(2)输入函数f、导数df、初始值x0和精度要求tol(3)调用函数计算根:[root, error] = newton(f, df, x0, tol)3. 不动点迭代法求方程根(1)编写不动点迭代法函数:function [root, error] = fixed_point(g, x0, tol)(2)输入函数g、初始值x0和精度要求tol(3)调用函数计算根:[root, error] = fixed_point(g, x0, tol)4. 弦截法求方程根(1)编写弦截法函数:function [root, error] = secant(f, x0, x1, tol)(2)输入函数f、初始值x0和x1,以及精度要求tol(3)调用函数计算根:[root, error] = secant(f, x0, x1, tol)五、实验结果与分析1. 二分法求方程根以方程f(x) = x^2 - 2 = 0为例,输入初始区间[a, b]为[1, 3],精度要求tol 为1e-6。