d9r6.frm

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

VERSION 5.00
Begin VB.Form Form1 
   Caption         =   "Form1"
   ClientHeight    =   3405
   ClientLeft      =   60
   ClientTop       =   345
   ClientWidth     =   5280
   LinkTopic       =   "Form1"
   ScaleHeight     =   3405
   ScaleWidth      =   5280
   StartUpPosition =   3  'Windows Default
   Begin VB.CommandButton Command1 
      Caption         =   "Command1"
      Height          =   375
      Left            =   3600
      TabIndex        =   0
      Top             =   2520
      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 D9R6
    'Driver for routine POWELL
    NP = 3
    NDIM = NP
    FTOL = 0.000001
    Dim P(3), XI(3, 3)
    For I = 1 To NP
        For J = 1 To NP
            XI(I, J) = 0#
        Next J
    Next I
    XI(1, 1) = 1#
    XI(2, 2) = 1#
    XI(3, 3) = 1#
    P(1) = 1.5: P(2) = 1.5: P(3) = 2.5
    Call POWELL(P(), XI(), NDIM, NP, FTOL, ITER, FRET)
    Print
    Print Tab(5); "Iterations:  "; Format$(ITER, "##")
    Print
    Print Tab(5); "Minimum found at: "
    Print
    For I = 1 To NP
        Print Tab(5 + (I - 1) * 12); Format$(P(I), "#.#####0");
    Next I
    Print Tab(5)
    Print Tab(5); "Minimum function value =  "; Format$(FRET, ".#####0")
    Print
    Print Tab(5); "True minimum of function is at: 1.0   2.0   3.0"
End Sub
Function FUNC(X)
    FUNC = F1DIM(X)
End Function
Function FUNC2(X(), N)
    AAA = (X(1) - 1#) ^ 2 + (X(2) - 2#) ^ 2 + (X(3) - 3#) ^ 2
    FUNC2 = 0.5 - BESSJ0(AAA)
End Function
Sub POWELL(P(), XI(), N, NP, FTOL, ITER, FRET)
    ITMAX = 200
    Dim PT(20), PTT(20), XIT(20)
    FRET = FUNC2(P(), N)
    For J = 1 To N
        PT(J) = P(J)
    Next J
    ITER = 0
    Do
      Do
        Do
          ITER = ITER + 1
          FP = FRET
          IBIG = 0
          DEL = 0#
          For I = 1 To N
              For J = 1 To N
                  XIT(J) = XI(J, I)
              Next J
              FPTT = FRET
              Call LINMIN(P(), XIT(), N, FRET)
              If Abs(FPTT - FRET) > DEL Then
                  DEL = Abs(FPTT - FRET)
                  IBIG = I
              End If
          Next I
          If 2# * Abs(FP - FRET) <= FTOL * (Abs(FP) + Abs(FRET)) Then
              Erase XIT, PTT, PT
              Exit Sub
          End If
          If ITER = ITMAX Then
              Print " POWELL exceeding maximum iterations"
              Exit Sub
          End If
          For J = 1 To N
              PTT(J) = 2# * P(J) - PT(J)
              XIT(J) = P(J) - PT(J)
              PT(J) = P(J)
          Next J
          FPTT = FUNC2(PTT(), N)
        Loop While FPTT >= FP
        DUM = FP - 2 * FRET + FPTT
        T = 2# * DUM * (FP - FRET - DEL) ^ 2 - DEL * (FP - FPTT) ^ 2
      Loop While T >= 0#
      Call LINMIN(P(), XIT(), N, FRET)
      For J = 1 To N
          XI(J, IBIG) = XIT(J)
      Next J
    Loop
End Sub