/* write.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 {
    char koment[80], title[80];
} titles_;

#define titles_1 titles_

struct {
    char elemnt[214];
} elemts_;

#define elemts_1 elemts_

struct {
    doublereal geo[258]	/* was [3][86] */;
} geom_;

#define geom_1 geom_

struct {
    integer natoms, labels[86], na[86], nb[86], nc[86];
} geokst_;

#define geokst_1 geokst_

struct {
    doublereal h[23220];
} hmatrx_;

#define hmatrx_1 hmatrx_

struct {
    doublereal f[23220], fb[23220];
} fokmat_;

#define fokmat_1 fokmat_

struct {
    doublereal c[46225], eigs[215], cbeta[46225], eigb[215];
} vector_;

#define vector_1 vector_

struct {
    doublereal p[23220], pa[23220], pb[23220];
} densty_;

#define densty_1 densty_

struct {
    integer ndep, locpar[258], idepfn[258], locdep[258];
} geosym_;

#define geosym_1 geosym_

struct {
    integer latom, lparam;
    doublereal react[200];
} path_;

#define path_1 path_

struct {
    integer nscf;
} numscf_;

#define numscf_1 numscf_

struct {
    real wj[219386], wk[219386];
} wmatrx_;

#define wmatrx_1 wmatrx_

struct {
    doublereal atheat;
} atheat_;

#define atheat_1 atheat_

struct {
    doublereal core[107];
} core_;

#define core_1 core_

struct {
    doublereal rxyz[23220], xdumy[43344];
} scrach_;

#define scrach_1 scrach_

struct {
    doublereal engyci[3], vectci[9], eci[6];
} cimats_;

#define cimats_1 cimats_

struct {
    integer iflepo, iiter;
} mesage_;

#define mesage_1 mesage_

struct {
    doublereal atmass[86];
} atmass_;

#define atmass_1 atmass_

struct {
    doublereal enuclr;
} enuclr_;

#define enuclr_1 enuclr_

struct {
    doublereal elect;
} elect_;

#define elect_1 elect_

struct {
    doublereal dxyz[6966]	/* was [3][2322] */;
} xyzgra_;

#define xyzgra_1 xyzgra_

struct {
    doublereal grad[258], gnorm;
} gradnt_;

#define gradnt_1 gradnt_

union {
    struct {
	integer numat, nat[86], nfirst[86], nmidle[86], nlast[86], norbs, 
		nelecs, nalpha, nbeta, nclose, nopen, ndumy;
	doublereal fract;
    } _1;
    struct {
	integer numat, nat[86], nfirst[86], nmidle[86], nlast[86];
    } _2;
} molkst_;

#define molkst_1 (molkst_._1)
#define molkst_2 (molkst_._2)

struct {
    integer nvar, loc[516]	/* was [2][258] */, idumy;
    doublereal xparam[258];
} geovar_;

#define geovar_1 geovar_

/* Table of constant values */

static integer c__1 = 1;
static doublereal c_b36 = 5.;
static integer c__6 = 6;
static doublereal c_b85 = .4999;
static doublereal c_b171 = 2.;
static doublereal c_b173 = 1.;
static integer c__3 = 3;
static logical c_true = TRUE_;
static integer c__0 = 0;
static integer c__2 = 2;

/* Subroutine */ int write_(time0, funct)
doublereal *time0, *funct;
{
    /* Initialized data */

    static char type[11*3+1] = "BOND       ANGLE      DIHEDRAL   ";
    static char calcn[5*2+1] = "     ALPHA";
    static char numbrs[1*11+1] = "0123456789 ";
    static logical first = TRUE_;
    static char flepo[58*13+1] = " 1SCF WAS SPECIFIED, SO BFGS WAS NOT USED \
                 GRADIENTS WERE INITIALLY ACCEPTABLY SMALL                 H\
ERBERTS TEST WAS SATISFIED IN BFGS                       THE LINE MINIMIZATI\
ON FAILED TWICE IN A ROW.   TAKE CARE! BFGS FAILED DUE TO COUNTS EXCEEDED. T\
AKE CARE!            PETERS TEST WAS SATISFIED IN BFGS OPTIMIZATION         \
   THIS MESSAGE SHOULD NEVER APPEAR, CONSULT A PROGRAMMER!!  GRADIENT TEST N\
OT PASSED, BUT FURTHER WORK NOT JUSTIFIED  A FAILURE HAS OCCURRED, TREAT RES\
ULTS WITH CAUTION!!      GEOMETRY OPTIMIZED USING NLLSQ. GRADIENT NORM MINIM\
IZED   GEOMETRY OPTIMIZED USING POWSQ. GRADIENT NORM MINIMIZED   CYCLES EXCE\
EDED, GRADIENT NOT FULLY MINIMIZED IN NLLSQ                                 \
                            ";
    static char iter[58*2+1] = " SCF FIELD WAS ACHIEVED                     \
                ++++----**** FAILED TO ACHIEVE SCF. ****----++++        ";

    /* Format strings */
    static char fmt_210[] = "(//,1x,17(a2,a1,a1))";

    /* System generated locals */
    integer i__1, i__2;
    doublereal d__1, d__2, d__3;
    cilist ci__1;
    olist o__1;
    cllist cl__1;
    alist al__1;

    /* Builtin functions */
    /* Subroutine */ int s_copy();
    integer i_indx();
    double sqrt();
    integer s_wsfe(), e_wsfe(), do_fio();
    /* Subroutine */ int s_stop();
    double d_sign(), d_int();
    integer f_open(), f_rew(), s_wsue(), do_uio(), e_wsue(), f_clos();

    /* Local variables */
    extern doublereal meci_();
    integer loc11, loc21, ksec, iuhf;
    extern doublereal spcg_();
    integer ncis, ivar, mvar, nopn, nmos;
    doublereal dumy[3], sumw;
    integer i, j, k, l;
    extern doublereal reada_();
    doublereal q[215];
    extern /* Subroutine */ int fdate_();
    char idate[24];
#define w ((doublereal *)&wmatrx_1)
    doublereal x;
    extern /* Subroutine */ int chrge_(), local_();
    doublereal coord[258]	/* was [3][86] */;
    extern /* Subroutine */ int bonds_(), deriv_();
    doublereal xiiii;
    extern /* Subroutine */ int geout_();
    char gtype[13];
    integer kfrst;
    doublereal q2[215];
    logical ci;
    integer ik;
    doublereal pi, degree;
    integer it, linear;
    logical excitd;
    char ielemt[2*20];
    integer nelemt[107];
    logical singlt, prtgra;
    char caltyp[7], grtype[14];
    doublereal xreact, gcoord;
    logical triplt;
    doublereal eionis;
    extern doublereal second_();
    extern /* Subroutine */ int timout_(), gmetry_();
    integer kchrge;
    doublereal xi;
    extern /* Subroutine */ int vecprt_(), matout_();
    extern doublereal dipole_();
    integer nfilld;
    extern /* Subroutine */ int mpcsyb_(), denrot_();
    doublereal sz, ss2;
    char spntyp[7];
    extern /* Subroutine */ int molval_(), enpart_();
    integer ichfor;
    extern /* Subroutine */ int getarg_(), mullik_(), mpcpop_();
    integer iwrite, na1;
    extern /* Subroutine */ int xyzint_(), symtry_();
    doublereal dip;
    integer igo;
    logical uhf;
    extern doublereal dot_();
    doublereal tim, sum;
    integer nzs, iel1[107], iel2[107];
    logical xyz;

    /* Fortran I/O blocks */
    static cilist io___59 = { 0, 10, 0, 0, 0 };
    static cilist io___60 = { 0, 10, 0, 0, 0 };
    static cilist io___64 = { 0, 0, 0, fmt_210, 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 */

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

/*   WRITE PRINTS OUT MOST OF THE RESULTS. */
/*         IT SHOULD NOT ALTER ANY PARAMETERS, SO THAT IT CAN BE CALLED */
/*         AT ANY CONVENIENT TIME. */

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

/* SUMMARY OF RESULTS (NOTE: THIS IS IN A SUBROUTINE SO IT */
/*          CAN BE USED BY THE PATH OPTION) */
    pi = 3.141592653589;
    s_copy(idate, " ", 24L, 1L);
    if (mesage_1.iflepo == 0) {
	mesage_1.iflepo = 7;
    }
/* Computing MIN */
    i__1 = i_indx(keywrd_1.keywrd, "UHF", 80L, 3L);
    iuhf = min(i__1,1) + 1;
    prtgra = i_indx(keywrd_1.keywrd, "GRADIENTS", 80L, 9L) != 0;
    linear = molkst_1.norbs * (molkst_1.norbs + 1) / 2;
    xyz = i_indx(keywrd_1.keywrd, " XYZ", 80L, 4L) != 0;
    singlt = i_indx(keywrd_1.keywrd, "SINGLET", 80L, 7L) != 0;
    triplt = i_indx(keywrd_1.keywrd, "TRIPLET", 80L, 7L) != 0;
    excitd = i_indx(keywrd_1.keywrd, "EXCITED", 80L, 7L) != 0;
    s_copy(spntyp, "GROUND ", 7L, 7L);
    if (singlt) {
	s_copy(spntyp, "SINGLET", 7L, 7L);
    }
    if (triplt) {
	s_copy(spntyp, "TRIPLET", 7L, 7L);
    }
    if (excitd) {
	s_copy(spntyp, "EXCITED", 7L, 7L);
    }
    ci = i_indx(keywrd_1.keywrd, "C.I.", 80L, 4L) != 0;
    if (i_indx(keywrd_1.keywrd, "MINDO", 80L, 5L) != 0) {
	s_copy(caltyp, "MINDO/3", 7L, 7L);
    } else if (i_indx(keywrd_1.keywrd, "AM1", 80L, 3L) != 0) {
	s_copy(caltyp, "  AM1  ", 7L, 7L);
    } else if (i_indx(keywrd_1.keywrd, "PM3", 80L, 3L) != 0) {
	s_copy(caltyp, "  PM3  ", 7L, 7L);
    } else {
	s_copy(caltyp, " MNDO  ", 7L, 7L);
    }
    uhf = iuhf == 2;
    fdate_(idate, 24L);
    degree = 57.29577951;
    if (geokst_1.na[0] == 99) {
	degree = 1.;
	s_copy(type, "           ", 11L, 11L);
	s_copy(type + 11, "           ", 11L, 11L);
	s_copy(type + 22, "           ", 11L, 11L);
    }
    gradnt_1.gnorm = 0.;
    if (geovar_1.nvar != 0) {
	gradnt_1.gnorm = sqrt(dot_(gradnt_1.grad, gradnt_1.grad, &
		geovar_1.nvar));
    }
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(/,' ----',15('-----'))";
    s_wsfe(&ci__1);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, keywrd_1.keywrd, 80L);
    do_fio(&c__1, titles_1.koment, 80L);
    do_fio(&c__1, titles_1.title, 80L);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(//4X,A58)";
    s_wsfe(&ci__1);
    do_fio(&c__1, flepo + (mesage_1.iflepo - 1) * 58, 58L);
    e_wsfe();
    mesage_1.iiter = max(1,mesage_1.iiter);
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(4X,A58)";
    s_wsfe(&ci__1);
    do_fio(&c__1, iter + (mesage_1.iiter - 1) * 58, 58L);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(//30X,A7,'  CALCULATION')";
    s_wsfe(&ci__1);
    do_fio(&c__1, caltyp, 7L);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(55X,'VERSION ',F5.2)";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&c_b36, (ftnlen)sizeof(doublereal));
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(55X,A24)";
    s_wsfe(&ci__1);
    do_fio(&c__1, idate, 24L);
    e_wsfe();
    if (mesage_1.iiter == 2) {

/*   RESULTS ARE MEANINGLESS. DON'T PRINT ANYTHING! */

	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(//,' FOR SOME REASON THE SCF CALCULATION FAILED.',/,\
' THE RESULTS WOULD BE MEANINGLESS, SO WILL NOT BE PRINTED.')";
	s_wsfe(&ci__1);
	e_wsfe();
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(' TRY TO FIND THE REASON FOR THE FAILURE BY USING \
','\"PL\".',/,' CHECK YOUR GEOMETRY AND ALSO TRY USING SHIFT OR PULAY. ')";
	s_wsfe(&ci__1);
	e_wsfe();
	geout_();
	s_stop("", 0L);
    }
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(////10X,'FINAL HEAT OF FORMATION =',F17.5,' KCAL')";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&(*funct), (ftnlen)sizeof(doublereal));
    e_wsfe();
    if (path_1.latom == 0) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(/)";
	s_wsfe(&ci__1);
	e_wsfe();
    }
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(    10X,'TOTAL ENERGY            =',F17.5,' EV')";
    s_wsfe(&ci__1);
    d__1 = elect_1.elect + enuclr_1.enuclr;
    do_fio(&c__1, (char *)&d__1, (ftnlen)sizeof(doublereal));
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(    10X,'ELECTRONIC ENERGY       =',F17.5,' EV')";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&elect_1.elect, (ftnlen)sizeof(doublereal));
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(    10X,'CORE-CORE REPULSION     =',F17.5,' EV')";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&enuclr_1.enuclr, (ftnlen)sizeof(doublereal));
    e_wsfe();
    if (path_1.latom == 0) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(1X)";
	s_wsfe(&ci__1);
	e_wsfe();
    }
    prtgra = prtgra || gradnt_1.gnorm > 2.;
    if (prtgra) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(    10X,'GRADIENT NORM           =',F17.5)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&gradnt_1.gnorm, (ftnlen)sizeof(doublereal));
	e_wsfe();
    }
    if (path_1.latom != 0) {

/*   WE NEED TO CALCULATE THE REACTION COORDINATE GRADIENT. */

	mvar = geovar_1.nvar;
	loc11 = geovar_1.loc[0];
	loc21 = geovar_1.loc[1];
	geovar_1.nvar = 1;
	geovar_1.loc[0] = path_1.latom;
	geovar_1.loc[1] = path_1.lparam;
	xreact = geom_1.geo[path_1.lparam + path_1.latom * 3 - 4];
	deriv_(geom_1.geo, &gcoord);
	geovar_1.nvar = mvar;
	geovar_1.loc[0] = loc11;
	geovar_1.loc[1] = loc21;
	s_copy(grtype, " KCAL/ANGSTROM", 14L, 14L);
	if (path_1.lparam == 1) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(    10X,'FOR REACTION COORDINATE =',F17.5       \
 ,' ANGSTROMS')";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, (char *)&xreact, (ftnlen)sizeof(doublereal));
	    e_wsfe();
	} else {
	    if (geokst_1.na[0] != 99) {
		s_copy(grtype, " KCAL/RADIAN  ", 14L, 14L);
	    }
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(    10X,'FOR REACTION COORDINATE =',F17.5       \
 ,' DEGREES')";
	    s_wsfe(&ci__1);
	    d__1 = xreact * degree;
	    do_fio(&c__1, (char *)&d__1, (ftnlen)sizeof(doublereal));
	    e_wsfe();
	}
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(    10X,'REACTION GRADIENT       =',F17.5,A14    )";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&gcoord, (ftnlen)sizeof(doublereal));
	do_fio(&c__1, grtype, 14L);
	e_wsfe();
    }
    if (molkst_1.nalpha > 0) {
/* Computing MAX */
	d__1 = vector_1.eigs[molkst_1.nalpha - 1], d__2 = vector_1.eigb[
		molkst_1.nbeta - 1];
	eionis = -max(d__1,d__2);
    } else if (molkst_1.nelecs == 1) {
	eionis = -vector_1.eigs[0];
    } else if (molkst_1.nelecs > 1) {
/* Computing MAX */
	d__1 = vector_1.eigs[molkst_1.nclose - 1], d__2 = vector_1.eigs[
		molkst_1.nopen - 1];
	eionis = -max(d__1,d__2);
    } else {
	eionis = 0.;
    }
    nopn = molkst_1.nopen - molkst_1.nclose;
/*   CORRECTION TO I.P. OF DOUBLETS */
    if (nopn == 1) {
	i = molkst_1.nclose * molkst_1.norbs + 1;
	xiiii = spcg_(&vector_1.c[i - 1], &vector_1.c[i - 1], &vector_1.c[i - 
		1], &vector_1.c[i - 1], w, wmatrx_1.wj);
	eionis += xiiii * .5;
    }
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(       10X,'IONIZATION POTENTIAL    =',F17.5)";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&eionis, (ftnlen)sizeof(doublereal));
    e_wsfe();
    if (uhf) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(      10X,'NO. OF ALPHA ELECTRONS  =',I11)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&molkst_1.nalpha, (ftnlen)sizeof(integer));
	e_wsfe();
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(      10X,'NO. OF BETA  ELECTRONS  =',I11)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&molkst_1.nbeta, (ftnlen)sizeof(integer));
	e_wsfe();
    } else {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(      10X,'NO. OF FILLED LEVELS    =',I11)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&molkst_1.nclose, (ftnlen)sizeof(integer));
	e_wsfe();
	if (nopn != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(   10X,'AND NO. OF OPEN LEVELS  =',I11)";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, (char *)&nopn, (ftnlen)sizeof(integer));
	    e_wsfe();
	}
    }
    sumw = 0.;
    i__1 = molkst_1.numat;
    for (i = 1; i <= i__1; ++i) {
/* L10: */
	sumw += atmass_1.atmass[i - 1];
    }
    if (sumw > .1) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(    10X,'MOLECULAR WEIGHT        =',F11.3)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&sumw, (ftnlen)sizeof(doublereal));
	e_wsfe();
    }
    if (path_1.latom == 0) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(/)";
	s_wsfe(&ci__1);
	e_wsfe();
    }
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(10X,'SCF CALCULATIONS  =   ',I14 )";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&numscf_1.nscf, (ftnlen)sizeof(integer));
    e_wsfe();
    tim = second_() - *time0;
    i = (integer) (tim * 1e-6);
    tim -= i * 1000000;
    timout_(&c__6, &tim);
    if (geosym_1.ndep != 0) {
	symtry_();
    }
    if (geokst_1.na[0] != 99) {
	i__1 = geokst_1.natoms;
	for (j = 1; j <= i__1; ++j) {
	    for (i = 1; i <= 3; ++i) {
		x = geom_1.geo[i + j * 3 - 4];
		switch ((int)i) {
		    case 1:  goto L40;
		    case 2:  goto L30;
		    case 3:  goto L20;
		}
L20:
		d__1 = x / (pi * 2.) + d_sign(&c_b85, &x) - 1e-4;
		x -= d_int(&d__1) * pi * 2.;
		geom_1.geo[j * 3 - 1] = x;
		goto L40;
L30:
		d__1 = x / (pi * 2.);
		x -= d_int(&d__1) * pi * 2.;
		if (x < 0.) {
		    x += pi * 2.;
		}
		if (x > pi) {
		    geom_1.geo[j * 3 - 1] += pi;
		    x = pi * 2. - x;
		}
		geom_1.geo[j * 3 - 2] = x;
L40:
/* L50: */
		;
	    }
	}
    }
    i__1 = geovar_1.nvar;
    for (i = 1; i <= i__1; ++i) {
/* L60: */
	geovar_1.xparam[i - 1] = geom_1.geo[geovar_1.loc[(i << 1) - 1] + 
		geovar_1.loc[(i << 1) - 2] * 3 - 4];
    }
    gmetry_(geom_1.geo, coord);
    if (prtgra) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(///7X,'FINAL  POINT  AND  DERIVATIVES',/)";
	s_wsfe(&ci__1);
	e_wsfe();
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "('   PARAMETER     ATOM    TYPE  '    ,'          VAL\
UE       GRADIENT')";
	s_wsfe(&ci__1);
	e_wsfe();
    }
    sum = .5;
    i__1 = molkst_1.numat;
    for (i = 1; i <= i__1; ++i) {
/* L70: */
	sum += core_1.core[molkst_1.nat[i - 1] - 1];
    }
    i = (integer) sum;
    kchrge = i - molkst_1.nclose - molkst_1.nopen - molkst_1.nalpha - 
	    molkst_1.nbeta;

/*    WRITE OUT THE GEOMETRIC VARIABLES */

    if (prtgra) {
	i__1 = geovar_1.nvar;
	for (i = 1; i <= i__1; ++i) {
	    j = geovar_1.loc[(i << 1) - 1];
	    k = geovar_1.loc[(i << 1) - 2];
	    l = geokst_1.labels[k - 1];
	    xi = geovar_1.xparam[i - 1];
	    if (j != 1) {
		xi *= degree;
	    }
	    if (j == 1 || geokst_1.na[0] == 99) {
		s_copy(gtype, "KCAL/ANGSTROM", 13L, 13L);
	    } else {
		s_copy(gtype, "KCAL/RADIAN  ", 13L, 13L);
	    }
/* L80: */
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(I7,I11,1X,A2,4X,A11,F13.6,F13.6,2X,A13)";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
	    do_fio(&c__1, (char *)&k, (ftnlen)sizeof(integer));
	    do_fio(&c__1, elemts_1.elemnt + (l - 1 << 1), 2L);
	    do_fio(&c__1, type + (j - 1) * 11, 11L);
	    do_fio(&c__1, (char *)&xi, (ftnlen)sizeof(doublereal));
	    do_fio(&c__1, (char *)&gradnt_1.grad[i - 1], (ftnlen)sizeof(
		    doublereal));
	    do_fio(&c__1, gtype, 13L);
	    e_wsfe();
	}
    }

/*     WRITE OUT THE GEOMETRY */

    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(///)";
    s_wsfe(&ci__1);
    e_wsfe();
    geout_();
    if (i_indx(keywrd_1.keywrd, "NOINTER", 80L, 7L) == 0) {

/*   WRITE OUT THE INTERATOMIC DISTANCES */

	l = 0;
	i__1 = molkst_1.numat;
	for (i = 1; i <= i__1; ++i) {
	    i__2 = i;
	    for (j = 1; j <= i__2; ++j) {
		++l;
/* L90: */
/* Computing 2nd power */
		d__1 = coord[i * 3 - 3] - coord[j * 3 - 3];
/* Computing 2nd power */
		d__2 = coord[i * 3 - 2] - coord[j * 3 - 2];
/* Computing 2nd power */
		d__3 = coord[i * 3 - 1] - coord[j * 3 - 1];
		scrach_1.rxyz[l - 1] = sqrt(d__1 * d__1 + d__2 * d__2 + d__3 *
			 d__3);
	    }
	}
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(//10X,'  INTERATOMIC DISTANCES')";
	s_wsfe(&ci__1);
	e_wsfe();
	vecprt_(scrach_1.rxyz, &molkst_1.numat);
    }
    i__2 = molkst_1.norbs;
    for (i = 1; i <= i__2; ++i) {
/* L100: */
	if (vector_1.eigs[i - 1] < -999. || vector_1.eigs[i - 1] > 1e3) {
	    vector_1.eigs[i - 1] = 0.;
	}
    }
    i__2 = molkst_1.norbs;
    for (i = 1; i <= i__2; ++i) {
/* L110: */
	if (vector_1.eigb[i - 1] < -999. || vector_1.eigb[i - 1] > 1e3) {
	    vector_1.eigs[i - 1] = 0.;
	}
    }
    if (molkst_1.norbs > 0) {
	if (i_indx(keywrd_1.keywrd, "VECT", 80L, 4L) != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(//10X,A5,' EIGENVECTORS  ')";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, calcn + (iuhf - 1) * 5, 5L);
	    e_wsfe();
	    matout_(vector_1.c, vector_1.eigs, &molkst_1.norbs, &
		    molkst_1.norbs, &molkst_1.norbs);
	    if (uhf) {
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(//10X,' BETA EIGENVECTORS  ')";
		s_wsfe(&ci__1);
		e_wsfe();
		matout_(vector_1.cbeta, vector_1.eigb, &molkst_1.norbs, &
			molkst_1.norbs, &molkst_1.norbs);
	    }
	} else {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(//10X,A5,'   EIGENVALUES',/)";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, calcn + (iuhf - 1) * 5, 5L);
	    e_wsfe();
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(8F10.5)";
	    s_wsfe(&ci__1);
	    i__2 = molkst_1.norbs;
	    for (i = 1; i <= i__2; ++i) {
		do_fio(&c__1, (char *)&vector_1.eigs[i - 1], (ftnlen)sizeof(
			doublereal));
	    }
	    e_wsfe();
	    if (uhf) {
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(//10X,' BETA EIGENVALUES ')";
		s_wsfe(&ci__1);
		e_wsfe();
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(8F10.5)";
		s_wsfe(&ci__1);
		i__2 = molkst_1.norbs;
		for (i = 1; i <= i__2; ++i) {
		    do_fio(&c__1, (char *)&vector_1.eigb[i - 1], (ftnlen)
			    sizeof(doublereal));
		}
		e_wsfe();
	    }
	}
    }
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(//13X,' NET ATOMIC CHARGES AND DIPOLE ','CONTRIBUTIONS',\
/)";
    s_wsfe(&ci__1);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(8X,' ATOM NO.   TYPE          CHARGE        ATOM','  ELE\
CTRON DENSITY')";
    s_wsfe(&ci__1);
    e_wsfe();
    chrge_(densty_1.p, q);
    i__2 = molkst_1.numat;
    for (i = 1; i <= i__2; ++i) {
	l = molkst_1.nat[i - 1];
	q2[i - 1] = core_1.core[l - 1] - q[i - 1];
/* L120: */
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(I12,9X,A2,4X,F13.4,F16.4)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
	do_fio(&c__1, elemts_1.elemnt + (l - 1 << 1), 2L);
	do_fio(&c__1, (char *)&q2[i - 1], (ftnlen)sizeof(doublereal));
	do_fio(&c__1, (char *)&q[i - 1], (ftnlen)sizeof(doublereal));
	e_wsfe();
    }
    dip = dipole_(densty_1.p, q2, coord, dumy);
    if (i_indx(keywrd_1.keywrd, "NOXYZ", 80L, 5L) == 0) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(//10X,'CARTESIAN COORDINATES ',/)";
	s_wsfe(&ci__1);
	e_wsfe();
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(4X,'NO.',7X,'ATOM',15X,'X',  9X,'Y',9X,'Z',/)";
	s_wsfe(&ci__1);
	e_wsfe();
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(I6,8X,A2,14X,3F10.4)";
	s_wsfe(&ci__1);
	i__2 = molkst_1.numat;
	for (i = 1; i <= i__2; ++i) {
	    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
	    do_fio(&c__1, elemts_1.elemnt + (molkst_1.nat[i - 1] - 1 << 1), 
		    2L);
	    for (j = 1; j <= 3; ++j) {
		do_fio(&c__1, (char *)&coord[j + i * 3 - 4], (ftnlen)sizeof(
			doublereal));
	    }
	}
	e_wsfe();

/*  The following call was added by Pets 5/29/85 to write an output fi
le */
/*      readable by sybyl */
/* Computing MAX */
	i__2 = max(molkst_1.nclose,molkst_1.nalpha);
	nfilld = max(i__2,molkst_1.nbeta);
	mpcsyb_(&molkst_1.numat, coord, q2, &c__1, vector_1.eigs, &nfilld, 
		funct, &eionis, &kchrge, &dip);
    }
    if (molkst_1.norbs > 0) {
	if (i_indx(keywrd_1.keywrd, "FOCK", 80L, 4L) != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(' FOCK MATRIX IS ')";
	    s_wsfe(&ci__1);
	    e_wsfe();
	    vecprt_(fokmat_1.f, &molkst_1.norbs);
	}
	if (i_indx(keywrd_1.keywrd, "DENSI", 80L, 5L) != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(//,20X,' DENSITY MATRIX IS ')";
	    s_wsfe(&ci__1);
	    e_wsfe();
	    vecprt_(densty_1.p, &molkst_1.norbs);
	} else {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(//10X,'ATOMIC ORBITAL ELECTRON POPULATIONS',/)";
	    s_wsfe(&ci__1);
	    e_wsfe();
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(8F10.5)";
	    s_wsfe(&ci__1);
	    i__2 = molkst_1.norbs;
	    for (i = 1; i <= i__2; ++i) {
		do_fio(&c__1, (char *)&densty_1.p[i * (i + 1) / 2 - 1], (
			ftnlen)sizeof(doublereal));
	    }
	    e_wsfe();
	}
	if (i_indx(keywrd_1.keywrd, " PI", 80L, 3L) != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(//10X,'SIGMA-PI BOND-ORDER MATRIX')";
	    s_wsfe(&ci__1);
	    e_wsfe();
	    denrot_();
	}
	if (uhf) {
	    sz = (i__2 = molkst_1.nalpha - molkst_1.nbeta, abs(i__2)) * .5;
	    ss2 = sz * sz;
	    l = 0;
	    i__2 = molkst_1.norbs;
	    for (i = 1; i <= i__2; ++i) {
		i__1 = i;
		for (j = 1; j <= i__1; ++j) {
		    ++l;
		    densty_1.pa[l - 1] -= densty_1.pb[l - 1];
/* L130: */
/* Computing 2nd power */
		    d__1 = densty_1.pa[l - 1];
		    ss2 += d__1 * d__1;
		}
/* L140: */
/* Computing 2nd power */
		d__1 = densty_1.pa[l - 1];
		ss2 -= d__1 * d__1 * .5;
	    }
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(//20X,'(SZ)    =',F10.6)";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, (char *)&sz, (ftnlen)sizeof(doublereal));
	    e_wsfe();
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(  20X,'(S**2)  =',F10.6)";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, (char *)&ss2, (ftnlen)sizeof(doublereal));
	    e_wsfe();
	    if (i_indx(keywrd_1.keywrd, "SPIN", 80L, 4L) != 0) {
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(//10X,'SPIN DENSITY MATRIX')";
		s_wsfe(&ci__1);
		e_wsfe();
		vecprt_(densty_1.pa, &molkst_1.norbs);
	    } else {
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(//10X,'ATOMIC ORBITAL SPIN POPULATIONS',/)";
		s_wsfe(&ci__1);
		e_wsfe();
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(8F10.5)";
		s_wsfe(&ci__1);
		i__2 = molkst_1.norbs;
		for (i = 1; i <= i__2; ++i) {
		    do_fio(&c__1, (char *)&densty_1.pa[i * (i + 1) / 2 - 1], (
			    ftnlen)sizeof(doublereal));
		}
		e_wsfe();
	    }
	    if (i_indx(keywrd_1.keywrd, "HYPERFINE", 80L, 9L) != 0) {

/*  WORK OUT THE HYPERFINE COUPLING CONSTANTS. */

		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(//10X,'    HYPERFINE COUPLING COEFFICIENTS',\
/)";
		s_wsfe(&ci__1);
		e_wsfe();
		j = (molkst_1.nalpha - 1) * molkst_1.norbs;
		i__2 = molkst_1.numat;
		for (k = 1; k <= i__2; ++k) {
		    i = molkst_1.nfirst[k - 1];
/* #          WRITE(6,'('' PA:'',F13.6,'' C('',I2,''+'',I3
,''):'', */
/* #     +F13.5)')PA((I*(I+1))/2),I,J,C(I+J) */
/* L150: */
/* Computing 2nd power */
		    d__1 = vector_1.c[i + j - 1];
		    q[k - 1] = densty_1.pa[i * (i + 1) / 2 - 1] * .3333333 + 
			    d__1 * d__1 * .66666666;
		}
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(5(2X,A2,I2,F9.5,1X))";
		s_wsfe(&ci__1);
		i__2 = molkst_1.numat;
		for (i = 1; i <= i__2; ++i) {
		    do_fio(&c__1, elemts_1.elemnt + (molkst_1.nat[i - 1] - 1 
			    << 1), 2L);
		    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
		    do_fio(&c__1, (char *)&q[i - 1], (ftnlen)sizeof(
			    doublereal));
		}
		e_wsfe();
	    }
	    i__2 = linear;
	    for (i = 1; i <= i__2; ++i) {
/* L160: */
		densty_1.pa[i - 1] = densty_1.p[i - 1] - densty_1.pb[i - 1];
	    }
	}
	if (i_indx(keywrd_1.keywrd, "BONDS", 80L, 5L) != 0) {
	    if (molkst_1.nbeta == 0) {
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(/10X,'BONDING CONTRIBUTION OF EACH M.O.',/)";
		s_wsfe(&ci__1);
		e_wsfe();
		molval_(vector_1.c, &molkst_1.norbs, densty_1.p, &
			molkst_1.norbs, &c_b171);
	    } else {
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(/10X,'BONDING CONTRIBUTION OF EACH ALPHA M.O\
.',/)";
		s_wsfe(&ci__1);
		e_wsfe();
		molval_(vector_1.c, &molkst_1.norbs, densty_1.p, &
			molkst_1.norbs, &c_b173);
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(/10X,'BONDING CONTRIBUTION OF EACH BETA  M.O\
.',/)";
		s_wsfe(&ci__1);
		e_wsfe();
		molval_(vector_1.c, &molkst_1.norbs, densty_1.p, &
			molkst_1.norbs, &c_b173);
	    }
	    bonds_(densty_1.p);
	}
	i = molkst_1.nclose + molkst_1.nalpha;
	if (i_indx(keywrd_1.keywrd, "LOCAL", 80L, 5L) != 0) {
	    local_(vector_1.c, &molkst_1.norbs, &i, vector_1.eigs);
	    if (molkst_1.nbeta != 0) {
		ci__1.cierr = 0;
		ci__1.ciunit = 6;
		ci__1.cifmt = "(//10X,' LOCALIZED BETA MOLECULAR ORBITALS')";
		s_wsfe(&ci__1);
		e_wsfe();
		local_(vector_1.cbeta, &molkst_1.norbs, &molkst_1.nbeta, 
			vector_1.eigb);
	    }
	}
	if (i_indx(keywrd_1.keywrd, "1ELE", 80L, 4L) != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(' FINAL ONE-ELECTRON MATRIX ')";
	    s_wsfe(&ci__1);
	    e_wsfe();
	    vecprt_(hmatrx_1.h, &molkst_1.norbs);
	}
	if (i_indx(keywrd_1.keywrd, "ENPART", 80L, 6L) != 0) {
	    enpart_(&uhf, hmatrx_1.h, densty_1.pa, densty_1.pb, densty_1.p, q,
		     coord);
	}
    }
    for (i = 1; i <= 107; ++i) {
/* L170: */
	nelemt[i - 1] = 0;
    }
    i__2 = molkst_1.numat;
    for (i = 1; i <= i__2; ++i) {
	igo = molkst_1.nat[i - 1];
	if (igo > 107) {
	    goto L180;
	}
	++nelemt[igo - 1];
L180:
	;
    }
    ichfor = 0;
    if (nelemt[5] == 0) {
	goto L190;
    }
    ichfor = 1;
    s_copy(ielemt, elemts_1.elemnt + 10, 2L, 2L);
    nzs = nelemt[5];
    if (nzs < 10) {
	if (nzs == 1) {
	    iel1[0] = 11;
	} else {
	    iel1[0] = nzs + 1;
	}
	iel2[0] = 11;
    } else {
	kfrst = nzs / 10;
	ksec = nzs - kfrst * 10;
	iel1[0] = kfrst + 1;
	iel2[0] = ksec + 1;
    }
L190:
    nelemt[5] = 0;
    for (i = 1; i <= 107; ++i) {
	if (nelemt[i - 1] == 0) {
	    goto L200;
	}
	++ichfor;
	s_copy(ielemt + (ichfor - 1 << 1), elemts_1.elemnt + (i - 1 << 1), 2L,
		 2L);
	nzs = nelemt[i - 1];
	if (nzs < 10) {
	    if (nzs == 1) {
		iel1[ichfor - 1] = 11;
	    } else {
		iel1[ichfor - 1] = nzs + 1;
	    }
	    iel2[ichfor - 1] = 11;
	} else {
	    kfrst = nzs / 10;
	    ksec = nzs - kfrst * 10;
	    iel1[ichfor - 1] = kfrst + 1;
	    iel2[ichfor - 1] = ksec + 1;
	}
L200:
	;
    }
    if (i_indx(keywrd_1.keywrd, "DENOUT", 80L, 6L) != 0) {
	getarg_(&c__3, argz_1.argz, 512L);
	o__1.oerr = 0;
	o__1.ounit = 10;
	o__1.ofnmlen = 512;
	o__1.ofnm = argz_1.argz;
	o__1.orl = 0;
	o__1.osta = "UNKNOWN";
	o__1.oacc = 0;
	o__1.ofm = "UNFORMATTED";
	o__1.oblnk = 0;
	f_open(&o__1);
	al__1.aerr = 0;
	al__1.aunit = 10;
	f_rew(&al__1);
	s_wsue(&io___59);
	i__2 = linear;
	for (i = 1; i <= i__2; ++i) {
	    do_uio(&c__1, (char *)&densty_1.pa[i - 1], (ftnlen)sizeof(
		    doublereal));
	}
	e_wsue();
	if (uhf) {
	    s_wsue(&io___60);
	    i__2 = linear;
	    for (i = 1; i <= i__2; ++i) {
		do_uio(&c__1, (char *)&densty_1.pb[i - 1], (ftnlen)sizeof(
			doublereal));
	    }
	    e_wsue();
	}
	cl__1.cerr = 0;
	cl__1.cunit = 10;
	cl__1.csta = 0;
	f_clos(&cl__1);
    }
    if ((ci || molkst_1.nopen != molkst_1.nclose || i_indx(keywrd_1.keywrd, 
	    "SIZE", 80L, 4L) != 0) && i_indx(keywrd_1.keywrd, "MECI", 80L, 4L)
	     + i_indx(keywrd_1.keywrd, "ESR", 80L, 3L) != 0) {
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(//10X,'MULTI-ELECTRON CONFIGURATION INTERACTION CALC\
ULATION',//)";
	s_wsfe(&ci__1);
	e_wsfe();
	nmos = 0;
	ncis = 0;
	if (i_indx(keywrd_1.keywrd, "C.I.=", 80L, 5L) != 0) {
	    i__2 = i_indx(keywrd_1.keywrd, "C.I.=", 80L, 5L) + 5;
	    nmos = (integer) reada_(keywrd_1.keywrd, &i__2, 80L);
	}

/*   SET UP C.I. PARAMETERS */
/*   NMOS IS NO. OF M.O.S USED IN C.I. */
/*   NCIS IS CHANGE IN SPIN, OR NUMBER OF STATES */

	if (nmos == 0) {
	    nmos = molkst_1.nopen - molkst_1.nclose;
	}
	if (ncis == 0) {
	    if (triplt || i_indx(keywrd_1.keywrd, "QUAR", 80L, 4L) != 0) {
		ncis = 1;
	    }
	    if (i_indx(keywrd_1.keywrd, "QUIN", 80L, 4L) + i_indx(
		    keywrd_1.keywrd, "SEXT", 80L, 4L) != 0) {
		ncis = 2;
	    }
	}
	x = meci_(vector_1.eigs, vector_1.c, vector_1.cbeta, vector_1.eigb, &
		molkst_1.norbs, &nmos, &ncis, &c_true);
    }
    if (i_indx(keywrd_1.keywrd, "MULLIK", 80L, 6L) + i_indx(keywrd_1.keywrd, 
	    "GRAPH", 80L, 5L) != 0) {
	if (i_indx(keywrd_1.keywrd, "MULLIK", 80L, 6L) != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(/10X,' MULLIKEN POPULATION ANALYSIS')";
	    s_wsfe(&ci__1);
	    e_wsfe();
	}
	mullik_(vector_1.c, vector_1.cbeta, &uhf, hmatrx_1.h, fokmat_1.f, &
		molkst_1.norbs, densty_1.p, scrach_1.rxyz);
	if (i_indx(keywrd_1.keywrd, "GRAPH", 80L, 5L) != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = 6;
	    ci__1.cifmt = "(/10X,' DATA FOR GRAPH WRITTEN TO DISK')";
	    s_wsfe(&ci__1);
	    e_wsfe();
	}
    }

/*  NOTE THAT THE DENSITY, H AND F MATRICES ARE CORRUPTED BY A */
/*  CALL TO MULLIK. */

/*   FOLLOWING SUBROUTINE CALL ADDED BY VIC L. TO CALCULATE MULLIKEN */
/*   POPULATIONS FOR EACH ATOM  (28JULY86) */

    if (i_indx(keywrd_1.keywrd, "MULLIK", 80L, 6L) == 0) {

/*               MULLIKEN ANALYSIS WAS NOT PERFORMED */

	mpcpop_(vector_1.c, &c__0);
    } else {

/*               WRITE OUT POP AND CHARGES */

	mpcpop_(vector_1.c, &c__1);
    }

/*  NOTE THAT THE DENSITY, H AND F MATRICES ARE CORRUPTED BY A */
/*  CALL TO MULLIK. */
    if (first) {

/*  FOLLOWING OPEN STATEMENT WAS CHANGED BY VIC L. TO CHANGE STATUS */
/*  FROM NEW TO UNKNOWN */

	getarg_(&c__1, argz_1.argz, 512L);
/*        OPEN(UNIT=12,FILE=argz,STATUS='NEW') */
	getarg_(&c__1, argz_1.argz, 512L);
	o__1.oerr = 0;
	o__1.ounit = 12;
	o__1.ofnmlen = 512;
	o__1.ofnm = argz_1.argz;
	o__1.orl = 0;
	o__1.osta = "UNKNOWN";
	o__1.oacc = 0;
	o__1.ofm = 0;
	o__1.oblnk = 0;
	f_open(&o__1);
	al__1.aerr = 0;
	al__1.aunit = 12;
	f_rew(&al__1);
	first = FALSE_;
    }
    iwrite = 12;
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(//20X,' SUMMARY OF ',A7,' CALCULATION',/)";
    s_wsfe(&ci__1);
    do_fio(&c__1, caltyp, 7L);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(60X,'VERSION ',F5.2)";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&c_b36, (ftnlen)sizeof(doublereal));
    e_wsfe();
    io___64.ciunit = iwrite;
    s_wsfe(&io___64);
    i__2 = ichfor;
    for (i = 1; i <= i__2; ++i) {
	do_fio(&c__1, ielemt + (i - 1 << 1), 2L);
	do_fio(&c__1, numbrs + (iel1[i - 1] - 1), 1L);
	do_fio(&c__1, numbrs + (iel2[i - 1] - 1), 1L);
    }
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(55X,A24)";
    s_wsfe(&ci__1);
    do_fio(&c__1, idate, 24L);
    e_wsfe();
    for (ik = 80; ik >= 3; --ik) {
/* L220: */
	if (*(unsigned char *)&titles_1.koment[ik - 1] != ' ') {
	    goto L230;
	}
    }
L230:
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, titles_1.koment, ik);
    e_wsfe();
    for (it = 80; it >= 3; --it) {
/* L240: */
	if (*(unsigned char *)&titles_1.title[it - 1] != ' ') {
	    goto L250;
	}
    }
L250:
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, titles_1.title, it);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(//4X,A58)";
    s_wsfe(&ci__1);
    do_fio(&c__1, flepo + (mesage_1.iflepo - 1) * 58, 58L);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(4X,A58)";
    s_wsfe(&ci__1);
    do_fio(&c__1, iter + (mesage_1.iiter - 1) * 58, 58L);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(//10X,'HEAT OF FORMATION       =',F17.6,' KCAL')";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&(*funct), (ftnlen)sizeof(doublereal));
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(  10X,'ELECTRONIC ENERGY       =',F17.6,' EV')";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&elect_1.elect, (ftnlen)sizeof(doublereal));
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(  10X,'CORE-CORE REPULSION     =',F17.6,' EV')";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&enuclr_1.enuclr, (ftnlen)sizeof(doublereal));
    e_wsfe();
    if (prtgra) {
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(  10X,'GRADIENT NORM           =',F17.6)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&gradnt_1.gnorm, (ftnlen)sizeof(doublereal));
	e_wsfe();
    }
    if (path_1.latom != 0) {
	s_copy(grtype, " KCAL/ANGSTROM", 14L, 14L);
	if (path_1.lparam == 1) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = iwrite;
	    ci__1.cifmt = "(    10X,'FOR REACTION COORDINATE =',F17.4       \
 ,' ANGSTROMS')";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, (char *)&xreact, (ftnlen)sizeof(doublereal));
	    e_wsfe();
	} else {
	    if (geokst_1.na[0] != 99) {
		s_copy(grtype, " KCAL/RADIAN  ", 14L, 14L);
	    }
	    ci__1.cierr = 0;
	    ci__1.ciunit = iwrite;
	    ci__1.cifmt = "(    10X,'FOR REACTION COORDINATE =',F17.4       \
 ,' DEGREES')";
	    s_wsfe(&ci__1);
	    d__1 = xreact * degree;
	    do_fio(&c__1, (char *)&d__1, (ftnlen)sizeof(doublereal));
	    e_wsfe();
	}
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(    10X,'REACTION GRADIENT       =',F17.6,A14    )";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&gcoord, (ftnlen)sizeof(doublereal));
	do_fio(&c__1, grtype, 14L);
	e_wsfe();
    }
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(  10X,'DIPOLE                  =',F16.5, ' DEBYE')";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&dip, (ftnlen)sizeof(doublereal));
    e_wsfe();
    if (uhf) {
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(  10X,'(SZ)                    =',F17.6)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&sz, (ftnlen)sizeof(doublereal));
	e_wsfe();
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(  10X,'(S**2)                  =',F17.6)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&ss2, (ftnlen)sizeof(doublereal));
	e_wsfe();
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(  10X,'NO. OF ALPHA ELECTRONS  =',I10)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&molkst_1.nalpha, (ftnlen)sizeof(integer));
	e_wsfe();
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(  10X,'NO. OF BETA  ELECTRONS  =',I10)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&molkst_1.nbeta, (ftnlen)sizeof(integer));
	e_wsfe();
    } else {
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(  10X,'NO. OF FILLED LEVELS    =',I10)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&molkst_1.nclose, (ftnlen)sizeof(integer));
	e_wsfe();
	nopn = molkst_1.nopen - molkst_1.nclose;
	if (nopn != 0) {
	    ci__1.cierr = 0;
	    ci__1.ciunit = iwrite;
	    ci__1.cifmt = "(  10X,'AND NO. OF OPEN LEVELS  =',I10)";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, (char *)&nopn, (ftnlen)sizeof(integer));
	    e_wsfe();
	}
    }
    if (ci) {
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(  10X,'CONFIGURATION INTERACTION WAS USED')";
	s_wsfe(&ci__1);
	e_wsfe();
    }
    if (kchrge != 0) {
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(  10X,'CHARGE ON SYSTEM        =',I10)";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&kchrge, (ftnlen)sizeof(integer));
	e_wsfe();
    }
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(  10X,'IONIZATION POTENTIAL    =',F17.6,' EV')";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&eionis, (ftnlen)sizeof(doublereal));
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(  10X,'MOLECULAR WEIGHT        =',F14.3)";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&sumw, (ftnlen)sizeof(doublereal));
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(  10X,'SCF CALCULATIONS        =',I10)";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&numscf_1.nscf, (ftnlen)sizeof(integer));
    e_wsfe();
    tim = second_() - *time0;
    timout_(&iwrite, &tim);
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(//10X,'FINAL GEOMETRY OBTAINED',36X,'CHARGE')";
    s_wsfe(&ci__1);
    e_wsfe();
    for (i = 80; i >= 3; --i) {
/* L260: */
	if (*(unsigned char *)&keywrd_1.keywrd[i - 1] != ' ') {
	    goto L270;
	}
    }
L270:
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, keywrd_1.keywrd, i);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, titles_1.koment, ik);
    e_wsfe();
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, titles_1.title, it);
    e_wsfe();
    na1 = geokst_1.na[0];
    if (xyz) {
	xyzint_(geom_1.geo, &geokst_1.natoms, geokst_1.na, geokst_1.nb, 
		geokst_1.nc, &c_b173, coord);
    }
    degree = 57.29577951;
    coord[1] = 0.;
    coord[2] = 0.;
    coord[0] = 0.;
    coord[4] = 0.;
    coord[5] = 0.;
    coord[8] = 0.;
    ivar = 1;
    geokst_1.na[0] = 0;
    l = 0;
    i__2 = geokst_1.natoms;
    for (i = 1; i <= i__2; ++i) {
	for (j = 1; j <= 3; ++j) {
	    if (! xyz) {
		coord[j + i * 3 - 4] = geom_1.geo[j + i * 3 - 4];
	    }
/* L280: */
	    iel1[j - 1] = 0;
	}
L290:
	if (geovar_1.loc[(ivar << 1) - 2] == i) {
	    iel1[geovar_1.loc[(ivar << 1) - 1] - 1] = 1;
	    ++ivar;
	    goto L290;
	}
	if (i < 4) {
	    iel1[2] = 0;
	    if (i < 3) {
		iel1[1] = 0;
		if (i < 2) {
		    iel1[0] = 0;
		}
	    }
	}
	if (i == path_1.latom) {
	    iel1[path_1.lparam - 1] = -1;
	}
	q[0] = coord[i * 3 - 3];
	q[1] = coord[i * 3 - 2] * degree;
	q[2] = coord[i * 3 - 1] * degree;
	if (geokst_1.labels[i - 1] != 107 && geokst_1.labels[i - 1] != 99) {
	    ++l;
	    ci__1.cierr = 0;
	    ci__1.ciunit = iwrite;
	    ci__1.cifmt = "(1X,A2,3(F12.6,I3),I5,2I5,F13.4)";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, elemts_1.elemnt + (geokst_1.labels[i - 1] - 1 << 1),
		     2L);
	    for (k = 1; k <= 3; ++k) {
		do_fio(&c__1, (char *)&q[k - 1], (ftnlen)sizeof(doublereal));
		do_fio(&c__1, (char *)&iel1[k - 1], (ftnlen)sizeof(integer));
	    }
	    do_fio(&c__1, (char *)&geokst_1.na[i - 1], (ftnlen)sizeof(integer)
		    );
	    do_fio(&c__1, (char *)&geokst_1.nb[i - 1], (ftnlen)sizeof(integer)
		    );
	    do_fio(&c__1, (char *)&geokst_1.nc[i - 1], (ftnlen)sizeof(integer)
		    );
	    do_fio(&c__1, (char *)&q2[l - 1], (ftnlen)sizeof(doublereal));
	    e_wsfe();
	} else {
	    ci__1.cierr = 0;
	    ci__1.ciunit = iwrite;
	    ci__1.cifmt = "(1X,A2,3(F12.6,I3),I5,2I5,F13.4)";
	    s_wsfe(&ci__1);
	    do_fio(&c__1, elemts_1.elemnt + (geokst_1.labels[i - 1] - 1 << 1),
		     2L);
	    for (k = 1; k <= 3; ++k) {
		do_fio(&c__1, (char *)&q[k - 1], (ftnlen)sizeof(doublereal));
		do_fio(&c__1, (char *)&iel1[k - 1], (ftnlen)sizeof(integer));
	    }
	    do_fio(&c__1, (char *)&geokst_1.na[i - 1], (ftnlen)sizeof(integer)
		    );
	    do_fio(&c__1, (char *)&geokst_1.nb[i - 1], (ftnlen)sizeof(integer)
		    );
	    do_fio(&c__1, (char *)&geokst_1.nc[i - 1], (ftnlen)sizeof(integer)
		    );
	    e_wsfe();
	}
/* L300: */
    }
    geokst_1.na[0] = na1;
    i = 0;
    x = 0.;
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(I3,3(F12.6,I3),I5,2I5)";
    s_wsfe(&ci__1);
    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
    do_fio(&c__1, (char *)&x, (ftnlen)sizeof(doublereal));
    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
    do_fio(&c__1, (char *)&x, (ftnlen)sizeof(doublereal));
    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
    do_fio(&c__1, (char *)&x, (ftnlen)sizeof(doublereal));
    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
    do_fio(&c__1, (char *)&i, (ftnlen)sizeof(integer));
    e_wsfe();
    i__2 = geosym_1.ndep;
    for (i = 1; i <= i__2; ++i) {
/* L310: */
	ci__1.cierr = 0;
	ci__1.ciunit = iwrite;
	ci__1.cifmt = "(3(I4,','))";
	s_wsfe(&ci__1);
	do_fio(&c__1, (char *)&geosym_1.locpar[i - 1], (ftnlen)sizeof(integer)
		);
	do_fio(&c__1, (char *)&geosym_1.idepfn[i - 1], (ftnlen)sizeof(integer)
		);
	do_fio(&c__1, (char *)&geosym_1.locdep[i - 1], (ftnlen)sizeof(integer)
		);
	e_wsfe();
    }
    ci__1.cierr = 0;
    ci__1.ciunit = iwrite;
    ci__1.cifmt = "(///)";
    s_wsfe(&ci__1);
    e_wsfe();
    numscf_1.nscf = 0;
    return 0;
} /* write_ */

#undef w


/* Subroutine */ int timout_(nout, tim)
integer *nout;
doublereal *tim;
{
    /* Initialized data */

    static doublereal hrspd = 24.;
    static doublereal minphr = 60.;
    static doublereal secpd = 86400.;
    static doublereal secpmi = 60.;

    /* Format strings */
    static char fmt_10[] = "(10x,\002COMPUTATION TIME = \002,i2,1x,\002DAY\
S\002,2x,i2,1x,\002HOURS\002,1x,i2,1x,\002MINUTES AND\002,1x,f7.3,1x,\002SEC\
ONDS\002)";
    static char fmt_20[] = "(10x,\002COMPUTATION TIME = \002,i2,1x,\002DA\
Y\002,2x,i2,1x,\002HOURS\002,1x,i2,1x,\002MINUTES AND\002,1x,f7.3,1x,\002SEC\
ONDS\002)";
    static char fmt_30[] = "(10x,\002COMPUTATION TIME = \002i2,1x,\002HOUR\
S\002,1x,i2,1x,\002MINUTES AND\002,1x,f7.3,1x,\002SECONDS\002)";
    static char fmt_40[] = "(10x,\002COMPUTATION TIME = \002,i2,1x,\002MINUT\
ES AND\002,1x,f7.3,1x,\002SECONDS\002)";
    static char fmt_50[] = "(10x,\002COMPUTATION TIME = \002,f7.3,1x,\002SEC\
ONDS\002)";

    /* Builtin functions */
    integer s_wsfe(), do_fio(), e_wsfe();

    /* Local variables */
    doublereal secs, days, mins;
    integer idays, imins;
    doublereal hours;
    integer ihours;

    /* Fortran I/O blocks */
    static cilist io___81 = { 0, 0, 0, fmt_10, 0 };
    static cilist io___82 = { 0, 0, 0, fmt_20, 0 };
    static cilist io___83 = { 0, 0, 0, fmt_30, 0 };
    static cilist io___84 = { 0, 0, 0, fmt_40, 0 };
    static cilist io___85 = { 0, 0, 0, fmt_50, 0 };



/*     CONVERT THE TIME FROM SECONDS TO DAYS, HOURS, MINUTES, AND SECONDS 
*/





    days = *tim / secpd;
    idays = (integer) days;
    hours = (days - (real) idays) * hrspd;
    ihours = (integer) hours;
    mins = (hours - (real) ihours) * minphr;
    imins = (integer) mins;
    secs = (mins - (real) imins) * secpmi;

    if (idays > 1) {
	io___81.ciunit = *nout;
	s_wsfe(&io___81);
	do_fio(&c__1, (char *)&idays, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&ihours, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&imins, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&secs, (ftnlen)sizeof(doublereal));
	e_wsfe();
    } else if (idays == 1) {
	io___82.ciunit = *nout;
	s_wsfe(&io___82);
	do_fio(&c__1, (char *)&idays, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&ihours, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&imins, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&secs, (ftnlen)sizeof(doublereal));
	e_wsfe();
    } else if (ihours > 0) {
	io___83.ciunit = *nout;
	s_wsfe(&io___83);
	do_fio(&c__1, (char *)&ihours, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&imins, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&secs, (ftnlen)sizeof(doublereal));
	e_wsfe();
    } else if (imins > 0) {
	io___84.ciunit = *nout;
	s_wsfe(&io___84);
	do_fio(&c__1, (char *)&imins, (ftnlen)sizeof(integer));
	do_fio(&c__1, (char *)&secs, (ftnlen)sizeof(doublereal));
	e_wsfe();
    } else {
	io___85.ciunit = *nout;
	s_wsfe(&io___85);
	do_fio(&c__1, (char *)&secs, (ftnlen)sizeof(doublereal));
	e_wsfe();
    }

} /* timout_ */

/* Subroutine */ int mpcpop_(c, icok)
doublereal *c;
integer *icok;
{
    /* Format strings */
    static char fmt_70[] = "(5x,i4,4x,f11.6,6x,f11.6)";
    static char fmt_80[] = "(2f12.6)";

    /* System generated locals */
    integer i__1, i__2;
    cilist ci__1;

    /* Builtin functions */
    integer s_wsfe(), do_fio(), e_wsfe();

    /* Local variables */
    doublereal chrg[86];
    integer i, j, k, if__, il;
    doublereal pop[86], sum;

    /* Fortran I/O blocks */
    static cilist io___94 = { 0, 6, 0, fmt_70, 0 };
    static cilist io___95 = { 1, 16, 0, fmt_80, 0 };



/* This subroutine calculates the total Mulliken populations on the */
/*   atoms by summing the diagonal elements from the  Mulliken */
/*   population analysis. */

/* 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 */

    /* Parameter adjustments */
    --c;

    /* Function Body */
    ci__1.cierr = 1;
    ci__1.ciunit = 16;
    ci__1.cifmt = "(I4,5X' MULLIKEN POPULATION AND CHARGE')";
    i__1 = s_wsfe(&ci__1);
    if (i__1 != 0) {
	goto L40;
    }
    i__1 = do_fio(&c__1, (char *)&(*icok), (ftnlen)sizeof(integer));
    if (i__1 != 0) {
	goto L40;
    }
    i__1 = e_wsfe();
    if (i__1 != 0) {
	goto L40;
    }

/* ICOK = 1 ==> PRINT POPULATIONS */
/* ICOK = 0 ==> KEYWORD mulliken = .f. */
/*         NO POPULATION ANALYSIS PERFORMED */

    if (*icok != 0) {
	i__1 = molkst_2.numat;
	for (i = 1; i <= i__1; ++i) {
	    if__ = molkst_2.nfirst[i - 1];
	    il = molkst_2.nlast[i - 1];
	    sum = (float)0.;
	    pop[i - 1] = (float)0.;
	    chrg[i - 1] = (float)0.;
	    i__2 = il;
	    for (j = if__; j <= i__2; ++j) {

/*    Diagonal element of mulliken matrix */

		sum += c[j * (j + 1) / 2];
/* L10: */
	    }
	    k = molkst_2.nat[i - 1];

/*    Mulliken population for i'th atom */

	    pop[i - 1] = sum;
	    chrg[i - 1] = core_1.core[k - 1] - pop[i - 1];
/* L20: */
	}
	ci__1.cierr = 0;
	ci__1.ciunit = 6;
	ci__1.cifmt = "(///10X,'MULLIKEN POPULATIONS AND CHARGES')";
	s_wsfe(&ci__1);
	e_wsfe();
	i__1 = molkst_2.numat;
	for (j = 1; j <= i__1; ++j) {
	    s_wsfe(&io___94);
	    do_fio(&c__1, (char *)&j, (ftnlen)sizeof(integer));
	    do_fio(&c__1, (char *)&pop[j - 1], (ftnlen)sizeof(doublereal));
	    do_fio(&c__1, (char *)&chrg[j - 1], (ftnlen)sizeof(doublereal));
	    e_wsfe();
	    i__2 = s_wsfe(&io___95);
	    if (i__2 != 0) {
		goto L40;
	    }
	    i__2 = do_fio(&c__1, (char *)&pop[j - 1], (ftnlen)sizeof(
		    doublereal));
	    if (i__2 != 0) {
		goto L40;
	    }
	    i__2 = do_fio(&c__1, (char *)&chrg[j - 1], (ftnlen)sizeof(
		    doublereal));
	    if (i__2 != 0) {
		goto L40;
	    }
	    i__2 = e_wsfe();
	    if (i__2 != 0) {
		goto L40;
	    }
/* L30: */
	}
    }
    return 0;
L40:
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, "Error writing SYBYL Mulliken population output", 46L);
    e_wsfe();
    return 0;
/* L50: */
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, "Error opening output file for SYBYL Mulliken output", 51L);
    e_wsfe();
    return 0;
/* L60: */
} /* mpcpop_ */


/* This subroutine writes out the optimized geometry and atomic charges */
/*   for a MOPAC run. */

/* Subroutine */ int mpcsyb_(numat, coord, chr, icok, eigs, nclose, funct, 
	eionis, kchrge, dip)
integer *numat;
doublereal *coord, *chr;
integer *icok;
doublereal *eigs;
integer *nclose;
doublereal *funct, *eionis;
integer *kchrge;
doublereal *dip;
{
    /* Format strings */
    static char fmt_20[] = "(4f12.6,2x,i4,2x,\002HOMOs,LUMOs,# of occupied M\
Os\002)";

    /* System generated locals */
    integer i__1, i__2;
    cilist ci__1;
    olist o__1;

    /* Builtin functions */
    integer f_open(), s_wsfe(), do_fio(), e_wsfe();

    /* Local variables */
    integer i, j, i1, i2;
    extern /* Subroutine */ int getarg_();

    /* Fortran I/O blocks */
    static cilist io___100 = { 1, 16, 0, fmt_20, 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 */

    /* Parameter adjustments */
    --chr;
    coord -= 4;
    --eigs;

    /* Function Body */
    getarg_(&c__2, argz_1.argz, 512L);
    o__1.oerr = 1;
    o__1.ounit = 16;
    o__1.ofnmlen = 512;
    o__1.ofnm = argz_1.argz;
    o__1.orl = 0;
    o__1.osta = 0;
    o__1.oacc = 0;
    o__1.ofm = 0;
    o__1.oblnk = 0;
    i__1 = f_open(&o__1);
    if (i__1 != 0) {
	goto L40;
    }
/*  Write out the charge flag and number of atoms */
    ci__1.cierr = 1;
    ci__1.ciunit = 16;
    ci__1.cifmt = "(2I4)";
    i__1 = s_wsfe(&ci__1);
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = do_fio(&c__1, (char *)&(*icok), (ftnlen)sizeof(integer));
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = do_fio(&c__1, (char *)&(*numat), (ftnlen)sizeof(integer));
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = e_wsfe();
    if (i__1 != 0) {
	goto L30;
    }
/*  Write out the coordinates and charges */
    i__1 = *numat;
    for (i = 1; i <= i__1; ++i) {
	ci__1.cierr = 1;
	ci__1.ciunit = 16;
	ci__1.cifmt = "(4F12.6)";
	i__2 = s_wsfe(&ci__1);
	if (i__2 != 0) {
	    goto L30;
	}
	for (j = 1; j <= 3; ++j) {
	    i__2 = do_fio(&c__1, (char *)&coord[j + i * 3], (ftnlen)sizeof(
		    doublereal));
	    if (i__2 != 0) {
		goto L30;
	    }
	}
	i__2 = do_fio(&c__1, (char *)&chr[i], (ftnlen)sizeof(doublereal));
	if (i__2 != 0) {
	    goto L30;
	}
	i__2 = e_wsfe();
	if (i__2 != 0) {
	    goto L30;
	}
/* L10: */
    }
/* Computing MAX */
    i__1 = 1, i__2 = *nclose - 1;
    i1 = max(i__1,i__2);
/* Computing MIN */
    i__1 = 215, i__2 = *nclose + 2;
    i2 = min(i__1,i__2);

/*  Write out the 2 highest and 2 lowest orbital energies */

    i__1 = s_wsfe(&io___100);
    if (i__1 != 0) {
	goto L30;
    }
    i__2 = i2;
    for (j = i1; j <= i__2; ++j) {
	i__1 = do_fio(&c__1, (char *)&eigs[j], (ftnlen)sizeof(doublereal));
	if (i__1 != 0) {
	    goto L30;
	}
    }
    i__1 = do_fio(&c__1, (char *)&(*nclose), (ftnlen)sizeof(integer));
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = e_wsfe();
    if (i__1 != 0) {
	goto L30;
    }

/*  Write out the Heat of Formation and Ionisation Potential */

    ci__1.cierr = 1;
    ci__1.ciunit = 16;
    ci__1.cifmt = "(2F12.6,4X,'HF and IP')";
    i__1 = s_wsfe(&ci__1);
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = do_fio(&c__1, (char *)&(*funct), (ftnlen)sizeof(doublereal));
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = do_fio(&c__1, (char *)&(*eionis), (ftnlen)sizeof(doublereal));
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = e_wsfe();
    if (i__1 != 0) {
	goto L30;
    }

/*  Write out the Dipole Moment */

    if (*kchrge != 0) {
	*dip = (float)0.;
    }
    ci__1.cierr = 1;
    ci__1.ciunit = 16;
    ci__1.cifmt = "(I4,F10.3,'  Charge,Dipole Moment')";
    i__1 = s_wsfe(&ci__1);
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = do_fio(&c__1, (char *)&(*kchrge), (ftnlen)sizeof(integer));
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = do_fio(&c__1, (char *)&(*dip), (ftnlen)sizeof(doublereal));
    if (i__1 != 0) {
	goto L30;
    }
    i__1 = e_wsfe();
    if (i__1 != 0) {
	goto L30;
    }
    return 0;
L30:
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, "Error writing SYBYL MOPAC output", 32L);
    e_wsfe();
    return 0;
L40:
    ci__1.cierr = 0;
    ci__1.ciunit = 6;
    ci__1.cifmt = "(A)";
    s_wsfe(&ci__1);
    do_fio(&c__1, "Error opening SYBYL MOPAC output", 32L);
    e_wsfe();
    return 0;
} /* mpcsyb_ */

