nitgm2 Subroutine

subroutine nitgm2(n, xcur, fcur, step, eta, f, jacv, rpar, ipar, ijacv, irpre, iksmax, iresup, ifdord, nfe, njve, nrpre, nli, kdmax, kdmaxp1, vv, rr, svbig, svsml, w, rwork, rsnrm, dinpr, dnorm, itrmks)

FCNRM is no longer an argument, as it needs to be calculated from the residual, which isn't known until a little way into this routine. It could be calculated by the calling routine, but that would mean doing the calculation twice.

step = vector of length n, (trial) step.

eta = relative residual reduction factor.

f = name of user-supplied subroutine for evaluating the function the zero of which is sought; this routine has the form

              subroutine f(n, xcur, fcur, rpar, ipar, itrmf)

       where xcur is the array containing the current x value, fcur 
       is f(xcur) on output, rpar and ipar are, respectively, real 
       and integer parameter/work arrays for use by the subroutine,
       and itrmf is an integer termination flag.  The meaning of
       itrmf is as follows:
         0 => normal termination; desired function value calculated.
         1 => failure to produce f(xcur).

jacv = name of user-supplied subroutine for evaluating Jv or P(inverse)v, where J is the Jacobian of f and P is a right preconditioning operator. If neither analytic Jv evaluations nor right preconditioning is used, this can be a dummy subroutine; if right preconditioning is used but not analytic Jv evaluations, this need only evaluate P(inverse)*v. The form is

        subroutine jacv(n, xcur, fcur, ijob, v, z, rpar, ipar, itrmjv)

        where xcur and fcur are vectors of length n containing the 
        current x and f values, ijob is an integer flag indicating 
        which product is desired, v is a vector of length n to be 
        multiplied, z is a vector of length n containing the desired 
        product on output, rpar and ipar are, respectively, real 
        and integer parameter/work arrays for use by the subroutine, 
        and itrmjv is an integer termination 
        flag. The meaning of ijob is as follows: 
          0 => z = J*v
          1 => z = P(inverse)*v 
        The meaning of itrmjv is as follows:
          0 => normal termination; desired product evaluated. 
          1 => failure to produce J*v.
          2 => failure to produce P(inverse)*v. 
        This subroutine is called only from nitjv, and is always 
        called with v .ne. z.

rpar = real parameter/work array passed to the f and jacv routines.

ipar = integer parameter/work array passed to the f and jacv routines.

ijacv = flag for determining method of J*v evaluation. 0 => finite-difference evaluation (default). 1 => analytic evaluation.

irpre = flag for right preconditioning. 0 => no right preconditioning 1 => right preconditioning

iksmax = maximum allowable number of GMRES iterations.

iresup = residual update flag; on GMRES restarts, the residual is updated as follows: 0 => linear combination (default) 1 => direct evaluation The first is cheap (one n-vector saxpy) but may lose accuracy with extreme residual reduction; the second retains accuracy better but costs one J*v product

ifdord = order of the finite-difference formula (sometimes) used on GMRES restarts when J*v products are evaluated using finite- differences. When ijacv = 0 on input to nitsol, ifdord is set to 1, 2, or 4 in nitsol; otherwise, it is irrelevant. When ijacv = 0 on input to this subroutine, the precise meaning is as follows:

        With GMRES, ifdord matters only if iresup = 1, in which case 
        it determines the order of the finite-difference formula used 
        in evaluating the initial residual at each GMRES restart 
        (default 2). If iresup = 1 and ijacv = 0 on input to this 
        subroutine, then ijacv is temporarily reset to -1 at each 
        restart below to force a finite-difference evaluation of order 
        ifdord. NOTE: This only affects initial residuals at restarts; 
        first-order differences are always used within each GMRES 
        cycle. Using higher-order differences at restarts only should 
        give the same accuracy as if higher-order differences were 
        used throughout; see K. Turner and H. F. Walker, "Efficient 
        high accuracy solutions with GMRES(m)," SIAM J. Sci. Stat. 
        Comput., 13 (1992), pp. 815--825.

nfe = number of function evaluations.

njve = number of J*v evaluations.

nrpre = number of P(inverse)*v evaluations.

nli = number of linear iterations.

kdmax = maximum Krylov subspace dimension; default 10.

kdmaxp1 = kdmax + 1.

vv = n x (kdmax+1) matrix for storage of Krylov basis in GMRES; on return, the residual vector is contained in the first column.

rr = kdmax x kdmax matrix for storage of triangular matrix in GMRES.

svbig = vector of length kdmax for storage of estimate of singular vector of rr with largest singular value.

svsml = vector of length kdmax for storage of estimate of singular vector of rr with smallest singular value.

w = vector of length kdmax, contains right-hand side of triangular system and least-squares residual norm in GMRES.

rwork = vector of length n, work array.

rsnrm = GMRES residual norm on return.

dinpr = inner-product routine, either user-supplied or blas ddot.

dnorm = norm routine, either user-supplied or blas dnrm2.

itrmks = termination flag; values have the following meanings: 0 => normal termination: acceptable step found. 1 => Jv failure in nitjv. 2 => P(inverse)v failure in nitjv. 3 => acceptable step not found in iksmax GMRES iterations. 4 => insufficient residual norm reduction over a cycle of kdmax steps (stagnation) before an acceptable step has been found. 5 => dangerous ill-conditioning detected before an acceptable step has been found.

         Note: On return, nitsol terminates if itrmks is 1 or 2. If  
         itrmks is 3, 4, or 5, nitsol may terminate or continue. In 
         this event, a meaningful inexact Newton step is returned, 
         even though the desired inexact Newton condition may not 
         hold, and a decision on termination/continuation is made 
         in nitdrv according to whether there is sufficient residual 
         norm reduction.

Subroutines required by this and all called routines:

user supplied: f, jacv

nitsol routines: nitjv, nitfd

lapack routine: dlaic1, dlamch

blas routines: daxpy, dcopy, dscal

user supplied or blas: dinpr, dnorm

explanation: In nitsol, dinpr and dnorm are set to either the blas 
ddot and dnrm2 routines or the user-supplied usrnpr and usrnrm 
routines.

This subroutine called by: nitdrv

Subroutines called by this subroutine: daxpy, dcopy, dscal, dlaic1, dlamch, dinpr, dnorm, nitjv

Common block:

<<<<<<<<<<<<<<<<<<<<< TAKEN FROM nitprint.h >>>>>>>>>>>>>>>>>>>>>>>>>

Include these declaratinos and common blocks directly so that the include files are not needed in isoft.

 Do not set the step to 0, as it contains an initial guess at the 
 solution.

C do 20 i = 1, n C step(i) = 0.d0 C 20 continue

 Because the initial step is not zero, the residual must include an
 evaluation of the left-hand-side..

C call dcopy(n, fcur, 1, vv(1,1), 1)

Arguments

 iplvl = 0 => no printout
       = 1 => iteration numbers and F-norms
       = 2 => ... + some stats, step norms, and linear model norms
       = 3 => ... + some Krylov solver and backtrack information
       = 4 => ... + more Krylov solver and backtrack information

 ipunit = printout unit number.