! Copyright 2014 College of William and Mary ! ! Licensed under the Apache License, Version 2.0 (the "License"); ! you may not use this file except in compliance with the License. ! You may obtain a copy of the License at ! ! http://www.apache.org/licenses/LICENSE-2.0 ! ! Unless required by applicable law or agreed to in writing, software ! distributed under the License is distributed on an "AS IS" BASIS, ! WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. ! See the License for the specific language governing permissions and ! limitations under the License. !=============================================================================== !=============================================================================== ! SCHISM GLOBAL DATA MODULE !=============================================================================== !=============================================================================== module schism_glbl implicit none public !Default scope is public !#ifdef USE_SINGLE ! integer,parameter :: rkind = 4 !#else integer,parameter :: rkind = 8 ! Default real datatype !iwp*: working precision of single/double real option for some code integer,parameter :: iwp= 8 ! for ICM !#endif #ifdef DOUBLE_REAL_OUT integer,parameter::out_rkind=8 #else integer,parameter ::out_rkind=4 !Default output real type #endif ! Some constants integer,parameter :: nthfiles=5 !# of type I (ASCII) .th files (for dimensioning) integer,parameter :: nthfiles2=5 !# of *3D.th files (for dimensioning) integer,parameter :: nthfiles3=3 !# of source/sink .th (for dimensioning) integer,parameter :: max_ncoutvar=500 !max # of vars in global nc outputs real(rkind),parameter :: small1=real(1.d-6,rkind) !small non-negative number real(rkind),parameter :: small2=small1*100._rkind !slightly larger number real(rkind),parameter :: pi=3.141592653589793_rkind real(rkind),parameter :: grav=9.81_rkind real(rkind),parameter :: omega_e=real(7.292d-5,rkind) !angular freq. of earth rotation !# of threads in openMP. Note that this is the # at the start of run (thus the 'default'), but !you can change # of threads in some loops (after that, it should revert to default) integer,save :: nthreads !In/out dirs character(len=1000),save :: in_dir,out_dir integer,save :: len_in_dir,len_out_dir ! For timing integer,parameter :: mxtimer=20 ! Max number of wallclock timers real(rkind),save :: wtimer(0:mxtimer,2) ! Array of wallclock timers ! (:,1)=execution time ! (:,2)=communication time ! timer_ns: used in target routines for profiling. Increase its dimension and use for your own purpose. ! 1: upwind/TVD transport real(rkind),save :: timer_ns(11) ! For debugging character(len=72) :: fdb ! Name of debugging file integer :: lfdb ! Length of debugging file name ! Error message string character(len=2000),save :: errmsg ! Constants used in UB closure character(len=2),save :: mid,stab real(rkind),save :: ubd0,ubd1,ubd2,ubd3,ubd4,ubd5, & ubs0,ubs1,ubs2,ubs4,ubs5,ubs6, & a2_cm03,schk,schpsi integer,parameter :: natrm=12 !# of _available_ tracer models at the moment (including T,S) integer,parameter :: mntracers=30 !max # of tracers, used only for dimensioning btrack arrays. Must >=ntracers !Parameters from param.nml integer,save :: ipre,ipre2,indvel,imm,ihot,ics,iwbl,iharind,nws,impose_net_flux,iwindoff, & &ibc,ibdef,ihorcon,nstep_wwm,icou_elfe_wwm, & &fwvor_advxy_stokes,fwvor_advz_stokes,fwvor_gradpress,fwvor_breaking, & &fwvor_streaming,cur_wwm,wafo_obcramp,iwind_form,irec_nu,itur,ihhat,inu_elev, & &inu_uv,ibcc_mean,iflux,iout_sta,nspool_sta,nhot,nhot_write, & &moitn0,mxitn0,nchi,ibtrack_test,nramp_elev,islip,ibtp,inunfl,shorewafo, & &inv_atm_bnd,ieos_type,ieos_pres,iupwind_mom,inter_mom,ishapiro,isav, & &nstep_ice,niter_shap,iunder_deep,flag_fib,ielm_transport,max_subcyc, & &itransport_only,meth_sink,iloadtide,nc_out,nu_sum_mult,iflux_out_format integer,save :: ntrs(natrm),nnu_pts(natrm),mnu_pts,lev_tr_source(natrm) integer,save,dimension(:),allocatable :: iof_hydro,iof_wwm,iof_gen,iof_age,iof_sed,iof_eco, & &iof_icm,iof_icm_core,iof_icm_silica,iof_icm_zb,iof_icm_ph,iof_icm_sav,iof_icm_veg, & &iof_cos,iof_fib,iof_sed2d,iof_ice,iof_ana,iof_marsh,iof_dvd,iadjust_mass_consv, & &lev_tr_source2(:) real(rkind),save :: dt,h0,drampbc,drampwind,drampwafo,dramp,dramp_ss,wtiminc,npstime,npstiminc, & &surf_time1,surf_time2,time_nu,step_nu,time_nu_tr,step_nu_tr,dzb_min,vdmax_pp1, & &vdmin_pp1,tdmin_pp1,vdmax_pp2,vdmin_pp2,tdmin_pp2, & &h1_pp,h2_pp,dtb_min,dtb_max,thetai,theta2,rtol0, & &vnh1,vnh2,vnf1,vnf2,rnday,btrack_nudge,hmin_man, & &prmsl_ref,hmin_radstress,eos_a,eos_b,eps1_tvd_imp,eps2_tvd_imp, & &xlsc0,rearth_pole,rearth_eq,hvis_coef0,disch_coef(10),hw_depth,hw_ratio, & &slr_rate,rho0,shw,gen_wsett,turbinj,turbinjds,alphaw,h1_bcc,h2_bcc,vclose_surf_frac, & &hmin_airsea_ex,hmin_salt_ex,shapiro0,loadtide_coef,h_massconsv,rinflation_icm ! Misc. variables shared between routines integer,save :: nz_r,ieqstate,kr_co, & &ihconsv,isconsv,i_hmin_airsea_ex,i_hmin_salt_ex,ihdif,ntracers, & &ihydraulics,irouse_test,iwbl_itmax,nettype,nfltype, & &ntetype,nsatype,ntrtype1(natrm),nettype2,nnode_et,nfltype2,nnode_fl, & &ntetype2,nsatype2,nnode_tr2(natrm),inu_tr(natrm), & &nvar_sta,nout_sta,ntip,nbfr,itr_met,if_source,mass_source,nsources,nsinks, & &max_flreg,irange_tr(2,natrm),nea_wwm,mnei_wwm,ne_wwm,neg_wwm, & &max_iadjust_mass_consv,nsteps_from_cold real(rkind),save :: q2min,tempmin,tempmax,saltmin,saltmax, & &vis_coe1,vis_coe2,h_bcc1,velmin_btrack,h_tvd,rmaxvel1,rmaxvel2, & &difnum_max_l2,wtime1,wtime2,cmiu0, & &cpsi2,rpub,rmub,rnub,cpsi1,psimin,eps_min,tip_dp,sav_di0,sav_h0,sav_nv0, & &dtb_min_transport,bounds_lon(2) ! logical,save :: lm2d !2D or 3D model logical,save :: lhas_quad=.false. !existence of quads logical,save :: lflbc !flag to indicate existence of ifltype/=0 logical,save :: lreadll=.false. !if true, xlon and ylat are already read in and can be used by any routine/modules integer,save :: i2tier_case !how 2-tier ghosts are constructed ! Variables for global output files integer, parameter :: nbyte=4 !# bytes for output record size ! integer, parameter :: mnout=200 !max. # of output files ! integer, parameter :: mirec=1109000000 !max. record # to prevent output ~> 4GB ! character(len=11), parameter :: fileopenformat='unformatted' ! integer,save :: iwrite ! character(len=48),save :: start_time,version,data_format='DataFormat v5.0' integer,save :: start_year = -9999 integer,save :: start_month = -9999 integer,save :: start_day = -9999 real(rkind),save :: start_hour = -9999._rkind real(rkind),save :: utc_start = -9999._rkind character(len=12),save :: ifile_char ! character(len=48),save,dimension(mnout) :: outfile !,variable_nm,variable_dim integer,save :: ihfskip,nrec,nspool,ifile,ifile_len, & &noutput,it_main,iths_main,id_out_var(2000),ncount_2dnode, & &ncount_2delem,ncount_2dside,ncount_3dnode,ncount_3delem,ncount_3dside,nsend_varout integer,save,allocatable :: srqst7(:) real(rkind),save :: time_stamp !simulation time in sec !Send var buffers real(4),save,allocatable :: varout_3dnode(:,:,:),varout_3delem(:,:,:),varout_3dside(:,:,:) real(4),save,allocatable :: varout_2dnode(:,:),varout_2delem(:,:),varout_2dside(:,:) character(len=48),save,allocatable :: outfile_ns(:) character(len=48),save :: a_48 character(len=16),save :: a_16 character(len= 8),save :: a_8 character(len= 4),save :: a_4 integer,save :: ncid_nu(natrm),ncid_tr3D(natrm),ncid_elev2D,ncid_uv3D, & &istack0_schout,ncid_source ! ADT for global-to-local linked-lists type :: llist_type integer :: rank ! Processor rank assignment integer :: id=0 ! Local index on processor "rank" type(llist_type),pointer :: next=>null() ! Next entry in linked-list end type llist_type ! ADT for inter-subdomain backtracking type :: bt_type integer :: rank ! Originating processor rank integer :: l0 ! Originating node or side or element index (1<=l0<=3) integer :: i0gb ! Originating node or side or element _global_ index integer :: isbndy ! Flag for originating node or side on the boundary (for Kriging) integer :: j0 ! Originating vertical level integer :: adv ! Original advection flag (0-2) ! integer :: ndt ! Number of backtracking sub-steps fro Euler tracking integer :: iegb ! Global index of current encompassing element integer :: jvrt ! Current vertical level real(rkind) :: dtbk ! target time step real(rkind) :: vis ! vis_coe real(rkind) :: rt ! time remaining (total inside dt) real(rkind) :: rt2 ! time remaining from left-over from previous subdomain real(rkind) :: ut,vt,wt ! Current backtracking sub-step velocity real(rkind) :: xt,yt,zt ! Current backtracking sub-step point real(rkind) :: sclr(4+mntracers) ! Backtracked values for tracers etc real(rkind) :: gcor0(3) ! global coord. of the starting pt (for ics=2) real(rkind) :: frame0(3,3) ! frame tensor at starting pt (for ics=2) end type bt_type integer,save :: bt_mpitype ! MPI datatype for inter-subdomain backtracking real(rkind),save :: s1_mxnbt ! scale used in mxnbt for btlist() real(rkind),save :: s2_mxnbt ! scale used in the routine inter_btrack() integer,save :: mxnbt ! Dimension of btlist() ! Vertical layer data !Vertical coord. types. 2: SZ; -1: Z (only significant in sflux_subs.F90); 1: !localized sigma integer,save :: ivcor integer,save :: nvrt ! Number of vertical layers integer,save :: kz ! Number of Z levels integer,save :: nsig ! Number of S levels real(rkind),save :: h_s,h_c,theta_b,theta_f,s_con1 !constants used in vgrid.in real(rkind),save,allocatable :: ztot(:) ! Z coord. of Z levels (local frame) real(rkind),save,allocatable :: sigma(:) ! sigma coordinates real(rkind),save,allocatable :: cs(:) ! function in S-coordinate real(rkind),save,allocatable :: dcs(:) ! derivative of cs() real(rkind),save,allocatable :: sigma_lcl(:,:) !localized sigma; ivcor=1 ! Element geometry data integer,save :: ne_global ! Global number of elements integer,save :: ne ! Local number of resident elements integer,save :: neg ! Local number of ghost elements integer,save :: nea ! Local number of elements in augmented subdomain (ne+neg) integer,save :: neg2 ! Local number of 2-tier ghost elements integer,save :: nea2 ! Local number of elements in 2-tier augmented subdomain (ne+neg+neg2) integer,save,allocatable :: ielg(:) ! Local-to-global element index table (augmented) type(llist_type),save,pointer :: iegl(:) ! Global-to-local element index table (augmented) integer,save,allocatable :: iegrpv(:) ! Global element to resident processor vector integer,save :: nx(3,2) ! Cyclic index of nodes in an element (kept only for modules) !nxq(i,j,k): i: offset; j: local index; k: elem. type (3 or 4) integer,save :: nxq(3,4,4) ! Cyclic index of nodes in an element (tri/quads) integer,save :: ixi_n(4),iet_n(4) !local coord. of 4 vertices of quads integer,save,allocatable :: i34(:) !elem. type (3 or 4) integer,save,allocatable :: elnode(:,:) ! Element-node tables integer,save,allocatable :: elnode_wwm(:,:) ! Element-node tables after splitting quads (for WWM) integer,save,allocatable :: iself(:,:) ! Index of node in element-node table integer,save,allocatable :: ic3(:,:) ! Element-side-element tables integer,save,allocatable :: elside(:,:) ! Element-side tables real(rkind),save,allocatable :: ssign(:,:) ! Sign associated with each side of an element real(rkind),save,allocatable :: area(:) ! Element areas real(rkind),save,allocatable :: radiel(:) ! Element equivalent radii ! Cartesian coordinates of element centers; see comments for xnd real(rkind),save,allocatable :: xctr(:),yctr(:),zctr(:) real(rkind),save,allocatable :: xlon_el(:),ylat_el(:) ! Element center lat/lon coordinates in _degrees_ real(rkind),save,allocatable :: dpe(:) ! Depth at element (min of all nodes) integer,save,allocatable :: kbe(:) ! Element bottom vertical indices integer,save,allocatable :: idry_e(:) ! wet/dry flag integer,save,allocatable :: idry_e_2t(:) ! wet/dry flag including 2-tier ghost zone integer,save,allocatable :: interpol(:) ! interpolation mode ! integer,save,allocatable :: lqk(:) ! interpolation for S,T in btrack integer,save,allocatable :: ie_kr(:) ! used in Kriging integer,save,allocatable :: krvel(:) ! used in Kriging real(rkind),save,allocatable :: ze(:,:) ! z-coord (local frame - vertical up) integer,save,allocatable :: itvd_e(:) !for TVD transport !Derivatives of shape function in an element !For ics=1, the global coordinates are used !For ics=2, element frame is used !dldxy(j,k,ie)=dL_{j}/dx_{k}, where j=1:i34(ie) (shape function index), k=1:2 for x|y; !ie is the local element index. For quads, the derivative is evaluated at !centroid real(rkind),save,allocatable :: dldxy(:,:,:) !Transformation tensor for element (ll) frame: eframe(i,j,ie) for ics=2 !where j is the axis id, i is the component id, ie is the local element id (aug.) !Undefined for ics=1 real(rkind),save,allocatable :: eframe(:,:,:) !x,y coordinates of each element node/side in the _element_ frame !xel(4,nea), xs_el(4,nea) real(rkind),save,allocatable :: xel(:,:),yel(:,:),xs_el(:,:),ys_el(:,:) !shape_c2(4,2,nea)- for quads only. The shape function values for the mid-pts of 2 diagnonals inside ! the quad formed by 4 sidecenters real(rkind),save,allocatable :: shape_c2(:,:,:) integer,save,allocatable :: iflux_e(:) !for computing fluxes integer,save,allocatable :: ielg2(:) ! Local-to-global element index table (2-tier augmented) integer,save,allocatable :: iegl2(:,:) ! Global-to-local element index table (2-tier augmented) ! Node geometry data integer,save :: mnei ! Max number of neighboring elements surrounding a node integer,save :: mnei_p ! Max number of neighboring nodes surrounding a node integer,save :: mnei_kr ! Max # of Kriging nodes integer,save :: np_global ! Global number of nodes integer,save :: np ! Local number of resident nodes integer,save :: npg ! Local number of ghost nodes integer,save :: npa ! Local number of nodes in augmented subdomain (np+npg) integer,save :: npg2 ! Local number of 2-tier ghost nodes integer,save :: npa2 ! Local number of nodes in 2-tier augmented subdomain (np+npg+npg2) integer,save,allocatable :: iplg(:) ! Local-to-global node index table (augmented) type(llist_type),save,pointer :: ipgl(:) ! Global-to-local node index table (augmented) !Node cartesian coordinates. They mean different things for ics=1 (plane projection) or ics=2 (sphere); !for ics=1, znd=0, and xnd,ynd are the Cartesian coord. in the projection plane; !for ics=2, the triplet are the coordinate in a global frame with origin at center of earth, !x-axis point to prime meridian, z-axis to the north pole real(rkind),save,allocatable :: xnd(:),ynd(:),znd(:) ! Node cartesian coordinates real(rkind),save,allocatable :: xlon(:),ylat(:) ! Node lat/lon coordinates in _radians_ real(rkind),save,allocatable :: xlon_gb(:),ylat_gb(:) ! global lat/lon coordinates in _degrees_ (for some options) ! Modified node lat/lon coordinates in _radians_ for ice model (via FESOM2). ! The nodes are at same physical location as (xlon, ylat) ! real(rkind),save,allocatable :: xlon2(:),ylat2(:) real(rkind),save,allocatable :: dp(:),dp00(:) ! Node depths ! integer,save,allocatable :: ibad(:) ! Reliable bndry elevation flag ! integer,save,allocatable :: nnegb(:),inegb(:,:) ! Global node-element tables integer,save,allocatable :: nne(:),indel(:,:) ! Node-element tables integer,save,allocatable :: nne_wwm(:) ! Node-element tables for WWM (temp use) integer,save,allocatable :: nnp(:),indnd(:,:) ! Node-node tables (exlcuding the node itself) integer,save,allocatable :: isbnd(:,:) ! local node to _global_ open bndry segment flags integer,save,allocatable :: ibnd_ext_int(:) ! interior (-1) /exterior (1) bnd node flag for an aug. node (0: not on bnd) ! real(rkind),save,allocatable :: edge_angle(:,:) !angles (orientation) at a bnd node of 2 adjacent sides ! integer,save,allocatable :: isbnd_global(:) ! Node to open bndry segment flags (global) integer,save,allocatable :: kbp(:),kbp00(:),kbp_e(:) ! Node surface & bottom vertical indices integer,save,allocatable :: idry(:) ! wet/dry flag ! integer,save,allocatable :: iback(:) ! back-up flag for abnormal cases in S-coord. real(rkind),save,allocatable :: hmod(:) ! constrained depth real(rkind),save,allocatable :: znl(:,:) ! z-coord in local Z-axis (vertical up) ! Transformation tensor for node (ll) frame: pframe(i,j,ip) for ics=2. ! where j is the axis id, i is the component id, ip is the local node id (aug.) ! For ics=1, this is not used real(rkind),save,allocatable :: pframe(:,:,:) integer,save,allocatable :: iplg2(:) ! Local-to-global node index table (2-tier augmented) integer,save,allocatable :: ipgl2(:,:) ! Global-to-local node index table (2-tier augmented) ! Side geometry data integer,save :: ns_global ! Global number of sides integer,save :: ns ! Local number of resident sides integer,save :: nsg ! Local number of ghost sides integer,save :: nsa ! Local number of sides in augmented subdomain (ns+nsg) integer,save :: nsg2 ! Local number of 2-tier ghost sides integer,save :: nsa2 ! Local number of sides in 2-tier augmented subdomain (ns+nsg+nsg2) integer,save,allocatable :: islg(:) ! Local-to-global side index table (augmented) type(llist_type),save,pointer :: isgl(:) ! Global-to-local side index table (augmented) integer,save,allocatable :: isdel(:,:) ! Side-element tables integer,save,allocatable :: isidenode(:,:) ! Side-node tables !Cartesian coordinates of side centers; see the comments for xnd real(rkind),save,allocatable :: xcj(:),ycj(:),zcj(:) real(rkind),save,allocatable :: dps(:) ! Depth at side centers real(rkind),save,allocatable :: distj(:) ! Side lengths ! Distance between adjacent elements of an internal side; used only in horizontal diffusion real(rkind),save,allocatable :: delj(:) integer,save,allocatable :: isbs(:) ! local side to _global_ open bndry segment mapping ! integer,save,allocatable :: isbs_global(:) ! Side to open bndry segment mapping (global) integer,save,allocatable :: kbs(:) ! Side bottom vertical indices integer,save,allocatable :: idry_s(:) ! wet/dry flag integer,save,allocatable :: isidenei2(:,:) !side neighborhood real(rkind),save,allocatable :: zs(:,:) ! z-coord. (local frame - vertical up) !Not used: Transformation tensor for side frame: sframe(i,j=1,isd) ! where j=1 is the axis id, i is the component id, isd is the local side id (aug.) ! In all cases, x-axis is from adjacent elem 1 to 2 (i.e. local normal dir). ! We only need j=1 ! real(rkind),save,allocatable :: sframe(:,:,:) !new37: sframe2 !2nd version of side frame, aligning with local lon/lat frame. Only defined if ics=2 real(rkind),save,allocatable :: sframe2(:,:,:) !cos/sin of side normals (for ics=2: these are done inside local lon/lat frame) real(rkind),save,allocatable :: snx(:),sny(:) real(rkind),save,allocatable :: isblock_sd(:,:) integer,save,allocatable :: islg2(:) ! Local-to-global side index table (2-tier augmented) integer,save,allocatable :: isgl2(:,:) ! Global-to-local side index table (2-tier augmented) ! Open boundary segment data integer,save :: nope_global ! Global number of open bndry segments integer,save :: neta_global ! Global number of open bndry nodes integer,save :: nope ! Local number of local open bndry segments integer,save :: neta ! Local number of local open bndry nodes integer,save :: mnond ! Max # nodes per open bndry segment integer,save :: mnond_global ! Max # nodes per open bndry segment (global) integer,save,allocatable :: iopelg(:) ! Local-to-global open bndry segment table integer,save,allocatable :: iopegl(:,:) ! Global-to-Local open bndry segment table integer,save,allocatable :: nond(:) ! Number of nodes in each open bndry segment integer,save,allocatable :: iond(:,:) ! Node list for each open bndry segment integer,save,allocatable :: nond_global(:) ! Number of nodes in each open global bndry segment integer,save,allocatable :: iond_global(:,:) ! Node list for each open bndry segment (global) real(rkind),save,allocatable :: cwidth(:) ! length of each global open bnd segment for imposing discharge real(rkind),save,allocatable :: th_dt(:,:) !time step for .th (ascii; each tracer has its own dt) real(rkind),save,allocatable :: th_time(:,:,:) !2 time stamps for reading .th (ascii) real(rkind),save :: th_dt2(nthfiles2) !time step for .th (binary) real(rkind),save :: th_time2(2,nthfiles2) !2 time stamps for reading .th (binary) integer,save :: irec_th(nthfiles2) !binary record # for reading .th (binary) real(rkind),save :: th_dt3(nthfiles3) !time step for source/sink .th (ascii) real(rkind),save :: th_time3(2,nthfiles3) !2 time stamps for reading source/sink .th (ascii) real(rkind),save,allocatable :: trth(:,:,:,:) !time series of b.c. for T,S, tracers logical,save,allocatable :: lelbc(:) integer,save,allocatable :: iettype(:),ifltype(:),itrtype(:,:), & &jspc(:) integer,save,allocatable :: isblock_nd(:,:),isblock_el(:),iq_block_lcl(:),iq_block(:) real(rkind),save,allocatable :: trobc(:,:),vobc1(:),vobc2(:),tamp(:), & &tnf(:),tfreq(:),tear(:),amig(:),ff(:),face(:),rloadtide(:,:,:), & &emo(:,:,:),efa(:,:,:),umo(:,:,:),ufa(:,:,:), & &vmo(:,:,:),vfa(:,:,:),eth(:,:), & &qthcon(:),uthnd(:,:,:),vthnd(:,:,:), & &ath(:,:,:,:),carea(:),clen(:),eta_mean(:),q_block(:),vnth_block(:,:), & &dir_block(:,:),q_block_lcl(:) real(4),save,allocatable :: ath2(:,:,:,:,:) !used to read *.nc for b.c. time series real(4),save,allocatable :: ath3(:,:,:,:) !used to read source/sink inputs ! Land boundary segment data integer,save :: nland_global ! Global number of land bndry segments integer,save :: nvel_global ! Global number of land bndry nodes integer,save :: nland ! Local number of local land bndry segments integer,save :: nvel ! Local number of local land bndry nodes integer,save :: mnlnd ! Max # nodes per land bndry segment integer,save :: mnlnd_global ! Max # nodes per land bndry segment (global) integer,save,allocatable :: nlnd_global(:) ! Number of nodes in each land bndry segment (global) integer,save,allocatable :: ilnd_global(:,:) ! Node list for each land bndry segment (global) integer,save,allocatable :: nlnd(:) ! Number of nodes in each land bndry segment integer,save,allocatable :: ilnd(:,:) ! Node list for each land bndry segment ! Dynamic quantities integer,save,allocatable :: ieg_source(:) !global elem. indices for volume/mass sources integer,save,allocatable :: ieg_sink(:) !global elem. indices for volume/mass sinks !tracer concentration @ prism center; used as temp. storage. tr_el(ntracers,nvrt,nea2) but last index usually !is valid up to nea only except for TVD real(rkind),save,allocatable,target :: tr_el(:,:,:) real(rkind),save,allocatable,target :: tr_nd(:,:,:) !tracer concentration @ node and whole levels !Store T,S offline mode ts_offline(4,nvrt,nea), where [1 3] are for T, [2 4] !for S real(rkind),save,allocatable :: ts_offline(:,:,:) real(rkind),save,allocatable :: bdy_frc(:,:,:) !body force at prism center Q_{i,k} real(rkind),save,allocatable :: flx_sf(:,:) !surface b.c. \kappa*dC/dz = flx_sf (at element center) real(rkind),save,allocatable :: flx_bt(:,:) !bottom b.c. real(rkind),save,allocatable :: hdif(:,:) !horizontal diffusivity real(rkind),save,allocatable :: tr_nd0(:,:,:) ! Initial tracer conc. at nodes real(rkind),save,allocatable :: rkai_num(:,:,:) !DVD (numerical mixing) [C^2]/sec real(rkind),save,allocatable :: eta1(:) ! Elevation at nodes at previous timestep real(rkind),save,allocatable,target :: eta2(:) ! Elevation at nodes at current timestep !cumsum of elev from 1st _cold_ start (regardles of ihot) until the cum time step (nsteps_from_cold) !exceeds certain amount. This is used to compute anomaly real(rkind),save,allocatable :: cumsum_eta(:) !Vertical velocity at element centers & whole levels, along local vertical direction (element frame) real(rkind),save,allocatable :: we(:,:) !Vertical velocity at element centers & whole levels, calculated using F.V.M. For hydrostatic !model, this is the same as we(); for non-hydrostatic model, this is only used in upwind transport ! real(rkind),save,allocatable :: we_fv(:,:) real(rkind),save,allocatable :: flux_adv_vface(:,:,:) !unmodified vertical fluxes (positive upward) real(rkind),save,allocatable :: total_mass_error(:) !(ntracers) Total mass error after advection step for mass correction !x & y-component of velocity at side centers & whole levels !For ics=1, these are defined in the global frame !For ics=2, these are defined in the ll frame (sframe2) real(rkind),save,allocatable :: su2(:,:),sv2(:,:) !velocity at nodes in an element. In ll if ics=2 !ufg(1:4,1:nvrt,1:nea) !real(rkind),save,allocatable :: ufg(:,:,:),vfg(:,:,:) real(rkind),save,allocatable :: prho(:,:) ! Density at whole levels and either nodes or elements real(rkind),save,allocatable :: q2(:,:) ! Turbulent kinetic energy at sides & half levels real(rkind),save,allocatable :: xl(:,:) ! Turbulent mixing length at sides & half levels real(rkind),save,allocatable :: xlmin2(:) ! min. turbulent mixing length fro non-surface layers !Velocity at nodes & whole levels at current timestep. If ics=2, it's in ll frame real(rkind),save,allocatable, target :: uu2(:,:),vv2(:,:),ww2(:,:) real(rkind),save,allocatable :: bdef(:) !bottom deformation real(rkind),save,allocatable :: bdef1(:) !bottom deformation real(rkind),save,allocatable :: bdef2(:) !bottom deformation real(rkind),save,allocatable :: dfh(:,:) !diffusivity real(rkind),save,allocatable :: dfv(:,:) !viscosity integer,save,allocatable :: itier_nd(:,:) !multi-tier neighborhood; used in Kriging real(rkind),save,allocatable :: akrmat_nd(:,:,:) ! Kriging matrix real(rkind),save,allocatable :: albedo(:) ! albedo(npa) real(rkind),save,allocatable :: z_r(:) ! z-cor. used in ts.ic real(rkind),save,allocatable :: tem1(:) ! T profile in ts.ic real(rkind),save,allocatable :: sal1(:) ! S profile in ts.ic real(rkind),save,allocatable :: cspline_ypp(:,:) ! 2nd derivative in cubic spline for tem1,sal1 real(rkind),save,allocatable :: rho_mean(:,:) ! mean density real(rkind),save,allocatable :: Cdp(:) ! drag at node real(rkind),save,allocatable :: rmanning(:) ! Manning's n at node ! real(rkind),save,allocatable :: shapiro_min(:) !min of Shapiro filter strength (used with some ishapiro options) real(rkind),save,allocatable,target :: windx(:),windy(:) !wind vector real(rkind),save,allocatable,target :: prec_rain(:),prec_snow(:) !precipitation of rain and snow real(rkind),save,allocatable,target :: sdbt(:,:,:),shapiro(:), & &windx1(:),windy1(:),windx2(:),windy2(:), & &surf_t1(:),surf_t2(:),surf_t(:), & !YC !WARNING: airt[12] are in C not K. The !original air T in sflux_air*.nc is in K but !get_wind() converts it to C &tau(:,:),tau_bot_node(:,:),windfactor(:),pr1(:),airt1(:), & &shum1(:),pr2(:),airt2(:),shum2(:),pr(:), & &sflux(:),srad(:),tauxz(:),tauyz(:),fluxsu(:), & &fluxlu(:),hradu(:),hradd(:),cori(:), & !chi(:) &Cd(:),rough_p(:),erho(:,:),hvis_coef(:,:), & &sparsem(:,:), & &tr_nudge(:,:),fun_lat(:,:), & &elev_nudge(:),uv_nudge(:),fluxprc(:),fluxevp(:), & &dav(:,:),elevmax(:),dav_max(:,:),dav_maxmag(:), & &etaic(:),diffmax(:),diffmin(:),dfq1(:,:),dfq2(:,:) !(2,npa). ocean-ice stress (junk if no ice) [m^2/s/s] real(rkind),save,allocatable :: tau_oi(:,:) !(npa). freshwater flux due to ice melting [kg/s/m/m]. >0: precip; <0: evap real(rkind),save,allocatable :: fresh_wa_flux(:) !(npa). net heat flux into the ocean surface [W/m/m]. >0: warm the ocean real(rkind),save,allocatable :: net_heat_flux(:) real(rkind),save,allocatable :: wind_rotate_angle(:) !in radians !(npa). shortwave radiation through ice in ice model [W/m/m]. real(rkind),save,allocatable :: srad_th_ice(:) !(npa). evap water flux in ice model [kg/s/m/m]. real(rkind),save,allocatable :: ice_evap(:) real(rkind),save,allocatable :: srad_o(:) logical,save,allocatable :: lhas_ice(:) logical,save :: lice_free_gb real(4),save,dimension(:,:,:),allocatable :: trnd_nu1,trnd_nu2,trnd_nu integer,save,allocatable :: iadv(:),iwater_type(:) !weno> integer,save,allocatable :: isbe(:) !(ne): bnd seg flags, isbe(ie)=1 if any node of element ie lies on bnd; isbe(ie)=0 otherwise logical,save,allocatable :: is_inter(:) !identifier of interface sides (between two ranks), for debugging only integer,save,allocatable :: iside_table(:) !a record of all interface sides within the current rank integer, save :: ip_weno !order of the polynomials used for weno stencils, see param.in.sample real(rkind),save :: courant_weno !Courant number for weno transport real(rkind),save :: epsilon1 !coefficient for 2nd order weno smoother real(rkind),save :: epsilon2 !1st coefficient for 3rd order weno smoother real(rkind),save :: epsilon3 !2nd coefficient for 3rd order weno smoother integer, save :: nquad !number of quad points used for 3rd order weno !levels of time discretization, mainly for testing purposes integer, save :: ntd_weno !(1) one-level, reduces to Euler; (2) not implemented yet; (3) 3rd-order Runge-Kutta temporal discretization (Shu and Osher, 1988) !Elad filter integer, save :: ielad_weno !switch for elad filter, not used at the moment integer, save :: i_prtnftl_weno !switch for printing invalid T/S to nonfatal_* real(rkind),save :: small_elad !criteria for ELAD, not used at the moment real(rkind),save,allocatable :: xqp(:,:),yqp(:,:) !quadrature point coordinates integer,save :: mnweno1 !maxium number of p1 polynomial integer,save,allocatable :: nweno1(:) !number of p1 polynomial !stencil of P1 polynomial (3 elements #,mnweno1 polynomials,ne) integer,save,allocatable :: isten1(:,:,:) !polynomial coefficients at quadrature points (3 poly. coeffcients, mnweno1 ploynomials, 2 quadrature points, 3 sides, ne elem.) real(rkind),save,allocatable :: wmat1(:,:,:,:,:) !coefficients for final p1 polynomial weight (3 components, xy direction,mnweno1 polynomials,ne) real(rkind),save,allocatable :: wts1(:,:,:,:) !stencil quality, check if at least 1 stencil is on one side of an element side logical,save,allocatable :: isten_qual1(:), isten_qual2(:) integer,save :: mnweno2 !maxium number of p2 polynomial integer,save,allocatable :: nweno2(:) !number of p2 polynomial !stencil of P2 polynomial (6 element #,mnweno2 polynomials,ne) integer,save,allocatable :: isten2(:,:,:) !polynomial coefficient at quadrature points (3 poly. coeffcients, mnweno2 ploynomials, 2 quadrature points, 3 sides, ne elem.) real(rkind),save,allocatable :: wmat2(:,:,:,:,:) !coefficient for calculating final p2 polynomial weight (6 directions,5 components,mnweno2 polynomials,ne) real(rkind),save,allocatable :: wts2(:,:,:,:) !coefficients used in smoother calculation (3,ne) real(rkind),save,allocatable :: fwts2(:,:) !diagnostic variables integer,save,allocatable :: ie_all_stencils1(:,:,:),n_all_stencils1(:) !(3,mnweno1,ne), (ne), keep a record of all p1 stencils real(rkind),save,allocatable :: det_all_stencils1(:,:) !matrix determinants of all p1 stencils integer,save,allocatable :: iremove1(:,:) !(3 x mnweno1, ne), keep a record of the p1 stencils removed due to small determinants integer,save,allocatable :: nremove1(:) !size: ne, number of p1 stencils removed at each element real(rkind),save,allocatable :: rremove1(:,:) !size: ne, determinants (under local coordinates) of p1 stencils removed at each element integer,save,allocatable :: ie_all_stencils2(:,:,:),n_all_stencils2(:) !(6,mnweno2,ne), (ne), keep a record of all p2 stencils real(rkind),save,allocatable :: det_all_stencils2(:,:) !matrix determinants of all p2 stencils integer,save,allocatable :: iremove2(:,:) !(6 x mnweno2, ne), keep a record of the p2 stencils removed due to small determinants integer,save,allocatable :: nremove2(:) !size: ne, number of p2 stencils removed at each element real(rkind),save,allocatable :: rremove2(:,:) !size: ne, determinants (under local coordinates) of p2 stencils removed at each element !0 in ! this case. !3) iwsett(itrc)=0: variable settling (or swimming) vel along vertical. Code ! will assume wsett(:,kbe|nvrt,:)=0 to accumulate mass there. real(rkind),save,allocatable :: wsett(:,:,:) integer,save,allocatable :: iwsett(:) !iwsett(ntracers) integer,save,allocatable :: inu_pts_gb(:,:) !Declarations for other modules ! ANALYSIS real(rkind),save,allocatable :: dtbe(:) ! vertical flux diversion closure fraction applied at surface ! real(rkind) :: vclose_surf_frac ! 1.0:flux applied at surface, 0.5:half at top half at bottom ! WWM !#ifdef USE_WWM integer,save :: msc2,mdc2 real(rkind),save,allocatable :: wwave_force(:,:,:), jpress(:), sbr(:,:), sbf(:,:), srol(:,:),sds(:,:), eps_w(:), eps_r(:),eps_br(:) real(rkind),save,allocatable :: stokes_hvel(:,:,:), stokes_wvel(:,:), stokes_hvel_side(:,:,:), stokes_wvel_side(:,:) real(rkind),save,allocatable :: roller_stokes_hvel(:,:,:), roller_stokes_hvel_side(:,:,:) !real(rkind),save,allocatable :: stokes_vel(:,:,:), stokes_w_nd(:,:), stokes_vel_sd(:,:,:) real,save,allocatable :: out_wwm(:,:), out_wwm_rol(:,:), out_wwm_windpar(:,:) real(rkind),save,allocatable :: tanbeta_x(:), tanbeta_y(:) ! MP from KM: bottom slope real(rkind),save,allocatable :: curx_wwm(:),cury_wwm(:) !BM: coupling current real(rkind),save,allocatable :: wafo_opbnd_ramp(:) !BM: ramp on wave forces at open boundary real(rkind),save,allocatable :: taub_wc(:) real(rkind),save,allocatable :: delta_wbl(:) real(rkind),save,allocatable :: wave_sbrtot(:) real(rkind),save,allocatable :: wave_sbftot(:) real(rkind),save,allocatable :: wave_sintot(:) real(rkind),save,allocatable :: wave_sdstot(:) !#endif ! TIMOR !#ifdef USE_TIMOR real(rkind), parameter :: rhosed = 2650. ! kg/m3 logical,save :: laddmud_d !switch on/off density effects logical,save :: laddmud_v !switch on/off viscosity effects real(rkind),save,allocatable :: vts(:,:) !vts(nvrt,npa) rheological viscosity [m^2/s] ! real(rkind),save,allocatable :: wsink(:,:,:) !wsink([tr],nvrt,npa); sink velocity>=0 [m/s] real(rkind),save,allocatable :: rhomud(:,:,:) ! rhomud([tr],nvrt,npa): Mud floc particle density [kg/m3] !#endif /*USE_TIMOR*/ !#ifdef USE_SED real(rkind),save,allocatable :: dave(:),total_sus_conc(:,:) INTEGER :: ddensed ! activation key for sediment density effects on water density !#endif !USE_AGE integer,save,allocatable :: nelem_age(:),ielem_age(:,:),level_age(:) !NAPZD !#ifdef USE_NAPZD real(rkind),save,allocatable :: Bio_bdef(:,:) ! biological deficit for NAPZD model !#endif !Harmonic analysis !#ifdef USE_HA integer,save :: MNHARF,NTSTEPS,ITMV,ITHAS,NHAGE,NHAGV,ITHAF,ICHA,NHAINC logical,save :: CHARMV real(rkind),save :: TIMEBEG,FMV real(rkind),save,allocatable :: XVELAV(:),YVELAV(:),XVELVA(:),YVELVA(:),ELAV(:),ELVA(:) !#endif !#ifdef USE_ECO !... MFR, Nov/2015 - Updated, ECO uses nws=2 !... MFR - other variables to atmospheric parameters (when nws=0 ... probably to clean later...) ! real(rkind),save,allocatable :: Pair(:), Tair(:), Hair(:), Uwind(:), Vwind(:), cloud(:) !... MFR - Tracer models real(rkind),save :: tr_tmp1 !#endif #ifdef USE_SIMPLE_WIND ! nws=5,6 option real(rkind),save,allocatable :: cf_x1(:),cf_x2(:),cf_y1(:),cf_y2(:),cf_denom(:) integer,save,allocatable :: cf_i(:), cf_j(:) #endif ! Marsh model integer,save,allocatable :: imarsh(:),ibarrier_m(:) ! SAV real(rkind),save,allocatable :: sav_alpha(:),sav_h(:),sav_nv(:),sav_di(:),sav_cd(:) !Tsinghua group:0825 REAL(rkind),save :: Cbeta,beta0,c_miu,Cv_max,ecol,ecol1,sigf,sigepsf,Ceps1,Ceps2,Ceps3,Acol,sig_s,fi_c,ksi_c,kpz !1013+kpz REAL(rkind),save,ALLOCATABLE :: Dpzz(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: Tpzz(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: Vpx(:,:) !x drift velocity at side centers & whole levels REAL(rkind),save,ALLOCATABLE :: Vpy(:,:) !y drift velocity at side centers & whole levels REAL(rkind),save,ALLOCATABLE :: Vpx2(:,:),Vpy2(:,:),Vpz2(:,:) !x,y,z drift velocity at nodes & whole levels 0927.1 1006 REAL(rkind),save,ALLOCATABLE :: TDxz(:,:),TDyz(:,:) !x,y diff stress tensor side centers & whole levels 1006 REAL(rkind),save,ALLOCATABLE :: Phai(:,:,:) !correction of settle term at nodes & whole levels 1006 REAL(rkind),save,ALLOCATABLE :: dfhm(:,:,:) !diff in transport Eq. at nodes & whole levels 1007 REAL(rkind),save,ALLOCATABLE :: Dpxz(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: Dpyz(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: taufp_t(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: ws(:,:) !Ws at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: epsf(:,:) !epsf at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: q2fp(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: q2p(:,:) !kp at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: q2f(:,:) !kp at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: kppian(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: miuft(:,:) !miuft at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: miuepsf(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: trndtot(:,:) !tot sed volumetric conc at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: miup(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: miup_t(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: miup_c(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: Kp_tc(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: Kp_t(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: Kp_c(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: Kft(:,:) !Kf at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: taup(:,:) !relax time at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: taup_c(:,:) !collide time at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: g0(:,:) !at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: SDav(:,:) !average SD50 at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: Srhoav(:,:) !average Srho at nodes & whole levels REAL(rkind),save,ALLOCATABLE :: kesit(:,:) !kesi_tau Srho at nodes & whole levels !0821 !Tsinghua group---------------------------------- !End module variables declaration contains subroutine release_gl(n,gl_ll) ! Free memory associated with global-to-local linked-list implicit none integer,intent(in) :: n type(llist_type),pointer :: gl_ll(:) integer i do i=1,n call release_llist(gl_ll(i)%next) enddo deallocate(gl_ll) nullify(gl_ll) end subroutine release_gl recursive subroutine release_llist(llentry) ! Free memory associated with linked-list implicit none type(llist_type),pointer :: llentry if(associated(llentry)) then call release_llist(llentry%next) deallocate(llentry) nullify(llentry) endif end subroutine release_llist end module schism_glbl !=============================================================================== !=============================================================================== ! END SCHISM GLOBAL DATA MODULE !=============================================================================== !===============================================================================