数值分析(王能超)matlab上机作业
数值分析上机作业(MATLAB)
将系数矩阵 A 分解为：A=L+U+D
Ax=b
⇔ (D + L +U)x = b ⇔ Dx = −(L + U )x + b ⇔ x = −D −1(L + U )x + D −1b x(k +1) = −D −1 (L + U ) x(k ) + D −1b
输入 A，b 和初始向量 x
迭代矩阵 BJ , BG
否
ρ(B) < 1?
按雅各比方法进行迭代
否
|| x (k+1) − x(k) ||< ε ?
按高斯-塞德尔法进行迭代
否
|| x(k+1) − x (k ) ||< ε ?
输出迭代结果
图 1 雅各布和高斯-赛德尔算法程序流程图
1.2 问题求解
按图 1 所示的程序流程，用 MATLAB 编写程序代码，具体见附录 1。解上述三个问题 如下
16
-0.72723528355328
0.80813484897616
0.25249261987171
17
-0.72729617968010
0.80805513082418
0.25253982509100
18
-0.72726173942623
0.80809395746552
0.25251408253388
0.80756312717373
8
-0.72715363032573
0.80789064377799
9
-0.72718652854079
数值分析上机题3
数值分析上机题目3实验一1．根据Matlab 语言特点，描述Jacobi 迭代法、Gauss-Seidel 迭代法和SOR 迭代法。
2．编写Jacobi 迭代法、Gauss-Seidel 迭代法和SOR 迭代法的M 文件。
3．给定2020⨯∈R A 为五对角矩阵⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣⎡------------------321412132141412132141412132141412132141213O O O O (1)选取不同的初始向量)0(x 及右端面项向量b ，给定迭代误差要求，分别用编写的Jacobi 迭代法和Gauss-Seidel 迭代法程序求解，观察得到的序列是否收敛？若收敛，通过迭代次数分析计算结果并得出你的结论。
(2)用编写的SOR 迭代法程序， 对于(1)所选取的初始向量)0(x 及右端面项向量b 进行求解，松驰系数ω取1<ω<2的不同值，在5)1()(10-+≤-k k x x 时停止迭代，通过迭代次数分析计算结果并得出你的结论。
实验11、 根据MATLAB 语言特点，描述Jacobi 迭代法，Gauss-Seidel 迭代法和SOR 迭代法。
2、 编写Jacobi 迭代法，Gauss-Seidel 迭代法和SOR 迭代法的M 文件。
Jacobi 迭代法function [x1,k]=GS_2(A,b)n=length(A);D=diag(diag(A));L=-tril(A,-1);U=-triu(A,1);x1=zeros(n,1);x0=3*ones(n,1);k=0; while norm(x1-x0,1)>10^(-7)&k<100 k=k+1;x0=x1;x1=D\((L+U)*x0+b);endk=kx=x1Gauss-Seidel迭代法function [x1,k]=GS_h(A,b)n=length(A);D=diag(diag(A));L=-tril(A,-1);U=-triu(A,1);x1=zeros(n,1);x0=3*ones(n,1);k=0; while norm(x1-x0,1)>10^(-7)&k<100 k=k+1;x0=x1;x1=(D-L)\U*x0-D\b;endk=kx=x1SOR迭代法function [x1,k]=SOR_h(A,b)n=length(A);D=diag(diag(A));L=-tril(A,-1);U=-triu(A,1);x1=zeros(n,1);x0=3*ones(n,1);k=0;w=0.96;while norm(x1-x0,1)>10^(-7)&k<100k=k+1;x0=x1;x1=(D-w*U)\(((1-w)*D+w*L)*x0+w*b);endk=kx=x13、采用Jacobi迭代法，Gauss-Seidel迭代法求解五对角矩阵clear,clcA=diag(3*ones(20,1))+diag((-0.5)*ones(19,1),-1)+diag((-0.5)*ones(19,1 ),1)+diag((-0.25)*ones(18,1),-2)+diag((-0.25)*ones(18,1),2);b=sum(A')';[x1,k1]=Jacob_h(A,b)[x2,k2]=GS_h(A,b)运行结果：两种方法都收敛，k1=27,k2=13。
数值计算上机实习题目（matlab编程）
数值计算上机实习题目（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进行数值计算的基本方法。
在具体实验中，我们将实现三种常见的数值分析算法：二分法、牛顿法和追赶法，分别应用于解决非线性方程、方程组和线性方程组的求解问题。
二、实验原理与方法1.二分法二分法是一种常见的求解非线性方程的数值方法。
根据函数在给定区间端点处的函数值的符号，不断缩小区间的长度，直到满足精度要求。
2.牛顿法牛顿法是求解方程的一种迭代方法，通过构造方程的泰勒展开式进行近似求解。
根据泰勒展式可以得到迭代公式，利用迭代公式不断逼近方程的解。
3.追赶法追赶法是用于求解三对角线性方程组的一种直接求解方法。
通过构造追赶矩阵，采用较为简便的向前追赶和向后追赶的方法进行计算。
本次实验中，我们选择了一组非线性方程、方程组和线性方程组进行求解。
具体的实验步骤如下：1.调用二分法函数，通过输入给定区间的上下界、截止误差和最大迭代次数，得到非线性方程的数值解。
2.调用牛顿法函数，通过输入初始迭代点、截止误差和最大迭代次数，得到方程组的数值解。
3.调用追赶法函数，通过输入追赶矩阵的三个向量与结果向量，得到线性方程组的数值解。
三、实验结果与分析在进行实验过程中，我们分别给定了不同的参数，通过调用相应的函数得到了实验结果。
下面是实验结果的汇总及分析。
1.非线性方程的数值解我们通过使用二分法对非线性方程进行求解，给定了区间的上下界、截止误差和最大迭代次数。
实验结果显示，根据给定的输入，我们得到了方程的数值解。
通过与解析解进行比较，可以发现二分法得到的数值解与解析解的误差在可接受范围内，说明二分法是有效的。
2.方程组的数值解我们通过使用牛顿法对方程组进行求解，给定了初始迭代点、截止误差和最大迭代次数。
实验结果显示，根据给定的输入，我们得到了方程组的数值解。
与解析解进行比较，同样可以发现牛顿法得到的数值解与解析解的误差在可接受范围内，说明牛顿法是有效的。
数值分析matlab上机实验报告
数值分析matlab上机实验报告matlab软件实验报告数学上机课实验报告matlab实验报告总结数值分析试卷篇一：《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) ); y(n+1) = b*x(n)+a*(y(n)-x(n) )2. 编程实现奥运5环图，允许用户输入环的直径。
3. 实现对输入任意长度向量元素的冒泡排序的升序排列。
不允许使用sort函数。
四、实验数据及结果分析题目一：①在Editor窗口编写函数代码如下：并将编写的函数文件用“draw.m”储存在指定地址；②在Command窗口输入如下命令：③得到图形结果如下：题目二：①在Editor窗口编写函数代码如下：并将编写的函数文件用“circle.m”储存在指定地址；②再次在Editor窗口编写代码：并将编写的函数文件用“Olympic.m”储存在指定地址；③在Command窗口输入如下指令（半径可任意输入）：④按回车执行，将在图形窗口获得五环旗：题目三：①在Editor窗口编写函数代码如下：并用.将编写的函数文件用“qipaofa.m”储存在指定地址；②在Command窗口输入一组乱序数值，则可以得到升序排序结果如下：五、总结及心得体会1. 要熟悉MATLAB编译软件的使用方法，明白有关语法，语句的基本用法，才可以在编写程序的时候游刃有余，不至于寸步难行。
数值分析作业-matlab上机作业
数值分析———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上机作业报告
一、给定向量x ≠0，计算初等反射阵H k 。
1．程序功能：给定向量x ≠0，计算初等反射阵H k 。
2．基本原理： 若()xx x R x ∈=,,, 的分量不全为零，则由12112212122()x (,,,)1()22n T T sign x e x x x x σσσρσσρ-⎧=⎪=+=+⎪⎪⎪==+⎨⎪⎪=-=-⎪⎪⎩u x u uu H I I uuu 确定的镜面反射阵H 使得y e Hx =-=σ；当(1)k n ≤<时，由21/2k ()T 1()()()k 1()()()(())(0,,0,,,,)1()()=()2()nk i i kk nk k k n k T k k Tk k k kk k T k k sign x x x x x x σσρσσσρ=+-⎧=⎪⎪⎪=+∈⎨⎪==+⎪⎪=-⎩∑u R u u u H I u u 有T 121(,,,,,0,,0)n k k k x x x σ-=-∈H x R 算法：(1)输入x ，若x 为零向量，则报错 (2)将x 规范化，{}x x x M ,,,max =如果M ＝0，则报错同时转出停机 否则n i M x x i i ,,2,1, =←(3)计算2x =σ，如果0<1x ，则σσ-= (4))(1x +=σσρ (5)计算1,(1)x σ==+u x u (6)1Tρ-=-H I uu (7)(M ,0,,0)σ=-y(8)按要求输出,结束3．变量说明：x －输入的n维向量；n －n维向量x的维数；M －M是向量x的无穷范数，即x中绝对值最大的一项的绝对值；p －Householder初等变换阵的系数ρ；u －Householder初等变换阵的向量Us －向量x的二范数；x －输入的n维向量；n －n维向量x的维数；p －Householder初等变换阵的系数ρ；u －Householder初等变换阵的向量Uk －数k，H*x=y，使得y的第k+1项到最后项全为零；4．程序代码：(1)function [p,u]=holder2(x)%HOLDER2 给定向量x≠0,计算Householder初等变换阵的p,u%程序功能：函数holder2给定向量x≠0,计算Householder初等变换阵的p,u；%输入：n维向量x；%输出：[p,u]。
数值分析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方法运行结果：代码截图：。
《数值分析》上机实验报告
数值分析上机实验报告《数值分析》上机实验报告1.用Newton 法求方程 X 7-X 4+14=0在（0.1，1.9）中的近似根（初始近似值取为区间端点，迭代6次或误差小于0.00001）。
1.1 理论依据：设函数在有限区间[a ，b]上二阶导数存在，且满足条件{}αϕ上的惟一解在区间平方收敛于方程所生的迭代序列迭代过程由则对任意初始近似值达到的一个中使是其中上不变号在区间],[0)(3,2,1,0,)(')()(],,[x |))(),((|,|,)(||)(|.4;0)(.3],[)(.20)()(.110......b a x f x k x f x f x x x Newton b a b f a f mir b a c x f ab c f x f b a x f b f x f k k k k k k ==-==∈≤-≠>+令)9.1()9.1(0)8(4233642)(0)16(71127)(0)9.1(,0)1.0(,1428)(3225333647>⋅''<-=-=''<-=-='<>+-=f f x x x x x f x x x x x f f f x x x f故以1.9为起点⎪⎩⎪⎨⎧='-=+9.1)()(01x x f x f x x k k k k 如此一次一次的迭代，逼近x 的真实根。
当前后两个的差<=ε时，就认为求出了近似的根。
本程序用Newton 法求代数方程(最高次数不大于10)在（a,b ）区间的根。
1.2 C语言程序原代码：#include<stdio.h>#include<math.h>main(){double x2,f,f1;double x1=1.9; //取初值为1.9do{x2=x1;f=pow(x2,7)-28*pow(x2,4)+14;f1=7*pow(x2,6)-4*28*pow(x2,3);x1=x2-f/f1;}while(fabs(x1-x2)>=0.00001||x1<0.1); //限制循环次数printf("计算结果：x=%f\n",x1);}1.3 运行结果：1.4 MATLAB上机程序function y=Newton(f,df,x0,eps,M)d=0;for k=1:Mif feval(df,x0)==0d=2;breakelsex1=x0-feval(f,x0)/feval(df,x0);ende=abs(x1-x0);x0=x1;if e<=eps&&abs(feval(f,x1))<=epsd=1;breakendendif d==1y=x1;elseif d==0y='迭代M次失败';elsey= '奇异'endfunction y=df(x)y=7*x^6-28*4*x^3;Endfunction y=f(x)y=x^7-28*x^4+14;End>> x0=1.9;>> eps=0.00001;>> M=100;>> x=Newton('f','df',x0,eps,M);>> vpa(x,7)1.5 问题讨论：1.使用此方法求方解，用误差来控制循环迭代次数，可以在误差允许的范围内得到比较理想的计算结果。
数值分析上机作业
第二次上机作业一. 任务：用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中，可以通过多项式逼近、三次样条、拉格朗日插值、最小二乘逼近等方法来实现数值逼近的计算。
例如，使用MATLAB的插值函数interp1可以实现一维函数的插值计算，使用MATLAB 的polyfit函数可以拟合一维数据集合的多项式曲线。
而对于二维函数和三维函数的逼近，可以使用MATLAB的interp2和interp3函数来实现。
数值微分和数值积分是数值分析中的基本操作之一、它们可以根据给定的函数计算函数的导数和积分。
在MATLAB中，使用diff函数可以计算一维函数的导数，使用trapz和quad函数可以计算一维函数的定积分和数值积分。
而对于二维函数和三维函数的微分和积分，可以使用MATLAB 的grad函数和integral2函数来实现。
此外，MATLAB还提供了很多其他的数学函数和工具，包括解微分方程、优化问题、曲线拟合和最小二乘等。
对于一些复杂的数学问题，可以通过使用MATLAB的符号计算工具箱来实现符号计算。
第四章数值分析上机作业
1.用两种方法解sin 5cos 20.31x x +=-，如何一次求出此方程的四个根。
解：方法一：在Ｍatlab 命令窗口输入命令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 ２.用三种方法解方程11852912312x x x x -+-=。
解：方法一：在Ｍatlab 命令窗口输入命令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 方法二：在Ｍatlab 命令窗口输入命令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方法三：在Ｍatlab命令窗口输入命令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*３ 用MATLAB 方法求函数11852()912312f x x x x x =-+--的导数'()f x 。
数值分析大作业(利用MATLAB软件)
实验报告课程名称：数值分析实验项目：曲线拟合/数值积分专业班级：姓名：学号：实验室号：实验组号：实验时间：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。
《MATLAB及其应用》上机作业.doc(修订于2009.11.19)
《MATLAB及其应用》上机作业学院名称：（四号宋体）专业班级：（四号宋体）学生姓名：（四号宋体）学生学号：（四号宋体）年月作业11.用MATLAB 可以识别的格式输入下面两个矩阵12342357135732391894A ⎛⎫⎪ ⎪ ⎪= ⎪ ⎪ ⎪⎝⎭,144367723355422675342189543i i B i +⎛⎫ ⎪+⎪= ⎪+ ⎪⎪⎝⎭再求出它们的乘积矩阵C ,并将C 矩阵的右下角23⨯子矩阵赋给D 矩阵。
赋值完成后,调用相应的命令查看MATLAB 工作空间的占有情况。
2.设矩阵16231351110897612414152⎛⎫⎪⎪ ⎪ ⎪⎝⎭,求A ,1A -,3A ,12A A -+,1'3A A --,并求矩阵A 的特征值和特征向量。
3.解下列矩阵方程:010100143100001201001010120X -⎛⎫⎛⎫⎛⎫ ⎪ ⎪ ⎪=- ⎪ ⎪ ⎪ ⎪ ⎪ ⎪-⎝⎭⎝⎭⎝⎭4.一球从100米高度自由落下,每次落地后反跳回原高度的一半,再落下。
求它在第10次落地时,共经过多少米?第10次反弹有多高?5.用MATLAB 语言实现下面的分段函数5,1(),25,y f x x ⎧⎪⎪==⎨⎪-⎪⎩101010x x x >≤<-6.分别用for 和while 循环编写程序,求出6323626312122222i i K ===++++++∑并考虑一种避免循环的简洁方法来进行求和,并比较各种算法的运行时间。
7.应用MA TLAB 语言及二分法编写求解一元方程32()1459700f x x x x =-+-=在区间[3，6]的实数解的算法，要求绝对误差不超过0.001。
8.二阶系统的单位阶跃响应为1cos )t y a ζζ-=+，在同一平面绘制ζ分别为0，0.3，0.5，0.707的单位阶跃响应曲线。
要求：(1)四条曲线的颜色分别为蓝、绿、红、黄，线型分别为“——”、“……”、“oooooo”、“++++++”； (2)添加横坐标轴和纵坐标轴名分别为“时间t”和“响应y”，并在平面图上添加标题“二阶系统曲线”和网格；(3)在右上角添加图例（即用对应的字符串区分图形上的线），并分别在对应的响应曲线的第一个峰值处标示“zeta ＝0”、“zeta ＝0.3”、“zeta ＝0.5”、“zeta ＝0.707”。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
引论习题题目：用二分法求方程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)求得结果如下（整理之后）：其中红色数据表示插值节点将以上结果绘制成图，如下所示:习题二题目编写一通用型复化辛普生公式，能够对任意长度的等间距离散数据进行积分运算。
要求：命令格式：y=simpson(x,h);y：积分结果x：输入数组h：间距程序：function [ y ] = simpson( x,h )% 一通用型复化辛普森求积公式% y:积分结果；x:输入数组；h:间距if nargin<2 disp('输入参数不够，请准确输入积分数据和步长');return;endn=length(x);% X=zeros(n+1);if n==1 y=x;return;else if mod(n,2)==1% X=zeros(n+1);A=zeros((n+1)/2);B=zeros((n-1)/2);A=x(1:2:(end-2));B=x(2:2:(end-1));y=2*h/6*(x(end)-x(1)+4*sum(B)+2*sum(A));% X(1:n)=x;else mod(n,2)==1X=zeros(n-1);X=x(1:(end-1));m=length(X);A=zeros((m+1)/2);B=zeros((m-1)/2);A=X(1:2:(end-2));B=X(2:2:(end-1));y=2*h/6*(X(end)-X(1)+4*sum(B)+2*sum(A))+(x(n)+x(n-1))/2*h;endendend运算结果：以计算10sin()xI dxx=⎰为例在commond window中运行结果如下：clear>> t=0:(1/8):1;>> x=sin(t)./t;>> h=1/8;>> [ y ] = simpson( x,h )y =0.9461以上积分结果与教材P65所给结果完全一致习题三题目分别用显式和隐式的二阶亚当姆斯方法求解初值问题y'=1-y，y(0)=0，令y(0.2)=0.181，取h=0.2计算y(1.0)。
程序：1、显式亚当姆斯方法：function [ X,Y ] = Adams_shown( x0,y0,y1,h,N )%显式二阶亚当姆斯求解微分方程%x0，y0,y1为输入初值，h为步长，N为计算点数（包括除x0,y0在内的输入初值）%X,Y分别为输出数组,f为y'=f(x,y)X=zeros(N+1,1);Y=zeros(N+1,1);%为X,Y预置空间X=linspace(x0,x0+N*h,N+1);X=X(:);Y(1)=y0;Y(2)=y1;for n=2:Nyy0=f(x0,y0);yy1=f(x0+h,y1);y2=y1+h/2*(3*yy1-yy0);yy2=f(x0+2*h,y2);Y(n+1)=y2;x0=x0+h;y0=y1;y1=y2;endend2、隐式亚当姆斯方法：function [ X,Y ] = Adams_hidden( x0,y0,h,N )%隐式二阶亚当姆斯求解微分方程%x0，y0为输入初值，h为步长，N为计算点数（不包括输入初值）%X,Y分别为输出数组,f为y'=f(x,y)X=zeros(N+1,1);Y=zeros(N+1,1);%为X,Y预置空间X=linspace(x0,x0+N*h,N+1);X=X(:);Y(1)=y0;for n=1:Nyy0=f(x0,y0);% y1=solve('y0+h/2*(f(x0+h,y1)+yy0)-y1','y1');y1=(y0+h/2*(1+yy0))/(1+h/2);yy1=f(x0+h,y1);Y(n+1)=y1;x0=x0+h;y0=y1;endend3、二阶亚当姆斯预报—校正方法：function [ X,Y ] = Adams_shown_hidden( x0,y0,y1,h,N )%二阶亚当姆斯预报-校正求解微分方程%x0，y0,y1为输入初值，h为步长，N为计算点数（包括除x0,y0在内的输入初值）%X,Y分别为输出数组,f为y'=f(x,y)X=zeros(N+1,1);Y=zeros(N+1,1);%为X,Y预置空间X=linspace(x0,x0+N*h,N+1);X=X(:);Y(1)=y0;Y(2)=y1;for n=2:Nyy0=f(x0,y0);yy1=f(x0+h,y1);y2=y1+h/2*(3*yy1-yy0);yy2=f(x0+2*h,y2);y2=y1+h/2*(yy2+yy1);yy2=f(x0+2*h,y2);Y(n+1)=y2;x0=x0+h;y0=y1;y1=y2;endend运算结果：定义函数y'=1-y如下：function [yy] = f(x,y)%定义了测试函数% x，y为输入变量yy=1-y;end在commond window中运行结果如下：>> clear>> x0=0;y0=0;y1=0.181;h=0.2;N=5;>> [ X1,Y1 ] = Adams_shown( x0,y0,y1,h,N );>> [ X2,Y2 ] = Adams_hidden( x0,y0,h,N );>> [ X3,Y3 ] = Adams_shown_hidden( x0,y0,y1,h,N );>> plot(X1,Y1,'r',X2,Y2,'b',X3,Y3,'g');>> legend('显式亚当姆斯','隐式亚当姆斯',,'二阶亚当姆斯校正-预报方法') 最后所得结果如下表所示：根据上述数据绘制成图，如下所示：习题四题目 用牛顿法求下列方程的根，要求计算结果小数点后有4位有效数字（1）30310,2x x x --== （2）20320,1x x x e x --+==程序：function root=NewtonRoot(f,x0,eps)%牛顿法求方程的根% x0为迭代初值%f 为方程表达式：f(x)=0%eps 为迭代精度%root 为满足精度要求的近似根if(nargin==2)eps=1.0e-4;endtol=1;fx=sym(f);dfx=diff(sym(f));%求导数while(tol>eps)fx0=subs(fx,findsym(fx),x0);dfx0=subs(dfx,findsym(dfx),x0); %求该点的导数值 root=x0-fx0/dfx0; %迭代的核心公式 tol=abs(root-x0);x0=root;endformat long g;root=eval(root);end运算结果：在commond window 中运行结果如下：（1）30310,2x x x --== >> clear>> root=NewtonRoot('x^3-3*x-1',2,1.0e-4)root =1.87938524483667（2）20320,1x x x e x --+== >> root=NewtonRoot('x^2-3*x-exp(x)+2',1,1.0e-4)root =0.257530285426349。