第六章 计算方法  非线性方程求解 多项式插值与曲线拟合 数值微分与数值积分 求常微分方程数值解命令.

Presentation on theme: "第六章 计算方法  非线性方程求解 多项式插值与曲线拟合 数值微分与数值积分 求常微分方程数值解命令."— Presentation transcript:

1 第六章 计算方法  非线性方程求解 多项式插值与曲线拟合 数值微分与数值积分 求常微分方程数值解命令

2 已知函数 f(x),求 x 使 f(x)=0 —— 非线性方程
例 x = tan x x y (k = 1,2,·······) f=inline('x-tan(x)') x=fzero(f,3.5) x = 求一元函数零点 fzero(fun, x0) X = fzero(fun, x0)寻找fun在x0附近的零点

3 r x 例 x 3 – 30 x =0 p=[ ]; roots(p) 半径r =10 cm的球体,密度 =0.638.浸入水深度 x=? ans = 15/18

4 求多项式零点: roots(p); 例2 求 x3 –8x2 +6x– 30=0的解 P=[ ]; r=roots(P) r = 7.7260 i i 例1.求方程 x3 – 4x + 5 = 0 的解 P=[ ]; roots(P) ans = i i

5 P10(t) f(t) f(x) 多项式插值原理 已知 f(x) 在点 xi 上的函数值 yi=f(xi), (i=0,1,2,···,n)
求多项式 P(x)=a0 + a1x +···+ anxn 满足: P(xk)= yk (k = 0,1,…,n) P10(t) f(t) f(x)

6 zi=interp2(x, y, z, xi, yi, ‘method’)
一维插值: yi = interp1(x, y, xi, ‘method ’) x = 0:10; y = sin(x); xi = 0:.25:10; yi = interp1(x,y,xi); plot(x,y,'o',xi,yi) method nearest 最近点插值 linear 线性插值 spline 样条插值 cubic 立方插值 二维插值 zi=interp2(x, y, z, xi, yi, ‘method’) 三维插值 vi = interp3(x,y,z,v,xi,yi,zi, ‘method’) 17/18

7 例1.正弦曲线数据插值试验

8 离散数据的多项式拟合 x x1 x2 ·········· xm f(x) y1 y2 ·········· ym 使得 已知数据表
求拟合函数:  (x) = a0 + a1x + ······+ anxn 使得 达到最小

9 P(x)=a0xn + a1xn-1 + ······ + an-1x + an
多项式拟合命令 P=polyfit(X,Y,n) 求出(最小二乘意义下)n次拟合多项式 P(x)=a0xn + a1xn-1 + ······ + an-1x + an 计算结果为系数: P=[ a0, a1, ······ , an-1, an ] 多项式求值命令: y1=polyval(P, x) 其中,P是n次多项式的系数,x是自变量的值,y1是多项式在x处的值 17/18

10 定积分与积分和式 Sn 5.2908 5.1044 4.9835 ······ 4.8999 右矩形和 h 1 0.5 0.2 ······
4/16

11 梯形数值积分公式 左矩形 梯形 右矩形 5/16

12 梯形数值积分命令 trapz() clear x=0:pi/100:pi; y=sin(x); trapz(x,y) ans = 1.9998 clear x=sort(rand(1,101)*pi); y=sin(x); trapz(x,y) ans = 1.9981 rand(1,101)产生101个均匀随机数,每个数都介于0-1之间

13 辛卜生求数值积分命令quad() clear fun=inline('1./(x.^3-2*x-5)') ezplot(fun,[0,2]) [q,n]=quad(fun,0,2) q = n = 53

14 二重积分dblquad()与三重积分 fun=inline('y*sin(x)+x*cos(y)') Q=dblquad(fun,pi,2*pi,0,pi) Q = [x,y]=meshgrid(pi:.1:2*pi,0:.1:pi); z=fun(x,y); mesh(x,y,z)

15 物理量在空间的分布情况其及随时间变化的规律。
研究具有实际应用背景的数学函数 物理量在空间的分布情况其及随时间变化的规律。 物理现象  物理定律  微分方程  数值解 建模型: 物理学家、化学家、生物学家、社会学家 2/18

16  y(t)=y0 e r t 马尔萨斯(Malthus)人口模型 某时刻人口数为y0, 在以后某一时间 t 人口数y 的变化率正比于y
Thomas Malthus ( ) MATLAB符号求解方法 y=dsolve('Dy=r*y','t') y =C1*exp(r*t) 弹性力学问题 单摆问题 天体力学二体问题 电磁场计算············ 令 r = a, y0 = exp( b )  y(t) = exp( at + b) a = ?, b = ?  人口数据拟合实验 3/18

17 MATLAB求常微分方程数值解 简例 f=inline('y-x.*y.^2'); [x,y]=ode23(f,[0,2],1)
f=inline('y-x.*y.^2'); [x,y]=ode23(f,[0,2],1) 4/18

18 ode23(),ode45() 求常微分方程初值问题数值解函数 解一阶常微分方程组初值问题命令格式
[T,y] = ode23(' F ',Tspan,y0) 其中, F是ODE文件, 表示 y' = F(t,y)中右端函数 Tspan = [t0 Tfinal] —— 求解区域; y0 —— 初始条件 注: 函数F(t,y) 必须返回列向量. 数值解 y 的每一行对应于列向量T中的每一时刻 5/18

19 捕食者和被捕食者模型 ( 0≤ t ≤ 15 ) y1(0)=20, y2(0)=20 数会急增，数学家将鱼分成捕食者和被捕食者，分别
生物学家研究发现，战争期间捕鱼量减少,食肉鱼总 数会急增，数学家将鱼分成捕食者和被捕食者，分别 以y2(t) , y1(t)表示它们在时刻t 的总数。假定被捕食者 (食用鱼)有充分的食物。 ( 0≤ t ≤ 15 ) y1(0)=20, y2(0)=20 6/18

20 function yp = lotka(t,y) % predator-prey model.
MATLAB求解: (1)定义函数; (2)调用ode23 函数文件,文件名:lotka.m function yp = lotka(t,y) % predator-prey model. yp=diag([ *y(2), *y(1)])*y; t0=0;tf=15; y0=[20 20]; [t,y]=ode23(‘lotka’,[t0,tf],y0); plot(t,y(:,1),t,y(:,2),'r') —— 捕食者y2 —— 被捕食者y1 7/18

21 洛伦兹模型与“蝴蝶效应” x(0)=0，y(0)=0，z(0)=0.01 function z=flo(t,y)
Y-X平面 x(0)=0，y(0)=0，z(0)=0.01  = 8/3  =  =28 function z=flo(t,y) A=[-8/3 0 y(2); ; -y(2) 28 -1]; z=A*y; y0=[0;0;0.01]; [x,y]=ode23(‘flo’,[0, 80],y0); plot(y(:,2),y(:,1)) 8/18