06 数值微积分和常微分方程求解.ppt

第6章数值微积分与常微分方程求解;;目录;6.1 数 值 微 分;6.1.2 数值微分的实现;【例6.1】设f(x)?=?sinx，用不同的方法求函数f(x)的数值导数，并在同一个坐标系中做出f?(x)的图像。 为确定计算数值导数的点，假设在[0，pi]区间内以π/24为步长求数值导数。下面用3种方法求f(x)在这些点的导数。首先用一个5次多项式p(x)拟合函数f(x)，并对p(x)求一般意义下的导数dp(x)，求出dp(x)在假设点的值；第2种方法用diff函数直接求f(x)在假设点的数值导数；第3种方法先求出导函数f?(x)?=?cosx，然后直接求f?(x)在假设点的导数。 x=0:pi/24:pi; %用5次多项式p拟合f(x),并对拟合多项式p求导数dp在假设点的函数值 p=polyfit(x,sin(x),5); dp=polyder(p); dpx=polyval(dp,x); dx=diff(sin([x,pi+pi/24]))/(pi/24); %直接对sin(x)求数值导数 gx=cos(x); %求函数f的导函数g在假设点的导数 plot(x,dpx,x,dx,o,x,gx,+); %作图 ;对于求矩阵的差分，即为求各列或各行向量的差分，从向量的差分值可以判断列或行向量的单调性、是否等间距以及是否有重复的元素。 【例6.2】生成一个5阶魔方矩阵，按列进行差分运算。 M=magic(5) M= 17 24 1 8 15 23 5 7 14 16 4 6 13 20 22 10 12 19 21 3 11 18 25 2 9 DM=diff(M) %计算M的一阶差分 DM= 6 ?19 6 6 1 ?19 1 6 6 6 6 6 6 1 ?19 1 6 6 ?19 6 可以看出，diff函数对矩阵的每一列都进行差分运算，因而结果矩阵的列数是不变的，只有行数减1。矩阵DM第3列值相同，表明原矩阵第3列是等间距的。;四猎硼初谤暮糊蟹曰蛊港啦谢尸请募咽昌骗睁颤象弘怕钓垢晕脓指头符嫉06 数值微积分和常微分方程求解06 数值微积分和常微分方程求解;6.2.2 定积分的数值求解实现;;;;;3．梯形积分法 在MATLAB中，对由表格形式定义的函数关系的求定积分问题用梯形积分函数trapz。该函数调用格式如下。 T ?=? trapz(Y) 若Y是一向量，则从1开始取单位步长，以Y的值为函数值计算积分值。若Y是一矩阵，则计算Y的每一列的积分。例如： trapz([1:5; 2:6]) ans= 12 16 T?=?trapz(X, Y)：向量X、Y定义函数关系Y?=?f(X)。X、Y是两个等长的向量：X?=?(x1，x2，…，xn)，Y?=?(y1，y2，…，yn)，并且x1x2…xn，积分区间是[x1，xn]。;;;【例6.7】计算二重定积分 （1）建立一个函数文件fxy.m： function f=fxy(x,y) global ki; ki=ki+1; %ki用于统计被积函数的调用次数 f=exp(-x.^2/2).*sin(x.^2+y); （2）调用dblquad函数求解： global ki; ki=0; I=dblquad(@fxy,-2,2,-1,1) Ki 如果使用inline函数，则命令如下： f=inline(exp(-x.^2/2).*sin(x.^2+y),x,y); I=dblquad(f,?2,2,?1,1);6.3 常微分方程的数值求解;;6.3.2 常微分方程数值求解的实现 MATLAB提供了多个求常微分方程数值解的函数，一般调用格式为 [t,y]=solver(fname,tspan,y0[,options]) 其中t和y分别给出时间向量和相应的状态向量。 solver为求常微分方程数值解的函数ode23、ode45、ode113、ode23t、ode15s、ode23s、ode23tb、ode15i之一，表6.1所示为各函数的采用方法和适用问题。 fname是定义f(t，y)的函数文件名，该函数文件必须返回一个列向