! module fields_module ! ! There are allocatable arrays and an array of field structures for ! 3-d and 4-d fields. Subdomains are allocated to the allocatable ! arrays, and structure data pointers point to the allocatable arrays, ! e.g. f4d(i_tn)%data(:,:,:,:) => tn(:,:,:,:), where the 4 dimensions ! are (pressure,longitude,latitude,2). The last dimension is for ! previous and current time steps. 3-d fields are the same as 4-d ! without the final time dimension. ! ! btf 6/11/04: Broke sub init_fields into subs init_f4d and init_f3d. ! This avoids error from pgf90 compiler (s.a. init.F). ! use params_module,only: nlonp1,nlonp4,nlat,nlatp1,nlevp1, | mxfsech,spval,nmlat,nmlonp1,nmlev,nmagphrlon,nmagphrlat implicit none save integer,parameter :: | nf4d = 52, ! number of 4-d fields | nf4d_hist = 50, ! number of 4-d fields on primary histories | nf3d = 84 ! number of 3-d fields real,parameter :: | field_initval = 0. ! initialization value for fields data integer,parameter :: | longname_len = 80, ! length of field long name | shortname_len = 16, ! length of field short name | units_len = 16 ! length of field units attribute integer :: levd0,levd1, lond0,lond1, latd0,latd1 ! subdomain dimensions ! ! 4-d fields structure type: ! Data pointers will point to allocatable target arrays. ! type fields_4d character(len=longname_len) :: long_name character(len=shortname_len) :: short_name character(len=units_len) :: units character(len=16) :: vcoord ! midpoints or interfaces logical :: | magnetic, ! true if field is on magnetic 3d grid | magnetos, ! true if field is on magnetospheric grid | prognostic, ! true if prognostic (diagnostic if false) | mpi ! flag used by some mpi routines real :: polesign real,pointer :: data(:,:,:,:) end type fields_4d type (fields_4d),save :: f4d(nf4d) ! ! Indices, to f4d fields, e.g. f4d(i_tn)%data is neutral temperature. integer :: | i_tn ,i_un ,i_vn ,i_o2 ,i_ox , | i_n4s ,i_noz ,i_co ,i_co2 ,i_h2o , | i_h2 ,i_hox ,i_op ,i_ch4 ,i_ar , | i_he ,i_nat ,i_o21d ,i_no2 ,i_no , | i_o3 ,i_o1 ,i_oh ,i_ho2 ,i_h , | i_n2d ,i_ti ,i_te ,i_ne ,i_o2p , | i_w ,i_z ,i_poten ,i_tn_nm ,i_un_nm , | i_vn_nm ,i_o2_nm ,i_ox_nm ,i_n4s_nm ,i_noz_nm , | i_co_nm ,i_co2_nm ,i_h2o_nm ,i_h2_nm ,i_hox_nm , | i_op_nm ,i_ch4_nm ,i_ar_nm ,i_he_nm ,i_nat_nm , | i_barm ,i_vc ! ! Allocatable target arrays for 4-d fields: ! Allocation will be: tn(levd0:levd1,lond0 :lond1 ,latd0 :latd1 ,2) ! which is the same as: tn(1 :nlev ,lon0-2:lon1+2,lat0-2:lat1+2,2) ! real,dimension(:,:,:,:),target,allocatable,save :: ! (k,i,j,2) | tn ,un ,vn ,o2 ,ox , | n4s ,noz ,co ,co2 ,h2o , | h2 ,hox ,op ,ch4 ,ar , | he ,nat ,o21d ,no2 ,no , | o3 ,o1 ,oh ,ho2 ,h , | n2d ,ti ,te ,ne ,o2p , | w ,z ,poten ,tn_nm ,un_nm , | vn_nm ,o2_nm ,ox_nm ,n4s_nm ,noz_nm , | co_nm ,co2_nm ,h2o_nm ,h2_nm ,hox_nm , | op_nm ,ch4_nm ,ar_nm ,he_nm ,nat_nm , | barm ,vc ! ! 3-d fields structure type: ! Data pointers will point to allocatable target arrays. ! type fields_3d character(len=longname_len) :: long_name character(len=shortname_len) :: short_name character(len=units_len) :: units character(len=16) :: vcoord ! midpoints or interfaces logical :: | magnetic, ! true if field is on magnetic 3d grid | magnetos, ! true if field is on magnetospheric grid | prognostic, ! true if prognostic (diagnostic if false) | mpi, ! flag used by some mpi routines | task0_only ! if true, field was defined at root task only real :: polesign real,pointer :: data(:,:,:) end type fields_3d type (fields_3d),save :: f3d(nf3d) ! ! 3-d fields and indices, e.g., reference f3d(kldt)%data integer :: | i_kldt ,i_kldu ,i_kldv ,i_kldo2 ,i_kldox , | i_kldn4s ,i_kldnoz ,i_kldco ,i_kldco2,i_kldh2o, | i_kldh2 ,i_kldhox ,i_kldch4 ,i_kldar ,i_kldhe , | i_kldnat ,i_xnmbar ,i_xnmbari ,i_xnmbarm ,i_cp , | i_kt ,i_km ,i_ui ,i_vi ,i_wi , | i_vo2 ,i_vo1 ,i_vn2 ,i_vo3 ,i_vco2, i_vno, | i_sco2 ,i_sco1 ,i_scn2 ,i_sco3 ,i_scco2, i_scno, | i_xiop2p ,i_xiop2d ,i_hplus ,i_nplus ,i_n2p , | i_nop ,i_lxx ,i_lyy ,i_lxy ,i_lyx , | i_qji_ti ,i_qji_tn ,i_cool_implicit,i_cool_explicit,i_hdt, | i_hdu ,i_hdv ,i_hdo2 ,i_hdox ,i_hdn4s, | i_hdnoz ,i_hdco ,i_hdco2 ,i_hdh2o ,i_hdh2 , | i_hdhox ,i_hdch4 ,i_hdar ,i_hdhe ,i_hdnat, | i_ped ,i_hall ,i_nocool ,i_n2o ,i_cl , | i_clo ,i_o1d ,i_o21s ,i_gwv ,i_gwu , | i_gwt ,i_h2o2 ,i_phoxic,i_difkk ,i_difkt, | i_difkv ,i_lam1 ! ! Allocatable target arrays 3-d: ! Allocation will be: tn(levd0:levd1,lond0 :lond1 ,latd0 :latd1 ) ! which is the same as: tn(1 :nlev ,lon0-2:lon1+2,lat0-2:lat1+2) ! real,dimension(:,:,:),target,allocatable,save :: ! (k,i,j) | kldt ,kldu ,kldv ,kldo2 ,kldox , | kldn4s ,kldnoz ,kldco ,kldco2 ,kldh2o , | kldh2 ,kldhox ,kldch4 ,kldar ,kldhe , | kldnat ,xnmbar ,xnmbari ,xnmbarm ,cp , | kt ,km ,ui ,vi ,wi , | vo2 ,vo1 ,vn2 ,vo3 ,vco2 , vno, | sco2 ,sco1 ,scn2 ,sco3 ,scco2 , scno, | xiop2p ,xiop2d ,hplus ,nplus ,n2p , | nop ,lxx ,lyy ,lxy ,lyx , | qji_ti ,qji_tn ,cool_implicit,cool_explicit,hdt, | hdu ,hdv ,hdo2 ,hdox ,hdn4s , | hdnoz ,hdco ,hdco2 ,hdh2o ,hdh2 , | hdhox ,hdch4 ,hdar ,hdhe ,hdnat , | ped ,hall ,nocool ,n2o ,cl , | clo ,o1d ,o21s ,gwu ,gwv , | gwt ,h2o2 ,phoxic ,difkk ,difkt , | difkv ,lam1 ! ! 2-d field type (used for 2d secondary history fields): type fields_2d character(len=longname_len) :: long_name character(len=shortname_len) :: short_name character(len=units_len) :: units logical :: | magnetic, ! true if field is on magnetic 3d grid | magnetos, ! true if field is on magnetospheric grid | prognostic, ! true if prognostic (diagnostic if false) | mpi, ! flag used by some mpi routines | task0_only ! if true, field was defined at root task only real :: polesign real,pointer :: data(:,:) ! (k,i) end type fields_2d ! ! Electric potential on geographic and magnetic grids: ! (full domains until dynamo is parallelized) real :: | dynpot(nlonp1,0:nlatp1,nlevp1), ! 3d electric potential geographic | phim3d(nmlonp1,nmlat,nmlev) ! 3d electric potential magnetic ! ! Secondary history field structures. ! These are initialized in sub init_fsech, and data is set in addfld.F. type fields_sech character(len=longname_len) :: long_name character(len=shortname_len) :: short_name character(len=units_len) :: units logical :: task0_only logical :: prognostic logical :: mag,geo integer :: ndims ! 2d or 3d character(len=8) :: dimnames(3) ! dims on history integer :: dimsizes(3) ! dim sizes on history real,pointer :: data(:,:,:) ! allocated by addfld end type fields_sech type(fields_sech) :: fsechist(mxfsech) ! ! Time indices for rightmost dimension of 4d data, itp for current timestep, ! itc for next timestep. Fields at the previous timestep (time n-1) are ! saved at both itp and itc (e.g., tn_nm, un_nm, etc). ! integer :: itc,itp ! ! If fakeflds is true, use fake dimensions for fields ! (for dry runs, testing, etc) ! logical,parameter :: fakeflds = .false. ! real,dimension(nlevp1,nlonp4+1,-2:nlat) :: | fnrh, ! eddy viscosity | fkmh ! M/T ! ! Lower boundary conditions (bottom interface level) for t,u,v from ! source history, and at current model time: real,allocatable,dimension(:,:),save :: | tlbc, ulbc, vlbc, ! subdomains (lond0:lond1,latd0:latd1) | tlbc_nm, ulbc_nm, vlbc_nm, ! subdomains (lond0:lond1,latd0:latd1) | tlbc_glb, ulbc_glb, vlbc_glb, ! global (nlonp4,nlat) | tlbc_nm_glb, ulbc_nm_glb, vlbc_nm_glb ! global (nlonp4,nlat) real,allocatable,dimension(:,:),save :: | t_bgrd, z_bgrd,u_bgrd, v_bgrd, ! background | u_prt, v_prt, u_prt_tn, v_prt_tn,u_prt_tnm, v_prt_tnm ! perturbations ! ! Full 3d grid with all primary history fields for writing to netcdf ! history files. This will be allocated only on the root task ! (see allocdata): real,allocatable,save :: foutput(:,:,:,:) ! (nlevp1,nlonp4,nlat,nf4d_hist) real,allocatable,save :: fzg(:,:,:) ! (nlevp1,nlonp4,nlat) ! contains !----------------------------------------------------------------------- ! subroutine init_4d(lon0,lon1,lat0,lat1,mytid,iprint) ! ! Set names, units, indices and pointers for f4d and f3d field structures, ! and allocate 3d and 4d field arrays. Also make other data allocations. ! implicit none ! ! Args: integer,intent(in) :: lon0,lon1,lat0,lat1,mytid,iprint ! ! Local: integer :: i,n,istat ! ! Fields are allocated at full task subdomain, including ghost cells: levd0 = 1 ; levd1 = nlevp1 lond0 = lon0-2 ; lond1 = lon1+2 latd0 = lat0-2 ; latd1 = lat1+2 !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - n = 1 ! ! 4-d fields (long and short names, units, index): ! (pointer definition must follow allocate statement) ! Note it is difficult to make a subroutine to do the allocations because ! allocatable arrays cannot be dummy arguments. (Once they are allocated, ! they can be passed to subroutines and declared as real, but subroutine ! dummy arguments cannot be declared as allocatable arrays.) ! allocate(tn(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NEUTRAL TEMPERATURE')") f4d(n)%short_name = "TN" f4d(n)%units = "K" f4d(n)%vcoord = "midpoints" f4d(n)%data => tn i_tn = n ; n = n+1 ! allocate(un(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NEUTRAL ZONAL WIND')") f4d(n)%short_name = "UN" f4d(n)%units = "CM/S" f4d(n)%vcoord = "midpoints" f4d(n)%data => un i_un = n ; n = n+1 ! allocate(vn(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NEUTRAL MERIDIONAL WIND')") f4d(n)%short_name = "VN" f4d(n)%units = "CM/S" f4d(n)%vcoord = "midpoints" f4d(n)%data => vn i_vn = n ; n = n+1 ! allocate(o2(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('MOLECULAR OXYGEN')") f4d(n)%short_name = ("O2") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => o2 i_o2 = n ; n = n+1 ! allocate(ox(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('OX GROUP (O+O3)')") f4d(n)%short_name = ("OX") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => ox i_ox = n ; n = n+1 ! allocate(n4s(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('N4S')") f4d(n)%short_name = ("N4S") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => n4s i_n4s = n ; n = n+1 ! allocate(noz(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NO GROUP (NO+NO2)')") f4d(n)%short_name = ("NOZ") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => noz i_noz = n ; n = n+1 ! allocate(co(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('CARBON MONOXIDE (CO)')") f4d(n)%short_name = ("CO") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => co i_co = n ; n = n+1 ! allocate(co2(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) co2 = 0. write(f4d(n)%long_name,"('CARBON DIOXIDE (CO2)')") f4d(n)%short_name = ("CO2") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => co2 i_co2 = n ; n = n+1 ! allocate(h2o(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('WATER VAPOR (H2O)')") f4d(n)%short_name = ("H2O") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => h2o i_h2o = n ; n = n+1 ! allocate(h2(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('MOLECULAR HYDROGEN (H2)')") f4d(n)%short_name = ("H2") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => h2 i_h2 = n ; n = n+1 allocate(hox(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('HYDROGEN GROUP (OH+HO2+H)')") f4d(n)%short_name = ("HOX") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => hox i_hox = n ; n = n+1 ! allocate(op(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('O+ ION')") f4d(n)%short_name = ("OP") ! 11/13/08 btf: O+ units should be cm^3, and at midpoints: f4d(n)%units = "CM^3" f4d(n)%vcoord = "midpoints" f4d(n)%data => op i_op = n ; n = n+1 ! allocate(ch4(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('METHANE (CH4)')") f4d(n)%short_name = ("CH4") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => ch4 i_ch4 = n ; n = n+1 ! allocate(ar(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('ARGON (AR)')") f4d(n)%short_name = ("AR") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => ar i_ar = n ; n = n+1 ! allocate(he(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('HELIUM (HE)')") f4d(n)%short_name = ("HE") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => he i_he = n ; n = n+1 ! allocate(nat(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('TOTAL SODIUM (NAT)')") f4d(n)%short_name = ("NAT") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => nat i_nat = n ; n = n+1 ! allocate(o21d(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('O21D')") f4d(n)%short_name = ("O21D") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => o21d i_o21d = n ; n = n+1 ! allocate(no2(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NITROGEN DIOXIDE( NO2)')") f4d(n)%short_name = ("NO2") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => no2 i_no2 = n ; n = n+1 ! allocate(no(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NITRIC OXIDE (NO)')") f4d(n)%short_name = ("NO") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => no i_no = n ; n = n+1 ! allocate(o3(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('OZONE (O3)')") f4d(n)%short_name = ("O3") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => o3 i_o3 = n ; n = n+1 ! allocate(o1(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('ATOMIC OXYGEN')") f4d(n)%short_name = ("O1") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => o1 i_o1 = n ; n = n+1 ! allocate(oh(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('OH')") f4d(n)%short_name = ("OH") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => oh i_oh = n ; n = n+1 ! allocate(ho2(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('HYDROGEN DIOXIDE (HO2)')") f4d(n)%short_name = ("HO2") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => ho2 i_ho2 = n ; n = n+1 ! allocate(h(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('HYDROGEN (H)')") f4d(n)%short_name = ("H") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => h i_h = n ; n = n+1 ! allocate(n2d(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('N2D')") f4d(n)%short_name = ("N2D") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => n2d i_n2d = n ; n = n+1 ! allocate(ti(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('ION TEMPERATURE')") f4d(n)%short_name = ("TI") f4d(n)%units = "DEG K" f4d(n)%vcoord = "midpoints" f4d(n)%data => ti i_ti = n ; n = n+1 ! allocate(te(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('ELECTRON TEMPERATURE')") f4d(n)%short_name = ("TE") f4d(n)%units = "DEG K" f4d(n)%vcoord = "midpoints" f4d(n)%data => te i_te = n ; n = n+1 ! allocate(ne(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('ELECTRON DENSITY')") f4d(n)%short_name = ("NE") f4d(n)%units = "CM^3" f4d(n)%vcoord = "interfaces" f4d(n)%data => ne i_ne = n ; n = n+1 ! allocate(o2p(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('O2+ ION')") f4d(n)%short_name = ("O2P") f4d(n)%units = "CM^3" f4d(n)%vcoord = "midpoints" f4d(n)%data => o2p i_o2p = n ; n = n+1 ! ! 5/16/06 btf: "W" is replaced by "OMEGA": ! ! allocate(w(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) ! write(f4d(n)%long_name,"('VERTICAL VELOCITY (PLUS UP)')") ! f4d(n)%short_name = ("W") ! f4d(n)%units = "CM/S" ! f4d(n)%data => w ! i_w = n ; n = n+1 ! allocate(w(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('VERTICAL MOTION')") f4d(n)%short_name = ("OMEGA") f4d(n)%units = "s-1" f4d(n)%vcoord = "interfaces" f4d(n)%data => w i_w = n ; n = n+1 ! allocate(z(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('GEOPOTENTIAL HEIGHT')") f4d(n)%short_name = ("Z") f4d(n)%units = "CM" f4d(n)%vcoord = "interfaces" f4d(n)%data => z i_z = n ; n = n+1 ! allocate(poten(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('ELECTRIC POTENTIAL')") f4d(n)%short_name = ("POTEN") f4d(n)%units = "VOLTS" f4d(n)%vcoord = "interfaces" f4d(n)%data => poten i_poten = n ; n = n+1 ! allocate(tn_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NEUTRAL TEMPERATURE (TIME N-1)')") f4d(n)%short_name = ("TN_NM") f4d(n)%units = "DEG K" f4d(n)%vcoord = "midpoints" f4d(n)%data => tn_nm i_tn_nm = n ; n = n+1 ! allocate(un_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NEUTRAL ZONAL WIND (TIME N-1)')") f4d(n)%short_name = ("UN_NM") f4d(n)%units = "CM/S" f4d(n)%vcoord = "midpoints" f4d(n)%data => un_nm i_un_nm = n ; n = n+1 ! allocate(vn_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NEUTRAL MERIDIONAL WIND (TIME N-1)')") f4d(n)%short_name = ("VN_NM") f4d(n)%units = "CM/S" f4d(n)%vcoord = "midpoints" f4d(n)%data => vn_nm i_vn_nm = n ; n = n+1 ! allocate(o2_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('MOLECULAR OXYGEN (TIME N-1)')") f4d(n)%short_name = ("O2_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => o2_nm i_o2_nm = n ; n = n+1 ! allocate(ox_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('OX GROUP (TIME N-1)')") f4d(n)%short_name = ("OX_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => ox_nm i_ox_nm = n ; n = n+1 ! allocate(n4s_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('N4S (TIME N-1)')") f4d(n)%short_name = ("N4S_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => n4s_nm i_n4s_nm = n ; n = n+1 ! allocate(noz_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NOZ (TIME N-1)')") f4d(n)%short_name = ("NOZ_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => noz_nm i_noz_nm = n ; n = n+1 ! allocate(co_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('CO (TIME N-1)')") f4d(n)%short_name = ("CO_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => co_nm i_co_nm = n ; n = n+1 ! allocate(co2_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('CO2 (TIME N-1)')") f4d(n)%short_name = ("CO2_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => co2_nm i_co2_nm = n ; n = n+1 ! allocate(h2o_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('H2O (TIME N-1)')") f4d(n)%short_name = ("H2O_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => h2o_nm i_h2o_nm = n ; n = n+1 ! allocate(h2_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('H2 (TIME N-1)')") f4d(n)%short_name = ("H2_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => h2_nm i_h2_nm = n ; n = n+1 ! allocate(hox_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('HOX (TIME N-1)')") f4d(n)%short_name = ("HOX_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => hox_nm i_hox_nm = n ; n = n+1 ! allocate(op_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('OP (TIME N-1)')") f4d(n)%short_name = ("OP_NM") f4d(n)%units = "CM^3" f4d(n)%vcoord = "midpoints" f4d(n)%data => op_nm i_op_nm = n ; n = n+1 ! allocate(ch4_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('CH4 (TIME N-1)')") f4d(n)%short_name = ("CH4_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => ch4_nm i_ch4_nm = n ; n = n+1 ! allocate(ar_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('AR (TIME N-1)')") f4d(n)%short_name = ("AR_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => ar_nm i_ar_nm = n ; n = n+1 ! allocate(he_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('HE (TIME N-1)')") f4d(n)%short_name = ("HE_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => he_nm i_he_nm = n ; n = n+1 ! allocate(nat_nm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('NAT (TIME N-1)')") f4d(n)%short_name = ("NAT_NM") f4d(n)%units = "MMR" f4d(n)%vcoord = "midpoints" f4d(n)%data => nat_nm i_nat_nm = n ; n = n+1 ! allocate(barm(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('MEAN MOLECULAR WEIGHT')") f4d(n)%short_name = ("BARM") f4d(n)%units = ' ' f4d(n)%vcoord = "midpoints" f4d(n)%data => barm i_barm = n ; n = n+1 ! allocate(vc(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) write(f4d(n)%long_name,"('COS(PHI)*VN')") f4d(n)%short_name = ("VC") f4d(n)%units = ' ' f4d(n)%vcoord = "midpoints" f4d(n)%data => vc i_vc = n ; n = n+1 ! ! Check number of fields allocated: if (n-1 /= nf4d) then write(6,"(/,'>>> init_4d: wrong number of 4-d fields?', | ' n-1=',i5,' nf4d=',i5)") n-1,nf4d endif ! ! itp,itc are indices to rightmost dimension of field data, for ! previous and current time steps. itp = 1 itc = 2 ! ! Set polesign for crossing the poles (see mk_polelats and mp_bndlats) f4d%polesign = 1. f4d(i_un)%polesign = -1. f4d(i_un_nm)%polesign = -1. f4d(i_vn)%polesign = -1. f4d(i_vn_nm)%polesign = -1. f4d(i_n2d)%polesign = 0. f4d(i_ne )%polesign = 0. f4d(i_o2p)%polesign = 0. ! ! All f4d fields are on geographic grid: f4d%magnetic = .false. f4d%magnetos = .false. ! ! f4d fields are "prognostic": f4d%prognostic = .true. ! ! Data init: do i=1,nf4d f4d(i)%data = field_initval enddo ! ! Report to stdout: if (iprint > 0) then write(6,"(/,'There are ',i3,' 4-d fields:')") nf4d do n=1,nf4d call print_f4d(f4d(n),n) enddo write(6,"('Initialized f4d%data to ',e12.4)") field_initval endif end subroutine init_4d !----------------------------------------------------------------------- ! subroutine init_3d(lon0,lon1,lat0,lat1,mytid,iprint) ! ! Set names, units, indices and pointers for f4d and f3d field structures, ! and allocate 3d and 4d field arrays. Also make other data allocations. ! implicit none ! ! Args: integer,intent(in) :: lon0,lon1,lat0,lat1,mytid,iprint ! ! Local: integer :: i,n,istat ! !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ! n = 1 ! ! 3-d fields (long and short names, units, index): ! allocate(kldt(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('TN HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDT" f3d(n)%units = " " f3d(n)%data => kldt i_kldt = n ; n = n+1 ! allocate(kldu(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('UN HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDU" f3d(n)%units = " " f3d(n)%data => kldu i_kldu = n ; n = n+1 ! allocate(kldv(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('VN HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDV" f3d(n)%units = " " f3d(n)%data => kldv i_kldv = n ; n = n+1 ! allocate(kldo2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O2 HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDO2" f3d(n)%units = " " f3d(n)%data => kldo2 i_kldo2 = n ; n = n+1 ! allocate(kldox(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('OX HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDOX" f3d(n)%units = " " f3d(n)%data => kldox i_kldox = n ; n = n+1 ! allocate(kldn4s(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('N4S HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDN4S" f3d(n)%units = " " f3d(n)%data => kldn4s i_kldn4s = n ; n = n+1 ! allocate(kldnoz(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('NOZ HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDNOZ" f3d(n)%units = " " f3d(n)%data => kldnoz i_kldnoz = n ; n = n+1 ! allocate(kldco(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CO HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDCO" f3d(n)%units = " " f3d(n)%data => kldco i_kldco = n ; n = n+1 ! allocate(kldco2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CO2 HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDCO2" f3d(n)%units = " " f3d(n)%data => kldco2 i_kldco2 = n ; n = n+1 ! allocate(kldh2o(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('H2O HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDH2O" f3d(n)%units = " " f3d(n)%data => kldh2o i_kldh2o = n ; n = n+1 ! allocate(kldh2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('H2 HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDH2" f3d(n)%units = " " f3d(n)%data => kldh2 i_kldh2 = n ; n = n+1 ! allocate(kldhox(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('HOX HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDHOX" f3d(n)%units = " " f3d(n)%data => kldhox i_kldhox = n ; n = n+1 ! allocate(kldch4(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CH4 HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDCH4" f3d(n)%units = " " f3d(n)%data => kldch4 i_kldch4 = n ; n = n+1 ! allocate(kldar(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('AR HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDAR" f3d(n)%units = " " f3d(n)%data => kldar i_kldar = n ; n = n+1 ! allocate(kldhe(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('HE HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDHE" f3d(n)%units = " " f3d(n)%data => kldhe i_kldhe = n ; n = n+1 ! allocate(kldnat(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('NAT HORIZONTAL DIFFUSION COEFF')") f3d(n)%short_name = "KLDNAT" f3d(n)%units = " " f3d(n)%data => kldnat i_kldnat = n ; n = n+1 ! allocate(xnmbar(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('p0*e(-z)/kT*barm')") f3d(n)%short_name = ("XNMBAR") f3d(n)%units = ' ' f3d(n)%data => xnmbar i_xnmbar = n ; n = n+1 ! allocate(xnmbari(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('p0*e(-z)/kT*barm at interfaces')") f3d(n)%short_name = ("XNMBARI") f3d(n)%units = ' ' f3d(n)%data => xnmbari i_xnmbari = n ; n = n+1 ! allocate(xnmbarm(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('p0*e(-z)/kT*barm at midpoints')") f3d(n)%short_name = ("XNMBARM") f3d(n)%units = ' ' f3d(n)%data => xnmbarm i_xnmbarm = n ; n = n+1 ! allocate(cp(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('SPECIFIC HEAT')") f3d(n)%short_name = "CP" f3d(n)%units = " " f3d(n)%data => cp i_cp = n ; n = n+1 ! allocate(kt(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('MOLECULAR THERMAL CONDUCTIVITY')") f3d(n)%short_name = "KT" f3d(n)%units = " " f3d(n)%data => kt i_kt = n ; n = n+1 ! allocate(km(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('MOLECULAR DIFFUSION')") f3d(n)%short_name = "KM" f3d(n)%units = " " f3d(n)%data => km i_km = n ; n = n+1 ! allocate(ui(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('ZONAL ION DRIFT VELOCITY')") f3d(n)%short_name = "UI" f3d(n)%units = " " f3d(n)%data => ui i_ui = n ; n = n+1 ! allocate(vi(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('MERIDIONAL ION DRIFT VELOCITY')") f3d(n)%short_name = "VI" f3d(n)%units = " " f3d(n)%data => vi i_vi = n ; n = n+1 ! allocate(wi(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('VERTICAL ION DRIFT VELOCITY')") f3d(n)%short_name = "WI" f3d(n)%units = " " f3d(n)%data => wi i_wi = n ; n = n+1 ! allocate(vo2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O2 VERTICAL LINE INTEGRAL')") f3d(n)%short_name = "VO2" f3d(n)%units = " " f3d(n)%data => vo2 i_vo2 = n ; n = n+1 ! allocate(vo1(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O1 VERTICAL LINE INTEGRAL')") f3d(n)%short_name = "VO1" f3d(n)%units = " " f3d(n)%data => vo1 i_vo1 = n ; n = n+1 ! allocate(vn2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('N2 VERTICAL LINE INTEGRAL')") f3d(n)%short_name = "VN2" f3d(n)%units = " " f3d(n)%data => vn2 i_vn2 = n ; n = n+1 ! allocate(vo3(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O3 VERTICAL LINE INTEGRAL')") f3d(n)%short_name = "VO3" f3d(n)%units = " " f3d(n)%data => vo3 i_vo3 = n ; n = n+1 ! allocate(vco2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CO2 VERTICAL LINE INTEGRAL')") f3d(n)%short_name = "CO2" f3d(n)%units = " " f3d(n)%data => vco2 i_vco2 = n ; n = n+1 ! allocate(vno(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('NO VERTICAL LINE INTEGRAL')") f3d(n)%short_name = "NO" f3d(n)%units = " " f3d(n)%data => vno i_no = n ; n = n+1 ! allocate(sco2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O2 SLANT-LINE INTEGRAL')") f3d(n)%short_name = "SCO2" f3d(n)%units = " " f3d(n)%data => sco2 i_sco2 = n ; n = n+1 ! allocate(sco1(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O1 SLANT-LINE INTEGRAL')") f3d(n)%short_name = "SCO1" f3d(n)%units = " " f3d(n)%data => sco1 i_sco1 = n ; n = n+1 ! allocate(scn2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('N2 SLANT-LINE INTEGRAL')") f3d(n)%short_name = "SCN2" f3d(n)%units = " " f3d(n)%data => scn2 i_scn2 = n ; n = n+1 ! allocate(sco3(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O3 SLANT-LINE INTEGRAL')") f3d(n)%short_name = "SCO3" f3d(n)%units = " " f3d(n)%data => sco3 i_sco3 = n ; n = n+1 ! allocate(scco2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CO2 SLANT-LINE INTEGRAL')") f3d(n)%short_name = "SCCO2" f3d(n)%units = " " f3d(n)%data => scco2 i_scco2 = n ; n = n+1 ! allocate(scno(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('NO SLANT-LINE INTEGRAL')") f3d(n)%short_name = "SCNO" f3d(n)%units = " " f3d(n)%data => scno i_scno = n ; n = n+1 ! allocate(xiop2p(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('XIOP2P')") f3d(n)%short_name = "XIOP2P" f3d(n)%units = " " f3d(n)%data => xiop2p i_xiop2p = n ; n = n+1 ! allocate(xiop2d(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('XIOP2D')") f3d(n)%short_name = "XIOP2D" f3d(n)%units = " " f3d(n)%data => xiop2d i_xiop2d = n ; n = n+1 ! allocate(hplus(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('HPLUS')") f3d(n)%short_name = "HPLUS" f3d(n)%units = " " f3d(n)%data => hplus i_hplus = n ; n = n+1 ! allocate(nplus(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('N+ ION')") f3d(n)%short_name = "NPLUS" f3d(n)%units = " " f3d(n)%data => nplus i_nplus = n ; n = n+1 ! allocate(n2p(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('N2+ ION')") f3d(n)%short_name = "N2P" f3d(n)%units = " " f3d(n)%data => n2p i_n2p = n ; n = n+1 ! allocate(nop(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('NO+ ION')") f3d(n)%short_name = "NOP" f3d(n)%units = " " f3d(n)%data => nop i_nop = n ; n = n+1 ! allocate(lxx(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('LAMDA ION DRAG XX')") f3d(n)%short_name = "LXX" f3d(n)%units = " " f3d(n)%data => lxx i_lxx = n ; n = n+1 ! allocate(lyy(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('LAMDA ION DRAG YY')") f3d(n)%short_name = "LYY" f3d(n)%units = " " f3d(n)%data => lyy i_lyy = n ; n = n+1 ! allocate(lxy(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('LAMDA ION DRAG XY')") f3d(n)%short_name = "LXY" f3d(n)%units = " " f3d(n)%data => lxy i_lxy = n ; n = n+1 ! allocate(lyx(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('LAMDA ION DRAG YX')") f3d(n)%short_name = "LYX" f3d(n)%units = " " f3d(n)%data => lyx i_lyx = n ; n = n+1 ! allocate(qji_ti(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('ION JOULE HEATING FOR TI')") f3d(n)%short_name = "QJI_TI" f3d(n)%units = " " f3d(n)%data => qji_ti i_qji_ti = n ; n = n+1 ! allocate(qji_tn(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('ION JOULE HEATING FOR TN')") f3d(n)%short_name = "QJI_TN" f3d(n)%units = " " f3d(n)%data => qji_tn i_qji_tn = n ; n = n+1 ! allocate(cool_implicit(levd0:levd1,lond0:lond1,latd0:latd1), | stat=istat) write(f3d(n)%long_name,"('IMPLICIT COOLING TERM')") f3d(n)%short_name = "cool_implicit" f3d(n)%units = " " f3d(n)%data => cool_implicit i_cool_implicit = n ; n = n+1 ! allocate(cool_explicit(levd0:levd1,lond0:lond1,latd0:latd1), | stat=istat) write(f3d(n)%long_name,"('EXPLICIT COOLING TERM')") f3d(n)%short_name = "cool_explicit" f3d(n)%units = " " f3d(n)%data => cool_explicit i_cool_explicit = n ; n = n+1 ! allocate(hdt(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('TN HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdt" f3d(n)%units = " " f3d(n)%data => hdt i_hdt = n ; n = n+1 ! allocate(hdu(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('UN HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdu" f3d(n)%units = " " f3d(n)%data => hdu i_hdu = n ; n = n+1 ! allocate(hdv(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('VN HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdv" f3d(n)%units = " " f3d(n)%data => hdv i_hdv = n ; n = n+1 ! allocate(hdo2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O2 HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdo2" f3d(n)%units = " " f3d(n)%data => hdo2 i_hdo2 = n ; n = n+1 ! allocate(hdox(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('OX HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdox" f3d(n)%units = " " f3d(n)%data => hdox i_hdox = n ; n = n+1 ! ! Minor sp (time-gcm only): allocate(hdn4s(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('N4S HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdn4s" f3d(n)%units = " " f3d(n)%data => hdn4s i_hdn4s = n ; n = n+1 ! allocate(hdnoz(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('NOZ HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdnoz" f3d(n)%units = " " f3d(n)%data => hdnoz i_hdnoz = n ; n = n+1 ! allocate(hdco(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CO HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdco" f3d(n)%units = " " f3d(n)%data => hdco i_hdco = n ; n = n+1 ! allocate(hdco2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CO2 HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdco2" f3d(n)%units = " " f3d(n)%data => hdco2 i_hdco2 = n ; n = n+1 ! allocate(hdh2o(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('H2O HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdh2o" f3d(n)%units = " " f3d(n)%data => hdh2o i_hdh2o = n ; n = n+1 ! allocate(hdh2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('H2 HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdh2" f3d(n)%units = " " f3d(n)%data => hdh2 i_hdh2 = n ; n = n+1 ! allocate(hdhox(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('HOX HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdhox" f3d(n)%units = " " f3d(n)%data => hdhox i_hdhox = n ; n = n+1 ! allocate(hdch4(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CH4 HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdch4" f3d(n)%units = " " f3d(n)%data => hdch4 i_hdch4 = n ; n = n+1 ! allocate(hdar(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('AR HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdar" f3d(n)%units = " " f3d(n)%data => hdar i_hdar = n ; n = n+1 ! allocate(hdhe(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('HE HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdhe" f3d(n)%units = " " f3d(n)%data => hdhe i_hdhe = n ; n = n+1 ! allocate(hdnat(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('NAT HORIZONTAL DIFFUSION')") f3d(n)%short_name = "hdnat" f3d(n)%units = " " f3d(n)%data => hdnat i_hdnat = n ; n = n+1 ! allocate(ped(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('PEDERSEN CONDUCTIVITY')") f3d(n)%short_name = "ped" f3d(n)%units = " " f3d(n)%data => ped i_ped = n ; n = n+1 ! allocate(hall(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('HALL CONDUCTIVITY')") f3d(n)%short_name = "hall" f3d(n)%units = " " f3d(n)%data => hall i_hall = n ; n = n+1 ! allocate(nocool(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('NO COOLING')") f3d(n)%short_name = "nocool" f3d(n)%units = " " f3d(n)%data => nocool i_nocool = n ; n = n+1 ! ! n2o, cl, and clo are defined from solomon-garcia imports ! (see solgar.F and cmpsolgar.F) ! allocate(n2o(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('N2O (SOLOMON-GARCIA)')") f3d(n)%short_name = "n2o" f3d(n)%units = " " f3d(n)%data => n2o i_n2o = n ; n = n+1 ! allocate(cl(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CL (SOLOMON-GARCIA)')") f3d(n)%short_name = "cl" f3d(n)%units = " " f3d(n)%data => cl i_cl = n ; n = n+1 ! allocate(clo(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('CLO (SOLOMON-GARCIA)')") f3d(n)%short_name = "clo" f3d(n)%units = " " f3d(n)%data => clo i_clo = n ; n = n+1 ! allocate(o1d(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O(1DELTA)')") f3d(n)%short_name = "o1d" f3d(n)%units = " " ! MMR? f3d(n)%data => o1d i_o1d = n ; n = n+1 ! allocate(o21s(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('O2(1SIGMA)')") f3d(n)%short_name = "o21s" f3d(n)%units = " " ! MMR? f3d(n)%data => o21s i_o21s = n ; n = n+1 ! allocate(gwu(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('GRAVITY WAVE ZONAL TENDENCY')") f3d(n)%short_name = "gwu" f3d(n)%units = "CM/S" f3d(n)%data => gwu i_gwu = n ; n = n+1 ! allocate(gwv(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('GRAVITY WAVE MERIDIONAL TENDENCY')") f3d(n)%short_name = "gwv" f3d(n)%units = "CM/S" f3d(n)%data => gwv i_gwv = n ; n = n+1 ! allocate(gwt(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('GRAVITY WAVE TEMPERATURE TENDENCY')") f3d(n)%short_name = "gwt" f3d(n)%units = "CM/S" f3d(n)%data => gwt i_gwt = n ; n = n+1 ! allocate(h2o2(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('H2O2')") f3d(n)%short_name = "h2o2" f3d(n)%units = "MMR" f3d(n)%data => h2o2 i_h2o2 = n ; n = n+1 ! allocate(phoxic(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('PHOXIC')") f3d(n)%short_name = "phoxic" f3d(n)%units = " " f3d(n)%data => phoxic i_phoxic = n ; n = n+1 ! allocate(difkk(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('EDDY VISCOCITY')") f3d(n)%short_name = "difkk" f3d(n)%units = " " f3d(n)%data => difkk i_difkk = n ; n = n+1 ! allocate(difkt(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('EDDY THERMAL DIFFUSION')") f3d(n)%short_name = "difkt" f3d(n)%units = " " f3d(n)%data => difkt i_difkt = n ; n = n+1 ! allocate(difkv(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('EDDY DIFFUSION')") f3d(n)%short_name = "difkv" f3d(n)%units = " " f3d(n)%data => difkv i_difkv = n ; n = n+1 ! allocate(lam1(levd0:levd1,lond0:lond1,latd0:latd1),stat=istat) write(f3d(n)%long_name,"('LAMDA ION DRAG 1')") f3d(n)%short_name = "LAM1" f3d(n)%units = " " f3d(n)%data => lam1 i_lam1 = n ; n = n+1 ! ! Check number of fields allocated: if (n-1 /= nf3d) then write(6,"(/,'>>> init_3d: wrong number of 3-d fields?', | ' n-1=',i5,' nf3d=',i5)") n-1,nf3d endif ! ! f3d fields are diagnostics: f3d%prognostic = .false. f3d%magnetic = .false. f3d%magnetos = .false. f3d%task0_only = .false. f3d%vcoord = 'midpoints' do n=1,nf3d if (trim(f3d(n)%short_name)=='ZG') f3d(n)%vcoord='interfaces' enddo ! ! Data init: do i=1,nf3d f3d(i)%data = field_initval enddo ! ! Report to stdout: if (iprint > 0) then write(6,"(/,'There are ',i3,' 3-d fields:')") nf3d do n=1,nf3d call print_f3d(f3d(n),n) enddo write(6,"('Initialized f3d%data to ',e12.4)") field_initval endif end subroutine init_3d !----------------------------------------------------------------------- subroutine init_fsech ! ! Initialize secondary history field structures. This does not include ! data, which is allocated and defined by sub addfld (addfld.F). ! use input_module,only: secflds use hist_module,only: nfsech implicit none ! ! Local: integer :: i,iprog,ier ! ! External: integer,external :: strloc ! ! Secondary history fields: fsechist(:)%long_name = ' ' fsechist(:)%short_name= ' ' fsechist(:)%units = ' ' fsechist(:)%prognostic = .false. ! 5/9/06 btf: Set default task0_only=true to cover case where addfld is ! called by task0 only (e.g., from serial dynamo) fsechist(:)%task0_only = .false. ! fsechist(:)%task0_only = .true. fsechist(:)%mag = .false. fsechist(:)%geo = .false. do i=1,nfsech if (len_trim(secflds(i)) > 0) then iprog = strloc(f4d%short_name,nf4d,secflds(i)) ! ! Is a prognostic: define fsech(i) from f4d(iprog): ! if (iprog > 0) then fsechist(i)%prognostic = .true. fsechist(i)%short_name = f4d(iprog)%short_name fsechist(i)%long_name = f4d(iprog)%long_name fsechist(i)%units = f4d(iprog)%units ! ! Is a diagnostic: define fsech(i)%name from input field name. ! Set units and long_name blank (will be optionally defined in ! user called sub addfld) ! else ! is diagnostic fsechist(i)%long_name = secflds(i) fsechist(i)%short_name = secflds(i) fsechist(i)%units = ' ' fsechist(i)%long_name = ' ' endif endif enddo ! i=1,nfsech end subroutine init_fsech !----------------------------------------------------------------------- subroutine init_lbc(lon0,lon1,lat0,lat1) implicit none ! ! Arg: integer,intent(in) :: lon0,lon1,lat0,lat1 ! ! Local: integer :: istat ! ! Allocate t,u,v lbc (t,u,v at bottom interface level) ! These will be read from source history, then set in dt.F and duv.F. ! ! Subdomains: allocate(tlbc(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of tlbc') allocate(ulbc(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of ulbc') allocate(vlbc(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of vlbc') write(6,"('init_lbc: allocated subdomains tlbc, ulbc, vlbc')") allocate(tlbc_nm(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of tlbc_nm') allocate(ulbc_nm(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of ulbc_nm') allocate(vlbc_nm(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of vlbc_nm') write(6,"('init_lbc_nm: allocated subdomains tlbc_nm, ulbc_nm,', | ' vlbc_nm')") ! ! Global domain (subdomains will be gathered into these arrays ! for writing to the history, see sub mp_gather2root_lbc in mpi.F): allocate(tlbc_glb(nlonp4,nlat),stat=istat) if (istat /= 0) call shutdown('bad allocate of tlbc_glb') allocate(ulbc_glb(nlonp4,nlat),stat=istat) if (istat /= 0) call shutdown('bad allocate of ulbc_glb') allocate(vlbc_glb(nlonp4,nlat),stat=istat) if (istat /= 0) call shutdown('bad allocate of vlbc_glb') write(6,"('init_lbc: allocated globals tlbc_glb, ulbc_glb,', | ' vlbc_glb: nlonp4=',i3,' nlat=',i3)") nlonp4,nlat allocate(tlbc_nm_glb(nlonp4,nlat),stat=istat) if (istat /= 0) call shutdown('bad allocate of tlbc_nm_glb') allocate(ulbc_nm_glb(nlonp4,nlat),stat=istat) if (istat /= 0) call shutdown('bad allocate of ulbc_nm_glb') allocate(vlbc_nm_glb(nlonp4,nlat),stat=istat) if (istat /= 0) call shutdown('bad allocate of vlbc_nm_glb') write(6,"('init_lbc: allocated globals tlbc_nm_glb, ulbc_nm_glb,', | ' vlbc_nm_glb: nlonp4=',i3,' nlat=',i3)") nlonp4,nlat tlbc = 0. ; ulbc = 0. ; vlbc = 0. tlbc_nm = 0. ; ulbc_nm = 0. ; vlbc_nm = 0. tlbc_glb = 0. ; vlbc_glb = 0. ; vlbc_glb = 0. tlbc_nm_glb = 0. ; vlbc_nm_glb = 0. ; vlbc_nm_glb = 0. ! ! Background and perturbed lbc (see lbc.F and uvbnd.F) allocate(t_bgrd(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of t_bgrd') allocate(z_bgrd(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of z_bgrd') allocate(u_bgrd(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of u_bgrd') allocate(v_bgrd(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of v_bgrd') t_bgrd = 0. ; z_bgrd = 0. ; u_bgrd = 0. ; v_bgrd = 0. allocate(u_prt(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of u_prt') allocate(v_prt(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of v_prt') allocate(u_prt_tn(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of u_prt_tn') allocate(v_prt_tn(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of v_prt_tn') allocate(u_prt_tnm(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of u_prt_tnm') allocate(v_prt_tnm(lond0:lond1,latd0:latd1),stat=istat) if (istat /= 0) call shutdown('bad allocate of v_prt_tnm') u_prt = 0. ; v_prt = 0. ; u_prt_tn = 0. ; v_prt_tn = 0. u_prt_tnm = 0. ; v_prt_tnm = 0. ! end subroutine init_lbc !----------------------------------------------------------------------- subroutine print_f4d(f,n) implicit none type(fields_4d),intent(in) :: f integer,intent(in) :: n ! if (.not.associated(f%data)) then write(6,"('Field ',i3,': ',a,' Short name: ',a,' Units: ',a)") | n,f%long_name(1:40),f%short_name(1:8),f%units(1:8) else write(6,"('Field ',i3,': ',a,' Short name: ',a,' Units: ',a, | ' data min,max=',2e12.4)") | n,f%long_name(1:40),f%short_name(1:8),f%units(1:8), | minval(f%data),maxval(f%data) endif end subroutine print_f4d !----------------------------------------------------------------------- subroutine print_f3d(f,n) implicit none type(fields_3d),intent(in) :: f integer,intent(in) :: n ! write(6,"('Field ',i3,': ',a,' Short name: ',a,' Units: ',a)") | n,f%long_name(1:40),f%short_name(1:8),f%units(1:8) end subroutine print_f3d !----------------------------------------------------------------------- subroutine print_f2d(f,n) implicit none type(fields_2d),intent(in) :: f integer,intent(in) :: n ! write(6,"('Field ',i3,': ',a,' Short name: ',a,' Units: ',a)") | n,f%long_name(1:40),f%short_name(1:8),f%units(1:8) end subroutine print_f2d !----------------------------------------------------------------------- subroutine fprint(msg,k0,k1,i0,i1,j0,j1,itx) implicit none character(len=*),intent(in) :: msg integer,intent(in) :: k0,k1,i0,i1,j0,j1,itx ! ! integer,parameter :: nf = 2 integer,parameter :: nf = 1 integer :: i,ii(nf) ! 4d allocation: ! allocate(tn(levd0:levd1,lond0:lond1,latd0:latd1,2),stat=istat) return ! temporarilly turn this off w/o having to remove the calls ! ii = (/i_co,i_co2/) ii = (/i_un/) do i=1,nf if (.not.associated(f4d(ii(i))%data)) then write(6,"('Field ',a,' (data not allocated) msg=''',a,'''')") | f4d(ii(i))%short_name,trim(msg) else write(6,"('Field ',a,' data min,max=',2e12.4, | ' msg=''',a,'''')") | f4d(ii(i))%short_name, | minval(f4d(ii(i))%data(k0:k1,i0:i1,j0:j1,itx)), | maxval(f4d(ii(i))%data(k0:k1,i0:i1,j0:j1,itx)), | trim(msg) endif enddo end subroutine fprint !----------------------------------------------------------------------- end module fields_module