Matlab 科学计算.ppt

1 Matlab与科学计算 计算机学院刘咏梅Email liuyongmei 2 第二章MATLAB数值计算基础 IntroductiontoMATLABMATLABBasicsofNumericalComputingMATLABProgramming branchingandLoops MATLABProgramming PredefinedFunctions MATLABGraphicandImageProcessingScientificComputingusingMATLAB 3 本章学习的目标 利用MATLAB进行插值和拟合利用MATLAB求积分和微分利用MATLAB求解线性方程组 4 插值和拟合 插值 利用函数f x 在某区间中若干个点的函数值 做适当特定的函数f x 使区间上其它点的值用f x 作为f x 的近似 拟合 已知某函数f x 若干离散值 通过调整f x 中若干待定系数 使f x 与已知点集的差别最小 5 插值和拟合示意图 插值示意图 拟合示意图 6 插值与拟合的应用背景 在生产实际和科学研究中 会遇到大量函数 其中相当一部分是通过测量和实验观测得到的 虽然其函数关系y f x 在某区间上客观存在 但却不知道具体的解析表达式 只能得到函数在区间 a b 上的一些离散点的函数值 导数值等 我们希望对这样的函数用一个比较简单的函数的表达式来近似地给出其整体描述 有时是虽然有明确表达式 但函数本身过于复杂而不便于进行数值计算 同样希望构造一个既能反映函数性质又便于计算的简单函数 来近似代替原来的函数 插值就是寻求近似函数的方法之一 7 计算方法中的插值法 多项式插值 设函数y f x 在区间 a b 上有定义 已知 a b 区间上n 1互异点x0 x1 xn及其函数值y0 y1 yn 若存在一个简单函数y p x 使其经过y f x 上的这n 1个已知点 x0 y0 x1 y1 xn yn 那么函数p x 称插值函数 点x0 x1 xn称为插值节点 点 x0 y0 x1 y1 xn yn 称为插值点 包含插值点的区间称为插值区间 求p x 的方法叫做插值法 f x 称为被插函数 若p x 表示为pn x a0 a1x a2x2 anxn则称pn x 为n次插值多项式 相应的插值法称为多项式插值 8 计算方法中的插值法 Lagrange插值法分段线性插值法 9 Lagrange插值法 因为pn x a0 a1x a2x2 anxn且pn xi yi i 0 1 2 n 插值问题实际上就是根据给定的n 1个节点上的函数值 确定n次多项式pn x 的n 1个系数 即由n 1个条件确定n 1个待定系数 也就是求解线性方程组的解 但是系数矩阵为Vandermonde矩阵 条件数常常很大 是病态问题 小的扰动就会产生大的偏移量 求精确解失去意义 10 Lagrange插值法 Lagrange插值公式如下 对任意k k 0 n Lagrange插值基函数 11 分段线性插值法 高次插值多项式可能产生极大的误差 这一现象称为荣格 Runge 现象 因此不宜用太多的点来做插值多项式 为了解决荣格现象 引入了分段线性插值 即通过n 1个插值点用折线段连接起来逼近原曲线 这也是计算机绘制图形的基本原理 12 荣格 Runge 现象 x 5 1 5 y 1 1 x 2 x0 5 0 1 5 y0 lagrange x y x0 y1 1 1 x0 2 plot x0 y0 r holdonplot x0 y1 b 13 一维插值函数interp1 YI interp1 X Y XI method interpolatestofindYI ThevaluesoftheunderlyingfunctionYatthepointsinthearrayXI XmustbeavectoroflengthN IfYisavector thenitmustalsohavelengthN andYIisthesamesizeasXI Method nearest linear default spline pchip cubic 14 插值实例 x 0 10 y sin x xi 0 25 10 yi interp1 x y xi plot x y o xi yi 15 插值实例 x 5 1 5 y 1 1 x 2 x0 5 0 1 5 y0 interp1 x y x0 plot x0 y0 m holdonplot x y bo 16 二维插值函数interp2 ZI interp2 X Y Z XI YI method interpolatestofindZI Thevaluesoftheunderlying2 DfunctionZatthepointsinmatricesXIandYI MatricesXandYspecifythepointsatwhichthedataZisgiven Method nearest linear default spline cubic 17 拟合 在实际中 往往需要从一组实验数据 xi yi 中找到自变量和因变量的函数关系y f x 由于观测数据往往不够准确 因此不要求y f x 经过所有的点 xi yi 而只要求在给定点xi上的误差det i f xi yi按某种要求达到最小 最小二乘拟合采用的策略是使误差的平方和达到最小 18 Matlab拟合函数 利用polyfit函数进行多项式拟合格式 p polyfit X Y N p 1 X N p 2 X N 1 p N X p N 1 实例 x 0 51 01 52 02 53 0 y 1 752 453 814 808 008 60 a polyfit x y 2 a 0 56140 82871 1560 19 x1 0 5 0 05 3 0 y1 a 1 x1 2 a 2 x1 a 3 plot x y holdon plot x1 y1 r 20 计算方法中的求积公式 微积分是理工科大学生踏入大学校门后首先要接触到的一门重要基础课程 实践证明 微积分运算也是实际工作中广泛应用的基本工具 它是科学计算最重要的组成部分 21 计算方法中的求积公式 Newton Cotes系列数值求积公式梯形公式复合梯形公式 22 Newton Cotes系列数值求积公式 一种等间隔闭式求积公式将积分区间 a b 划为n等分 步长h b a n 选取等距离节点xk a kh 构造下面的插值型求积公式 Cotes系数 23 Cotes系数 当n 1时 c0 1 2 c1 1 2当n 2时 c0 1 6 c1 2 3 c2 1 6当n 3时 c0 1 8 c1 3 8 c2 3 8 c3 1 8 梯形公式 Simpson公式 24 梯形公式 25 复合梯形求积公式 复合梯形求积公式的由来 为了提高精度 常把区间 a b 分成若干个小区间 在每个小区间上使用次数不高的Newton Cotes求积公式 如梯形公式和Simpson公式 再将所有的子区间上的积分值累加起来 复合梯形求积公式 f x 在区间 a b 上的积分值为 26 MATLAB求积函数trapz trapz trapezoidalnumericalintegration 梯形数值积分 Form Z TRAPZ Y computesanapproximationoftheintegralofYviathetrapezoidalmethod withunitspacing Tocomputetheintegralforspacingdifferentfromone multiplyZbythespacingincrement Comments Forvectors TRAPZ Y istheintegralofY Formatrices TRAPZ Y isarowvectorwiththeintegralovereachcolumn 27 trapz实例 Z TRAPZ X Y DIM orTRAPZ Y DIM integratesacrossdimensionDIMofY Example IfY 012345 thentrapz Y 1 is 1 52 53 5 andtrapz Y 2 is 28 28 trapz实例 x 0 2 pi 100 2 pi y sin x z trapz y z 0 x 0 pi 100 pi y sin x z trapz y pi 100z 1 9998 29 思考 如何利用trapz求下列积分的值 30 解线性方程组的应用背景 求解线性方程组在科学计算中最常遇到 数值方法如三次样条 有限元法等都要涉及到求解线性方程组 求解线性方程组的实际应用也是非常广泛 天气预报 石油勘探 地质测绘 电路分析 分子结构预测等 31 解线性方程组的理论方法 克莱姆法则 Cramer sRule 若线性方程组 的系数行列式D 0 则线性方程组 有唯一解 其解为xj Dj D 其中Dj是把D中第j列元素对应地换成常数项 j 1 n 而其余各列保持不变所得到的行列式 32 解线性方程组的理论方法 如果按克莱姆法则求解 假设计算机每秒可计算1亿次乘法 共n 1个行列式 每个行列式计算需要n n 1 次乘法 故共需 n2 1 n 次乘法 当问题规模n 10时的计算量耗时3 59秒 而n 100时已经非常大了 不能接受 因此说 理论上比较完善的方法在实际中却用处不大 33 利用计算机来解线性方程组的方法 解线性方程组的数值方法有两类 直接法和迭代法 直接法是指通过有限次的四则运算求解方程组的解的方法 但由于舍入误差的存在 此类方法只能得到近似解 只能适用于低阶方程 迭代法是先给出一个初始近似值 然后按一定的法则 逐步求得更加精确的近似解的方法 高阶方程组的数值解常用迭代法 上述两种方法均包含许多种算法 有的解法是针对一般性方程组 有的则是针对一些有特点的方程组 34 MATLAB求解线性方程组 AX B的解可以用X A B求解 XA B的解可以用X A B求解 如果A是m n的矩阵 当m n时可以得到唯一解 在MATLAB中 只需用左除或右除就可以解得 虽然表面上是一个简单的符号 而它内部却包含着许多自适应算法 对不同类型的方程组用不同的方法 35 高斯消去法 高斯消去法是通过一系列初等行变换来实现将增广矩阵简化为一个上三角阵 向量b在这一过程中也将被修改 再由回代过程解得向量x 类似于我们中学学过的消元法 36 高斯消去法实例 用Gauss消去法解方程组解消元得 37 再由回代过程 得解x1 x2 x3 1 38 第六章作业 第六章需要提交实验报告实验目的掌握利用Matlab进行插值和拟合的方法掌握利用Matlab求积分和微分的方法掌握利用Matlab求解线性方程组的方