(高等教育)大学数学实验基础知识整理(清华大学)

Contents差分方程和数值微分实验41.1 差分方程的基本定义41.2 一阶线性常系数差分方程41.3高阶线性常系数差分方程41.4 线性常系数差分方程组51.5 非线性差分方程52 数值微分5插值与数值积分61 插值与拟合61.1 插值与拟合的基本概念61.2 三种插值方法62 数值积分82.1 数值积分的基本思路82.2 三种常用数值积分方法8常微分方程数值解10常微分方程的初值问题102.初值问题的数值解法102.1 欧拉方法102.2 龙格-库塔方法11常微分方程组和高阶方程初值问题的数值方法112.3 龙格-库塔方法的MATLAB实现122.4 算法的收敛性、稳定性分析12刚性现象与刚性方程13线性代数方程组数值解法13线性代数方程组的一般形式和解法132.求解线性代数方程组的直接法132.1 高斯消元法132.2 LU分解142.3 解的误差分析P95143.求解线性代数方程组的迭代法153.1 雅可比迭代法153.2 高斯-赛德尔迭代法153.3 迭代法的收敛性和收敛速度153.4 超松弛迭代164.超定线性代数方程组的最小二乘解164.1 超定线性方程组的概念164.2 最小二乘准则164.3 最小二乘解164.4 基函数的选取17MATLAB实现17非线性方程求解171 非线性方程(组)的定义及特点172 非线性方程的基本解法182.1 图形法和二分法182.2 迭代法182.3 牛顿法193 非线性方程组的牛顿法、拟牛顿法194 用MATLAB工具箱解非线性方程(组)204.1 fzero的基本用法204.2 fsolve的基本用法214.3 roots的基本用法22无约束优化231无约束优化的基本原理、解法231.1 无约束优化的一般形式231.2 最优性条件231.3 下降法的基本思想231.4 用MATLAB优化工具箱解无约束优化问题232非线性最小二乘拟合的基本原理、解法252.1 非线性最小二乘拟合问题252.2 非线性最小二乘拟合问题的解法252.3 用MATLAB优化工具箱解非线性最小二乘拟合问题26约束优化2711.线性规划的基本原理、解法281.1 线性规划的图解法281.2 线性规划的标准形281.3基本可行解281.4 线性规划的基本性质281.5 单纯形法的基本思路281.6 线性规划解的几种可能291.7 用MATLAB优化工具包解线性规划292.非线性规划的基本原理、解法312.1 非线性规划的一般形式312.2 可行方向与下降方向312.3 最优解的必要条件312.4 二次规划的一般形式322.5 二次规划的有效集方法322.6 用MATLAB优化工具包解二次规划332.7 非线性规划的解法342.8 用MATLAB优化工具包解非线性规划34数据的统计与分析361 统计的基本概念362 频数表和直方图373 统计量374 统计中几个重要的概率分布384.1 分布函数、密度函数和分位数384.2 统计中几个重要的概率分布384.3 MATLAB统计工具箱(ToolboxStats)中的概率分布 P246395 正态总体统计量的分布396. 用随机模拟计算数值积分406.1两种方法406.2重积分的计算406.3MATLAB实现40统计推断401、参数估计40概述401.1 点估计411.2 点估计的评价标准411.3 总体均值的区间估计421.4 总体方差的区间估计441.5 参数估计的MATLAB实现442、假设检验45概述452.1 均值的假设检验452.2 方差（或标准差）的假设检验462.3 两总体的假设检验462.4 0-1分布总体均值的假设检验472.5 总体分布正态性检验472.6 假设检验与Matlab命令汇总49差分方程和数值微分实验1.1 差分方程的基本定义差分方程是在离散时段上描述现实世界中变化过程的数学模型。现实中的问题通常是连续变化的，但我们常常只能在离散的时间点上对其进行观测和描述。为了表述这一类的数学模型，我们引入了差分方程的方法。1.2 一阶线性常系数差分方程 一阶线性常系数差分方程的一般形式 差分方程的平衡点 差分方程的解 平衡点稳定的条件 1.3高阶线性常系数差分方程高阶线性常系数差分方程的一般形式 特征方程 特征根 平衡点 差分方程的解 平衡点稳定的条件所有特征值的模均小于1 (用roots(c)-c:多项式的系数（降幂）P125)1.4 线性常系数差分方程组当我们研究的对象是若干变量构成的一个向量的离散动态过程时，就需要引入差分方程组来描述，详见前面对一阶或高阶线性常系数差分方程的描述。 平衡点X=Ax+b稳定条件：A的所有特征根小于1（eig）1.5 非线性差分方程2 数值微分数值微分是用离散方法近似地计算函数y=f(x)在某点x=a的导数值。常用公式有： 前差公式 后差公式 中点公式 三点公式 插值与数值积分1 插值与拟合1.1 插值与拟合的基本概念 插值与插值函数：已知由 (可能未知或非常复杂)产生的一批离散数据 ，且 个互异插值节点 ，在插值区间内寻找一个相对简单的函数 ，使其满足下列插值条件： 再利用已求得的 计算任一非插值节点 的近似值 ，这就是插值。其中 称为插值函数， 称为被插函数。 最小二乘拟合： 已知一批离散的数据 ， 互不相同，寻求一个拟合函数 ，使 与 的误差平方和在最小二乘意义下最小。在最小二乘意义下确定的 称为最小二乘拟合函数。1.2 三种插值方法1）Lagrange插值法 a待定系数法： 假设插值多项式 ，利用待定系数法即可求得满足插值条件 的插值函数。关键在于确定待定系数 。 b利用基函数的构造方法 首先构造 个满足条件： 的 次插值基函数 ，再将其线性组合即可得如下的Lagrange插值多项式： 其中 cLagrange插值余项 注：上述两种构造方法所得的Lagrange插值多项式是一样的，即满足插值条件 的Lagrange插值多项式是唯一的。Lagrange插值会发生Runge现象。2）分段线性插值 作分段线性插值的目的在于克服Lagrange插值方法可能发生的不收敛性缺点。所谓分段线性插值就是利用每两个相邻插值节点作线性插值，即可得如下分段线性插值函数： 其中 特点：插值函数序列 具有一致收敛性，克服了高次Lagrange插值方法的缺点，故可通过增加插值节点的方法提高其插值精度。但存在于节点处不光滑、插值精度低的缺点。 3）三次样条插值 三次样条插值的目的在于克服Lagrange插值的不收敛性和提高分段线性插值函数在节点处的光滑性。所谓三次样条插值方法就是在满足下列条件： a b 在每个子区间 上是三次多项式的三次样条函数中寻找满足如下插值条件： 一及形如 等边界条件的插值函数 的方法。 特点：三次样条插值函数序列 一致收敛于被插函数，因此可通过增加节点的方法提高插值的精度。4）插值方法的Matlab实现 a对于Lagrange插值必须自编程序 b低次插值的Matlab命令分段线性插值： y=interp1(x0, y0, x)，其中输入离散数据x0、y0、x，输出对应x的插值y。三次样条插值： y=interp1(x0, y0, spline) 或 y=spline(x0, y0, x) 其中，x0、y0、x和y的意义同上。2 数值积分2.1 数值积分的基本思路2.2 三种常用数值积分方法1) 梯形公式2) 辛普森公式3) 高斯求积公式Gauss-Lobatto公式 P604）数值积分的Matlab实现trapz(x)用梯形公式计算(h=1)，输入数组x为各区间端点的函数值。trapz(x,y)用梯形公式计算，输入x,y为同长度的数组，输出y对x的积分(步长可不相等)。quad(fun,a,b,tol)用自适应辛普森公式计算，输入被积函数fun可以自定义如exp(-x.2)，也可以是fun.m命名的函数M文件，积分区间(a,b)，绝对误差tol，输出积分值。quadl(fun,a,b,tol)用自适应的Gauss-Lobatto公式计算，其余同上。常微分方程数值解常微分方程的初值问题2.初值问题的数值解法2.1 欧拉方法欧拉方法的基本思想向前欧拉公式向后欧拉公式改进的欧拉公式精度归纳：向前1阶 向后1阶 梯形2阶 改进欧拉2阶O（hp+1）p阶精度2.2 龙格-库塔方法龙格-库塔方法的基本思想龙格-库塔方法一般形式经典的龙格-库塔方法常微分方程组和高阶方程初值问题的数值方法P7374高阶方程，需要先降阶化为一阶常微分方程组2.3 龙格-库塔方法的MATLAB实现2.4 算法的收敛性、稳定性分析收敛性分析P81稳定性分析P81向后欧拉公式无条件稳定刚性现象与刚性方程精度慢稳态解的特征根决定步长快稳态解快慢稳态解衰减速度（两个特征根）相差悬殊刚性现象刚性方程求解ode23s,ode15s线性代数方程组数值解法线性代数方程组的一般形式和解法2.求解线性代数方程组的直接法2.1 高斯消元法高斯消元法列主元消去法2.2 LU分解LU分解和Cholesky分解求解三对角线性方程组的追赶法2.3 解的误差分析P95病态是解的固有性质，与解法无关。向量范数和矩阵范数P96 相容性条件3.求解线性代数方程组的迭代法3.1 雅可比迭代法3.2 高斯-赛德尔迭代法高斯赛德尔收敛快于雅可比3.3 迭代法的收敛性和收敛速度迭代公式收敛B的谱半径（B）1。谱半径不超过任一种范数（B）2)次代数方程(a0xn+a1xn-1+an=0)和超越方程(包含超越函数(如sinx, lnx)的方程) 通称为非线性方程。方程中的未知数也称为变量或变元，只含一个未知数的方程（即一元方程或单变量方程）可以记作，该方程的解也称为方程的根（或函数的零点）。n次代数方程有且只有n个根(包括复根、重根)； 5次以上的代数方程无求根公式；超越方程有无根，有几个根通常难以判断。这里仅讨论方程的实根。包含n个未知数的m个方程称为方程组，可以记作 ，其中是一个向量，是一个向量值函数。当中至少有一个非线性函数时，称为非线性方程组。多数情况下，方程组中包含的方程的个数等于未知数的个数（即m=n) 。求解非线性方程（组）的一般方法是迭代法，会出现分岔混沌现象。2 非线性方程的基本解法2.1 图形法和二分法解方程的第一步通常是确定根的近似位置或大致范围。有两种方法：图形法和二分法。图形法是利用MATLAB的函数图形功能作f(x)的图形,观察f(x)与x轴的交点，确定根的个数和范围。二分法是基于连续函数的零点存在定理，通过试探，确定函数值异号的区间后，可以用简单的二分法将区间缩小，具体步骤如下：取的中点, 若, 则即是根。否则，如, 令; 如, 令。 在内至少有一个根，且 。 再取 的中点, 如此进行下去，包含根的区间的长度每次缩小一半 (n=1, 2, )，n足够大时即可达到满意的精度。图形法和二分法都可提供迭代法的初始迭代点。2.2 迭代法迭代法的基本思想是将原方程改写成等价形式, 选择适当的初值, 按照迭代公式 计算，若迭代序列收敛到, 则满足，称为迭代函数的不动点，即为原方程f(x)=0的根。迭代法的关键在于如何构造迭代函数，使迭代序列以较快速度收敛。迭代法是否收敛取决于曲线的斜率。（P118）关于迭代法的收敛性，理论上有如下的所谓局部收敛性定理： 设在的一个邻域内连续、可微，且 则对于该邻域内的任意初值，序列xn收敛于。对迭代序列，记，若 ， 为一个正数，其中|表示某种范数（对实数可以认为就是绝对值），则称序列为阶收敛。特别地，1阶收敛称线性收敛，二阶收敛称平方收敛；若p=1, c=0，则称为超线性收敛。P越大收敛越快。利用在的泰勒展开：可得，从而可知 若，则为1阶收敛（线性收敛）；若，则为阶收敛。2.3 牛顿法将在作泰勒展开，去掉2阶及2阶以上项（即线性化）后得。设，令上面的，用代替右端的，就得到迭代公式。对应的迭代函数为，其几何意义是过点的曲线的切线与轴的交点即为(点击看图1)，称为牛顿切线法。由 知，若是的单根，即，则，这时牛顿切线法2阶收敛。当是 的重根时，牛顿切线法只是1阶收敛，并且重数越高收敛越慢。用差商代替，迭代公式变为 ，其几何意义是用割线代替了原来的切线(点击看图2)，称为割线法（或称弦截法）其收敛速度比切线法稍慢（对于单根其收敛阶数是1.618），且需要两个初值x0, x1开始迭代。3 非线性方程组的牛顿法、拟牛顿法将求解非线性方程的牛顿切线法推广到解方程组F(x)=0，其中。设是第步近似解，在作泰勒展开，线性化后用代替可得，其中为F的雅可比 (Jacobi) 矩阵 在处的值。若可逆，则可得求解方程组F(x)=0的牛顿迭代公式。实际计算中，在计算过程的第k步，通常是先计算和，再解线性方程组得到后，令即可。牛顿迭代公式是超线性收敛的（即收敛阶不小于1），稍加条件就至少是平方收敛的。当函数F比较复杂时计算雅可比矩阵很不方便，所以希望能用较简单的矩阵近似，即。仿照割线法中用差商代替的作法，使满足和计算之。这种方法称为拟牛顿法。至于如何确定，又有不同的构造方法，例如DFP、BFGS等。4 用MATLAB工具箱解非线性方程(组)4.1 fzero的基本用法fzero命令用于求单变量方程的根，所采用的算法主要是二分法、割线法和逆二次插值法等的混合方法。其最简单的调用方式为x= fzero(f,x0)函数简单可以用句柄形式：fzero(inline(x3-2*x-5),0)%初值取0或fzero(inline(x3-2*x-5),1,3)%有根区间取1,3或fzero(x)x3-2*x-5,0)fzero求的可能只是变号点而不是零点：连续函数近似零点，不连续函数，间断点；连续没变号找不到其最一般的调用方式为x,fv,ef,out = fzero(f,x0,opt,P1,P2,.)输出参数 输入参数 注意事项（控件）4.1.1 命令的输出参数其中fzero命令输出参数的含义为：x：变号点的近似值 fv：x点所对应的函数值 ef：程序停止时的状态 l 1：找到异号点l -1：没有找到异号点 Out：包含以下数据的一个结构变量 l Iterations： 迭代次数 l funcCount： 函数被调用的次数l algorithm： 实际使用的算法4.1.2 命令的输入参数其中fzero命令输出参数的含义为：1. f 函数名（必须输入的参数） 2. x0迭代初值（或有根区间）（必须输入的参数） 3. opt控制参数的结构变量，设定（或显示）控制参数的命令为Optimset(参见约束优化实验)，用户不指定或指定为时将采用缺省值。对fzero命令可选择的参数只有Display和TolX（含义见约束优化实验) 4. P1,P2,.是传给f函数的参数（如果需要的话） 4.1.3 命令注意事项1. 对简单函数f(x)可直接用MATLAB提供的inline函数输入（inline函数返回一个字符串表示的函数的句柄） 2. fzero实际上求得的不一定是函数的零点，而只是函数值发生变号的点。对于连续函数，这个点就是近似零点；但对于不连续的函数，这个点很可能只是一个间断点（且在该点两边，函数值异号）。4.2 fsolve的基本用法fsolve命令用于非线性方程组的求解（当然也可以用于方程求根，但效果一般不如fzero程序），最一般的调用方式是：x,fv,ef,out,jac = fsolve(F,x0,opt,P1,P2, . )输出参数 输入参数（控件）4.2.1 命令的输出参数其中fsolve命令输出参数的含义为：1. x： 方程组的解 2. fv：解所对应的向量函数值 3. ef：程序停止时的状态 l 1：收敛l -1：不收敛l 0：达到了迭代或函数调用的最大次数4. out：包含以下数据的一个结构变量 l Iterations： 迭代次数 l funcCount： 函数被调用的次数l algorithm： 实际使用的算法l firstorderopt：结果处梯度向量的1-范数5. jac：x点所对应的雅可比矩阵 4.2.2 命令的输入参数其中输入参数的含义为：1. f 函数名（必须输入的参数） 2. x0迭代初值（必须输入的参数） 3. opt控制参数的结构变量，设定（或显示）控制参数的命令为Optimset，有以下一些用法： Optimset /显示控制参数optimset optfun /显示程序optfun的控制参数opt=optimset /控制参数设为（即缺省值opt=optimset(optfun)/ 设定为程序optfun的控制参数缺省值Opt=optimset(par1,val1,par2,val2,.)Opt=optimset(oldopts,par1,val1,.)opt=optimset(oldopts,newopts)对fsolve命令可选择的常用参数有Diagnostics 是否显示诊断信息（ on 或off）Display 显示信息的级别（off ， iter ， final，notify）LargeScale 是否采用大规模算法（ on 或off）缺省值为offMaxIter 最大迭代次数TolFun 函数计算的误差限TolX 决策变量的误差限 Jacobian 目标函数是否采用分析Jacobi矩阵（on ，off）MaxFunEvals 目标函数最大调用次数4. P1,P2,.是传给f函数的参数（如果需要的话）4.3 roots的基本用法Roots命令用于一元多项式求根。两个相关的命令是：l r=roots(c) 输入多项式的系数c（按降幂排列），输出r为的全部根（包括复根）；l c=poly(r) 输入的全部根r，输出c为多项式的系数（按降幂排列）。无约束优化1无约束优化的基本原理、解法1.1 无约束优化的一般形式无约束非线性规划问题为其最优解通常都是局部最优解，寻找全局最优解需要对局部最优解进行比较以后得到（如果能够求出所有局部最优解的话）。1.2 最优性条件是最优解的必要条件为 ；充分条件为，且正定。1.3 下降法的基本思想在迭代的第k步，确定一个搜索方向和一个步长，使沿此方向、按此步长走一步到达下一点时，函数值下降。其基本步骤为1） 选初始解；2） 对于第次迭代解，确定搜索方向并在此方向确定搜索步长令，使0：收敛 0：函数调用次数或迭代次数达到最大值（该值在options中指定） 1 g = . % gradient of the function if nargout 2 H = . % Hessian of the function end options：包含算法控制参数的结构设定（或显示）控制参数的命令为Optimset，有以下一些用法：Optimset /显示控制参数optimset optfun /显示程序optfun的控制参数opt=optimset /控制参数设为（即缺省值opt=optimset(optfun)/ 设定为程序optfun的控制参数缺省值Opt=optimset(par1,val1,par2,val2,.)Opt=optimset(oldopts,par1,val1,.)opt=optimset(oldopts,newopts)可以设定的参数比较多，对fminunc，常用的有以下一些参数：Diagnostics 是否显示诊断信息（ on 或off）Display 显示信息的级别（off ， iter ， final，notify）LargeScale 是否采用大规模算法（ on 或off）缺省值为onMaxIter 最大迭代次数TolFun 函数计算的误差限TolX 决策变量的误差限 GradObj 目标函数是否采用分析梯度（on ，off）Hessian 目标函数是否采用分析Hess矩阵（on ，off）MaxFunEvals 目标函数最大调用次数HessUpdate 拟牛顿法修改方法（bfgs(缺省值)，dfp,gillmurray,steepdesc） LineSearchType 线搜索方法(cubicpoly,quadcubic(缺省值)1.4.3 注意事项 fminunc中输出变量、输入参数不一定写全，可以缺省。采用缺省值效果一般会很好。 当中间某个输入参数缺省时，需用占据其位置。 当函数高度非线性或严重不连续时，用程序fminsearch代替fminunc。2非线性最小二乘拟合的基本原理、解法2.1 非线性最小二乘拟合问题有一组数据 要拟合一个已知函数y=f(x, t), x=(x1,x2,xm), , x为待定系数。记误差， ，拟合误差定义为的平方和，于是问题表示为如下的优化模型：当对（的某些分量）是非线性函数时，称非线性最小二乘拟合。2.2 非线性最小二乘拟合问题的解法1）高斯牛顿法：下降方向由解出，称高斯牛顿法（GN法）。其收敛速度依赖于对的线性程度，及误差的大小，适用于对近似于线性及较小的情况。当对高度非线性或很大时，高斯牛顿法可能不收敛。2） LM方法将高斯牛顿法中修正为 其中是单位阵，0是在每次迭代中修正的参数。这个方法称为Levenbery-Marquardt法（LM法）。当很小时此法接近高斯牛顿法，当很大时接近于梯度法，所以每次迭代的下降方向在牛顿方向与负梯度方向之间。2.3 用MATLAB优化工具箱解非线性最小二乘拟合问题l 非线性最小二乘问题min s.t. v1xv2求解程序名为lsqnonlin，其最简单的调用格式为：x=lsqnonlin(F,x0, v1,v2)其最复杂的调用格式为：x,norm,res,ef,out,lam,jac = lsqnonlin(F,x0,v1,v2,opt,P1,P2, . )l 非线性拟合问题min s.t. v1xv2求解程序名为lsqcurvefit，其最简单的调用格式为：x=lsqcurvefit(F, x0,t,y,v1,v2)其最复杂的调用格式为：x,norm,res,ef,out,lam,jac = lsqcurvefit(F,x0,t,y,v1,v2,opt,P1,P2,.)输出参数 输入参数 注意事项2.3.1 程序lsqnonlin和lsqcurvefit的输出参数其中输出变量的含义为：1） x ： 最优解2） norm ： 误差的平方和，即res*res3）res: 误差向量4） ef ： 程序结束时的状态指示： 0：收敛 0：函数调用次数或迭代次数达到最大值（该值在options中指定） 1 % two output argumentsJ = . % Jacobian of the function evaluated at xend t,y: 拟合数据 v1,v2: 上下界 options：包含算法控制参数的结构设定（或显示）控制参数的命令为Optimset，有以下一些用法：Optimset /显示控制参数optimset optfun /显示程序optfun的控制参数opt=optimset /控制参数设为（即缺省值opt=optimset(optfun)/ 设定为程序optfun的控制参数缺省值Opt=optimset(par1,val1,par2,val2,.)Opt=optimset(oldopts,par1,val1,.)opt=optimset(oldopts,newopts)可以设定的参数比较多，对lsqnonlin和lsqcurvefit，常用的有以下一些参数：Diagnostics 是