FV3DYCORE  Version1.0.0
fv_arrays_mod::fv_flags_type Type Reference

Public Attributes

character(len=80) grid_name = 'Gnomonic'
 
character(len=120) grid_file = 'Inline'
 
integer grid_type = 0
 -1: read from file; 0: ED Gnomonic 0: the "true" equal-distance Gnomonic grid 1: the traditional equal-distance Gnomonic grid 2: the equal-angular Gnomonic grid 3: the lat-lon grid – to be implemented 4: double periodic boundary condition on Cartesian grid 5: channel flow on Cartesian grid More...
 
integer hord_mt = 9
 Horizontal advection scheme for momentum fluxes. A complete list of kord options is given in the corresponding table in Appendix A of the FV3 technical document. The default value is 9, which uses the third-order piecewise-parabolic method with the monotonicity constraint of Huynh, which is less diffusive than other constraints. For hydrostatic simulation, 8 (the L04 monotonicity constraint) is recommended; for nonhydrostatic simulation, the completely unlimited (“linear” or non-monotone) PPM scheme is recommended. If no monotonicity constraint is applied, enabling the flux damping (do_vort_damp = .true.) is highly recommended to control grid-scale noise. It is also recommended that hord_mt, hord_vt, hord_tm, and hord_dp use the same value, to ensure consistent transport of all dynamical fields, unless a positivity constraint on mass advection (hord_dp) is desired. More...
 
integer kord_mt = 8
 Vertical remapping scheme for the winds. 8 by default; 9 is recommended as the safest option, although 10, and 11 can also be useful. See corresponding table in Appendix A of the FV3 technical document for a complete list of kord options. More...
 
integer kord_wz = 8
 Vertical remapping scheme for vertical velocity in nonhydrostatic simulations. 8 by default; 9 recommended. It is also recommended to use the same value for 'kord_wz' as for 'kord_mt'. More...
 
integer hord_vt = 9
 Vorticity & w transport options: More...
 
integer hord_tm = 9
 Horizontal advection scheme for potential temperature and layer thickness in nonhydrostatic simulations. 9 by default. More...
 
integer hord_dp = 9
 Horizontal advection scheme for mass. A positivity constraint may be warranted for hord_dp but not strictly necessary. 9 by default. More...
 
integer kord_tm = -8
 Vertical remapping scheme for temperature. If positive (not recommended), then vertical remapping is performed on total energy instead of temperature (see 'remap_t'). The default value is -8. More...
 
integer hord_tr = 12
 Tracer transport options: More...
 
integer kord_tr = 8
 The vertical remapping scheme for tracers. The default is 8. 9 or 11 recommended. It is often recommended to use the same value for 'kord_tr' as for 'kord_tm'. More...
 
real scale_z = 0.
 diff_z = scale_z**2 * 0.25 More...
 
real w_max = 75.
 max w (m/s) threshold for hydostatiic adjustment More...
 
real z_min = 0.05
 min ratio of dz_nonhydrostatic/dz_hydrostatic More...
 
real lim_fac = 1.0
 linear scheme limiting factor, 1: hord = 5, 3: hord = 6 More...
 
integer nord = 1
 Order of divergence damping: 0 for second-order; 1 for fourth-order (default); 2 for sixth-order; 3 for eighth-order. Sixth-order may yield a better solution for low resolutions (one degree or coarser) by virtue of it being more scale-selective and will not damp moderately well-resolved disturbances as much as does lower-order damping. More...
 
integer nord_tr = 0
 Order of tracer damping; values mean the same as for 'nord'. The default value is 0. More...
 
real dddmp = 0.0
 Dimensionless coefficient for the second-order Smagorinsky-type divergence damping. The default is value is 0.0. 0.2 (the Smagorinsky constant) is recommended if ICs are noisy. More...
 
real d2_bg = 0.0
 Coefficient for background second-order divergence damping. This option remains active even if nord is nonzero. The default value is 0.0. The proper range is 0 to 0.02. More...
 
real d4_bg = 0.16
 Dimensionless coefficient for background higher-order divergence damping. 0.0 by default. If no second-order divergence damping is used, then values between 0.1 and 0.16 are recommended. Requires 'nord' > 0. Note that the scaling for 'd4_bg' differs from that of 'd2_bg'; 'nord' >= 1 and 'd4_bg' = 0.16 will be less diffusive than 'nord' = 0 and 'd2_bg' = 0.02. More...
 
real vtdm4 = 0.0
 Coefficient for background other-variable damping. The value of 'vtdm4' should be less than that of 'd4_bg'. A good first guess for 'vtdm4' is about one-third the value of d4_bg. Requires 'do_vort_damp' to be .true. Disabled for values less than 1.e-3. Other- variable damping should not be used if a monotonic horizontal advection scheme is used. The default value is 0.0. More...
 
real trdm2 = 0.0
 Coefficient for del-2 tracer damping. More...
 
real d2_bg_k1 = 4.
 Strength of second-order diffusion in the top sponge layer. Value must be specified. This value, and d2_bg_k2, will be changed appropriately in the model (depending on the height of model top), so the actual damping may be very reduced. See atmos_cubed_sphere/model/dyncore.F90 for details. Recommended range is 0. to 0.2. Note that since diffusion is converted to heat if d_con > 0 larger amounts of sponge-layer diffusion may be less stable. More...
 
real d2_bg_k2 = 2.
 Strength of second-order diffusion in the second sponge layer from the model top. This value must be specified, and should be less than 'd2_bg_k1'. More...
 
real d2_divg_max_k1 = 0.15
 d2_divg max value (k=1) More...
 
real d2_divg_max_k2 = 0.08
 d2_divg max value (k=2) More...
 
real damp_k_k1 = 0.2
 damp_k value (k=1) More...
 
real damp_k_k2 = 0.12
 damp_k value (k=2) More...
 
integer n_zs_filter = 0
 Additional (after the fact) terrain filter (to further smooth the terrain after cold start) More...
 
integer nord_zs_filter = 4
 Order of the topography filter applied to n_zs_filter. Set to 2 to get a second-order filter, or 4 to get a fourth-order filter; other values do no filtering. 0 by default. This should not be set to a non-zero value on multiple successive simulations; the filter is applied every time the model restarts. This option is useful for testing the terrain filter, and SHOULD NOT BE USED FOR REGULAR RUNS.
use del-2 (2) OR del-4 (4) More...
 
logical full_zs_filter = .false.
 Whether to apply the on-line topography filter during initialization. Only active if get_nggps_ic = .true. This is so topography filtering can be performed on the initial conditions output by the pre-processing tools, which currently do not support topography filter- ing for some configurations (such as the nested grid); this also allows the user to easily test changes to the topography filtering on the simulation. Note that for all other initialization methods (if external_ic = .true.) the on-line topography filter will be applied automatically during the initialization of the topography. The default value is .false. More...
 
logical rf_fast = .false.
 Option controlling whether to apply Rayleigh damping (for tau > 0) on the dynamic/acoustic timestep rather than on the physics timestep. This can help stabilize the model by applying the damping more weakly more frequently, so the instantaneous amount of damping (and thereby heat added) is reduced. The default is .false., which applies the Rayleigh drag every physics timestep. More...
 
logical consv_am = .false.
 Whether to enable Angular momentum fixer. The default is .false. More...
 
logical do_sat_adj = .false.
 
logical do_f3d = .false.
 
logical no_dycore = .false.
 Disables execution of the dynamical core, only running the initialization, diagnostic, and I/O routines, and any physics that may be enabled. Essentially turns the model into a column physics model. The default is .false. More...
 
logical convert_ke = .false.
 If .true., adds energy dissipated through mechanical damping to heat throughout the entire depth of the domain; if .false. (default) this is only done in the sponge layer at the top of the domain. This option is only enabled if d_con > 1.e-5. More...
 
logical do_vort_damp = .false.
 Whether to apply flux damping (of strength governed by 'vtdm4') to the fluxes of vorticity, air mass, and nonhydrostatic vertical velocity (there is no dynamically correct way to add explicit diffusion to the tracer fluxes). The form is the same as is used for the divergence damping, including the same order (from 'nord') damping, unless 'nord' = 0, in which case this damping is fourth-order, or if 'nord' = 3,in which case this damping is sixth-order (instead of eighth-order). We recommend enabling this damping when the linear or non-monotonic horizontal advection schemes are enabled, but is unnecessary and not recommended when using monotonic advection. The default is .false. More...
 
logical use_old_omega = .true.
 
real beta = 0.0
 Parameter specifying fraction of time-off-centering for backwards evaluation of the pressure gradient force. The default is 0.0, which produces a fully backwards evaluation of the pressure gradient force that is entirely evaluated using the updated (time n+1) dynamical fields. A value of 0.5 will equally weight the PGF determined at times n and n+1, but may not be stable; values larger than 0.45 are not recommended. A value of 0.4 is recommended for most hydrostatic simulations, which allows an improved representation of inertia-gravity waves in the tropics. In non-hydrostatic simulations using the semi-implicit solver (a_imp > 0.5) the values of 'a_imp' and 'beta' should add to 1, so that the time-centering is consistent between the PGF and the nonhydrostatic solver. The proper range is 0 to 0.45. More...
 
integer n_sponge = 1
 Controls the number of layers at the upper boundary on which the 2Dx filter is applied. This does not control the sponge layer. The default value is 0. More...
 
real d_ext = 0.02
 Coefficient for external (barotropic) mode damping. Proper range is 0 to 0.02. A value of 0.01 or 0.02 may help improve the models maximum stable time step in low-resolution (2-degree or lower) simulations; otherwise a value of 0 is recommended. The default value is 0.02. More...
 
integer nwat = 3
 Number of water species to be included in condensate and water vapor loading. The masses of the first nwat tracer species will be added to the dry air mass, so that p is the mass of dry air, water vapor, and the included condensate species. The value used depends on the microphysics in the physics package you are using. For GFS physics with only a single condensate species, set to 2. For schemes with prognostic cloud water and cloud ice, such as GFDL AM2/AM3/AM4 Rotsteyn-Klein or Morrison-Gettlean microphysics, set to 3. For warm-rain (Kessler) microphysics set to 4 (with an inactive ice tracer), which only handles three species but uses 4 to avoid interference with the R-K physics. For schemes such as WSM5 or Ferrier that have prognostic rain and snow but not hail, set to 5 (not yet implemented). For six-category schemes that also have prognostic hail or graupel, such as the GFDL, Thompson, or WSM6 microphysics, set to 6. A value of 0 turns off condensate loading. The default value is 3. More...
 
logical warm_start = .true.
 Whether to start from restart files, instead of cold-starting the model. True by default; if this is set to .true. and restart files cannot be found the model will stop. More...
 
logical inline_q = .false.
 Whether to compute tracer transport in-line with the rest of the dynamics instead of sub-cycling, so that tracer transport is done at the same time and on the same time step as is p and potential temperature. False by default; if true, q_split and z_tracer are ignored. More...
 
logical adiabatic = .false.
 Run without physics (full or idealized). More...
 
real shift_fac = 18.
 Westward zonal rotation (or shift) of cubed-sphere grid from its natural orientation with cube face centers at 0, 90, 180, and 270 degrees longitude. The shift, in degrees, is 180/shift_fac. This shift does not move the poles. By default this is set to 18, shifting the grid westward 180/18=10 degrees, so that the edges of the cube do not run through the mountains of Japan; all standard CM2.x, AM3, CM3, and HiRAM simulations use this orientation of the grid. Requires do_schmidt = .false. More...
 
logical do_schmidt = .false.
 Whether to enable grid stretching and rotation using stretch_fac, target_lat, and target_lon. The default value is .false. More...
 
real(kind=r_gridstretch_fac = 1.
 Stretching factor for the Schmidt transformation. This is the factor by which tile 6 of the cubed sphere will be shrunk, with the grid size shrinking accordingly. The default value is 1, which performs no grid stretching. Requires do_schmidt =.true. THE MODEL WILL CRASH IF stretch_fac IS SET TO ZERO. Values of up to 40 have been found useful and stable for short-term cloud-scale integrations. More...
 
real(kind=r_gridtarget_lat = -90.
 Latitude (in degrees) to which the center of tile 6 will be rotated; if stretching is done with stretch_fac the center of the high-resolution part of the grid will be at this latitude. -90 by default, which does no grid rotation (the Schmidt transformation rotates the south pole to the appropriate target). Requires do_schmidt = .true. More...
 
real(kind=r_gridtarget_lon = 0.
 Longitude to which the center of tile 6 will be rotated. 0 by default. Requires do_schmidt = .true. More...
 
logical reset_eta = .false.
 
real p_fac = 0.05
 Safety factor for minimum nonhydrostatic pressures, which will be limited so the full pressure is no less than p_fac times the hydrostatic pressure. This is only of concern in mid-top or high-top models with very low pressures near the model top, and has no effect in most simulations. The pressure limiting activates only when model is in danger of blowup due to unphysical negative total pressures. Only used if 'hydrostatic' = .false.and the semi-implicit solver is used. The proper range is 0 to 0.25. The default value is 0.05. More...
 
real a_imp = 0.75
 Controls behavior of the non-hydrostatic solver. Values > 0.5 enable the semi-implicit solver, in which the value of 'a_imp' controls the time-off-centering: use a_imp = 1.0 for a fully backward time stepping. For consistency, the sum of 'beta' and 'a_imp' should be 1 when the semi-implicit solver is used. The semi-implicit algorithm is substantially more efficient except at very high (km-scale) resolutions with an acoustic time step of a few seconds or less. Proper values are 0, or between 0.5
and 1. The default value is 0.75. Only used if 'hydrostatic' = .false. More...
 
integer n_split = 0
 The number of small dynamics (acoustic) time steps between vertical remapping. 0 by default, in which case the model produces a good first guess by examining the resolution, dt_atmos, and k_split. More...
 
real fac_n_spl = 1.0
 factor multiplying n_split up tp forecast hour fhouri More...
 
real fhouri = 0.0
 forecast hour up to which the number of small dynamics (acoustic) time steps are nint(n_split*fac_n_spl) More...
 
integer m_split = 0
 Number of time splits for Riemann solver. More...
 
integer k_split = 1
 Number of vertical remappings per dt_atmos (physics timestep). 1 by default. More...
 
logical use_logp = .false.
 Enables a variant of the Lin pressure-gradient force algorithm, which uses the logarithm of pressure instead of the Exner function (as in [lin1997explicit]). This yields more accurate results for regions that are nearly isothermal. Ignored if 'hydrostatic' = .true. The default is .false. More...
 
integer q_split = 0
 number of time steps for sub-cycled tracer advection. The default value is 0 (recommended), in which case the model determines the number of time steps from the global maximum wind speed at each call to the tracer advection. More...
 
integer print_freq = 0
 number of hours between print out of max/min and air/tracer mass diagnostics to standard output. 0 by default, which never prints out any output; set to -1 to see output after every dt_at-mos. Computing these diagnostics requires some computational overhead More...
 
logical write_3d_diags = .true.
 whether to write out three-dimensional dynamical diagnostic fields (those defined in fv_diagnostics.F90). This is useful for runs with multiple grids if you only want very large 3D diagnostics written out for (say) a nested grid, and not for the global grid. False by default. More...
 
integer npx
 Number of grid corners in the x-direction on one tile of the domain; so one more than the number of grid cells across a tile. On the cubed sphere this is one more than the number of cells across a cube face. Must be set. More...
 
integer npy
 Number of grid corners in the y-direction on one tile of the domain. This value should be identical to npx on a cubed-sphere grid; doubly periodic or nested grids do not have this restriction. Must be set. More...
 
integer npz
 Number of vertical levels. Each choice of npz comes with a pre-defined set of hybrid sigma-pressure levels and model top (see fv_eta.F90). Must be set. More...
 
integer npz_rst = 0
 If using a restart file with a different number of vertical levels, set npz_rst to be the number of levels in your restart file. The model will then remap the restart file data to the vertical coordinates specified by npz. 0 by default; if 0 or equal to npz no remapping is done. More...
 
integer ncnst = 0
 Number of tracer species advected by fv_tracer in the dynamical core.
Typically this is set automatically by reading in values from field_table, but ncnst can be set to a smaller value so only the first ncnst tracers listed in field_table are not advected. 0 by default, which will use the value from field_table. More...
 
integer pnats = 0
 The number of tracers not to advect by the dynamical core. Unlike dnats, these tracers are not seen by the dynamical core. The last pnats entries in field_table are not advected. The default value is 0. More...
 
integer dnats = 0
 The number of tracers which are not to be advected by the dynamical core, but still passed into the dynamical core; the last dnats+pnats tracers in field_table are not advected. 0 by default. More...
 
integer ntiles = 1
 Number of tiles on the domain. For the cubed sphere, this should be 6, one tile for each face of the cubed sphere; normally for most other domains (including nested grids) this should be set to 1. Must be set. More...
 
integer ndims = 2
 
integer nf_omega = 1
 Number of times to apply second-order smoothing to the diagnosed omega. When 0 the filter is disabled. 1 by default. More...
 
integer fv_sg_adj = -1
 Timescale (in seconds) at which to remove two-delta-z instability when the local (between two adjacent levels) Richardson number is less than 1. This is achieved by local mixing, which conserves mass, momentum, and total energy. Values of 0 or smaller disable this feature. If n_sponge < 0 then the mixing is applied only to the top n_sponge layers of the domain. Set to -1 (inactive) by default. The proper range is 0 to 3600. More...
 
integer na_init = 0
 Number of forward-backward dynamics steps used to initialize adiabatic solver. This is useful for spinning up the nonhydrostatic state from the hydrostatic GFS analyses. 0 by default. Recommended to set this to a non-zero value (1 or 2 is typically sufficient) when initializing from GFS or ECMWF analyses. More...
 
logical nudge_dz = .false.
 During the adiabatic initialization (na_init > 0), if set to .true., delz is nudged back to the value specified in the initial conditions, instead of nudging the temperature back to the initial value. Nudging delz is simpler (faster), doesn’t require consideration of the virtual temperature effect, and may be more stable. .false.by default. More...
 
real p_ref = 1.E5
 Surface pressure used to construct a horizontally-uniform reference vertical pressure profile, used in some simple physics packages in the solo_core and in the Rayleigh damping. This should not be confused with the actual, horizontally-varying pressure levels used for all other dynamical calculations. The default value is 1.e5. CHANGING THIS VALUE IS STRONGLY DISCOURAGED. More...
 
real dry_mass = 98290.
 If adjust_dry_mass is .true., sets the global dry air mass, measured in the globally-averaged surface pressure (Pascals) by adding or removing mass from the lowest layer of the atmosphere as needed. The default value is 98290. (Pa). More...
 
integer nt_prog = 0
 
integer nt_phys = 0
 
real tau_h2o = 0.
 Time-scale (days) for simple methane chemistry to act as a source of water in the stratosphere. Can be useful if the stratosphere dries out too quickly; consider a value between 60 and 120 days if this is the case. The default value is 0., which disables the methane chemistry. Values less than zero apply the chemistry above 100 mb; else applied above 30 mb. Requires 'adiabatic' to be .false. More...
 
real delt_max = 1.
 Maximum allowed magnitude of the dissipative heating rate, K s−1; larger magnitudes are clipped to this amount. This can help avoid instability that can occur due to strong heating when d_con > 0. A value of 0.008 (a rate equivalent to about 800 K/day) is sufficient to stabilize the model at 3-km resolution. Set to 1. by default, which effectively disables this limitation. More...
 
real d_con = 0.
 Fraction of kinetic energy lost to explicit damping to be converted to heat. Acts as a dissipative heating mechanism in the dynamical core. The default is 0. Proper range is 0 to 1. Note that this is a local, physically correct, energy fixer. More...
 
real ke_bg = 0.
 background KE production (m^2/s^3) over a small step Use this to conserve total energy if consv_te=0 More...
 
real consv_te = 0.
 Fraction of total energy lost during the adiabatic integration between calls of the physics, to be added back globally as heat; essentially the strength of the energy fixer in the physics. Note that this is a global energy fixer and cannot add back energy locally. The default algorithm increments the potential temperature so the pressure gradients are unchanged. The default value is 0. Proper range is 0 to 1. 1 will restore the energy completely to its original value before entering the physics; a value of 0.7 roughly causes the energy fixer to compensate for the amount of energy changed by the physics in GFDL HiRAM or AM3. More...
 
real tau = 0.
 Time scale (in days) for Rayleigh friction applied to horizontal and vertical winds; lost kinetic energy is converted to heat, except on nested grids. The default value is 0.0, which disables damping. Larger values yield less damping. For models with tops at 1 mb or lower values between 10 and 30 are useful for preventing overly-strong polar night jets; for higher-top hydrostatic models values between 5 and 15 should be considered; and for non-hydrostatic models values of 10 or less should be
considered, with smaller values for higher-resolution. More...
 
real rf_cutoff = 30.E2
 Pressure below which no Rayleigh damping is applied if tau > 0. More...
 
logical filter_phys = .false.
 
logical dwind_2d = .false.
 Whether to use a simpler & faster algorithm for interpolating the A-grid (cell-centered) wind tendencies computed from the physics to the D-grid. Typically, the A-grid wind tendencies are first converted in 3D cartesian coordinates and then interpolated before converting back to 2D local coordinates. When this option enabled, a much simpler but less accurate 2D interpolation is used. False by default. More...
 
logical breed_vortex_inline = .false.
 Whether to bogus tropical cyclones into the model, which are specified by an external file. Options are set in fv_nwp_nudge_nml. False by default. More...
 
logical range_warn = .false.
 Checks whether the values of the prognostic variables are within a reasonable range at the end of a dynamics time step, and prints a warning if not. The default is .false.; adds computational, overhead so we only recommend using this when debugging. More...
 
logical fill = .false.
 Fills in negative tracer values by taking positive tracers from the cells above and below. This option is useful when the physical parameterizations produced negatives. The default is .false. More...
 
logical fill_dp = .false.
 Like 'fill' except for p, the hydrostatic pressure thickness. When the filling occurs a diagnostic message is printed out, which is helpful for diagnosing where the problem may be occurring. Typically, a crash is inevitable if the pressure filling is needed; thus, this option is often better for debugging than as a safety valve. The default is .false. More...
 
logical fill_wz = .false.
 
logical check_negative = .false.
 Whether to print the most negative global value of microphysical tracers. More...
 
logical non_ortho = .true.
 
logical moist_phys = .true.
 Run with moist physics. More...
 
logical do_held_suarez = .false.
 Whether to use Held-Suarez forcing. Requires adiabatic to be false. The default is .false.; this option has no effect if not running solo_core. More...
 
logical do_reed_physics = .false.
 
logical reed_cond_only = .false.
 
logical reproduce_sum = .true.
 uses an exactly-reproducible global sum operation performed when computing the global energy for consv_te. This is used because the FMS routine mpp_sum() is not bit-wise reproducible due to its handling of floating-point arithmetic, and so can return different answers for (say) different processor layouts. The default is .true. More...
 
logical adjust_dry_mass = .false.
 Whether to adjust the global dry-air mass to the value set by dry_mass. This is only done in an initialization step, particularly when using an initial condition from an external dataset, interpolated from another resolution (either horizontal or vertical), or when changing the topography, so that the global mass of the atmosphere matches some estimate of observed value. False by default. It is recommended to only set this to .true. when initializing the model. More...
 
logical fv_debug = .false.
 Whether to turn on additional diagnostics in fv_dynamics. The default is .false. More...
 
logical srf_init = .false.
 
logical mountain = .true.
 Takes topography into account when initializing the model. Set this to .true. to apply the terrain filter (if n_zs_filter = 2 or 4) upon startup; also set to True when cold starting so that the topography can be initialized. Only set this to .false. if you wish to cold-start without any topography; this value is ignored for the aquaplanet test_case = 14. The default is .true. It is highly recommended TO NOT ALTER this value unless you know what you are doing. More...
 
logical old_divg_damp = .false.
 parameter to revert damping parameters back to values defined in a previous revision old_values: d2_bg_k1 = 6. d2_bg_k2 = 4. d2_divg_max_k1 = 0.02 d2_divg_max_k2 = 0.01 damp_k_k1 = 0. damp_k_k2 = 0. current_values: d2_bg_k1 = 4. d2_bg_k2 = 2. d2_divg_max_k1 = 0.15 d2_divg_max_k2 = 0.08 damp_k_k1 = 0.2 damp_k_k2 = 0.12 More...
 
logical remap_t = .true.
 Whether the vertical remapping is performed on (virtual) temperature instead of (virtual) potential temperature. Since typically potential temperature increases exponentially from layer to layer near the top boundary, the cubic-spline interpolation in the vertical remapping will have difficulty with the exponential profile. Temperature does not have this problem and will often yield a more accurate result. The default is .true. More...
 
logical z_tracer = .false.
 Whether to transport sub-cycled tracers layer-by-layer, each with its own computed sub-cycling time step (if q_split = 0). This may improve efficiency for very large numbers of tracers. The default value is .false.; currently not implemented. More...
 
logical fv_land = .false.
 Whether to create terrain deviation and land fraction for output to mg_drag restart files, for use in mg_drag and in the land model. The default is .false; .true. is recommended when, and only when, initializing the model, since the mg_drag files created provide a much more accurate terrain representation for the mountain gravity wave drag parameterization and for the land surface roughness than either computes internally. This has no effect on the representation of the terrain in the dynamics. More...
 
logical nudge = .false.
 Whether to use the nudging towards the state in some externally-supplied file (such as from reanalysis or another simulation). Further nudging options are set in fv_nwp_nudge_nml. The default is .false. More...
 
logical nudge_ic = .false.
 Same as nudge, but works in adiabatic solo_core simulations to nudge the field to a single external analysis file. The default is .false. More...
 
logical ncep_ic = .false.
 If external_ic = .true., this variable says whether the file in res_latlon_dynamics is an NCEP analysis or reanalysis file. This option zeros out all tracer fields except specific humidity. The default is .false. More...
 
logical nggps_ic = .false.
 If external_ic = .true., reads initial conditions from horizontally-interpolated output from chgres. The default is .false. Additional options are available through external_ic_nml. More...
 
logical ecmwf_ic = .false.
 If external_ic = .true., reads initial conditions from ECMWF analyses. The default is .false. More...
 
logical gfs_phil = .false.
 if .T., compute geopotential inside of GFS physics More...
 
logical agrid_vel_rst = .false.
 Whether to write the unstaggered latitude-longitude winds (ua and va) to the restart files. This is useful for data assimilation cycling systems which do not handle staggered winds. The default is .false. More...
 
logical use_new_ncep = .false.
 
logical use_ncep_phy = .false.
 
logical fv_diag_ic = .false.
 
logical external_ic = .false.
 Whether to initialize the models state using the data in an externally specified file, given in res_latlon_dynamics. By default this file is assumed to be a legacy lat-lon FV core restart file; set either ncep_ic or fv_diag_ic to .true.to override this behavior. The default is .false. Note that external_ic = .true. will cause the model to re-initialize the dynamical fields from the input dataset regardless of whether warm_start is set. More...
 
logical external_eta = .false.
 If .true., reads the interface coefficients ak and bk from either the restart file (if restarting) or from the external initial condition file (if nggps_ic or ecwmf_ic are .true.). This overrides the hard-coded levels in fv_eta. The default is .false. More...
 
logical read_increment = .false.
 read in analysis increment and add to restart More...
 
logical do_skeb = .false.
 save dissipation estimate More...
 
integer skeb_npass = 11
 Filter dissipation estimate "skeb_npass" times. More...
 
character(len=128) res_latlon_dynamics = 'INPUT/fv_rst.res.nc'
 If external_ic =.true.gives the filename of the input IC file. The default is 'INPUT/fv_rst.res.nc'. More...
 
character(len=128) res_latlon_tracers = 'INPUT/atmos_tracers.res.nc'
 If external_ic =.true.and both ncep_ic and fv_diag_ic are.false., this variable gives the filename of the initial conditions for the tracers, assumed to be a legacy lat-lon FV core restart file. The default is 'INPUT/atmos_tracers.res.nc'. More...
 
logical hydrostatic = .true.
 Whether to use the hydrostatic or nonhydrostatic solver. The default is .true. More...
 
logical phys_hydrostatic = .true.
 Option to enable hydrostatic application of heating from the physics in a nonhydrostatic simulation: heating is applied in hydrostatic balance, causing the entire atmospheric column to expand instantaneously. If .false., heating from the physics is applied simply as a temperature tendency. The default value is .true.; ignored if hydrostatic = .true. More...
 
logical use_hydro_pressure = .false.
 Whether to compute hydrostatic pressure for input to the physics. Currently only enabled for the fvGFS model. Ignored in hydrostatic simulations. The default is .false. More...
 
logical do_uni_zfull = .false.
 Whether to compute z_full (the height of each modellayer, as opposed to z_half, the height of each model interface) as the midpoint of the layer, as is done for the nonhydrostatic solver, instead of the height of the location where p = p the mean pressure in the layer. This option is not available for fvGFS or the solo_core. The default is .false. More...
 
logical hybrid_z = .false.
 Whether to use a hybrid-height coordinate, instead of the usual sigma-p coordinate. The default value is .false. (Not currently maintained.) More...
 
logical make_nh = .false.
 Whether to re-initialize the nonhydrostatic state, by recomputing dz from hydrostatic balance and setting w to 0. The default is false. More...
 
logical make_hybrid_z = .false.
 Converts the vertical coordinate to a hybrid-height coordinate, instead of the usual sigma-p coordinate. Requires hybrid_z = .true. The default value is .false. More...
 
logical nudge_qv = .false.
 During the adiabatic initialization (na_init > 0), if set to .true., the water vapor profile is nudged to an analytic fit to the HALOE climatology. This is to improve the water vapor concentrations in GFS initial conditions, especially in the stratosphere, where values can be several times higher than observed. This nudging is unnecessary for other ICs, especially the ECMWF initial conditions. The default is .false. More...
 
real add_noise = -1.
 Amplitude of random thermal noise (in K) to add upon startup. Useful for perturbing initial conditions. -1 by default; disabled if 0 or negative. More...
 
integer a2b_ord = 4
 Order of interpolation used by the pressure gradient force to interpolate cell-centered (A-grid) values to the grid corners. The default value is 4 (recommended), which uses fourth-order interpolation; otherwise second-order interpolation is used. More...
 
integer c2l_ord = 4
 Order of interpolation from the solvers native D-grid winds to latitude-longitude A-grid winds, which are used as input to the physics routines and for writing to history files. The default value is 4 (recommended); fourth-order interpolation is used unless c2l_ord = 2. More...
 
real(kind=r_griddx_const = 1000.
 Specifies the (uniform) grid-cell-width in the x-direction on a doubly-periodic grid (grid_type = 4) in meters. The default value is 1000. More...
 
real(kind=r_griddy_const = 1000.
 Specifies the (uniform) grid-cell-width in the y-direction on a doubly-periodic grid (grid_type = 4) in meters. The default value is 1000. More...
 
real(kind=r_griddeglat = 15.
 Latitude (in degrees) used to compute the uniform f-plane Coriolis parameter for doubly-periodic simulations (grid_type = 4). The default value is 15. More...
 
real(kind=r_griddeglon_start = -30.
 
real(kind=r_griddeglon_stop = 30.
 boundaries of latlon patch More...
 
real(kind=r_griddeglat_start = -30.
 
real(kind=r_griddeglat_stop = 30.
 
logical regional = .false.
 Default setting for the regional domain. More...
 
integer bc_update_interval = 3
 Default setting for interval (hours) between external regional BC data files. More...
 
integer, pointer grid_number
 Convenience pointers. More...
 
logical adj_mass_vmr = .false.
 

Detailed Description

Definition at line 245 of file fv_arrays.F90.

Member Data Documentation

◆ a2b_ord

integer fv_arrays_mod::fv_flags_type::a2b_ord = 4

Order of interpolation used by the pressure gradient force to interpolate cell-centered (A-grid) values to the grid corners. The default value is 4 (recommended), which uses fourth-order interpolation; otherwise second-order interpolation is used.

Definition at line 977 of file fv_arrays.F90.

◆ a_imp

real fv_arrays_mod::fv_flags_type::a_imp = 0.75

Controls behavior of the non-hydrostatic solver. Values > 0.5 enable the semi-implicit solver, in which the value of 'a_imp' controls the time-off-centering: use a_imp = 1.0 for a fully backward time stepping. For consistency, the sum of 'beta' and 'a_imp' should be 1 when the semi-implicit solver is used. The semi-implicit algorithm is substantially more efficient except at very high (km-scale) resolutions with an acoustic time step of a few seconds or less. Proper values are 0, or between 0.5
and 1. The default value is 0.75. Only used if 'hydrostatic' = .false.

Definition at line 589 of file fv_arrays.F90.

◆ add_noise

real fv_arrays_mod::fv_flags_type::add_noise = -1.

Amplitude of random thermal noise (in K) to add upon startup. Useful for perturbing initial conditions. -1 by default; disabled if 0 or negative.

Definition at line 973 of file fv_arrays.F90.

◆ adiabatic

logical fv_arrays_mod::fv_flags_type::adiabatic = .false.

Run without physics (full or idealized).

Definition at line 531 of file fv_arrays.F90.

◆ adj_mass_vmr

logical fv_arrays_mod::fv_flags_type::adj_mass_vmr = .false.

Definition at line 1013 of file fv_arrays.F90.

◆ adjust_dry_mass

logical fv_arrays_mod::fv_flags_type::adjust_dry_mass = .false.

Whether to adjust the global dry-air mass to the value set by dry_mass. This is only done in an initialization step, particularly when using an initial condition from an external dataset, interpolated from another resolution (either horizontal or vertical), or when changing the topography, so that the global mass of the atmosphere matches some estimate of observed value. False by default. It is recommended to only set this to .true. when initializing the model.

Definition at line 819 of file fv_arrays.F90.

◆ agrid_vel_rst

logical fv_arrays_mod::fv_flags_type::agrid_vel_rst = .false.

Whether to write the unstaggered latitude-longitude winds (ua and va) to the restart files. This is useful for data assimilation cycling systems which do not handle staggered winds. The default is .false.

Definition at line 894 of file fv_arrays.F90.

◆ bc_update_interval

integer fv_arrays_mod::fv_flags_type::bc_update_interval = 3

Default setting for interval (hours) between external regional BC data files.

Definition at line 1007 of file fv_arrays.F90.

◆ beta

real fv_arrays_mod::fv_flags_type::beta = 0.0

Parameter specifying fraction of time-off-centering for backwards evaluation of the pressure gradient force. The default is 0.0, which produces a fully backwards evaluation of the pressure gradient force that is entirely evaluated using the updated (time n+1) dynamical fields. A value of 0.5 will equally weight the PGF determined at times n and n+1, but may not be stable; values larger than 0.45 are not recommended. A value of 0.4 is recommended for most hydrostatic simulations, which allows an improved representation of inertia-gravity waves in the tropics. In non-hydrostatic simulations using the semi-implicit solver (a_imp > 0.5) the values of 'a_imp' and 'beta' should add to 1, so that the time-centering is consistent between the PGF and the nonhydrostatic solver. The proper range is 0 to 0.45.

Definition at line 446 of file fv_arrays.F90.

◆ breed_vortex_inline

logical fv_arrays_mod::fv_flags_type::breed_vortex_inline = .false.

Whether to bogus tropical cyclones into the model, which are specified by an external file. Options are set in fv_nwp_nudge_nml. False by default.

Definition at line 782 of file fv_arrays.F90.

◆ c2l_ord

integer fv_arrays_mod::fv_flags_type::c2l_ord = 4

Order of interpolation from the solvers native D-grid winds to latitude-longitude A-grid winds, which are used as input to the physics routines and for writing to history files. The default value is 4 (recommended); fourth-order interpolation is used unless c2l_ord = 2.

Definition at line 982 of file fv_arrays.F90.

◆ check_negative

logical fv_arrays_mod::fv_flags_type::check_negative = .false.

Whether to print the most negative global value of microphysical tracers.

Definition at line 804 of file fv_arrays.F90.

◆ consv_am

logical fv_arrays_mod::fv_flags_type::consv_am = .false.

Whether to enable Angular momentum fixer. The default is .false.

Definition at line 416 of file fv_arrays.F90.

◆ consv_te

real fv_arrays_mod::fv_flags_type::consv_te = 0.

Fraction of total energy lost during the adiabatic integration between calls of the physics, to be added back globally as heat; essentially the strength of the energy fixer in the physics. Note that this is a global energy fixer and cannot add back energy locally. The default algorithm increments the potential temperature so the pressure gradients are unchanged. The default value is 0. Proper range is 0 to 1. 1 will restore the energy completely to its original value before entering the physics; a value of 0.7 roughly causes the energy fixer to compensate for the amount of energy changed by the physics in GFDL HiRAM or AM3.

Definition at line 751 of file fv_arrays.F90.

◆ convert_ke

logical fv_arrays_mod::fv_flags_type::convert_ke = .false.

If .true., adds energy dissipated through mechanical damping to heat throughout the entire depth of the domain; if .false. (default) this is only done in the sponge layer at the top of the domain. This option is only enabled if d_con > 1.e-5.

Definition at line 426 of file fv_arrays.F90.

◆ d2_bg

real fv_arrays_mod::fv_flags_type::d2_bg = 0.0

Coefficient for background second-order divergence damping. This option remains active even if nord is nonzero. The default value is 0.0. The proper range is 0 to 0.02.

Definition at line 346 of file fv_arrays.F90.

◆ d2_bg_k1

real fv_arrays_mod::fv_flags_type::d2_bg_k1 = 4.

Strength of second-order diffusion in the top sponge layer. Value must be specified. This value, and d2_bg_k2, will be changed appropriately in the model (depending on the height of model top), so the actual damping may be very reduced. See atmos_cubed_sphere/model/dyncore.F90 for details. Recommended range is 0. to 0.2. Note that since diffusion is converted to heat if d_con > 0 larger amounts of sponge-layer diffusion may be less stable.

Definition at line 364 of file fv_arrays.F90.

◆ d2_bg_k2

real fv_arrays_mod::fv_flags_type::d2_bg_k2 = 2.

Strength of second-order diffusion in the second sponge layer from the model top. This value must be specified, and should be less than 'd2_bg_k1'.

Definition at line 373 of file fv_arrays.F90.

◆ d2_divg_max_k1

real fv_arrays_mod::fv_flags_type::d2_divg_max_k1 = 0.15

d2_divg max value (k=1)

Definition at line 377 of file fv_arrays.F90.

◆ d2_divg_max_k2

real fv_arrays_mod::fv_flags_type::d2_divg_max_k2 = 0.08

d2_divg max value (k=2)

Definition at line 378 of file fv_arrays.F90.

◆ d4_bg

real fv_arrays_mod::fv_flags_type::d4_bg = 0.16

Dimensionless coefficient for background higher-order divergence damping. 0.0 by default. If no second-order divergence damping is used, then values between 0.1 and 0.16 are recommended. Requires 'nord' > 0. Note that the scaling for 'd4_bg' differs from that of 'd2_bg'; 'nord' >= 1 and 'd4_bg' = 0.16 will be less diffusive than 'nord' = 0 and 'd2_bg' = 0.02.

Definition at line 350 of file fv_arrays.F90.

◆ d_con

real fv_arrays_mod::fv_flags_type::d_con = 0.

Fraction of kinetic energy lost to explicit damping to be converted to heat. Acts as a dissipative heating mechanism in the dynamical core. The default is 0. Proper range is 0 to 1. Note that this is a local, physically correct, energy fixer.

Definition at line 743 of file fv_arrays.F90.

◆ d_ext

real fv_arrays_mod::fv_flags_type::d_ext = 0.02

Coefficient for external (barotropic) mode damping. Proper range is 0 to 0.02. A value of 0.01 or 0.02 may help improve the models maximum stable time step in low-resolution (2-degree or lower) simulations; otherwise a value of 0 is recommended. The default value is 0.02.

Definition at line 502 of file fv_arrays.F90.

◆ damp_k_k1

real fv_arrays_mod::fv_flags_type::damp_k_k1 = 0.2

damp_k value (k=1)

Definition at line 379 of file fv_arrays.F90.

◆ damp_k_k2

real fv_arrays_mod::fv_flags_type::damp_k_k2 = 0.12

damp_k value (k=2)

Definition at line 380 of file fv_arrays.F90.

◆ dddmp

real fv_arrays_mod::fv_flags_type::dddmp = 0.0

Dimensionless coefficient for the second-order Smagorinsky-type divergence damping. The default is value is 0.0. 0.2 (the Smagorinsky constant) is recommended if ICs are noisy.

Definition at line 342 of file fv_arrays.F90.

◆ deglat

real(kind=r_grid) fv_arrays_mod::fv_flags_type::deglat = 15.

Latitude (in degrees) used to compute the uniform f-plane Coriolis parameter for doubly-periodic simulations (grid_type = 4). The default value is 15.

Definition at line 997 of file fv_arrays.F90.

◆ deglat_start

real(kind=r_grid) fv_arrays_mod::fv_flags_type::deglat_start = -30.

Definition at line 1002 of file fv_arrays.F90.

◆ deglat_stop

real(kind=r_grid) fv_arrays_mod::fv_flags_type::deglat_stop = 30.

Definition at line 1002 of file fv_arrays.F90.

◆ deglon_start

real(kind=r_grid) fv_arrays_mod::fv_flags_type::deglon_start = -30.

Definition at line 1002 of file fv_arrays.F90.

◆ deglon_stop

real(kind=r_grid) fv_arrays_mod::fv_flags_type::deglon_stop = 30.

boundaries of latlon patch

Definition at line 1002 of file fv_arrays.F90.

◆ delt_max

real fv_arrays_mod::fv_flags_type::delt_max = 1.

Maximum allowed magnitude of the dissipative heating rate, K s−1; larger magnitudes are clipped to this amount. This can help avoid instability that can occur due to strong heating when d_con > 0. A value of 0.008 (a rate equivalent to about 800 K/day) is sufficient to stabilize the model at 3-km resolution. Set to 1. by default, which effectively disables this limitation.

Definition at line 736 of file fv_arrays.F90.

◆ dnats

integer fv_arrays_mod::fv_flags_type::dnats = 0

The number of tracers which are not to be advected by the dynamical core, but still passed into the dynamical core; the last dnats+pnats tracers in field_table are not advected. 0 by default.

Definition at line 679 of file fv_arrays.F90.

◆ do_f3d

logical fv_arrays_mod::fv_flags_type::do_f3d = .false.

Definition at line 420 of file fv_arrays.F90.

◆ do_held_suarez

logical fv_arrays_mod::fv_flags_type::do_held_suarez = .false.

Whether to use Held-Suarez forcing. Requires adiabatic to be false. The default is .false.; this option has no effect if not running solo_core.

Definition at line 807 of file fv_arrays.F90.

◆ do_reed_physics

logical fv_arrays_mod::fv_flags_type::do_reed_physics = .false.

Definition at line 810 of file fv_arrays.F90.

◆ do_sat_adj

logical fv_arrays_mod::fv_flags_type::do_sat_adj = .false.

Definition at line 419 of file fv_arrays.F90.

◆ do_schmidt

logical fv_arrays_mod::fv_flags_type::do_schmidt = .false.

Whether to enable grid stretching and rotation using stretch_fac, target_lat, and target_lon. The default value is .false.

Definition at line 546 of file fv_arrays.F90.

◆ do_skeb

logical fv_arrays_mod::fv_flags_type::do_skeb = .false.

save dissipation estimate

Definition at line 918 of file fv_arrays.F90.

◆ do_uni_zfull

logical fv_arrays_mod::fv_flags_type::do_uni_zfull = .false.

Whether to compute z_full (the height of each modellayer, as opposed to z_half, the height of each model interface) as the midpoint of the layer, as is done for the nonhydrostatic solver, instead of the height of the location where p = p the mean pressure in the layer. This option is not available for fvGFS or the solo_core. The default is .false.

Definition at line 946 of file fv_arrays.F90.

◆ do_vort_damp

logical fv_arrays_mod::fv_flags_type::do_vort_damp = .false.

Whether to apply flux damping (of strength governed by 'vtdm4') to the fluxes of vorticity, air mass, and nonhydrostatic vertical velocity (there is no dynamically correct way to add explicit diffusion to the tracer fluxes). The form is the same as is used for the divergence damping, including the same order (from 'nord') damping, unless 'nord' = 0, in which case this damping is fourth-order, or if 'nord' = 3,in which case this damping is sixth-order (instead of eighth-order). We recommend enabling this damping when the linear or non-monotonic horizontal advection schemes are enabled, but is unnecessary and not recommended when using monotonic advection. The default is .false.

Definition at line 432 of file fv_arrays.F90.

◆ dry_mass

real fv_arrays_mod::fv_flags_type::dry_mass = 98290.

If adjust_dry_mass is .true., sets the global dry air mass, measured in the globally-averaged surface pressure (Pascals) by adding or removing mass from the lowest layer of the atmosphere as needed. The default value is 98290. (Pa).

Definition at line 720 of file fv_arrays.F90.

◆ dwind_2d

logical fv_arrays_mod::fv_flags_type::dwind_2d = .false.

Whether to use a simpler & faster algorithm for interpolating the A-grid (cell-centered) wind tendencies computed from the physics to the D-grid. Typically, the A-grid wind tendencies are first converted in 3D cartesian coordinates and then interpolated before converting back to 2D local coordinates. When this option enabled, a much simpler but less accurate 2D interpolation is used. False by default.

Definition at line 774 of file fv_arrays.F90.

◆ dx_const

real(kind=r_grid) fv_arrays_mod::fv_flags_type::dx_const = 1000.

Specifies the (uniform) grid-cell-width in the x-direction on a doubly-periodic grid (grid_type = 4) in meters. The default value is 1000.

Definition at line 989 of file fv_arrays.F90.

◆ dy_const

real(kind=r_grid) fv_arrays_mod::fv_flags_type::dy_const = 1000.

Specifies the (uniform) grid-cell-width in the y-direction on a doubly-periodic grid (grid_type = 4) in meters. The default value is 1000.

Definition at line 993 of file fv_arrays.F90.

◆ ecmwf_ic

logical fv_arrays_mod::fv_flags_type::ecmwf_ic = .false.

If external_ic = .true., reads initial conditions from ECMWF analyses. The default is .false.

Definition at line 889 of file fv_arrays.F90.

◆ external_eta

logical fv_arrays_mod::fv_flags_type::external_eta = .false.

If .true., reads the interface coefficients ak and bk from either the restart file (if restarting) or from the external initial condition file (if nggps_ic or ecwmf_ic are .true.). This overrides the hard-coded levels in fv_eta. The default is .false.

Definition at line 911 of file fv_arrays.F90.

◆ external_ic

logical fv_arrays_mod::fv_flags_type::external_ic = .false.

Whether to initialize the models state using the data in an externally specified file, given in res_latlon_dynamics. By default this file is assumed to be a legacy lat-lon FV core restart file; set either ncep_ic or fv_diag_ic to .true.to override this behavior. The default is .false. Note that external_ic = .true. will cause the model to re-initialize the dynamical fields from the input dataset regardless of whether warm_start is set.

Definition at line 903 of file fv_arrays.F90.

◆ fac_n_spl

real fv_arrays_mod::fv_flags_type::fac_n_spl = 1.0

factor multiplying n_split up tp forecast hour fhouri

Definition at line 605 of file fv_arrays.F90.

◆ fhouri

real fv_arrays_mod::fv_flags_type::fhouri = 0.0

forecast hour up to which the number of small dynamics (acoustic) time steps are nint(n_split*fac_n_spl)

Definition at line 606 of file fv_arrays.F90.

◆ fill

logical fv_arrays_mod::fv_flags_type::fill = .false.

Fills in negative tracer values by taking positive tracers from the cells above and below. This option is useful when the physical parameterizations produced negatives. The default is .false.

Definition at line 792 of file fv_arrays.F90.

◆ fill_dp

logical fv_arrays_mod::fv_flags_type::fill_dp = .false.

Like 'fill' except for p, the hydrostatic pressure thickness. When the filling occurs a diagnostic message is printed out, which is helpful for diagnosing where the problem may be occurring. Typically, a crash is inevitable if the pressure filling is needed; thus, this option is often better for debugging than as a safety valve. The default is .false.

Definition at line 795 of file fv_arrays.F90.

◆ fill_wz

logical fv_arrays_mod::fv_flags_type::fill_wz = .false.

Definition at line 803 of file fv_arrays.F90.

◆ filter_phys

logical fv_arrays_mod::fv_flags_type::filter_phys = .false.

Definition at line 773 of file fv_arrays.F90.

◆ full_zs_filter

logical fv_arrays_mod::fv_flags_type::full_zs_filter = .false.

Whether to apply the on-line topography filter during initialization. Only active if get_nggps_ic = .true. This is so topography filtering can be performed on the initial conditions output by the pre-processing tools, which currently do not support topography filter- ing for some configurations (such as the nested grid); this also allows the user to easily test changes to the topography filtering on the simulation. Note that for all other initialization methods (if external_ic = .true.) the on-line topography filter will be applied automatically during the initialization of the topography. The default value is .false.

Definition at line 399 of file fv_arrays.F90.

◆ fv_debug

logical fv_arrays_mod::fv_flags_type::fv_debug = .false.

Whether to turn on additional diagnostics in fv_dynamics. The default is .false.

Definition at line 827 of file fv_arrays.F90.

◆ fv_diag_ic

logical fv_arrays_mod::fv_flags_type::fv_diag_ic = .false.

Definition at line 901 of file fv_arrays.F90.

◆ fv_land

logical fv_arrays_mod::fv_flags_type::fv_land = .false.

Whether to create terrain deviation and land fraction for output to mg_drag restart files, for use in mg_drag and in the land model. The default is .false; .true. is recommended when, and only when, initializing the model, since the mg_drag files created provide a much more accurate terrain representation for the mountain gravity wave drag parameterization and for the land surface roughness than either computes internally. This has no effect on the representation of the terrain in the dynamics.

Definition at line 861 of file fv_arrays.F90.

◆ fv_sg_adj

integer fv_arrays_mod::fv_flags_type::fv_sg_adj = -1

Timescale (in seconds) at which to remove two-delta-z instability when the local (between two adjacent levels) Richardson number is less than 1. This is achieved by local mixing, which conserves mass, momentum, and total energy. Values of 0 or smaller disable this feature. If n_sponge < 0 then the mixing is applied only to the top n_sponge layers of the domain. Set to -1 (inactive) by default. The proper range is 0 to 3600.

Definition at line 692 of file fv_arrays.F90.

◆ gfs_phil

logical fv_arrays_mod::fv_flags_type::gfs_phil = .false.

if .T., compute geopotential inside of GFS physics

Definition at line 892 of file fv_arrays.F90.

◆ grid_file

character(len=120) fv_arrays_mod::fv_flags_type::grid_file = 'Inline'

Definition at line 262 of file fv_arrays.F90.

◆ grid_name

character(len=80) fv_arrays_mod::fv_flags_type::grid_name = 'Gnomonic'

Definition at line 261 of file fv_arrays.F90.

◆ grid_number

integer, pointer fv_arrays_mod::fv_flags_type::grid_number

Convenience pointers.

Definition at line 1010 of file fv_arrays.F90.

◆ grid_type

integer fv_arrays_mod::fv_flags_type::grid_type = 0

-1: read from file; 0: ED Gnomonic 0: the "true" equal-distance Gnomonic grid 1: the traditional equal-distance Gnomonic grid 2: the equal-angular Gnomonic grid 3: the lat-lon grid – to be implemented 4: double periodic boundary condition on Cartesian grid 5: channel flow on Cartesian grid

Definition at line 263 of file fv_arrays.F90.

◆ hord_dp

integer fv_arrays_mod::fv_flags_type::hord_dp = 9

Horizontal advection scheme for mass. A positivity constraint may be warranted for hord_dp but not strictly necessary. 9 by default.

Definition at line 308 of file fv_arrays.F90.

◆ hord_mt

integer fv_arrays_mod::fv_flags_type::hord_mt = 9

Horizontal advection scheme for momentum fluxes. A complete list of kord options is given in the corresponding table in Appendix A of the FV3 technical document. The default value is 9, which uses the third-order piecewise-parabolic method with the monotonicity constraint of Huynh, which is less diffusive than other constraints. For hydrostatic simulation, 8 (the L04 monotonicity constraint) is recommended; for nonhydrostatic simulation, the completely unlimited (“linear” or non-monotone) PPM scheme is recommended. If no monotonicity constraint is applied, enabling the flux damping (do_vort_damp = .true.) is highly recommended to control grid-scale noise. It is also recommended that hord_mt, hord_vt, hord_tm, and hord_dp use the same value, to ensure consistent transport of all dynamical fields, unless a positivity constraint on mass advection (hord_dp) is desired.

Definition at line 273 of file fv_arrays.F90.

◆ hord_tm

integer fv_arrays_mod::fv_flags_type::hord_tm = 9

Horizontal advection scheme for potential temperature and layer thickness in nonhydrostatic simulations. 9 by default.

Definition at line 305 of file fv_arrays.F90.

◆ hord_tr

integer fv_arrays_mod::fv_flags_type::hord_tr = 12

Tracer transport options:

Horizontal advection scheme for tracers. The default is 12. This value can differ from the other hord options since tracers are subcycled (if inline_q == .false.) and require positive-definite advection to control the appearance of non-physical negative masses. 8 (fastest) or 10 (least diffusive) are typically recommended.

Definition at line 317 of file fv_arrays.F90.

◆ hord_vt

integer fv_arrays_mod::fv_flags_type::hord_vt = 9

Vorticity & w transport options:

Horizontal advection scheme for absolute vorticity and for vertical velocity in nonhydrostatic simulations. 9 by default.

Definition at line 300 of file fv_arrays.F90.

◆ hybrid_z

logical fv_arrays_mod::fv_flags_type::hybrid_z = .false.

Whether to use a hybrid-height coordinate, instead of the usual sigma-p coordinate. The default value is .false. (Not currently maintained.)

Definition at line 953 of file fv_arrays.F90.

◆ hydrostatic

logical fv_arrays_mod::fv_flags_type::hydrostatic = .true.

Whether to use the hydrostatic or nonhydrostatic solver. The default is .true.

Definition at line 933 of file fv_arrays.F90.

◆ inline_q

logical fv_arrays_mod::fv_flags_type::inline_q = .false.

Whether to compute tracer transport in-line with the rest of the dynamics instead of sub-cycling, so that tracer transport is done at the same time and on the same time step as is p and potential temperature. False by default; if true, q_split and z_tracer are ignored.

Definition at line 527 of file fv_arrays.F90.

◆ k_split

integer fv_arrays_mod::fv_flags_type::k_split = 1

Number of vertical remappings per dt_atmos (physics timestep). 1 by default.

Definition at line 610 of file fv_arrays.F90.

◆ ke_bg

real fv_arrays_mod::fv_flags_type::ke_bg = 0.

background KE production (m^2/s^3) over a small step Use this to conserve total energy if consv_te=0

Definition at line 748 of file fv_arrays.F90.

◆ kord_mt

integer fv_arrays_mod::fv_flags_type::kord_mt = 8

Vertical remapping scheme for the winds. 8 by default; 9 is recommended as the safest option, although 10, and 11 can also be useful. See corresponding table in Appendix A of the FV3 technical document for a complete list of kord options.

Definition at line 290 of file fv_arrays.F90.

◆ kord_tm

integer fv_arrays_mod::fv_flags_type::kord_tm = -8

Vertical remapping scheme for temperature. If positive (not recommended), then vertical remapping is performed on total energy instead of temperature (see 'remap_t'). The default value is -8.

Definition at line 312 of file fv_arrays.F90.

◆ kord_tr

integer fv_arrays_mod::fv_flags_type::kord_tr = 8

The vertical remapping scheme for tracers. The default is 8. 9 or 11 recommended. It is often recommended to use the same value for 'kord_tr' as for 'kord_tm'.

Definition at line 324 of file fv_arrays.F90.

◆ kord_wz

integer fv_arrays_mod::fv_flags_type::kord_wz = 8

Vertical remapping scheme for vertical velocity in nonhydrostatic simulations. 8 by default; 9 recommended. It is also recommended to use the same value for 'kord_wz' as for 'kord_mt'.

Definition at line 295 of file fv_arrays.F90.

◆ lim_fac

real fv_arrays_mod::fv_flags_type::lim_fac = 1.0

linear scheme limiting factor, 1: hord = 5, 3: hord = 6

Definition at line 331 of file fv_arrays.F90.

◆ m_split

integer fv_arrays_mod::fv_flags_type::m_split = 0

Number of time splits for Riemann solver.

Definition at line 609 of file fv_arrays.F90.

◆ make_hybrid_z

logical fv_arrays_mod::fv_flags_type::make_hybrid_z = .false.

Converts the vertical coordinate to a hybrid-height coordinate, instead of the usual sigma-p coordinate. Requires hybrid_z = .true. The default value is .false.

Definition at line 961 of file fv_arrays.F90.

◆ make_nh

logical fv_arrays_mod::fv_flags_type::make_nh = .false.

Whether to re-initialize the nonhydrostatic state, by recomputing dz from hydrostatic balance and setting w to 0. The default is false.

Definition at line 957 of file fv_arrays.F90.

◆ moist_phys

logical fv_arrays_mod::fv_flags_type::moist_phys = .true.

Run with moist physics.

Definition at line 806 of file fv_arrays.F90.

◆ mountain

logical fv_arrays_mod::fv_flags_type::mountain = .true.

Takes topography into account when initializing the model. Set this to .true. to apply the terrain filter (if n_zs_filter = 2 or 4) upon startup; also set to True when cold starting so that the topography can be initialized. Only set this to .false. if you wish to cold-start without any topography; this value is ignored for the aquaplanet test_case = 14. The default is .true. It is highly recommended TO NOT ALTER this value unless you know what you are doing.

Definition at line 830 of file fv_arrays.F90.

◆ n_split

integer fv_arrays_mod::fv_flags_type::n_split = 0

The number of small dynamics (acoustic) time steps between vertical remapping. 0 by default, in which case the model produces a good first guess by examining the resolution, dt_atmos, and k_split.

Definition at line 600 of file fv_arrays.F90.

◆ n_sponge

integer fv_arrays_mod::fv_flags_type::n_sponge = 1

Controls the number of layers at the upper boundary on which the 2Dx filter is applied. This does not control the sponge layer. The default value is 0.

Definition at line 499 of file fv_arrays.F90.

◆ n_zs_filter

integer fv_arrays_mod::fv_flags_type::n_zs_filter = 0

Additional (after the fact) terrain filter (to further smooth the terrain after cold start)

Number of times to apply a diffusive filter to the topography upon startup, if mountain is True and the model is not being cold-started. This is applied every time the model is warm-started, so if you want to smooth the topography make sure this is set to 0 after the first simulation. If initializing the model from cold-start the topography is already being filtered by an amount appropriate for the model resolution. 0 by default.

Definition at line 383 of file fv_arrays.F90.

◆ na_init

integer fv_arrays_mod::fv_flags_type::na_init = 0

Number of forward-backward dynamics steps used to initialize adiabatic solver. This is useful for spinning up the nonhydrostatic state from the hydrostatic GFS analyses. 0 by default. Recommended to set this to a non-zero value (1 or 2 is typically sufficient) when initializing from GFS or ECMWF analyses.

Definition at line 700 of file fv_arrays.F90.

◆ ncep_ic

logical fv_arrays_mod::fv_flags_type::ncep_ic = .false.

If external_ic = .true., this variable says whether the file in res_latlon_dynamics is an NCEP analysis or reanalysis file. This option zeros out all tracer fields except specific humidity. The default is .false.

Definition at line 880 of file fv_arrays.F90.

◆ ncnst

integer fv_arrays_mod::fv_flags_type::ncnst = 0

Number of tracer species advected by fv_tracer in the dynamical core.
Typically this is set automatically by reading in values from field_table, but ncnst can be set to a smaller value so only the first ncnst tracers listed in field_table are not advected. 0 by default, which will use the value from field_table.

Definition at line 668 of file fv_arrays.F90.

◆ ndims

integer fv_arrays_mod::fv_flags_type::ndims = 2

Definition at line 688 of file fv_arrays.F90.

◆ nf_omega

integer fv_arrays_mod::fv_flags_type::nf_omega = 1

Number of times to apply second-order smoothing to the diagnosed omega. When 0 the filter is disabled. 1 by default.

Definition at line 689 of file fv_arrays.F90.

◆ nggps_ic

logical fv_arrays_mod::fv_flags_type::nggps_ic = .false.

If external_ic = .true., reads initial conditions from horizontally-interpolated output from chgres. The default is .false. Additional options are available through external_ic_nml.

Definition at line 885 of file fv_arrays.F90.

◆ no_dycore

logical fv_arrays_mod::fv_flags_type::no_dycore = .false.

Disables execution of the dynamical core, only running the initialization, diagnostic, and I/O routines, and any physics that may be enabled. Essentially turns the model into a column physics model. The default is .false.

Definition at line 421 of file fv_arrays.F90.

◆ non_ortho

logical fv_arrays_mod::fv_flags_type::non_ortho = .true.

Definition at line 805 of file fv_arrays.F90.

◆ nord

integer fv_arrays_mod::fv_flags_type::nord = 1

Order of divergence damping: 0 for second-order; 1 for fourth-order (default); 2 for sixth-order; 3 for eighth-order. Sixth-order may yield a better solution for low resolutions (one degree or coarser) by virtue of it being more scale-selective and will not damp moderately well-resolved disturbances as much as does lower-order damping.

Definition at line 333 of file fv_arrays.F90.

◆ nord_tr

integer fv_arrays_mod::fv_flags_type::nord_tr = 0

Order of tracer damping; values mean the same as for 'nord'. The default value is 0.

Definition at line 339 of file fv_arrays.F90.

◆ nord_zs_filter

integer fv_arrays_mod::fv_flags_type::nord_zs_filter = 4

Order of the topography filter applied to n_zs_filter. Set to 2 to get a second-order filter, or 4 to get a fourth-order filter; other values do no filtering. 0 by default. This should not be set to a non-zero value on multiple successive simulations; the filter is applied every time the model restarts. This option is useful for testing the terrain filter, and SHOULD NOT BE USED FOR REGULAR RUNS.
use del-2 (2) OR del-4 (4)

Definition at line 391 of file fv_arrays.F90.

◆ npx

integer fv_arrays_mod::fv_flags_type::npx

Number of grid corners in the x-direction on one tile of the domain; so one more than the number of grid cells across a tile. On the cubed sphere this is one more than the number of cells across a cube face. Must be set.

Definition at line 651 of file fv_arrays.F90.

◆ npy

integer fv_arrays_mod::fv_flags_type::npy

Number of grid corners in the y-direction on one tile of the domain. This value should be identical to npx on a cubed-sphere grid; doubly periodic or nested grids do not have this restriction. Must be set.

Definition at line 655 of file fv_arrays.F90.

◆ npz

integer fv_arrays_mod::fv_flags_type::npz

Number of vertical levels. Each choice of npz comes with a pre-defined set of hybrid sigma-pressure levels and model top (see fv_eta.F90). Must be set.

Definition at line 659 of file fv_arrays.F90.

◆ npz_rst

integer fv_arrays_mod::fv_flags_type::npz_rst = 0

If using a restart file with a different number of vertical levels, set npz_rst to be the number of levels in your restart file. The model will then remap the restart file data to the vertical coordinates specified by npz. 0 by default; if 0 or equal to npz no remapping is done.

Definition at line 663 of file fv_arrays.F90.

◆ nt_phys

integer fv_arrays_mod::fv_flags_type::nt_phys = 0

Definition at line 726 of file fv_arrays.F90.

◆ nt_prog

integer fv_arrays_mod::fv_flags_type::nt_prog = 0

Definition at line 725 of file fv_arrays.F90.

◆ ntiles

integer fv_arrays_mod::fv_flags_type::ntiles = 1

Number of tiles on the domain. For the cubed sphere, this should be 6, one tile for each face of the cubed sphere; normally for most other domains (including nested grids) this should be set to 1. Must be set.

Definition at line 683 of file fv_arrays.F90.

◆ nudge

logical fv_arrays_mod::fv_flags_type::nudge = .false.

Whether to use the nudging towards the state in some externally-supplied file (such as from reanalysis or another simulation). Further nudging options are set in fv_nwp_nudge_nml. The default is .false.

Definition at line 872 of file fv_arrays.F90.

◆ nudge_dz

logical fv_arrays_mod::fv_flags_type::nudge_dz = .false.

During the adiabatic initialization (na_init > 0), if set to .true., delz is nudged back to the value specified in the initial conditions, instead of nudging the temperature back to the initial value. Nudging delz is simpler (faster), doesn’t require consideration of the virtual temperature effect, and may be more stable. .false.by default.

Definition at line 707 of file fv_arrays.F90.

◆ nudge_ic

logical fv_arrays_mod::fv_flags_type::nudge_ic = .false.

Same as nudge, but works in adiabatic solo_core simulations to nudge the field to a single external analysis file. The default is .false.

Definition at line 876 of file fv_arrays.F90.

◆ nudge_qv

logical fv_arrays_mod::fv_flags_type::nudge_qv = .false.

During the adiabatic initialization (na_init > 0), if set to .true., the water vapor profile is nudged to an analytic fit to the HALOE climatology. This is to improve the water vapor concentrations in GFS initial conditions, especially in the stratosphere, where values can be several times higher than observed. This nudging is unnecessary for other ICs, especially the ECMWF initial conditions. The default is .false.

Definition at line 965 of file fv_arrays.F90.

◆ nwat

integer fv_arrays_mod::fv_flags_type::nwat = 3

Number of water species to be included in condensate and water vapor loading. The masses of the first nwat tracer species will be added to the dry air mass, so that p is the mass of dry air, water vapor, and the included condensate species. The value used depends on the microphysics in the physics package you are using. For GFS physics with only a single condensate species, set to 2. For schemes with prognostic cloud water and cloud ice, such as GFDL AM2/AM3/AM4 Rotsteyn-Klein or Morrison-Gettlean microphysics, set to 3. For warm-rain (Kessler) microphysics set to 4 (with an inactive ice tracer), which only handles three species but uses 4 to avoid interference with the R-K physics. For schemes such as WSM5 or Ferrier that have prognostic rain and snow but not hail, set to 5 (not yet implemented). For six-category schemes that also have prognostic hail or graupel, such as the GFDL, Thompson, or WSM6 microphysics, set to 6. A value of 0 turns off condensate loading. The default value is 3.

Definition at line 507 of file fv_arrays.F90.

◆ old_divg_damp

logical fv_arrays_mod::fv_flags_type::old_divg_damp = .false.

parameter to revert damping parameters back to values defined in a previous revision old_values: d2_bg_k1 = 6. d2_bg_k2 = 4. d2_divg_max_k1 = 0.02 d2_divg_max_k2 = 0.01 damp_k_k1 = 0. damp_k_k2 = 0. current_values: d2_bg_k1 = 4. d2_bg_k2 = 2. d2_divg_max_k1 = 0.15 d2_divg_max_k2 = 0.08 damp_k_k1 = 0.2 damp_k_k2 = 0.12

Definition at line 837 of file fv_arrays.F90.

◆ p_fac

real fv_arrays_mod::fv_flags_type::p_fac = 0.05

Safety factor for minimum nonhydrostatic pressures, which will be limited so the full pressure is no less than p_fac times the hydrostatic pressure. This is only of concern in mid-top or high-top models with very low pressures near the model top, and has no effect in most simulations. The pressure limiting activates only when model is in danger of blowup due to unphysical negative total pressures. Only used if 'hydrostatic' = .false.and the semi-implicit solver is used. The proper range is 0 to 0.25. The default value is 0.05.

Definition at line 579 of file fv_arrays.F90.

◆ p_ref

real fv_arrays_mod::fv_flags_type::p_ref = 1.E5

Surface pressure used to construct a horizontally-uniform reference vertical pressure profile, used in some simple physics packages in the solo_core and in the Rayleigh damping. This should not be confused with the actual, horizontally-varying pressure levels used for all other dynamical calculations. The default value is 1.e5. CHANGING THIS VALUE IS STRONGLY DISCOURAGED.

Definition at line 713 of file fv_arrays.F90.

◆ phys_hydrostatic

logical fv_arrays_mod::fv_flags_type::phys_hydrostatic = .true.

Option to enable hydrostatic application of heating from the physics in a nonhydrostatic simulation: heating is applied in hydrostatic balance, causing the entire atmospheric column to expand instantaneously. If .false., heating from the physics is applied simply as a temperature tendency. The default value is .true.; ignored if hydrostatic = .true.

Definition at line 936 of file fv_arrays.F90.

◆ pnats

integer fv_arrays_mod::fv_flags_type::pnats = 0

The number of tracers not to advect by the dynamical core. Unlike dnats, these tracers are not seen by the dynamical core. The last pnats entries in field_table are not advected. The default value is 0.

Definition at line 674 of file fv_arrays.F90.

◆ print_freq

integer fv_arrays_mod::fv_flags_type::print_freq = 0

number of hours between print out of max/min and air/tracer mass diagnostics to standard output. 0 by default, which never prints out any output; set to -1 to see output after every dt_at-mos. Computing these diagnostics requires some computational overhead

Definition at line 638 of file fv_arrays.F90.

◆ q_split

integer fv_arrays_mod::fv_flags_type::q_split = 0

number of time steps for sub-cycled tracer advection. The default value is 0 (recommended), in which case the model determines the number of time steps from the global maximum wind speed at each call to the tracer advection.

Definition at line 633 of file fv_arrays.F90.

◆ range_warn

logical fv_arrays_mod::fv_flags_type::range_warn = .false.

Checks whether the values of the prognostic variables are within a reasonable range at the end of a dynamics time step, and prints a warning if not. The default is .false.; adds computational, overhead so we only recommend using this when debugging.

Definition at line 786 of file fv_arrays.F90.

◆ read_increment

logical fv_arrays_mod::fv_flags_type::read_increment = .false.

read in analysis increment and add to restart

Definition at line 916 of file fv_arrays.F90.

◆ reed_cond_only

logical fv_arrays_mod::fv_flags_type::reed_cond_only = .false.

Definition at line 811 of file fv_arrays.F90.

◆ regional

logical fv_arrays_mod::fv_flags_type::regional = .false.

Default setting for the regional domain.

Definition at line 1005 of file fv_arrays.F90.

◆ remap_t

logical fv_arrays_mod::fv_flags_type::remap_t = .true.

Whether the vertical remapping is performed on (virtual) temperature instead of (virtual) potential temperature. Since typically potential temperature increases exponentially from layer to layer near the top boundary, the cubic-spline interpolation in the vertical remapping will have difficulty with the exponential profile. Temperature does not have this problem and will often yield a more accurate result. The default is .true.

Definition at line 848 of file fv_arrays.F90.

◆ reproduce_sum

logical fv_arrays_mod::fv_flags_type::reproduce_sum = .true.

uses an exactly-reproducible global sum operation performed when computing the global energy for consv_te. This is used because the FMS routine mpp_sum() is not bit-wise reproducible due to its handling of floating-point arithmetic, and so can return different answers for (say) different processor layouts. The default is .true.

Definition at line 812 of file fv_arrays.F90.

◆ res_latlon_dynamics

character(len=128) fv_arrays_mod::fv_flags_type::res_latlon_dynamics = 'INPUT/fv_rst.res.nc'

If external_ic =.true.gives the filename of the input IC file. The default is 'INPUT/fv_rst.res.nc'.

Definition at line 921 of file fv_arrays.F90.

◆ res_latlon_tracers

character(len=128) fv_arrays_mod::fv_flags_type::res_latlon_tracers = 'INPUT/atmos_tracers.res.nc'

If external_ic =.true.and both ncep_ic and fv_diag_ic are.false., this variable gives the filename of the initial conditions for the tracers, assumed to be a legacy lat-lon FV core restart file. The default is 'INPUT/atmos_tracers.res.nc'.

Definition at line 923 of file fv_arrays.F90.

◆ reset_eta

logical fv_arrays_mod::fv_flags_type::reset_eta = .false.

Definition at line 578 of file fv_arrays.F90.

◆ rf_cutoff

real fv_arrays_mod::fv_flags_type::rf_cutoff = 30.E2

Pressure below which no Rayleigh damping is applied if tau > 0.

Definition at line 771 of file fv_arrays.F90.

◆ rf_fast

logical fv_arrays_mod::fv_flags_type::rf_fast = .false.

Option controlling whether to apply Rayleigh damping (for tau > 0) on the dynamic/acoustic timestep rather than on the physics timestep. This can help stabilize the model by applying the damping more weakly more frequently, so the instantaneous amount of damping (and thereby heat added) is reduced. The default is .false., which applies the Rayleigh drag every physics timestep.

Definition at line 409 of file fv_arrays.F90.

◆ scale_z

real fv_arrays_mod::fv_flags_type::scale_z = 0.

diff_z = scale_z**2 * 0.25

Definition at line 328 of file fv_arrays.F90.

◆ shift_fac

real fv_arrays_mod::fv_flags_type::shift_fac = 18.

Westward zonal rotation (or shift) of cubed-sphere grid from its natural orientation with cube face centers at 0, 90, 180, and 270 degrees longitude. The shift, in degrees, is 180/shift_fac. This shift does not move the poles. By default this is set to 18, shifting the grid westward 180/18=10 degrees, so that the edges of the cube do not run through the mountains of Japan; all standard CM2.x, AM3, CM3, and HiRAM simulations use this orientation of the grid. Requires do_schmidt = .false.

Definition at line 536 of file fv_arrays.F90.

◆ skeb_npass

integer fv_arrays_mod::fv_flags_type::skeb_npass = 11

Filter dissipation estimate "skeb_npass" times.

Definition at line 919 of file fv_arrays.F90.

◆ srf_init

logical fv_arrays_mod::fv_flags_type::srf_init = .false.

Definition at line 829 of file fv_arrays.F90.

◆ stretch_fac

real(kind=r_grid) fv_arrays_mod::fv_flags_type::stretch_fac = 1.

Stretching factor for the Schmidt transformation. This is the factor by which tile 6 of the cubed sphere will be shrunk, with the grid size shrinking accordingly. The default value is 1, which performs no grid stretching. Requires do_schmidt =.true. THE MODEL WILL CRASH IF stretch_fac IS SET TO ZERO. Values of up to 40 have been found useful and stable for short-term cloud-scale integrations.

Definition at line 549 of file fv_arrays.F90.

◆ target_lat

real(kind=r_grid) fv_arrays_mod::fv_flags_type::target_lat = -90.

Latitude (in degrees) to which the center of tile 6 will be rotated; if stretching is done with stretch_fac the center of the high-resolution part of the grid will be at this latitude. -90 by default, which does no grid rotation (the Schmidt transformation rotates the south pole to the appropriate target). Requires do_schmidt = .true.

Definition at line 558 of file fv_arrays.F90.

◆ target_lon

real(kind=r_grid) fv_arrays_mod::fv_flags_type::target_lon = 0.

Longitude to which the center of tile 6 will be rotated. 0 by default. Requires do_schmidt = .true.

Definition at line 565 of file fv_arrays.F90.

◆ tau

real fv_arrays_mod::fv_flags_type::tau = 0.

Time scale (in days) for Rayleigh friction applied to horizontal and vertical winds; lost kinetic energy is converted to heat, except on nested grids. The default value is 0.0, which disables damping. Larger values yield less damping. For models with tops at 1 mb or lower values between 10 and 30 are useful for preventing overly-strong polar night jets; for higher-top hydrostatic models values between 5 and 15 should be considered; and for non-hydrostatic models values of 10 or less should be
considered, with smaller values for higher-resolution.

Definition at line 762 of file fv_arrays.F90.

◆ tau_h2o

real fv_arrays_mod::fv_flags_type::tau_h2o = 0.

Time-scale (days) for simple methane chemistry to act as a source of water in the stratosphere. Can be useful if the stratosphere dries out too quickly; consider a value between 60 and 120 days if this is the case. The default value is 0., which disables the methane chemistry. Values less than zero apply the chemistry above 100 mb; else applied above 30 mb. Requires 'adiabatic' to be .false.

Definition at line 727 of file fv_arrays.F90.

◆ trdm2

real fv_arrays_mod::fv_flags_type::trdm2 = 0.0

Coefficient for del-2 tracer damping.

Definition at line 363 of file fv_arrays.F90.

◆ use_hydro_pressure

logical fv_arrays_mod::fv_flags_type::use_hydro_pressure = .false.

Whether to compute hydrostatic pressure for input to the physics. Currently only enabled for the fvGFS model. Ignored in hydrostatic simulations. The default is .false.

Definition at line 942 of file fv_arrays.F90.

◆ use_logp

logical fv_arrays_mod::fv_flags_type::use_logp = .false.

Enables a variant of the Lin pressure-gradient force algorithm, which uses the logarithm of pressure instead of the Exner function (as in [lin1997explicit]). This yields more accurate results for regions that are nearly isothermal. Ignored if 'hydrostatic' = .true. The default is .false.

Definition at line 613 of file fv_arrays.F90.

◆ use_ncep_phy

logical fv_arrays_mod::fv_flags_type::use_ncep_phy = .false.

Definition at line 900 of file fv_arrays.F90.

◆ use_new_ncep

logical fv_arrays_mod::fv_flags_type::use_new_ncep = .false.

Definition at line 899 of file fv_arrays.F90.

◆ use_old_omega

logical fv_arrays_mod::fv_flags_type::use_old_omega = .true.

Definition at line 444 of file fv_arrays.F90.

◆ vtdm4

real fv_arrays_mod::fv_flags_type::vtdm4 = 0.0

Coefficient for background other-variable damping. The value of 'vtdm4' should be less than that of 'd4_bg'. A good first guess for 'vtdm4' is about one-third the value of d4_bg. Requires 'do_vort_damp' to be .true. Disabled for values less than 1.e-3. Other- variable damping should not be used if a monotonic horizontal advection scheme is used. The default value is 0.0.

Definition at line 356 of file fv_arrays.F90.

◆ w_max

real fv_arrays_mod::fv_flags_type::w_max = 75.

max w (m/s) threshold for hydostatiic adjustment

Definition at line 329 of file fv_arrays.F90.

◆ warm_start

logical fv_arrays_mod::fv_flags_type::warm_start = .true.

Whether to start from restart files, instead of cold-starting the model. True by default; if this is set to .true. and restart files cannot be found the model will stop.

Definition at line 523 of file fv_arrays.F90.

◆ write_3d_diags

logical fv_arrays_mod::fv_flags_type::write_3d_diags = .true.

whether to write out three-dimensional dynamical diagnostic fields (those defined in fv_diagnostics.F90). This is useful for runs with multiple grids if you only want very large 3D diagnostics written out for (say) a nested grid, and not for the global grid. False by default.

Definition at line 643 of file fv_arrays.F90.

◆ z_min

real fv_arrays_mod::fv_flags_type::z_min = 0.05

min ratio of dz_nonhydrostatic/dz_hydrostatic

Definition at line 330 of file fv_arrays.F90.

◆ z_tracer

logical fv_arrays_mod::fv_flags_type::z_tracer = .false.

Whether to transport sub-cycled tracers layer-by-layer, each with its own computed sub-cycling time step (if q_split = 0). This may improve efficiency for very large numbers of tracers. The default value is .false.; currently not implemented.

Definition at line 856 of file fv_arrays.F90.


The documentation for this type was generated from the following file: