subroutine dynamics(nstep,istep) ! ! Call dynamics layer of the model for current time step. ! Routines that do not have to do message-passing are generally called ! from inside a latitude scan in this routine. Routines that do have ! to make mpi calls are called once from this routine, and latitude ! scan(s) are inside the called routine. ! ! Include the option of setting the lower boundary by using GSWM output ! this affects subroutine dt and duv ! ! am 6/06: change call to qjoule_ti to include ion and neutral density ! factors, include vertical component ! use fields_module use cons_module,only: check_exp use mpi_module,only: lon0,lon1,lat0,lat1,mytidi use chemrates_module,only: chemrates_tdep,disn2p use lbc,only: bndry_comp,add_uvbnd use chapman_module,only: chapman use qrj_module,only: qrj,qtotal,qop,qo2p,qn2p,qnp,qnop use heatnirco2_module,only: heatnearirco2 use input_module,only: iaurora=>aurora,solar_protons, | step,step_compqrj,ncep_reanalysis,ecmwf_ncfile use aurora_module,only: aurora use aurora_module,only: cusp2d,drzl2d,alfa2d,nflx2d,ec,ed use oplus_module,only: oplus use mgw_module,only: gwsource,mgw use comp_meta_module,only: comp_meta use hdif_module,only: hdif_bndlons,hdif3 use hox_module,only: comp_hox,minor_hox use n4s_module,only: comp_n4s,minor_n4s use noz_module,only: comp_noz,minor_noz use ch4_module,only: comp_ch4,minor_ch4 use co2_module,only: comp_co2,minor_co2 use co_module,only: comp_co,minor_co use h2o_module,only: comp_h2o,minor_h2o use h2_module,only: comp_h2,minor_h2 use ar_module,only: comp_ar,minor_ar use ox_module,only: comp_ox use qjchem_module,only: qjchem use qjion_module,only: qjion use radcool_module,only: radcool use dt_module,only: dt use duv_module,only: duv use timing_module,only: timer,timing use uv_bndry,only: uvbnd,uvbgrd use addfld_module,only: addfld use diags_module,only: mkdiag_CUSP,mkdiag_DRIZZLE,mkdiag_ALFA, | mkdiag_NFLUX ! implicit none ! ! Args: integer,intent(in) :: nstep,istep ! ! Local: integer :: lat,nlats,ier integer :: i0,i1,nk logical,parameter :: debug=.false. ! if true, write to stdout ! ! External: logical,external :: time2print ! ! If diags is .true., selected fields are saved to secondary histories ! (only fields specified by SECFLDS are actually saved, see below) logical,parameter :: diags=.false. ! if true, put fields on sec hist ! ! The following fields will be saved to secondary histories if diags=.true.: ! (subsets of this list can be specified by SECFLDS) ! !SECFLDS = 'OMEGA', 'UIVEL', 'VIVEL', 'WIVEL', 'CP', 'KT', 'KM', ! 'SCO2', 'SCO1', 'SCN2', 'QTOTAL', 'QOP', 'QO2P', 'QN2P', 'QNP', ! 'QNOP', 'OPLUS', 'XIOP2P', 'XIOP2D', 'HPLUS', 'NPLUS', 'ELD_N2P', ! 'ELD_NOP', 'ELD_O2P', 'ELD_NE', 'LXX', 'LYY', 'LXY', 'LYX', ! 'SIGMAPED', 'SIGMAHAL', 'GWU', 'GWV', 'GWT', 'DIFKK', 'DIFKT', ! 'DIFKV', 'CMP_N2D', 'CMP_O21D', 'CMP_O21S', 'CMP_O1D', 'CMP_HOX', ! 'CMP_OH', 'CMP_HO2', 'CMP_H', 'CMP_N4S', 'CMP_NOZ', 'CMP_NO', ! 'CMP_NO2', 'CMP_CH4', 'CMP_CO2', 'CMP_CO', 'CMP_H2O', 'CMP_H2', ! 'QCHEM','QIC','QPHOTO','QJI_TI', 'TE_UPD', 'TI_UPD', 'NOCOOL', ! 'QJI_TN', 'TN_UPD', 'U_UPD', 'V_UPD', 'UPD_O2', 'UPD_OX', ! 'UPD_O1', 'UPD_O3' ! real :: time0,time1,time1_qrj,time1_gw,time1_cmpminor, | time1_cmpmajor real,dimension(nlevp1,lon0:lon1) :: | qchem,qic,qphoto ! diagnostics time0=0.; time1=0.; time1_qrj=0.; time1_gw=0. ; time1_cmpminor=0. ; time1_cmpmajor=0. i0=lon0 ; i1=lon1 ; nk=nlevp1 nlats = lat1-lat0+1 ! ! Calculate omega for vertical velocity: ! (latitude scan is inside swdot) ! call swdot( | un(levd0,lond0,latd0,itp), ! un input | vc(levd0,lond0,latd0,itp), ! vc input | z (levd0,lond0,latd0,itc), ! input for WN diag | w (levd0,lond0,latd0,itc), ! omega output | 1,nlevp1,lon0,lon1,lat0,lat1) if (debug) write(6,"('dynamics after swdot: istep=',i5)") | istep if (diags) then do lat=lat0,lat1 call addfld('OMEGA' ,' ',' ',w(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Calculate ExB ion drift velocities from electric field: ! call ionvel( | z (levd0,lond0,latd0,itc), ! geopotential input | ui(levd0,lond0,latd0), ! U ion velocity output | vi(levd0,lond0,latd0), ! V ion velocity output | wi(levd0,lond0,latd0), ! W ion velocity output | 1,nlevp1,lon0,lon1,lat0,lat1) if (debug) write(6,"('dynamics after ionvel: istep=',i5)") | istep if (diags) then do lat=lat0,lat1 call addfld('UIVEL' ,' ',' ',ui(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('VIVEL' ,' ',' ',vi(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('WIVEL' ,' ',' ',wi(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Calculate specific heat and molecular viscosity. ! f4d itp fields are input, f3d fields are output. ! do lat=lat0,lat1 call cpktkm( | tn(levd0,lond0,lat,itp), ! 4d input | o2(levd0,lond0,lat,itp), ! 4d input | ox(levd0,lond0,lat,itp), ! 4d input | n2(levd0,lond0,lat), ! 3d input | he(levd0,lond0,lat,itp), ! 4d input | cp(levd0,lond0,lat), ! 3d output specific heat | kt(levd0,lond0,lat), ! 3d output | km(levd0,lond0,lat), ! 3d output | 1,nlevp1,lon0,lon1,lat) enddo ! lat=lat0,lat1 if (debug) write(6,"('dynamics after cpktkm: istep=',i5)") | istep if (diags) then do lat=lat0,lat1 call addfld('CP' ,' ',' ',cp(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('KT' ,' ',' ',kt(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('KM' ,' ',' ',km(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Calculate column densities and line integrals of o2, o, and n2. ! Prognostics z,tn,o2,o1,ox,o3,co2,barm are input, diagnostics ! vo2,vo1,vn2,vo3,vco2,vno, sco2,sco1,scn2,sco3,scco2,scno are ! output. ! do lat=lat0,lat1 call chapman( | z (levd0,lond0,lat,itc), ! updated Z from addiag | tn (levd0,lond0,lat,itp), ! 4d in | o2 (levd0,lond0,lat,itp), ! 4d in | o1 (levd0,lond0,lat,itp), ! 4d in | ox (levd0,lond0,lat,itp), ! 4d in | n2 (levd0,lond0,lat), ! 3d in | o3 (levd0,lond0,lat,itp), ! 4d in | co2 (levd0,lond0,lat,itp), ! 4d in | no (levd0,lond0,lat,itp), ! 4d in | barm(levd0,lond0,lat,itp), ! 4d in | vo2 (levd0,lond0,lat), ! 3d out | vo1 (levd0,lond0,lat), ! 3d out | vn2 (levd0,lond0,lat), ! 3d out | vo3 (levd0,lond0,lat), ! 3d out | vco2(levd0,lond0,lat), ! 3d out | vno (levd0,lond0,lat), ! 3d out | sco2(levd0,lond0,lat), ! 3d out | sco1(levd0,lond0,lat), ! 3d out | scn2(levd0,lond0,lat), ! 3d out | sco3(levd0,lond0,lat), ! 3d out | scco2(levd0,lond0,lat), ! 3d out | scno (levd0,lond0,lat), ! 3d out | 1,nlevp1,lon0,lon1,lat) enddo ! lat=lat0,lat1 if (debug) write(6,"('dynamics after chapman: istep=',i5)") | istep if (diags) then do lat=lat0,lat1 call addfld('SCO2' ,' ',' ',sco2(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('SCO1' ,' ',' ',sco1(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('SCN2' ,' ',' ',scn2(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Calculate temperature dependent reaction rates (chemrates_module): ! do lat=lat0,lat1 call chemrates_tdep( | tn (levd0,lond0,lat,itp), | te (levd0,lond0,lat,itp), | ti (levd0,lond0,lat,itp), | sco2(levd0,lond0,lat), | vo2 (levd0,lond0,lat), | 1,nlevp1,lon0,lon1,lat) enddo if (debug) write(6,"('dynamics after chemrates: istep=',i5)") | istep ! ! Latitude scan: time1_qrj = 0. do lat=lat0,lat1 ! ! Get lower boundary for major composition. call bndry_comp( | barm(levd0,lond0,lat,itp), | 1,nlevp1,lon0,lon1,lat) if (debug) | write(6,"('dynamics after bndry_comp: lat=',i3,' istep=',i5)") | lat,istep ! ! Calculate ionization and heating rates. call timer(time0,time1,'QRJ',0,0) ! start qrj timing for current lat ! ! 10/10/06 btf: new qrj from Liying: ! 7/23/07 btf: check step_compqrj: if (istep==1.or.mod(istep,step_compqrj/step)==0) then call qrj( | sco2(levd0,lond0,lat), | sco1(levd0,lond0,lat), | scn2(levd0,lond0,lat), | sco3(levd0,lond0,lat), | scno(levd0,lond0,lat), | vo2 (levd0,lond0,lat), | tn (levd0,lond0,lat,itp), | no (levd0,lond0,lat,itp), | o2 (levd0,lond0,lat,itp), | ox (levd0,lond0,lat,itp), | he (levd0,lond0,lat,itp), | o3 (levd0,lond0,lat,itp), | n4s (levd0,lond0,lat,itp), | o21d(levd0,lond0,lat,itp), | xnmbari(levd0,lond0,lat), | cp(levd0,lond0,lat), | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after qrj: lat=',i3,' istep=', | i5)") lat,istep call timer(time0,time1,'QRJ',1,0) ! end qrj timing for current lat time1_qrj = time1_qrj+time1 if (diags) then call addfld('QTOTAL' ,' ',' ',qtotal(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('QOP' ,' ',' ',qop(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('QO2P' ,' ',' ',qo2p(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('QN2P' ,' ',' ',qn2p(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('QNP' ,' ',' ',QNP(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('QNOP' ,' ',' ',qnop(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) endif endif ! step_compqrj ! ! Near infrared co2 heating. All args are input. Output heatnirco2 ! and heatkday are in heatnirco2_module (heatnirco2.F): ! (note heatnirco2_init is called from init.F for initialization) ! call heatnearirco2( | o2 (levd0,lond0,lat,itp), | o1 (levd0,lond0,lat,itp), ! note this is not ox | n2 (levd0,lond0,lat), | co2 (levd0,lond0,lat,itp), | barm (levd0,lond0,lat,itp), | vco2 (levd0,lond0,lat), | scco2(levd0,lond0,lat), | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after heatnearirco2: lat=',i3, | ' istep=',i5)") lat,istep ! ! subroutine qinight(tn,o2,o1,xnmbari,1,nk,lon0,lon1,lat) ! Night-time background ionization: ! call qinite( | tn (levd0,lond0,lat,itp), | o2 (levd0,lond0,lat,itp), | ox (levd0,lond0,lat,itp), | he (levd0,lond0,lat,itp), | no (levd0,lond0,lat,itp), | xnmbari(levd0,lond0,lat), | vo2 (levd0,lond0,lat), | vo1 (levd0,lond0,lat), | vn2 (levd0,lond0,lat), | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after qinite: lat=',i3, | ' istep=',i5)") lat,istep ! ! Ionization from xrays: ! 10/10/06 btf: xray is not called with new qrj. ! call xray( ! | tn (levd0,lond0,lat,itp), ! | o2 (levd0,lond0,lat,itp), ! | ox (levd0,lond0,lat,itp), ! | n2 (levd0,lond0,lat), ! | barm(levd0,lond0,lat,itp), ! | sco2(levd0,lond0,lat), ! | sco1(levd0,lond0,lat), ! | scn2(levd0,lond0,lat), ! | 1,nlevp1,lon0,lon1,lat) ! if (debug) write(6,"('dynamics after xray: lat=',i3, ! | ' istep=',i5)") lat,istep ! ! Add aurora particle precip and ionization: if (iaurora > 0) then call aurora( | tn (levd0,lond0,lat,itp), | o2 (levd0,lond0,lat,itp), | ox (levd0,lond0,lat,itp), | n2 (levd0,lond0,lat), | barm(levd0,lond0,lat,itp), | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after aurora: lat=',i3, | ' istep=',i5)") lat,istep endif ! ! End lat loop for bndry_comp, qrj, heatnirco2, qinite, xray, aurora enddo ! lat=lat0,lat1 ! ! Save 2d auroral diagnostics (see diags.F): ! if (iaurora > 0) then call mkdiag_CUSP('CUSP',cusp2d(lon0:lon1,lat0:lat1),ec, | lon0,lon1,lat0,lat1) call mkdiag_DRIZZLE('DRIZZLE',drzl2d(lon0:lon1,lat0:lat1),ed, | lon0,lon1,lat0,lat1) call mkdiag_ALFA('ALFA',alfa2d(lon0:lon1,lat0:lat1), | lon0,lon1,lat0,lat1) call mkdiag_NFLUX('NFLUX',nflx2d(lon0:lon1,lat0:lat1), | lon0,lon1,lat0,lat1) endif ! ! Report timing from above latitude loop: if (timing%level >= 2.and.time2print(nstep,istep)) | write(6,"('Time in QRJ = ', | f12.2,' Dynamics: step ',i5)") time1_qrj,istep ! ! Advance O+: ! Oplus makes several calls to mp_bndlons_f3d, mp_periodic_f3d, ! mp_gatherlons_f3d, and mp_scatterlons_f3. ! Full task subdomains (including ghost cells) are passed. ! There are several latitude loops inside oplus. ! call timer(time0,time1,'OPLUS',0,0) ! start oplus timing call oplus( | tn (levd0,lond0,latd0,itp), | te (levd0,lond0,latd0,itp), | ti (levd0,lond0,latd0,itp), | o2 (levd0,lond0,latd0,itp), | o1 (levd0,lond0,latd0,itp), | ox (levd0,lond0,latd0,itp), | n2 (levd0,lond0,latd0), | he (levd0,lond0,latd0,itp), | n2d (levd0,lond0,latd0,itp), | ne (levd0,lond0,latd0,itp), | un (levd0,lond0,latd0,itp), | vn (levd0,lond0,latd0,itp), | w (levd0,lond0,latd0,itp), | h (levd0,lond0,latd0,itp), | h2 (levd0,lond0,latd0,itp), | h2o (levd0,lond0,latd0,itp), | co2 (levd0,lond0,latd0,itp), | barm (levd0,lond0,latd0,itp), | xnmbar(levd0,lond0,latd0), ! output | ui (levd0,lond0,latd0), ! U ion velocity from ionvel | vi (levd0,lond0,latd0), ! V ion velocity from ionvel | wi (levd0,lond0,latd0), ! W ion velocity from ionvel | op (levd0,lond0,latd0,itp), ! O+ from previous step | op_nm (levd0,lond0,latd0,itp), ! O+ at time n-1 | op (levd0,lond0,latd0,itc), ! output new O+ | op_nm (levd0,lond0,latd0,itc), ! output new O+ time n-1 | xiop2p, ! output O+(2p) | xiop2d, ! output O+(2d) | hplus, ! output H+ | nplus, ! output N+ | 1,nlevp1,lon0,lon1,lat0,lat1) ! if (debug) write(6,"('dynamics after oplus')") call timer(time0,time1,'OPLUS',1,0) ! end oplus timing if (timing%level >= 2.and.time2print(nstep,istep)) | write(6,"('Time in OPLUS = ', | f12.2,' Dynamics: step ',i5)") time1,istep if (diags) then do lat=lat0,lat1 call addfld('OPLUS' ,' ',' ',op(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('XIOP2P',' ',' ',xiop2p(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('XIOP2D',' ',' ',xiop2d(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('HPLUS' ,' ',' ',hplus(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('NPLUS' ,' ',' ',nplus(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Electron density: do lat=lat0,lat1 call elden( | tn (levd0,lond0,lat,itp), | barm (levd0,lond0,lat,itp), | op (levd0,lond0,lat,itp), | o2 (levd0,lond0,lat,itp), | o1 (levd0,lond0,lat,itp), | ox (levd0,lond0,lat,itp), | n2 (levd0,lond0,lat), | n2d (levd0,lond0,lat,itp), | no (levd0,lond0,lat,itp), | n4s (levd0,lond0,lat,itp), | co2 (levd0,lond0,lat,itp), | xiop2p(levd0,lond0,lat), ! O+(2P) from oplus | xiop2d(levd0,lond0,lat), ! O+(2D) from oplus | z (levd0,lond0,lat,itc), ! updated Z from addiag | nplus (levd0,lond0,lat), ! N+ from oplus | n2p (levd0,lond0,lat), ! N2+ output | nop (levd0,lond0,lat), ! NO+ output | o2p (levd0,lond0,lat,itc), ! O2+ output | ne (levd0,lond0,lat,itc), ! electron density output | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after elden: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('ELD_N2P',' ',' ',n2p(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('ELD_NOP',' ',' ',nop(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('ELD_O2P',' ',' ',o2p(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('ELD_NE' ,' ',' ',ne (:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) endif enddo ! lat=lat0,lat1 ! ! Calculate ion drag coefficients and conductivities. ! time1_gw = 0. do lat=lat0,lat1 call lamdas( | tn (levd0,lond0,lat,itp), | barm(levd0,lond0,lat,itp), | o2 (levd0,lond0,lat,itp), | ox (levd0,lond0,lat,itp), | n2 (levd0,lond0,lat), | ti (levd0,lond0,lat,itp), | te (levd0,lond0,lat,itp), | o2p (levd0,lond0,lat,itp), | op (levd0,lond0,lat,itp), | nop (levd0,lond0,lat), | lxx (levd0,lond0,lat), ! LXX output | lyy (levd0,lond0,lat), ! LYY output | lxy (levd0,lond0,lat), ! LXY output | lyx (levd0,lond0,lat), ! LYX output | lam1(levd0,lond0,lat), ! LAM1 output | ped (levd0,lond0,lat), ! pedersen conductivity output | hall(levd0,lond0,lat), ! hall conductivity output | 1,nlevp1,lon0,lon1,lat) ! if (debug) write(6,"('dynamics after lamdas: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('LXX',' ',' ',lxx(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('LYY',' ',' ',lyy(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('LXY',' ',' ',lxy(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('LYX',' ',' ',lyx(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('SIGMAPED',' ',' ',ped(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('SIGMAHAL',' ',' ',hall(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) endif ! ! Add qn2p to pdn2 (formerly qtef): call qtieff(1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after qtieff: lat=',i3, | ' istep=',i5)") lat,istep ! ! Thermal and momentum inputs due to gravity waves: call timer(time0,time1,'GW',0,0) ! start gravity wave timing call gwsource( | z(levd0,lond0,lat,itp), ! geopotential input | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after gwsource: lat=',i3, | ' istep=',i5)") lat,istep ! ! Gravity wave parameterization: mgw calls waccm/cam1 gw module (gw_drag.F). call mgw( | tn (levd0,lond0,lat,itp), | un (levd0,lond0,lat,itp), | vn (levd0,lond0,lat,itp), | barm(levd0,lond0,lat,itp), | z (levd0,lond0,lat,itp), | cp (levd0,lond0,lat), ! specific heat | gwu (levd0,lond0,lat), ! U tendency output | gwv (levd0,lond0,lat), ! V tendency output | gwt (levd0,lond0,lat), ! T tendency output | difkk(levd0,lond0,lat), ! eddy viscosity output | difkt(levd0,lond0,lat), ! eddy thermal diffusion output | difkv(levd0,lond0,lat), ! eddy diffusion output | 1,nlevp1,lon0,lon1,lat) call timer(time0,time1,'GW',1,0) ! end gravity wave timing time1_gw = time1_gw+time1 if (debug) write(6,"('dynamics after mgw: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('GWU',' ',' ',gwu(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('GWV',' ',' ',gwv(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('GWT',' ',' ',gwt(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('DIFKK',' ',' ',difkk(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('DIFKT',' ',' ',difkt(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('DIFKV',' ',' ',difkv(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) endif ! ! Define 3-d n2o, cl, clo from solomon-garcia imports (see solgar.F) call comp_solgar( | n2o(levd0,lond0,lat), | cl (levd0,lond0,lat), | clo(levd0,lond0,lat), | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_solgar: lat=',i3, | ' istep=',i5)") lat,istep ! ! Define meta-species n(2d), o(1d), o2(1s), o2(1d), and h2o2: ! call comp_meta( | barm(levd0,lond0,lat,itp), | xnmbar(levd0,lond0,lat), ! from oplus | o2 (levd0,lond0,lat,itp), | o1 (levd0,lond0,lat,itp), | ox (levd0,lond0,lat,itp), | o3 (levd0,lond0,lat,itp), | h2o (levd0,lond0,lat,itp), | h2 (levd0,lond0,lat,itp), | ch4 (levd0,lond0,lat,itp), | co (levd0,lond0,lat,itp), | co2 (levd0,lond0,lat,itp), | no (levd0,lond0,lat,itp), | n4s (levd0,lond0,lat,itp), | ne (levd0,lond0,lat,itp), | n2o (levd0,lond0,lat), ! n2o from comp_solgar | op (levd0,lond0,lat,itp), ! O+ from oplus | o2p (levd0,lond0,lat,itp), ! O2+ from elden | n2p (levd0,lond0,lat), ! N2+ from elden | nop (levd0,lond0,lat), ! NO+ from elden | n2d (levd0,lond0,lat,itc), ! n2d output | o1d (levd0,lond0,lat), ! o(1delta) output | o21s(levd0,lond0,lat), ! o2(1sigma) output | o21d(levd0,lond0,lat,itc), ! o2(1delta) output | h2o2(levd0,lond0,lat), ! h2o2 output | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_meta: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('CMP_N2D',' ',' ',n2d(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('CMP_O21D',' ',' ',o21d(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('CMP_O21S',' ',' ',o21s(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('CMP_O1D',' ',' ',o1d(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) endif ! diags ! ! HOX production from D-region ion chemistry (only if solar_protons > 0): ! (default solar_protons == 0) ! call hoxpion( | tn (levd0,lond0,lat,itp), | o2 (levd0,lond0,lat,itp), | ox (levd0,lond0,lat,itp), | o1 (levd0,lond0,lat,itp), | n2 (levd0,lond0,lat), | no (levd0,lond0,lat,itp), | h2o(levd0,lond0,lat,itp), | co2(levd0,lond0,lat,itp), | ne (levd0,lond0,lat,itp), | barm (levd0,lond0,lat,itp), | xnmbar(levd0,lond0,lat), | phoxic(levd0,lond0,lat), ! phoxic output | solar_protons,1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after hoxpion: lat=',i3, | ' istep=',i5)") lat,istep if (diags) | call addfld('PHOXIC',' ',' ',phoxic(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 ! ! Report timings from above latitude loop: if (timing%level >= 2.and.time2print(nstep,istep)) | write(6,"('Time in GW = ', | f12.2,' Dynamics: step ',i5)") time1_gw,istep ! ! Make boundary longitudes available for inputs to hdif3: ! call hdif_bndlons(kldt,kldu,kldv,kldo2,kldox,kldn4s,kldnoz, | kldco,kldco2,kldh2o,kldh2,kldhox,kldch4,kldar,kldhe,kldnat, | 1,nlevp1,lon0,lon1,lat0,lat1) ! ! Calculate horizontal diffusion terms using coefficients from ! hdif2 (hdif2 was called from advance). ! do lat=lat0,lat1 call hdif3( | cp(levd0,lond0,lat), ! specific heat input | kldt ,kldu ,kldv ,kldo2 ,kldox ,kldn4s ,kldnoz ,kldco , | kldco2 ,kldh2o ,kldh2 ,kldhox ,kldch4 ,kldar ,kldhe ,kldnat, | hdt (levd0,lond0,lat), ! 2d tn output at current lat | hdu (levd0,lond0,lat), ! 2d un output at current lat | hdv (levd0,lond0,lat), ! 2d vn output at current lat | hdo2 (levd0,lond0,lat), ! 2d o2 output at current lat | hdox (levd0,lond0,lat), ! 2d ox output at current lat | hdn4s(levd0,lond0,lat), ! 2d n4s output at current lat | hdnoz(levd0,lond0,lat), ! 2d noz output at current lat | hdco (levd0,lond0,lat), ! 2d co output at current lat | hdco2(levd0,lond0,lat), ! 2d co2 output at current lat | hdh2o(levd0,lond0,lat), ! 2d h2o output at current lat | hdh2 (levd0,lond0,lat), ! 2d h2 output at current lat | hdhox(levd0,lond0,lat), ! 2d hox output at current lat | hdch4(levd0,lond0,lat), ! 2d ch4 output at current lat | hdar (levd0,lond0,lat), ! 2d ar output at current lat | hdhe (levd0,lond0,lat), ! 2d he output at current lat | hdnat(levd0,lond0,lat), ! 2d nat output at current lat | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after hdif3: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! ! Advance hox: ! CMPMINOR timing will include hox, n4s, noz, ch4, co2, co, h2o, h2 ! call timer(time0,time1,'CMPMINOR',0,0) ! start comp minor timing do lat=lat0,lat1 call comp_hox( | tn (levd0,lond0,lat,itp), | o2 (levd0,lond0,lat,itp), | ox (levd0,lond0,lat,itp), | op (levd0,lond0,lat,itp), | h2o (levd0,lond0,lat,itp), | oh (levd0,lond0,lat,itp), | ho2 (levd0,lond0,lat,itp), | h2 (levd0,lond0,lat,itp), | h (levd0,lond0,lat,itp), | hox (levd0,lond0,lat,itp), | o1 (levd0,lond0,lat,itp), | o3 (levd0,lond0,lat,itp), | no (levd0,lond0,lat,itp), | co (levd0,lond0,lat,itp), | n4s (levd0,lond0,lat,itp), | ch4 (levd0,lond0,lat,itp), | z (levd0,lond0,lat,itp), | phoxic(levd0,lond0,lat), ! from hoxpion | o1d (levd0,lond0,lat), ! o(1delta) from comp_meta | h2o2 (levd0,lond0,lat), ! h2o2 from comp_meta | hplus (levd0,lond0,lat), ! H+ from oplus | xnmbar(levd0,lond0,lat), | barm(levd0,lond0,lat,itp), | hox_nm(:,lon0-2:lon1+2,lat-2:lat+2,itp), ! hox at prev time n-1 (5 lats) | hox_nm (levd0,lond0,lat,itc), ! hox at time n-1 output | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_hox: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! ! Complete hox and partition into oh, ho2, and h. ! Fields hox, hox_nm, oh, ho2, and h are output by minor_hox. ! Minor_hox calls minor2, which includes fft filtering. ! Latitude loops are internal to minor_hox and minor2. ! call minor_hox( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | hox (levd0,lond0,latd0,itp), ! hox from previous step | hox_nm (levd0,lond0,latd0,itp), ! hox at time n-1 | hdhox (levd0,lond0,latd0), ! hox diffusion from hdif3 | hox (levd0,lond0,latd0,itc), ! hox output | hox_nm (levd0,lond0,latd0,itc), ! hox output at time n-1 | oh (levd0,lond0,latd0,itc), ! oh output | ho2 (levd0,lond0,latd0,itc), ! ho2 output | h (levd0,lond0,latd0,itc), ! h output | 1,nlevp1,lon0,lon1,lat0,lat1) ! if (debug) write(6,"('dynamics after minor_hox', | ' istep=',i5)") istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 if (diags) then do lat=lat0,lat1 call addfld('CMP_HOX',' ',' ',hox(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('CMP_OH',' ',' ' ,oh(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('CMP_HO2',' ',' ',ho2(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('CMP_H',' ',' ' ,h(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 endif ! diags ! ! Advance n4s: ! call timer(time0,time1,'CMPMINOR',0,0) ! resume comp minor timing do lat=lat0,lat1 call comp_n4s( | tn (levd0,lond0,lat,itp), ! tn (deg K) | o2 (levd0,lond0,lat,itp), ! o2 (mmr) | ox (levd0,lond0,lat,itp), ! ox (mmr) | o1 (levd0,lond0,lat,itp), ! o1 (mmr) | n2 (levd0,lond0,lat), ! n2 (mmr) | no (levd0,lond0,lat,itp), ! no (mmr) | no2 (levd0,lond0,lat,itp), ! no2 (mmr) | oh (levd0,lond0,lat,itp), ! oh (mmr) | hox (levd0,lond0,lat,itp), ! hox (mmr) | ne (levd0,lond0,lat,itp), ! ne (cm3) | n2d (levd0,lond0,lat,itc), ! updated n2d from comp_meta | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! n*mbar | o2p (levd0,lond0,lat,itp), ! o2+ | op (levd0,lond0,lat,itp), ! o+ | n2p (levd0,lond0,lat), ! n2+ | nplus (levd0,lond0,lat), ! n+ | nop (levd0,lond0,lat), ! no+ | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_n4s: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! call minor_n4s( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | n4s (levd0,lond0,latd0,itp), ! n4s from previous step | n4s_nm (levd0,lond0,latd0,itp), ! n4s at time n-1 | hdn4s (levd0,lond0,latd0) , ! n4s diffusion from hdif3 | n4s (levd0,lond0,latd0,itc), ! n4s output | n4s_nm (levd0,lond0,latd0,itc), ! n4s output at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) ! if (debug) write(6,"('dynamics after minor_n4s: istep=',i5)") | istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 if (diags) then do lat=lat0,lat1 call addfld('CMP_N4S',' ',' ',n4s(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 endif ! diags ! ! Advance noz: ! call timer(time0,time1,'CMPMINOR',0,0) ! resume comp minor timing do lat=lat0,lat1 call comp_noz( | tn (levd0,lond0,lat,itp), ! tn (deg K) | o2 (levd0,lond0,lat,itp), ! o2 (mmr) | o1 (levd0,lond0,lat,itp), ! o1 (mmr) | o3 (levd0,lond0,lat,itp), ! o3 (mmr) | no (levd0,lond0,lat,itp), ! no (mmr) | no2 (levd0,lond0,lat,itp), ! no2 (mmr) | n2o (levd0,lond0,lat) , ! n2o from comp_solgar | n4s (levd0,lond0,lat,itp), ! n4s (mmr) | oh (levd0,lond0,lat,itp), ! oh (mmr) | hox (levd0,lond0,lat,itp), ! hox (mmr) | ne (levd0,lond0,lat,itp), ! ne (cm3) | n2d (levd0,lond0,lat,itc), ! updated n2d from comp_meta | o1d (levd0,lond0,lat), ! o(1delta) from comp_meta | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! n*mbar | o2p (levd0,lond0,lat,itp), ! o2+ | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_noz: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! call minor_noz( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | noz (levd0,lond0,latd0,itp), ! noz from previous step | noz_nm (levd0,lond0,latd0,itp), ! noz at time n-1 | hdnoz (levd0,lond0,latd0) , ! noz diffusion from hdif3 | noz (levd0,lond0,latd0,itc), ! noz output | noz_nm (levd0,lond0,latd0,itc), ! noz output at time n-1 | no (levd0,lond0,latd0,itc), ! no output | no2 (levd0,lond0,latd0,itc), ! no2 output | 1,nlevp1,lon0,lon1,lat0,lat1) ! if (debug) write(6,"('dynamics after minor_noz', | ' istep=',i5)") istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 if (diags) then do lat=lat0,lat1 call addfld('CMP_NOZ',' ',' ',noz(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('CMP_NO' ,' ',' ',no(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('CMP_NO2',' ',' ',no2(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 endif ! diags ! ! Advance ch4: ! call timer(time0,time1,'CMPMINOR',0,0) ! resume comp minor timing if (istep==1.or.mod(istep,step_compqrj/step)==0) then do lat=lat0,lat1 call comp_ch4( | o1 (levd0,lond0,lat,itp), ! o1 | ch4 (levd0,lond0,lat,itp), ! ch4 | hox (levd0,lond0,lat,itp), ! hox | o1d (levd0,lond0,lat), ! o(1delta) from comp_meta | cl (levd0,lond0,lat), ! chlorine from comp_solgar | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! xnmbar | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_ch4: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! call minor_ch4( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | ch4 (levd0,lond0,latd0,itp), ! ch4 from previous step | ch4_nm (levd0,lond0,latd0,itp), ! ch4 at time n-1 | hdch4 (levd0,lond0,latd0) , ! ch4 diffusion from hdif3 | ch4 (levd0,lond0,latd0,itc), ! ch4 output | ch4_nm (levd0,lond0,latd0,itc), ! ch4 output at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) ! if (debug) write(6,"('dynamics after minor_ch4: istep=',i5)") | istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 endif if (diags) then do lat=lat0,lat1 call addfld('CMP_CH4',' ',' ',ch4(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 endif ! diags ! ! Advance co2: ! call timer(time0,time1,'CMPMINOR',0,0) ! resume comp minor timing do lat=lat0,lat1 call comp_co2( | o1 (levd0,lond0,lat,itp), ! o1 | co (levd0,lond0,lat,itp), ! co | hox (levd0,lond0,lat,itp), ! hox | o1d (levd0,lond0,lat), ! o1d | op (levd0,lond0,lat,itp), ! o+ | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! xnmbar | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_co2: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! call minor_co2( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | co2 (levd0,lond0,latd0,itp), ! co2 from previous step | co2_nm (levd0,lond0,latd0,itp), ! co2 at time n-1 | hdco2 (levd0,lond0,latd0) , ! co2 diffusion from hdif3 | co2 (levd0,lond0,latd0,itc), ! co2 output | co2_nm (levd0,lond0,latd0,itc), ! co2 output at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) if (debug) write(6,"('dynamics after minor_co2: istep=',i5)") | istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 if (diags) then do lat=lat0,lat1 call addfld('CMP_CO2',' ',' ',co2(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 endif ! diags ! ! Advance co: ! call timer(time0,time1,'CMPMINOR',0,0) ! resume comp minor timing do lat=lat0,lat1 call comp_co( | o1 (levd0,lond0,lat,itp), ! o1 | co (levd0,lond0,lat,itp), ! co | co2 (levd0,lond0,lat,itp), ! co2 | hox (levd0,lond0,lat,itp), ! hox | o1d (levd0,lond0,lat), ! o1d | op (levd0,lond0,lat,itp), ! op | ch4 (levd0,lond0,lat,itp), ! ch4 | cl (levd0,lond0,lat), ! from solgar | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! xnmbar | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_co: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! call minor_co( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | co (levd0,lond0,latd0,itp), ! co from previous step | co_nm (levd0,lond0,latd0,itp), ! co at time n-1 | hdco (levd0,lond0,latd0) , ! co diffusion from hdif3 | co (levd0,lond0,latd0,itc), ! co output | co_nm (levd0,lond0,latd0,itc), ! co output at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) ! if (debug) write(6,"('dynamics after minor_co: istep=',i5)") | istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 if (diags) then do lat=lat0,lat1 call addfld('CMP_CO',' ',' ',co(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 endif ! diags ! ! Advance h2o: ! call timer(time0,time1,'CMPMINOR',0,0) ! resume comp minor timing do lat=lat0,lat1 call comp_h2o( | o1 (levd0,lond0,lat,itp), ! o1 | h2 (levd0,lond0,lat,itp), ! h2 | h2o (levd0,lond0,lat,itp), ! h2o | h2o2 (levd0,lond0,lat), ! h2o2 | hox (levd0,lond0,lat,itp), ! hox | o1d (levd0,lond0,lat), ! o1d | op (levd0,lond0,lat,itp), ! o+ | ch4 (levd0,lond0,lat,itp), ! methane | cl (levd0,lond0,lat), ! from solgar | phoxic(levd0,lond0,lat), ! from hoxpion | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! xnmbar | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_h2o: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! call minor_h2o( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | h2o (levd0,lond0,latd0,itp), ! h2o from previous step | h2o_nm (levd0,lond0,latd0,itp), ! h2o at time n-1 | hdh2o (levd0,lond0,latd0) , ! h2o diffusion from hdif3 | h2o (levd0,lond0,latd0,itc), ! h2o output | h2o_nm (levd0,lond0,latd0,itc), ! h2o output at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) if (debug) write(6,"('dynamics after minor_h2o: istep=',i5)") | istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 if (diags) then do lat=lat0,lat1 call addfld('CMP_H2O',' ',' ',h2o(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 endif ! diags ! ! Advance h2: ! call timer(time0,time1,'CMPMINOR',0,0) ! resume comp minor timing do lat=lat0,lat1 call comp_h2( | o1 (levd0,lond0,lat,itp), ! o1 | h2 (levd0,lond0,lat,itp), ! h2 | h2o (levd0,lond0,lat,itp), ! h2o | hox (levd0,lond0,lat,itp), ! hox | o1d (levd0,lond0,lat), ! o1d | op (levd0,lond0,lat,itp), ! op | ch4 (levd0,lond0,lat,itp), ! methane | barm (levd0,lond0,lat,itp), ! barm | xnmbar(levd0,lond0,lat), ! xnmbar | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_h2: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! call minor_h2( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | h2 (levd0,lond0,latd0,itp), ! h2 from previous step | h2_nm (levd0,lond0,latd0,itp), ! h2 at time n-1 | hdh2 (levd0,lond0,latd0) , ! h2 diffusion from hdif3 | h2 (levd0,lond0,latd0,itc), ! h2 output | h2_nm (levd0,lond0,latd0,itc), ! h2 output at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) ! if (debug) write(6,"('dynamics after minor_h2: istep=',i5)") istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 if (timing%level >= 2.and.time2print(nstep,istep)) | write(6,"('Time in CMPMINOR = ',f12.2, | ' Dynamics: step ',i5)") time1_cmpminor,istep if (diags) then do lat=lat0,lat1 call addfld('CMP_H2',' ',' ',h2(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Advance argon: ! call timer(time0,time1,'CMPMINOR',0,0) ! resume comp minor timing do lat=lat0,lat1 call comp_ar( | ar (levd0,lond0,lat,itp), ! argon from previous step | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_ar: lat=',i3, | ' istep=',i5)") lat,istep enddo ! lat=lat0,lat1 ! call minor_ar( | tn (levd0,lond0,latd0,itp), ! tn (deg K) | o2 (levd0,lond0,latd0,itp), ! o2 (mmr) | ox (levd0,lond0,latd0,itp), ! ox (mmr) | n2 (levd0,lond0,latd0), ! n2 (mmr) | he (levd0,lond0,latd0,itp), ! he (mmr) | w (levd0,lond0,latd0,itc), ! w (from swdot) | difkk (levd0,lond0,latd0) , ! eddy viscosity from mgw | ar (levd0,lond0,latd0,itp), ! ar from previous step | ar_nm (levd0,lond0,latd0,itp), ! ar at time n-1 | hdar (levd0,lond0,latd0) , ! ar diffusion from hdif3 | ar (levd0,lond0,latd0,itc), ! ar output | ar_nm (levd0,lond0,latd0,itc), ! ar output at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) ! if (debug) write(6,"('dynamics after minor_ar: istep=',i5)") istep call timer(time0,time1,'CMPMINOR',1,0) ! suspend comp minor timing time1_cmpminor = time1_cmpminor+time1 if (timing%level >= 2.and.time2print(nstep,istep)) | write(6,"('Time in CMPMINOR = ',f12.2, | ' Dynamics: step ',i5)") time1_cmpminor,istep if (diags) then do lat=lat0,lat1 ! write(6,"('dynamics after minor_ar: lat=',i4,' ar min,max=', ! | 2e12.4)") lat,minval(ar(:,i0:i1,lat,itc)), ! | maxval(ar(:,i0:i1,lat,itc)) call addfld('CMP_AR',' ',' ',ar(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Calculate matrix coefficients for partitioning of OX: ! do lat=lat0,lat1 call comp_ox( | tn (levd0,lond0,lat,itp), ! tn (deg K) | o2 (levd0,lond0,lat,itp), ! o2 (mmr) | ox (levd0,lond0,lat,itp), ! ox (mmr) | n2 (levd0,lond0,lat), ! n2 (mmr) | op (levd0,lond0,lat,itc), ! o+ (mmr) | h2o (levd0,lond0,lat,itc), ! h2o (mmr) | h2o2 (levd0,lond0,lat), ! h2o2(mmr) | oh (levd0,lond0,lat,itc), ! oh (mmr) | ho2 (levd0,lond0,lat,itc), ! ho2 (mmr) | h2 (levd0,lond0,lat,itc), ! h2 (mmr) | h (levd0,lond0,lat,itc), ! h (mmr) | o1 (levd0,lond0,lat,itp), ! o (mmr) | o1d (levd0,lond0,lat), ! o(1d) (mmr) | o3 (levd0,lond0,lat,itp), ! o3 (mmr) | no (levd0,lond0,lat,itc), ! no (mmr) | no2 (levd0,lond0,lat,itc), ! no2 (mmr) | n4s (levd0,lond0,lat,itc), ! n4s (mmr) | n2d (levd0,lond0,lat,itc), ! n(2d) (mmr) | o2p (levd0,lond0,lat,itp), ! o2+ | n2p (levd0,lond0,lat), ! n2+ | nplus (levd0,lond0,lat), ! n+ | hplus (levd0,lond0,lat), ! h+ | nop (levd0,lond0,lat), ! no+ | ne (levd0,lond0,lat,itc), ! ne (cm3) | co (levd0,lond0,lat,itc), ! co (mmr) | ch4 (levd0,lond0,lat,itc), ! ch4 (mmr) | cl (levd0,lond0,lat), ! cl (solgar) | clo (levd0,lond0,lat), ! clo (solgar) | xiop2p(levd0,lond0,lat), ! o+(2p) | xiop2d(levd0,lond0,lat), ! o+(2d) | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! p0e(-z)... | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after comp_ox: lat=',i3, | ' istep=',i5)") lat,istep ! ! Heating and o2 dissociation due to odd nitrogen and ion chemistry: ! call qjchem( | tn (levd0,lond0,lat,itp), ! tn (deg K) | o2 (levd0,lond0,lat,itp), ! o2 (mmr) | ox (levd0,lond0,lat,itp), ! ox (mmr) | he (levd0,lond0,lat,itp), ! he (mmr) | no (levd0,lond0,lat,itc), ! no (mmr) | no2 (levd0,lond0,lat,itc), ! no2 (mmr) | n4s (levd0,lond0,lat,itc), ! n4s (mmr) | n2d (levd0,lond0,lat,itc), ! n(2d) (mmr) | ne (levd0,lond0,lat,itp), ! ne (cm3) | o3 (levd0,lond0,lat,itp), ! o3 (mmr) | o1 (levd0,lond0,lat,itp), ! o (mmr) | o1d (levd0,lond0,lat), ! o(1d) (mmr) | co (levd0,lond0,lat,itc), ! co (mmr) | co2 (levd0,lond0,lat,itc), ! co2 (mmr) | h2o (levd0,lond0,lat,itc), ! h2o (mmr) | h2 (levd0,lond0,lat,itc), ! h2 (mmr) | h (levd0,lond0,lat,itc), ! h (mmr) | oh (levd0,lond0,lat,itc), ! oh (mmr) | ho2 (levd0,lond0,lat,itc), ! ho2 (mmr) | h2o2 (levd0,lond0,lat), ! h2o2 | ch4 (levd0,lond0,lat,itc), ! ch4 (mmr) | cl (levd0,lond0,lat), ! cl (solgar) | clo (levd0,lond0,lat), ! clo (solgar) | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! p0e(-z)*barm/(boltz*T) | cp (levd0,lond0,lat), ! specific heat | qchem (levd0,lon0), ! chem heating output (local for diag) | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after qjchem: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('QCHEM',' ','DEG/DAY',qchem, | 'lev',1,nk,'lon',i0,i1,lat) endif ! ! Ion heating: ! call qjion( | tn (levd0,lond0,lat,itp), ! tn (deg K) | o2 (levd0,lond0,lat,itp), ! o2 (mmr) | ox (levd0,lond0,lat,itp), ! ox (mmr) | n2 (levd0,lond0,lat), ! n2 (mmr) | o2p (levd0,lond0,lat,itp), ! o2+ | op (levd0,lond0,lat,itp), ! o+ | n2p (levd0,lond0,lat), ! n2+ | nplus (levd0,lond0,lat), ! n+ | nop (levd0,lond0,lat), ! no+ | xiop2p(levd0,lond0,lat), ! o+(2p) from oplus | xiop2d(levd0,lond0,lat), ! o+(2d) from oplus | n4s (levd0,lond0,lat,itp), ! n(4s) | n2d (levd0,lond0,lat,itp), ! n(2d) | no (levd0,lond0,lat,itp), ! no | ne (levd0,lond0,lat,itp), ! ne | o1 (levd0,lond0,lat,itp), ! o1 | barm (levd0,lond0,lat,itp), ! mbar | cp (levd0,lond0,lat), ! specific heat | qic (levd0,lon0), ! ion chem heating output | qphoto(levd0,lon0), ! photo electron heating output | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after qjion: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('QIC' ,' ',' ',qic , | 'lev',1,nk,'lon',i0,i1,lat) call addfld('QPHOTO',' ',' ',qphoto, | 'lev',1,nk,'lon',i0,i1,lat) endif ! ! Joule heating for ion temperature: ! call qjoule_ti( | un (levd0,lond0,lat,itp), ! zonal neutral velocity input | vn (levd0,lond0,lat,itp), ! meridional velocity input | w (levd0,lond0,lat,itp), ! vertical velocity input | ui (levd0,lond0,lat), ! zonal ion velocity input (ExB) | vi (levd0,lond0,lat), ! meridional ion velocity input (ExB) | wi (levd0,lond0,lat), ! vertical ion velocity input (ExB) | lam1(levd0,lond0,lat), ! lambda_1 ion drag coefficient [1/s] | op (levd0,lond0,lat,itp), ! O+ 1/cm^3 | nplus(levd0,lond0,lat), ! N+ 1/cm^3 | n2p(levd0,lond0,lat), ! N2+ 1/cm^3 | nop(levd0,lond0,lat), ! NO+ 1/cm^3 | o2p(levd0,lond0,lat,itp), ! O2+ 1/cm^3 | barm(levd0,lond0,lat,itp),! mean molecular weight | tn (levd0,lond0,lat,itp),! neutral temperature | qji_ti(levd0,lond0,lat), ! ion joule heating output | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after qjoule_ti: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('QJI_TI',' ',' ',qji_ti(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) endif ! ! Electron and Ion temperatures: ! call settei( | tn (levd0,lond0,lat,itp), ! tn (deg K) | o2 (levd0,lond0,lat,itp), ! o2 (mmr) | ox (levd0,lond0,lat,itp), ! ox (mmr) | n2 (levd0,lond0,lat), ! n2 (mmr) | ne (levd0,lond0,lat,itp), ! ne | op (levd0,lond0,lat,itp), ! o+ | o2p (levd0,lond0,lat,itp), ! o2+ | nop (levd0,lond0,lat), ! no+ | barm (levd0,lond0,lat,itp), ! mbar | cp (levd0,lond0,lat), ! specific heat | qji_ti(levd0,lond0,lat), ! ion joule heating from qjoule_ti | te (levd0,lond0,lat,itp), ! te from previous step | ti (levd0,lond0,lat,itp), ! ti from previous step | te (levd0,lond0,lat,itc), ! updated electron temp (output) | ti (levd0,lond0,lat,itc), ! updated ion temp (output) | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after settei: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('TE_UPD',' ',' ',te(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('TI_UPD',' ',' ',ti(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) endif ! ! Radiative cooling (co2, o3, no, o(3p)): ! call radcool( | tn (levd0,lond0,lat,itp), ! tn (deg K) | o2 (levd0,lond0,lat,itp), ! o2 (mmr) | ox (levd0,lond0,lat,itp), ! ox (mmr) | o1 (levd0,lond0,lat,itp), ! o (mmr) | n2 (levd0,lond0,lat), ! n2 (mmr) | o3 (levd0,lond0,lat,itp), ! o3 (mmr) | no (levd0,lond0,lat,itp), ! no (mmr) | co2 (levd0,lond0,lat,itp), ! co2(mmr) | barm (levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! p0e(-z)*barm/(boltz*T) | cp (levd0,lond0,lat), ! specific heat | cool_implicit(levd0,lond0,lat), ! implicit cooling output | cool_explicit(levd0,lond0,lat), ! explicit cooling output | nocool(levd0,lond0,lat), ! nocool output | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after radcool: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('NOCOOL',' ',' ',nocool(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) endif ! ! Joule heating for tn: ! call qjoule_tn( | un (levd0,lond0,lat,itp), ! zonal neutral velocity input | vn (levd0,lond0,lat,itp), ! meridional velocity input | w (levd0,lond0,lat,itp), ! vertical velocity input | ui (levd0,lond0,lat), ! zonal ion velocity input | vi (levd0,lond0,lat), ! meridional ion velocity input | wi (levd0,lond0,lat), ! vertical ion velocity input | lam1(levd0,lond0,lat), ! lambda_1 ion drag coefficient | barm(levd0,lond0,lat,itp),! mean molecular weight | tn (levd0,lond0,lat,itp),! neutral temperature | qji_tn(levd0,lond0,lat), ! ion joule heating output | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after qjoule_tn: lat=',i3, | ' istep=',i5)") lat,istep if (diags) then call addfld('QJI_TN',' ',' ',qji_tn(:,i0:i1,lat), | 'lev',1,nk,'lon',i0,i1,lat) endif enddo ! lat=lat0,lat1 ! ! Calculate ubnd,vbnd if not using ecmwf or ncep_rean data files: if (len_trim(ecmwf_ncfile) == 0 .and. | len_trim(ncep_reanalysis) == 0) then ! ! calculation of Un,Vn at LBC (background and perturbations) if(istep < 3) then call uvbnd( | un (levd0,lond0,latd0,itp), | un_nm(levd0,lond0,latd0,itp), | vn (levd0,lond0,latd0,itp), | vn_nm(levd0,lond0,latd0,itp), | w (levd0,lond0,latd0,itp), | z (levd0,lond0,latd0,itp), | barm (levd0,lond0,latd0,itp), | hdu (levd0,lond0,latd0), ! horiz diffusion of U from hdif3 | hdv (levd0,lond0,latd0), ! horiz diffusion of V from hdif3 | 1,nlevp1,lon0,lon1,lat0,lat1,istep) if (debug) write(6,"('dynamics after uvbnd: istep=',i5)") | istep ! else ! am 8/09 calculation of Un,Vn at LBC (background only) call uvbgrd( | un (levd0,lond0,latd0,itp), | un_nm(levd0,lond0,latd0,itp), | vn (levd0,lond0,latd0,itp), | vn_nm(levd0,lond0,latd0,itp), | w (levd0,lond0,latd0,itp), | z (levd0,lond0,latd0,itp), | barm (levd0,lond0,latd0,itp), | hdu (levd0,lond0,latd0), ! horiz diffusion of U from hdif3 | hdv (levd0,lond0,latd0), ! horiz diffusion of V from hdif3 | 1,nlevp1,lon0,lon1,lat0,lat1,istep) if (debug) write(6,"('dynamics after uvbgrd: istep=',i5)") | istep endif ! ! Update lbc module data u_lbc,v_lbc w/ ubnd,vbnd output of uvbnd: ! add background and perturbations together if this applies call add_uvbnd endif ! ! Advance neutral temperature: ! call timer(time0,time1,'DT',0,0) ! start dt timing call dt( | tn (levd0,lond0,latd0,itp), | tn_nm (levd0,lond0,latd0,itp), | un (levd0,lond0,latd0,itp), | vn (levd0,lond0,latd0,itp), | o2 (levd0,lond0,latd0,itp), | o1 (levd0,lond0,latd0,itp), | ox (levd0,lond0,latd0,itp), | n2 (levd0,lond0,latd0), | barm (levd0,lond0,latd0,itp), ! mean molecular weight | cp (levd0,lond0,latd0), ! specific heat (cpktkm.F) | kt (levd0,lond0,latd0), ! molecular diffusion (cpktkm.F) | km (levd0,lond0,latd0), ! molecular viscosity (cpktkm.F) | hdt (levd0,lond0,latd0), ! horizontal diffusion (hdif.F) | qji_tn(levd0,lond0,latd0), ! joule heating (sub qjoule_tn) | gwt (levd0,lond0,latd0), ! T grav wave tendency from mgw | difkt (levd0,lond0,latd0), ! eddy thermal diffusion | cool_implicit(levd0,lond0,latd0), ! implicit cooling (newton.F) | cool_explicit(levd0,lond0,latd0), ! explicit cooling (newton.F) | w (levd0,lond0,latd0,itc), ! updated W (swdot.F) | tn (levd0,lond0,latd0,itc), ! output updated tn | tn_nm (levd0,lond0,latd0,itc), ! output updated tn at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) if (debug) write(6,"('dynamics after dt: istep=',i5)") istep call timer(time0,time1,'DT',1,0) ! end dt timing if (timing%level >= 2.and.time2print(nstep,istep)) | write(6,"('Time in DT = ',f12.2, | ' Dynamics: step ',i5)") time1,istep if (diags) then do lat=lat0,lat1 call addfld('TN_UPD',' ',' ',tn(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo ! lat=lat0,lat1 endif ! ! Advance neutral velocities: ! call timer(time0,time1,'DUV',0,0) ! start duv timing call duv( | tn (levd0,lond0,latd0,itp), ! neutral temperature | tn (levd0,lond0,latd0,itc), ! updated neutral temperature (dt.F) | tn_nm (levd0,lond0,latd0,itp), ! tn at time n-1 | un (levd0,lond0,latd0,itp), ! zonal velocity | vn (levd0,lond0,latd0,itp), ! meridional velocity | un_nm (levd0,lond0,latd0,itp), ! zonal velocity at time n-1 | vn_nm (levd0,lond0,latd0,itp), ! meridional velocity at time n-1 | w (levd0,lond0,latd0,itc), ! updated vertical velocity (swdot.F) | o2 (levd0,lond0,latd0,itp), ! | ox (levd0,lond0,latd0,itp), ! | n2 (levd0,lond0,latd0), ! molecular nitrogen | barm (levd0,lond0,latd0,itp), ! mean molecular weight | z (levd0,lond0,latd0,itp), ! geopotential height | hdu (levd0,lond0,latd0), ! horizontal diffusion of U (hdif3.F) | hdv (levd0,lond0,latd0), ! horizontal diffusion of V (hdif3.F) | gwu (levd0,lond0,latd0), ! U grav wave tendency from mgw | gwv (levd0,lond0,latd0), ! V grav wave tendency from mgw | difkv (levd0,lond0,latd0), ! eddy diffusion | ui (levd0,lond0,latd0), ! zonal ion velocity input | vi (levd0,lond0,latd0), ! meridional ion velocity input | lxx (levd0,lond0,latd0), ! xx ion drag coefficient | lyy (levd0,lond0,latd0), ! yy ion drag coefficient | lxy (levd0,lond0,latd0), ! xy ion drag coefficient | lyx (levd0,lond0,latd0), ! yx ion drag coefficient | km (levd0,lond0,latd0), ! molecular viscosity (cpktkm.F) | un (levd0,lond0,latd0,itc), ! output updated un | un_nm (levd0,lond0,latd0,itc), ! output updated un at time n-1 | vn (levd0,lond0,latd0,itc), ! output updated vn | vn_nm (levd0,lond0,latd0,itc), ! output updated vn at time n-1 | 1,nlevp1,lon0,lon1,lat0,lat1) if (debug) write(6,"('dynamics after duv: istep=',i5)") istep ! call mpi_finalize(ier) ! stop 'dynamics after duv' call timer(time0,time1,'DUV',1,0) ! end duv timing if (timing%level >= 2.and.time2print(nstep,istep)) | write(6,"('Time in DUV = ',f12.2, | ' Dynamics: step ',i5)") time1,istep if (diags) then do lat=lat0,lat1 call addfld('U_UPD',' ',' ',un(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('V_UPD',' ',' ',vn(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Advance O2, OX, HE: ! subroutine comp_major(tn,o2,o2_nm,ox,ox_nm,he,he_nm,n2,un,vn,w, ! | difkk,hdo2,hdox,hdhe,o2_upd,o2nm_upd,ox_upd,oxnm_upd,he_upd, ! | henm_upd,lev0,lev1,lon0,lon1,lat0,lat1) ! time1_cmpmajor = 0. call timer(time0,time1,'CMPMAJOR',0,0) ! start comp major timing call comp_major( | tn (levd0,lond0,latd0,itp), ! neutral temperature | o2 (levd0,lond0,latd0,itp), ! O2 (mmr) | o2_nm (levd0,lond0,latd0,itp), ! O2 (mmr) at time n-1 | ox (levd0,lond0,latd0,itp), ! OX (mmr) | ox_nm (levd0,lond0,latd0,itp), ! OX (mmr) at time n-1 | he (levd0,lond0,latd0,itp), ! He (mmr) | he_nm (levd0,lond0,latd0,itp), ! He (mmr) at time n-1 | n2 (levd0,lond0,latd0), ! N2 (mmr) input | un (levd0,lond0,latd0,itp), ! zonal velocity | vn (levd0,lond0,latd0,itp), ! meridional velocity | w (levd0,lond0,latd0,itp), ! vertical velocity | difkk (levd0,lond0,latd0), ! eddy viscosity from mgw | hdo2 (levd0,lond0,latd0), ! O2 horizontal diffusion (hdif3) | hdox (levd0,lond0,latd0), ! O horizontal diffusion (hdif3) | hdhe (levd0,lond0,latd0), ! He horizontal diffusion (hdif3) | o2 (levd0,lond0,latd0,itc), ! output: O2 updated for current step | o2_nm (levd0,lond0,latd0,itc), ! output: O2 updated for previous step | ox (levd0,lond0,latd0,itc), ! output: OX updated for current step | ox_nm (levd0,lond0,latd0,itc), ! output: OX updated for previous step | he (levd0,lond0,latd0,itc), ! output: He updated for current step | he_nm (levd0,lond0,latd0,itc), ! output: He updated for previous step | 1,nlevp1,lon0,lon1,lat0,lat1) call timer(time0,time1,'CMPMAJOR',1,0) ! suspend comp major timing time1_cmpmajor = time1_cmpmajor+time1 if (debug) write(6,"('dynamics after comp_major: istep=',i5)") if (diags) then do lat=lat0,lat1 call addfld('UPD_O2',' ',' ',o2(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('UPD_OX',' ',' ',ox(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) enddo endif ! ! Partition ox into o and o3 (this is included in CMPMAJOR timing). ! do lat=lat0,lat1 call timer(time0,time1,'CMPMAJOR',0,0) ! resume comp major timing call compart( | o2 (levd0,lond0,lat,itc), ! updated o2 from comp_major | ox (levd0,lond0,lat,itc), ! updated ox from comp_major | oh (levd0,lond0,lat,itc), ! updated oh from minor_hox | ho2 (levd0,lond0,lat,itc), ! updated ho2 from minor_hox | h (levd0,lond0,lat,itc), ! updated h from minor_hox | no (levd0,lond0,lat,itc), ! updated no from minor_noz | no2 (levd0,lond0,lat,itc), ! updated no2 from minor_noz | o1d (levd0,lond0,lat), ! o(1delta) from comp_meta | cl (levd0,lond0,lat), ! chlorine from comp_solgar | o1 (levd0,lond0,lat,itp), ! current o1 | o3 (levd0,lond0,lat,itp), ! current o3 | barm(levd0,lond0,lat,itp), ! mbar | xnmbar(levd0,lond0,lat), ! | o1 (levd0,lond0,lat,itc), ! o1 output | o3 (levd0,lond0,lat,itc), ! o3 output | 1,nlevp1,lon0,lon1,lat) if (debug) write(6,"('dynamics after compart: lat=',i3, | ' istep=',i5)") lat,istep call timer(time0,time1,'CMPMAJOR',1,0) ! suspend comp major timing time1_cmpmajor = time1_cmpmajor+time1 if (diags) then call addfld('UPD_O1',' ',' ',o1(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) call addfld('UPD_O3',' ',' ',o3(:,i0:i1,lat,itc), | 'lev',1,nk,'lon',i0,i1,lat) endif enddo ! lat=lat0,lat1 if (timing%level >= 2.and.time2print(nstep,istep)) | write(6,"('Time in CMPMAJOR = ',f12.2, | ' Dynamics: step ',i5)") time1_cmpmajor,istep end subroutine dynamics