d8r2.frm

常用的数值算法的VB程序

VERSION 5.00
Begin VB.Form Form1 
   Caption         =   "Form1"
   ClientHeight    =   6180
   ClientLeft      =   1860
   ClientTop       =   945
   ClientWidth     =   4680
   LinkTopic       =   "Form1"
   ScaleHeight     =   6180
   ScaleWidth      =   4680
   Begin VB.CommandButton Command1 
      Caption         =   "Command1"
      Height          =   375
      Left            =   3240
      TabIndex        =   0
      Top             =   5520
      Width           =   1215
   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 D8R2
    'Driver for routine EIGSRT
    NP = 3
    Dim D(3), V(3, 3), A(3, 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#
    Call JACOBI(A(), NP, D(), V(), NROT)
    Print Tab(5); "Unsorted Eigenvectors:"
    For I = 1 To NP
        Print Tab(5)
        Print Tab(5); " Eigenvalue"; I; " = ";
        Print Format$(D(I), ".#####0")
        Print Tab(5); " Eigenvector:"
        For J = 1 To NP
            Print Format$(V(J, I), "#.######"),
        Next J
    Next I
    Print Tab(5)
    Print Tab(5); "********** sorting **********"
    Print
    Call EIGSRT(D(), V(), NP)
    Print Tab(5); " Sorted Eigenvectors:"
    For I = 1 To NP
        Print Tab(5)
        Print Tab(5); " Eigenvalue"; I; " = ";
        Print Format$(D(I), ".#####0")
        Print Tab(5); " Eigenvector:"
        For J = 1 To NP
            Print Format$(V(J, I), ".#####0"),
        Next J
    Next I
End Sub
Sub EIGSRT(D(), V(), N)
    For I = 1 To N - 1
        K = I
        P = D(I)
        For J = I + 1 To N
            If D(J) >= P Then
                K = J
                P = D(J)
            End If
        Next J
        If K <> I Then
            D(K) = D(I)
            D(I) = P
            For J = 1 To N
                P = V(J, I)
                V(J, I) = V(J, K)
                V(J, K) = P
            Next J
        End If
    Next I
End Sub
Sub JACOBI(A(), N, D(), V(), NROT)
    Dim B(100), Z(100)
    For IP = 1 To N
      For IQ = 1 To N
        V(IP, IQ) = 0#
      Next IQ
      V(IP, IP) = 1#
    Next IP
    For IP = 1 To N
      B(IP) = A(IP, IP)
      D(IP) = B(IP)
      Z(IP) = 0#
    Next IP
    NROT = 0
    For I = 1 To 50
      SM = 0#
      For IP = 1 To N - 1
        For IQ = IP + 1 To N
          SM = SM + Abs(A(IP, IQ))
        Next IQ
      Next IP
      If SM = 0# Then Exit Sub
        If I < 4 Then
          TRESH = 0.2 * SM / N ^ 2
        Else
          TRESH = 0#
        End If
        For IP = 1 To N - 1
          For IQ = IP + 1 To N
            G = 100# * Abs(A(IP, IQ))
            SSS = Abs(D(IP)) + G
            DDD = Abs(D(IQ)) + G
            If (I > 4) And (SSS = Abs(D(IP))) And (DDD = Abs(D(IQ))) Then
               A(IP, IQ) = 0#
            ElseIf Abs(A(IP, IQ)) > TRESH Then
               H = D(IQ) - D(IP)
               If Abs(H) + G = Abs(H) Then
               T = A(IP, IQ) / H
            Else
              THETA = 0.5 * H / A(IP, IQ)
              T = 1# / (Abs(THETA) + Sqr(1# + THETA ^ 2))
              If THETA < 0# Then T = -T
            End If
              C = 1# / Sqr(1# + T ^ 2)
              S = T * C
              TAU = S / (1# + C)
              H = T * A(IP, IQ)
              Z(IP) = Z(IP) - H
              Z(IQ) = Z(IQ) + H
              D(IP) = D(IP) - H
              D(IQ) = D(IQ) + H
              A(IP, IQ) = 0#
              For J = 1 To IP - 1
                G = A(J, IP)
                H = A(J, IQ)
                A(J, IP) = G - S * (H + G * TAU)
                A(J, IQ) = H + S * (G - H * TAU)
              Next J
              For J = IP + 1 To IQ - 1
                G = A(IP, J)
                H = A(J, IQ)
                A(IP, J) = G - S * (H + G * TAU)
                A(J, IQ) = H + S * (G - H * TAU)
              Next J
              For J = IQ + 1 To N
                G = A(IP, J)
                H = A(IQ, J)
                A(IP, J) = G - S * (H + G * TAU)
                A(IQ, J) = H + S * (G - H * TAU)
              Next J
              For J = 1 To N
                G = V(J, IP)
                H = V(J, IQ)
                V(J, IP) = G - S * (H + G * TAU)
                V(J, IQ) = H + S * (G - H * TAU)
              Next J
              NROT = NROT + 1
            End If
          Next IQ
        Next IP
        For IP = 1 To N
              B(IP) = B(IP) + Z(IP)
              D(IP) = B(IP)
              Z(IP) = 0#
        Next IP
      Next I
      Print " 50 iterations should never happen"
End Sub