数模竞赛-插值与拟合

插值方法

待定系数法插值

略。

拉格朗日插值

例题1：

已知未知函数$y=f(x)$的$6$个观测点，求插值函数$y=\hat{f}(x)$，并求$x=1.5,2.6$处函数估计值。

代码：

1
2
3
4
5
6
7
import numpy,scipy
x0=numpy.arange(1,7)
y0=numpy.array([16,18,21,17,15,12])
p=scipy.interpolate.lagrange(x0,y0)
print(numpy.round(p,4)) #高次幂到低次幂系数
yh=numpy.polyval(p,[1.5,2.6])
print(numpy.round(yh,4))

牛顿插值

略。

分段线性插值

龙格震荡现象：

例题2：在区间$[-5,5]$上，用$n+1$个等距节点作多项式$P_n(x)$，使得它在节点处的值与函数$y=\dfrac{1}{1+x^2}$在对应节点处的值相等，考查$n=6,8,10$时，多项式次数与逼近误差的关系。

代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
import numpy,pylab,scipy
yx=lambda x:1/(1+x**2)
def fun(n):
    x=numpy.linspace(-5,5,n+1)
    p=scipy.interpolate.lagrange(x,yx(x))
    return p
x0=numpy.linspace(-5,5,100)
pylab.rc('text',usetex=True)
pylab.rc('font',size=15)
N=[6,8,10]
s=['--*b','-.','-p']
for k in range(len(N)):
    p=fun(N[k])
    pylab.plot(x0,numpy.polyval(p,x0),s[k])
pylab.plot(x0,yx(x0))
pylab.legend(['$n=6$','$n=8$','$n=10$','$\\frac{1}{1+x^2}$'])
pylab.show()

分段线性插值略。

三次样条插值

略。

二维插值

略。

代码总结

例题3：数据如下，画出曲线，求出$x=0$处曲线斜率和$13\leqslant x\leqslant15$范围内$y$的最小值。

代码：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
import numpy,scipy,pylab
x0=[0,3,5,7,9,11,12,13,14,15]
y0=[0,1.2,1.7,2.0,2.1,2.0,1.8,1.2,1.0,1.6]
x=numpy.linspace(0,15,151)
pylab.rc('font',family='SimHei')
pylab.rc('axes',unicode_minus=False) #用来正常显示负号
pylab.subplots_adjust(wspace=0.5) #各子图水平间距

pylab.subplot(131)
yx1=scipy.interpolate.interp1d(x0,y0)
y1=yx1(x)
pylab.plot(x,y1)
pylab.title('分段线性插值')

pylab.subplot(132)
p2=scipy.interpolate.lagrange(x0,y0)
y2=numpy.polyval(p2,x)
pylab.plot(x,y2)
pylab.title('拉格朗日插值')

pylab.subplot(133)
yx3=scipy.interpolate.interp1d(x0,y0,'cubic')
y3=yx3(x)
pylab.plot(x,y3)
pylab.title('三次样条插值')

pylab.show()

dx=numpy.diff(x)
dy=numpy.diff(y3)
dyx=dy/dx
dyx0=dyx[0] #x=0
xt=x[130:] #>=13
yt=y3[130:] #>=13
ymin=min(yt)
xmin=[xt[ind]for ind,v in enumerate(yt) if v==ymin]
print(dyx0) #x=0处斜率 
print(xmin,ymin)