d1r2.frm

矩阵特征值的求解过程之二

VERSION 5.00
Begin VB.Form Form1 
   Caption         =   "Form1"
   ClientHeight    =   3810
   ClientLeft      =   60
   ClientTop       =   345
   ClientWidth     =   4680
   LinkTopic       =   "Form1"
   ScaleHeight     =   3810
   ScaleWidth      =   4680
   StartUpPosition =   3  'Windows Default
   Begin VB.CommandButton Command1 
      Caption         =   "Command1"
      Height          =   375
      Left            =   3000
      TabIndex        =   0
      Top             =   3240
      Width           =   1335
   End
End
Attribute VB_Name = "Form1"
Attribute VB_GlobalNameSpace = False
Attribute VB_Creatable = False
Attribute VB_PredeclaredId = True
Attribute VB_Exposed = False
Private Sub Command1_Click()
    'program D1R2
    '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 I = 1 To N
        For J = 1 To N
            A1(I, J) = A(I, J)
        Next J
    Next I
    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 J = 1 To N
            B(L) = B(L) + A(L, J) * X(J)
        Next J
    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
Sub LUBKSB(A(), N, INDX(), B())
    II = 0
    For I = 1 To N
        LL = INDX(I)
        Sum = B(LL)
        B(LL) = B(I)
        If II <> 0 Then
            For J = II To I - 1
                Sum = Sum - A(I, J) * B(J)
            Next J
        ElseIf Sum <> 0# Then
            II = I
        End If
        B(I) = Sum
    Next I
    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(J)
            Next J
        End If
        B(I) = Sum / A(I, I)
    Next I
End Sub
Sub LUDCMP(A(), N, INDX(), D)
    NMAX = 100
    TINY = 1E-20
    Dim VV(100)
    D = 1#
    For I = 1 To N
        AAMAX = 0#
        For J = 1 To N
            If Abs(A(I, J)) > AAMAX Then AAMAX = Abs(A(I, J))
        Next J
        If AAMAX = 0# Then Print "Singular matrix."
        VV(I) = 1# / AAMAX
    Next I
    For J = 1 To N
        If J > 1 Then
            For I = 1 To J - 1
                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
            Next I
        End If
        AAMAX = 0#
        For I = J To N
            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
            DUM = VV(I) * Abs(Sum)
            If DUM >= AAMAX Then
                IMAX = I
                AAMAX = DUM
            End If
        Next I
        If J <> IMAX Then
            For K = 1 To N
                DUM = A(IMAX, K)
                A(IMAX, K) = A(J, K)
                A(J, K) = DUM
            Next K
            D = -D
            VV(IMAX) = VV(J)
        End If
        INDX(J) = IMAX
        If J <> N Then
            If A(J, J) = 0# Then A(J, J) = TINY
            DUM = 1# / A(J, J)
            For I = J + 1 To N
                A(I, J) = A(I, J) * DUM
            Next I
        End If
    Next J
    If A(N, N) = 0# Then A(N, N) = TINY
End Sub