Matlab课件——第5章 数值计算
合集下载
相关主题
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
• method用于指定插值的方法，有‘nearest’, ‘linear’, ‘spline’,
‘pchip’, ‘cubic’, ‘v5cubic’等方法。各种插值方法在速度 上,Nearest最快,然后是Linear再到Cubic,最慢的是 Spline.但是精度和曲线的平滑度恰好相反,Nearest甚至不
• 用法： x=lsqcurvefit(fun,x0,xdata,ydata) x=lsqcurvefit(fun,x0,xdata,ydata,lb,ub) x=lsqcurvefit(fun,x0,xdata,ydata,lb,ub,options) [x,resnorm]=lsqcurvefit(fun,x0,xdata,ydata,...) [x,resnorm,residual,exitflag]=lsqcurvefit(fun,x0,xdata ,ydata,...)
例子：
>> P=1:6, poly2str(P,'x'), x=[2 3;4 5], polyval(P,x) P= 1 2 3 4 5 6 ans = x^5 + 2 x^4 + 3 x^3 + 4 x^2 + 5 x + 6 x= 2 3 4 5 ans = 120 543 1818 4881
Matlab提供了下列关于多项式的函数：
1.多项式表示法
MATLAB采用行向量表示多项式系数，多 项式系数按降幂排列。
poly2str函数将多项式系数向量转换为完 整形式。
>> poly2str([1 0 2],'s') ans = s^2 + 2
2.多项式求值
polyval函数计算多项式的值，其用法为： Y = polyval(P,X),
>> P=[3 0.5 7 -0.09],X=polyder(P),poly2str(X,'x') P= 3.0000 0.5000 7.0000 -0.0900 X= 9 1 7 ans = 9 x^2 + x + 7
（2）多项式的积分
polyint函数计算多项式的不定积分， polyint(P,C)
第五章 数值计算
目 录
5.1 多项式运算
5.2 插 值 运 算 5.3 数 据 分 析 5.4 功 能 函 数
5.5 微分方程组数值解
习 题
5.1 多项式运算
1.多项式Matlab表示法 2.多项式求值
3.多项式乘法和除法
4.多项式的微积分
5.多项式的根和由根创建多项式
6.多项式部分分式展开 7.多项式曲线拟合 8.多项式曲线拟合图形用户接口
插值函数
5.2.1 一维插值
一维插值就是对函数 y=f(x) 进行插值，一维 插值的原理如下图所示。
interp1函数实现一维插值，其用法如下：
yi=interp1(x,y,xi)，x、y是已知数据集且具有相同长度的向
量；
yi=interp1(y,xi)，默认x为1:n，其中n为向量y的长度； yi=interp1(x,y,xi,method)
p = poly(A)，输入变量A是方阵，返回值为A 的特征多项式的系数向量。
6.多项式部分分式展开
• 部分分式展开又叫部分分式分解，是将有理函
数分解成许多次数较低有理函数和的形式，来降 低分子或分母多项式的次数，
• 例如：
• 分解后的分式需满足以下条件：
–分式的分母需为不可约多项式；
–分式的分子多项式次数需比其分母多项式次数要低。
4.多项式的微积分
（1）多项式的微分
polyder函数计算多项式的微分 polyder(P)，返回多项式P微分的系数向量；
polyder(A,B)，返回多项式A*B微分的系数向量；
[Q,D]=polyder(B,A)，返回多项式B/A微分的系
数向量。
实例
• 如何用MATLAB对一个已知的函数： y=3*x^3+0.5x^2+7*x-0.09进行求导，并 分别作出求导前和求导后的相应曲线。
连续~~系统默认的是Linear。
内插运算与外插运算
(1)只对已知数据点集内部的点进行的插值运算称 为内插，可比较准确的估测插值点上的函数值；
(2)当插值点落在已知数据集的外部时的插值称为 外插，要估计外插函数值很难。 y=interp1([x,]y,xi,[method],['extrap'],[extrapval]) []代表可选。 • MATLAB对已知数据集外部点上函数值的预测都 返回NaN，但可通过为interp1函数添加'extrap'参 数指明用于外插，MATLAB的外插结果偏差较大。
B(x) R(1) R(2) R(n) ---- = -------- + -------- + ... + -------- + K(x) A(x) x - P(1) x - P(2) x - P(n)
7.多项式曲线拟合
已知在点集上的函数值,构造一个解析函数(其图形 为一曲线)使在原离散点上尽可能接近给定的值，
Y = polyval(P,X),
把矩阵或向量X中的每个元素逐个代入多项式中进行
计算；
Y = polyvalm(P,X)，把矩阵X作为整体代入多项式 中进行计算，X必须为方阵。
例子：
>> P=1:6, poly2str(P,'x'), x=[2 3;4 5], polyval(P,x) P= 1 2 3 4 5 6 ans = x^5 + 2 x^4 + 3 x^3 + 4 x^2 + 5 x + 6 x= 2 3 4 5 ans = 120 543 1818 4881
>> P=[3 0.5 7 -0.09], x=roots(P) P= 3.0000 0.5000 7.0000 -0.0900 x= -0.0898 + 1.5256i -0.0898 - 1.5256i 0.0128
（2）由根创建多项式
poly函数实现由根创建多项式，其用法如下：
p = poly(r)，输入变量r是多项式的所有根， 返回值为多项式的系数向量，与roots函数是 两个可逆的过程；
3.多项式乘法和除法
deconv函数进行反卷积(deconvolution)和多项式
除法运算， [Q,R] = deconv(B,A)
对向量A和B进行反卷积运算，返回结果为向量Q和余量 R; 如果A和B是多项式的系数向量，它们的反卷积相当于两 个多项式相除，A是被除数，B是除数，Q是商，R是余 数。
ห้องสมุดไป่ตู้
3.多项式乘法和除法
conv函数进行卷积(convolution)和多项 式乘法运算，
C = conv(A, B)
对向量A和B进行卷积，返回结果为长度为 length(A)+length(B)-1的向量; 如果A和B是多项式的系数向量，它们的卷积相 当于两个多项式相乘。
向量卷积
所谓两个向量卷积，就是多项式乘法。 比如:p=[1 2 3], q=[1 1] 是两个向量， p和q的卷积如下： • 把p的元素作为一个多项式的系数按降幂排列，写 出对应的多项式：x^2+2x+3; • 把q的元素也作为多项式的系数按降幂排列，写出 对应的多项式：x+1； • 对两个多项式相乘， (x^2+2x+3)×(x+1)=x^3+3x+5x^2+3 • 取系数，所以p和q卷积的结果就是[1 3 5 3]。
P为多项式系数向量；
C为不定积分常数项，为标量。
polyint(P)，假设 C=0.
返回值为多项式不定积分的系数向量。
5.多项式的根和由根创建多项式
（1）多项式的根
roots函数用于求多项式的根，其用法如下： r = roots(C) 返回以C向量为系数的多项式的所有根r。
实例
• 求方程 3*x^3+0.5x^2+7*x-0.09=0的根。
plot(xdata,ydata,'ro',xdata,yfit) legend('ydata', 'yfit')
5.2 插 值 运 算
插值运算是根据数据的分布规律，找到一个函数表达 式可以连接已知的各点，并用此函数表达式预测两点
之间任意位置上的函数值。
5.2.1 一维插值
5.2.2 二维插值
插值是根据已知输入/输出数据集和当前输 入估计输出值。MATLAB提供大量的插值函 数，如下表所示。
这一过程称为曲线拟合。
polyfit函数采用最小二乘法对给定数据进行多项式 拟合，其用法如下： • p=polyfit(x,y,n) • [p,s]= polyfit(x,y,n) • 说明：x,y为数据点，n为多项式阶数，返回p为 幂次从高到低的多项式系数向量p。x必须是单调 的。矩阵s用于生成预测值的误差估计。
8.多项式曲线拟合图形用户接口
曲线拟合的图形用户接口可通过图形窗口的
【Tools】菜单中【Basic Fitting】选项启动。
>> x=0:0.1:10; >> y=0.25*x+20*sin(x); >> plot(x,y,'ro')
最小二乘法曲线拟合
• 利用lsqcurvefit函数进行任意表达式曲线拟合 • 原理： min sum {(fun(x,xdata)-ydata).^2}
实例
xdata = 0:0.1:10; ydata = 3*sin(xdata)+6; fun = inline('x(1)*sin(xdata)+x(2)','x','xdata'); x = lsqcurvefit(fun,[0 0], xdata, ydata) yfit=feval(fun,x,xdata);
• P——多项式系数行向量
• X——代入多项式的值
其意义为：
当P为长度为n+1的行向量 [p1,p2, ...,pn,pn+1] 时，
Y = p1xn + p2xn-1 + ... + pnx + pn+1
例子：
>> P=1:6, poly2str(P,'x'), x=5, polyval(P,x) P= 1 2 3 4 5 6 ans = x^5 + 2 x^4 + 3 x^3 + 4 x^2 + 5 x + 6 x= 5 ans = 4881
实 例 1
给出如下数据的拟合曲线，x=[0.5,1.0,1.5,2.0,2.5,3.0]， y=[1.75,2.45,3.81,4.80,7.00,8.60] MATLAB程序如下： x=[0.5,1.0,1.5,2.0,2.5,3.0]; y=[1.75,2.45,3.81,4.80,7.00,8.60]; p=polyfit(x,y,2) x1=0.5:0.05:3.0; y1=polyval(p,x1); plot(x,y,'*r',x1,y1,'-b') 计算结果为： p =0.5614 0.8287 1.1560 即所得多项式为y=0.5614x^2+0.08287x+1.15560
问题：已知相渗曲线（如图所示），原油的密度为0.9，粘 度为9，地层水的密度为1，粘度为0.45，计算一定饱和度 时的含水率。
实 例
Sw 0.290 0.356 0.449 0.496 0.543 0.590 0.684 0.730 0.780 Kro Krw 1 0.000 0.63 0.017 0.35 0.069 0.27 0.113 0.2 0.165 0.133 0.227 0.05 0.372 0.01 0.456 0 0.546
6.多项式部分分式展开
residue函数将多项式之比按部分分式展开， 其用法如下：
[R,P,K] = residue(B,A) ，求多项式B/A的 部分分式展开(分解) [B,A]=residue(R,P,K)，从部分分式得到多 项式向量。
6.多项式部分分式展开
[R,P,K] = residue(B,A)
例子：
>> P=1:6, poly2str(P,'x'), x=[2 3;4 5], polyvalm(P,X) P= 1 2 3 4 5 6 ans = x^5 + 2 x^4 + 3 x^3 + 4 x^2 + 5 x + 6 x= 2 3 4 5 >> x^5 + 2*x^4 + 3*x^3 + 4*x^2 + 5*x + 6 ans = 8256 10881 ans = 8256 10887 14508 19137 14514 19137
‘pchip’, ‘cubic’, ‘v5cubic’等方法。各种插值方法在速度 上,Nearest最快,然后是Linear再到Cubic,最慢的是 Spline.但是精度和曲线的平滑度恰好相反,Nearest甚至不
• 用法： x=lsqcurvefit(fun,x0,xdata,ydata) x=lsqcurvefit(fun,x0,xdata,ydata,lb,ub) x=lsqcurvefit(fun,x0,xdata,ydata,lb,ub,options) [x,resnorm]=lsqcurvefit(fun,x0,xdata,ydata,...) [x,resnorm,residual,exitflag]=lsqcurvefit(fun,x0,xdata ,ydata,...)
例子：
>> P=1:6, poly2str(P,'x'), x=[2 3;4 5], polyval(P,x) P= 1 2 3 4 5 6 ans = x^5 + 2 x^4 + 3 x^3 + 4 x^2 + 5 x + 6 x= 2 3 4 5 ans = 120 543 1818 4881
Matlab提供了下列关于多项式的函数：
1.多项式表示法
MATLAB采用行向量表示多项式系数，多 项式系数按降幂排列。
poly2str函数将多项式系数向量转换为完 整形式。
>> poly2str([1 0 2],'s') ans = s^2 + 2
2.多项式求值
polyval函数计算多项式的值，其用法为： Y = polyval(P,X),
>> P=[3 0.5 7 -0.09],X=polyder(P),poly2str(X,'x') P= 3.0000 0.5000 7.0000 -0.0900 X= 9 1 7 ans = 9 x^2 + x + 7
（2）多项式的积分
polyint函数计算多项式的不定积分， polyint(P,C)
第五章 数值计算
目 录
5.1 多项式运算
5.2 插 值 运 算 5.3 数 据 分 析 5.4 功 能 函 数
5.5 微分方程组数值解
习 题
5.1 多项式运算
1.多项式Matlab表示法 2.多项式求值
3.多项式乘法和除法
4.多项式的微积分
5.多项式的根和由根创建多项式
6.多项式部分分式展开 7.多项式曲线拟合 8.多项式曲线拟合图形用户接口
插值函数
5.2.1 一维插值
一维插值就是对函数 y=f(x) 进行插值，一维 插值的原理如下图所示。
interp1函数实现一维插值，其用法如下：
yi=interp1(x,y,xi)，x、y是已知数据集且具有相同长度的向
量；
yi=interp1(y,xi)，默认x为1:n，其中n为向量y的长度； yi=interp1(x,y,xi,method)
p = poly(A)，输入变量A是方阵，返回值为A 的特征多项式的系数向量。
6.多项式部分分式展开
• 部分分式展开又叫部分分式分解，是将有理函
数分解成许多次数较低有理函数和的形式，来降 低分子或分母多项式的次数，
• 例如：
• 分解后的分式需满足以下条件：
–分式的分母需为不可约多项式；
–分式的分子多项式次数需比其分母多项式次数要低。
4.多项式的微积分
（1）多项式的微分
polyder函数计算多项式的微分 polyder(P)，返回多项式P微分的系数向量；
polyder(A,B)，返回多项式A*B微分的系数向量；
[Q,D]=polyder(B,A)，返回多项式B/A微分的系
数向量。
实例
• 如何用MATLAB对一个已知的函数： y=3*x^3+0.5x^2+7*x-0.09进行求导，并 分别作出求导前和求导后的相应曲线。
连续~~系统默认的是Linear。
内插运算与外插运算
(1)只对已知数据点集内部的点进行的插值运算称 为内插，可比较准确的估测插值点上的函数值；
(2)当插值点落在已知数据集的外部时的插值称为 外插，要估计外插函数值很难。 y=interp1([x,]y,xi,[method],['extrap'],[extrapval]) []代表可选。 • MATLAB对已知数据集外部点上函数值的预测都 返回NaN，但可通过为interp1函数添加'extrap'参 数指明用于外插，MATLAB的外插结果偏差较大。
B(x) R(1) R(2) R(n) ---- = -------- + -------- + ... + -------- + K(x) A(x) x - P(1) x - P(2) x - P(n)
7.多项式曲线拟合
已知在点集上的函数值,构造一个解析函数(其图形 为一曲线)使在原离散点上尽可能接近给定的值，
Y = polyval(P,X),
把矩阵或向量X中的每个元素逐个代入多项式中进行
计算；
Y = polyvalm(P,X)，把矩阵X作为整体代入多项式 中进行计算，X必须为方阵。
例子：
>> P=1:6, poly2str(P,'x'), x=[2 3;4 5], polyval(P,x) P= 1 2 3 4 5 6 ans = x^5 + 2 x^4 + 3 x^3 + 4 x^2 + 5 x + 6 x= 2 3 4 5 ans = 120 543 1818 4881
>> P=[3 0.5 7 -0.09], x=roots(P) P= 3.0000 0.5000 7.0000 -0.0900 x= -0.0898 + 1.5256i -0.0898 - 1.5256i 0.0128
（2）由根创建多项式
poly函数实现由根创建多项式，其用法如下：
p = poly(r)，输入变量r是多项式的所有根， 返回值为多项式的系数向量，与roots函数是 两个可逆的过程；
3.多项式乘法和除法
deconv函数进行反卷积(deconvolution)和多项式
除法运算， [Q,R] = deconv(B,A)
对向量A和B进行反卷积运算，返回结果为向量Q和余量 R; 如果A和B是多项式的系数向量，它们的反卷积相当于两 个多项式相除，A是被除数，B是除数，Q是商，R是余 数。
ห้องสมุดไป่ตู้
3.多项式乘法和除法
conv函数进行卷积(convolution)和多项 式乘法运算，
C = conv(A, B)
对向量A和B进行卷积，返回结果为长度为 length(A)+length(B)-1的向量; 如果A和B是多项式的系数向量，它们的卷积相 当于两个多项式相乘。
向量卷积
所谓两个向量卷积，就是多项式乘法。 比如:p=[1 2 3], q=[1 1] 是两个向量， p和q的卷积如下： • 把p的元素作为一个多项式的系数按降幂排列，写 出对应的多项式：x^2+2x+3; • 把q的元素也作为多项式的系数按降幂排列，写出 对应的多项式：x+1； • 对两个多项式相乘， (x^2+2x+3)×(x+1)=x^3+3x+5x^2+3 • 取系数，所以p和q卷积的结果就是[1 3 5 3]。
P为多项式系数向量；
C为不定积分常数项，为标量。
polyint(P)，假设 C=0.
返回值为多项式不定积分的系数向量。
5.多项式的根和由根创建多项式
（1）多项式的根
roots函数用于求多项式的根，其用法如下： r = roots(C) 返回以C向量为系数的多项式的所有根r。
实例
• 求方程 3*x^3+0.5x^2+7*x-0.09=0的根。
plot(xdata,ydata,'ro',xdata,yfit) legend('ydata', 'yfit')
5.2 插 值 运 算
插值运算是根据数据的分布规律，找到一个函数表达 式可以连接已知的各点，并用此函数表达式预测两点
之间任意位置上的函数值。
5.2.1 一维插值
5.2.2 二维插值
插值是根据已知输入/输出数据集和当前输 入估计输出值。MATLAB提供大量的插值函 数，如下表所示。
这一过程称为曲线拟合。
polyfit函数采用最小二乘法对给定数据进行多项式 拟合，其用法如下： • p=polyfit(x,y,n) • [p,s]= polyfit(x,y,n) • 说明：x,y为数据点，n为多项式阶数，返回p为 幂次从高到低的多项式系数向量p。x必须是单调 的。矩阵s用于生成预测值的误差估计。
8.多项式曲线拟合图形用户接口
曲线拟合的图形用户接口可通过图形窗口的
【Tools】菜单中【Basic Fitting】选项启动。
>> x=0:0.1:10; >> y=0.25*x+20*sin(x); >> plot(x,y,'ro')
最小二乘法曲线拟合
• 利用lsqcurvefit函数进行任意表达式曲线拟合 • 原理： min sum {(fun(x,xdata)-ydata).^2}
实例
xdata = 0:0.1:10; ydata = 3*sin(xdata)+6; fun = inline('x(1)*sin(xdata)+x(2)','x','xdata'); x = lsqcurvefit(fun,[0 0], xdata, ydata) yfit=feval(fun,x,xdata);
• P——多项式系数行向量
• X——代入多项式的值
其意义为：
当P为长度为n+1的行向量 [p1,p2, ...,pn,pn+1] 时，
Y = p1xn + p2xn-1 + ... + pnx + pn+1
例子：
>> P=1:6, poly2str(P,'x'), x=5, polyval(P,x) P= 1 2 3 4 5 6 ans = x^5 + 2 x^4 + 3 x^3 + 4 x^2 + 5 x + 6 x= 5 ans = 4881
实 例 1
给出如下数据的拟合曲线，x=[0.5,1.0,1.5,2.0,2.5,3.0]， y=[1.75,2.45,3.81,4.80,7.00,8.60] MATLAB程序如下： x=[0.5,1.0,1.5,2.0,2.5,3.0]; y=[1.75,2.45,3.81,4.80,7.00,8.60]; p=polyfit(x,y,2) x1=0.5:0.05:3.0; y1=polyval(p,x1); plot(x,y,'*r',x1,y1,'-b') 计算结果为： p =0.5614 0.8287 1.1560 即所得多项式为y=0.5614x^2+0.08287x+1.15560
问题：已知相渗曲线（如图所示），原油的密度为0.9，粘 度为9，地层水的密度为1，粘度为0.45，计算一定饱和度 时的含水率。
实 例
Sw 0.290 0.356 0.449 0.496 0.543 0.590 0.684 0.730 0.780 Kro Krw 1 0.000 0.63 0.017 0.35 0.069 0.27 0.113 0.2 0.165 0.133 0.227 0.05 0.372 0.01 0.456 0 0.546
6.多项式部分分式展开
residue函数将多项式之比按部分分式展开， 其用法如下：
[R,P,K] = residue(B,A) ，求多项式B/A的 部分分式展开(分解) [B,A]=residue(R,P,K)，从部分分式得到多 项式向量。
6.多项式部分分式展开
[R,P,K] = residue(B,A)
例子：
>> P=1:6, poly2str(P,'x'), x=[2 3;4 5], polyvalm(P,X) P= 1 2 3 4 5 6 ans = x^5 + 2 x^4 + 3 x^3 + 4 x^2 + 5 x + 6 x= 2 3 4 5 >> x^5 + 2*x^4 + 3*x^3 + 4*x^2 + 5*x + 6 ans = 8256 10881 ans = 8256 10887 14508 19137 14514 19137