MATLAB应用第五章MATLAB数值计算●——MATLAB强大的数值计算能力，使其成为在科学计算方面的首先解题工具。
本章主要内容●5.1 数据处理与多项式计算●5.2 数值微积分●5.3 离散傅里叶变换●5.4 线性方程组求解●5.5 非线性方程组求解●5.6 常微分方程求解●5.7 稀疏矩阵5.1 数据处理与多项式计算●一、数据统计与分析●二、多项式计算●三、曲线拟合●四、数据插值一、数据统计与分析●1. 求矩阵最大元素和最小元素●max●min●语法：C = max(A)C = max(A,B)C = max(A,[],dim)[C,I] = max(...)●例求矩阵A的每行及每列的最大和最小元素，并求整个矩阵的最大和最小元素。
●2.求矩阵的平均值和中值●mean●median●语法：M = mean(A)M = mean(A,dim)●3. 矩阵元素求和与求积●sum●prod●语法：●B = sum(A)●B = sum(A, dim)●4. 标准方差●std●语法:s = std(X)s = std(X,flag)s = std(X,flag,dim)●说明：●对于具有N个元素的数据序列（x1,x2,…,x N)，标准方差计算公式如下：●或●其中●flag取0或1，其取值决定计算标准方差时所使用的公式。
●flag=0，按S1计算；●flag=1，按S2计算。
●5. 求元素的累加和与累乘积●cumsum●cumprod●语法：●B = cumsum(A)●B = cumsum(A,dim)●参考P142●6. 相关系数●corrcoef●可求出两组数据的相关系数。
●例5.1数据统计函数使用举例●（1 ）某城市的3 个路口24 小时的车流量数据。
>>load count（2 ）P144 例6.5●7. 元素排序●sort(x)●返回一个对X 中的元素按升序排列的新向量。
二、多项式计算●1.多项式的表示●2.多项式的四则运算●3.多项式的导函数●4.多项式求值●5.多项式求根●1.多项式的表示●在MATLAB 中，n 次多项式用一个长度为n+1 的行向量表示，缺少的幂次项系数为0.●如果n 次多项式表示为：●则在MATLAB中,p(x)表示为向量形式:●2. 多项式的四则运算●(1) 多项式的加减运算●(2) 多项式的乘法运算( 卷积)●conv(P1,P2 )●——求多项式P1 ，P2 的乘积●(3) 多项式除法( 去卷积)●[Q,r]=deconv(P1,P2)●——P1 除以P2 ，Q 为商式，r 为余式。
●deconv 是conv 的逆函数，即P1=conv(p2,Q)+r●例5.2 设a(s)=s2+2s+3,b(s)=4s2+5s+6，则求c(s)=a(s)*b(s).●源代码见examp502.m●2.多项式的导函数●polyder●语法：●p=polyder(P)▪求多项式P 的导函数●p=polyder(P,Q)▪求多项式P·Q 的导函数●[p,q]=polyder(P,Q)▪求多项式P/Q 的导函数，导函数的分子存入p ，分母存入q●3. 多项式求值●polyval●代数多项式求值●polyvalm●矩阵多项式求值●语法：●polyval(P,x)●polyvalm(P,x)●例5.3以多项式x4+8x3-10为例，取一个2×2矩阵为自变量分别用polyval和polyvalm计算该多项式的值。
●见examp503.m●4. 多项式求根——roots●roots●语法：●x=roots(P)●若已知多项式的全部根，则可用poly 函数建立多项式。
●语法：●P=poly(x)●例5.4 已知(1) 计算f(x)=0 的全部根。
(2) 由方程f(x)=0 的根构造一个多项式g(x) ，并与f(x)进行对比。
源代码见examp504.m●5. 曲线拟合●曲线拟合是进行数据分析时经常遇到的问题，它是指根据一组或多组测量数据找出数学上可以描述此数据走向的一条曲线的过程。
评价一条曲线是否准确地描述了测量数据的最通用方法，是看测量数据点与该曲线上对应点之间的平方误差是否达到最小，这种曲线拟合的方法称为最小二乘曲线拟合。
●在MATLAB 中用polyfit 实现●语法：●[P,S]=polyfit(X,Y,m)●X,Y 是两个等长的向量●P 是一个长度为m+1 的向量，●P 的元素为多项式系数●S 为采样点的误差向量●例5.5用一个3次多项式在区间[0，2π]内逼近函数sin x。
●源代码见examp505.m6.2 数据插值●为了通过离散的数据点来获得更为丰富的信息，就必须对数据进行插值，插值本身也是数值计算中常用到的处理方法，MATLAB对此作了充分的考虑。
●测量得n各点(x1,y1)、(x2,y2)、…、(x n,y n)，这些点反映一个函数关系y=f(x),●然而，f(x)的解析式不知道。
●数据插值的任务就是根据上述条件构造一个函数y=g(x),使得对于x i(i=1,2,…,n),有g(x i)=f(x i)，且在相邻的采样点（x i,x i+1)(i=1,2,…,n-1)，g(x)光滑过渡。
●如果被插值函数f(x)是光滑的，且采样点足够密，一般在采样区间内，f(x)与g(x)比较接近。
●MATLAB提供了一维、二维、三维、N维数据插值函数interp1、interp2、interp3、interpn，以及三次样条插值函数spline等。
●1.一维数据插值——interp1●调用格式：●Y1=interp1(X,Y,X1,method)●说明：●X ，Y 是两个等长的已知向量，分别描述采样点和样本值；●X1 是一个向量或标量，描述欲插值的点。
●Y1 是一个与X1 等长的插值结果。
●method是插值方法，有四种：●…li near‟：线性插值●…nearest‟：最近点插值●…cubic‟：三次多项式插值●…spline‟：3次样条插值●注意：●X1 的取值范围不能超出X 的给定范围，否则会出现“NaN”错误。
●例5.6 在[-2，2]区间上分别用4中插值方法来绘制y=x2的图形。
●源代码见examp506.m●2. 二维数据插值——interp2●调用格式:●Z1=interp2(X,Y,Z,X1,Y1,method)●说明：●X,Y 分别描述两个参数的采样点，Z 是与采样点对应的函数值。
●X1 ，Y1 是两个向量或标量，描述欲插值的点。
●Z1 是根据相应的插值方法得到的插值结果。
●method 的取值与一维插值函数相同●例5.7 某实验对一根长10米的钢轨进行热源的温度传播测试。
用x表示测量点(米)，用h表示测量时间(秒)，用T表示测得各点的温度(℃)，测量结果如下表所示。
●试用用3次多项式插值求出在一分钟内每隔10秒、钢轨每隔0.5米处的温度。
二、数值微积分●1. 数值微分●(1) 数值差分与差商●任意函数f(x)在x点的导数可定义为：●当步长h（h>0）充分小时，则有●2. 数值微分的实现●(1) 多项式求导求数值微分●用多项式或样条函数g(x) 对f(x) 进行逼近（插值或拟合），然后用g(x) 在点x 处的导数作为f(x) 在x 处的导数。
●例5.8用5阶多项式拟合函数cos(x)，并利用多项式的求导来求π处的一阶和二阶导数。
（源代码见examp508.m）●2. 用diff计算差分求数值微分●语法：●Y = diff(X)●Y = diff(X,n)●Y = diff(X,n,dim)●例5.9 设：用不同的方法求函数f(x)的数值导数，并在同一坐标系中画出f‟(x)的图像。
解题思路：方法1：用一个5次多项式p(x)拟合函数f(x),对p(x)求导。
方法2：直接求f(x)在假设点的数值导数方法3：求出f‟(x),然后直接求f‟(x)在假设点的导数源代码见examp509.m●2.数值积分的实现●(1) 被积函数是一个解析式●使用quad和quadl函数来求定积分●语法：▪quad(filename,a,b,tol,trace)▪quadl(filename,a,b,tol,trace)●例5.10 用两种不同的方法求：先建立一个函数文件ex.m：function ex=ex(x)ex=exp(-x.^2);然后在MATLAB 命令窗口，输入命令: format longI=quad('ex',0,1)I=quadl('ex',0,1)使用内联函数求解，命令如下：g=inline(…exp(-x.^2)‟);I=quadl(g,0,1)●(2) 被积函数有一个表格定义●使用trapz 求数值积分●语法：●trapz(X,Y)●(3)二重积分数值求解●dblquard●设被积函数为：●语法：●I=dblquad(f,a,b,c,d,tol,trace)6.3 离散傅里叶变换（DFT）●一. 离散傅里叶变换算法简述●在某时间片等距地抽取N 个抽样时间t m处的样本值f(t m) ，这里m=0,1,2,…,N-1, 称向量F(k)(k=0,1,2,…,N-1) 为f(m) 的一个离散傅里叶变换，其中：●由于MATLAB不允许有零下标，将上式m下标均移1，得：●由f(m)求F(k)的过程，称为求f(m)的离散傅里叶变换，又称F(k)为f(m)的离散频谱。
●反之，由F(k)逆求f(m)的过程，称为离散傅里叶逆变换。
二、离散傅里叶变换的实现●一维离散傅里叶变换函数——fft●语法：●Y = fft(X)●Y = fft(X,n)●Y = fft(X,[],dim)●Y = fft(X,n,dim)●例5.11 给定数学函数x(t)=12sin(2π×10t+π/4)+5cos(2π×40t)取N=128，试对t从0~1秒采样，用FFT作快速傅立叶变换，绘制相应的振幅-频率图。
5.4 线性方程组求解●直接解法●X=A\B●X=inv(A)*B5.5 非线性方程与最优化问题求解●略5.6 常微分方程的数值求解●求解常微分方程的函数：●ode45●ode23●一、微分方程的求解过程●1. 将微分方程变换成一阶微分方程，及表示成右函数形式。
（这是利用龙格-库塔法求解微分方程的前提）●例如：有微分方程●若令●则可得到一阶微分方程●相应，可确定y1(0),y2(0),…,y n(0)●2. 将一阶微分方程组编写成M文件，形成一个函数，设为funf(t,y)●function dy=funf(t,y)●dy=[y(2);y(3);…;f(t,y(1),y(2),…,y(n-1)]●注：实际编写程序时不能用省略号代替，应根据系统方程书写；有时方程中随不出现t ，但为了方便ode45 等函数调用，必须将t 用作自变量●3.利用MATLAB提供的函数求解微分方程,以ode45为例：●[t,y]=ode45(…funf‟,tspan,y0)●例设有初值问题：试求其数值解，并与精确解相比较(精确解为)。