/* nllsq.f -- translated by f2c (version 19950314).
   You must link the resulting object file with the libraries:
	-lF77 -lI77 -lm   (in that order)
*/

#include "f2c.h"

/* Common Block Declarations */

struct {
    char argz[512];
} argz_;

#define argz_1 argz_

struct {
    char keywrd[80];
} keywrd_;

#define keywrd_1 keywrd_

struct {
    integer iflepo, iiter;
} mesage_;

#define mesage_1 mesage_

struct {
    doublereal time0;
} time_;

#define time_1 time_

struct {
    integer ncount;
} nllsqi_;

#define nllsqi_1 nllsqi_

struct {
    integer nscf;
} numscf_;

#define numscf_1 numscf_

struct {
    integer last;
} last_;

#define last_1 last_

struct {
    doublereal dddum[6], efslst[258], q[66564]	/* was [258][258] */, r[66564]
	    	/* was [258][258] */, xlast[258];
    integer iiium[7], idumy[132077];
} nllcom_;

#define nllcom_1 nllcom_

/* Table of constant values */

static integer c__1 = 1;
static integer c__0 = 0;
static logical c_true = TRUE_;
static integer c__2 = 2;

/* Subroutine */ int nllsq_(doublereal *x, integer *n)
{

    /* Format strings */
    static char fmt_640[] = "(\002  RESTART FILE WRITTEN,  TIME LEFT:\002,f9"
	    ".1,\002 GRAD.:\002,f10.3,\002 HEAT:\002,g14.7)";
    static char fmt_650[] = "(\002 CYCLE:\002,i5,\002 TIME:\002,f6.1,\002 TI"
	    "ME LEFT:\002,f9.1,\002 GRAD.:\002,f10.3,\002 HEAT:\002,g14.7)";
    static char fmt_810[] = "(\0020TEST ON X SATISFIED, NUMBER OF FUNCTION C"
	    "ALLS = \002,i5)";
    static char fmt_820[] = "(\0020TEST ON SSQ SATISFIED, NUMBER OF FUNCTION"
	    " CALLS = \002,i5)";
    static char fmt_720[] = "(\002 \002,5x,\002ATTEMPT TO GO DOWNHILL IS UNS"
	    "UCCESSFUL AFTER\002,i5,5x,\002ORTHOGONAL SEARCHES\002)";

    /* System generated locals */
    integer i__1, i__2, i__3;
    doublereal d__1, d__2;
    cllist cl__1;

    /* Builtin functions */
    integer i_indx(char *, char *, ftnlen, ftnlen), s_wsfe(cilist *), do_fio(
	    integer *, char *, ftnlen), e_wsfe(void), s_cmp(char *, char *, 
	    ftnlen, ftnlen), f_clos(cllist *);
    double sqrt(doublereal), d_sign(doublereal *, doublereal *);
    /* Subroutine */ int s_stop(char *, ftnlen);

    /* Local variables */
    static doublereal escf;
    static integer jend;
#define icyc ((integer *)&nllcom_1 + 267301)
    static integer nrem, ierr;
    static doublereal temp;
#define irst ((integer *)&nllcom_1 + 267302)
#define jrst ((integer *)&nllcom_1 + 267303)
    static integer ixso;
    static doublereal work, ttmp;
    static integer nrst, nsst;
    static doublereal tolx, time1, time2, tols1, tols2, tols5, tols6;
    static integer i, j, k;
    extern doublereal reada_(char *, integer *, ftnlen);
    static integer m;
    static doublereal p[258], t, y[258], tleft;
    static integer ifrtl;
    static doublereal ytail, efsss;
    extern /* Subroutine */ int geout_(void);
    static doublereal const__, tlast, tdump;
    static integer i1, j1;
    static char ch[1];
    static integer ii, jj;
    static logical middle;
#define pn ((doublereal *)&nllcom_1 + 3)
    static doublereal yn;
    extern doublereal second_(void);
    extern /* Subroutine */ int compfg_(doublereal *, logical *, doublereal *,
	     logical *, doublereal *, logical *), locmin_(integer *, 
	    doublereal *, integer *, doublereal *, doublereal *, doublereal *,
	     doublereal *, integer *, doublereal *);
    static doublereal tcycle;
    static logical resfil;
    static integer irepet;
    extern /* Subroutine */ int parsav_(integer *, integer *, integer *);
    static integer mxcycl;
    static doublereal pnlast;
    static integer iprint, np1, np2;
    static doublereal pn2, ymaxst, ssqlst;
#define alf ((doublereal *)&nllcom_1 + 1)
    static doublereal bet, efs[258], cos__;
    extern doublereal dot_(doublereal *, doublereal *, integer *);
    static doublereal eps, tim, sin__;
    static integer nto;
    static doublereal tmp, prt;
#define ssq ((doublereal *)&nllcom_1 + 2)
    static logical scf1;
    static doublereal tol2;

    /* Fortran I/O blocks */
    static cilist io___15 = { 0, 6, 0, "(/,A)", 0 };
    static cilist io___18 = { 0, 6, 0, "(/10X,'NUMBER OF CYCLES TO BE RUN ='"
	    ",I5)", 0 };
    static cilist io___61 = { 0, 6, 0, "(' SYSTEM DOES NOT APPEAR TO BE OPTI"
	    "MIZABLE.',/     ,' THIS CAN HAPPEN IF (A) IT WAS OPTIMIZED TO BE"
	    "GIN WITH',/     ,' OR                 (B) IT IS NEITHER A GROUND"
	    " NOR A',        ' TRANSITION STATE')", 0 };
    static cilist io___74 = { 0, 6, 0, fmt_640, 0 };
    static cilist io___75 = { 0, 6, 0, fmt_650, 0 };
    static cilist io___76 = { 0, 6, 0, fmt_810, 0 };
    static cilist io___77 = { 0, 6, 0, fmt_820, 0 };
    static cilist io___78 = { 0, 6, 0, fmt_720, 0 };


/* COMDECK SIZES */
/************************************************************************
****/
/*  THIS FILE CONTAINS ALL THE ARRAY SIZES FOR USE IN MOPAC.              
**/
/*                                                                        
**/
/*    THERE ARE ONLY  PARAMETERS THAT THE PROGRAMMER NEED SET:            
**/
/*    MAXHEV = MAXIMUM NUMBER OF HEAVY ATOMS (HEAVY: NON-HYDROGEN ATOMS)  
**/
/*    MAXLIT = MAXIMUM NUMBER OF HYDROGEN ATOMS.                          
**/
/*    MAXTIM = DEFAULT TIME FOR A JOB. (SECONDS)                          
**/
/*    MAXDMP = DEFAULT TIME FOR AUTOMATIC RESTART FILE GENERATION (SECS)  
**/
/*                                                                        
**/
/*                                                                        
**/
/************************************************************************
****/
/*                                                                        
**/
/*  THE FOLLOWING CODE DOES NOT NEED TO BE ALTERED BY THE PROGRAMMER      
**/
/*                                                                        
**/
/************************************************************************
****/
/*                                                                        
**/
/*   ALL OTHER PARAMETERS ARE DERIVED FUNCTIONS OF THESE TWO PARAMETERS   
**/
/*                                                                        
**/
/*     NAME                   DEFINITION                                  
**/
/*    NUMATM         MAXIMUM NUMBER OF ATOMS ALLOWED.                     
**/
/*    MAXORB         MAXIMUM NUMBER OF ORBITALS ALLOWED.                  
**/
/*    MAXPAR         MAXIMUM NUMBER OF PARAMETERS FOR OPTIMISATION.       
**/
/*    N2ELEC         MAXIMUM NUMBER OF TWO ELECTRON INTEGRALS ALLOWED.    
**/
/*    MPACK          AREA OF LOWER HALF TRIANGLE OF DENSITY MATRIX.       
**/
/*    MORB2          SQUARE OF THE MAXIMUM NUMBER OF ORBITALS ALLOWED.    
**/
/*    MAXHES         AREA OF HESSIAN MATRIX                               
**/
/************************************************************************
****/
/************************************************************************
****/
/*  FOR SHORT VERSION USE LINE WITH NMECI=1, FOR LONG VERSION USE LINE    
**/
/*  WITH NMECI=10                                                         
**/
/************************************************************************
****/
/*     PARAMETER (NMECI=1,   NPULAY=1) */
/************************************************************************
****/
/* DECK MOPAC */
/* next line added for Unix implementation for command line arguments */

/* ***********************************************************************
 */

/*  NLLSQ IS A NON-DERIVATIVE, NONLINEAR LEAST-SQUARES MINIMIZER. IT USES 
*/
/*        BARTEL'S PROCEDURE TO MINIMIZE A FUNCTION WHICH IS A SUM OF */
/*        SQUARES. */

/*    ON INPUT N    = NUMBER OF UNKNOWNS */
/*             X    = PARAMETERS OF FUNCTION TO BE MINIMIZED. */

/*    ON EXIT  X    = OPTIMIZED PARAMETERS. */

/*    THE FUNCTION TO BE MINIMIZED IS "COMPFG". COMPFG MUST HAVE THE */
/*    CALLING SEQUENCE */
/*                  CALL COMPFG(XPARAM,.TRUE.,ESCF,.TRUE.,EFS,.TRUE.) */
/*                  SSQ=DOT(EFS,EFS,N) */
/*    WHERE   EFS  IS A VECTOR WHICH  COMPFG  FILLS WITH THE N INDIVIDUAL 
*/
/*                 COMPONENTS OF THE ERROR FUNCTION AT THE POINT X */
/*            SSQ IS THE VALUE OF THE SUM OF THE  EFS  SQUARED. */
/*    IN THIS FORMULATION OF NLLSQ M AND N ARE THE SAME. */
/*    THE PRECISE DEFINITIONS OF THESE TWO QUANTITIES IS: */

/*     N = NUMBER OF PARAMETERS TO BE OPTIMIZED. */
/*     M = NUMBER OF REFERENCE FUNCTIONS. M MUST BE GREATER THEN, OR */
/*         EQUAL TO, N */
/* ***********************************************************************
 */
/*     Q = ORTHOGONAL MATRIX   (M BY M) */
/*     R = RIGHT-TRIANGULAR MATRIX   (M BY N) */
/*     MXCNT(1) = MAX ALLOW OVERALL FUN EVALS */
/*     MXCNT(2) = MAX ALLOW NO OF FNC EVALS PER LIN SEARCH */
/*     TOLS1 = RELATIVE TOLERANCE ON X OVERALL */
/*     TOLS2 = ABSOLUTE TOLERANCE ON X OVERALL */
/*     TOLS5 = RELATIVE TOLERANCE ON X FOR LINEAR SEARCHES */
/*     TOLS6 = ABSOLUTE TOLERANCE ON X FOR LINEAR SEARCHES */
/*     IPRINT = PRINT SWITCH */
/*     NRST = NUMBER OF CYCLES BETWEEN SIDESTEPS */
/*     ********** */
    /* Parameter adjustments */
    --x;

    /* Function Body */
    middle = i_indx(keywrd_1.keywrd, "RESTART", 80L, 7L) != 0;
    scf1 = i_indx(keywrd_1.keywrd, "1SCF", 80L, 4L) != 0;
    mesage_1.iflepo = 10;
/* * */
    m = *n;
/* * */
    tol2 = .4;
    if (i_indx(keywrd_1.keywrd, "GNORM", 80L, 5L) != 0) {
	i__1 = i_indx(keywrd_1.keywrd, "GNORM", 80L, 5L);
	tol2 = reada_(keywrd_1.keywrd, &i__1, 80L);
	if (tol2 < .001 && i_indx(keywrd_1.keywrd, "LET", 80L, 3L) == 0) {
	    s_wsfe(&io___15);
	    do_fio(&c__1, "  GNORM HAS BEEN SET TOO LOW, RESET TO 0   .001", 
		    47L);
	    e_wsfe();
	    tol2 = .001;
	}
    }
    last_1.last = 0;
    mxcycl = 100;
    i = i_indx(keywrd_1.keywrd, "CYCLES", 80L, 6L);
    if (i != 0) {
	mxcycl = (integer) reada_(keywrd_1.keywrd, &i, 80L);
	s_wsfe(&io___18);
	do_fio(&c__1, (char *)&mxcycl, (ftnlen)sizeof(integer));
	e_wsfe();
    }
    tols1 = 1e-12;
    tols2 = 1e-10;
    tols5 = 1e-6;
    tols6 = .001;
    iprint = -1;
    nrst = 4;
    ymaxst = 1.;
    tleft = 3600.;
    i = i_indx(keywrd_1.keywrd, " T=", 80L, 3L);
    if (i != 0) {
	tim = reada_(keywrd_1.keywrd, &i, 80L);
	for (j = i + 3; j <= 80; ++j) {
	    i__1 = j;
	    if (s_cmp(keywrd_1.keywrd + i__1, " ", j + 1 - i__1, 1L) == 0) {
		*(unsigned char *)ch = *(unsigned char *)&keywrd_1.keywrd[j - 
			1];
		if (*(unsigned char *)ch == 'M') {
		    tim *= 60;
		}
		if (*(unsigned char *)ch == 'H') {
		    tim *= 3600;
		}
		if (*(unsigned char *)ch == 'D') {
		    tim *= 86400;
		}
		goto L20;
	    }
/* L10: */
	}
L20:
	tleft = tim;
    }
    tlast = tleft;
    tdump = 3600.;
    i = i_indx(keywrd_1.keywrd, " DUMP", 80L, 5L);
    if (i != 0) {
	tdump = reada_(keywrd_1.keywrd, &i, 80L);
	for (j = i + 6; j <= 80; ++j) {
	    i__1 = j;
	    if (s_cmp(keywrd_1.keywrd + i__1, " ", j + 1 - i__1, 1L) == 0) {
		*(unsigned char *)ch = *(unsigned char *)&keywrd_1.keywrd[j - 
			1];
		if (*(unsigned char *)ch == 'M') {
		    tdump *= 60;
		}
		if (*(unsigned char *)ch == 'H') {
		    tdump *= 3600;
		}
		if (*(unsigned char *)ch == 'D') {
		    tdump *= 86400;
		}
		goto L40;
	    }
/* L30: */
	}
L40:
	;
    }
    resfil = FALSE_;
    tleft = tleft - second_() + time_1.time0;
/*     ********** */
/*     SET UP COUNTERS AND SWITCHES */
/*     ********** */
    nto = *n / 6;
    ifrtl = 0;
    nrem = *n - nto * 6;
    irepet = 0;
    nsst = 0;
    if (ixso == 0) {
	ixso = *n;
    }
    np1 = *n + 1;
    np2 = *n + 2;
    *icyc = 0;
    *irst = 0;
    *jrst = 1;
    eps = tols5;
    t = tols6;
/*     ********** */
/*     GET STARTING-POINT FUNCTION VALUE */
/*     SET UP ESTIMATE OF INITIAL LINE STEP */
/*     ********** */
    if (middle) {
	parsav_(&c__0, n, &m);
	numscf_1.nscf = nllcom_1.iiium[0];
	cl__1.cerr = 0;
	cl__1.cunit = 13;
	cl__1.csta = 0;
	f_clos(&cl__1);
	nllsqi_1.ncount = nllcom_1.iiium[4];
	mxcycl += *icyc;
	i__1 = *n;
	for (i = 1; i <= i__1; ++i) {
/* L50: */
	    x[i] = nllcom_1.xlast[i - 1];
	}
	time1 = second_();
	goto L100;
    }
    compfg_(&x[1], &c_true, &escf, &c_true, nllcom_1.efslst, &c_true);
    if (scf1) {
	goto L920;
    }
    *ssq = dot_(nllcom_1.efslst, nllcom_1.efslst, n);
    nllsqi_1.ncount = 1;
/* L60: */
    i__1 = m;
    for (i = 1; i <= i__1; ++i) {
	i__2 = *n;
	for (j = 1; j <= i__2; ++j) {
	    nllcom_1.r[i + j * 258 - 259] = 0.;
	    if (i == j) {
		nllcom_1.r[i + j * 258 - 259] = 1.;
	    }
/* L70: */
	}
	i__2 = m;
	for (j = i; j <= i__2; ++j) {
	    nllcom_1.q[i + j * 258 - 259] = 0.;
	    nllcom_1.q[j + i * 258 - 259] = 0.;
	    if (i == j) {
		nllcom_1.q[i + i * 258 - 259] = 1.;
	    }
/* L80: */
	}
    }
    temp = 0.;
    i__2 = *n;
    for (i = 1; i <= i__2; ++i) {
/* L90: */
/* Computing 2nd power */
	d__1 = x[i];
	temp += d__1 * d__1;
    }
    *alf = (eps * sqrt(temp) + t) * 100.;
/*     ********** */
/*     MAIN LOOP */
/*     ********** */
    time1 = second_();
L100:
/*     ********** */
/*     UPDATE COUNTERS AND TEST FOR PRINTING THIS CYCLE */
/*     ********** */
    ++ifrtl;
    ++(*icyc);
    ++(*irst);
/*     ********** */
/*     SET  PRT,  THE LEVENBERG-MARQUARDT PARAMETER. */
/*     ********** */
    prt = sqrt(*ssq);
/*     ********** */
/*     IF A SIDESTEP IS TO BE TAKEN, GO TO 31 */
/*     ********** */
    if (*irst >= nrst) {
	goto L250;
    }
/*     ********** */
/*     SOLVE THE SYSTEM    Q*R*P = -EFSLST    IN THE LEAST-SQUARES SENSE 
*/
/*     ********** */
    nsst = 0;
    i__2 = m;
    for (i = 1; i <= i__2; ++i) {
	temp = 0.;
	i__1 = m;
	for (j = 1; j <= i__1; ++j) {
/* L110: */
	    temp -= nllcom_1.q[j + i * 258 - 259] * nllcom_1.efslst[j - 1];
	}
/* L120: */
	efs[i - 1] = temp;
    }
    i__2 = *n;
    for (j = 1; j <= i__2; ++j) {
	jj = np1 - j;
	i__1 = j;
	for (i = 1; i <= i__1; ++i) {
	    ii = np2 - i;
/* L130: */
	    nllcom_1.r[ii + jj * 258 - 259] = nllcom_1.r[i + j * 258 - 259];
	}
    }
    i__1 = *n;
    for (i = 1; i <= i__1; ++i) {
	i1 = i + 1;
	y[i - 1] = prt;
	efsss = 0.;
	if (i >= *n) {
	    goto L150;
	}
	i__2 = *n;
	for (j = i1; j <= i__2; ++j) {
/* L140: */
	    y[j - 1] = 0.;
	}
L150:
	i__2 = *n;
	for (j = i; j <= i__2; ++j) {
	    ii = np2 - j;
	    jj = np1 - j;
	    if ((d__1 = y[j - 1], abs(d__1)) < (d__2 = nllcom_1.r[ii + jj * 
		    258 - 259], abs(d__2))) {
		goto L160;
	    }
/* Computing 2nd power */
	    d__1 = nllcom_1.r[ii + jj * 258 - 259] / y[j - 1];
	    temp = y[j - 1] * sqrt(d__1 * d__1 + 1.);
	    goto L170;
L160:
/* Computing 2nd power */
	    d__1 = y[j - 1] / nllcom_1.r[ii + jj * 258 - 259];
	    temp = nllcom_1.r[ii + jj * 258 - 259] * sqrt(d__1 * d__1 + 1.);
L170:
	    sin__ = nllcom_1.r[ii + jj * 258 - 259] / temp;
	    cos__ = y[j - 1] / temp;
	    nllcom_1.r[ii + jj * 258 - 259] = temp;
	    temp = efs[j - 1];
	    efs[j - 1] = sin__ * temp + cos__ * efsss;
	    efsss = sin__ * efsss - cos__ * temp;
	    if (j >= *n) {
		goto L200;
	    }
	    j1 = j + 1;
	    i__3 = *n;
	    for (k = j1; k <= i__3; ++k) {
		jj = np1 - k;
		temp = nllcom_1.r[ii + jj * 258 - 259];
		nllcom_1.r[ii + jj * 258 - 259] = sin__ * temp + cos__ * y[k 
			- 1];
/* L180: */
		y[k - 1] = sin__ * y[k - 1] - cos__ * temp;
	    }
/* L190: */
	}
L200:
	;
    }
    p[*n - 1] = efs[*n - 1] / nllcom_1.r[1];
    i = *n;
L210:
    --i;
    if (i <= 0) {
	goto L240;
    } else {
	goto L220;
    }
L220:
    temp = efs[i - 1];
    k = i + 1;
    ii = np2 - i;
    i__1 = *n;
    for (j = k; j <= i__1; ++j) {
	jj = np1 - j;
/* L230: */
	temp -= nllcom_1.r[ii + jj * 258 - 259] * p[j - 1];
    }
    jj = np1 - i;
    p[i - 1] = temp / nllcom_1.r[ii + jj * 258 - 259];
    goto L210;
L240:
    goto L270;
/*     ********** */
/*     SIDESTEP SECTION */
/*     ********** */
L250:
    ++(*jrst);
    ++nsst;
    if (nsst >= ixso) {
	goto L710;
    }
    if (*jrst > *n) {
	*jrst = 2;
    }
    *irst = 0;
/*     ********** */
/*     PRODUCTION OF A VECTOR ORTHOGONAL TO THE LAST P-VECTOR */
/*     ********** */
    work = *pn * (abs(p[0]) + *pn);
    temp = p[*jrst - 1];
    p[0] = temp * (p[0] + d_sign(pn, p));
    i__1 = *n;
    for (i = 2; i <= i__1; ++i) {
/* L260: */
	p[i - 1] = temp * p[i - 1];
    }
    p[*jrst - 1] -= work;
/*     ********** */
/*     COMPUTE NORM AND NORM-SQUARE OF THE P-VECTOR */
/*     ********** */
L270:
    pnlast = *pn;
    *pn = 0.;
    pn2 = 0.;
    i__1 = *n;
    for (i = 1; i <= i__1; ++i) {
	*pn += (d__1 = p[i - 1], abs(d__1));
/* L280: */
/* Computing 2nd power */
	d__1 = p[i - 1];
	pn2 += d__1 * d__1;
    }
    if (*pn < 1e-20) {
	s_wsfe(&io___61);
	e_wsfe();
	geout_();
	s_stop("", 0L);
    }
    if (pn2 < 1e-20) {
	pn2 = 1e-20;
    }
    *pn = sqrt(pn2);
    if (*alf > 1e20) {
	*alf = 1e20;
    }
    if (*icyc > 1) {
	*alf = *alf * 1e-20 * pnlast / *pn;
	if (*alf > 1e10) {
	    *alf = 1e10;
	}
	*alf *= 1e20;
    }
    ttmp = *alf * *pn;
    if (ttmp < 1e-4) {
	*alf = .001 / *pn;
    }
/*     ********** */
/*     PRINTING SECTION */
/*     ********** */
/* #      WRITE(6,501)TLEFT,ICYC,SSQ */
    i__1 = *n;
    for (i = 1; i <= i__1; ++i) {
	efs[i - 1] = x[i];
/* L290: */
    }
/*     ********** */
/*     PERFORM LINE-MINIMIZATION FROM POINT X IN DIRECTION P OR -P */
/*     ********** */
    ssqlst = *ssq;
    i__1 = *n;
    for (i = 1; i <= i__1; ++i) {
	efs[i - 1] = 0.;
/* L300: */
	nllcom_1.xlast[i - 1] = x[i];
    }
    locmin_(&m, &x[1], n, p, ssq, alf, efs, &ierr, &escf);
    if (ssqlst < *ssq) {
	if (ierr == 0) {
	    *ssq = ssqlst;
	}
	i__1 = *n;
	for (i = 1; i <= i__1; ++i) {
/* L310: */
	    x[i] = nllcom_1.xlast[i - 1];
	}
	*irst = nrst;
	*pn = pnlast;
	time2 = time1;
	time1 = second_();
	tcycle = time1 - time2;
	tleft -= tcycle;
	if (tleft > tcycle * 2) {
	    goto L100;
	}
	goto L670;
    }
    irepet = 0;
/*     ********** */
/*     PRODUCE THE VECTOR   R*P */
/*     ********** */
    i__1 = *n;
    for (i = 1; i <= i__1; ++i) {
	temp = 0.;
	i__2 = *n;
	for (j = i; j <= i__2; ++j) {
/* L320: */
	    temp += nllcom_1.r[i + j * 258 - 259] * p[j - 1];
	}
/* L330: */
	y[i - 1] = temp;
    }
/*     ********** */
/*     PRODUCE THE VECTOR ... */
/*                  Y  =    (EFS-EFSLST-ALF*Q*R*P)/(ALF*(NORMSQUARE(P)) */
/*     COMPUTE NORM OF THIS VECTOR AS WELL */
/*     ********** */
    work = *alf * pn2;
    yn = 0.;
    i__1 = m;
    for (i = 1; i <= i__1; ++i) {
	temp = 0.;
	i__2 = *n;
	for (j = 1; j <= i__2; ++j) {
/* L340: */
	    temp += nllcom_1.q[i + j * 258 - 259] * y[j - 1];
	}
	temp = efs[i - 1] - nllcom_1.efslst[i - 1] - *alf * temp;
	nllcom_1.efslst[i - 1] = efs[i - 1];
/* Computing 2nd power */
	d__1 = temp;
	yn += d__1 * d__1;
/* L350: */
	efs[i - 1] = temp / work;
    }
    yn = sqrt(yn) / work;
/*     ********** */
/*     THE BROYDEN UPDATE   NEW MATRIX = OLD MATRIX + Y*(P-TRANS) */
/*     HAS BEEN FORMED.  IT IS NOW NECESSARY TO UPDATE THE  QR DECOMP. */
/*     FIRST LET    Y = (Q-TRANS)*Y. */
/*     ********** */
    i__1 = m;
    for (i = 1; i <= i__1; ++i) {
	temp = 0.;
	i__2 = m;
	for (j = 1; j <= i__2; ++j) {
/* L360: */
	    temp += nllcom_1.q[j + i * 258 - 259] * efs[j - 1];
	}
/* L370: */
	y[i - 1] = temp;
    }
/*     ********** */
/*     REDUCE THE VECTOR Y TO A MULTIPLE OF THE FIRST UNIT VECTOR USING */
/*     A HOUSEHOLDER TRANSFORMATION FOR COMPONENTS N+1 THROUGH M AND */
/*     ELEMENTARY ROTATIONS FOR THE FIRST N+1 COMPONENTS.  APPLY ALL */
/*     TRANSFORMATIONS TRANSPOSED ON THE RIGHT TO THE MATRIX Q, AND */
/*     APPLY THE ROTATIONS ON THE LEFT TO THE MATRIX R. */
/*     THIS GIVES    (Q*(V-TRANS))*((V*R) + (V*Y)*(P-TRANS)),    WHERE */
/*     V IS THE COMPOSITE OF THE TRANSFORMATIONS.  THE MATRIX */
/*     ((V*R) + (V*Y)*(P-TRANS))    IS UPPER HESSENBERG. */
/*     ********** */
    if (m <= np1) {
	goto L430;
    }

/* THE NEXT THREE LINES WERE INSERTED TO TRY TO GET ROUND OVERFLOW BUGS. 
*/

    const__ = 1e-12;
    i__1 = m;
    for (i = np1; i <= i__1; ++i) {
/* L380: */
/* Computing MAX */
	d__2 = (d__1 = y[np1 - 1], abs(d__1));
	const__ = max(d__2,const__);
    }
    ytail = 0.;
    i__1 = m;
    for (i = np1; i <= i__1; ++i) {
/* L390: */
/* Computing 2nd power */
	d__1 = y[i - 1] / const__;
	ytail += d__1 * d__1;
    }
    ytail = sqrt(ytail) * const__;
    bet = 1e25 / ytail / (ytail + (d__1 = y[np1 - 1], abs(d__1)));
    d__2 = ytail + (d__1 = y[np1 - 1], abs(d__1));
    y[np1 - 1] = d_sign(&d__2, &y[np1 - 1]);
    i__1 = m;
    for (i = 1; i <= i__1; ++i) {
	tmp = 0.;
	i__2 = m;
	for (j = np1; j <= i__2; ++j) {
/* L400: */
	    tmp += nllcom_1.q[i + j * 258 - 259] * y[j - 1] * 1e-25;
	}
	tmp = bet * tmp;
	i__2 = m;
	for (j = np1; j <= i__2; ++j) {
/* L410: */
	    nllcom_1.q[i + j * 258 - 259] -= tmp * y[j - 1];
	}
/* L420: */
    }
    y[np1 - 1] = ytail;
    i = np1;
    goto L440;
L430:
    i = m;
L440:
L450:
    j = i;
    --i;
    if (i <= 0) {
	goto L520;
    } else {
	goto L460;
    }
L460:
    if (y[j - 1] != 0.) {
	goto L470;
    } else {
	goto L450;
    }
L470:
    if ((d__1 = y[i - 1], abs(d__1)) < (d__2 = y[j - 1], abs(d__2))) {
	goto L480;
    }
/* Computing 2nd power */
    d__2 = y[j - 1] / y[i - 1];
    temp = (d__1 = y[i - 1], abs(d__1)) * sqrt(d__2 * d__2 + 1.);
    goto L490;
L480:
/* Computing 2nd power */
    d__2 = y[i - 1] / y[j - 1];
    temp = (d__1 = y[j - 1], abs(d__1)) * sqrt(d__2 * d__2 + 1.);
L490:
    cos__ = y[i - 1] / temp;
    sin__ = y[j - 1] / temp;
    y[i - 1] = temp;
    i__1 = m;
    for (k = 1; k <= i__1; ++k) {
	temp = cos__ * nllcom_1.q[k + i * 258 - 259] + sin__ * nllcom_1.q[k + 
		j * 258 - 259];
	work = -sin__ * nllcom_1.q[k + i * 258 - 259] + cos__ * nllcom_1.q[k 
		+ j * 258 - 259];
	nllcom_1.q[k + i * 258 - 259] = temp;
/* L500: */
	nllcom_1.q[k + j * 258 - 259] = work;
    }
    if (i > *n) {
	goto L450;
    }
    nllcom_1.r[j + i * 258 - 259] = -sin__ * nllcom_1.r[i + i * 258 - 259];
    nllcom_1.r[i + i * 258 - 259] = cos__ * nllcom_1.r[i + i * 258 - 259];
    if (j > *n) {
	goto L450;
    }
    i__1 = *n;
    for (k = j; k <= i__1; ++k) {
	temp = cos__ * nllcom_1.r[i + k * 258 - 259] + sin__ * nllcom_1.r[j + 
		k * 258 - 259];
	work = -sin__ * nllcom_1.r[i + k * 258 - 259] + cos__ * nllcom_1.r[j 
		+ k * 258 - 259];
	nllcom_1.r[i + k * 258 - 259] = temp;
/* L510: */
	nllcom_1.r[j + k * 258 - 259] = work;
    }
    goto L450;
L520:
/*     ********** */
/*     REDUCE THE UPPER-HESSENBERG MATRIX TO UPPER-TRIANGULAR FORM */
/*     USING ELEMENTARY ROTATIONS.  APPLY THE SAME ROTATIONS, TRANSPOSED, 
*/
/*     ON THE RIGHT TO THE MATRIX  Q. */
/*     ********** */
    i__1 = *n;
    for (k = 1; k <= i__1; ++k) {
/* L530: */
	nllcom_1.r[k * 258 - 258] += yn * p[k - 1];
    }
    jend = np1;
    if (m == *n) {
	jend = *n;
    }
    i__1 = jend;
    for (j = 2; j <= i__1; ++j) {
	i = j - 1;
	if (nllcom_1.r[j + i * 258 - 259] != 0.) {
	    goto L540;
	} else {
	    goto L600;
	}
L540:
	if ((d__1 = nllcom_1.r[i + i * 258 - 259], abs(d__1)) < (d__2 = 
		nllcom_1.r[j + i * 258 - 259], abs(d__2))) {
	    goto L550;
	}
/* Computing 2nd power */
	d__2 = nllcom_1.r[j + i * 258 - 259] / nllcom_1.r[i + i * 258 - 259];
	temp = (d__1 = nllcom_1.r[i + i * 258 - 259], abs(d__1)) * sqrt(d__2 *
		 d__2 + 1.);
	goto L560;
L550:
/* Computing 2nd power */
	d__2 = nllcom_1.r[i + i * 258 - 259] / nllcom_1.r[j + i * 258 - 259];
	temp = (d__1 = nllcom_1.r[j + i * 258 - 259], abs(d__1)) * sqrt(d__2 *
		 d__2 + 1.);
L560:
	cos__ = nllcom_1.r[i + i * 258 - 259] / temp;
	sin__ = nllcom_1.r[j + i * 258 - 259] / temp;
	nllcom_1.r[i + i * 258 - 259] = temp;
	if (j > *n) {
	    goto L580;
	}
	i__2 = *n;
	for (k = j; k <= i__2; ++k) {
	    temp = cos__ * nllcom_1.r[i + k * 258 - 259] + sin__ * nllcom_1.r[
		    j + k * 258 - 259];
	    work = -sin__ * nllcom_1.r[i + k * 258 - 259] + cos__ * 
		    nllcom_1.r[j + k * 258 - 259];
	    nllcom_1.r[i + k * 258 - 259] = temp;
/* L570: */
	    nllcom_1.r[j + k * 258 - 259] = work;
	}
L580:
	i__2 = m;
	for (k = 1; k <= i__2; ++k) {
	    temp = cos__ * nllcom_1.q[k + i * 258 - 259] + sin__ * nllcom_1.q[
		    k + j * 258 - 259];
	    work = -sin__ * nllcom_1.q[k + i * 258 - 259] + cos__ * 
		    nllcom_1.q[k + j * 258 - 259];
	    nllcom_1.q[k + i * 258 - 259] = temp;
/* L590: */
	    nllcom_1.q[k + j * 258 - 259] = work;
	}
L600:
	;
    }
/*     ********** */
/*     CHECK THE STOPPING CRITERIA */
/*     ********** */
    temp = 0.;
    i__1 = *n;
    for (i = 1; i <= i__1; ++i) {
/* L610: */
/* Computing 2nd power */
	d__1 = x[i];
	temp += d__1 * d__1;
    }
    tolx = tols1 * sqrt(temp) + tols2;
    if (sqrt(*alf * pn2) <= tolx) {
	goto L690;
    }
    if (*ssq >= *n * 2.) {
	goto L630;
    }
    i__1 = *n;
    for (i = 1; i <= i__1; ++i) {
/* ***** */
/*     The stopping criterion is that no individual gradient be */
/*         greater than TOL2 */
/* ***** */
	if ((d__1 = nllcom_1.efslst[i - 1], abs(d__1)) >= tol2) {
	    goto L630;
	}
/* L620: */
    }
/* #      WRITE(6,730) SSQ */
    goto L700;
L630:
    if (*icyc >= mxcycl) {
	mesage_1.iflepo = 12;
	goto L920;
    }
    time2 = time1;
    time1 = second_();
    tcycle = time1 - time2;
    tleft -= tcycle;
    if (resfil) {
	s_wsfe(&io___74);
	do_fio(&c__1, (char *)&tleft, (ftnlen)sizeof(doublereal));
	d__1 = sqrt(*ssq);
	do_fio(&c__1, (char *)&d__1, (ftnlen)sizeof(doublereal));
	do_fio(&c__1, (char *)&escf, (ftnlen)sizeof(doublereal));
	e_wsfe();
	resfil = FALSE_;
    } else {
	s_wsfe(&io___75);
	do_fio(&c__1, (char *)&(*icyc), (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&tcycle, (ftnlen)sizeof(doublereal));
	do_fio(&c__1, (char *)&tleft, (ftnlen)sizeof(doublereal));
	d__1 = sqrt(*ssq);
	do_fio(&c__1, (char *)&d__1, (ftnlen)sizeof(doublereal));
	do_fio(&c__1, (char *)&escf, (ftnlen)sizeof(doublereal));
	e_wsfe();
    }
    if (tlast - tleft > tdump) {
	tlast = tleft;
	resfil = TRUE_;
	i__1 = *n;
	for (i = 1; i <= i__1; ++i) {
/* L660: */
	    nllcom_1.xlast[i - 1] = x[i];
	}
	nllcom_1.iiium[0] = numscf_1.nscf;
	parsav_(&c__2, n, &m);
    }
    if (tleft > tcycle * 2) {
	goto L100;
    }
L670:
    nllcom_1.iiium[4] = nllsqi_1.ncount;
    i__1 = *n;
    for (i = 1; i <= i__1; ++i) {
/* L680: */
	nllcom_1.xlast[i - 1] = x[i];
    }
    nllcom_1.iiium[0] = numscf_1.nscf;
    parsav_(&c__1, n, &m);
    s_stop("", 0L);
L690:
    s_wsfe(&io___76);
    do_fio(&c__1, (char *)&nllsqi_1.ncount, (ftnlen)sizeof(integer));
    e_wsfe();
    goto L920;
L700:
    s_wsfe(&io___77);
    do_fio(&c__1, (char *)&nllsqi_1.ncount, (ftnlen)sizeof(integer));
    e_wsfe();
    goto L920;
L710:
    s_wsfe(&io___78);
    do_fio(&c__1, (char *)&ixso, (ftnlen)sizeof(integer));
    e_wsfe();
    goto L920;
/* #  730 FORMAT(1H ,'FINAL GRADIENT =',F15.7) */
/* L740: */
/* L750: */
/* L760: */
/* L770: */
/* L780: */
/* L790: */
/* L800: */
/* L830: */
/* L840: */
/* L850: */
/* L860: */
/* L870: */
/* L880: */
/* L890: */
/* L900: */
/* L910: */
L920:
    last_1.last = 1;
    compfg_(&x[1], &c_true, &escf, &c_true, nllcom_1.efslst, &c_true);
    return 0;
} /* nllsq_ */

#undef ssq
#undef alf
#undef pn
#undef jrst
#undef irst
#undef icyc