!=============================================================================
!
! $Id: uci_control.F90,v 1.2 2007/05/11 01:29:39 kouatch Exp $
!
! CODE DEVELOPER
! Deanne Proctor, et. al., LLNL (Original code from Michael Prather, UCI)
! ddproctor@llnl.gov
!
! FILE
! uci_control.F
!
! ROUTINES
! Uci_Control
! Set_Vert_Uci
! Wind
! Exitc
! Jshift1
! Ad_To_Mass
! Jshiftn
! Stt_To_Const
! Get_Diag_1d_Dump_Flag
!
!=============================================================================
!-----------------------------------------------------------------------------
!
! ROUTINE
! Uci_Control
!
! DESCRIPTION
! This routine is the wrapper for UCI 3D transport model.
!
! This transport driver is based on the GISS GCM.
! This version is meant to be used for both trop and strat:
! GISS 4x5 / 8x10 x 9 layer GCM history tapes with 5-day averages,
! GISS 7.8x10 x 23 layer GCM history tapes with monthly convection.
!
! ARGUMENTS
! const_ncar_infile : const input file name
! met_uci_infile : met input file name (UCI format)
!c pr_diag_1d : should UCI 1D diagnostic output be written?
!c pr_diag_2d : should UCI 2D diagnostic output be written?
!c time_for_diag_init : time for UCI diagnostics initialization?
! gmi_lm : UCI vertical grid dimension
! num_time_steps : number of time steps
! nymd : year/month/day (YYYYMMDD)
! nhms : hour/min/sec (HHMMSS)
! gmi_sec : total Gmimod seconds (s)
! tdt : model time step (s)
! uci_lat_out : UCI latitudes (deg)
! press : midpoint pressures by level (mb)
! prslay : interface pressures by level (mb)
! mass : total mass of the atmosphere within each grid box
! (kg)
! const : species concentration, known at zone centers
! (mixing ratio)
!
!-----------------------------------------------------------------------------
subroutine Uci_Control &,53
!c & (const_uci_infile, met_uci_infile, pr_diag_1d, pr_diag_2d,
!c & time_for_diag_init, gmi_lm, num_time_steps, nymd, nhms,
& (const_uci_infile, met_uci_infile, gmi_lm, num_time_steps, &
& nymd, nhms, gmi_sec, tdt, uci_lat_out, press, prslay, &
& mass, const)
!xx!? implicit none
# include "netcdf.inc"
# include "gmi_dims.h"
# include "uci.h"
! ----------------------
! Argument declarations.
! ----------------------
character*(*) :: const_uci_infile
character*(*) :: met_uci_infile
!c logical :: pr_diag_1d, pr_diag_2d
!c logical :: time_for_diag_init
integer :: gmi_lm
integer :: num_time_steps
integer :: nymd, nhms
real*8 :: gmi_sec
real*8 :: tdt
real*8 :: uci_lat_out(ju1_gl:j2_gl)
real*8 :: press (k1:k2)
real*8 :: prslay(k1-1:k2)
real*8 :: mass (i1:i2, ju1:j2, k1:k2)
real*8 :: const (i1:i2, ju1:j2, k1:k2, num_species)
! ----------------------
! Variable declarations.
! ----------------------
character*78 :: titver
character*128 :: uci_read_file
character*128 :: uci_write_file
character*8 :: tlabl(10)
logical :: dump_diag
logical, save :: first_time = .true.
integer :: cdfid
integer :: dnum_time_steps
integer :: dnymd, dnhms
!xx integer :: idt
integer :: ierr_getlun
integer :: ijx
integer :: il, ij, ik, ic
integer :: gmi_rec
real*8 :: gmi_day
real*8 :: mw(num_species)
real*8 :: ad_jsh (i1:i2, ju1:j2, k1:k2)
real*8 :: const_jsh(i1:i2, ju1:j2, k1:k2, num_species)
! ====
save
! ====
data tlabl / 'INITIAL ', 'WINDS+P ', 'DYNAMICS', 'ALL-CONV', &
& '________', 'SOURCES ', 'CHEMSTRY', 'DIFF-MIX', &
& 'DIAGNOST', 'SAVETAPE' /
data titver &
& / ' CTM (GISS/Harvard/UCIrvine) p-1.0 1/95(Prather)' /
! ----------------
! Begin execution.
! ----------------
if (pr_diag) then
if (first_time) then
Write (6,*) 'Uci_Control called for the first time.'
else
Write (6,*) 'Uci_Control called, ', ntau * 3600.0d0
end if
end if
if (.not. first_time) then
! =======
goto 20
! =======
end if
!rt
!rt
!rt
!rt first time through UCI
!rt
!rt
!rt
!xx ymd=start_ymd
!xx hms=start_hms
!xx
KYY = 93
KMM = 7
KDD = 1
! CALL IDATE(KMM,KDD,KYY)
! CALL SYSTEM('date')
!---------------------input files--(___.ctm)----------------------------
!#1 CTM INITIALIZE/RESTART FILE
!xx OPEN(1,FILE='init.uci',status='OLD',form='UNFORMATTED')
!#2 ALTITUDE/LAND-FRACT FOR GCM GRID
!xx OPEN(2,FILE='landalt.uci',STATUS='OLD',form='FORMATTED')
!#3 GCM DATA TAPE (winds, etc.)
!rt
!rt OPEN(3,FILE='wind.ctm',STATUS='OLD',FORM='UNFORMATTED')
call opencdf
(met_uci_infile,cdfid)
!rt
!rt
!xx
call FI_Getlun (lun8, ierr_getlun)
!xx OPEN(8,FILE='input.uci',STATUS='OLD',form='FORMATTED')
OPEN(unit=lun8,FILE='input.uci',STATUS='OLD',form='FORMATTED')
!xx
!#9 TRACER SETUP DATA CONC/SOURCE/CHEM
call FI_Getlun (lun9, ierr_getlun)
!xx OPEN(9,FILE=const_uci_infile,STATUS='OLD',form='FORMATTED')
OPEN(unit=lun9,FILE=const_uci_infile, &
& STATUS='OLD',form='FORMATTED')
!--------------------output files--(ctm.___)----------------------------
!#10 (FORMATTED DIAGNOSTIC OUTPUT, 80 col)
call FI_Getlun (lun10, ierr_getlun)
!xx OPEN(10,FILE='uci.pch',form='FORMATTED')
OPEN(unit=lun10,FILE='uci.pch',form='FORMATTED')
!#11 (CTM DAILY RESTART FILE, UNFORMATTED)
!xx OPEN(11,FILE='uci.sav',form='UNFORMATTED')
!#12 (UNFORMATTED DIAGNOSTICS)
call FI_Getlun (lun12, ierr_getlun)
!xx OPEN(12,FILE='uci.unf',form='UNFORMATTED')
OPEN(unit=lun12,FILE='uci.unf',form='UNFORMATTED')
!xx
call FI_Getlun (lun13, ierr_getlun)
!xx OPEN(13,FILE='output.uci',form='FORMATTED')
OPEN(unit=lun13,FILE='output.uci',form='FORMATTED')
!xx
!-----------------------------------------------------------------------
!xx
WRITE (lun13,909)
WRITE (lun13,900) TITVER,TITCOM,KDD,KMM,KYY
WRITE (lun13,909)
!xx
!rt
!rt
!rt
!rt name UCI NetCDF restart file
!rt
UCI_read_file = "UCI_read"
UCI_write_file = "UCI_write"
!rt
!rt
!rt
!xx CALL INPUT
CALL INPUT
(tdt,UCI_read_file)
!xx
ijx = 1
do ij = ju1_gl, j2_gl
uci_lat_out(ij) = uci_lat(ijx)
ijx = ijx + 1
end do
!xx
!rt
!xx
!
CALL INPTR
!
!
CALL GRDSET
WRITE(lun10,981) (XLABEL(K),K=1,15),NTRACE,IM,JM,LM,I0,J0,IYEAR
WRITE(lun12) (XLABEL(K),K=1,15),NTRACE,IM,JM,LM,I0,J0,IYEAR
!xx CLOSE (1)
!xx CLOSE (2)
CLOSE (lun9)
KDACC = 0
KWACC = 0
DO K=1,KACCZ
SAVACC(K) = 0.0
ENDDO
!-----------------------------------------------------------------------
!-------Begin CTM calculation at time = NTAU hrs relative to IYEAR-----
!-------NTAU = ABSOLUTE TIME (I*4), ITAU = RELATIVE TIME (0 - <8760 HRS)
!-------TAU (R*4) = NTAU, calculation not resolved finer than 1 hour
NTAU = (TAU+0.1)
NTAUE = (TAUE+0.1)
IF (NTAU .GE. NTAUE) GOTO 99
ITAU = MOD(NTAU,8760)
!-----------------------------------------------------------------------
!xx CALL WIND (0,ITAU,ITAUW)
CALL WIND
(0,ITAU,ITAUW,cdfid)
!-----------------------------------------------------------------------
IF (NDYN*(NREAD/NDYN) .NE. NREAD ) GOTO 97
NTAU = (NTAU/NREAD)*NREAD
CALL DAILY
LPRT = ND70.GT.0
WRITE (lun13,901) TLABL(1),NTAU,JDATE,JMONTH,JYEAR
NTAU00 = NTAU
!rt
IF (LCONT) THEN
NTAU00 = TAU0
TAU0 = TAU
ENDIF
!
if (LCONT) Write (lun13,905) NTAU,NTAU00,TAU0,TAU
if (.not.LCONT) Write (lun13,906) NTAU,NTAU00,tAU0,TAU
!rt
!
NTSAVE = -9999
!-----------------------------------------------------------------------
!
CALL AIRSET
!
CALL TRACE0
!
!-----------------------------------------------------------------------
!rt
!rt IF (LCONT) THEN
!rt CALL DIAG0(9)
!rt ELSE
CALL DIAG0
(1)
! should really pass mw in from Gmimod? !!
do k=1,num_species
mw(k)=28.97/tcvv(k)
end do
!xx call Get_Diag_1d_Prnt_Flag (pr_diag_1d)
!xx call Get_Diag_2d_Prnt_Flag (pr_diag_2d)
!c if(pr_diag_1D) then
! need j shifts in const and ad being fed to diag0_gmi
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (-1,const_jsh,ad_jsh)
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (1,const_jsh,ad_jsh)
!c endif
!c if(pr_diag_2d) then
! note there's no shift for Gmi_Diag_Save so should
! rename const_jsh,ad_jsh when get rid of 1d diagnostics
!c call stt_to_const(stt,const_jsh)
!c call ad_to_mass(ad,ad_jsh)
!c call Gmi_Diag_Save
!c & (1,nymd,nhms,gmi_sec,mw,ad_jsh,const_jsh)
!c endif
!rt ENDIF
!rt
MONTH0 = 0
!rt
!rt
!rt
! set month_co2
! to first of year relative to IYEAR input value
! month_co2 = 0, 12, 24, ... ,180
!
if (run_type .eq. " CO2") then
month_co2 = (JYEAR - IYEAR)*12
endif
!rt
!rt
!rt
DTWIND = 3600.D0 *REAL(NDYN)
DTCONV = 3600.D0 *REAL(NCONV)
!xx
! =================
call Set_Vert_Uci
&
! =================
& (k1, k2, gmi_lm, ltm, ptop, sige, plev, press, prslay)
first_time = .false.
!xx
!-----------------------------------------------------------------------
!------------------------MAIN LOOP--------------------------------------
!------------increment by NDYN hours as basic step----------------------
!-----------------------------------------------------------------------
! ========
20 continue
! ========
!rt
!rt
!rt
if (edit_type .eq. "debug") then
Write (lun13,914) NDYN,ITAU,ITAUW,NTAU
endif
!rt
!rt
!rt
!
!---------do not read new winds if NDYN < NREAD and keep using same winds
IF (ITAU+NDYN.GT.ITAUW) THEN
CALL WIND
(1,ITAU,ITAUW,cdfid)
CALL DYN0
(1)
! in principle only have to fill mass array when you want to
! have pr_ascii1_c true, but...
! (also assumes only change mass when get new windfield+dyn0
call jshift1
(ad,mass)
! could be using mass below instead of ad_jsh?
CALL DIAG0
(9)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (9,const_jsh,ad_jsh)
!c endif
IF (LPRT) WRITE(lun13,901) TLABL(1),NTAU,JDATE,JMONTH,JYEAR
ENDIF
!----------CHECK FOR DAILY / MONTHLY HOUSEKEEPING
IF (MOD(NTAU,24).EQ.0) THEN
CALL DAILY
ENDIF
IF (MONTH .NE. MONTH0) THEN
MONTH0 = MONTH
!rt
!rt
!rt
if ((run_type .eq. " N2O") .or. &
& (run_type .eq. " F11")) then
!rt
CALL STRATL
!rt
endif
!rt
!rt
!rt
if (run_type .eq. " CO2") then
month_co2 = month_co2 + 1
!rt
if (edit_type .eq. "debug") then
Write (lun13,916) month0,month_co2,tparm_co2(1,10,month_co2,1)
endif
endif
!rt
!rt
!rt
ENDIF
!rt
!rt
!rt
if (run_type .eq. " CO2") then
!rt
CALL STRATL_co2
!rt
endif
!rt
!rt
!rt
if (run_type .eq. " SF6") then
!rt
CALL STRATL_sf6
!rt
endif
!rt
!rt
!rt
!-----------------------------------------------------------------------
!-----Begin sequence of operations on the tracer distribution:
!----- currently the order is SOURCE/CHEM/ADVECT/PFILT/DIFF/CONV/DIAG
!----- N.B. the order may affect diagnostics, esp. for "rapid change"
!----- species with large source/sink terms.
!-----This loop calculates changes in tracer from NTAU to NTAU+NDYN
!-----------------------------------------------------------------------
!xx call Get_Diag_Init_Flag(time_for_diag_init)
!c if(time_for_diag_init) then
!c call stt_to_const(stt,const_jsh)
!c call ad_to_mass(ad,ad_jsh)
!c call Gmi_Diag_Save
!c & (1,nymd,nhms,gmi_sec,mw,ad_jsh,const_jsh)
!c endif
CALL DIAG0
(9)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (9,const_jsh,ad_jsh)
!c endif
!----------SOURCE TERMS:------------------------------------------------
IF (MOD(NTAU,NSRCE).EQ.0) THEN
!rt
!rt
!rt
if ((run_type .eq. " N2O") .or. &
& (run_type .eq. " F11")) then
!rt
CALL SOURCE
!rt
endif
!rt
!rt
!rt
CALL DIAG0
(6)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (3,const_jsh,ad_jsh)
!c endif
!c if(pr_diag_2d) then
!c call stt_to_const(stt,const_jsh)
!c call ad_to_mass(ad,ad_jsh)
!c call Gmi_Diag_Save
!c & (6,nymd,nhms,gmi_sec,mw,ad_jsh,const_jsh)
!c endif
IF (LPRT) WRITE (lun13,901) TLABL(6),NTAU,JDATE,JMONTH,JYEAR
ENDIF
!----------CHEMISTRY:---------------------------------------------------
IF (MOD(NTAU,NCHEM).EQ.0) THEN
!
CALL CHEM
!
CALL DIAG0
(7)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (4,const_jsh,ad_jsh)
!c endif
!c if(pr_diag_2d) then
!c call stt_to_const(stt,const_jsh)
!c call ad_to_mass(ad,ad_jsh)
!c call Gmi_Diag_Save
!c & (3,nymd,nhms,gmi_sec,mw,ad_jsh,const_jsh)
!c endif
IF (LPRT) WRITE (lun13,901) TLABL(7),NTAU,JDATE,JMONTH,JYEAR
ENDIF
!----------ADVECT TRACERS BY THE WINDS: U-V-W pattern based on 8x10 winds
IF (MOD(NTAU,NDYN).EQ.0) THEN
DTALFA = DTWIND/4.
DTBETA = DTWIND/3.
DTGAMA = DTWIND/2.
!
CALL DYN2U
(DTALFA)
!
CALL DYN2V
(DTBETA)
!
CALL DYN2W
(DTGAMA)
!
CALL DYN2U
(DTALFA)
!
CALL DYN2V
(DTBETA)
!
CALL DYN2U
(DTALFA)
!
CALL DYN2W
(DTGAMA)
!
CALL DYN2V
(DTBETA)
!
CALL DYN2U
(DTALFA)
!
CALL DIAG0
(3)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (5,const_jsh,ad_jsh)
!c endif
!c if(pr_diag_2d) then
!c call stt_to_const(stt,const_jsh)
!c call ad_to_mass(ad,ad_jsh)
!c call Gmi_Diag_Save
!c & (5,nymd,nhms,gmi_sec,mw,ad_jsh,const_jsh,)
!c endif
!
!
! WRITE(lun13,'(11HDYN2. LD42=,I3)') LD42
! DO J = 14,10,-1
! WRITE(lun13,'(1P10E8.1)') (AIJ(I,J,1+ND42),I=1,10)
! ENDDO
IF (LPRT) WRITE (lun13,901) TLABL(3),NTAU,JDATE,JMONTH,JYEAR
ENDIF
!----------DIFFUSION OF TRACER:-----------------------------------------
IF (MOD(NTAU,NDYN).EQ.0) THEN
CALL DIFFU
CALL DIFFV
CALL DIFFW
CALL DIAG0
(8)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (6,const_jsh,ad_jsh)
!c endif
!c if(pr_diag_2d) then
!c call stt_to_const(stt,const_jsh)
!c call ad_to_mass(ad,ad_jsh)
!c call Gmi_Diag_Save
!c & (7,nymd,nhms,gmi_sec,mw,ad_jsh,const_jsh)
!c endif
IF (LPRT) WRITE (lun13,901) TLABL(8),NTAU,JDATE,JMONTH,JYEAR
ENDIF
!----------CONVECTION: includes boundary layer, moist cumulus, large-scale
IF (MOD(NTAU,NCONV).EQ.0) THEN
CALL CNVWET
(DTCONV)
CALL DIAG0
(4)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (7,const_jsh,ad_jsh)
!c endif
!c if(pr_diag_2d) then
!c call stt_to_const(stt,const_jsh)
!c call ad_to_mass(ad,ad_jsh)
!c call Gmi_Diag_Save
!c & (8,nymd,nhms,gmi_sec,mw,ad_jsh,const_jsh)
!c endif
CALL CNVDRY
(DTCONV)
CALL DIAG0
(5)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (2,const_jsh,ad_jsh)
!c endif
!c if(pr_diag_2d) then
!c call stt_to_const(stt,const_jsh)
!c call ad_to_mass(ad,ad_jsh)
!c call Gmi_Diag_Save
!c & (8,nymd,nhms,gmi_sec,mw,ad_jsh,const_jsh)
!c endif
IF (LPRT) WRITE (lun13,901) TLABL(4),NTAU,JDATE,JMONTH,JYEAR
! WRITE(lun13,'(11HCONV. LD42=,I3)') LD42
! DO J = 14,10,-1
! WRITE(lun13,'(1P10E8.1)') (AIJ(I,J,1+ND42),I=1,10)
! ENDDO
ENDIF
!-----------------------------------------------------------------------
!----------COMPLETED A TIME STEP: update NTAU by NDYN (basic time step)
!-----------------------------------------------------------------------
NTAU = NTAU+NDYN
!rt
!rt
!rt
ITAU = MOD(NTAU-1,8760)+1
!rt
!rt
!rt
TAU = NTAU
IDAY = 1+NTAU/24
TOFDAY = TAU-(IDAY-1)*24
!----------DIAGNOSTICS: ACCUMULATE SOME DIAGNOSTICS AT END OF TIME STEP
IF (MOD(NTAU,NDIAG).EQ.0) THEN
CALL DIAG1
IF (LPRT) WRITE (lun13,902) NTAU
ENDIF
IF (MOD(NTAU,NINST).EQ.0) CALL DIAG6
!----------DIAGNOSTICS: print / zero-out accumulating arrays:-----------
IF (MOD(NTAU,24).EQ.0) THEN
CALL DAILY
! njday is table (input) with value 3 on days when dump wanted
JDO = NJDAY(JDAY)
LZERO = MOD(JDO,2)
LPRNT = MOD(JDO/2,2)
!----------WRITE DIAGNOSTICS to UNIT=6 (printout) and UNIT=10 (pch)-----
IF (LPRNT.EQ.1) THEN
CALL DIAG0
(9)
!c if(pr_diag_1d) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (9,const_jsh,ad_jsh)
!c endif
CALL DIAG0
(99)
CALL DIAG2
CALL DIAG4
CALL DIAG3
CALL DIAG5
!---------may need: ENDFILE 10 / BACKSPACE 10 to clear buffer if crash
ENDIF
!---------ZERO OUT AND RESET ACCUMLATING ARRAYS-------------------------
IF (LZERO.EQ.1) THEN
TAU0 = TAU
NTAU0 = NTAU
IDAY0 = IDAY
TOFDY0= TOFDAY
JDATE0= JDATE
JMNTH0= JMONTH
JYEAR0= JYEAR
KDACC = 0
KWACC = 0
DO 50 K=1,KACCZ
50 SAVACC(K) = 0.0
CALL DIAG0
(1)
ENDIF
ENDIF
! GMIMOD DIAGNOSTICS DUMP?
!c if(pr_diag_1d) then
! since Gmimod ticks clock then checks time, in order to sync
! diagnostics need to get dump flag for incremented ymd
!c dnum_time_steps=num_time_steps
!c dgmi_sec=gmi_sec
!c dnymd=nymd
!c dnhms=nhms
!xx idt=Nint(tdt)
!xx call Gmi_Tick(dnum_time_steps,dgmi_sec,dnymd,dnhms,idt,tdt)
!c call Gmi_Tick
!c & (dnum_time_steps,dgmi_sec,dnymd,dnhms,tdt)
!c call Get_Diag_1d_Dump_Flag (dnymd,dump_diag)
!c if(dump_diag) then
!c call jshiftn(stt,const_jsh)
!c call jshift1(ad,ad_jsh)
!c call diag0_gmi (99,const_jsh,ad_jsh)
!c call diag0_gmi (-1,const_jsh,ad_jsh)
!c call diag0_gmi ( 1,const_jsh,ad_jsh)
!c endif
!c endif
! END GMIMOD DIAGNOSTICS DUMP
!----------SAVE COMPLETE CTM CALCULATION: for restart or re-init.-------
!rt IF (MOD(NTAU,NWRITE).EQ.0) THEN
!xx REWIND 11
!xx WRITE(11) TAU,STT,SUT,SVT,SWT,AD,SAVLAB,SAVDAT,SAVACC,SAVTRA
!xx IF (LSOM) WRITE(11) TAU,SUU,SVV,SWW,SUV,SUW,SVW
!rt
!rt
!rt
!rt Write NetCDF restart file every NWRITE hours
!rt
if (mod(NTAU,NWRITE).eq.0) then
gmi_rec=1
gmi_day=NTAU/24
N=1
!rt
call Write_Mpra2
(UCI_write_file, gmi_rec, gmi_day, &
& STT, ad, SUT, SVT, SWT, SUU, SVV, SWW, SUV, SUW, SVW, &
& IM, JM, LM, N, &
& NTAU00, TAU)
!rt
if (LCONT) Write (lun13,907) NTAU,NTAU00,TAU0,TAU
if (.not.LCONT) Write (lun13,908) NTAU,NTAU00,tAU0,TAU
!rt
!rt
!rt
NTSAVE = NTAU
IF (LPRT) WRITE (lun13,901) TLABL(10),NTSAVE,JDATE,JMONTH,JYEAR
!xx REWIND 11
ENDIF
!-----------------------------------------------------------------------
!----------TEST FOR TERMINATION OF RUN
!-----------------------------------------------------------------------
!xx
!xx IF (NTAU.LT.NTAUE) GOTO 20
! ========
goto 888
! ========
!xx
!-----------------------------------------------------------------------
!----------END OF MAIN LOOP, CLOSE DOWN AND SAVE RESULTS----------------
!-----------------------------------------------------------------------
!xx IF (NTAU.NE.NTSAVE) THEN
!xx Write(lun13,*) ' begin to write ctm.sav file'
!xx REWIND 11
!xx WRITE(11) TAU,STT,SUT,SVT,SWT,AD,SAVLAB,SAVDAT,SAVACC,SAVTRA
!xx IF (LSOM) WRITE(11) TAU,SUU,SVV,SWW,SUV,SUW,SVW
!xx REWIND 11
!xx WRITE (lun13,901) TLABL(10),NTAU,JDATE,JMONTH,JYEAR
!xx ENDIF
!xx GOTO 99
!-----------------------------------------------------------------------
!---------PROBLEM WITH MONTHLY FILE-------------------------------------
!xx96 WRITE (lun13,961) NTAU,IMON
!xx WRITE (lun13,961) NTAU
!xx GOTO 99
!---------INVALID TIME STEPS--------------------------------------------
97 WRITE(lun13,971) NDYN,NREAD
GOTO 99
!-------------------------FINAL-----------------------------------------
99 CONTINUE
WRITE (lun13,909)
WRITE (lun13,900) TITVER,TITCOM,KDD,KMM,KYY
WRITE (lun13,910) NTAU00,NTAU,JDATE,JMONTH,JYEAR
! CALL SYSTEM('date')
WRITE (lun13,909)
STOP
!xx
888 continue
!rt
!rt
!rt
do 890 ic = 1, num_species
do 890 ik = k1, k2
do 890 ij = ju1, j2
do 890 il = i1, i2
const(il,ij,ik,ic) = stt(il,ij+1,ik,ic)
const(il,ij,ik,ic) = TCVV(ic) * const(il,ij,ik,ic)
const(il,ij,ik,ic) = const(il,ij,ik,ic)/ad(il,ij+1,ik)
!rt
890 continue
!rt
!rt
!rt
return
!xx
!-----------------------------------------------------------------------
900 FORMAT(1X,A80/' COMMON BLOCKS: ',A40,' dd/mm/yy=',3I3)
901 FORMAT(1X,A8,' TAU=',I10,I6,A4,I5)
902 FORMAT(' DIAGNOSTICS ACCUMULATED TAU=',I10,10I4)
!c904 FORMAT(' READ 5-DAY/MONTHLY STATS: ',3I10)
909 FORMAT(1X,32('****'))
910 FORMAT (' CTM finish: NTAU0/NTAU=',2I10,I6,A5,I5)
!c911 FORMAT(' END OF GCM WINDS AT TAUW=',F10.1)
!c912 FORMAT(' REWOUND GCM WIND FILE',2I10)
!c921 FORMAT(' INCORRECT FILE??? NTAU,NTAUW,NTAUC,NREAD=',4I15)
!c931 FORMAT(' =====ERROR ON WIND/CTM FILES')
!c951 FORMAT(' =====UNEXPECTED EOF ON WINDS',F10.1)
!c961 FORMAT(' =====ERROR ON MONTHLY FILES',2I10)
961 FORMAT(' =====ERROR ON MONTHLY FILES',I10)
971 FORMAT(' =====ERROR ON NDYN/NREAD:',2I10)
981 FORMAT(15A4,6I2,I8)
!rt
!rt
!rt
905 format (1x,"restart",2x,"NTAU=",i8,1x,"NTAU00=",i8, &
& 1x,"TAU0=",e12.5,1x,"TAU=",e12.5)
906 format (1x," start",2x,"NTAU=",i8,1x,"NTAU00=",i8, &
& 1x,"TAU0=",e12.5,1x,"TAU=",e12.5)
907 format (1x,"restart (write)",2x,"NTAU=",i8,1x,"NTAU00=",i8, &
& 1x,"TAU0=",e12.5,1x,"TAU=",e12.5)
908 format (1x," start (write)",2x,"NTAU=",i8,1x,"NTAU00=",i8, &
& 1x,"TAU0=",e12.5,1x,"TAU=",e12.5)
914 format (1x,"MPRA - CONTROL",1x,"NDYN=",i2,1x,"ITAU=",i8, &
& 1x,"ITAUW=",i8,1x,"NTAU=",i8)
916 format (1x,"CALL STRATL",1x,"month0=",i3,1x,"month_co2=",i3, &
& 1x,"tparm_co2=",e12.5)
!rt
!rt
!rt
END
!-----------------------------------------------------------------------------
!
! ROUTINE
! Set_Uci_Vert
!
! DESCRIPTION
! This routine sets Gmimod press & prslay for 2D diagnostics (this is an
! extra layer of stuff since 2D diagnostics were originally written for
! UCI, but....).
!
! ARGUMENTS
! k1, k2 : Gmimod vertical grid range
! gmi_lm : UCI vertical grid dimension
! ltm : top level still defined by sigma coordinates at lower boundary
! ptop : pressure at trop/strat interface (mb)
! sige : UCI sigma,
! sige = (p-ptop)/(psf-ptop) in strat,
! sige = p/ptop in trop?
! plev_uci : midpoint pressures by level (1) &
! interface pressures by level (2) (mb)
! press : midpoint pressures by level (mb)
! prslay : interface pressures by level (mb)
!
!-----------------------------------------------------------------------------
subroutine Set_Vert_Uci & 1
& (k1, k2, gmi_lm, ltm, ptop, sige, plev_uci, press, prslay)
implicit none
!c# include "gmi_diag_constants_uci.h"
! ----------------------
! Argument declarations.
! ----------------------
integer k1, k2
integer gmi_lm
integer ltm
!c? For tmp fix.
real*8 ptop
real*8 sige(*)
real*8 plev_uci(22,2)
real*8 press (k1:k2)
real*8 prslay(k1-1:k2)
! ----------------------
! Variable declarations.
! ----------------------
integer ik
! ----------------
! Begin execution.
! ----------------
press(k1:k2) = plev_uci(1:k2,1)
do ik = 1, ltm+1
prslay(ik-1) = sige(ik) * (PSF - ptop) + ptop
end do
do ik = ltm+2, gmi_lm+1
prslay(ik-1) = sige(ik) * ptop
end do
return
end
!xx SUBROUTINE WIND(ISW,ITAU,ITAUW)
SUBROUTINE WIND(ISW,ITAU,ITAUW,cdfid) 2,9
!
!----ISW=0 = locate pressure fields for initialization of dry-air mass
!----ISW>0 = read in next wind, pressure, temperature, q fields.
!----ITAU = starting time (hours, 0-8759) of initial calculation,
!---- used to position file (ONLY used for init: ISW=0)
!----ITAUW = hours (1-8760) of last read wind record, i.e.,
!---- pressure at ITAUW, winds & q averaged over ITAUW-NREAD to ITAUW
!-----------------------------------------------------------------------
# include "uci.h"
!xx
integer cdfid
!xx
!rt
!rt
!rt
integer ihyper
integer number
!rt
!rt
!rt
!xx$$
save
!xx$$
!-----------------------------------------------------------------------
IF (ISW .GT. 0) GOTO 41
!------------locate first wind/press/etc field--------------------------
! NTBEG/NTEND = START/END HOURS OF WIND TAPE
! NAVG = # HOURS FOR AVERAGE STATS ***NOT USED NOW***
! NREAD = # HOURS FOR WINDS (4 OR 8 HOURS USUALLY)
! II0/IIM = GRID OFFSET/SIZE IN THE W -> E DIRECTION
! JJ0/JJM = S -> N
! LL0/LLM = VERTICAL
!*********THE OPTION (II0,JJ0,LL0) IS NO LONGER ALLOWED,****************
!*********MUST HAVE GLOBAL GRID OF WINDS: II0=JJ0=LL0=0****************
!-----------------------------------------------------------------------
!rt READ(3) NTBEG,NTEND,NAVG,NREAD,I00,IIM,J00,JJM,L00,LLM
!rt
!rt
!rt
call met_readh
(cdfid,number, &
& NTBEG,NTEND,NAVG,NREAD, &
& I00,IIM,J00,JJM,L00,LLM)
!rt
!rt
!rt
I00 = 0
J00 = 0
!---------LREWND = .FALSE. IF A FULL YEAR OF HISTORY DATA IS NOT AVAILABLE!
IF (NTEND-NTBEG .NE. 8760) LREWND = .FALSE.
WRITE(lun13,982) LHEADR,LREWND,NTBEG,NTEND,NAVG,NREAD
!xx
II0 = 0
JJ0 = 0
LL0 = 0
!xx
WRITE(lun13,983) II0,IIM,JJ0,JJM,LL0,LLM
!-------------CHECK WINDS & GRID SIZE-----------------------------------
IF (IIM.LT.IM .OR. JJM.LT.JM) GOTO 91
IF (IIM.NE.IMWIND .OR. JJM.NE.JMWIND .OR. LLM.NE.LMWIND) GOTO 91
!---------wind fields: TAUWB+NREAD = TAUW, winds = avg over NREAD hours,
!--------- pressure at end of intervval = TAUW
!---------locate wind field such that TAUW .LE. ITAU and ITAU .LT. TAUW + NREAD
!---------map the location of ITAU into the range NTBEG to NTEND:-------
ITAU1 = (ITAU /NREAD) *NREAD
!rt
!rt
!rt
if (edit_type .eq. "debug") then
Write (lun13,906) number
Write (lun13,900) ITAU,NREAD,ITAU1
endif
!rt
!rt
!rt
ITAU1 = MOD(ITAU1, 8760)
!rt
!rt
!rt
if (edit_type .eq. "debug") then
Write (lun13,901) ITAU1
endif
!rt
!rt
!rt
IF (ITAU1 .LT. NTBEG) ITAU1 = ITAU1 +((NTBEG -1) /8760 +1) *8760
!rt
!rt
!rt
if (edit_type .eq. "debug") then
Write (lun13,902) ITAU1,NTBEG
endif
!rt
!rt
!rt
IF (ITAU1.GT.NTEND .OR. ITAU1.EQ.NTBEG) GOTO 91
NSKIP = (ITAU1 -NTBEG) /NREAD - 1
!rt
!rt
!rt
if (edit_type .eq. "debug") then
Write (lun13,903) NSKIP
endif
!rt
!rt
!rt
!rt DO ISKIP = 1, NSKIP
!rt READ(3,ERR=93,END=95) TAUWB,TAUW
!rt ENDDO
!rt READ(3,ERR=93,END=95) TAUWB,TAUW,U,V,P,CWET,CDRY,BLH,EFFK,TAUX
!rt
!rt
!rt
ihyper=NSKIP+1
if (edit_type .eq. "debug") then
Write (lun13,905) ihyper
endif
!rt
call met_readv
(cdfid,ihyper, &
& IWPAR,JWPAR,LWPAR,IIPAR,JJPAR, &
& TAUWB,TAUW,TAUX,U,V,P,CWET,CDRY,BLH,EFFK)
!rt
IF (TAUW .NE. TAUX) GOTO 93
NTAUW = (TAUW +0.01)
ITAUW = MOD(NTAUW,8760)
!rt
if (edit_type .eq. "debug") then
Write (lun13,904) ITAU,ITAUW,NTAUW,TAUW
endif
!rt
!rt
!rt
GOTO 99
!---------read in next wind/press/etc record----------------------------
!rt41 READ(3,ERR=93,END=61) TAUWB,TAUW,U,V,P,CWET,CDRY,BLH,EFFK,TAUX
!rt
!rt
!rt
41 continue
ihyper=ihyper+1
if (ihyper.gt.number) goto 61
!rt
if (edit_type .eq. "debug") then
Write (lun13,905) ihyper
endif
!rt
call met_readv
(cdfid,ihyper, &
& IWPAR,JWPAR,LWPAR,IIPAR,JJPAR, &
& TAUWB,TAUW,TAUX,U,V,P,CWET,CDRY,BLH,EFFK)
!rt
IF (TAUW .NE. TAUX) GOTO 93
NTAUW = (TAUW+0.1)
ITAUW = MOD(NTAUW-1,8760) + 1
!rt
if (edit_type .eq. "debug") then
Write (lun13,904) ITAU,ITAUW,NTAUW,TAUW
endif
!rt
!rt
!rt
GOTO 99
!---------end of wind file, REWIND & recycle----------------------------
!rt61 IF (LREWND) THEN
!rt REWIND 3
!rt READ(3,ERR=93,END=95) NTAU1,NTAU2
!rt WRITE (lun13,'(A,2I10)') ' REWOUND GCM WIND FILE',NTAU1,NTAU2
!rt READ(3,ERR=93,END=95) TAUWB,TAUW,U,V,P,CWET,CDRY,BLH,EFFK,TAUX
!rt
61 if (LREWND) then
!rt
!rt
!rt
call met_readh
(cdfid,number, &
& NTBEG,NTEND,NAVG,NREAD, &
& I00,IIM,J00,JJM,L00,LLM)
!rt
WRITE (lun13,'(A,2I10)') ' REWOUND GCM WIND FILE',NTBEG,NTEND
!rt
ihyper=1
!rt
call met_readv
(cdfid,ihyper, &
& IWPAR,JWPAR,LWPAR,IIPAR,JJPAR, &
& TAUWB,TAUW,TAUX,U,V,P,CWET,CDRY,BLH,EFFK)
!rt
!rt
!rt
IF (TAUW .NE. TAUX) GOTO 93
NTAUW = (TAUW+0.1)
ITAUW = MOD(NTAUW-1,8760) + 1
ELSE
CALL EXITC
('=====NOT ALLOWED TO REIWND WIND FILE====')
ENDIF
GOTO 99
!----------------------------------------------------------------------------
!---------error on initial wind files------------------------------------
91 CALL EXITC
('****ERROR ILLOGICAL WIND TAPE/TRACER RUN')
!---------error on unformatted input------------------------------------
93 CALL EXITC
(' =====ERROR ON WIND/CTM FILES===== ')
!---------unexpected end of wind file-----------------------------------
95 CALL EXITC
(' =====UNEXPECTED EOF ON WINDS===== ')
982 FORMAT(' GCM DATA: LHEADR/LREWND/NTBEG/NTEND/NAVG/NREAD',2L2,4I10)
983 FORMAT(' GCM DATA: II0/IIM/JJ0/JJM/LL0/LLM',6I5)
!rt
!rt
!rt
900 format (1x,'ITAU=',i8,1x,'NREAD=',i8,1x,'ITAU1=',i8)
901 format (1x,'ITAU1=',i8)
902 format (1x,'ITAU1=',i8,1x,'NTBEG=',i8)
903 format (1x,'NSKIP=',i8)
904 format (1x,'WIND',18x,'ITAU=',i8,1x,'ITAUW=',i8, &
& 1x,'NTAUW=',i8,1x,'TAUW=',e13.6)
905 format (1x,'WIND',1x,'ihyper=',i8)
906 format (1x,'number=',i8)
!rt
!rt
!rt
99 RETURN
END
SUBROUTINE EXITC (MESSAG) 4
!xx IMPLICIT NONE
!xx
# include "uci.h"
!xx
CHARACTER*40 MESSAG
!xx$$
save
!xx$$
WRITE(lun13,'(A)') MESSAG
STOP
END
subroutine jshift1(inarr,outarr) 1
! shift mass/ad array in j by one index
#include "gmi_dims.h"
#include "uci.h"
real*8 outarr(i1:i2, ju1:j2, k1:k2)
real*8 inarr(IIPAR,JJPAR,LLPAR)
integer il, ij, ik
do 888 ik=k1,k2
do 887 ij=ju1,j2
do 886 il=i1,i2
outarr(il,ij,ik)=inarr(il,ij+1,ik)
886 continue
887 continue
888 continue
return
end
subroutine ad_to_mass(inarr,outarr)
! convert ad array to mass Gmimod form array
#include "gmi_dims.h"
#include "uci.h"
real*8 outarr(i1:i2, ju1:j2, k1:k2)
real*8 inarr(IIPAR,JJPAR,LLPAR)
integer il, ij, ik
do 888 ik=k1,k2
do 887 ij=1,JJPAR
do 886 il=i1,i2
outarr(il,ij,ik)=inarr(il,ij,ik)
886 continue
887 continue
888 continue
return
end
subroutine jshiftn(inarr,outarr)
! shift stt to const and convert to mixing ratio for diag0_gmi
#include "gmi_dims.h"
#include "uci.h"
real*8 outarr(i1:i2, ju1:j2, k1:k2, 1:num_species)
real*8 inarr(IIPAR,JJPAR,LLPAR,NNPAR)
integer il, ij, ik, ic
do 889 ic=1,num_species
do 888 ik=k1,k2
do 887 ij=ju1,j2
do 886 il=i1,i2
outarr(il,ij,ik,ic)=tcvv(ic)*inarr(il,ij+1,ik,ic)/ &
& ad(il,ij+1,ik)
886 continue
887 continue
888 continue
889 continue
return
end
subroutine stt_to_const(inarr,outarr)
! routine to convert stt to const form for diagnostics, no shift
! for gmi_diag_save
#include "gmi_dims.h"
#include "uci.h"
real*8 outarr(i1:i2, ju1:j2, k1:k2, 1:num_species)
real*8 inarr(IIPAR,JJPAR,LLPAR,NNPAR)
integer il, ij, ik, ic
do 889 ic=1,num_species
do 888 ik=k1,k2
do 887 ij=1,JJPAR ! FIX ME
do 886 il=i1,i2
outarr(il,ij,ik,ic)=tcvv(ic)*inarr(il,ij,ik,ic)/ &
& ad(il,ij,ik)
886 continue
887 continue
888 continue
889 continue
return
end
!-----------------------------------------------------------------------------
!
! ROUTINE
! Get_Diag_1d_Dump_Flag
!
! DESCRIPTION
! This routine determines whether or not it is time for the UCI 1D
! diagnostics to be dumped.
!
! ARGUMEMTS
! dynmd : year/month/day (YYYYMMDD) at end of current time step
! dump_diag : time to dump diagnostics?
!
!-----------------------------------------------------------------------------
!c subroutine Get_Diag_1d_Dump_Flag (dnymd, dump_diag)
!c implicit none
!c# include "gmi_control.h"
!c# include "gmi_gen_io.h"
! ----------------------
! Argument declarations.
! ----------------------
!c integer :: dnymd
!c logical :: dump_diag
! ----------------------
! Variable declarations.
! ----------------------
!c logical, save :: first_diag_uci = .true.
!c integer :: idumday
!c integer :: idumyear
!c integer :: this_month_diag_uci
!c integer, save :: save_month_diag_uci = 0
! ----------------
! Begin execution.
! ----------------
!c dump_diag = .false.
! ===============
!c if (pr_diag_1d) then
! ===============
!c if (pr_diag_uci_step_interval < 0) then
! -------------------------------------------
! UCI diagnostic output at monthly intervals.
! -------------------------------------------
!c if (first_diag_uci) then
!c dump_diag = .true.
!c call Split_Time_Flds
!c & (dnymd, idumyear, save_month_diag_uci, idumday)
!c first_diag_uci = .false.
!c else
!c call Split_Time_Flds
!c & (dnymd, idumyear, this_month_diag_uci, idumday)
!c if (this_month_diag_uci /= save_month_diag_uci) then
!c dump_diag = .true.
!c save_month_diag_uci = this_month_diag_uci
!c end if
!c end if
! ----------------------------------------------------------
! The +1 below is because the clock has not been ticked yet.
! ----------------------------------------------------------
!c else if (Mod (num_time_steps+1,
!c & pr_diag_uci_step_interval) == 0) then
! --------------------------------------------
! UCI diagnostic output at specified interval.
! --------------------------------------------
!c dump_diag = .true.
!c end if
!c end if
!c return
!c end