/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:31:16 */ /*FOR_C Options SET: ftn=u io=c no=p op=aimnv s=dbov str=l x=f - prototypes */ #include #include "fcrt.h" #include "dnlagu.h" #include /* PARAMETER translations */ #define D 27 #define J 70 #define NEXTV 47 #define NFCALL 6 #define NFGCAL 7 #define R 61 #define REGD 67 #define REGD0 82 #define TOOBIG 2 #define VNEED 4 /* end of PARAMETER translations */ void /*FUNCTION*/ dnlagu( long n, long p, double x[], void (*dcalcr)(long,long,double[],long*,double[]), void (*dcalcj)(long,long,double[],long*,double[]), long iv[], long liv, long lv, double v[]) { long int d1, dr1, iv1, n1, n2, nf, r1, rd1; /* OFFSET Vectors w/subscript range: 1 to dimension */ long *const Iv = &iv[0] - 1; double *const V = &v[0] - 1; double *const X = &x[0] - 1; /* end of OFFSET VECTORS */ /* Copyright (c) 1996 California Institute of Technology, Pasadena, CA. * ALL RIGHTS RESERVED. * Based on Government Sponsored Research NAS7-03001. *>> 1996-04-27 DNLAGU Krogh Changes to get desired C prototypes. *>> 1994-10-20 DNLAGU Krogh Changes to use M77CON *>> 1990-07-02 DNLAGU C. L. Lawson, JPL *>> 1990-01-31 C. L. Lawson, JPL * * *** VERSION OF NL2SOL THAT CALLS DRN2G *** * * *** PARAMETERS *** * */ /*/ */ /* *** PARAMETER USAGE *** * * N....... TOTAL NUMBER OF RESIDUALS. * P....... NUMBER OF PARAMETERS (COMPONENTS OF X) BEING ESTIMATED. * X....... PARAMETER VECTOR BEING ESTIMATED (INPUT = INITIAL GUESS, * OUTPUT = BEST VALUE FOUND). * DCALCR... SUBROUTINE FOR COMPUTING RESIDUAL VECTOR. * DCALCJ... SUBROUTINE FOR COMPUTING JACOBIAN MATRIX = MATRIX OF FIRST * PARTIALS OF THE RESIDUAL VECTOR. * IV...... INTEGER VALUES ARRAY. * LIV..... LENGTH OF IV (SEE DISCUSSION BELOW). * LV...... LENGTH OF V (SEE DISCUSSION BELOW). * V....... FLOATING-POINT VALUES ARRAY. * * * *** DISCUSSION *** * * NOTE... NL2SOL (MENTIONED BELOW) IS A CODE FOR SOLVING * NONLINEAR LEAST-SQUARES PROBLEMS. IT IS DESCRIBED IN * ACM TRANS. MATH. SOFTWARE, VOL. 9, PP. 369-383 (AN ADAPTIVE * NONLINEAR LEAST-SQUARES ALGORITHM, BY J.E. DENNIS, D.M. GAY, * AND R.E. WELSCH). * * LIV GIVES THE LENGTH OF IV. IT MUST BE AT LEAST 82+P. IF NOT, * THEN DNLAGU RETURNS WITH IV(1) = 15. WHEN DNLAGU RETURNS, THE * MINIMUM ACCEPTABLE VALUE OF LIV IS STORED IN IV(LASTIV) = IV(44), * (PROVIDED THAT LIV .GE. 44). * * LV GIVES THE LENGTH OF V. THE MINIMUM VALUE FOR LV IS * LV0 = 105 + P*(N + 2*P + 17) + 2*N. IF LV IS SMALLER THAN THIS, * THEN DNLAGU RETURNS WITH IV(1) = 16. WHEN DNLAGU RETURNS, THE * MINIMUM ACCEPTABLE VALUE OF LV IS STORED IN IV(LASTV) = IV(45) * (PROVIDED LIV .GE. 45). * * RETURN CODES AND CONVERGENCE TOLERANCES ARE THE SAME AS FOR * NL2SOL, WITH SOME SMALL EXTENSIONS... IV(1) = 15 MEANS LIV WAS * TOO SMALL. IV(1) = 16 MEANS LV WAS TOO SMALL. * * THERE ARE TWO NEW V INPUT COMPONENTS... V(LMAXS) = V(36) AND * V(SCTOL) = V(37) SERVE WHERE V(LMAX0) AND V(RFCTOL) FORMERLY DID * IN THE SINGULAR CONVERGENCE TEST -- SEE THE NL2SOL DOCUMENTATION. * * *** DEFAULT VALUES *** * * DEFAULT VALUES ARE PROVIDED BY SUBROUTINE DIVSET, RATHER THAN * DFAULT. THE CALLING SEQUENCE IS... * CALL DIVSET(1, IV, LIV, LV, V) * THE FIRST PARAMETER IS AN INTEGER 1. IF LIV AND LV ARE LARGE * ENOUGH FOR DIVSET, THEN DIVSET SETS IV(1) TO 12. OTHERWISE IT * SETS IV(1) TO 15 OR 16. CALLING DNLAGU WITH IV(1) = 0 CAUSES ALL * DEFAULT VALUES TO BE USED FOR THE INPUT COMPONENTS OF IV AND V. * IF YOU FIRST CALL DIVSET, THEN SET IV(1) TO 13 AND CALL DNLAGU, * THEN STORAGE ALLOCATION ONLY WILL BE PERFORMED. IN PARTICULAR, * IV(D) = IV(27), IV(J) = IV(70), AND IV(R) = IV(61) WILL BE SET * TO THE FIRST SUBSCRIPT IN V OF THE SCALE VECTOR, THE JACOBIAN * MATRIX, AND THE RESIDUAL VECTOR RESPECTIVELY, PROVIDED LIV AND LV * ARE LARGE ENOUGH. IF SO, THEN DNLAGU RETURNS WITH IV(1) = 14. * WHEN CALLED WITH IV(1) = 14, DNLAGU ASSUMES THAT STORAGE HAS * BEEN ALLOCATED, AND IT BEGINS THE MINIMIZATION ALGORITHM. * * *** SCALE VECTOR *** * * ONE DIFFERENCE WITH NL2SOL IS THAT THE SCALE VECTOR D IS * STORED IN V, STARTING AT A DIFFERENT SUBSCRIPT. THE STARTING * SUBSCRIPT VALUE IS STILL STORED IN IV(D) = IV(27). THE * DISCUSSION OF DEFAULT VALUES ABOVE TELLS HOW TO HAVE IV(D) SET * BEFORE THE ALGORITHM IS STARTED. * * *** REGRESSION DIAGNOSTICS *** * * IF IV(RDREQ) SO DICTATES, THEN ESTIMATES ARE COMPUTED OF THE * INFLUENCE EACH RESIDUAL COMPONENT HAS ON THE FINAL PARAMETER * ESTIMATE X. THE GENERAL IDEA IS THAT ONE MAY WISH TO EXAMINE * RESIDUAL COMPONENTS (AND THE DATA BEHIND THEM) FOR WHICH THE * INFLUENCE ESTIMATE IS SIGNIFICANTLY LARGER THAN MOST OF THE OTHER * INFLUENCE ESTIMATES. THESE ESTIMATES, HEREAFTER CALLED * REGRESSION DIAGNOSTICS, ARE ONLY COMPUTED IF IV(RDREQ) = 2 OR 3. * IN THIS CASE, FOR I = 1(1)N, * SQRT( G(I)**T * H(I)**-1 * G(I) ) * IS COMPUTED AND STORED IN V, STARTING AT V(IV(REGD)), WHERE * RDREQ = 57 AND REGD = 67. HERE G(I) STANDS FOR THE GRADIENT * RESULTING WHEN THE I-TH OBSERVATION IS DELETED AND H(I) STANDS * FOR AN APPROXIMATION TO THE CORRESPONDING HESSIAN AT X, THE SOLU- * TION CORRESPONDING TO ALL OBSERVATIONS. (THIS APPROXIMATION IS * OBTAINED BY SUBTRACTING THE FIRST-ORDER CONTRIBUTION OF THE I-TH * OBSERVATION TO THE HESSIAN FROM A FINITE-DIFFERENCE HESSIAN * APPROXIMATION. IF H IS INDEFINITE, THEN IV(REGD) IS SET TO -1. * IF H(I) IS INDEFINITE, THEN -1 IS RETURNED AS THE DIAGNOSTIC FOR * OBSERVATION I. IF NO DIAGNOSTICS ARE COMPUTED, PERHAPS BECAUSE * OF A FAILURE TO CONVERGE, THEN IV(REGD) = 0 IS RETURNED.) * PRINTING OF THE REGRESSION DIAGNOSTICS IS CONTROLLED BY * IV(COVPRT) = IV(14)... IF IV(COVPRT) = 3, THEN BOTH THE * COVARIANCE MATRIX AND THE REGRESSION DIAGNOSTICS ARE PRINTED. * IV(COVPRT) = 2 CAUSES ONLY THE REGRESSION DIAGNOSTICS TO BE * PRINTED, IV(COVPRT) = 1 CAUSES ONLY THE COVARIANCE MATRIX TO BE * PRINTED, AND IV(COVPRT) = 0 CAUSES NEITHER TO BE PRINTED. * * RDREQ = 57 AND REGD = 67. * */ /* *** GENERAL *** * * CODED BY DAVID M. GAY. * *+++++++++++++++++++++++++++ DECLARATIONS +++++++++++++++++++++++++++ * * *** EXTERNAL SUBROUTINES *** * */ /*--D replaces "?": ?NLAGU, ?RN2G, ?IVSET, ?N2RDP, ?CALCR, ?CALCJ * * DIVSET.... PROVIDES DEFAULT IV AND V INPUT COMPONENTS. * DRN2G ... CARRIES OUT OPTIMIZATION ITERATIONS. * DN2RDP... PRINTS REGRESSION DIAGNOSTICS. * * *** NO INTRINSIC FUNCTIONS *** * * *** LOCAL VARIABLES *** * */ /* *** IV COMPONENTS *** * */ /* -------------------------------- BODY ------------------------------ * */ if (Iv[1] == 0) divset( 1, iv, liv, lv, v ); iv1 = Iv[1]; if (iv1 == 14) goto L_10; if (iv1 > 2 && iv1 < 12) goto L_10; if (iv1 == 12) Iv[1] = 13; if (Iv[1] == 13) Iv[VNEED] += p + n*(p + 2); drn2g( x, v, iv, liv, lv, n, n, &n1, &n2, p, v, v, v, x ); if (Iv[1] != 14) goto L_999; /* *** STORAGE ALLOCATION *** * */ Iv[D] = Iv[NEXTV]; Iv[R] = Iv[D] + p; Iv[REGD0] = Iv[R] + n; Iv[J] = Iv[REGD0] + n; Iv[NEXTV] = Iv[J] + n*p; if (iv1 == 13) goto L_999; L_10: d1 = Iv[D]; dr1 = Iv[J]; r1 = Iv[R]; rd1 = Iv[REGD0]; L_20: drn2g( &V[d1], &V[dr1], iv, liv, lv, n, n, &n1, &n2, p, &V[r1], &V[rd1], v, x ); switch (IARITHIF(Iv[1] - 2)) { case -1: goto L_30; case 0: goto L_50; case 1: goto L_60; } /* *** NEW FUNCTION VALUE (R VALUE) NEEDED *** * */ L_30: nf = Iv[NFCALL]; (*dcalcr)( n, p, x, &nf, &V[r1] ); if (nf > 0) goto L_40; /* CALL DCALCR(N, P, X, NF, V(R1)) */ Iv[TOOBIG] = 1; goto L_20; L_40: if (Iv[1] > 0) goto L_20; /* *** COMPUTE DR = GRADIENT OF R COMPONENTS *** * */ L_50: ; (*dcalcj)( n, p, x, &Iv[NFGCAL], &V[dr1] ); if (Iv[NFGCAL] == 0) Iv[TOOBIG] = 1; /* CALL DCALCJ(N, P, X, IV(NFGCAL), V(DR1)) */ goto L_20; /* *** INDICATE WHETHER THE REGRESSION DIAGNOSTIC ARRAY WAS COMPUTED * *** AND PRINT IT IF SO REQUESTED... * */ L_60: if (Iv[REGD] > 0) Iv[REGD] = rd1; dn2rdp( iv, liv, n, &V[rd1] ); L_999: return; /* *** LAST CARD OF DNLAGU FOLLOWS *** */ } /* end of function */