Attribute VB_Name = "ModNewtonMethod"
Option Explicit
Public Sub NewtonMethod(canvas As PictureBox, AAA As String, nEq As Integer, stEqs() As String, VvarsSt() As String, VvarsI() As Double, NTry As Integer, Tol As Double, XAns() As Double)
'C***********************************************************************
'C                                                                      *
'C NEWTON            'S METHOD FOR SYSTEMS ALGORITHM 10.1               *
'C                                                                      *
'C***********************************************************************
'C
'C
'C
'C     TO APPROXIMATE THE SOLUTION OF THE NONLINEAR SYSTEM F(X)=0 GIVEN
'C     AN INITIAL APPROXIMATION X:
'C
'C     INPUT:   NUMBER n OF EQUATIONS AND UNKNOWNS; INITIAL APPROXIMATION
'C              X=(X(1),...,X(n)); TOLERANCE TOL; MAXIMUM NUMBER OF
'C              ITERATIONS N.
'C
'C     OUTPUT:  APPROXIMATE SOLUTION X=(X(1),...,X(n)) OR A MESSAGE
'C              THAT THE NUMBER OF ITERATIONS WAS EXCEEDED.
'C
'C INITIALIZATION
      Dim AA() As Double, yY() As Double, Xx() As Double
      ReDim AA(1 To nEq, 1 To nEq + 1), yY(1 To nEq), Xx(1 To nEq)
      Dim i As Integer, K As Integer, j As Integer
      Dim R
'C     USE AA FOR J;  NN FOR MAXIMUM NUMBER OF ITERATIONS N
    
      Dim OK As Boolean
      'OPEN(UNIT=5,FILE='CON',ACCESS='SEQUENTIAL')
      'OPEN(UNIT=6,FILE='CON',ACCESS='SEQUENTIAL')
      'WRITE(6,*) 'This is the Newton Method for Nonlinear Systems.'
      'WRITE(6,*) 'Has the function F been defined and has the'
      'WRITE(6,*) 'Jacobian of partial derivatives been defined'
      'WRITE(6,*) 'at lines 105 and 115 respectively.'
      ''WRITE(6,*) 'Enter Y or N '
      'WRITE(6,*) ' '
      'READ(5,*)  AAA
      If (LCase(AAA) = "y") Then
         OK = False
101      If (OK) Then GoTo 11
            'WRITE(6,*) 'Input the number n of equations.'
            'WRITE(6,*) ' '
            'EAD(5,*) NEq
            If (nEq >= 2) Then
               OK = True
            Else
               MsgBox ("NEq must be greater than 1.")
            End If
         GoTo 101
11       OK = False
14       If OK Then GoTo 15
            'WRITE(6,*) 'Input tolerance '
            'WRITE(6,*) ' '
            'READ(5,*) TOL
            If (Tol <= 0#) Then
               MsgBox ("Tolerance must be positive ")
            Else
               OK = True
            End If
         GoTo 14
15       OK = False
12       If OK Then GoTo 13
           ' WRITE(6,*) 'Input the maximum number of iterations.'
            'WRITE(6,*) ' '
            'READ(5,*) Ntry
            If (NTry <= 0) Then
              MsgBox ("NTry Must be positive integer ")
            Else
              OK = True
            End If
         GoTo 12
13       For i = 1 To nEq 'DO 16 I=1,NEq
            'WRITE(6,*) 'Input initial approximation X(',I,')'
            'WRITE(6,*) ' '
            'READ(5,*) X(I)
            Xx(i) = VvarsI(i)
            canvas.Print "i=" & i & "x(" & i & ")=" & Xx(i)
16       Next i ' CONTINUE
      Else
        MsgBox ("program will end so that the functions F(X) and J(X) can be created ")
         OK = False
      End If
      
      If Not OK Then GoTo 400
     
      canvas.Print "NEWTON METHOD FOR NONLINEAR SYSTEMS"
      'Pi = 4 * ATAN(1#)
'C     STEP 1
      K = 1
'C     STEP 2
100   If (K > NTry) Then GoTo 200
'C STEP3
'C          COMPUTE J(X)

           'AA(1, 1) = 3#
           'AA(1, 2) = x(3) * Sin(x(2) * x(3))
           'AA(1, 3) = x(2) * Sin(x(2) * x(3))
           
           'AA(2, 1) = 2 * x(1)
           'AA(2, 2) = -162 * (x(2) + 0.1)
           'AA(2, 3) = Cos(x(3))
           
           'AA(3, 1) = -x(2) * Exp(-x(1) * x(2))
           'AA(3, 2) = -x(1) * Exp(-x(1) * x(2))
           'AA(3, 3) = 20
           
           
           For i = 1 To nEq
               ' MsgBox ("steq= " & stEqs(i) & " ;xvst= " & VvarsSt(i) & " ; xv= " & xX(i))
                For j = 1 To nEq
                AA(i, j) = PartialDif(glVbScriptA(i - 1), stEqs(i), VvarsSt, Xx, j)
                'canvas.Print "AA(" & I & "," & J & ")= " & AA(I, J)
                Next j
           Next i
           
'C COMPUTE - F(X)
           'AA(1, 4) = -(3 * x(1) - Cos(x(2) * x(3)) - 0.5)
           'AA(2, 4) = -(x(1) ^ 2 - 81 * (x(2) + 0.1) ^ 2 + Sin(x(3)) + 1.06)
           'AA(3, 4) = -(Exp(-x(1) * x(2)) + 20 * x(3) + (10 * PI - 3) / 3)
           For i = 1 To nEq
           AA(i, nEq + 1) = -FunctionValues(glVbScriptA(i - 1), stEqs(i), VvarsSt, Xx)
           'canvas.Print "AA(" & I & "," & (nEq + 1) & ")= " & AA(I, nEq + 1)
           Next i
           
'C          STEP 4
'C          SOLVES THE N X N LINEAR SYSTEM J(X) Y = -F(X)
           Call LIN(nEq, nEq + 1, AA, yY)
'C          STEPS 5 AND 6
'C          R = INFINITY NORM OF Y

           R = 0
           For i = 1 To nEq 'DO 10 I=1,N
                If (Abs(yY(i)) > R) Then R = Abs(yY(i))
                ' canvas.Print "«ek=" & K & "; x(" & I & ")=" & xX(I)
10              Xx(i) = Xx(i) + yY(i)
                'canvas.Print "«ák=" & K & "; x(" & I & ")=" & xX(I)
           Next i
           canvas.Print "IterNo=" & K & "Error=" & R
           'Form1.Caption = "IterNo=" & K & " ;Error=" & R
           For i = 1 To nEq
           canvas.Print "x(" & i & ")= " & Xx(i)
           ReDim Preserve XAns(1 To i)
           XAns(i) = Xx(i)
           Next i
'C          STEP 6
           If (R < Tol) Then
'C PROCESS Is COMPLETE
                '?MsgBox ("SUCCESS WITHIN TOLERANCE 1.0E-05=" & R)
                GoTo 400
           End If
'C          STEP 7
           K = K + 1
      GoTo 100
      Exit Sub
'C     STEP 8
'C DIVERGENCE
200    MsgBox ("DIVERGENCE - STOPPED AFTER ITERNN" & NTry)
400
      End Sub
           Sub LIN(n As Integer, M As Integer, A() As Double, x() As Double)

         
'C          LIN SOLVES THE LINEAR SYSTEM: J(X) Y = -F(X) USING PARTIAL
'C PIVOTING
           ReDim Preserve A(1 To n, 1 To M), x(1 To n)
           Dim i As Integer
           Dim Zero As Double, L As Integer, K As Integer, IR As Integer, IA As Integer, j As Integer, JA As Integer
           Dim y  As Double, C As Double
           Zero = 1E-20
           K = n - 1
           For i = 1 To K 'DO 10 I=1,K
                y = Abs(A(i, i))
                IR = i
                IA = i + 1
                For j = IA To n ' DO 20 J=IA,N
                     If (Abs(A(j, i)) > y) Then
                          IR = j
                          y = Abs(A(j, i))
                     End If
20             Next j ' CONTINUE

                If (y <= Zero) Then
                     MsgBox ("LINEAR SYSTEM HAS NO SOLUTION") 'WRITE(6,1)
                     Exit Sub
                End If
                
                If (IR <> i) Then
                     For j = i To M 'DO 30 J=I,M
                          C = A(i, j)
                          A(i, j) = A(IR, j)
30                   A(IR, j) = C
                     Next j
                End If
                
                For j = IA To n ' DO 40 J=IA,N
                     C = A(j, i) / A(i, i)
                     If (Abs(C) <= Zero) Then C = 0
                     For L = i To M 'DO 40 L=I,M
40                   A(j, L) = A(j, L) - C * A(i, L)
                     Next L
                Next j
                
10         Next i 'CONTINUE

           If (Abs(A(n, n)) < 1E-20) Then
                 MsgBox ("LINEAR SYSTEM HAS NO SOLUTION") 'WRITE(6,1)WRITE(6,1)
                Exit Sub
           End If
      x(n) = A(n, M) / A(n, n)
      
      For i = 1 To K 'DO 50 I=1,K
           j = n - i
           JA = j + 1
           C = A(j, M)
           For L = JA To n 'DO 60 L=JA,N
60        C = C - A(j, L) * x(L)
         Next L
50    x(j) = C / A(j, j)
      Next i
     
'1     FORMAT(1X,'LINEAR SYSTEM HAS NO SOLUTION')
      
 
  End Sub
'-------------------------------------------

