      SUBROUTINE ADVNCE
C     ****                   ADVANCES THE MODEL IN TIME
      include "params.h"
      include "blnk.h"
      include "vscr.h"
      include "buff.h"
      include "cons.h"
      include "index.h"
      include "strt.h"
      include "unit.h"
      include "phys.h"
C     ****
C     ****     DIMENSIONS FOR AMIE FIELDS
C     ****
      COMMON/AMIE/IAMIE,IUTS1,IUTS2,NSPSEC,NTIMS,NTST,SECTGCM,CRAD(2),
     1 PHIDA(2),SECUTA(2),HPSHA(2),HPNHA(2),CPSHA(2),CPNHA(2),
     2 CUSPSTA(2),CUSPNTA(2),CUSPSLA(2),CUSPNLA(2),
     3 POTKVA1(73,36),VIXYZA1(73,36,3),EKEVA1(73,36),EFLXA1(73,36),
     4 POTKVA2(73,36),VIXYZA2(73,36,3),EKEVA2(73,36),EFLXA2(73,36)
      WRITE(6,5) ITER
    5 FORMAT(1H0,17X,*ADVANCE ENTERED AT ITER = *I6)
C     ****                   INITIALIZE BUFFER INDICES
      MODULO=8
      CALL INIT
C     ****                   TIME LOOP
      TIME = SECOND(X)
      TIMET = TIMEF(X)
  100 CONTINUE
      ITER=ITER+1
      JTOTS=0.
      JTOTN=0.
C     ****
C     ****     FIND UT TIME BETWEEN 0-24 HR (MODEL DAY = UT DAY, 11/89)
C     ****
      SECS = AMOD(FLOAT(ITER)*C(4),86400.)
      IF (IADVDA .EQ. 1) THEN
       IDAYIT = ITER*IFIX(C(4))/86400
       IF (IDAYIT*86400 .EQ. ITER*IFIX(C(4))) THEN
C       ****
C       ****     ADVANCE DAY
C       ****
        IDAPR = IIDAY
        IYRPR = IYEAR
        IIDAY = IIDAY + 1
C       ****
C       ****     LPYR = 1(0)  IF IS (NOT) A LEAP YEAR
C       ****
        IYR4 = IYEAR/4
        IYR100 = IYEAR/100
        LPYR = 0
        IF (IYR4*4 .EQ. IYEAR .AND. IYR100*100 .NE. IYEAR) LPYR=1
        IENDA = 365 + LPYR
        IF (IIDAY .GT. IENDA) THEN
         IYEAR = IYEAR + 1
         IIDAY = IIDAY - IENDA
        ENDIF                          !   FOR PAST YEAR'S END
C       ****
C       ****     RECALCULATE SUN'S DECLINATION
C       ****
        PI = 3.14159265358979
        DELTA = ATAN(TAN(23.5*PI/180.)*SIN(2.*PI*FLOAT(IIDAY-80)/365.))
        C(95) = SIN(DELTA)
        C(96) = COS(DELTA)
        WRITE(6,"(1X,'ADVNCE:  ADVANCING DAY (PREVIOUS,PRESENT)=',4I5)")
     1   IDAPR,IYRPR,IIDAY,IYEAR
       ENDIF                           !   FOR ADVANCING DAY
      ENDIF                           !   FOR IADVDA FLAG ON
C     ****
C     ****     CALCULATE AURORAL PARAMETERS FOR THIS TIME STEP
C     ****
C     ****     FIND UT TIME  (MODEL DAY = UT DAY, 11/89)
      SEC = FLOAT(ITER)*C(4)
      SECTGCM = FLOAT(ITER)*C(4)
      CALL SUN
      IF (IAMIE.EQ.1) THEN
       IF (NTST.LT.NTIMS .AND. SEC.GT.SECUTA(2)) THEN
        NTST = NTST + 1
        SECUTA(1) = SECUTA(2)
        HPNHA(1) = HPNHA(2)
        HPSHA(1) = HPSHA(2)
        CPNHA(1) = CPNHA(2)
        CPSHA(1) = CPSHA(2)
        CUSPSTA(1) = CUSPSTA(2)
        CUSPSLA(1) = CUSPSLA(2)
        CUSPNTA(1) = CUSPNTA(2)
        CUSPNLA(1) = CUSPNLA(2)
        DO 120 JA=1,36
          DO 120 IA=1,73
            POTKVA1(IA,JA) = POTKVA2(IA,JA)
            EKEVA1(IA,JA) = EKEVA2(IA,JA)
            EFLXA1(IA,JA) = EFLXA2(IA,JA)
  120   CONTINUE
        DO 121 KA=1,3
          DO 121 JA=1,36
            DO 121 IA=1,73
              VIXYZA1(IA,JA,KA) = VIXYZA2(IA,JA,KA)
  121   CONTINUE
        READ (IUAMI) SECUTA(2),HPSHA(2),HPNHA(2),CPSHA(2),CPNHA(2),
     1    CUSPSTA(2),CUSPNTA(2),CUSPSLA(2),CUSPNLA(2),POTKVA2,VIXYZA2,
     2    EKEVA2,EFLXA2
       ENDIF                    !  FOR READING AMIE
      ENDIF                     !  FOR IAMIE=1
      CALL TAIL (SEC)
      IHIS=IHIS-1
      ISAV=ISAV-1
C     ****                   FOR THE FIRST FORWARD TIME STEP OF A
C     ****                   FORECAST START, HALVE DT.
      IF(MODEFC)150,130,140
  130 C(4) = C(111)/2.
      C(6)   = 2.*C(4)
      C(7)   = 1./C(6)
      MODEFC=1
      GO TO 150
C     ****                   RESTORE DT TO ORIGINAL VALUE
  140 C(4) = C(111)
      C(6)   = 2.*C(4)
      C(7)   = 1./C(6)
      MODEFC=-1
  150 CONTINUE
C     ****                   BEGIN SWEEP OF MESH FROM SOUTH TO NORTH
C     ****                   BRING IN LINE 1
      CALL VDRIFT
      CALL FLIP
      CALL UNLCRD(IUNIN,NJIN,LBDSKI,NSIN)
      CALL UNLCCK(IUNIN,NSIN,LEN)
      IF(LEN.EQ.LBDSKI.AND.NSIN.EQ.0) GO TO 155
      NUM=155
      GO TO 4000
C     ****                   BRING IN LINES 2 AND 3
C     ****                   AND REARRANGE LINES 1 AND 2
C     ****                   TO FORM LINES 0 AND -1
  155 DO 270 JP2=1,2
      CALL FLIP
      CALL UNLCRD(IUNIN,NJIN,LBDSKI,NSIN)
C     ****                   ADD DIAGNOSTIC QUANTITIES ONTO LINE JP2
      J=JP2-2
      CALL ADDIAG
C     ****                   REARRANGE LINE JP2 TO GET LINE -JP2+1
      NJPUT=NJP1
      IF(JP2.EQ.2) NJPUT=NJM1
      IF(NHEMI.NE.0) GO TO 200
C     ****                   CONSTRUCT LINES -1,0 FOR GLOBAL CASE
      ILIM=IMAXH+2
      DO 195 N=1,NFLDS
C     ****                   CHECK IF FIELD IS NEEDED ACROSS POLE
      IF(KPOLE(N)) 160,195,160
C     ****                   TRANSFER Z TYPES FIRST ITERATION ONLY
  160 IF(IFRST.EQ.1 .AND. KPOLE(N).EQ.2) GO TO 195
C     ****                   DETERMINE LENGTH OF FIELD
      IF(KFLDS(N)-2) 165,167,170
  165 KLIM=1
      LEN=LEN1
      GO TO 175
  167 KLIM=KMAX
      LEN=LEN2
      GO TO 175
  170 KLIM=KMAXP1
      LEN=LEN3
C     ****                   SET POINTERS TO INPUT AND OUTPUT BUFFERS
  175 NXK   = NJP2+NDEXA(N+1)
      NYK   = NJPUT+NDEXA(N+1)
C     ****                   CHECK IF FIELD IS TO HAVE SIGN CHANGED
      IF(KPOLE(N)) 180,195,185
C     ****                   CHANGE SIGN
  180 DO 184 K=1,KLIM
      ID=IMAXH+1
      DO 182 I=1,ILIM
      F(ID+2,NYK)=-F(I+2,NXK)
      F(I   ,NYK)=-F(ID ,NXK)
      ID=ID+1
  182 CONTINUE
      NXK=NXK+1
      NYK=NYK+1
  184 CONTINUE
      GO TO 195
C     ****                   DO NOT CHANGE SIGN
  185 DO 188 K=1,KLIM
      ID=IMAXH+1
      DO 187 I=1,ILIM
      F(ID+2,NYK)=F(I+2,NXK)
      F(I   ,NYK)=F(ID ,NXK)
      ID=ID+1
  187 CONTINUE
      NXK=NXK+1
      NYK=NYK+1
  188 CONTINUE
  195 CONTINUE
      GO TO 267
C     ****                   CONSTRUCT LINES -1,0 FOR HEMISPHERIC CASE
  200 DO 265 N=1,NFLDS
C     ****                   CHECK IF FIELD IS NEEDED ACROSS EQUATOR
      IF(KEQTR(N)) 205,265,205
C     ****                   TRANSFER Z TYPES FIRST ITERATION ONLY
  205 IF(IFRST.EQ.1 .AND. KEQTR(N).EQ.2) GO TO 265
C     ****                   DETERMINE LENGTH OF FIELD
      IF(KFLDS(N)-2) 210,215,220
  210 KLIM=1
      LEN=LEN1
      GO TO 225
  215 KLIM=KMAX
      LEN=LEN2
      GO TO 225
  220 KLIM=KMAXP1
      LEN=LEN3
C     ****                   SET POINTERS TO FIELDS IN SCM BUFFERS
  225 NJP2K = NJP2+NDEXA(N+1)
      NJPUTK= NJPUT+NDEXA(N+1)
C     ****                   CHECK IF FIELD IS TO HAVE SIGN CHANGED
      IF(KEQTR(N)) 230,265,260
C     ****                   CHANGE SIGN AND STORE IN PROPER SCM BUFFER
  230 DO 235 I=1,LEN
  235 F(I,NJPUTK)=-1.*F(I,NJP2K)
      GO TO 265
C     ****                   DO NOT CHANGE SIGN AND STORE IN PROPER
C     ****                   SCM BUFFER.
  260 DO 263 I=1,LEN
  263 F(I,NJPUTK)=F(I,NJP2K)
  265 CONTINUE
  267 CONTINUE
C     ****                   CHECK FOR COMPLETION OF BUFFERING
      CALL UNLCCK(IUNIN,NSIN,LEN)
      IF(LEN.EQ.LBDSKI.AND.NSIN.EQ.0) GO TO 270
      NUM=270
      GO TO 4000
  270 CONTINUE
C     ****                   CHANGE BUFFER POINTERS SO THAT LINE J IS -2
C     ****                   SO THAT INITIALIZATION CALLS TO HDIF
C     ****                   CAN BE MADE.
      DO 280 J=1,5
  280 CALL FLIP
C     ****                   ITITIALIZE HORIZONTAL DIFFUSION QUANTITIES
      DO 290 JJ=1,3
      CALL FLIP
      J=JJ-3
  290 CALL HDIF
C     ****                   GENERAL LINE LOOP (LATITUDINAL VARIATION)
      DO 3000 J=1,JMAX
      RACS=1./(C(51)*CS(J))
      CALL FLIP
C     ****                   START IN NEXT LINE NORTH
      IF(J.GE.JMAXM2) GO TO 310
      CALL UNLCRD(IUNIN,NJIN,LBDSKI,NSIN)
C     ****                   START OUT PREVIOUS LINE SOUTH
  310 IF(J.EQ.1) GO TO 315
      CALL UNLCWT(IUNOT,NJOT,LBDSKI,NSOT)
C     ****                   ADD DIAGNOSTIC QUANTITIES ONTO LINE J+2
  315 CONTINUE
      IF(J.GT.JMAXM2)GO TO 405
      CALL ADDIAG
  405 CONTINUE
      IF(J.EQ.JMAX) IFRST=1
C     ******************************************************************
C     ****
C     ****                   DYNAMICS
C     ****
C     ******************************************************************
C     ****
C     ****    CALL TO BNDCMP SHOULD LEAD CALLS IN DYNAMICS SECTION OF
C     ****      ADVNCE TO ENSURE THAT MATRICES B(ZIMXP,2,2) AND VECTORS
C     ****      FB(ZIMXP,2) ARE AVAILABLE WHEN NEEDED
C     ****
      CALL BNDCMP
C     ****
C     ****    CALCULATE W AT N*DT
C     ****
      CALL SWDOT
C     ****
C     ****     CALCULATE ION DRIFT VELOCITY
C     ****
      CALL VDRIFT2
C     ****     CALCULATE CP, KT, AND KM
      CALL CPKTKM
C     ****      CLEAR NQO2P, NQOP, NQN2P, NQNOP, NQNP
      NQPK=NQO2P
      DO 1 I=1,5*LEN3
      F(I,NQPK)=0.
    1 CONTINUE
C     ****     CALCULATE COLUMN NUMBER DENSITIES
      CALL CHAPMN
C     ****     CALCULATE SOLAR HEATING
C     ****     CALCULATE REACTION RATES
      CALL RATES
      CALL QRJ
C     ****     CALCULATE BACKGROUND (NIGHT-TIME) IONIZATION
      CALL QINITE
C     ****     ADD IONIZATION DUE TO SOLAR X-RAYS
      CALL XRAY
C     ****     GENERATE AURORAL FIELDS, UI, VI, WI.  FIELDS RETURNED
C     ****     IN T1 THRU T7
      CALL HEELIS
C     ****     CALCULATE AURORAL ADDITIONS TO IONIZATION RATES
      CALL ORORA
C     ****     SOLVE FOR N(O+)
C     ****     SET FLUX AT UPPER BOUNDARY
      CALL OPFLUX
      CALL OPLUS
C     ****     COMPUTE ELECTRON DENSITY, N(E), AND ION NUMBER DENSITIES,
C     ****     N(O2+), N(N2+), N(NO+), N(N+)
      CALL ELDEN
C     ****     ION-DRAG PARAMETERS, LXX, LYY, LXY (IN DIPOLE COORDINATE
C     ****     SYSTEM)
      CALL LAMDAS
C     ****     ROTATE LAMDAS TO GEOGRAPHIC SYSTEM
      CALL ROTATE
C     ****     ELECTRON PRODUCTION OF N4S AND N2D
      CALL QTIEFF
C     ****     ODD NITROGEN CHEMISTRY
      CALL CMPN2D
      CALL CMPN4S
      CALL CMPNO
      CALL CMPO2O
C     ****     HEATING AND O2 DISSOCIATION DUE TO ODD NITROGEN AND ION
C     ****     CHEMISTRY
      CALL QJNNO
      CALL QJION
C     ****     CALCULATE ELECTRON AND ION TEMPERATURES
      CALL SETTEI
C     ****     CALCULATE IMPLICIT AND EXPLICIT COOLING TERMS IN NWTI
C     ****     AND NWTE
      CALL NEWTON
C     ****     ADD IN O3P COOLING
      CALL NEWTO3P
C     ****     CALCULATE FLF, FPH, NQDH, NKMH, NPSDH, NPSDH2
C     ****     AT (N-1)*DT
      CALL HDIF
C     ****     ADVANCE T, U, V BY ONE TIME STEP
      CALL DT
      CALL DUV
C     ****     COMPOSITION OF MAJOR SPECIES
      CALL COMP
C     ******************************************************************
C     ****
C     ****                   END OF DYNAMICS
C     ****
C     ******************************************************************
C     ****                   CHECK IF ON LINE JMAXM2 OR JMAXM1
      IF(J.LT.JMAXM2) GO TO 2980
      IF(J.EQ.JMAX) GO TO 2990
C     ****                   REARRANGE LINE NJGET TO GET LINE NJIN
      NJGET=NJP2
      IF(J.EQ.JMAXM1) NJGET=NJ
      ILIM=IMAXH+2
      DO 2975 N=1,NFLDS
C     ****                   CHECK IF FIELD IS NEEDED ACROSS POLE
      IF(KPOLE(N)) 2916,2975,2916
C     ****                   TRANSFER Z TYPES FIRST ITERATION ONLY
 2916 IF(IFRST.EQ.1.AND.KPOLE(N).EQ.2) GO TO 2975
C     ****                   DETERMINE LENGTH OF FIELD
      IF(KFLDS(N)-2) 2920,2925,2930
 2920 KLIM=1
      LEN=LEN1
      GO TO 2935
 2925 KLIM=KMAX
      LEN=LEN2
      GO TO 2935
 2930 KLIM=KMAXP1
      LEN=LEN3
C     ****                   BRING IN FIELD FROM SCM BUFFER
 2935 NXK   = NJGET+NDEXA(N+1)
      NYK   = NJIN+NDEXA(N+1)
C     ****                   CHECK IF FIELD IS TO HAVE SIGN CHANGED
      IF(KPOLE(N)) 2940,2975,2955
C     ****                   CHANGE SIGN
 2940 DO 2950 K=1,KLIM
      ID=IMAXH+1
      DO 2945 I=1,ILIM
      F(ID+2,NYK)=-F(I+2,NXK)
      F(I   ,NYK)=-F(ID ,NXK)
      ID=ID+1
 2945 CONTINUE
      NXK=NXK+1
      NYK=NYK+1
 2950 CONTINUE
      GO TO 2975
C     ****                   DO NOT CHANGE SIGN
 2955 DO 2965 K=1,KLIM
      ID=IMAXH+1
      DO 2960 I=1,ILIM
      F(ID+2,NYK)=F(I+2,NXK)
      F(I   ,NYK)=F(ID ,NXK)
      ID=ID+1
 2960 CONTINUE
      NXK=NXK+1
      NYK=NYK+1
 2965 CONTINUE
C     ****                   TRANSFER REARRANGED FIELD
C     ****                   TO CORRECT LOCATION IN SCM BUFFER
 2975 CONTINUE
      GO TO 2990
C     ****                   CHECK FOR COMPLETION OF BUFFERING
 2980 CALL UNLCCK(IUNIN,NSIN,LEN)
      IF(LEN.EQ.LBDSKI.AND.NSIN.EQ.0) GO TO 2990
      NUM=2990
      GO TO 4000
 2990 CONTINUE
      IF(J.EQ.1)GO TO 3000
      CALL UNLCCK(IUNOT,NSOT,LEN)
      IF(LEN.EQ.LBDSKI.AND.NSOT.EQ.0) GO TO 3000
      NUM=3000
      GO TO 4000
 3000 CONTINUE
C     ****                   TRANSFER LINE JMAX TO DISK
      CALL FLIP
      CALL UNLCWT(IUNOT,NJOT,LBDSKI,NSOT)
      CALL UNLCCK(IUNOT,NSOT,LEN)
      IF(LEN.EQ.LBDSKI.AND.NSOT.EQ.0) GO TO 3100
      NUM=3100
      GO TO 4000
 3100 ITEMP=IUNIN
      IUNIN=IUNOT
      IUNOT=ITEMP
      REWIND IUNIN
      REWIND IUNOT
c     JTOTN=JTOTN*C(81)*C(55)*C(3)*C(53)*C(1)*C(2)
c     IF(NHEMI.EQ.0)JTOTS=JTOTS*C(81)*C(55)*C(3)*C(53)*C(1)*C(2)
C     ****
C     ****     SET UP AND SOLVE DYNAMO EQUATIONS
C     ****
      CALL DYN
C     ****                   CHECK FOR EXIT
C     ****                   CHECK FOR SAVE TAPE WRITE
 3130 IF(ISAV.EQ.0) GO TO 3150
C     ****                   CHECK FOR HISTORY TAPE WRITE
      IF(IHIS.NE.0) GO TO 3140
      IF(ITER.GE.NSTP)ISAV=0
      GO TO 3150
C     ****                   CHECK FOR END OF RUN
 3140 IF(ITER.LT.NSTP) GO TO 100
      ISAV=0
 3150 TIME = SECOND(X)-TIME
      TIMET = TIMEF(X)-TIMET
      WRITE(6,3170) TIMET
 3170 FORMAT(1H ,17X,*TIME FROM TIMEF CALL = * E12.5)
      WRITE(6,3160) ITER,TIME
 3160 FORMAT(1H ,17X,*ADVANCE EXITED  AT ITER = *I6,10X,
     A       *TIME IN INNER LOOP = *E12.5,* SECONDS*)
      RETURN
 4000 WRITE(6,4005) ITER
 4005 FORMAT(*1ERROR EXIT AT ITER =*I6,5X,*RESTART*)
      WRITE(6,4010)NUM,LEN,NSIN,NSOT,NSIN1,NSIN2,J
 4010 FORMAT(*0ADVANCE- I/O PROBLEM AT NUM=*I5,2X,*LEN=*I5,2X,*NSIN=*I5,
     A2X,*NSOT=*I5,2X,*NSIN1=*I5,2X,*NSIN2=*I5,2X,*J=*I2)
      CALL EXIT
      END
C
