数学模型讲稿

1、数学建模与数学实验数学建模与数学实验 MATLAB介绍介绍 MATLABMATLAB介绍介绍 一、变 量 二、函 数 三、 矩阵（数值）运算 四、 MATLABMATLAB编程 六、 符号运算 五、 优化求解 七、 多项式运算 数学运算符号及标点符号数学运算符号及标点符号 +加法运算，适用于两个数或两个同阶矩阵相加. 减法运算 *乘法运算 .*点乘运算 /除法运算 ./点除运算 乘幂运算 .点乘幂运算 反斜杠表示左除. （1）MATLAB的每条命令后，若为逗号逗号或或无标点无标点符号， 则显示命令的结果；若命令后为分号分号，则禁止显示结果. （2）“%” 后面所有文字为注释. （3） “.”表

2、示续行. 函 数名 称函 数名 称 sin(x)正弦函数asin(x)反正弦函数 cos(x)余弦函数acos(x)反余弦函数 tan(x)正切函数atan(x)反正切函数 abs(x) 绝对值max(x) 最大值 min(x) 最小值sum(x)元素的总和 sqrt(x) 开平方exp(x)以 e为底的指数 log(x)自然对数 )(log10 x 以 10 为底的对数 sign(x)符号函数fix(x)取整 数学函数数学函数 MATLAB的内部函数是有限的，有时为了研究某一个函数的各种性态，需要为的内部函数是有限的，有时为了研究某一个函数的各种性态，需要为 MATLAB定义新函数，为此必须

3、编写函数文件定义新函数，为此必须编写函数文件. 函数文件是文件名后缀为函数文件是文件名后缀为M的文件的文件 ，这类文件的，这类文件的第一行必须是一特殊字符第一行必须是一特殊字符function开始开始，格式为：，格式为： function 因变量名因变量名=函数名（自变量名）函数名（自变量名） 函数值的获得必须通过具体的运算实现，并赋给因变量函数值的获得必须通过具体的运算实现，并赋给因变量. M文件文件 M文件建立方法：文件建立方法：1. 在在Matlab中，点中，点:File-New-M-file 2. 在编辑窗口中输入程序内容在编辑窗口中输入程序内容 3. 点：点：File-Save，存盘

4、，存盘，M文件名必须文件名必须 与函数名一致。与函数名一致。 Matlab的应用程序也以的应用程序也以M文件保存。文件保存。 例：定义函数例：定义函数 f(x1,x2)=100(x2-x12)2+(1-x1)2 function f=fun(x) f=100*(x(2)-x(1)2)2+(1-x(1)2 1.建立建立M文件：文件：fun.m MATLAB(fu n) 2. 可以直接使用函数可以直接使用函数fun.m 例如：计算例如：计算 f(1,2), 只需在只需在Matlab命令窗口键入命令：命令窗口键入命令： x=1 2 fun(x) 数据的输入输出 1. input函数：用于向计算机输入

5、一个参数。 调用格式： A=input(提示信息，选项)； 注：s选项，则允许用户输入一个字符串。 例如想输入一个人的姓名，可采用命令 xm=input(Whats your name:,s) 【例2】 求一元二次方程a2 +bx+c=0的根。 二数据的输入输出 a=input(a=?); b=input(b=?); c=input(c=?); d=b*b-4*a*c; x=(-b+sqrt(d)/(2*a),(-b-sqrt(d)/(2*a) 将该程序以aa.m文件存盘，然后运行aa.m文件。 2 pause函数：暂停程序的执行。 调用格式： pause(延迟秒数) 注：如果省略延迟时间，直

6、接使用pause，则将暂 停程序，直到用户按任一键后程序继续执行。 二数据的输入输出 3 disp函数：函数：命令窗口输出函数。 调用格式： disp(输出项) 注：输出项为字符串或矩阵。 例如 A=Hello,MATLAB; disp(A) 输出为： Hello,MATLAB 三、三、 矩矩 阵阵 逗号逗号或或空格空格用于分隔某一行的元素，用于分隔某一行的元素，分号分号用于区分不用于区分不 同的行同的行. 除了分号，在输入矩阵时，按除了分号，在输入矩阵时，按Enter键也表示开键也表示开 始一新行始一新行. 输入矩阵时，严格要求所有行有相同的列输入矩阵时，严格要求所有行有相同的列. 例例 m

7、=1 2 3 4 ；5 6 7 8；9 10 11 12 p=1 1 1 1 2 2 2 2 3 3 3 3 1、矩阵的建立、矩阵的建立 特殊矩阵特殊矩阵的建立：. MATLAB(m atrix1) d=eye(m，n) 产生一个产生一个m行、行、n列的单位矩阵列的单位矩阵 c=ones(m，n) 产生一个产生一个m行、行、n列的元素列的元素 全为全为1的矩阵的矩阵 b=zeros(m，n) 产生一个产生一个m行、行、n列的零矩阵列的零矩阵 a= 产生一个空矩阵，当对一项操作无结产生一个空矩阵，当对一项操作无结 果时，返回空矩阵，空矩阵的大小为零果时，返回空矩阵，空矩阵的大小为零. 关系与逻辑

8、运算关系与逻辑运算 1、关系操作符、关系操作符 关关系系操操作作符符说说明明 小小于于 大大于于 =大大于于或或等等于于 = =等等于于 =不不等等于于 1、for循环：循环：允许一组命令以固定的和预定的次数重复允许一组命令以固定的和预定的次数重复 for x=array commands end 在在for和和end语句之间的命令串语句之间的命令串commands按数组（按数组（array）中的每一列执行一次）中的每一列执行一次. 在每一次迭代中，在每一次迭代中，x被指定为数组的下一列，即在第被指定为数组的下一列，即在第n次循环中，次循环中，x=array(：，：，n) 控制流控制流 MAT

9、LAB提供三种决策或控制流结构：提供三种决策或控制流结构： for循环、循环、while循环、循环、if-else-end结构结构. 这些结构经常包含大量的这些结构经常包含大量的MATLAB命令，故经常出现在命令，故经常出现在MATLAB程序中，而程序中，而 不是直接加在不是直接加在MATLAB提示符下提示符下. 例例 对对n=1,2,10,求求xn= 的值的值 10 sin n MATLAB(for 1) while expression commands end 只要在表达式只要在表达式(expression)里的所有元素为真，就执行里的所有元素为真，就执行while和和end语句之间的命

10、语句之间的命 令串令串commands. 2、While循环循环 与与for循环以固定次数求一组命令相反，循环以固定次数求一组命令相反，while循环以不定的次数求一组语句的值循环以不定的次数求一组语句的值. MATLAB(w hile1) 例例 设银行年利率为设银行年利率为11.25%。将。将10000元钱存入银行，问多长时间会连本元钱存入银行，问多长时间会连本 带利翻一番？带利翻一番？ 3、If-Else-End结构结构 （1）有一个选择的一般形式是：）有一个选择的一般形式是： if expression commands end 如果在表达式如果在表达式(expression)里的所有元

11、素为真里的所有元素为真 ，就执行，就执行if和和end语句之间的命令串语句之间的命令串commands. MATLAB(fu n1) ) 1(),2( , 12 11 )( 2 ff xx xx xf求设例 先建立先建立M文件文件fun1.m定义函数定义函数f（x），再在），再在Matlab命令窗口命令窗口 输入输入fun1(2),fun1(-1)即可。即可。 2) 有三个或更多的选择的一般形式是：有三个或更多的选择的一般形式是： if （expression1） commands1 else if （expression2） commands2 else if （expression3） c

12、ommands3 else if else end 或switch I case otherwise end ) 1(),5 . 0(),2( , 0 x 1x02 1x1x f(x) 3 2 fff x x求设例 先建立先建立M文件文件fun2.m定义函数定义函数f（x），再在），再在Matlab命令窗口命令窗口 输入输入fun2(2),fun2(0.5), fun2(-1)即可。即可。 MATLAB(fun2) 046 62 975 1352 4321 432 421 4321 xxxx xxx xxx xxxx 程序：A=2,1,-5,1; 1,-5,0,7; 0,2,1,-1; 1,6

13、,-1,-4; b=13,-9,6,0; x=Ab; ans = -66.5556 25.6667 -18.7778 26.5556 Eg 101 Eg 101 求解线性方程组求解线性方程组 Eg 102 Eg 102 求解非线性方程求解非线性方程 附近的根在2052)( 0 3 xxxxf Z = 2.0946 （1）建立函数文件eg102.m function y = eg102(x) y = x.3-2*x-5; （2）用fzero函数求根 Z=fzero(eg102,2) Eg 103 Eg 103 求解非线性方程组求解非线性方程组 附近的数值解在)5 . 0 , 5 . 0( 0si

14、n3 . 0cos6 . 0 0cos3 . 0sin6 . 0 yxy yxx （1）建立函数文件eg103.m function q=eg103(p) x=p(1); y=p(2); q(1)=x-0.6*sin(x)-0.3*cos(y); q(2)=y-0.6*cos(x)+0.3*sin(y); （2）调用fsolve函数求方程的根 X=fsolve(eg103,0.5,0.5); Q=eg103(X) X = 0.6354 0.3734 Q = 1.0e-009 * 0.2375 0.2957 Eg 201 Eg 201 微分方程解析解微分方程解析解 用dsolve 求解： u=d

15、solve(Du=1+u2,t) 结果为： u = tan(t+C1) 的通解求解 1 2 u dt du （1） （2） 的特解求解 15)0(, 0)0( 0294 2 2 yy y dx dy dx yd 用dsolve 求解： Y= dsolve(D2y+4*Dy+29*y=0,y(0)=0,Dy(0)=15,x) 结果为： Y = 3*exp(-2*x)*sin(5*x) Eg 202 Eg 202 微分方程数值解微分方程数值解 上的数值解在求解12, 0, 1)0(, 1)0(, 0)0( 51. 0 321 321 213 312 321 yyy yyy yyy yyy yyy

16、（1）建立函数文件eg202.m function dy = eg202(t,y) dy = zeros(3,1); % a column vector dy(1) = y(2) * y(3); dy(2) = -y(1) * y(3); dy(3) = -0.51 * y(1) * y(2); （2）调用ode45函数求方程的根,（to=0,tf=12） T,Y = ode45(eg202,0 12,0 1 1); plot(T,Y(:,1),-,T,Y(:,2),-.,T,Y(:,3),.); 024681012 -1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0

17、.8 1 Eg 301 Eg 301 二维图形的画法（一）二维图形的画法（一） （1）在 上画出sin(x)与cos(x)的图形 2 , 0 用plot命令： x=linspace(0,2*pi,30); y=sin(x); z=cos(x); plot(x,y,r,x,z,b*); 01234567 -1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 （2）符号函数（显函数、隐函数、参数方程）画图 画出显函数 在 -3,3 上的图形 3 xy 用ezplot命令： ezplot(x3,-3,3); -3-2-10123 -25 -20 -15 -10 -5

18、0 5 10 15 20 25 x x3 画出隐函数 在 -2,0.5，0,2 上的图形 0)sin(xye x 用ezplot命令： ezplot(exp(x)+sin(x*y),-2,0.5,0,2); -2-1.5-1-0.500.5 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 x y exp(x)+sin(x y) = 0 Eg 301 Eg 301 二维图形的画法（二）二维图形的画法（二） 画出参数方程 在 上的图形 tytx 33 sin,cos 用ezplot命令： ezplot(cos(t).3,sin(t).3,-pi,pi); , -1-0.

19、500.51 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 x y x = cos(t).3, y = sin(t).3 做法一：用fplot命令： fplot(exp(2*x)+sin(3*x.2),-1,2); 做法二：先建立函数文件eg301.m function y=eg301(x) y=exp(2*x)+sin(3*x.2); 在用fplot命令，或ezplot命令： fplot(eg301,-1,2); 或ezplot(eg301,-1,2); 画出函数 在 -1,2 上的图形)3sin( 22 xey x -1-0.500.511.52 0 10 2

20、0 30 40 50 60 Eg 302 Eg 302 三维图形的画法（一）三维图形的画法（一） 画出参数方程 x=sin(t),y=cos(t),z=t 在 上的图形 t=0:pi/50:10*pi; plot3(sin(t),cos(t),t,b); 10, 0 （1）空间三维曲线图形 plot3 -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 0 10 20 30 40 （2）空间三维曲面图形surf,mesh,meshz 用surf画出 在 -3,3*1,5 上的图形 2 )(yxz x=-3:0.1:3; y=1:0.1:5; xx,yy=meshgrid(x,y);

21、 zz=(xx+yy).2; surf(xx,yy,zz); %将数据点所表示的曲面画出来。-4 -2 0 2 4 1 2 3 4 5 0 20 40 60 80 Eg 302 Eg 302 三维图形的画法（二）三维图形的画法（二） 用mesh画出 在 -3,3*1,5 上的图形 2 )(yxz x=-3:0.1:3; y=1:0.1:5; xx,yy=meshgrid(x,y); zz=(xx+yy).2; mesh(xx,yy,zz); %将数据点用网格面表示出来。 -4 -2 0 2 4 1 2 3 4 5 0 20 40 60 80 用meshz画出 山峰的图形 xx,yy=meshg

22、rid(-3:0.125:3); zz=peaks(xx,yy); meshz(xx,yy,zz); %在网格周围画一个参考平面 -4 -2 0 2 4 -4 -2 0 2 4 -10 -5 0 5 10 Eg 401 Eg 401 用用matlabmatlab优化工具箱解线性规划（一）优化工具箱解线性规划（一） 求解：使问题 321 436xxxz在约束条件下的最小值 321 321 20 ,500 ,30 120 xxx xxx 用linprog求解： c=6,3,4; A=0,1,0; b=50; Aeq=1,1,1; beq=120; vlb=30,0,20; vub=; v,fval

23、=linprog(c,A,b,Aeq,beq,vlb,vub) 结果为：Optimization terminated successfully. v = 30.0000 50.0000 40.0000 fval = 490.0000 vubXvlb beqXAeq bAX tscXz . . min vubXvlb beqXAeq bAX tscXz . . min c = -5 -4 -6 A = 1 -1 1; 3 2 4; 3 2 0; b = 20; 42; 30; lb = zeros(3,1); x,fval = linprog(c,A,b,lb) 结果为：Optimizatio

24、n terminated successfully. v = 0 15.0000 3.0000 fval = -78.0000 321 21 321 321 0 ,0 ,0 3023 42423 20 xxx xx xxx xxx 例2使问题 321 645)(xxxxf 在约束条件下的最小值 Eg 401 Eg 401 用用matlabmatlab优化工具箱解线性规划优化工具箱解线性规划 （二）（二） 例：某厂每日8h的产量不低于1800件。为了进行质量控制，计 划聘请两种不同水平的检验员。一级检验员标准为：速度25件/h, 正确率98，计时工资4元/h，二级检验员的标准为：速度15件 /h

25、，正确率95，计时工资3元/h，检验员每检错一次，工厂要 损失2元。为使总检验费用最省，该工厂应聘请一级、二级检验 员各几名？ 解：设需要一级x1人，二级x2人， 则应付工资为84x183x2 因检验员出错造成的损失为： （8252x18155x2)2= 故：目标函数为 约束条件为 21 2432xx 21 128xx 212121 3640)128()2432( minxxxxxxz 0, 0 1800*15*8 1800*25*8 1800*15*8*25*8 21 2 1 21 xx x x xx ;150 ; 90 4535 . . 3640 min 21 21 21 xx xx ts

26、 xxz即求解： 用linprog求解： c=40,36; A=-5,-3; b=-45; Aeq=; beq=; vlb=zeros(2,1); vub=9;15; v,fval=linprog(c,A,b,Aeq,beq,vlb,vub) 结果为：Optimization terminated successfully. v = 9.0000 0.0000 fval = 360 即：只需聘用9个一级检验员。 vubXvlb beqXAeq bAX tscXz . . min Eg 402 Eg 402 用用matlabmatlab优化工具箱解非线性优化工具箱解非线性 规划规划 vubXvl

27、b XCeq XG beqXAeq bAX tsXFz 0)( 0)( . . )( min求解： ; 010 05 . 1 0 . . ) 12424()( min 21 2121 21 221 2 2 2 1 1 xx xxxx xx ts xxxxxexf x （1）先建立文件eg402.m定义目标函数： function f=eg402(x) f=exp(x(1)*(4*x(1)2+2*x(2)2+4*x(1)*x(2)+2*x(2)+1; （2）再建立文件eg4021.m定义非线性约束条件： function G,Ceq=eg4021(x) G=1.5+x(1)*x(2)-x(1)-

28、x(2);-x(1)*x(2)-10; Ceq=x(1)+x(2); （3）再用fmincon求解： x0=-1;1; A=; b=; Aeq=1,1; beq=0; vlb=; vub=; x,fval=fmincon(eg402,x0,A,b,Aeq,beq,vlb,vub,eg4021) 结果为：Optimization terminated successfully. x = -1.2250 1.2250 fval = 1.8951 vubXvlb XCeq XG beqXAeq bAX tsXFz 0)( 0)( . . )( min Eg 403 Eg 403 用用matlabma

29、tlab优化工具箱解二次规优化工具箱解二次规 划划 VUBXVLB beqXAeq bAX tsXHXXz T . . c 2 1 min T (非线性中特殊的一种:目标函数为二次函数) 例：求解问题 0, 0 22 2 . . 2262),( min 21 21 21 2 221 2 12121 xx xx xx ts xxxxxxxxf 用quadprog求解： H=1,-1;-1,2; c=-2;-6; A=1,1;-1,2; b=2;2; Aeq=; beq=; VLB=0;0; VUB=; x,z=quadprog(H,c,A,b,Aeq,beq,VLB,VUB) 结果为：Optim

30、ization terminated successfully. x = 0.6667 1.3333 z = -8.2222 用样本的信息来推断总体的特征叫统计推断. . 统计分析的主要目的是由样本推断总体。故统计学的主体是 统计推断 统计分析-统计推断 参数 点估计 估计 区间估计 统 u u检验 计 参数 t t检验 推 假设 F F检验 断 检验 秩和检验 非参 中位数检验 卡方检验 一元线性回归分析的主要任务是： 1、用试验值（样本值）对 0 、 1 和 作点估计； 2、对回归系数 0 、 1 作假设检验； 3、在 x= 0 x处对 y作预测，对 y作区间估计. 一般地，称由xy 10

31、 确定的模型为一元线性回归模型一元线性回归模型， 记为 2 10 , 0 DE xy 固定的未知参数 0 、 1 称为回归系数，自变量 x 也称为回归变量. 1、 2、求回归系数的点估计和区间估计、并检验回归模型：求回归系数的点估计和区间估计、并检验回归模型： b, bint,r,rint,stats=regress(Y,X,alpha) 回归系数的区间估计 残差 用于检验回归模型的统计量， 有三个数值：相关系数r2、 F值、与F对应的概率p 置信区间 显著性水平 （缺省时为0.05） 相关系数 r2越接近 1，说明回归方程越显著； F F1-（k，n-k-1）时拒绝 H0，F 越大，说明回归

32、方程越显著； 与 F 对应的概率 p时拒绝 H0，回归模型成立. Eg 501 Eg 501 用用matlabmatlab统计工具箱作回归分析统计工具箱作回归分析 （一元线性回归）（一元线性回归） 例 测16名成年女子的身高与腿长所得数据如下： 身高143145146147149150153154155156157158159160162164 腿长8885889192939395969897969899100102 以身高x为横坐标，以腿长y为纵坐标将这些数据点（xi，yi）在平面 直角坐标系上标出. 140145150155160165 84 86 88 90 92 94 96 98 10

33、0 102 散点图 xy 10 1 1、输入数据 x=143,145,146,147,149,150,153,154,155,156,157,158,159, 160,162,164; X=ones(16,1),x; Y=88 85 88 91 92 93 93 95 96 98 97 96 98 99 100 102; 2 2、回归分析及检验： b,bint,r,rint,stats=regress(Y,X) 140145150155160165 84 86 88 90 92 94 96 98 100 102 246810121416 -5 -4 -3 -2 -1 0 1 2 3 4 Res

34、idual Case Order Plot Residuals Case Number 参数参数参数估计值参数估计值置信区间置信区间 0-16.0730-33.7071 1.5612 10.71940.6047 0.8340 R2= 0.9282 F= 180.9531 p= 0.0000 tt xy071940730.16 R R2 20.92890.9289，拟合度较高，回归模型建立成功 Eg 502 Eg 502 用用matlabmatlab作多元线性回归分析作多元线性回归分析 问问 题题 建立投资额模型，研究某地区实际投资额与国民生产总值 ( GNP ) ( GNP ) 及物 价指数

35、( PI ) ( PI ) 的关系 2.06883073.0424.5201.00001185.9195.010 1.95142954.7474.9190.96011077.6166.49 1.78422631.7401.9180.9145 992.7144.28 1.63422417.8423.0170.8679 944.0149.37 1.50422163.9386.6160.8254 873.4133.36 1.40051918.3324.1150.7906 799.0122.85 1.32341718.0257.9140.7676 756.0125.74 1.25791549.2206

36、.1130.7436 691.1113.53 1.15081434.2228.7120.7277 637.797.42 1.05751326.4 229.8110.7167 596.7 90.91 物价物价 指数指数 国民生国民生 产总值产总值 投资额投资额年份年份 序号序号 物价物价 指数指数 国民生产国民生产 总值总值 投资额投资额年份年份 序号序号 根据对未来GNPGNP及PIPI的估计，预测未来投资额 该地区该地区连续连续20年的统计数据年的统计数据 基本回归模型基本回归模型 投资额与投资额与 GNP及物价指数间均有很强的线性关系及物价指数间均有很强的线性关系 tttt xxy 221

37、10 t 年份，年份， yt 投资额，投资额，x1t GNP, x2t 物价指数物价指数 0, 1, 2 回归系数回归系数 x1t yt x2t yt t 对对t相互相互独立的零均值正态随机变量独立的零均值正态随机变量 x1= 596.7 637.7 691.1 756.0 799.0 873.4 944.0 992.7 1077.6 1185.9. 1326.4 1434.2 1549.2 1718.0 1918.3 2163.9 2417.8 2631.7 2954.7 3073.0 ; x2=0.7167 0.7277 0.7436 0.7676 0.7906 0.8254 0.8679

38、 0.9145 0.9601 1.0. 1.0575 1.1508 1.2579 1.3234 1.4005 1.5042 1.6342 1.7842 1.9514 2.0688 ; Y=90.9 97.4 113.5 125.7 122.8 133.3 149.3 144.2 166.4 195.0. 229.8 228.7 206.1 257.9 324.1 386.6 423 401.9 474.9 424.5 ; X=ones(20,1) x1 x2; b,bint,r,rint,stats=regress(Y,X) 用用matlab统计工具箱中统计工具箱中regress函数函数求解过

39、程：求解过程： 2468101214161820 -40 -30 -20 -10 0 10 20 30 40 Residual Case Order Plot Residuals Case Number 0 1000 2000 3000 4000 0.5 1 1.5 2 2.5 -3000 -2000 -1000 0 1000 2000 基本回归模型的结果与分析基本回归模型的结果与分析 ttt xxy 21 479.8596185. 0725.322 MATLAB 统计工具箱统计工具箱 参数参数参数估计值参数估计值置信区间置信区间 0322.7250224.3386 421.1114 10.6

40、1850.4773 0.7596 2-859.4790-1121.4757 -597.4823 R2= 0.9908 F= 919.8529 p=0.0000 剩余标准差剩余标准差 s=12.7164 没有考虑时间序列数据的没有考虑时间序列数据的滞后性影响滞后性影响 R20.9908，拟合度高，拟合度高模型优点模型优点 模型缺点模型缺点 可能忽视了随机误差存在可能忽视了随机误差存在自相关自相关；如果；如果 存在自相关性，用此模型会有不良后果存在自相关性，用此模型会有不良后果 一维插值的定义一维插值的定义 已知已知 n+1个节点个节点, 1 , 0(),(njyx jj 其中其中 j x 互不相

41、同，不妨设互不相同，不妨设), 10 bxxxa n 求任一插值点求任一插值点)( * j xx 处的插值处的插值. * y 0 x 1 x n x 0 y 1 y 节点可视为由节点可视为由 )(xgy 产生产生,， g 表达式复杂表达式复杂,， 或无封闭形式或无封闭形式,， 或未知或未知.。 * x * y Eg 601 Eg 601 用用matlabmatlab解插值问题解插值问题 构造一个构造一个(相对简单的相对简单的)函数函数),(xfy 通过全部节点通过全部节点, 即即 ), 1 ,0()(njyxf jj 再用再用 )(xf 计算插值，即计算插值，即).( * xfy 0 x 1

42、x n x 0 y 1 y * x * y 用用MATLABMATLAB作插值计算作插值计算 一维插值函数：一维插值函数： yi=interp1(x，y，xi，method) 插值方法插值方法被插值点被插值点插值节点插值节点 xixi处的插处的插 值结果值结果 nearest ：最邻近插值：最邻近插值linear ： 线性插值；线性插值； spline ： 三次样条插值；三次样条插值； cubic ： 立方插值。立方插值。 缺省时：缺省时： 分段线性插值。分段线性插值。 注意：所有的插值方法都要求注意：所有的插值方法都要求x x是单调的，并且是单调的，并且xi不不 能够超过能够超过x的范围。的

43、范围。 例例1 1：在：在1-121-12的的1111小时内，每隔小时内，每隔1 1小时测量一次小时测量一次 温度，测得的温度依次为：温度，测得的温度依次为：5 5，8 8，9 9，1515，2525，2929， 3131，3030，2222，2525，2727，2424。试估计每隔。试估计每隔1/101/10小时的小时的 温度值。温度值。 hours=1:12; temps=5 8 9 15 25 29 31 30 22 25 27 24; h=1:0.1:12; t=interp1(hours,temps,h,spline); % (直接输出数据将是很多的) plot(hours,temp

44、s,+,h,t,hours,temps,r:) %作图 xlabel(Hour),ylabel(Degrees Celsius) 024681012 5 10 15 20 25 30 35 Hour Degrees Celsius 例例 2 山区地貌：山区地貌： 在某山区测得一些地点的高程如下表。平面区域为在某山区测得一些地点的高程如下表。平面区域为 1200=x=4000,1200=y=3600) 试作出该山区的地貌图和等高线图，并对几种插值方法进行比较。试作出该山区的地貌图和等高线图，并对几种插值方法进行比较。 X Y 1200 1600 2000 2400 2800 3200 3600

45、4000 1200 1130 1250 1280 1230 1040 900 500 700 1600 1320 1450 1420 1400 1300 700 900 850 2000 1390 1500 1500 1400 900 1100 1060 950 2400 1500 1200 1100 1350 1450 1200 1150 1010 2800 1500 1200 1100 1550 1600 1550 1380 1070 3200 1500 1550 1600 1550 1600 1600 1600 1550 3600 1480 1500 1550 1510 1430 130

46、0 1200 980 通过此例对最近邻点插值、双线性插值方法和双三次插值方法的插 值效果进行比较。 0 5000 10000 0 5000 0 1000 2000 X Y Z 0 5000 10000 0 5000 0 1000 2000 X Y Z 0 5000 10000 0 5000 0 1000 2000 X Y Z 0 5000 10000 0 5000 0 1000 2000 X Y Z 05000 0 500 1000 1500 2000 2500 3000 3500 4000 4500 05000 0 500 1000 1500 2000 2500 3000 3500 4000

47、 4500 05000 0 500 1000 1500 2000 2500 3000 3500 4000 4500 1.原始数据山区表面 2.最近邻点插值山区表面 3.双线性插值山区表面 4.双三次插值山区表面 5.三种插值法的等高线 拟合与插值的关系拟合与插值的关系 函数插值与曲线拟合都是要根据一组数据构造一个函数作为近似，由于近似的函数插值与曲线拟合都是要根据一组数据构造一个函数作为近似，由于近似的 要求不同，二者的数学方法上是完全不同的。要求不同，二者的数学方法上是完全不同的。 实例：实例：下面数据是某次实验所得，希望得到X和 f之间的关系？ x124791 21 31 51 7 f1

48、.53 .96 .611 .71 5 .61 8 .81 9 .62 0 .62 1 .1 问题：问题：给定一批数据点，需确定满足特定要求的曲线或曲面 解决方案：解决方案： 若不要求曲线（面）通过所有数据点，而是要求它反映对象整体的变化趋势，这 就是数据拟合数据拟合，又称曲线拟合或曲面拟合。 若要求所求曲线（面）通过所给所有数据点，就是插值问题插值问题； Eg 701 Eg 701 用用matlabmatlab解拟合问题解拟合问题 最邻近插值、线性插值、样条插值与曲线拟合结果：最邻近插值、线性插值、样条插值与曲线拟合结果： 024681012141618 0 5 10 15 20 25 已知数

49、据点 spline 三次多项式插值 024681012141618 0 5 10 15 20 25 已知数据点 linest 三次多项式插值 024681012141618 0 5 10 15 20 25 已知数据点 nearest 三次多项式插值 MATLAB的符号计算的符号计算 所谓符号计算是指在运算时,无须事先对变量赋值,而将所得到结果以标准的 符号形式来表示。 MathWorks公司以Maple的内核作为符号计算引擎（Engine），依赖 Maple已有的函数库，开发了实现符号计算的两个工具箱：基本符号工具 箱和扩展符号工具箱。 一、符号计算基础 一、符号计算基础 （一） 定义符号变量

50、 参与符号运算的对象可以是符号变量、符号表达式或符号矩阵。符号变量要 先定义，后引用。可以用sym函数、syms函数将运算量定义为符号型数据。 引用符号运算函数时，用户可以指定函数执行过程中的变量参数；若用户没 有指定变量参数，则使用findsym函数默认的变量作为函数的变量参数。 1、sym函数 sym函数的主要功能是创建符号变量，以便进行符号运 算，也可以用于创建符号表达式或符号矩阵。用sym函 数创建符号变量的一般格式为： x = sym(x) 其目的是将x创建为符号变量，以x作为输出变量名。 每次调用该函数,可以定义一个符号变量。 一、符号计算基础 （一） 定义符号变量 【例1】作符号

51、计算： a,b,x,y均为符号运算量。在符号运算前，应先将a,b,x,y定义为符号运算 量 1 5 byax byax 一、符号计算基础 （一） 定义符号变量 a=sym(a); %定义a为符号运算量，输出变量名为a b=sym(b); x=sym(x); y=sym(y”); x,y=solve(a*x-b*y-1,a*x+b*y-5,x,y) %以a,b为符号常数，x,y为 符号变量 即可得到方程组的解： x =3/a y =2/b 一、符号计算基础 （一） 定义符号变量 【例2】已知一复数表达式 z=x+i*y, 试求其共轭复数,并求该表达式与其共轭复数 乘积的多项式。 为了使乘积表达式

52、x2+y2非负，这里，把变量x和y定义为实数。 x=sym(x,real); y=sym(y,real); 一、符号计算基础 （一） 定义符号变量 z=x+i*y; %定义复数表达式 conj(z); %求共轭复数 expand(z*conj(z) %求表达式与其共轭复数乘积的多项式 ans = x2+y2 若要去掉x的属性,可以使用下面语句 x = sym(x,unreal) 将x创建为纯格式的符号变量。 一、符号计算基础 （一） 定义符号变量 2、syms函数 syms函数的功能与sym函数类似。syms函数可以在一个语句中同时定义多 个符号变量，其一般格式为： syms arg1 arg

53、2 argN 用于将rg1, arg2,argN等符号创建为符号型数据。 一、符号计算基础 （一） 定义符号变量 （二）默认符号变量 在数学表达式中，一般习惯于使用排在字母表中前面的字母作为变量的系数， 而用排在后面的字母表示变量。例如： f=ax2+bx+c 表达式中的a,b,c通常被认为是常数，用作变量的系数；而将x看作自变量。 一、符号计算基础 例如，数学表达式 f=xn g=sin(at+b) 根据数学式中表示自变量的习惯，默认a,b,c为符号常 数，x为符号变量。 若在MATLAB中表示上述表达式，首先用syms 函数 定义a，b，n，t，x为符号对象。在进行导数运算时， 由于没有指

54、定符号变量，则系统采用数学习惯来确定 表达式中的自变量，默认a,b,c为符号常数，x，t为符 号变量。 即 ： 对函数f求导为：df/dx 对函数g求导为：dg/dt 一、符号计算基础 （二）默认符号变量 为了了解函数引用过程中使用的符号变量个数及变 量名，可以用findsym函数查询默认的变量。该函 数的引用格式为： findsym（f,n） 说明：f为用户定义的符号函数， n为正整数，表示查询变量的个数。 n=i,表示查询i个系统默认变量。n值省略时表示 查询符号函数中全部系统默认变量。 一、符号计算基础 （二）默认符号变量 【例3 】查询符号函数 f=xn g=sin(at+b) 中的系

55、统默认变量。 syms a b n t x %定义符号变量 f=xn; %给定符号函数 g=sin(a*t+b); findsym(f,1) %在f函数中查询1个系统默认变量 ans= x 表示f函数中查询的1个系统默认变量为x。 一、符号计算基础 （二）默认符号变量 （三） 符号表达式 符号表达式由符号变量、函数、算术运算符等组成。 符号表达式的书写格式与数值表达式相同。例如,数 学表达式 其符号表达式为： 1+sqr(5*x)/2 注意，在定义表达式前应先将表达式中的字符x定义 为符号变量。 2 51x 一、符号计算基础 （四） 生成符号函数 将表达式中的自变量定义为符号变量后，赋值给符号

56、函数名，即可生成符号 函数。例如有一数学表达式： 2 22 ),( c byax yxf 一、符号计算基础 其用符号表达式生成符号函数fxy的过程为： syms a b c x y %定义符号运算量 fxy=(a*x2+b*y2)/c2 %生成符号函数 生成符号函数fxy后，即可用于微积分等符号计算。 一、符号计算基础 （四） 生成符号函数 【例4】定义一个符号函数 fxy=(a*x2+b*y2)/c2 ，分别求该 函数对x、y的导数和对x的积分。 syms a b c x y %定义符号变量 fxy=(a*x2+b*y2)/c2； %生成符号函数 diff(fxy，x) %符号函数fxy对x

57、求导数 ans =2*a*x/c2 diff(fxy， y) %符号函数fxy对y求导数 ans =2*b*y/c2 %符号函数fxy对x求积分 int(fxy， x) ans =1/c2*(1/3*a*x3+b*y2*x) 一、符号计算基础 （四） 生成符号函数 二、微积分 （一） 微积分函数 1.求极限 函数limit用于求符号函数f的极限。系统可以根据用户要求，计算变量从不 同方向趋近于指定值的极限值。该函数的格式及功能： 二、微积分 limit(f,x,a)：求符号函数f（x）的极限值。即计算当变量x趋近于常 数a时，f（x）函数的极限值。 limit(f,a)：求符号函数f（x）的极

58、限值。由于没有指定符号函数f （x）的自变量，则使用该格式时，符号函数f（x）的变量为函数 findsym(f)确定的默认自变量，既变量x趋近于a。 limit(f)：求符号函数f（x）的极限值。符号函数f（x）的变量为函 数findsym(f)确定的默认变量；没有指定变量的目标值时，系统默 认变量趋近于0，即a=0的情况。 limit(f,x,a,right)：求符号函数f的极限值。right表示变量x从右边 趋近于a。 limit(f,x,a,left)：求符号函数f的极限值。left表示变量x从左边趋 近于a。 二、微积分 【例5】求极限 syms x； %定义符号变量 f=(x*(ex

59、p(sin(x)+1)-2*(exp(tan(x)-1)/sin(x)3； %确定符 号表达式 w=limit(f) %求函数的极限 w = -1/2 x eex tgxx x 3 sin 0 sin ) 1(2) 1( lim 二、微积分 2. 微分函数 diff函数用于对符号表达式s求微分。该函数的一般 引用格式为： diff(s,v,n) 二、微积分 说明： 应用diff（s）没有指定微分变量和微分阶数，则系统按findsym函数指示的默 认变量对符号表达式s求一阶微分。 应用diff（s，v）或diff（s，sym（v） 格式，表示以v为自变量，对符 号表达式s求一阶微分。 应用dif

60、f（s，n）格式，表示对符号表达式s求n阶微分，n为正整数。 应用diff（s，v，n）diff（s，n，v） 格式，表示以v为自变量，对符号表 达式s求n阶微分。 【例6】求导数： x = sym(x); %定义符号变量 t = sym(t); diff(sin(x2) %求导运算 ans = 2*cos(x2)*x dx xd 2 sin 二、微积分 3积分函数 积分函数int（s ，v，a，b)可以对被积函数或符号表达式s求积分。其引用格式 为： int（s ，v，a，b) 说明： 应用int（s）格式，表示没有指定积分变量和积分阶数时，系统 按findsym函数指示的默认变量对被积函数