PROGRAM drivert3dedyn use read_ncfile ,only: & get_geogrid,read_waccm,read_tgcm, & ! subroutines ntime,iday,year,ut, & ! time and date tn,un,vn,omega,zpot,barm,ped,hall ! 3d input fields use namelist_T3DEdyn, only: get_t3dedyn_namelist,nfiles,input_ncfiles,output_ncfiles,timestepsec use fileInput, only: nFilesIO,inNCFiles,outNCFiles,secPerFile,secRunBeg,secRunEnd ! ! Get namelist input ! CALL get_T3DEdynNamelist() ! ! Get input and output file names ! CALL get_IOFileNames() ! ! Read first and last input files to get first time and number of times for this run ! CALL get_InTimeWACCM() ! ! Find number of input time steps for time range of input files ! nTimeSteps = FIX( ( secRunEnd - secRunBeg ) / timestepsec ) runTimeSec = secRunBeg DO iTimeStep = 1, nTimeSteps ! ! Get input GCM data ! CALL get_GCMData(iTimeStep,timestepsec,runTimeSec) ! TIME3D module ! Input: TimeStep: time step; Z: altitude in meter; Tn: neutral temperature in K; ! cO,cO2,cN2,cH,cHe,cN: O, O2, N2, H, He, and N densities in m-3 (in now cersion, [He] is set to be 0 or 1.0D-25) ! Ws,We,Wu: southward, eastward, and upward wind velocities in m/s; ! Qop,QHep,QHp,Qep: production rates of O+,He+,H+, and photoelectron in M-3*s-1 (in now cersion, QHep,QHp are set to be 0 or 1.0D-25) ! OpL,HepL,HpL: loss rates of O+,He+,and H+ in s-1 (in now cersion, HepL,HpL are set to be 0 or 1.0D-25) ! cO2p,cNOp,cN2p: O2+, NO+, and N2+ density in m-3; (above are output from GCM) ! ed1,ed2: electric field index; (output from dynamo) ! output: Op: O+ density in m-3; Te,Ti: electron and ion temperature in K; Vi: anti-field velocity of O+; (input into GCM) ! Ped,Hall: Pedersen and Hall conductivities in s/m; (input into dynamo and GCM) ! zigm11,zigm22,zigm2,zigmc,rim1,rim2; dynamo index in S or S*m/s*Be3 (input into dynamo); the output of [subroutine fieldline_integrals] in dynamo CALL TIME3Dmodule(timestepsec,Z,Tn,cO,cO2,cN2,cH,cHe,cN,Ws,We, Wu,Qop,QHep,OpL,Qep,QHp,HepL,HpL,cO2p,cNOp, cN2p,ed1,ed2,Op,Te,Ti,Vi,Ped,Hall,zigm11, zigm22,zigm2,zigmc,rim1,rim2,t) runTimeSec = runTimeSec + timestepsec ENDDO use geogrid_t3d ,only: nlev,nlat,nlon use maggrid ,only: nmlat,nmlon,nmlonp1 use edynamo_driver,only: edyn_init,edyn_driver INCLUDE 'EIG.h' COMMON/input/ F107,F107A,AP,Uth,TimeStep,Hindex(6),iyd,IYR,IMO,IDA COMMON/GRIDI/ Calti(NLi),Cli(NLi),SIi(NLi),MLATi(lmi+1,NLi),Mloni( & NFi),altH(NHN),ClH(NHN),SIH(NHN),MlatH(NHN,2,NLi), & MlonH(NHN) real*8,allocatable,dimension(:,:,:) :: Z,Tn,cO,cO2,cN2,cH,cHe,cN, & Ws,We,Wu,Qop,QHep,OpL,Qep,QHp,HepL,HpL, & cO2p,cNOp,cN2p,Op,Te,Ti,Vi,Ped,Hall REAL*8 ed1(nmlat,nmlonp1),ed2(nmlat,nmlonp1),zigm11(nmlat,nmlonp1 & ),zigm22(nmlat,nmlonp1),zigm2(nmlat,nmlonp1),zigmc(nmlat, & nmlonp1),rim1(nmlat,nmlonp1),rim2(nmlat,nmlonp1) ! ! Set use_conduct = 'edynamo' to have edynamo do fieldline ! integrations, ! or set use_conduct = 'time3d' to have TIME3D do fieldline ! integrations. ! In the latter case, time3d will pass 2d conductances to edynamo. ! ! character(len=16) :: use_conduct='edynamo ' character(len=16) :: use_conduct='time3d ' write(6,"('time3d: use_conduct=',a)") use_conduct CALL PreGrid() !Set parameters for geogrid, remove when coupling with Dynamo allocate(Z(nlat,nlon,Nlev)) allocate(Tn(nlat,nlon,Nlev)) allocate(cO(nlat,nlon,Nlev)) allocate(cO2(nlat,nlon,Nlev)) allocate(cN2(nlat,nlon,Nlev)) allocate(cH(nlat,nlon,Nlev)) allocate(cHe(nlat,nlon,Nlev)) allocate(cN(nlat,nlon,Nlev)) allocate(Ws(nlat,nlon,Nlev)) allocate(We(nlat,nlon,Nlev)) allocate(Wu(nlat,nlon,Nlev)) allocate(Qop(nlat,nlon,Nlev)) allocate(QHep(nlat,nlon,Nlev)) allocate(OpL(nlat,nlon,Nlev)) allocate(Qep(nlat,nlon,Nlev)) allocate(QHp(nlat,nlon,Nlev)) allocate(HepL(nlat,nlon,Nlev)) allocate(HpL(nlat,nlon,Nlev)) allocate(cO2p(nlat,nlon,Nlev)) allocate(cNOp(nlat,nlon,Nlev)) allocate(cN2p(nlat,nlon,Nlev)) allocate(Op(nlat,nlon,Nlev)) allocate(Te(nlat,nlon,Nlev)) allocate(Ti(nlat,nlon,Nlev)) allocate(Vi(nlat,nlon,Nlev)) allocate(Ped(nlat,nlon,Nlev)) allocate(Hall(nlat,nlon,Nlev)) IYR=2001 !year iyd=1 !doy of year ! ! Timestep in seconds (default is 300, or 5 minutes): ! TimeStep=3.0D2 ! time step TimeStep=6.0D2 CALL GCMpara() !Set parameters for GCM, remove when coupling Cwith GCM model call edyn_init CALL TIME3DParaGrid(IYR,iyd,HB) ! Set TIME3D Parameters and CGrid System Do ifsn=3,3 C--------------------Parameters for GCM and dynamo, remove when coupling Cwith GCM and dynamo model-------------------- F107=65.0D0 F107A=F107 AP=4.0D0 iyd=int(30.5D0*ifsn)-9 !87 IYR=2001 UTh=0.0D0 Hindex(1)=0.0D0 Hindex(2)=-0.5D0 Hindex(3)=3.0D2 Hindex(4:5)=1.0D1 Hindex(6)=2.0D1 ih=3600.0D0/Timestep Nnum=ih*24 write(*,*)ifsn,F107,iyd CALL DoytoMonDay(IYR,iyd,IMO,IDA) CALL EUVFluxCrossSection(F107,F107A,iyd) CALL Prealt(Z) C----------------------------------------------------------------------------- CALL TIME3DInitialCondition(iyd,F107,F107A,AP,Uth) !TIME3D CInitilze CALL TIME3D2GCM(Z,Op,Te,Ti,Vi,Ped,Hall) ! output Op,Te,Ti,Vi for CGCM code ed1=0.0D0 ed2=0.0D0 ! DO iday=1,3 DO iday=1,1 DO i=1,Nnum UTh=DBLE(i)*Timestep/3.6D3 IF(UTh.GT.2.4d1)UTh=UTh-2.4d1 ! Instead of GCM code, remove when coupling with GCM ! input: Op: O+ density in m-3; Te,Ti: electron and ion temperature in ! K ! output: Z: altitude in meter; Tn: neutral temperature in K; ! cO,cO2,cN2,cH,cHe,cN: O, O2, N2, H, He, and N densities in m-3 ! (in now cersion, [He] is set to be 0 or 1.0D-25) ! Ws,We,Wu: southward, eastward, and upward wind velocities in ! m/s; ! Qop,QHep,QHp,Qep: production rates of O+,He+,H+, and ! photoelectron in M-3*s-1 (in now cersion, QHep,QHp are set to ! be 0 or 1.0D-25) ! OpL,HepL,HpL: loss rates of O+,He+,and H+ in s-1 (in now ! cersion, HepL,HpL are set to be 0 or 1.0D-25) ! cO2p,cNOp,cN2p: O2+, NO+, and N2+ density in m-3 CALL GCMsim(Z,Tn,cO,cO2,cN2,cH,cHe,cN,Ws,We,Wu,Qop,QHep,OpL, & Qep,QHp,HepL,HpL,cO2p,cNOp,cN2p,Op,Te,Ti) ! ! TIME3D module ! Input: TimeStep: time step; Z: altitude in meter; Tn: neutral ! temperature in K; ! cO,cO2,cN2,cH,cHe,cN: O, O2, N2, H, He, and N densities in m-3 ! (in now cersion, [He] is set to be 0 or 1.0D-25) ! Ws,We,Wu: southward, eastward, and upward wind velocities in ! m/s; ! Qop,QHep,QHp,Qep: production rates of O+,He+,H+, and ! photoelectron in M-3*s-1 (in now cersion, QHep,QHp are set to ! be 0 or 1.0D-25) ! OpL,HepL,HpL: loss rates of O+,He+,and H+ in s-1 (in now ! cersion, HepL,HpL are set to be 0 or 1.0D-25) ! cO2p,cNOp,cN2p: O2+, NO+, and N2+ density in m-3; (above are ! output from GCM) ! ed1,ed2: electric field index; (output from dynamo) ! output: Op: O+ density in m-3; Te,Ti: electron and ion temperature in ! K; Vi: anti-field velocity of O+; (input into GCM) ! Ped,Hall: Pedersen and Hall conductivities in s/m; (input into ! dynamo and GCM) ! zigm11,zigm22,zigm2,zigmc,rim1,rim2; dynamo index in S or ! S*m/s*Be3 (input into dynamo); the output of [subroutine ! fieldline_integral s] in dynamo CALL TIME3Dmodule(TimeStep,Z,Tn,cO,cO2,cN2,cH,cHe,cN,Ws,We, & Wu,Qop,QHep,OpL,Qep,QHp,HepL,HpL,cO2p,cNOp & ,cN2p,ed1,ed2,Op,Te,Ti,Vi,Ped,Hall,zigm11, & zigm22,zigm2,zigmc,rim1,rim2,t) !