d8r4.frm

常用的数值算法的VB程序

VERSION 5.00
Begin VB.Form Form1 
   Caption         =   "Form1"
   ClientHeight    =   5385
   ClientLeft      =   60
   ClientTop       =   345
   ClientWidth     =   4680
   LinkTopic       =   "Form1"
   ScaleHeight     =   5385
   ScaleWidth      =   4680
   StartUpPosition =   3  'Windows Default
   Begin VB.CommandButton Command1 
      Caption         =   "Command1"
      Height          =   375
      Left            =   3240
      TabIndex        =   0
      Top             =   4800
      Width           =   1095
   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 D8R4
    'Driver for routine TQLI
    NP = 3
    TINY = 0.000001
    Dim A(3, 3), C(3, 3), D(3), E(3), F(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#
    For I = 1 To NP
        For J = 1 To NP
            C(I, J) = A(I, J)
        Next J
    Next I
    Call TRED2(C(), NP, D(), E())
    Call TQLI(D(), E(), NP, C())
    Print
    Print Tab(5); " Eigenvectors for a real symmetric matrix"
    For I = 1 To NP
        For J = 1 To NP
            F(J) = 0#
            For K = 1 To NP
                F(J) = F(J) + A(J, K) * C(K, I)
            Next K
        Next J
        Print
        Print Tab(5); " Eigenvalue "; I; "= ";
        Print Format$(D(I), "0.000000")
        Print
        Print Tab(5); " Vector     Mtrx*vect.    Ratio"
        For J = 1 To NP
            If Abs(C(J, I)) < TINY Then
              Print Tab(5); Format$(C(J, I), "#.000000"),
              Print Format$(F(J), "#.000000"),
              Print "div. by 0"
            Else
              Print Tab(5); Format$(C(J, I), "#.000000"),
              Print Format$(F(J), "#.000000"),
              Print Format$(F(J) / C(J, I), "0.000000")
            End If
        Next J
    Next I
End Sub
Sub TQLI(D(), E(), N, Z())
    If N > 1 Then
        For I = 2 To N
            E(I - 1) = E(I)
        Next I
        E(N) = 0#
        For L = 1 To N
            ITER = 0
1           For M = L To N - 1
                DD = Abs(D(M)) + Abs(D(M + 1))
                If Abs(E(M)) + DD = DD Then GoTo 2
            Next M
            M = N
2           If M <> L Then
                If ITER = 30 Then Print " too many iterations "
                ITER = ITER + 1
                G = (D(L + 1) - D(L)) / (2# * E(L))
                R = Sqr(G ^ 2 + 1#)
                G = D(M) - D(L) + E(L) / (G + Abs(R) * Sgn(G))
                S = 1#
                C = 1#
                P = 0#
                For I = M - 1 To L Step -1
                    F = S * E(I)
                    B = C * E(I)
                    If Abs(F) >= Abs(G) Then
                        C = G / F
                        R = Sqr(C ^ 2 + 1#)
                        E(I + 1) = F * R
                        S = 1# / R
                        C = C * S
                    Else
                        S = F / G
                        R = Sqr(S ^ 2 + 1#)
                        E(I + 1) = G * R
                        C = 1# / R
                        S = S * C
                    End If
                    G = D(I + 1) - P
                    R = (D(I) - G) * S + 2# * C * B
                    P = S * R
                    D(I + 1) = G + P
                    G = C * R - B
                    'Omit lines from here ...
                    For K = 1 To N
                        F = Z(K, I + 1)
                        Z(K, I + 1) = S * Z(K, I) + C * F
                        Z(K, I) = C * Z(K, I) - S * F
                    Next K
                   'to here when finding only eigenvalues.
                Next I
                D(L) = D(L) - P
                E(L) = G
                E(M) = 0#
                GoTo 1
            End If
       Next L
    End If
End Sub
   Sub TRED2(A(), N, D(), E())
      If N > 1 Then
          For I = N To 2 Step -1
              L = I - 1
              H = 0#
              SCALE1 = 0#
              If L > 1 Then
                  For K = 1 To L
                      SCALE1 = SCALE1 + Abs(A(I, K))
                  Next K
                  If SCALE1 = 0# Then
                      E(I) = A(I, L)
                  Else
                      For K = 1 To L
                          A(I, K) = A(I, K) / SCALE1
                          H = H + A(I, K) ^ 2
                      Next K
                      F = A(I, L)
                      G = -Sqr(H) * Sgn(F)
                      E(I) = SCALE1 * G
                      H = H - F * G
                      A(I, L) = F - G
                      F = 0#
                      For J = 1 To L
                          A(J, I) = A(I, J) / H
                          G = 0#
                          For K = 1 To J
                              G = G + A(J, K) * A(I, K)
                          Next K
                          If L > J Then
                              For K = J + 1 To L
                                  G = G + A(K, J) * A(I, K)
                              Next K
                          End If
                          E(J) = G / H
                          F = F + E(J) * A(I, J)
                      Next J
                      HH = F / (H + H)
                      For J = 1 To L
                          F = A(I, J)
                          G = E(J) - HH * F
                          E(J) = G
                          For K = 1 To J
                              A(J, K) = A(J, K) - F * E(K) - G * A(I, K)
                          Next K
                      Next J
                  End If
              Else
                  E(I) = A(I, L)
              End If
              D(I) = H
          Next I
      End If
      'Omit following line if finding only eigenvalues.
      D(1) = 0#
      E(1) = 0#
      For I = 1 To N
      'Delete lines from here ...
          L = I - 1
          If D(I) <> 0# Then
              For J = 1 To L
                  G = 0#
                  For K = 1 To L
                      G = G + A(I, K) * A(K, J)
                  Next K
                  For K = 1 To L
                      A(K, J) = A(K, J) - G * A(K, I)
                  Next K
              Next J
          End If
       '... to here when finding only eibenvalues.
          D(I) = A(I, I)
      'Also delete lines from here ...
          A(I, I) = 1#
          If L >= 1 Then
              For J = 1 To L
                  A(I, J) = 0#
                  A(J, I) = 0#
              Next J
          End If
      '... to here when finding only eigenvalues.
      Next I
   End Sub