差商matlab编程,Matlab数值计算差商与插值均差定义 若已知函数f(x)在点x0,x1,...xn处的值f(x0),f(x1),...f(xn).如果i≠j，则
⼀阶均差f[xj,xj+1]=f(xj+1)−f(xj)xj+1−xj(j=0,1,...n−1)
⼆阶均差f[xj,xj+1,xj+2]=f[xj+1,xj+2]−f[xj,xj+1]xj+2−xj(j=0,1,...n−2)
n阶均差f[x0,x1,...,xn]=f[x1,...,xn]−f[x0,...,xn−1]xn−x0
例 由函数表求各阶均差
x
-2
-1
1
3
y
-56
-16
-2
-2
4
解：按公式计算⼀阶差商、⼆阶差商、三阶差商如下
x
f(x)
⼀阶差商
⼆阶差商
三阶差商
-
2
-56
-1
-16
40
-2
14
-13
1
-2
-7
2
3
4
3
1
2
Matlab代码
clear
clc
x=[-2 -1 0 1 3]
y=[-56 -16 -2 -2 4]
deltx=diff(x);
delty=diff(y);
firstorder=delty./deltx %⼀阶
for i=1:length(x)-2
delt2x(i)=x(i+2)-x(i);
end
delt2y=diff(firstorder); secondorder=delt2y./delt2x %⼆阶for i=1:length(x)-3
delt3x(i)=x(i+3)-x(i);
end
delt3y=diff(secondorder); thirdorder=delt3y./delt3x %三阶for i=1:length(x)-4diff函数
delt4x(i)=x(i+4)-x(i);
end
delt4y=diff(thirdorder);
fourorder=delt4y./delt4x %四阶
结果
x =
-2 -1 0 1 3
y =
-56 -16 -2 -2 4
firstorder =
40 14 0 3
secondorder =
-
13 -7 1
thirdorder =
2 2
fourorder =
这⾥⽤到了diff，就再次介绍⼀下差分函数
补充：差分函数diff
diff(X) X为向量时(⾏列均可)，计算相邻两数的差[X(2)-X(1) X(3)-X(2) … X(n)-X(n-1)]
diff(X) X为矩阵时，计算矩阵的2~n⾏与1~n-1⾏的差，[X(2:n,:) - X(1:n-1,:)]
diff(X,N) 对上⾯函数diff(X)的扩充，这⾥的N指定N阶差分，⼆阶差分是对⼀阶差分的结果再做差分运算
DIFF(X,N,DIM) 对上⾯函数diff(X,N)的扩充，DIM取1或2，取1时按⾏差分，与上⾯结果⼀样，取2时按列差分把上⾯的命令⽤字符串改造了⼀下，不过太难看懂了,no zuo no die
eval()函数的功能就是将括号内的字符串视为语句并运⾏，简单记为字符串转语句
num2str()函数的功能就是将括号内的数字转换为字符串，简单记为数字转字符串
clear
clc
x=[-2 -1 0 1 3]
y=[-56 -16 -2 -2 4]
deltx=diff(x);
delty=diff(y);
order1=delty./deltx %⼀阶
for j=2:4
str=['for i=1:length(x)-',num2str(j),char(10)];
str=[str,'delt',num2str(j),'x(i)=x(i+',num2str(j),')-x(i);',char(10)];
str=[str,'end',char(10)];
str=[str,'delt',num2str(j),'y=diff(order',num2str(j-1),');',char(10)];
str=[str,'order',num2str(j),'=delt',num2str(j),'y./delt',num2str(j),'x',char(10)];
eval(str)
x =
-2 -1 0 1 3
y =
-56 -16 -2 -2 4
order1 =
40 14 0 3
order2 =
-13 -7 1
order3 =
2 2
order4 =
⽜顿插值
⽜顿插值公式及其余项
当n=1时：
差商N1(x)=f(x0)+f[x0,x1](x−x0)
余项R1(x)=f[x,x0,x1](x−x0)(x−x1)
当n=2时：
差商N2(x)=f(x0)+f[x0,x1](x−x0)+f[x0,x1,x2](x−x0)(x−x1)余项R2(x)=f[x,x0,x1,x2](x−x0)(x−x1)(x−x2)
n阶：
差商Nn(x)=f(x0)+f[x0,x1](x−x0)+f[x0,x1,x2](x−x0)(x−x1) +f[x0,x1,...,xn](x−x0)(x−x1)...(x−xn−1)
余项Rn(x)=f[x,x0,x1,...,xn](x−x0)(x−x1)...(x−xn)
例 已知x=1,4,9的平⽅根为1,2,3，利⽤⽜顿基本差商公式7√的值解:
xi
xi−−√
f[xi,xi+1]
f[xi,xi+1,xi+2]
1
1
2−14−1=0.33333
9
3
3−29−4=0.2
0.2−0.333339−1=−0.016
从⽽得⼆阶⽜顿基本差商公式为
P2(x)=1+0.33333(x−1)−0.01667(x−1)(x−4)
因此计算得7√的近似值为P2(7)=2.69992
clear
clc
x=[1 4 9]
y=[1 2 3]
deltx=diff(x);
delty=diff(y);
order1=delty./deltx %⼀阶
for i=1:length(x)-2
delt2x(i)=x(i+2)-x(i);
end
delt2y=diff(order1);
order2=delt2y./delt2x %⼆阶
%%⽜顿插值需要的值是y(1)、order1(1)、order2(1)、x(1)、x(2) y(1),order1(1),order2(1),x(1),x(2) %%构造多项式P=[0 0 y(1)]+[0 order1(1)*poly(x(1))]+order2(1)*poly(x([1:2]))
%%求值 polyval(P,7)
结果
x =
1 4 9
y =
1 2 3
order1 =
0.3333 0.2000
order2 =
-0.0167

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。