《暑期培训Matlab之工具箱1教材》由会员分享，可在线阅读，更多相关《暑期培训Matlab之工具箱1教材（64页珍藏版）》请在金锄头文库上搜索。
1、 微 分 方 程 数学建模 微分方程的解析解 求微分方程（组）解析解的命令: dsolve(方程1,方程2,方程n,初始条件,自变量) 结 果：u = tg(t- c) 解 输入命令: y=dsolve(D2y+4*Dy+29*y=0,y(0)=0,Dy(0)=15,x ) 结 果 为 : y =3e-2xsin（5x） 解 输入命令 ： x,y,z=dsolve(Dx=2*x-3*y+3*z, Dy=4*x-5*y+3*z,Dz=4*x-4*y+2*z, t)； x=simple(x) % 将x化简 y=simple(y) z=simple(z) 结 果 为：x = (c1-c2+c3+c2
2、e -3t-c3e-3t)e2t y = -c1e-4t+c2e-4t+c2e-3t-c3e-3t+c1-c2+c3)e2t z = (-c1e-4t+c2e-4t+c1-c2+c3)e2t 微分方程的数值解 常微分方程数值解的定义 在生产和科研中所处理的微分方程往往很复杂且大多 得不出一般解。而在实际上对初值问题，一般是要求得 到解在若干个点上满足规定精确度的近似值，或者得到 一个满足精确度要求的便于计算的表达式。 因此，研究常微分方程的数值解法是十分必要的 。 返 回 用MATLAB软件求常微分方程的数值解 t，x=solver(f,ts,x0,options) ode45 ode23 o
3、de11 3ode1 5sode 23s 由待解 方程写 成的M 文件名 ts=t0,t f，t0、 tf为自变 量的初值 和终值 函数 的初 值 ode23：组合的2/3阶龙格库塔费尔贝格算法 ode45：运用组合的4/5阶龙格库塔费尔贝格算法 自变 量值 函数 值 用于设定误差限(缺省时设定相对误差10-3, 绝对误差10-6), 命令为：options=odeset（reltol,rt,abstol,at）, rt，at：分别为设定的相对误差和绝对误差 1在解含n个未知数的方程组时，x0和x均为n维向量， M文件中的待解方程组应以x的分量形式写出 2使用MATLAB软件求数值解时，高阶微
4、分方程必 须等价地变换成一阶微分方程组 注意: 解: 令 y1=x，y2=y1 1建立M文件vdp1000m如下： function dy=vdp1000(t,y) dy=zeros(2,1); dy(1)=y(2); dy(2)=1000*(1-y(1)2)*y(2)-y(1); 2取t0=0，tf=3000，输入命令： T,Y=ode15s(vdp1000,0 3000,2 0); plot(T,Y(:,1),-) 3结果如图To MATLAB（ff4） 解 1建立M文件rigidm如下： function dy=rigid(t,y) dy=zeros(3,1); dy(1)=y(2)*y
5、(3); dy(2)=-y(1)*y(3); dy(3)=-051*y(1)*y(2); 2取t0=0，tf=12，输入命令： T,Y=ode45(rigid,0 12,0 1 1); plot(T,Y(:,1),-,T,Y(:,2),*,T,Y(:,3),+) 3结果如图 To MATLAB（ff5） 图中，y1的图形为实线，y2的图形为“*”线，y3的图形为“+”线 返 回 插 值 数学建模 一维插值的定义 已知 n+1个节点 其中 互不相同，不妨设 求任一插值点处的插值 节点可视为由 产生, 表达式复杂, 或无封闭形式, 或未知. 构造一个(相对简单的)函数通过全部节点, 即 再用计算插
6、值，即 用MATLAB作插值计算 一维插值函数： yi=interp1(x，y，xi，method) 插值方法被插值点 插值节点 xi处的插 值结果 nearest 最邻近插值； linear 线性插值； spline 三次样条插值； cubic 立方插值； 缺省时 分段线性插值 注意：所有的插值方法 都要求x是单调的，并且xi不 能够超过x的范围 例：从1点12点的11小时内，每隔1小时测量一次温度， 测得的温度的数值依次为：5，8，9，15，25，29，31，30 ，22，25，27，24试估计每隔1/10小时的温度值 hours=1:12; temps=5 8 9 15 25 29 31
7、 30 22 25 27 24; h=1:0.1:12; t=interp1(hours,temps,h,spline); (直接输出数据将是很多的) plot(hours,temps,+,h,t,hours,temps,r:) %作图 xlabel(Hour),ylabel(Degrees Celsius) 二维插值的定义 x y O 第一种（网格节点）： 已知 mn个节点 其中互不相同，不妨设 构造一个二元函数通过全部已知节点,即 再用计算插值，即 第二种（散乱节点）： y x O 已知n个节点 其中互不相同， 构造一个二元函数通过全部已知节点,即 再用计算插值，即 要求x0,y0单调；x
8、，y可取为矩阵，或x 取行向量，y取为列向量，x,y的值分别不能超 出x0,y0的范围 z=interp2(x0,y0,z0,x,y,method) 被插值点 插值方法 用MATLAB作网格节点数据的插值 插值 节点 被插值点 的函数值 nearest 最邻近插值； linear 双线性插值； cubic 双三次插值； 缺省时 双线性插值. 例：测得平板表面35网格点处的温度分别为： 82 81 80 82 84 79 63 61 65 81 84 84 82 85 86 试作出平板表面的温度分布曲面z=f(x,y)的图形 输入以下命令： x=1:5; y=1:3; temps=82 81 8
9、0 82 84;79 63 61 65 81;84 84 82 85 86; mesh(x,y,temps) 1.先在三维坐标画出原始数据，画出粗糙的温度分布曲线图. 再输入以下命令: xi=1:0.2:5; yi=1:0.2:3; zi=interp2(x,y,temps,xi,yi,cubic) ; mesh(xi,yi,zi) 画出插值后的温度分布曲面图. 2以平滑数据,在 x、y方向上每隔0.2个单位的地方进行插值. 插值函数griddata格式为: cz =griddata（x，y，z，cx，cy，method） 用MATLAB作散点数据的插值计算 要求cx取行向量，cy取为列向量
10、被插值点 插值方法 插值 节点 被插值点 的函数值 nearest最邻近插值 linear 双线性插值 cubic 双三次插值 v4- MATLAB提供的插值方法 缺省时, 双线性插值 例 在某海域测得一些点(x,y)处的水深z由下表 给出，船的吃水深度为5英尺，在矩形区域（75， 200）（-50，150）里的哪些地方船要避免进入 4.作出水深小于5的海域范围,即z=5的等高线. 3 作海底曲面图 数学建模 拟 合 拟 合 问 题 引 例 1 温度t(C) 20.5 32.7 51.0 73.0 95.7 电阻R() 765 826 873 942 1032 已知热敏电阻数据 ： 求60C时
11、的电阻R 设 R=at+b a,b为待定系数 拟 合 问 题 引 例 2 t (h) 0.25 0.5 1 1.5 2 3 4 6 8 c (g/ml) 19.21 18.15 15.36 14.10 12.89 9.32 7.45 5.24 3.01 已知一室模型快速静脉注射下的血药浓度数据(t=0注射300mg) 求血药浓度随时间的变化规律c(t). 作半对数坐标系(semilogy)下的图形 曲 线 拟 合 问 题 的 提 法 已知一组（二维）数据，即平面上 n个点（xi,yi) i=1,n, 寻求一个函数（曲线）y=f(x), 使 f(x) 在某种准则下与所有 数据点最为接近，即曲线拟
12、合得最好 + + + + + + + + + x y y=f(x) (xi,yi) i i 为点（xi,yi) 与曲线 y=f(x) 的距离 拟合与插值的关系 函数插值与曲线拟合都是要根据一组数据构造一个函数作 为近似，由于近似的要求不同，二者在数学方法上是完全不同 的 实例：下面数据是某次实验所得，希望得到X和 f之间的关系？ MATLAB(cn) 问题：给定一批数据点，需确定满足特定要求的曲线或曲面 解决方案： 若不要求曲线（面）通过所有数据点，而是要求它反映对象 整体的变化趋势，这就是数据拟合，又称曲线拟合或曲面拟合 若要求所求曲线（面）通过所给所有数据点，就是插值问题； 最临近插值、线
13、性插值、样条插值与曲线拟合结果： 线性最小二乘拟合 f(x)=a1r1(x)+ +amrm(x)中 函数r1(x),rm(x)的选取 1. 通过机理分析建立数学模型来确定 f(x)； + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + f=a1+a2x f=a1+a2x+a3x2 f=a1+a2x+a3x2 f=a1+a2/xf=aebx f=ae-bx 2. 将数据 (xi,yi) i=1, ,n 作图，通过直观判断确定 f(x)： 用MATLAB解拟合问题 1.线性最小二乘拟合 2.非线性最小二乘拟合 用MATLAB作线性
14、最小二乘拟合 1. 作多项式f(x)=a1xm+ +amx+am+1拟合,可利用已有程序: a=polyfit(x,y,m) 2. 对超定方程组 可得最小二乘意义下的解 ，用 3.多项式在x处的值y可用以下命令计算： y=polyval（a，x） 输出拟合多项式系数 a=a1, ,am , am+1 (数组)） 输入同长度 的数组x，y 拟合多项 式次数 即要求 出二次多项式: 中 的使得: 例 对下面一组数据作二次多项式拟合 1）输入以下命令： x=0:0.1:1; y=-0.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2; R=(x.2) x ones(11,1)； A=Ry MATLAB(zxec1) 解法1用解超定方程的方法 2）计算结果： = -9.8108 20.1293 -0.0317 1）输入以下命令： x=0:0.1:1; y=-0.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2; A=polyfit(x,y,2) z=polyval(A,x); plot(x,y,k+,x,z,r) %作出数据点和拟合曲线的图形 2）计算结果： = -9.8108 20.1293 -0.
