?? interpmodule.bas
字號:
s(j) = s(j) / (2# + (1# - alpha) * dy(j - 1))
h0 = h1
Next j
' 一階導數值
dy(n) = (3# * (y(n) - y(n - 1)) / h1 + ddy(n) * h1 / 2# - s(n - 1)) / (2# + dy(n - 1))
For j = n - 1 To 1 Step -1
dy(j) = dy(j) * dy(j + 1) + s(j)
Next j
For j = 1 To n - 1
s(j) = x(j + 1) - x(j)
Next j
' 二階導數值
For j = 1 To n - 1
h1 = s(j) * s(j)
ddy(j) = 6# * (y(j + 1) - y(j)) / h1 - 2# * (2# * dy(j) + dy(j + 1)) / s(j)
Next j
h1 = s(n - 1) * s(n - 1)
ddy(n) = 6# * (y(n - 1) - y(n)) / h1 + 2# * (2# * dy(n) + dy(n - 1)) / s(n - 1)
g = 0#
For i = 1 To n - 1
h1 = 0.5 * s(i) * (y(i) + y(i + 1))
h1 = h1 - s(i) * s(i) * s(i) * (ddy(i) + ddy(i + 1)) / 24#
g = g + h1
Next i
' 插值
For j = 1 To m
If (t(j) >= x(n)) Then
i = n - 2
Else
i = 0
While (t(j) > x(i + 1))
i = i + 1
Wend
End If
h1 = (x(i + 1) - t(j)) / s(i)
h0 = h1 * h1
z(j) = (3# * h0 - 2# * h0 * h1) * y(i)
z(j) = z(j) + s(i) * (h0 - h0 * h1) * dy(i)
dz(j) = 6# * (h0 - h1) * y(i) / s(i)
dz(j) = dz(j) + (3# * h0 - 2# * h1) * dy(i)
ddz(j) = (6# - 12# * h1) * y(i) / (s(i) * s(i))
ddz(j) = ddz(j) + (2# - 6# * h1) * dy(i) / s(i)
h1 = (t(j) - x(i)) / s(i)
h0 = h1 * h1
z(j) = z(j) + (3# * h0 - 2# * h0 * h1) * y(i + 1)
z(j) = z(j) - s(i) * (h0 - h0 * h1) * dy(i + 1)
dz(j) = dz(j) - 6# * (h0 - h1) * y(i + 1) / s(i)
dz(j) = dz(j) + (3# * h0 - 2# * h1) * dy(i + 1)
ddz(j) = ddz(j) + (6# - 12# * h1) * y(i + 1) / (s(i) * s(i))
ddz(j) = ddz(j) - (2# - 6# * h1) * dy(i + 1) / s(i)
Next j
' 返回積分值
INSpline2 = g
End Function
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
' 模塊名：InterpModule.bas
' 函數名：INSpline3
' 功能： 實現第三種邊界條件的三次樣條函數插值、微商與積分
' 參數： n - Integer型變量，給定結點的點數
' x - Double型一維數組，長度為n，存放給定的n個結點的值x(i)，要求x(1)<x(2)<...<x(n)
' y - Double型一維數組，長度為n，存放給定的n個等距結點的函數值y(i)，y(i) = f(x(i)), i=1,2,...,n
' 要求y(1)=y(n)
' dy - Double型一維數組，長度為n，返回時，存放n個給定點處的一階導數值y'(i)，i=1,2,…n
' ddy - Double型一維數組，長度為n，返回時，存放n個給定點處的二階導數值y''(i)，i=1,2,…n
' m - Integer型變量，指定插值點的個數
' t - Double型一維數組，長度為m，存放m個指定的插值點的值。要求x(1)<t(j)<x(n), j=1,2,…,m-1
' z - Double型一維數組，長度為m，存放m個指定的插值點處的函數值
' dz - Double型一維數組，長度為m，存放m個指定的插值點處的一階導數值
' ddz - Double型一維數組，長度為m，存放m個指定的插值點處的二階導數值
' 返回值：Double型，指定函數的x(1)到x(n)的定積分值
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function INSpline3(n As Integer, x() As Double, y() As Double, dy() As Double, ddy() As Double, m As Integer, t() As Double, z() As Double, dz() As Double, ddz() As Double) As Double
Dim i As Integer, j As Integer
Dim y0 As Double, y1 As Double, h0 As Double, h1 As Double, alpha As Double, beta As Double, u As Double, g As Double
ReDim s(n) As Double
' 初值
h0 = x(n) - x(n - 1)
y0 = y(n) - y(n - 1)
dy(1) = 0#
ddy(1) = 0#
ddy(n) = 0#
s(1) = 1#
s(n) = 1#
For j = 2 To n
h1 = h0
y1 = y0
h0 = x(j) - x(j - 1)
y0 = y(j) - y(j - 1)
alpha = h1 / (h1 + h0)
beta = 3# * ((1# - alpha) * y1 / h1 + alpha * y0 / h0)
If (j < n) Then
u = 2# + (1# - alpha) * dy(j - 1)
dy(j) = -alpha / u
s(j) = (alpha - 1#) * s(j - 1) / u
ddy(j) = (beta - (1# - alpha) * ddy(j - 1)) / u
End If
Next j
' 二階導數值
For j = n - 1 To 2 Step -1
s(j) = dy(j) * s(j + 1) + s(j)
ddy(j) = dy(j) * ddy(j + 1) + ddy(j)
Next j
' 一階導數值
dy(n - 1) = (beta - alpha * ddy(2) - (1# - alpha) * ddy(n - 1)) / (alpha * s(2) + (1# - alpha) * s(n - 1) + 2#)
For j = 3 To n
dy(j - 2) = s(j - 1) * dy(n - 1) + ddy(j - 1)
Next j
dy(n) = dy(1)
For j = 1 To n - 1
s(j) = x(j + 1) - x(j)
Next j
For j = 1 To n - 1
h1 = s(j) * s(j)
ddy(j) = 6# * (y(j + 1) - y(j)) / h1 - 2# * (2# * dy(j) + dy(j + 1)) / s(j)
Next j
h1 = s(n - 1) * s(n - 1)
ddy(n) = 6# * (y(n - 1) - y(n)) / h1 + 2# * (2# * dy(n) + dy(n - 1)) / s(n - 1)
g = 0#
For i = 1 To n - 1
h1 = 0.5 * s(i) * (y(i) + y(i + 1))
h1 = h1 - s(i) * s(i) * s(i) * (ddy(i) + ddy(i + 1)) / 24#
g = g + h1
Next i
' 插值
For j = 1 To m
h0 = t(j)
While (h0 >= x(n))
h0 = h0 - (x(n) - x(1))
Wend
While (h0 < x(1))
h0 = h0 + (x(n) - x(1))
Wend
i = 1
While (h0 > x(i + 1))
i = i + 1
Wend
u = h0
h1 = (x(i + 1) - u) / s(i)
h0 = h1 * h1
z(j) = (3# * h0 - 2# * h0 * h1) * y(i)
z(j) = z(j) + s(i) * (h0 - h0 * h1) * dy(i)
dz(j) = 6# * (h0 - h1) * y(i) / s(i)
dz(j) = dz(j) + (3# * h0 - 2# * h1) * dy(i)
ddz(j) = (6# - 12# * h1) * y(i) / (s(i) * s(i))
ddz(j) = ddz(j) + (2# - 6# * h1) * dy(i) / s(i)
h1 = (u - x(i)) / s(i)
h0 = h1 * h1
z(j) = z(j) + (3# * h0 - 2# * h0 * h1) * y(i + 1)
z(j) = z(j) - s(i) * (h0 - h0 * h1) * dy(i + 1)
dz(j) = dz(j) - 6# * (h0 - h1) * y(i + 1) / s(i)
dz(j) = dz(j) + (3# * h0 - 2# * h1) * dy(i + 1)
ddz(j) = ddz(j) + (6# - 12# * h1) * y(i + 1) / (s(i) * s(i))
ddz(j) = ddz(j) - (2# - 6# * h1) * dy(i + 1) / s(i)
Next j
' 返回積分值
INSpline3 = g
End Function
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
' 模塊名：InterpModule.bas
' 函數名：INTqip
' 功能： 實現二元三點插值
' 參數： n - Integer型變量，x方向上給定結點的點數
' m - Integer型變量，y方向上給定結點的點數
' x - Double型一維數組，長度為n，存放給定n x m 個結點x方向上的n個值x(i)，要求x(1)<x(2)<...<x(n)
' y - Double型一維數組，長度為m，存放給定n x m 個結點y方向上的m個值y(i)，要求y(1)<y(2)<...<y(m)
' z - Double型n x m二維數組，存放給定的n x m個結點的函數值z(i,j)，z(i,j) = f(x(i), y(j)), i=1,2,...,n, j=1,2,...,m
' u - Double型變量，存放插值點x坐標
' v - Double型變量，存放插值點y坐標
' 返回值：Double型，指定函數值f(u, v)
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function INTqip(n As Integer, m As Integer, x() As Double, y() As Double, z() As Double, u As Double, v As Double) As Double
Dim nn As Integer, mm As Integer, ip As Integer, iq As Integer, i As Integer, j As Integer, k As Integer, l As Integer
Dim b(3) As Double, h As Double, w As Double
nn = 3
If (n <= 3) Then
ip = 0
nn = n
Else
If (u <= x(2)) Then
ip = 0
Else
If (u >= x(n - 1)) Then
ip = n - 3
Else
i = 1
j = n
While (((i - j) <> 1) And ((i - j) <> -1))
l = (i + j) / 2
If (u < x(l)) Then
j = l
Else
i = l
End If
Wend
If (Abs(u - x(i)) < Abs(u - x(j))) Then
ip = i - 2
Else
ip = i - 1
End If
End If
End If
End If
mm = 3
If (m <= 3) Then
iq = 0
mm = m
Else
If (v <= y(2)) Then
iq = 0
Else
If (v >= y(m - 1)) Then
iq = m - 3
Else
i = 1
j = m
While (((i - j) <> 1) And ((i - j) <> -1))
l = (i + j) / 2
If (v < y(l)) Then
j = l
Else
i = l
End If
Wend
If (Abs(v - y(i)) < Abs(v - y(j))) Then
iq = i - 2
Else
iq = i - 1
End If
End If
End If
End If
For i = 1 To nn
b(i) = 0#
For j = 1 To mm
h = z(ip + i, iq + j)
For k = 1 To mm
If (k <> j) Then
' 二元拉格朗日公式
h = h * (v - y(iq + k)) / (y(iq + j) - y(iq + k))
End If
Next k
b(i) = b(i) + h
Next j
Next i
w = 0#
For i = 1 To nn
h = b(i)
For j = 1 To nn
If (j <> i) Then
' 二元拉格朗日公式
h = h * (u - x(ip + j)) / (x(ip + i) - x(ip + j))
End If
Next j
w = w + h
Next i
' 返回函數值
INTqip = w
End Function
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
' 模塊名：InterpModule.bas
' 函數名：INLagrn2
' 功能： 實現二元全區間插值
' 參數： n - Integer型變量，x方向上給定結點的點數
' m - Integer型變量，y方向上給定結點的點數
' x - Double型一維數組，長度為n，存放給定n x m 個結點x方向上的n個值x(i)，要求x(1)<x(2)<...<x(n)
' y - Double型一維數組，長度為m，存放給定n x m 個結點y方向上的m個值y(i)，要求y(1)<y(2)<...<y(m)
' z - Double型n x m二維數組，存放給定的n x m個結點的函數值z(i,j)，z(i,j) = f(x(i), y(j)), i=1,2,...,n, j=1,2,...,m
' u - Double型變量，存放插值點x坐標
' v - Double型變量，存放插值點y坐標
' 返回值：Double型，指定函數值f(u, v)
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function INLagrn2(n As Integer, m As Integer, x() As Double, y() As Double, z() As Double, u As Double, v As Double) As Double
Dim ip As Integer, ipp As Integer, iq As Integer, iqq As Integer, i As Integer, j As Integer, k As Integer, l As Integer
Dim b(10) As Double, h As Double, w As Double
If (u <= x(1)) Then
ip = 1
ipp = 4
Else
If (u >= x(n)) Then
ip = n - 3
ipp = n
Else
i = 1
j = n
While (((i - j) <> 1) And ((i - j) <> -1))
l = (i + j) / 2
If (u < x(l)) Then
j = l
Else
i = l
End If
Wend
ip = i - 3
ipp = i + 4
End If
End If
If (ip < 1) Then
ip = 1
End If
If (ipp > n) Then
ipp = n
End If
If (v <= y(1)) Then
iq = 1
iqq = 4
Else
If (v >= y(m)) Then
iq = m - 3
iqq = m
Else
i = 1
j = m
While (((i - j) <> 1) And ((i - j) <> -1))
l = (i + j) / 2
If (v < y(l)) Then
j = l
Else
i = l
End If
Wend
iq = i - 3
iqq = i + 4
End If
End If
If (iq < 1) Then
iq = 1
End If
If (iqq > m) Then
iqq = m
End If
For i = ip To ipp
b(i - ip + 1) = 0#
For j = iq To iqq
h = z(i, j)
For k = iq To iqq
' 二元拉格朗日公式
If (k <> j) Then h = h * (v - y(k)) / (y(j) - y(k))
Next k
b(i - ip + 1) = b(i - ip + 1) + h
Next j
Next i
w = 0#
For i = ip To ipp
h = b(i - ip + 1)
For j = ip To ipp
' 二元拉格朗日公式
If (j <> i) Then h = h * (u - x(j)) / (x(i) - x(j))
Next j
w = w + h
Next i
' 返回函數值
INLagrn2 = w
End Function
?? 快捷鍵說明
復制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -