lufjf.frm

解线性方程组的一种行主元分解法

VERSION 5.00
Begin VB.Form Form1 
   Caption         =   "Form1"
   ClientHeight    =   3195
   ClientLeft      =   60
   ClientTop       =   345
   ClientWidth     =   4680
   LinkTopic       =   "Form1"
   ScaleHeight     =   3195
   ScaleWidth      =   4680
   StartUpPosition =   3  '窗口缺省
   Begin VB.CommandButton Command1 
      Caption         =   "计算"
      Height          =   375
      Left            =   3240
      TabIndex        =   0
      Top             =   1320
      Width           =   1215
   End
End
Attribute VB_Name = "Form1"
Attribute VB_GlobalNameSpace = False
Attribute VB_Creatable = False
Attribute VB_PredeclaredId = True
Attribute VB_Exposed = False
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
''过程LUDCMP是将系数矩阵A分解为上三角矩阵U和下三角矩陈L
''N：A的阶数，整型变量，为输入参数
''A：实型数组，在LUDCMP中，输入时按列存放实方陈，输出时，对角线以下部分存放单位下三角矩阵L，对角线及以上部分存放上三角矩阵U。
''   在LUBKSB中，将LUDCMP中输出的结果A作为输入。
''INDX（）：整型数组，在子过程LUDCMP中为输出参数，用于记录置换矩阵，在子过程LUBKSB中为输入参数，输入子过程LUDCMP的输出结果。
''D：1或-1，为输出参数，依赖于行交换次数为偶（+1）还是奇（-1）。
''B（）：实型数组，输入、输出参数，输入实向量了，输出时，方程组的解X存储在数组B（）中。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub LUDCMP(A(), N, INDX(), D)
NMAX = 100
TINY = 1E-20
Dim VV(100)
D = 1#
For j = 1 To N
AAMAX = 0#
For i = 1 To N
If Abs(A(i, j)) > AAMAX Then AAMAX = Abs(A(i, j))
Next i
If AAMAX = 0# Then Print "singular matrix."
VV(j) = 1# / AAMAX
Next j
For i = 1 To N
If i > 1 Then
For j = 1 To i - 1
Sum = A(i, j)
If j > 1 Then
For K = 1 To j - 1
Sum = Sum - A(i, K) * A(K, j)
Next K
A(i, j) = Sum
End If
Next j
End If
AAMAX = 0#
For j = i To N
Sum = A(i, j)
If i > 1 Then
For K = 1 To i - 1
Sum = Sum - A(i, K) * A(K, j)
Next K
A(i, j) = Sum
End If
DUM = VV(j) * Abs(Sum)
If DUM >= AAMAX Then
IMAX = j
AAMAX = DUM
End If
Next j
If i <> IMAX Then
For K = 1 To N
DUM = A(IMAX, K)
A(IMAX, K) = A(i, K)
A(i, K) = DUM
Next K
D = -D
VV(IMAX) = VV(i)
End If
INDX(i) = IMAX
If i <> N Then
If A(i, i) = 0# Then A(i, i) = TINY
DUM = 1# / A(i, i)
For j = i + 1 To N
A(i, j) = A(i, j) * DUM
Next j
End If
Next i
If A(N, N) = 0# Then A(N, N) = TINY
End Sub
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
''子过程LUBKSB利用上面LUDCMP的分解结果求得线性方程组AX=B的解。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub LUBKSB(A(), N, INDX(), B())
II = 0
For j = 1 To N
LL = INDX(j)
Sum = B(LL)
B(LL) = B(j)
If II <> 0 Then
For i = II To j - 1
Sum = Sum - A(i, j) * B(i)
Next i
ElseIf Sum <> 0# Then
II = j
End If
B(i) = Sum
Next j
For i = N To 1 Step -1
Sum = B(i)
If i < N Then
For j = i + 1 To N
Sum = Sum - A(i, j) * B(i)
Next j
End If
B(i) = Sum / A(i, i)
Next i
End Sub

Private Sub Command1_Click()
'PROGRAM DIR2
'DRIVER PROGRAM FOR ROUTINE LUBKSB,LUDCMP
N = 3  ''系数矩阵维数
Dim A(3, 3), B(3), A1(3, 3), INDX(3), X(3)
'输入已知的方程组的系数矩阵
A(1, 1) = 1#: A(1, 2) = 2#: A(1, 3) = 3#
A(2, 1) = 2#: A(2, 2) = 2#: A(2, 3) = 3#
A(3, 1) = 3#: A(3, 2) = 3#: A(3, 3) = 3#
''输入已知方程组的右端向量
B(1) = 1#
B(2) = 2#
B(3) = 3#
Print
Print Tab(5); "已知的方程组的右端向量："
Print Tab(14); Format$(B(1), "##.##")
Print Tab(14); Format$(B(2), "##.##")
Print Tab(14); Format$(B(3), "##.##")
For j = 1 To N
For i = 1 To N
A1(i, j) = A(i, j)
Next i
Next j
Call LUDCMP(A1(), N, INDX(), P)
For K = 1 To N
For L = 1 To N
X(L) = B(L)
Next L
Next K
Call LUBKSB(A1(), N, INDX(), X())
Print
Print Tab(5); "计算出的方程组的解："
Print Tab(14); Format$(X(1), "#.####E+00")
Print Tab(14); Format$(X(2), "#.####E+00")
Print Tab(14); Format$(X(3), "#.####E+00")
''将计算出的B乘以系数矩阵，以验证计算结果正确
For L = 1 To N
B(L) = 0#
For i = 1 To N
B(L) = B(L) + A(L, i) * X(i)
Next i
Next L
Print
Print Tab(5); "计算出的解乘以系数矩阵验证的结果："
Print Tab(14); Format$(B(1), "##.##")
Print Tab(14); Format$(B(2), "##.##")
Print Tab(14); Format$(B(3), "##.##")
End Sub