!-----------------------------------------------------------------------
! SCHISM model parameter inputs.
! Format rules for param.in:
! (1) Lines beginning with "!" are comments; blank lines are ignored;
! (2) one line for each parameter in the format: keywords= value;
!     keywords are case sensitive; spaces allowed between keywords and "=" and value;
!     comments starting with "!"  allowed after value;
! (3) value is an integer, double, or 2-char string (no single quote needed); for double, any of the format is acceptable:
!     40 40. 4.e1
!     Use of decimal point in integers is OK but discouraged.
! (4) If duplicate entries are present, the first entry gets used.
!-----------------------------------------------------------------------

!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
! Model configuration parameters 
!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
! Starting time
  start_year = 2011 !int
  start_month = 7 !int
  start_day = 27 !int
  start_hour = 0 !double
  utc_start = 0 !double
!-----------------------------------------------------------------------
! Coordinate option: 1: Cartesian; 2: lon/lat (hgrid.gr3=hgrid.ll in this case,
! and orientation of element is outward of earth)
!-----------------------------------------------------------------------
  ics = 2 !Coordinate option

!-----------------------------------------------------------------------
! Pre-processing option. Useful for checking grid violations.
!-----------------------------------------------------------------------
  ipre = 0 !Pre-processor flag (1: on; 0: off)

!-----------------------------------------------------------------------
! Equation of State type used
! ieos_type=0: UNICEF (nonlinear); =1: linear function of T ONLY, i.e. 
! \rho=eos_b+eos_a*T, where eos_a<=0 in kg/m^3/C
!-----------------------------------------------------------------------
  ieos_type = 0
  ieos_pres = 0
!  eos_a = -0.1 !needed if ieos_type=1; should be <=0 
!  eos_b = 1001. !needed if ieos_type=1

!-----------------------------------------------------------------------
! If WWM is used, set coupling/decoupling flag. Not used if USE_WWM is distabled in Makefile
!       0: decoupled so 2 models will run independently; 
!       1: full coupled (elevation, vel, and wind are all passed to WWM); 
!       2: elevation and currents in wwm, no wave force in selfe;  
!       3: no elevation and no currents in wwm, wave force in selfe;
!       4: elevation but no currents in wwm, wave force in selfe;
!       5: elevation but no currents in wwm, no wave force in selfe;
!       6: no elevation but currents in wwm, wave force in selfe;
!       7: no elevation but currents in wwm, no wave force in selfe;
! Note that all these parameters must be present in this file (even though not used).
!-----------------------------------------------------------------------
  icou_elfe_wwm = 0 
  nstep_wwm = 1 !call WWM every this many time steps. If /=1, consider using quasi-steady mode in WWM
  iwbl = 0 !1: modified Grant-Madsen formulation for wave boundary layer; used only if icou_elfe_wwm/=0; if icou_elfe_wwm=0, set iwbl=0
  msc2 = 24 !same as MSC in .nml ... for consitency check between SCHISM and WWM-II
  mdc2 = 30 !same as MDC in .nml
  hmin_radstress = 1. !min. total water depth used only in radiation stress calculation [m]
  nrampwafo = 1       !ramp-up option for the wave forces (1: on; 0: off)
  drampwafo = 1.      !needed if nrampwafo=1; ramp-up period in days
  turbinj = 0.15      !% of depth-induced wave breaking energy injected in turbulence (default: 0.15 (15%), as proposed by Feddersen, 2012)

!-----------------------------------------------------------------------
! Hydrological option (0 for normal simulations)
!-----------------------------------------------------------------------
  ihydlg = 0

!-----------------------------------------------------------------------
! Non-hydrostatic model switch (0: hydrostatic model; 1: non-hydrostatic model)
! For non-hydrostatic model, use ihydro_region to indicate if  hydrostatic 
! region needs to be specified in hydro_region.gr3 (depth=1 in hydrostatic region)
!-----------------------------------------------------------------------

!-----------------------------------------------------------------------
! Define # of tracers in each module (if enabled)
!-----------------------------------------------------------------------
  ntracer_gen = 2 !user defined module
  ntracer_age = 4 !age calculation. Must be =2*N where N is # of age tracers
  sed_class = 5 !SED3D
  eco_class = 27 !EcoSim: must be between [25,60]

!-----------------------------------------------------------------------
! Bed deformation option (0: off; 1: vertical deformation only; 2: 3D bed deformation). 
! If imm=1, bdef.gr3 is needed; if imm=2, user needs to update depth info etc
! in the code (not working for ics=2 yet).
!-----------------------------------------------------------------------
  imm = 0
! ibdef = 10 !needed if imm=1; # of steps used in deformation

!-----------------------------------------------------------------------
! Reference latitude for beta-plane approximation when ncor=1 (not used if ics=2)
!-----------------------------------------------------------------------
  cpp_lon = -77.07  !lon - not really used
  cpp_lat = 24 !lat

!-----------------------------------------------------------------------
! Baroclinic/barotropic option. If ibcc=0 (baroclinic model), itransport is not used.
!-----------------------------------------------------------------------
  ibcc = 0 !Baroclinic option
  itransport = 0 
  nrampbc = 0 !ramp-up flag for baroclinic force
  drampbc = 4. !not used if nrampbc=0

!-----------------------------------------------------------------------
! Baroclinicity calculation in off/nearshore. The 'below-bottom' gradient
! is zeroed out if h>=h2_bcc (i.e. like Z) or uses const extrap
! (i.e. like terrain-following) if h<=h1_bcc(<h2_bcc) (and linear
! transition in between based on local depth)
!-----------------------------------------------------------------------
  h1_bcc = 50. ![m]
  h2_bcc = 100. ![m]

!-----------------------------------------------------------------------
! Hotstart option. 0: cold start; 1: hotstart with time reset to 0; 2: 
! continue from the step in hotstart.nc
!-----------------------------------------------------------------------
  ihot = 1

!-----------------------------------------------------------------------
! Hydraulic model option. If ihydraulics/=0, hydraulics.in 
! is required. This option cannot be used with non-hydrostatic model.
!-----------------------------------------------------------------------
  ihydraulics = 0

!-----------------------------------------------------------------------
! Point sources/sinks option (0: no; 1: on). If =1, needs source_sink.in,
! vsource,th, vsink.th, and msource.th
!-----------------------------------------------------------------------
  if_source = 1
  nramp_ss = 0 !needed if if_source=1; ramp-up flag for source/sinks
  dramp_ss = 1 !needed if if_source=1; ramp-up period in days

!-----------------------------------------------------------------------
! Earth's radii at pole and equator (to define an ellipsoid)
!-----------------------------------------------------------------------
  rearth_pole = 6378206.4
  rearth_eq = 6378206.4

!-----------------------------------------------------------------------
! Method for momentum advection. 0: ELM; 1: upwind (not quite working yet)
!-----------------------------------------------------------------------
  iupwind_mom = 0

!-----------------------------------------------------------------------
! Methods for computing velocity at nodes. 
! If indvel=0, conformal linear shape function is used; if indvel=1, averaging method is used.
! For indvel=0, some stabilization method is needed (see below). 
!-----------------------------------------------------------------------
  indvel = 0
 
!-----------------------------------------------------------------------
! 2 stabilization methods, mostly for indvel=0.
! (1) Horizontal viscosity option. ihorcon=0: no viscosity is used; =1: Lapacian;
! =2: bi-harmonic. If ihorcon=1, horizontal viscosity _coefficient_ (<=1/8, related
! to diffusion number) is given in hvis_coef0, and the diffusion # 
! is problem dependent; [0.001-1/8] seems to work well.
! If ihorcon=2, diffusion number is given by hvis_coef0 (<=0.025).
! If indvel=1, no hvis is needed. 
! (2) Shapiro filter (see below)
!-----------------------------------------------------------------------
  ihorcon = 2
  hvis_coef0 = 0.025 !const. diffusion # if ihorcon/=0; <=0.025 for ihorcon=2, <=0.125 for ihorcon=1
!  cdh = 0.01 !needed only if ihorcon/=0; land friction coefficient - not active yet

!-----------------------------------------------------------------------
! 2nd stabilization method via Shapiro filter. This should normally be used 
! if indvel=ihorcon=0. To transition between eddying/non-eddying regimes, use
! indvel=0, ihorcon/=0, and ishapiro=-1 (shapiro.gr3).
!-----------------------------------------------------------------------
  ishapiro = -1 !on/off flag
! shapiro = 0.5 !Shapiro filter strength, needed only if ishapiro=1; max is 0.5 
  niter_shap = 1 !needed of ishapiro/=0 - # of iterations with Shapiro filter. Suggested: 1

!-----------------------------------------------------------------------
! Horizontal diffusivity option. if ihdif=1, horizontal diffusivity is given in hdif.gr3
!-----------------------------------------------------------------------
  ihdif = 0 

!-----------------------------------------------------------------------
! Bottom friction. 
!           bfric=0: drag coefficients specified in drag.gr3; bfric=-1: Manning's 
!           formulation (even for 3D prisms). 
!           bfric=1: bottom roughness (in meters) specified in rough.gr3 (and in this case, negative
!           or 0 depths in rough.gr3 indicate time-independent Cd, not roughness!).
!           Cd is calculated using the log law, when dzb>=dzb_min; when dzb<dzb_min,
!           Cd=Cdmax*exp[dzb_decay*(1-dzb/dzb_min)], where Cdmax=Cd(dzb=dzb_min),
!           and dzb_decay (<=0) is a decay const specified below.
!           If iwbl=1, bfric must =1.
!-----------------------------------------------------------------------
  bfric = 1 !nchi in code
  dzb_min = 0.5 !needed if bfric==1; min. bottom boundary layer thickness [m].
  dzb_decay = 0.  !needed if bfric=1; a decay const. [-]
  hmin_man = 1. !needed if bfric=-1: min. depth in Manning's formulation [m]

!-----------------------------------------------------------------------
! Coriolis. If ncor=-1, specify "latitude" (in degrees); if ncor=0,
! specify Coriolis parameter in "coriolis"; if ncor=1, model uses
! lat/lon in hgrid.ll for beta-plane approximation if ics=1, and in this case,
! the latitude specified in CPP projection ('cpp_lat') is used. If ncor=1 and ics=2,
! Coriolis is calculated from local latitude, and 'cpp_lat' is not used.
!-----------------------------------------------------------------------
  ncor = 1 !must be 1 if ics=2
  !latitude = 46 !if ncor=-1
  !coriolis = 1.e-4 !if ncor=0

!-----------------------------------------------------------------------
! Elevation initial condition flag. If ic_elev=1, elev.ic (in *.gr3 format) is needed
! to specify the initial elevations; otherwise elevation is initialized to 0 everywhere 
! (cold start only)
!-----------------------------------------------------------------------
  ic_elev = 0

!-----------------------------------------------------------------------
! Elevation boundary condition ramp-up flag. =0: ramp up from 0; =1: ramp up from
! elev. values read in from elev.ic or hotstart.nc - if neither is present, from 0.
! This flag is mainly used to start the simulation from non-zero elev.
! The ramp-up period is same as 'dramp' below.
!-----------------------------------------------------------------------
  nramp_elev = 0

!-----------------------------------------------------------------------
! Optional inverse barometric effects on the elev. b.c.
! If inv_atm_bnd=1, the elev.'s at boundary are corrected by the difference
! between the actual atmos. pressure and a reference pressure (prmsl_ref below)
!-----------------------------------------------------------------------
  inv_atm_bnd = 0 !0: off; 1: on
  prmsl_ref = 101325. !reference atmos. pressure on bnd [Pa]

!-----------------------------------------------------------------------
! Initial condition for T,S. This value only matters for ihot=0 (cold start).
! If ic_*=1, the initial T,S field is read in from temp.ic and salt.ic (horizontally varying).
! If ic_*=2, the initial T,S field is read in from ts.ic (vertical varying).
! If ihot=0 && ic_*=2 || ibcc_mean=1, ts.ic is used for removing mean density profile.
!-----------------------------------------------------------------------
  ic_TEM = 0
  ic_SAL = 0 !must be same as ic_TEM

! initial conditions for other tracers.
! 1: needs inputs [MOD]_hvar_[1,2,...].ic ('1...' is tracer id); format of each file is similar to salt.ic;
!    i.e. horizontally varying i.c. is used for each tracer.
! 2: needs [MOD]_vvar_[1,2,...].ic. Format of each file (for each tracer in tis MOD) is similar to ts.ic
!    (i.e. level #, z-coord., tracer value). Verically varying i.c. is used for each tracer.
! 0: model sets own i.c. (EcoSim; TIMOR)
  ic_GEN = 1 !user defined module
  ic_AGE = 1 !Age
  ic_SED = 1 !SED3D
  ic_ECO = 1 !EcoSim
  ic_ICM = 1 !ICM
  ic_COS = 1 !CoSINE
  ic_FIB = 1 !FIB

!-----------------------------------------------------------------------
! Settling vel [m/s] for GEN module (positive downward)
!-----------------------------------------------------------------------
  gen_wsett=1.e-4

!-----------------------------------------------------------------------
! Mean T,S profile option. If ibcc_mean=1 (or ihot=0 and ic_TEM=2), mean profile
! is read in from ts.ic, and will be removed when calculating baroclinic force.
! No ts.ic is needed if ibcc_mean=0.
!-----------------------------------------------------------------------
  ibcc_mean = 0 

!-----------------------------------------------------------------------
! Max. horizontal velocity magnitude, used mainly to prevent problem in 
! bulk aerodynamic module
!-----------------------------------------------------------------------
  rmaxvel = 3.5

!-----------------------------------------------------------------------
!  Following 2 parameters control backtracking
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!  min. vel for invoking btrack and for abnormal exit in quicksearch
!-----------------------------------------------------------------------
  velmin_btrack = 1.e-4
!-----------------------------------------------------------------------
! Nudging factors for starting side/node - add noise to avoid underflow
! The starting location is nudged to: old*(1-btrack_nudge)+btrack_nudge*centroid
! Suggested value: btrack_nudge=9.013e-3
!-----------------------------------------------------------------------
  btrack_nudge= 9.013e-3 

!-----------------------------------------------------------------------
! Wetting and drying. If ihhat=1, \hat{H} is made non-negative to enhance
! robustness near wetting and drying; if ihhat=0, no retriction is imposed for
! this quantity. 
! inunfl=0 is used for normal cases and inunfl=1 is used for more accurate wetting
! and drying if grid resolution is sufficiently fine.
!-----------------------------------------------------------------------
  ihhat = 1 
  inunfl = 0
  h0 = 1.e-6 !min. water depth for wetting/drying
  shorewafo = 0

!-----------------------------------------------------------------------
! Implicitness factor (0.5<thetai<=1).
!-----------------------------------------------------------------------
  thetai = 1

!-----------------------------------------------------------------------
! Run time and ramp option
!-----------------------------------------------------------------------
  rnday = 50 !total run time in days
  nramp = 0 !ramp-up option (1: on; 0: off)
  dramp = 1. !needed if nramp=1; ramp-up period in days
  dt = 150. !Time step in sec

!-----------------------------------------------------------------------
! Solver option. JCG is used presently.
!-----------------------------------------------------------------------
  slvr_output_spool = 50 !output spool for solver info
  mxitn = 1500 !max. iteration allowed
  tolerance = 1.e-12 !error tolerance

!-----------------------------------------------------------------------
! Advection (ELM) option. If nadv=1, backtracking is done using Euler method; 
! nadv=2, using 2nd order Runge-Kutta; if nadv=0, advection in momentum 
! is turned off/on in adv.gr3 (the depths=0,1, or 2 also control methods 
! in backtracking as above). dtb_max/min are the max/min steps allowed -
! actual step is calculated adaptively based on local gradient.
!-----------------------------------------------------------------------
  nadv = 1
  dtb_max = 30. !in sec
  dtb_min = 10.

!-----------------------------------------------------------------------
! If inter_mom=0, linear interpolation is used for velocity at foot of char. line.
! If inter_mom=1 or -1, Kriging is used, and the choice of covariance function is
! specified in 'kr_co'. If inter_mom=1, Kriging is applied to whole domain;
! if inter_mom=-1, the regions where Kriging is used is specified in krvel.gr3 
! (depth=0: no kriging; depth=1: with kriging). 
!-----------------------------------------------------------------------
  inter_mom = 0
  kr_co = 2 !not used if inter_mom=0

!-----------------------------------------------------------------------
! Transport method. 
! If itr_met=1, upwind method is used. If itr_met=2 or 3, 2nd-order TVD method is used
! on an element/prism if the total depth (at all 3 nodes of the elem.)>=h_tvd and the flag in
! tvd.prop = 1 for the elem. (tvd.prop is required in this case); 
! otherwise upwind is used for efficiency. 
! itr_met=3: implicit TVD in the vertical dimension.
!-----------------------------------------------------------------------
  itr_met = 4
  h_tvd = 5. !used only if itr_met>=2; cut-off depth (m) 
  !If itr_met=3 or 4, need the following 2 tolerances of convergence. The convergence
  !is achieved when sqrt[\sum_i(T_i^s+1-T_i^s)^2]<=eps1_tvd_imp*sqrt[\sum_i(T_i^s)^2]+eps2_tvd_imp
  eps1_tvd_imp = 1.e-4 !suggested value is 1.e-4, but for large suspended load, need to use a smaller value (e.g. 1.e-9)
  eps2_tvd_imp = 1.e-14  

  !if itr_met = 4, the following parameters are needed
  !if itr_met=4 and ipre=1, diagnostic outputs are generated for weno accuracy and stencil quality, 
  !  see subroutine weno_diag in src/Hydro/misc_subs.F90 for details
  ip_weno = 2   !order of accuracy: 0- upwind; 1- linear polynomial, 2nd order; 2- quadratic polynomial, 3rd order
  courant_weno=0.5 !Courant number for weno transport
  nquad = 2  !number of quad points on each side, nquad= 1 or 2
  ntd_weno = 1 !order of temporal discretization: (1) Euler (default); (3): 3rd-order Runge-Kutta (only for benchmarking)
  epsilon1 = 1.e-3   !coefficient for 2nd order weno smoother
  epsilon2 = 1.e-5  !1st coefficient for 3rd order weno smoother (larger values are more prone to numerical dispersion
                     !, 1.e-10 should be fairly safe, recommended values: 1.e-8 ~ 1.e-6 (from the real applications so far)
  i_prtnftl_weno = 1 !option for writing nonfatal errors on invalid temp. or salinity for density: (0) off; (1) on.

  !Inactive at the moment:
  epsilon3 = 1.e-25  !2nd coefficient for 3rd order weno smoother (inactive at the moment)
  !Elad filter has not been implemented yet; preliminary tests showed it might not be necessary
  ielad_weno = 0      !ielad, if ielad=1, use ELAD method to suppress dispersion (inactive at the moment)
  small_elad = 1.e-4  !small (inactive at the moment)

!-----------------------------------------------------------------------
! Atmos. option. If nws=0, no atmos. forcing is applied. If nws=1, atmos.
! variables are read in from wind.th. If nws=2, atmos. variables are
! read in from sflux_ files.
! If nws=4, ascii format is used for wind and atmos. pressure at each node (see source code).
! If nws=5 (or 6), netcdf format for wind and pressure are read in from UVP.nc (UVP_direct.nc).
!   nws=6 assumes the interpolation has been done onto the hgrid.ll, while nws=5 uses structured-grid
!   as atmos. input.
! If nws>0, 'iwindoff' can be used to scale wind speed (with windfactor.gr3).
!
! Stress calculation:
! If nws=1 or >=4, or nws=2 and ihconsv=0, or nws=2 and iwind_form=-1, 
! the stress is calculated from Pond & Pichard formulation 
! If nws=2, ihconsv=1 and iwind_form=0, the stress is calculated from heat exchange
! routine; 
! If WWM is enabled and icou_elfe_wwm>0 and iwind_form=-2, stress is calculated by WWM;
! otherwise the formulations above are used.
!-----------------------------------------------------------------------
  nws = 2 
  wtiminc = 150. !time step for atmos. forcing
  nrampwind = 0 !ramp-up option for atmos. forcing
  drampwind = 2. !needed if nrampwind/=0; ramp-up period in days
  iwindoff = 0 !needed only if nws/=0; '1': needs windfactor.gr3
  iwind_form = -1 !change to -2 after adding wave model 
  !If impose_net_flux/=0 and nws=2, read in net _surface_ heat flux as var 'dlwrf' 
  !(Downward Long Wave) in sflux_rad (solar radiation is still used separately),
  !and if PREC_EVAP is on, also read in net P-E as 'prate' (Surface Precipitation Rate) in sflux_prc.
  impose_net_flux = 0

!-----------------------------------------------------------------------
! Heat and salt exchange. isconsv=1 needs ihconsv=1; ihconsv=1 needs nws=2.
! If isconsv=1, need to compile with precip/evap module turned on.
!-----------------------------------------------------------------------
  ihconsv = 1 !heat exchange option
  isconsv = 0 !evaporation/precipitation model

!-----------------------------------------------------------------------
! Turbulence closure.
!-----------------------------------------------------------------------
  itur = 3 
!  dfv0 = 1.e-6 !needed if itur=0
!  dfh0 = 1.e-6 !needed if itur=0
  turb_met = KL !needed if itur=3
  turb_stab = KC !needed if itur=3. Use 'GA' if turb_met='MY'; otherwise use 'KC'
  xlsc0 = 0.1 !needed if itur=3. Scale for surface & bottom mixing length (>0)

!-----------------------------------------------------------------------
! Sponge layer for elevation and vel.
! If inu_elev=0, no relaxation is applied to elev.
! If inu_elev=1, relax. constants are specified in elev_nudge.gr3
! and applied to eta=0 (thus a depth=0 means no relaxation).
! Similarly for inu_uv (with input uv_nudge.gr3)
!-----------------------------------------------------------------------
  inu_elev = 0
  inu_uv = 0

!-----------------------------------------------------------------------
! Nudging options for tracers. If inu_[MOD]=0, no nudging is used. If inu_[MOD]=1,
! nudge to initial condition according to relaxation constants specified
! in [MOD]_nudge.gr3. If inu_[MOD]=2, nudge to values in [MOD]_nu,in
! (with step 'step_nu_tr') according to [MOD]_nudge.gr3. 
! The final relaxation = horizontal relax (specified in [MOD]_nudge.gr3) times dt.
! [MOD] are tracer model names.
!-----------------------------------------------------------------------
  inu_TEM = 2
  inu_SAL = 2
  inu_GEN = 0 !user defined
  inu_AGE = 0 !Age
  inu_SED = 0 !SED3D
  inu_ECO = 0 !EcoSim
  inu_ICM = 0 !ICM
  inu_COS = 0 !CoSINE
  inu_FIB = 0 !FIB

  step_nu_tr = 86400. !time step [sec] in all [MOD]_nu.in (for inu_[MOD]=2)

  vnh1 = 400 !vertical nudging depth 1
  vnf1 = 0. !vertical relax \in [0,1]
  vnh2 = 401 !vertical nudging depth 2 (must >vnh1)
  vnf2 = 0. !vertical relax

!-----------------------------------------------------------------------
! Cut-off depth for cubic spline interpolation near bottom when computing horizontal gradients
! e.g. using hgrad_nodes() (baroc. force, radiation stress, and gradients of qnon and qhat in non-hydro model). 
! If depth > depth_zsigma ('deep'),
! a min. (e.g. max bottom z-cor for the element) is imposed in the spline and so a more
! conservative method is used without extrapolation beyond bottom; 
! otherwise constant extrapolation below bottom is used.
!-----------------------------------------------------------------------
  depth_zsigma = 100. !h_bcc1

!-----------------------------------------------------------------------
! Dimensioning parameters for inter-subdomain btrack. 
! If error occurs like 'bktrk_subs: overflow' or 'MAIN: nbtrk > mxnbt'
! gradually increasing these will solve the problem
!-----------------------------------------------------------------------
  s1_mxnbt = 0.5
  s2_mxnbt = 3.5

!-----------------------------------------------------------------------
! Global output options.
!-----------------------------------------------------------------------
  nspool = 24 !output step spool
  ihfskip = 2880 !stack spool; every ihfskip steps will be put into 1_*, 2_*, etc...

  elev = 1 !0: off; 1: on - elev. [m]
  pres = 0 !air pressure [Pa]
  airt = 0 !air temperature [C]
  shum = 0 !Specific humidity [-]
  srad = 0 !solar (shortwave) radiation [W/m/m]
  flsu = 0 !sensible flux (positive upward) [W/m/m] 
  fllu = 0 !latent heat flux (positive upward) [W/m/m]
  radu = 0 !upward longwave radiation (positive upward) [W/m/m]
  radd = 0 !downward longwave radiation (positive downward) [W/m/m]
  flux = 0 !total flux=-flsu-fllu-(radu-radd) [W/m/m]
  evap = 0 !evaporation rate [kg/m/m/s]
  prcp = 0 !precipitation rate [kg/m/m/s]
  bdrc = 0 !Bottom drag coefficient [-]
  wind = 0 !wind speed [m/s]
  wist = 0 !wind stress [m^2/s/s]
  dahv = 1 !depth-averaged vel. [m/s]
  vert = 1 !vertical velocity [m/s]
  temp = 1 !water temperature [C]
  salt = 1 !water salinity [PSU]
  conc = 0 !water density [kg/m^3]
  tdff = 0 !eddy diffusivity [m^2/s]
  vdff = 1 !eddy viscosity [m^2/s]
  kine = 0 !turbulent kinetic energy
  mixl = 0 !turbulent mixing length [m]
  zcor = 1 !z-coordinates [m]
  qnon = 0 !non-hydrostatic pressure
  hvel = 1 !horizontal vel. [m/s]

!-----------------------------------------------------------------------
! Outputs from WWM (USE_WWM must be on in Makefile)
!-----------------------------------------------------------------------
  wwm_1 = 0 !sig. height (m)
  wwm_2 = 0 !Mean average period (sec) - TM01
  wwm_3 = 0 !Zero down crossing period for comparison with buoy (s) - TM02
  wwm_4 = 0 !Average period of wave runup/overtopping - TM10
  wwm_5 = 0 !Mean wave number (1/m)
  wwm_6 = 0 !Mean wave length (m)
  wwm_7 = 0 !Mean average energy transport direction (Average direction; deg)
  wwm_8 = 0 !Mean directional spreading (deg)
  wwm_9 = 0 !Discrete peak period (sec)
  wwm_10 = 0 !Continuous peak period (Tp) based on higher order moments (sec) 
  wwm_11 = 0 !Peak phase vel. (m/s)
  wwm_12 = 0 !Peak n-factor [-]
  wwm_13 = 0 !Peak group vel. (m/s)
  wwm_14 = 0 !Peak wave number (1/m)
  wwm_15 = 0 !Peak wave length (m)
  wwm_16 = 0 !Peak (dominant) wave direction (degr) ... some buoys record this 
  wwm_17 = 0 !Peak directional spreading (deg) ... some buoys record this  
  wwm_18 = 0 !Discrete peak direction (deg) ... some buoys record this  
  wwm_19 = 0 !Orbital vel. (m/s)
  wwm_20 = 0 !RMS orbital vel. (m/s)
  wwm_21 = 0 !Bottom excursion period (sec)
  wwm_22 = 0 !bottom wave period (sec)
  wwm_23 = 0 !Ursell number based on peak period 
  wwm_24 = 0 !none 
  wwm_25 = 0 !Etot energy vector (m^2)
  wwm_26 = 0 !none 

!-----------------------------------------------------------------------
! Outputs for user-defined tracer module
!-----------------------------------------------------------------------
  GEN_1 = 0 
  GEN_2 = 0

!-----------------------------------------------------------------------
! Outputs for (age)
!-----------------------------------------------------------------------
  AGE_1 = 0 !indices from "1" to "ntr/2"; [days]
  AGE_2 = 0

!-----------------------------------------------------------------------
! Specific outputs in SED3D (USE_SED must be on in Makefile;
! otherwise these are not needed)
!-----------------------------------------------------------------------
  SED_1 = 0 !conc. of 1st class (one output need by each class) [g/L]
  SED_2 = 0 
  SED_bfrac_1 = 0 ! Bed fraction 1st tracer (one output need by each class) [-]
  SED_bfrac_2 = 0
  SED_qbdl_1 = 0 ! Bedload transport rate vector (kg.m-1.s-1) for 1st tracer (one output need by tracer)
  SED_qbdl_2 = 0

  SED_depth = 0 !bottom depth _change_ from init. condition (m)
  SED_bedd50  = 0 ! Bed median grain size in the active layer (mm)
  SED_bstress = 0 ! Bottom shear stress (Pa)
  SED_brough  = 0 ! Bottom roughness lenghth (mm)

!-----------------------------------------------------------------------
! EcoSim outputs 
!-----------------------------------------------------------------------
  ECO_1 = 0
  ECO_2 = 0

!-----------------------------------------------------------------------
! ICM outputs 
!-----------------------------------------------------------------------
  ICM_1 = 1 !Zoo. #1
  ICM_2 = 1 !Zoo. #2
  ICM_3 = 1 !phyto #1
  ICM_4 = 1 !phyto #2
  ICM_5 = 1 !phyto #3
  ICM_6 = 1 !RPOC
  ICM_7 = 1 !LPOC
  ICM_8 = 1 !DOC
  ICM_9 = 1 !RPON
  ICM_10 = 1 !LPON
  ICM_11 = 1 !DON
  ICM_12 = 1 !NH4
  ICM_13 = 1 !NO3
  ICM_14 = 1 !RPOP
  ICM_15 = 1 !LPOP
  ICM_16 = 1 !DOP
  ICM_17 = 1 !PO4t
  ICM_18 = 1 !Si- biogenic
  ICM_19 = 1 !available Si
  ICM_20 = 1 !COD: Chemical oxygen demand
  ICM_21 = 1 !DO
  ICM_Chl= 1 !Chlorophyll
! iPh=1 then following
  ICM_22 = 1 !TIC
  ICM_23 = 1 !ALK
  ICM_24 = 1 !CA
  ICM_25 = 1 !CACO3
  ICM_pH = 1 !PH values

!-----------------------------------------------------------------------
! CoSINE outputs
!-----------------------------------------------------------------------
  COS_1 = 1 !
  COS_2 = 1 !
  COS_3 = 1 !
  COS_4 = 1 !
  COS_5 = 1 !
  COS_6 = 1 !
  COS_7 = 1 !
  COS_8 = 1 !
  COS_9 = 1 !
  COS_10 = 1
  COS_11 = 1
  COS_12 = 1
  COS_13 = 1 !

!-----------------------------------------------------------------------
! Fecal indicating bacteria module
!-----------------------------------------------------------------------
  FIB_1 = 1

!-----------------------------------------------------------------------
! Specific outputs in SED2D (USE_SED2D must be on in Makefile;
! otherwise these are not needed)
!-----------------------------------------------------------------------
  SED2D_depth  = 0 !bottom depth _change_ from init. condition (m)
  SED2D_cdsed  = 0 !drag coefficient used in transport formulae
  SED2D_cflsed = 0 !Courant number (b.qtot.dt / h.dx)
  SED2D_d50    = 0 !Top layer d50 (m)
  SED2D_qtot   = 0 !total transport rate vector (kg/m/s)
  SED2D_qsus   = 0 !suspended tranport rate vector (kg/m/s)
  SED2D_qbdl   = 0 !bedload transport rate vector (kg/m/s)
  SED2D_qav    = 0 !time averaged total transport rate vector (kg/m/s)
  SED2D_dpdxy  = 0 !bottom slope vector (m/m); negative uphill
 
!-----------------------------------------------------------------------
! Ice module outputs (if USE_ICE is on)
!-----------------------------------------------------------------------
  ICE_tracer_1 = 0 !ice volume [m]
  ICE_tracer_2 = 0 !ice concentration [-]
  ICE_tracer_3 = 0 !snow volume [m]
  ICE_velocity = 0 !ice advective velcoity vector [m/s]
  ICE_strain_rate = 0 !strain rate @ elem [1/sec]
  ICE_net_heat_flux = 0 !net heat flux to ocean (>0 warm up SST) [W/m/m]
  ICE_fresh_water = 0 !net fresh water flux to ocean (>0 freshens up SSS) [kg/s/m/m]
  ICE_top_T = 0 !ice temperature [C] at air-ice interface

!-----------------------------------------------------------------------
! Native outputs section. Some of these need corresponding cpp flags
! to be on in order to be active.
!-----------------------------------------------------------------------
  hvel_side = 0 !horizontal vel vector defined at side [m/s]
  vert_elem = 0 !vertical vel. at centroids [m/s]
  temp_elem = 0 !T at prism centers [C]
  salt_elem = 0 !S at prism centers [PSU]
  bthk_elem = 0 ! total bed thickness (m) {module: SED}
  bage_elem = 0 ! total bed age over all layers (sec) {module: SED}
  z0st_elem = 0 ! Sediment transport roughness length (m) {module: SED}
  z0cr_elem = 0 !current-ripples roughness length (m) {module: SED}
  z0sw_elem = 0 !sand-waves roughness length (m) {module: SED}
  z0wr_elem = 0 !wave-ripples roughness length (m) {module: SED}
  SED2D_z0eq_elem = 0 !Total roughness length (m) {module: SED2D}
  SED2D_z0cr_elem = 0 !current-ripples roughness length (m) {module: SED2D}
  SED2D_z0sw_elem = 0 !sand-waves roughness length (m) {module: SED2D}
  SED2D_z0wr_elem = 0 !wave-ripples roughness length (m) {module: SED2D}
  bpgr_side = 1 ! Barotropic pressure gradient force vector (m.s-2) at side centers 
  wave_force_side = 1 ! Wave force vector (m.s-2) computed by wwm at side centers and whole levels {module: WWM}
  ICM_SED_BENDOC_elem = 0 !ICM bed sediment flux arrays: SED_BENDOC (output name:ICM_SED_BENDOC),[gC/(m^2 day)]
  ICM_SED_BENNH4_elem = 0 !ICM bed sediment flux arrays: SED_BENNH4 (output name:ICM_SED_BENNH4),[gC/(m^2 day)]
  ICM_SED_BENNO3_elem = 0 !ICM bed sediment flux arrays: SED_BENNO3 (output name:ICM_SED_BENNO3),[gC/(m^2 day)]
  ICM_SED_BENPO4_elem = 0 !ICM bed sediment flux arrays: SED_BENPO4 (output name:ICM_SED_BENPO4),[gC/(m^2 day)]
  ICM_SED_BENCOD_elem = 0 !ICM bed sediment flux arrays: SED_BENCOD (output name:ICM_SED_BENCOD),[gC/(m^2 day)]
  ICM_SED_BENDO_elem = 0 !ICM bed sediment flux arrays: SED_BENDO (output name:ICM_SED_BENDO),[gC/(m^2 day)]
  ICM_SED_BENSA_elem = 0 !ICM bed sediment flux arrays: SED_BENSA (output name:ICM_SED_BENSA),[gC/(m^2 day)]
  ICM_lfsav = 0 !ICM SAV leaf biomass [gC/m^3] (k=1 is surface)
  ICM_stsav = 0 !ICM SAV stem biomass [gC/m^3]
  ICM_rtsav = 0 !ICM SAV root biomass [gC/m^3]
  ICM_tlfsav = 0 !ICM SAV total leaf biomass [gC/m^2]
  ICM_tstsav = 0 !ICM SAV total stem biomass [gC/m^2]
  ICM_trtsav = 0 !ICM SAV total root biomass [gC/m^2]
  ICM_hcansav = 0 !ICM SAV canopy height [m]

  mrsh_elem = 0 !marsh flags (USE_MARSH on)

!-----------------------------------------------------------------------
! Analysis module outputs
!-----------------------------------------------------------------------
  ANA_air_pres_grad_x = 1 !x-component of \nabla pres /\rho_0 [m/s/s]
  ANA_air_pres_grad_y = 1 !y-component of \nabla pres /\rho_0 [m/s/s]
  ANA_tide_pot_grad_x = 1  !\alpha*g*\nabla \Psi [m/s/s]
  ANA_tide_pot_grad_y = 1  !\alpha*g*\nabla \Psi [m/s/s]
  ANA_hor_viscosity_x = 1 !\nabla \cdot (\mu \nabla u) [m/s/s]
  ANA_hor_viscosity_y = 1 !\nabla \cdot (\mu \nabla u) [m/s/s]
  ANA_bclinic_force_x = 1 !-g/rho0* \int_z^\eta dr_dxy dz  [m/s/s]
  ANA_bclinic_force_y = 1 !-g/rho0* \int_z^\eta dr_dxy dz  [m/s/s]
  ANA_vert_viscosity_x = 1 !d (\nu du/dz)/dz [m/s/s] - no vegetation effects
  ANA_vert_viscosity_y = 1 !d (\nu du/dz)/dz [m/s/s] - no vegetation effects
  ANA_mom_advection_x = 1 !(u \cdot \nabla) u [m/s/s]
  ANA_mom_advection_y = 1 !(u \cdot \nabla) u [m/s/s]
  ANA_Richardson  = 1       !gradient Richardson number [-]
  ANA_transport_min_dt_elem = 1 !min time step at each element over all 
                                !subcycles in horizontal transport solver [s]

!-----------------------------------------------------------------------
! Station output option. If iout_sta/=0, need output skip (nspool_sta) and
! a station.in. If ics=2, the cordinates in station.in must be in lon., lat,
! and z (positive upward; not used for 2D variables). 
!-----------------------------------------------------------------------
  iout_sta = 1
  nspool_sta = 12 !needed if iout_sta/=0

!-----------------------------------------------------------------------
! Flag for harmonic analysis for elevation. If used , need to turn on cpp flags
! in Makefile first. Otherwise set it to 0.
!-----------------------------------------------------------------------
  iharind = 0

!-----------------------------------------------------------------------
! Option for hotstart outputs
!-----------------------------------------------------------------------
  hotout = 1 !1: output *_hotstart every 'hotout_write' steps
  hotout_write = 2880

!-----------------------------------------------------------------------
! Conservation check option. If consv_check=1, some fluxes are computed
! in regions specified in fluxflag.prop (regional number from -1 to an arbitrary integer).
!-----------------------------------------------------------------------
  consv_check = 1


!-----------------------------------------------------------------------
! Test flags for debugging. These flags should be turned off normally.
!-----------------------------------------------------------------------
! Williamson test #5 (zonal flow over an isolated mount); if
! on, ics must =2
!-----------------------------------------------------------------------
  izonal5 = 0 !"0" - no test; otherwise on

!-----------------------------------------------------------------------
! Rotating Gausshill test with stratified T,S (1: on; 0: off)
! Surface T,S read in from *.ic; code generates stratification
!-----------------------------------------------------------------------
  ibtrack_test = 0

!-----------------------------------------------------------------------
! Rouse profile test (1: on; 0: off)
! If on, must turn on USE_TIMOR
!-----------------------------------------------------------------------
  irouse_test = 0

!-----------------------------------------------------------------------
! Flag to choose FIB model for bacteria decay (used with USE_FIB)
! flag_fib = 1 - Constant decay rate (/day) in .gr3 format
!                (kkfib_1.gr3 and kkfib_2.gr3)
! flag_fib = 2 - Decay rate computed from Canteras et al., 1995
! flag_fib = 3 - Decay rate computed from Servais et al., 2007
!----------------------------------------------------------------------
  flag_fib = 1

! hw_depth = 100.
! hw_ratio = 0.5


!-----------------------------------------------------------------------
! Specific heat of water (C_p) in J/kg/K
!-----------------------------------------------------------------------
  shw = 4184.d0


!----------------------------------------------------------------------
! Submerged aquatic vegetation model
! If isav=1, need 3 extra inputs: (1) sav_D.gr3 (depth is stem diameter in meters);
! (2) sav_N.gr3 (depth is # of stems per m^2);
! (3) sav_h.gr3 (height of canopy in meters). 
! If one of these depths=0 at a node, the code will set all to 0. 
!----------------------------------------------------------------------
  isav = 0 !on/off flag
  sav_cd = 1.13 !only needed if isav=1. Drag coefficient


!-----------------------------------------------------------------------
! Reference water density for Boussinesq approximation
!-----------------------------------------------------------------------
  rho0 = 1000.d0 !kg/m^3
