lagrange.frm

lagrange插值的一种方便算法

VERSION 5.00
Begin VB.Form Form1 
   Caption         =   "Form1"
   ClientHeight    =   5565
   ClientLeft      =   60
   ClientTop       =   450
   ClientWidth     =   6735
   LinkTopic       =   "Form1"
   ScaleHeight     =   5565
   ScaleWidth      =   6735
   StartUpPosition =   3  '窗口缺省
   Begin VB.CommandButton Command1 
      Caption         =   "插值计算"
      Height          =   375
      Left            =   5760
      TabIndex        =   0
      Top             =   120
      Width           =   975
   End
End
Attribute VB_Name = "Form1"
Attribute VB_GlobalNameSpace = False
Attribute VB_Creatable = False
Attribute VB_PredeclaredId = True
Attribute VB_Exposed = False
'过程POLINT
'N:节点个数,输入参数,整型变量
'XA():N个元素的一维实型数组,输入参数,存放插值节点
'YA():N个元素的一维实型数组,输入参数,存放函数值
'X:插值自变量,输入参数,实型参量
'Y:所求值,输出参数,实型变量
'DY:误差估计,输出参数,实型变量
Sub POLINT(XA(), YA(), N, X, Y, DY)
Dim C(10), D(10)
NS = 1
DIF = Abs(X - XA(1))
For I = 1 To N
DIFT = Abs(X - XA(I))
If DIFT < DIF Then
NS = I
DIF = DIFT
End If
C(I) = YA(I)
D(I) = YA(I)
Next I
Y = YA(NS)
NS = NS - 1
For M = 1 To N - 1
For I = 1 To N - M
HO = XA(I) - X
HP = XA(I + M) - X
W = C(I + 1) - D(I)
DEN = HO - HP
If DEN = 0# Then
Print "PAUSE"
Exit Sub
End If
DEN = W / DEN
D(I) = HP * DEN
C(I) = HO * DEN
Next I
If 2 * NS < N - M Then
DY = C(NS + 1)
Else
DY = D(NS)
NS = NS - 1
End If
Y = Y + DY
Next M
End Sub
'我们所取的例子是f(x)=sinx和f(x)=exp(x).sinx的区间是(0,pi),exp(x)的区间是(0,1.0).
'验证程序中的N是任意给定的(N<10)当N增长时,我们可以观察到结果在不断的改进.验证程序D2R1如下:
Private Sub Command1_Click()
'PROGRAM D2R1
'Drive for routine POLINT
N = 9 '插值节点个数
PI = 3.1415926
Dim XA(10), YA(10)
Print Tab(5); "Generation of interpolation tables"
Print Tab(5); "...sin(x)  0<x<PI"
Print Tab(5); "...EXP(X)  0<X<1"
Print Tab(5); "How many entries go in these tables? (Note:N<10)"
Print Tab(5); N
For NFUNC = 1 To 2
If NFUNC = 1 Then
Print Tab(5); "sine function from 0 to PI"
For I = 1 To N
XA(I) = I * PI / N
YA(I) = Sin(XA(I))
Next I
ElseIf NFUNC = 2 Then
Print Tab(5); "exponential function from 0 to 1"
For I = 1 To N
XA(I) = I * 1# / N
YA(I) = Exp(XA(I))
Next I
Else
Exit Sub
End If
Print Tab(5); "X              F(X)           interpolated   error"
For I = 1 To 10
If NFUNC = 1 Then
X = (-0.05 + I / 10#) * PI
F = Sin(X)
ElseIf NFUNC = 2 Then
X = (-0.05 + I / 10#)
F = Exp(X)
End If
Call POLINT(XA(), YA(), N, X, Y, DY)
Print Tab(5); Format$(X, "0.##0000");
Print Tab(20); Format$(F, "0.####00");
Print Tab(35); Format$(Y, "0.####00");
Print Tab(50); Format$(DY, ".##00E+00")
Next I
Next NFUNC
End Sub