From e218478b9aca181aca9e910299b2107099aa6b52 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Tue, 30 Nov 2021 13:02:27 -0500 Subject: [PATCH] Expanded comments describing src/core code (#13) Added comments describing all of the real variables in the src/core directory that did not previously have them, and corrected multiple spelling error in this same directory. This PR also updates the description of the MOM6 code structure in MOM.F90 and adds dimensional descriptions to the commends describing a number of variables. Descriptions of the units of some variables were also added or corrected in a number of the src/tracer files. Changes are largely restricted to comments describing varibles, although a few unused variables were removed. One comment added to MOM_dynamics_split_RK2.F90 notes a possible memory leak, due to a missing deallocate call, but the actual correction of this leak will be delayed to a later, much more targeted commit. All answers and output are bitwise identical. --- src/core/MOM.F90 | 128 +++++++++++------ src/core/MOM_CoriolisAdv.F90 | 42 +++--- src/core/MOM_PressureForce_FV.F90 | 8 +- src/core/MOM_PressureForce_Montgomery.F90 | 4 +- src/core/MOM_barotropic.F90 | 130 ++++++++++-------- src/core/MOM_continuity_PPM.F90 | 94 +++++++------ src/core/MOM_density_integrals.F90 | 14 +- src/core/MOM_dynamics_split_RK2.F90 | 107 +++++++------- src/core/MOM_dynamics_unsplit.F90 | 12 +- src/core/MOM_dynamics_unsplit_RK2.F90 | 15 +- src/core/MOM_forcing_type.F90 | 65 ++++----- src/core/MOM_grid.F90 | 18 +-- src/core/MOM_interface_heights.F90 | 2 +- src/core/MOM_isopycnal_slopes.F90 | 10 +- src/core/MOM_open_boundary.F90 | 109 +++++++-------- src/core/MOM_verticalGrid.F90 | 2 +- src/tracer/MOM_CFC_cap.F90 | 2 +- src/tracer/MOM_OCMIP2_CFC.F90 | 2 +- src/tracer/MOM_generic_tracer.F90 | 2 +- src/tracer/MOM_lateral_boundary_diffusion.F90 | 81 ++++++----- src/tracer/MOM_neutral_diffusion.F90 | 7 +- src/tracer/MOM_offline_aux.F90 | 4 +- src/tracer/MOM_tracer_diabatic.F90 | 51 +++---- src/tracer/MOM_tracer_registry.F90 | 2 +- src/tracer/advection_test_tracer.F90 | 2 +- src/tracer/boundary_impulse_tracer.F90 | 2 +- src/tracer/dye_example.F90 | 2 +- src/tracer/ideal_age_example.F90 | 2 +- src/tracer/oil_tracer.F90 | 4 +- 29 files changed, 498 insertions(+), 425 deletions(-) diff --git a/src/core/MOM.F90 b/src/core/MOM.F90 index f3a2dff337..90c3d7bef2 100644 --- a/src/core/MOM.F90 +++ b/src/core/MOM.F90 @@ -168,7 +168,7 @@ module MOM !>@{ 3-d state field diagnostic IDs integer :: id_u = -1, id_v = -1, id_h = -1 !>@} - !> 2-d state field diagnotic ID + !> 2-d state field diagnostic ID integer :: id_ssh_inst = -1 end type MOM_diag_IDs @@ -249,7 +249,7 @@ module MOM real :: dt_therm !< thermodynamics time step [T ~> s] logical :: thermo_spans_coupling !< If true, thermodynamic and tracer time !! steps can span multiple coupled time steps. - integer :: nstep_tot = 0 !< The total number of dynamic timesteps tcaaken + integer :: nstep_tot = 0 !< The total number of dynamic timesteps taken !! so far in this run segment logical :: count_calls = .false. !< If true, count the calls to step_MOM, rather than the !! number of dynamics steps in nstep_tot @@ -330,7 +330,7 @@ module MOM real :: bad_val_col_thick !< Minimum column thickness before triggering bad value message [Z ~> m] logical :: answers_2018 !< If true, use expressions for the surface properties that recover !! the answers from the end of 2018. Otherwise, use more appropriate - !! expressions that differ at roundoff for non-Boussinsq cases. + !! expressions that differ at roundoff for non-Boussinesq cases. logical :: use_particles !< Turns on the particles package character(len=10) :: particle_type !< Particle types include: surface(default), profiling and sail drone. @@ -493,7 +493,8 @@ subroutine step_MOM(forces_in, fluxes_in, sfc_state, Time_start, time_int_in, CS real :: dt_therm ! a limited and quantized version of CS%dt_therm [T ~> s] real :: dt_therm_here ! a further limited value of dt_therm [T ~> s] - real :: wt_end, wt_beg + real :: wt_end, wt_beg ! Fractional weights of the future pressure at the end + ! and beginning of the current time step [nondim] real :: bbl_time_int ! The amount of time over which the calculated BBL ! properties will apply, for use in diagnostics, or 0 ! if it is not to be calculated anew [T ~> s]. @@ -1517,12 +1518,12 @@ subroutine step_offline(forces, fluxes, sfc_state, Time_start, time_interval, CS ! 3D pointers real, dimension(:,:,:), pointer :: & - uhtr => NULL(), vhtr => NULL(), & - eatr => NULL(), ebtr => NULL(), & - h_end => NULL() + uhtr => NULL(), & ! Accumulated zonal thickness fluxes to advect tracers [H L2 ~> m3 or kg] + vhtr => NULL(), & ! Accumulated meridional thickness fluxes to advect tracers [H L2 ~> m3 or kg] + eatr => NULL(), & ! Layer entrainment rates across the interface above [H ~> m or kg m-2] + ebtr => NULL(), & ! Layer entrainment rates across the interface below [H ~> m or kg m-2] + h_end => NULL() ! Layer thicknesses at the end of a step [H ~> m or kg m-2] - ! 2D Array for diagnostics - real, dimension(SZI_(CS%G),SZJ_(CS%G)) :: eta_pre, eta_end type(time_type) :: Time_end ! End time of a segment, as a time type ! Grid-related pointer assignments @@ -1719,9 +1720,13 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, & integer :: turns ! Number of grid quarter-turns ! Initial state on the input index map - real, allocatable, dimension(:,:,:) :: u_in, v_in, h_in - real, allocatable, dimension(:,:), target :: frac_shelf_in - real, allocatable, dimension(:,:,:), target :: T_in, S_in + real, allocatable :: u_in(:,:,:) ! Initial zonal velocities [L T-1 ~> m s-1] + real, allocatable :: v_in(:,:,:) ! Initial meridional velocities [L T-1 ~> m s-1] + real, allocatable :: h_in(:,:,:) ! Initial layer thicknesses [H ~> m or kg m-2] + real, allocatable, target :: frac_shelf_in(:,:) ! Initial fraction of the total cell area occupied + ! by an ice shelf [nondim] + real, allocatable, target :: T_in(:,:,:) ! Initial temperatures [degC] + real, allocatable, target :: S_in(:,:,:) ! Initial salinities [ppt] type(ocean_OBC_type), pointer :: OBC_in => NULL() type(sponge_CS), pointer :: sponge_in_CSp => NULL() type(ALE_sponge_CS), pointer :: ALE_sponge_in_CSp => NULL() @@ -1754,7 +1759,7 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, & logical :: bulkmixedlayer ! If true, a refined bulk mixed layer scheme is used ! with nkml sublayers and nkbl buffer layer. - logical :: use_temperature ! If true, temp and saln used as state variables. + logical :: use_temperature ! If true, temperature and salinity used as state variables. logical :: use_frazil ! If true, liquid seawater freezes if temp below freezing, ! with accumulated heat deficit returned to surface ocean. logical :: bound_salinity ! If true, salt is added to keep salinity above @@ -1781,7 +1786,9 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, & integer :: nkml, nkbl, verbosity, write_geom integer :: dynamics_stencil ! The computational stencil for the calculations ! in the dynamic core. - real :: conv2watt, conv2salt + real :: conv2watt ! A conversion factor from temperature fluxes to heat + ! fluxes [J m-2 H-1 degC-1 ~> J m-3 degC-1 or J kg-1 degC-1] + real :: conv2salt ! A conversion factor for salt fluxes [m H-1 ~> 1] or [kg m-2 H-1 ~> 1] real :: RL2_T2_rescale, Z_rescale, QRZ_rescale ! Unit conversion factors character(len=48) :: flux_units, S_flux_units @@ -2178,7 +2185,7 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, & CS%G => G_in endif - ! TODO: It is unlikey that test_grid_copy and rotate_index would work at the + ! TODO: It is unlikely that test_grid_copy and rotate_index would work at the ! same time. It may be possible to enable both but for now we prevent it. if (test_grid_copy .and. CS%rotate_index) & call MOM_error(FATAL, "Grid cannot be copied during index rotation.") @@ -2218,7 +2225,7 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, & call rotate_hor_index(HI_in, turns, HI) ! NOTE: If indices are rotated, then G and G_in must both be initialized separately, and ! the dynamic grid must be created to handle the grid rotation. G%domain has already been - ! initialzed above. + ! initialized above. call MOM_grid_init(G, param_file, US, HI, bathymetry_at_vel=bathy_at_vel) call create_dyn_horgrid(dG, HI, bathymetry_at_vel=bathy_at_vel) call clone_MOM_domain(G%Domain, dG%Domain) @@ -3142,7 +3149,7 @@ subroutine extract_surface_state(CS, sfc_state_in) !! calculation of properties of the uppermost ocean [nondim] or [Z H-1 ~> 1 or m3 kg-1] ! After the ANSWERS_2018 flag has been obsoleted, H_rescale will be 1. real :: delT(SZI_(CS%G)) !< Depth integral of T-T_freeze [Z degC ~> m degC] - logical :: use_temperature !< If true, temp and saln used as state variables. + logical :: use_temperature !< If true, temperature and salinity are used as state variables. integer :: i, j, k, is, ie, js, je, nz, numberOfErrors, ig, jg integer :: isd, ied, jsd, jed integer :: iscB, iecB, jscB, jecB, isdB, iedB, jsdB, jedB @@ -3528,10 +3535,18 @@ end subroutine extract_surface_state !> Rotate initialization fields from input to rotated arrays. subroutine rotate_initial_state(u_in, v_in, h_in, T_in, S_in, & use_temperature, turns, u, v, h, T, S) - real, dimension(:,:,:), intent(in) :: u_in, v_in, h_in, T_in, S_in - logical, intent(in) :: use_temperature - integer, intent(in) :: turns - real, dimension(:,:,:), intent(out) :: u, v, h, T, S + real, dimension(:,:,:), intent(in) :: u_in !< Zonal velocity on the initial grid [L T-1 ~> m s-1] + real, dimension(:,:,:), intent(in) :: v_in !< Meridional velocity on the initial grid [L T-1 ~> m s-1] + real, dimension(:,:,:), intent(in) :: h_in !< Layer thickness on the initial grid [H ~> m or kg m-2] + real, dimension(:,:,:), intent(in) :: T_in !< Temperature on the initial grid [degC] + real, dimension(:,:,:), intent(in) :: S_in !< Salinity on the initial grid [ppt] + logical, intent(in) :: use_temperature !< If true, temperature and salinity are active + integer, intent(in) :: turns !< The number quarter-turns to apply + real, dimension(:,:,:), intent(out) :: u !< Zonal velocity on the rotated grid [L T-1 ~> m s-1] + real, dimension(:,:,:), intent(out) :: v !< Meridional velocity on the rotated grid [L T-1 ~> m s-1] + real, dimension(:,:,:), intent(out) :: h !< Layer thickness on the rotated grid [H ~> m or kg m-2] + real, dimension(:,:,:), intent(out) :: T !< Temperature on the rotated grid [degC] + real, dimension(:,:,:), intent(out) :: S !< Salinity on the rotated grid [ppt] call rotate_vector(u_in, v_in, turns, u, v) call rotate_array(h_in, turns, h) @@ -3852,29 +3867,48 @@ end subroutine MOM_end !! !! \verbatim !! ../MOM +!! |-- ac !! |-- config_src -!! | |-- coupled_driver -!! | |-- dynamic -!! | `-- solo_driver -!! |-- examples -!! | |-- CM2G +!! | |-- drivers +!! | ! |-- FMS_cap +!! | ! |-- ice_solo_driver +!! | ! |-- mct_cap +!! | ! |-- nuopc_cap +!! | ! |-- solo_driver +!! | ! `-- unit_drivers +!! | |-- external +!! | ! |-- drifters +!! | ! |-- GFDL_ocean_BGC +!! | ! `-- ODA_hooks +!! | |-- infra +!! | ! |-- FMS1 +!! | ! `-- FMS2 +!! | `-- memory +!! | ! |-- dynamic_nonsymmetric +!! | ! `-- dynamic_symmetric +!! |-- docs +!! |-- pkg +!! | |-- CVMix-src !! | |-- ... -!! | `-- torus_advection_test +!! | `-- MOM6_DA_hooks !! `-- src +!! |-- ALE !! |-- core !! |-- diagnostics !! |-- equation_of_state !! |-- framework !! |-- ice_shelf !! |-- initialization +!! |-- ocean_data_assim !! |-- parameterizations +!! | |-- CVMix !! | |-- lateral !! | `-- vertical !! |-- tracer !! `-- user !! \endverbatim !! -!! Rather than describing each file here, each directory contents +!! Rather than describing each file here, selected directory contents !! will be described to give a broad overview of the MOM code !! structure. !! @@ -3883,27 +3917,35 @@ end subroutine MOM_end !! Only one or two of these directories are used in compiling any, !! particular run. !! -!! * config_src/coupled_driver: +!! * config_src/drivers/FMS-cap: !! The files here are used to couple MOM as a component in a larger !! run driven by the FMS coupler. This includes code that converts !! various forcing fields into the code structures and flux and unit !! conventions used by MOM, and converts the MOM surface fields !! back to the forms used by other FMS components. !! -!! * config_src/dynamic: -!! The only file here is the version of MOM_memory.h that is used -!! for dynamic memory configurations of MOM. +!! * config_src/drivers/nuopc-cap: +!! The files here are used to couple MOM as a component in a larger +!! run driven by the NUOPC coupler. This includes code that converts +!! various forcing fields into the code structures and flux and unit +!! conventions used by MOM, and converts the MOM surface fields +!! back to the forms used by other NUOPC components. !! -!! * config_src/solo_driver: +!! * config_src/drivers/solo_driver: !! The files here are include the _main driver that is used when !! MOM is configured as an ocean-only model, as well as the files !! that specify the surface forcing in this configuration. !! -!! The directories under examples provide a large number of working -!! configurations of MOM, along with reference solutions for several -!! different compilers on GFDL's latest large computer. The versions -!! of MOM_memory.h in these directories need not be used if dynamic -!! memory allocation is desired, and the answers should be unchanged. +!! * config_src/external: +!! The files here are mostly just stubs, so that MOM6 can compile +!! with calls to the public interfaces external packages, but +!! without actually requiring those packages themselves. In more +!! elaborate configurations, would be linked to the actual code for +!! those external packages rather than these simple stubs. +!! +!! * config_src/memory/dynamic-symmetric: +!! The only file here is the version of MOM_memory.h that is used +!! for dynamic memory configurations of MOM. !! !! The directories under src contain most of the MOM files. These !! files are used in every configuration using MOM. @@ -3916,7 +3958,7 @@ end subroutine MOM_end !! subroutine argument lists. !! !! * src/diagnostics: -!! The files here calculate various diagnostics that are anciliary +!! The files here calculate various diagnostics that are ancilliary !! to the model itself. While most of these diagnostics do not !! directly affect the model's solution, there are some, like the !! calculation of the deformation radius, that are used in some @@ -3973,13 +4015,19 @@ end subroutine MOM_end !! to build an appropriate makefile, and path_names should be edited !! to reflect the actual location of the desired source code. !! +!! The separate MOM-examples git repository provides a large number +!! of working configurations of MOM, along with reference solutions for several +!! different compilers on GFDL's latest large computer. The versions +!! of MOM_memory.h in these directories need not be used if dynamic +!! memory allocation is desired, and the answers should be unchanged. +!! !! !! There are 3 publicly visible subroutines in this file (MOM.F90). !! * step_MOM steps MOM over a specified interval of time. !! * MOM_initialize calls initialize and does other initialization !! that does not warrant user modification. !! * extract_surface_state determines the surface (bulk mixed layer -!! if traditional isoycnal vertical coordinate) properties of the +!! if traditional isopycnal vertical coordinate) properties of the !! current model state and packages pointers to these fields into an !! exported structure. !! diff --git a/src/core/MOM_CoriolisAdv.F90 b/src/core/MOM_CoriolisAdv.F90 index 93dda759c2..953d64c1f0 100644 --- a/src/core/MOM_CoriolisAdv.F90 +++ b/src/core/MOM_CoriolisAdv.F90 @@ -163,7 +163,7 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS, pbv) real, dimension(SZI_(G),SZJB_(G)) :: & hArea_v, & ! The cell area weighted thickness interpolated to v points ! times the effective areas [H L2 ~> m3 or kg]. - KEy, & ! The meridonal gradient of Kinetic energy per unit mass [L T-2 ~> m s-2], + KEy, & ! The meridional gradient of Kinetic energy per unit mass [L T-2 ~> m s-2], ! KEy = d/dy KE. vh_center ! Transport based on arithmetic mean h at v-points [H L2 T-1 ~> m3 s-1 or kg s-1] real, dimension(SZI_(G),SZJ_(G)) :: & @@ -204,17 +204,17 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS, pbv) real :: uhc, vhc ! Centered estimates of uh and vh [H L2 T-1 ~> m3 s-1 or kg s-1]. real :: uhm, vhm ! The input estimates of uh and vh [H L2 T-1 ~> m3 s-1 or kg s-1]. - real :: c1, c2, c3, slope ! Nondimensional parameters for the Coriolis limiter scheme. + real :: c1, c2, c3, slope ! Nondimensional parameters for the Coriolis limiter scheme [nondim] - real :: Fe_m2 ! Nondimensional temporary variables asssociated with - real :: rat_lin ! the ARAKAWA_LAMB_BLEND scheme. + real :: Fe_m2 ! Temporary variable associated with the ARAKAWA_LAMB_BLEND scheme [nondim] + real :: rat_lin ! Temporary variable associated with the ARAKAWA_LAMB_BLEND scheme [nondim] real :: rat_m1 ! The ratio of the maximum neighboring inverse thickness - ! to the minimum inverse thickness minus 1. rat_m1 >= 0. + ! to the minimum inverse thickness minus 1 [nondim]. rat_m1 >= 0. real :: AL_wt ! The relative weight of the Arakawa & Lamb scheme to the - ! Arakawa & Hsu scheme, nondimensional between 0 and 1. + ! Arakawa & Hsu scheme [nondim], between 0 and 1. real :: Sad_wt ! The relative weight of the Sadourny energy scheme to - ! the other two with the ARAKAWA_LAMB_BLEND scheme, - ! nondimensional between 0 and 1. + ! the other two with the ARAKAWA_LAMB_BLEND scheme [nondim], + ! between 0 and 1. real :: Heff1, Heff2 ! Temporary effective H at U or V points [H ~> m or kg m-2]. real :: Heff3, Heff4 ! Temporary effective H at U or V points [H ~> m or kg m-2]. @@ -232,7 +232,7 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS, pbv) ! hf_gKEu, hf_gKEv, & ! accel. due to KE gradient x fract. thickness [L T-2 ~> m s-2]. ! hf_rvxu, hf_rvxv ! accel. due to RV x fract. thickness [L T-2 ~> m s-2]. ! 3D diagnostics hf_gKEu etc. are commented because there is no clarity on proper remapping grid option. - ! The code is retained for degugging purposes in the future. + ! The code is retained for debugging purposes in the future. ! Diagnostics for thickness multiplied momentum budget terms real, allocatable, dimension(:,:,:) :: h_gKEu, h_gKEv ! h x gKEu, h x gKEv [H L T-2 ~> m2 s-2]. @@ -676,7 +676,7 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS, pbv) endif enddo ; enddo endif - ! Add in the additonal terms with Arakawa & Lamb. + ! Add in the additional terms with Arakawa & Lamb. if ((CS%Coriolis_Scheme == ARAKAWA_LAMB81) .or. & (CS%Coriolis_Scheme == AL_BLEND)) then ; do j=js,je ; do I=Isq,Ieq CAu(I,j,k) = CAu(I,j,k) + & @@ -876,7 +876,7 @@ subroutine CorAdCalc(u, v, h, uh, vh, CAu, CAv, OBC, AD, G, GV, US, CS, pbv) ! Diagnostics for terms multiplied by fractional thicknesses ! 3D diagnostics hf_gKEu etc. are commented because there is no clarity on proper remapping grid option. - ! The code is retained for degugging purposes in the future. + ! The code is retained for debugging purposes in the future. !if (CS%id_hf_gKEu > 0) then ! allocate(hf_gKEu(G%IsdB:G%IedB,G%jsd:G%jed,GV%ke)) ! do k=1,nz ; do j=js,je ; do I=Isq,Ieq @@ -1025,7 +1025,7 @@ end subroutine CorAdCalc !> Calculates the acceleration due to the gradient of kinetic energy. subroutine gradKE(u, v, h, KE, KEx, KEy, k, OBC, G, GV, US, CS) - type(ocean_grid_type), intent(in) :: G !< Ocen grid structure + type(ocean_grid_type), intent(in) :: G !< Ocean grid structure type(verticalGrid_type), intent(in) :: GV !< Vertical grid structure real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), intent(in) :: u !< Zonal velocity [L T-1 ~> m s-1] real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), intent(in) :: v !< Meridional velocity [L T-1 ~> m s-1] @@ -1061,7 +1061,7 @@ subroutine gradKE(u, v, h, KE, KEx, KEy, k, OBC, G, GV, US, CS) G%areaCv(i,J-1)*(v(i,J-1,k)*v(i,J-1,k)) ) )*0.25*G%IareaT(i,j) enddo ; enddo elseif (CS%KE_Scheme == KE_SIMPLE_GUDONOV) then - ! The following discretization of KE is based on the one-dimensinal Gudonov + ! The following discretization of KE is based on the one-dimensional Gudonov ! scheme which does not take into account any geometric factors do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 up = 0.5*( u(I-1,j,k) + ABS( u(I-1,j,k) ) ) ; up2 = up*up @@ -1071,7 +1071,7 @@ subroutine gradKE(u, v, h, KE, KEx, KEy, k, OBC, G, GV, US, CS) KE(i,j) = ( max(up2,um2) + max(vp2,vm2) ) *0.5 enddo ; enddo elseif (CS%KE_Scheme == KE_GUDONOV) then - ! The following discretization of KE is based on the one-dimensinal Gudonov + ! The following discretization of KE is based on the one-dimensional Gudonov ! scheme but has been adapted to take horizontal grid factors into account do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 up = 0.5*( u(I-1,j,k) + ABS( u(I-1,j,k) ) ) ; up2a = up*up*G%areaCu(I-1,j) @@ -1108,16 +1108,16 @@ subroutine gradKE(u, v, h, KE, KEx, KEy, k, OBC, G, GV, US, CS) end subroutine gradKE -!> Initializes the control structure for coriolisadv_cs +!> Initializes the control structure for MOM_CoriolisAdv subroutine CoriolisAdv_init(Time, G, GV, US, param_file, diag, AD, CS) type(time_type), target, intent(in) :: Time !< Current model time - type(ocean_grid_type), intent(in) :: G !< Ocean grid structure - type(verticalGrid_type), intent(in) :: GV !< Vertical grid structure - type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type + type(ocean_grid_type), intent(in) :: G !< Ocean grid structure + type(verticalGrid_type), intent(in) :: GV !< Vertical grid structure + type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(param_file_type), intent(in) :: param_file !< Runtime parameter handles type(diag_ctrl), target, intent(inout) :: diag !< Diagnostics control structure - type(accel_diag_ptrs), target, intent(inout) :: AD !< Strorage for acceleration diagnostics - type(CoriolisAdv_CS), intent(inout) :: CS !< Control structure fro MOM_CoriolisAdv + type(accel_diag_ptrs), target, intent(inout) :: AD !< Storage for acceleration diagnostics + type(CoriolisAdv_CS), intent(inout) :: CS !< Control structure for MOM_CoriolisAdv ! Local variables ! This include declares and sets the variable "version". #include "version_variable.h" @@ -1405,7 +1405,7 @@ end subroutine CoriolisAdv_init !> Destructor for coriolisadv_cs subroutine CoriolisAdv_end(CS) - type(CoriolisAdv_CS), intent(inout) :: CS !< Control structure fro MOM_CoriolisAdv + type(CoriolisAdv_CS), intent(inout) :: CS !< Control structure for MOM_CoriolisAdv end subroutine CoriolisAdv_end !> \namespace mom_coriolisadv diff --git a/src/core/MOM_PressureForce_FV.F90 b/src/core/MOM_PressureForce_FV.F90 index 8436b92f80..5ead019717 100644 --- a/src/core/MOM_PressureForce_FV.F90 +++ b/src/core/MOM_PressureForce_FV.F90 @@ -199,7 +199,7 @@ subroutine PressureForce_FV_nonBouss(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_ if (use_EOS) then ! With a bulk mixed layer, replace the T & S of any layers that are - ! lighter than the the buffer layer with the properties of the buffer + ! lighter than the buffer layer with the properties of the buffer ! layer. These layers will be massless anyway, and it avoids any ! formal calculations with hydrostatically unstable profiles. if (nkmb>0) then @@ -230,7 +230,7 @@ subroutine PressureForce_FV_nonBouss(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_ ! If regridding is activated, do a linear reconstruction of salinity ! and temperature across each layer. The subscripts 't' and 'b' refer ! to top and bottom values within each layer (these are the only degrees - ! of freedeom needed to know the linear profile). + ! of freedom needed to know the linear profile). if ( use_ALE ) then if ( CS%Recon_Scheme == 1 ) then call TS_PLM_edge_values(ALE_CSp, S_t, S_b, T_t, T_b, G, GV, tv, h, CS%boundary_extrap) @@ -595,7 +595,7 @@ subroutine PressureForce_FV_Bouss(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_atm if (use_EOS) then ! With a bulk mixed layer, replace the T & S of any layers that are -! lighter than the the buffer layer with the properties of the buffer +! lighter than the buffer layer with the properties of the buffer ! layer. These layers will be massless anyway, and it avoids any ! formal calculations with hydrostatically unstable profiles. @@ -654,7 +654,7 @@ subroutine PressureForce_FV_Bouss(h, tv, PFu, PFv, G, GV, US, CS, ALE_CSp, p_atm ! If regridding is activated, do a linear reconstruction of salinity ! and temperature across each layer. The subscripts 't' and 'b' refer ! to top and bottom values within each layer (these are the only degrees - ! of freedeom needed to know the linear profile). + ! of freedom needed to know the linear profile). if ( use_ALE ) then if ( CS%Recon_Scheme == 1 ) then call TS_PLM_edge_values(ALE_CSp, S_t, S_b, T_t, T_b, G, GV, tv, h, CS%boundary_extrap) diff --git a/src/core/MOM_PressureForce_Montgomery.F90 b/src/core/MOM_PressureForce_Montgomery.F90 index 6a0831eca9..a827fb12d0 100644 --- a/src/core/MOM_PressureForce_Montgomery.F90 +++ b/src/core/MOM_PressureForce_Montgomery.F90 @@ -219,7 +219,7 @@ subroutine PressureForce_Mont_nonBouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pb ! Calculate in-situ specific volumes (alpha_star). ! With a bulk mixed layer, replace the T & S of any layers that are - ! lighter than the the buffer layer with the properties of the buffer + ! lighter than the buffer layer with the properties of the buffer ! layer. These layers will be massless anyway, and it avoids any ! formal calculations with hydrostatically unstable profiles. if (nkmb>0) then @@ -475,7 +475,7 @@ subroutine PressureForce_Mont_Bouss(h, tv, PFu, PFv, G, GV, US, CS, p_atm, pbce, ! Calculate in-situ densities (rho_star). ! With a bulk mixed layer, replace the T & S of any layers that are -! lighter than the the buffer layer with the properties of the buffer +! lighter than the buffer layer with the properties of the buffer ! layer. These layers will be massless anyway, and it avoids any ! formal calculations with hydrostatically unstable profiles. diff --git a/src/core/MOM_barotropic.F90 b/src/core/MOM_barotropic.F90 index ca1a7d20e5..a7e8194a84 100644 --- a/src/core/MOM_barotropic.F90 +++ b/src/core/MOM_barotropic.F90 @@ -1,4 +1,4 @@ -!> Baropotric solver +!> Barotropic solver module MOM_barotropic ! This file is part of MOM6. See LICENSE.md for the license. @@ -98,7 +98,7 @@ module MOM_barotropic type(group_pass_type) :: pass_eta_outer !< Structure for group halo pass end type BT_OBC_type -!> The barotropic stepping control stucture +!> The barotropic stepping control structure type, public :: barotropic_CS ; private real ALLOCABLE_, dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) :: frhatu !< The fraction of the total column thickness interpolated to u grid points in each layer [nondim]. @@ -132,8 +132,8 @@ module MOM_barotropic !< A limit on the rate at which eta_cor can be applied while avoiding instability !! [H T-1 ~> m s-1 or kg m-2 s-1]. This is only used if CS%bound_BT_corr is true. real ALLOCABLE_, dimension(NIMEMW_,NJMEMW_) :: & - ua_polarity, & !< Test vector components for checking grid polarity. - va_polarity, & !< Test vector components for checking grid polarity. + ua_polarity, & !< Test vector components for checking grid polarity [nondim] + va_polarity, & !< Test vector components for checking grid polarity [nondim] bathyT !< A copy of bathyT (ocean bottom depth) with wide halos [Z ~> m] real ALLOCABLE_, dimension(NIMEMW_,NJMEMW_) :: IareaT !< This is a copy of G%IareaT with wide halos, but will @@ -149,15 +149,15 @@ module MOM_barotropic real ALLOCABLE_, dimension(NIMEMBW_,NJMEMBW_) :: & q_D !< f / D at PV points [Z-1 T-1 ~> m-1 s-1]. - real, allocatable :: frhatu1(:,:,:) !< Predictor step values of frhatu stored for diagnostics. - real, allocatable :: frhatv1(:,:,:) !< Predictor step values of frhatv stored for diagnostics. + real, allocatable :: frhatu1(:,:,:) !< Predictor step values of frhatu stored for diagnostics [nondim] + real, allocatable :: frhatv1(:,:,:) !< Predictor step values of frhatv stored for diagnostics [nondim] type(BT_OBC_type) :: BT_OBC !< A structure with all of this modules fields !! for applying open boundary conditions. real :: dtbt !< The barotropic time step [T ~> s]. real :: dtbt_fraction !< The fraction of the maximum time-step that - !! should used. The default is 0.98. + !! should used [nondim]. The default is 0.98. real :: dtbt_max !< The maximum stable barotropic time step [T ~> s]. real :: dt_bt_filter !< The time-scale over which the barotropic mode solutions are !! filtered [T ~> s] if positive, or as a fraction of DT if @@ -166,7 +166,7 @@ module MOM_barotropic integer :: nstep_last = 0 !< The number of barotropic timesteps per baroclinic !! time step the last time btstep was called. real :: bebt !< A nondimensional number, from 0 to 1, that - !! determines the gravity wave time stepping scheme. + !! determines the gravity wave time stepping scheme [nondim]. !! 0.0 gives a forward-backward scheme, while 1.0 !! give backward Euler. In practice, bebt should be !! of order 0.2 or greater. @@ -209,7 +209,7 @@ module MOM_barotropic !! barotropic step when calculating the surface stress contribution to !! the barotropic acclerations. Otherwise use the depth based on bathyT. real :: BT_Coriolis_scale !< A factor by which the barotropic Coriolis acceleration anomaly - !! terms are scaled. + !! terms are scaled [nondim]. logical :: answers_2018 !< If true, use expressions for the barotropic solver that recover !! the answers from the end of 2018. Otherwise, use more efficient !! or general expressions. @@ -228,7 +228,7 @@ module MOM_barotropic logical :: tidal_sal_bug !< If true, the tidal self-attraction and loading anomaly in the !! barotropic solver has the wrong sign, replicating a long-standing !! bug. - real :: G_extra !< A nondimensional factor by which gtot is enhanced. + real :: G_extra !< A nondimensional factor by which gtot is enhanced [nondim]. integer :: hvel_scheme !< An integer indicating how the thicknesses at !! velocity points are calculated. Valid values are !! given by the parameters defined below: @@ -255,10 +255,10 @@ module MOM_barotropic !! truncated to maxvel [L T-1 ~> m s-1]. real :: CFL_trunc !< If clip_velocity is true, velocity components will !! be truncated when they are large enough that the - !! corresponding CFL number exceeds this value, nondim. + !! corresponding CFL number exceeds this value [nondim]. real :: maxCFL_BT_cont !< The maximum permitted CFL number associated with the !! barotropic accelerations from the summed velocities - !! times the time-derivatives of thicknesses. The + !! times the time-derivatives of thicknesses [nondim]. The !! default is 0.1, and there will probably be real !! problems if this were set close to 1. logical :: BT_cont_bounds !< If true, use the BT_cont_type variables to set limits @@ -321,7 +321,7 @@ module MOM_barotropic end type barotropic_CS -!> A desciption of the functional dependence of transport at a u-point +!> A description of the functional dependence of transport at a u-point type, private :: local_BT_cont_u_type real :: FA_u_EE !< The effective open face area for zonal barotropic transport !! drawing from locations far to the east [H L ~> m2 or kg m-1]. @@ -347,7 +347,7 @@ module MOM_barotropic !! time-integrated transport with INTEGRAL_BT_CONTINUITY [H L2 ~> m3 or kg]. end type local_BT_cont_u_type -!> A desciption of the functional dependence of transport at a v-point +!> A description of the functional dependence of transport at a v-point type, private :: local_BT_cont_v_type real :: FA_v_NN !< The effective open face area for meridional barotropic transport !! drawing from locations far to the north [H L ~> m2 or kg m-1]. @@ -451,13 +451,13 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, real, dimension(SZI_(G),SZJB_(G)), intent(out) :: vhbtav !< the barotropic meridional volume or mass !! fluxes averaged through the barotropic steps !! [H L2 T-1 ~> m3 s-1 or kg s-1]. - type(barotropic_CS), intent(inout) :: CS !< Barotropic control struct + type(barotropic_CS), intent(inout) :: CS !< Barotropic control structure real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), intent(in) :: visc_rem_u !< Both the fraction of the momentum !! originally in a layer that remains after a time-step of !! viscosity, and the fraction of a time-step's worth of a !! barotropic acceleration that a layer experiences after - !! viscosity is applied, in the zonal direction. Nondimensional - !! between 0 (at the bottom) and 1 (far above the bottom). + !! viscosity is applied, in the zonal direction [nondim]. + !! Visc_rem_u is between 0 (at the bottom) and 1 (far above). real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), intent(in) :: visc_rem_v !< Ditto for meridional direction [nondim]. type(accel_diag_ptrs), pointer :: ADp !< Acceleration diagnostic pointers type(ocean_OBC_type), pointer :: OBC !< The open boundary condition structure. @@ -489,19 +489,19 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, real :: wt_u(SZIB_(G),SZJ_(G),SZK_(GV)) ! wt_u and wt_v are the real :: wt_v(SZI_(G),SZJB_(G),SZK_(GV)) ! normalized weights to ! be used in calculating barotropic velocities, possibly with - ! sums less than one due to viscous losses. Nondimensional. + ! sums less than one due to viscous losses [nondim] real, dimension(SZIB_(G),SZJ_(G)) :: & - av_rem_u, & ! The weighted average of visc_rem_u, nondimensional. - tmp_u, & ! A temporary array at u points. + av_rem_u, & ! The weighted average of visc_rem_u [nondim] + tmp_u, & ! A temporary array at u points [L T-2 ~> m s-2] or [nondim] ubt_st, & ! The zonal barotropic velocity at the start of timestep [L T-1 ~> m s-1]. ubt_dt ! The zonal barotropic velocity tendency [L T-2 ~> m s-2]. real, dimension(SZI_(G),SZJB_(G)) :: & - av_rem_v, & ! The weighted average of visc_rem_v, nondimensional. - tmp_v, & ! A temporary array at v points. + av_rem_v, & ! The weighted average of visc_rem_v [nondim] + tmp_v, & ! A temporary array at v points [L T-2 ~> m s-2] or [nondim] vbt_st, & ! The meridional barotropic velocity at the start of timestep [L T-1 ~> m s-1]. vbt_dt ! The meridional barotropic velocity tendency [L T-2 ~> m s-2]. real, dimension(SZI_(G),SZJ_(G)) :: & - tmp_h, & ! A temporary array at h points. + tmp_h, & ! A temporary array at h points [nondim] e_anom ! The anomaly in the sea surface height or column mass ! averaged between the beginning and end of the time step, ! relative to eta_PF, with SAL effects included [H ~> m or kg m-2]. @@ -512,8 +512,8 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, real, dimension(SZIBW_(CS),SZJW_(CS)) :: & ubt, & ! The zonal barotropic velocity [L T-1 ~> m s-1]. bt_rem_u, & ! The fraction of the barotropic zonal velocity that remains - ! after a time step, the remainder being lost to bottom drag. - ! bt_rem_u is a nondimensional number between 0 and 1. + ! after a time step, the remainder being lost to bottom drag [nondim]. + ! bt_rem_u is between 0 and 1. BT_force_u, & ! The vertical average of all of the u-accelerations that are ! not explicitly included in the barotropic equation [L T-2 ~> m s-2]. u_accel_bt, & ! The difference between the zonal acceleration from the @@ -530,8 +530,8 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, uhbt_int, & ! The running time integral of uhbt over the time steps [H L2 ~> m3]. ubt_wtd, & ! A weighted sum used to find the filtered final ubt [L T-1 ~> m s-1]. ubt_trans, & ! The latest value of ubt used for a transport [L T-1 ~> m s-1]. - azon, bzon, & ! _zon & _mer are the values of the Coriolis force which - czon, dzon, & ! are applied to the neighboring values of vbtav & ubtav, + azon, bzon, & ! _zon and _mer are the values of the Coriolis force which + czon, dzon, & ! are applied to the neighboring values of vbtav and ubtav, amer, bmer, & ! respectively to get the barotropic inertial rotation cmer, dmer, & ! [T-1 ~> s-1]. Cor_u, & ! The zonal Coriolis acceleration [L T-2 ~> m s-2]. @@ -548,7 +548,7 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, vbt, & ! The meridional barotropic velocity [L T-1 ~> m s-1]. bt_rem_v, & ! The fraction of the barotropic meridional velocity that ! remains after a time step, the rest being lost to bottom - ! drag. bt_rem_v is a nondimensional number between 0 and 1. + ! drag [nondim]. bt_rem_v is between 0 and 1. BT_force_v, & ! The vertical average of all of the v-accelerations that are ! not explicitly included in the barotropic equation [L T-2 ~> m s-2]. v_accel_bt, & ! The difference between the meridional acceleration from the @@ -635,9 +635,9 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, ! of the reference geopotential with the sea surface height [nondim]. ! This is typically ~0.09 or less. real :: dgeo_de ! The constant of proportionality between geopotential and - ! sea surface height [nondim]. It is a nondimensional number of - ! order 1. For stability, this may be made larger - ! than the physical problem would suggest. + ! sea surface height [nondim]. It is of order 1, but for + ! stability this may be made larger than the physical + ! problem would suggest. real :: Instep ! The inverse of the number of barotropic time steps to take [nondim]. real :: wt_end ! The weighting of the final value of eta_PF [nondim] integer :: nstep ! The number of barotropic time steps to take. @@ -673,9 +673,23 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce, ! in roundoff and can be neglected [H ~> m or kg m-2]. real :: Idtbt ! The inverse of the barotropic time step [T-1 ~> s-1] - real, allocatable, dimension(:) :: wt_vel, wt_eta, wt_accel, wt_trans, wt_accel2 - real :: sum_wt_vel, sum_wt_eta, sum_wt_accel, sum_wt_trans - real :: I_sum_wt_vel, I_sum_wt_eta, I_sum_wt_accel, I_sum_wt_trans + real, allocatable :: wt_vel(:) ! The raw or relative weights of each of the barotropic timesteps + ! in determining the average velocities [nondim] + real, allocatable :: wt_eta(:) ! The raw or relative weights of each of the barotropic timesteps + ! in determining the average the average of eta [nondim] + real, allocatable :: wt_accel(:) ! The raw or relative weights of each of the barotropic timesteps + ! in determining the average accelerations [nondim] + real, allocatable :: wt_trans(:) ! The raw or relative weights of each of the barotropic timesteps + ! in determining the average transports [nondim] + real, allocatable :: wt_accel2(:) ! A potentially un-normalized copy of wt_accel [nondim] + real :: sum_wt_vel ! The sum of the raw weights used to find average velocities [nondim] + real :: sum_wt_eta ! The sum of the raw weights used to find average the average of eta [nondim] + real :: sum_wt_accel ! The sum of the raw weights used to find average accelerations [nondim] + real :: sum_wt_trans ! The sum of the raw weights used to find average transports [nondim] + real :: I_sum_wt_vel ! The inverse of the sum of the raw weights used to find average velocities [nondim] + real :: I_sum_wt_eta ! The inverse of the sum of the raw weights used to find the average of eta [nondim] + real :: I_sum_wt_accel ! The inverse of the sum of the raw weights used to find average accelerations [nondim] + real :: I_sum_wt_trans ! The inverse of the sum of the raw weights used to find average transports [nondim] real :: dt_filt ! The half-width of the barotropic filter [T ~> s]. real :: trans_wt1, trans_wt2 ! The weights used to compute ubt_trans and vbt_trans integer :: nfilter @@ -2914,13 +2928,9 @@ subroutine apply_velocity_OBCs(OBC, ubt, vbt, uhbt, vhbt, ubt_trans, vbt_trans, real :: v_inlet ! The meridional inflow velocity [L T-1 ~> m s-1] real :: uhbt_int_new ! The updated time-integrated zonal transport [H L2 ~> m3] real :: vhbt_int_new ! The updated time-integrated meridional transport [H L2 ~> m3] - real :: h_in ! The inflow thickess [H ~> m or kg m-2]. - real :: cff, Cx, Cy, tau - real :: dhdt, dhdx, dhdy + real :: h_in ! The inflow thickness [H ~> m or kg m-2]. real :: Idtbt ! The inverse of the barotropic time step [T-1 ~> s-1] integer :: i, j, is, ie, js, je - real, dimension(SZIB_(G),SZJB_(G)) :: grad - real, parameter :: eps = 1.0e-20 is = G%isc-halo ; ie = G%iec+halo ; js = G%jsc-halo ; je = G%jec+halo if (.not.(BT_OBC%apply_u_OBCs .or. BT_OBC%apply_v_OBCs)) return @@ -3262,7 +3272,7 @@ subroutine btcalc(h, G, GV, CS, h_u, h_v, may_use_default, OBC) type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure. real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), & intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2]. - type(barotropic_CS), intent(inout) :: CS !< Barotropic control struct + type(barotropic_CS), intent(inout) :: CS !< Barotropic control structure real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), & optional, intent(in) :: h_u !< The specified thicknesses at u-points [H ~> m or kg m-2]. real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), & @@ -3284,7 +3294,7 @@ subroutine btcalc(h, G, GV, CS, h_u, h_v, may_use_default, OBC) real :: h_harm ! The harmonic mean thicknesses [H ~> m or kg m-2]. real :: h_neglect ! A thickness that is so small it is usually lost ! in roundoff and can be neglected [H ~> m or kg m-2]. - real :: wt_arith ! The nondimensional weight for the arithmetic mean thickness. + real :: wt_arith ! The weight for the arithmetic mean thickness [nondim]. ! The harmonic mean uses a weight of (1 - wt_arith). real :: Rh ! A ratio of summed thicknesses, nondim. real :: e_u(SZIB_(G),SZK_(GV)+1) ! The interface heights at u-velocity and @@ -3606,7 +3616,7 @@ function uhbt_to_ubt(uhbt, BTC) result(ubt) real :: vsr ! Temporary variable used in the limiting the velocity [nondim]. real, parameter :: vs1 = 1.25 ! Nondimensional parameters used in limiting real, parameter :: vs2 = 2.0 ! the velocity, starting at vs1, with the - ! maximum increase of vs2, both nondim. + ! maximum increase of vs2, both [nondim]. integer :: itt, max_itt = 20 ! Find the value of ubt that gives uhbt. @@ -3741,7 +3751,7 @@ function vhbt_to_vbt(vhbt, BTC) result(vbt) real :: vsr ! Temporary variable used in the limiting the velocity [nondim]. real, parameter :: vs1 = 1.25 ! Nondimensional parameters used in limiting real, parameter :: vs2 = 2.0 ! the velocity, starting at vs1, with the - ! maximum increase of vs2, both nondim. + ! maximum increase of vs2, both [nondim]. integer :: itt, max_itt = 20 ! Find the value of vbt that gives vhbt. @@ -3937,7 +3947,7 @@ end subroutine set_local_BT_cont_types !> Adjust_local_BT_cont_types expands the range of velocities with a cubic curve -!! translating velocities into transports to match the inital values of velocities and +!! translating velocities into transports to match the initial values of velocities and !! summed transports when the velocities are larger than the first guesses of the cubic !! transition velocities used to set up the local_BT_cont types. subroutine adjust_local_BT_cont_types(ubt, uhbt, vbt, vhbt, BTCL_u, BTCL_v, & @@ -3964,10 +3974,6 @@ subroutine adjust_local_BT_cont_types(ubt, uhbt, vbt, vhbt, BTCL_u, BTCL_v, & !! provided if INTEGRAL_BT_CONTINUITY is true. ! Local variables - real, dimension(SZIBW_(MS),SZJW_(MS)) :: & - u_polarity, uBT_EE, uBT_WW, FA_u_EE, FA_u_E0, FA_u_W0, FA_u_WW - real, dimension(SZIW_(MS),SZJBW_(MS)) :: & - v_polarity, vBT_NN, vBT_SS, FA_v_NN, FA_v_N0, FA_v_S0, FA_v_SS real :: dt ! The baroclinic timestep [T ~> s] or 1.0 [nondim] real, parameter :: C1_3 = 1.0/3.0 integer :: i, j, is, ie, js, je, hs @@ -4072,9 +4078,9 @@ end subroutine BT_cont_to_face_areas !> Swap the values of two real variables subroutine swap(a,b) - real, intent(inout) :: a !< The first variable to be swapped. - real, intent(inout) :: b !< The second variable to be swapped. - real :: tmp + real, intent(inout) :: a !< The first variable to be swapped [arbitrary units] + real, intent(inout) :: b !< The second variable to be swapped [arbitrary units] + real :: tmp ! A temporary variable [arbitrary units] tmp = a ; a = b ; b = tmp end subroutine swap @@ -4089,7 +4095,7 @@ subroutine find_face_areas(Datu, Datv, G, GV, US, CS, MS, halo, eta, add_max) type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure. type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure. type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type - type(barotropic_CS), intent(in) :: CS !< Barotropic control struct + type(barotropic_CS), intent(in) :: CS !< Barotropic control structure integer, intent(in) :: halo !< The halo size to use, default = 1. real, dimension(MS%isdw:MS%iedw,MS%jsdw:MS%jedw), & optional, intent(in) :: eta !< The barotropic free surface height anomaly @@ -4183,7 +4189,7 @@ subroutine bt_mass_source(h, eta, set_cor, G, GV, CS) !! fluxes (and update the slowly varying part of eta_cor) !! (.true.) or whether to incrementally update the !! corrective fluxes. - type(barotropic_CS), intent(inout) :: CS !< Barotropic control struct + type(barotropic_CS), intent(inout) :: CS !< Barotropic control structure ! Local variables real :: h_tot(SZI_(G)) ! The sum of the layer thicknesses [H ~> m or kg m-2]. @@ -4249,8 +4255,8 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS, type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters. type(diag_ctrl), target, intent(inout) :: diag !< A structure that is used to regulate diagnostic !! output. - type(barotropic_CS), intent(inout) :: CS !< Barotropic control struct - type(MOM_restart_CS), intent(in) :: restart_CS !< MOM restart control struct + type(barotropic_CS), intent(inout) :: CS !< Barotropic control structure + type(MOM_restart_CS), intent(in) :: restart_CS !< MOM restart control structure logical, intent(out) :: calc_dtbt !< If true, the barotropic time step must !! be recalculated before stepping. type(BT_cont_type), pointer :: BT_cont !< A structure with elements that describe the @@ -4259,8 +4265,8 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS, type(tidal_forcing_CS), target, optional :: tides_CSp !< A pointer to the control structure of the !! tide module. -! This include declares and sets the variable "version". -#include "version_variable.h" + ! This include declares and sets the variable "version". +# include "version_variable.h" ! Local variables character(len=40) :: mdl = "MOM_barotropic" ! This module's name. real :: Datu(SZIBS_(G),SZJ_(G)) ! Zonal open face area [H L ~> m2 or kg m-1]. @@ -4286,7 +4292,9 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS, real :: det_de ! The partial derivative due to self-attraction and loading of the reference ! geopotential with the sea surface height when tides are enabled. ! This is typically ~0.09 or less. - real, allocatable, dimension(:,:) :: lin_drag_h + real, allocatable :: lin_drag_h(:,:) ! A spatially varying linear drag coefficient at tracer points + ! that acts on the barotropic flow [Z T-1 ~> m s-1]. + type(memory_size_type) :: MS type(group_pass_type) :: pass_static_data, pass_q_D_Cor type(group_pass_type) :: pass_bt_hbt_btav, pass_a_polarity @@ -4958,7 +4966,7 @@ end subroutine barotropic_init !> Copies ubtav and vbtav from private type into arrays subroutine barotropic_get_tav(CS, ubtav, vbtav, G, US) - type(barotropic_CS), intent(in) :: CS !< Barotropic control struct + type(barotropic_CS), intent(in) :: CS !< Barotropic control structure type(ocean_grid_type), intent(in) :: G !< Grid structure real, dimension(SZIB_(G),SZJ_(G)), intent(inout) :: ubtav !< Zonal barotropic velocity averaged !! over a baroclinic timestep [L T-1 ~> m s-1] @@ -5007,9 +5015,9 @@ end subroutine barotropic_end subroutine register_barotropic_restarts(HI, GV, param_file, CS, restart_CS) type(hor_index_type), intent(in) :: HI !< A horizontal index type structure. type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters. - type(barotropic_CS), intent(inout) :: CS !< Barotropic control struct + type(barotropic_CS), intent(inout) :: CS !< Barotropic control structure type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure. - type(MOM_restart_CS), intent(inout) :: restart_CS !< MOM restart control struct + type(MOM_restart_CS), intent(inout) :: restart_CS !< MOM restart control structure ! Local variables type(vardesc) :: vd(3) diff --git a/src/core/MOM_continuity_PPM.F90 b/src/core/MOM_continuity_PPM.F90 index eca822ac6b..90e1086bc5 100644 --- a/src/core/MOM_continuity_PPM.F90 +++ b/src/core/MOM_continuity_PPM.F90 @@ -106,15 +106,15 @@ subroutine continuity_PPM(u, v, hin, h, uh, vh, dt, G, GV, US, CS, OBC, pbv, uhb !< The fraction of zonal momentum originally !! in a layer that remains after a time-step of viscosity, and the !! fraction of a time-step's worth of a barotropic acceleration that - !! a layer experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). + !! a layer experiences after viscosity is applied [nondim]. + !! Visc_rem_u is between 0 (at the bottom) and 1 (far above the bottom). real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), & optional, intent(in) :: visc_rem_v !< The fraction of meridional momentum originally !! in a layer that remains after a time-step of viscosity, and the !! fraction of a time-step's worth of a barotropic acceleration that - !! a layer experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). + !! a layer experiences after viscosity is applied [nondim]. + !! Visc_rem_v is between 0 (at the bottom) and 1 (far above the bottom). real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), & optional, intent(out) :: u_cor !< The zonal velocities that give uhbt as the depth-integrated transport [L T-1 ~> m s-1]. @@ -239,11 +239,11 @@ subroutine zonal_mass_flux(u, h_in, uh, dt, G, GV, US, CS, LB, OBC, por_face_are optional, intent(in) :: visc_rem_u !< The fraction of zonal momentum originally in a layer that remains after a !! time-step of viscosity, and the fraction of a time-step's worth of a barotropic - !! acceleration that a layer experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). + !! acceleration that a layer experiences after viscosity is applied [nondim]. + !! Visc_rem_u is between 0 (at the bottom) and 1 (far above the bottom). real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), & optional, intent(out) :: u_cor - !< The zonal velocitiess (u with a barotropic correction) + !< The zonal velocities (u with a barotropic correction) !! that give uhbt as the depth-integrated transport, m s-1. type(BT_cont_type), optional, pointer :: BT_cont !< A structure with elements that describe the !! effective open face areas as a function of barotropic flow. @@ -254,13 +254,13 @@ subroutine zonal_mass_flux(u, h_in, uh, dt, G, GV, US, CS, LB, OBC, por_face_are real, dimension(SZIB_(G)) :: & du, & ! Corrective barotropic change in the velocity [L T-1 ~> m s-1]. du_min_CFL, & ! Min/max limits on du correction - du_max_CFL, & ! to avoid CFL violations + du_max_CFL, & ! to avoid CFL violations [L T-1 ~> m s-1] duhdu_tot_0, & ! Summed partial derivative of uh with u [H L ~> m2 or kg m-1]. uh_tot_0, & ! Summed transport with no barotropic correction [H L2 T-1 ~> m3 s-1 or kg s-1]. - visc_rem_max ! The column maximum of visc_rem. + visc_rem_max ! The column maximum of visc_rem [nondim]. logical, dimension(SZIB_(G)) :: do_I real, dimension(SZIB_(G),SZK_(GV)) :: & - visc_rem ! A 2-D copy of visc_rem_u or an array of 1's. + visc_rem ! A 2-D copy of visc_rem_u or an array of 1's [nondim]. real, dimension(SZIB_(G)) :: FAuI ! A list of sums of zonal face areas [H L ~> m2 or kg m-1]. real :: FA_u ! A sum of zonal face areas [H L ~> m2 or kg m-1]. real :: I_vrm ! 1.0 / visc_rem_max, nondim. @@ -533,8 +533,8 @@ subroutine zonal_flux_layer(u, h, h_L, h_R, uh, duhdu, visc_rem, dt, G, US, j, & real, dimension(SZIB_(G)), intent(in) :: visc_rem !< Both the fraction of the !! momentum originally in a layer that remains after a time-step !! of viscosity, and the fraction of a time-step's worth of a barotropic - !! acceleration that a layer experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). + !! acceleration that a layer experiences after viscosity is applied [nondim]. + !! Visc_rem is between 0 (at the bottom) and 1 (far above the bottom). real, dimension(SZI_(G)), intent(in) :: h !< Layer thickness [H ~> m or kg m-2]. real, dimension(SZI_(G)), intent(in) :: h_L !< Left thickness [H ~> m or kg m-2]. real, dimension(SZI_(G)), intent(in) :: h_R !< Right thickness [H ~> m or kg m-2]. @@ -635,8 +635,8 @@ subroutine zonal_face_thickness(u, h, h_L, h_R, h_u, dt, G, GV, US, LB, vol_CFL, optional, intent(in) :: visc_rem_u !< Both the fraction of the momentum originally in a layer that remains after !! a time-step of viscosity, and the fraction of a time-step's worth of a - !! barotropic acceleration that a layer experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). + !! barotropic acceleration that a layer experiences after viscosity is applied [nondim]. + !! Visc_rem_u is between 0 (at the bottom) and 1 (far above the bottom). ! Local variables real :: CFL ! The CFL number based on the local velocity and grid spacing [nondim] @@ -665,7 +665,7 @@ subroutine zonal_face_thickness(u, h, h_L, h_R, h_u, dt, G, GV, US, LB, vol_CFL, 3.0*curv_3*(CFL - 1.0)) else h_avg = 0.5 * (h_L(i+1,j,k) + h_R(i,j,k)) - ! The choice to use the arithmetic mean here is somewhat arbitrariy, but + ! The choice to use the arithmetic mean here is somewhat arbitrarily, but ! it should be noted that h_L(i+1,j,k) and h_R(i,j,k) are usually the same. h_marg = 0.5 * (h_L(i+1,j,k) + h_R(i,j,k)) ! h_marg = (2.0 * h_L(i+1,j,k) * h_R(i,j,k)) / & @@ -733,8 +733,8 @@ subroutine zonal_flux_adjust(u, h_in, h_L, h_R, uhbt, uh_tot_0, duhdu_tot_0, & real, dimension(SZIB_(G),SZK_(GV)), intent(in) :: visc_rem !< Both the fraction of the !! momentum originally in a layer that remains after a time-step of viscosity, and !! the fraction of a time-step's worth of a barotropic acceleration that a layer - !! experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). + !! experiences after viscosity is applied [nondim]. + !! Visc_rem is between 0 (at the bottom) and 1 (far above the bottom). real, dimension(SZIB_(G)), optional, intent(in) :: uhbt !< The summed volume flux !! through zonal faces [H L2 T-1 ~> m3 s-1 or kg s-1]. @@ -898,9 +898,9 @@ subroutine set_zonal_BT_cont(u, h_in, h_L, h_R, BT_cont, uh_tot_0, duhdu_tot_0, real, dimension(SZIB_(G),SZK_(GV)), intent(in) :: visc_rem !< Both the fraction of the !! momentum originally in a layer that remains after a time-step of viscosity, and !! the fraction of a time-step's worth of a barotropic acceleration that a layer - !! experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). - real, dimension(SZIB_(G)), intent(in) :: visc_rem_max !< Maximum allowable visc_rem. + !! experiences after viscosity is applied [nondim]. + !! Visc_rem is between 0 (at the bottom) and 1 (far above the bottom). + real, dimension(SZIB_(G)), intent(in) :: visc_rem_max !< Maximum allowable visc_rem [nondim]. integer, intent(in) :: j !< Spatial index. integer, intent(in) :: ish !< Start of index range. integer, intent(in) :: ieh !< End of index range. @@ -929,7 +929,7 @@ subroutine set_zonal_BT_cont(u, h_in, h_L, h_R, BT_cont, uh_tot_0, duhdu_tot_0, real :: FA_0 ! The effective face area with 0 barotropic transport [L H ~> m2 or kg m]. real :: FA_avg ! The average effective face area [L H ~> m2 or kg m], nominally given by ! the realized transport divided by the barotropic velocity. - real :: visc_rem_lim ! The larger of visc_rem and min_visc_rem [nondim] This + real :: visc_rem_lim ! The larger of visc_rem and min_visc_rem [nondim]. This ! limiting is necessary to keep the inverse of visc_rem ! from leading to large CFL numbers. real :: min_visc_rem ! The smallest permitted value for visc_rem that is used @@ -1059,11 +1059,11 @@ subroutine meridional_mass_flux(v, h_in, vh, dt, G, GV, US, CS, LB, OBC, por_fac optional, intent(in) :: visc_rem_v !< Both the fraction of the momentum !! originally in a layer that remains after a time-step of viscosity, !! and the fraction of a time-step's worth of a barotropic acceleration - !! that a layer experiences after viscosity is applied. Nondimensional between - !! 0 (at the bottom) and 1 (far above the bottom). + !! that a layer experiences after viscosity is applied [nondim]. + !! Visc_rem_v is between 0 (at the bottom) and 1 (far above the bottom). real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), & optional, intent(out) :: v_cor - !< The meridional velocitiess (v with a barotropic correction) + !< The meridional velocities (v with a barotropic correction) !! that give vhbt as the depth-integrated transport [L T-1 ~> m s-1]. type(BT_cont_type), optional, pointer :: BT_cont !< A structure with elements that describe !! the effective open face areas as a function of barotropic flow. @@ -1349,8 +1349,8 @@ subroutine merid_flux_layer(v, h, h_L, h_R, vh, dvhdv, visc_rem, dt, G, US, J, & real, dimension(SZI_(G)), intent(in) :: visc_rem !< Both the fraction of the !! momentum originally in a layer that remains after a time-step !! of viscosity, and the fraction of a time-step's worth of a barotropic - !! acceleration that a layer experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). + !! acceleration that a layer experiences after viscosity is applied [nondim]. + !! Visc_rem is between 0 (at the bottom) and 1 (far above the bottom). real, dimension(SZI_(G),SZJ_(G)), intent(in) :: h !< Layer thickness used to calculate fluxes, !! [H ~> m or kg m-2]. real, dimension(SZI_(G),SZJ_(G)), intent(in) :: h_L !< Left thickness in the reconstruction @@ -1456,8 +1456,8 @@ subroutine merid_face_thickness(v, h, h_L, h_R, h_v, dt, G, GV, US, LB, vol_CFL, real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), optional, intent(in) :: visc_rem_v !< Both the fraction !! of the momentum originally in a layer that remains after a time-step of !! viscosity, and the fraction of a time-step's worth of a barotropic - !! acceleration that a layer experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). + !! acceleration that a layer experiences after viscosity is applied [nondim]. + !! Visc_rem_v is between 0 (at the bottom) and 1 (far above the bottom). ! Local variables real :: CFL ! The CFL number based on the local velocity and grid spacing [nondim] @@ -1487,7 +1487,7 @@ subroutine merid_face_thickness(v, h, h_L, h_R, h_v, dt, G, GV, US, LB, vol_CFL, 3.0*curv_3*(CFL - 1.0)) else h_avg = 0.5 * (h_L(i,j+1,k) + h_R(i,j,k)) - ! The choice to use the arithmetic mean here is somewhat arbitrariy, but + ! The choice to use the arithmetic mean here is somewhat arbitrarily, but ! it should be noted that h_L(i+1,j,k) and h_R(i,j,k) are usually the same. h_marg = 0.5 * (h_L(i,j+1,k) + h_R(i,j,k)) ! h_marg = (2.0 * h_L(i,j+1,k) * h_R(i,j,k)) / & @@ -1556,8 +1556,8 @@ subroutine meridional_flux_adjust(v, h_in, h_L, h_R, vhbt, vh_tot_0, dvhdv_tot_0 !< Both the fraction of the momentum originally !! in a layer that remains after a time-step of viscosity, and the !! fraction of a time-step's worth of a barotropic acceleration that - !! a layer experiences after viscosity is applied. Non-dimensional - !! between 0 (at the bottom) and 1 (far above the bottom). + !! a layer experiences after viscosity is applied [nondim]. + !! Visc_rem is between 0 (at the bottom) and 1 (far above the bottom). real, dimension(SZI_(G)), & optional, intent(in) :: vhbt !< The summed volume flux through meridional faces !! [H L2 T-1 ~> m3 s-1 or kg s-1]. @@ -1719,9 +1719,9 @@ subroutine set_merid_BT_cont(v, h_in, h_L, h_R, BT_cont, vh_tot_0, dvhdv_tot_0, real, dimension(SZI_(G),SZK_(GV)), intent(in) :: visc_rem !< Both the fraction of the !! momentum originally in a layer that remains after a time-step !! of viscosity, and the fraction of a time-step's worth of a barotropic - !! acceleration that a layer experiences after viscosity is applied. - !! Non-dimensional between 0 (at the bottom) and 1 (far above the bottom). - real, dimension(SZI_(G)), intent(in) :: visc_rem_max !< Maximum allowable visc_rem. + !! acceleration that a layer experiences after viscosity is applied [nondim]. + !! Visc_rem is between 0 (at the bottom) and 1 (far above the bottom). + real, dimension(SZI_(G)), intent(in) :: visc_rem_max !< Maximum allowable visc_rem [nondim] integer, intent(in) :: j !< Spatial index. integer, intent(in) :: ish !< Start of index range. integer, intent(in) :: ieh !< End of index range. @@ -1755,7 +1755,7 @@ subroutine set_merid_BT_cont(v, h_in, h_L, h_R, BT_cont, vh_tot_0, dvhdv_tot_0, ! from leading to large CFL numbers. real :: min_visc_rem ! The smallest permitted value for visc_rem that is used ! in finding the barotropic velocity that changes the - ! flow direction. This is necessary to keep the inverse + ! flow direction [nondim]. This is necessary to keep the inverse ! of visc_rem from leading to large CFL numbers. real :: CFL_min ! A minimal increment in the CFL to try to ensure that the ! flow is truly upwind [nondim] @@ -1876,8 +1876,9 @@ subroutine PPM_reconstruction_x(h_in, h_L, h_R, G, LB, h_min, monotonic, simple_ ! Local variables with useful mnemonic names. real, dimension(SZI_(G),SZJ_(G)) :: slp ! The slopes. real, parameter :: oneSixth = 1./6. - real :: h_ip1, h_im1 - real :: dMx, dMn + real :: h_ip1, h_im1 ! Neighboring thicknesses or sensibly extrapolated values [H ~> m or kg m-2] + real :: dMx, dMn ! The difference between the local thickness and the maximum (dMx) or + ! minimum (dMn) of the surrounding values [H ~> m or kg m-2] character(len=256) :: mesg integer :: i, j, isl, iel, jsl, jel, n, stencil logical :: local_open_BC @@ -2011,8 +2012,9 @@ subroutine PPM_reconstruction_y(h_in, h_L, h_R, G, LB, h_min, monotonic, simple_ ! Local variables with useful mnemonic names. real, dimension(SZI_(G),SZJ_(G)) :: slp ! The slopes. real, parameter :: oneSixth = 1./6. - real :: h_jp1, h_jm1 - real :: dMx, dMn + real :: h_jp1, h_jm1 ! Neighboring thicknesses or sensibly extrapolated values [H ~> m or kg m-2] + real :: dMx, dMn ! The difference between the local thickness and the maximum (dMx) or + ! minimum (dMn) of the surrounding values [H ~> m or kg m-2] character(len=256) :: mesg integer :: i, j, isl, iel, jsl, jel, n, stencil logical :: local_open_BC @@ -2139,8 +2141,9 @@ subroutine PPM_limit_pos(h_in, h_L, h_R, h_min, G, iis, iie, jis, jie) integer, intent(in) :: jie !< End of j index range. ! Local variables - real :: curv, dh, scale - character(len=256) :: mesg + real :: curv ! The grid-normalized curvature of the three thicknesses [H ~> m or kg m-2] + real :: dh ! The difference between the edge thicknesses [H ~> m or kg m-2] + real :: scale ! A scaling factor to reduce the curvature of the fit [nondim] integer :: i,j do j=jis,jie ; do i=iis,iie @@ -2180,9 +2183,12 @@ subroutine PPM_limit_CW84(h_in, h_L, h_R, G, iis, iie, jis, jie) integer, intent(in) :: jie !< End of j index range. ! Local variables - real :: h_i, RLdiff, RLdiff2, RLmean, FunFac - character(len=256) :: mesg - integer :: i,j + real :: h_i ! A copy of the cell-average layer thickness [H ~> m or kg m-2] + real :: RLdiff ! The difference between the input edge values [H ~> m or kg m-2] + real :: RLdiff2 ! The squared difference between the input edge values [H2 ~> m2 or kg2 m-4] + real :: RLmean ! The average of the input edge thicknesses [H ~> m or kg m-2] + real :: FunFac ! A curious product of the thickness slope and curvature [H2 ~> m2 or kg2 m-4] + integer :: i, j do j=jis,jie ; do i=iis,iie ! This limiter monotonizes the parabola following diff --git a/src/core/MOM_density_integrals.F90 b/src/core/MOM_density_integrals.F90 index 3bdca94af3..8f26918253 100644 --- a/src/core/MOM_density_integrals.F90 +++ b/src/core/MOM_density_integrals.F90 @@ -11,7 +11,6 @@ module MOM_density_integrals use MOM_EOS, only : calculate_spec_vol use MOM_EOS, only : calculate_specific_vol_derivs use MOM_error_handler, only : MOM_error, FATAL, WARNING, MOM_mesg -use MOM_file_parser, only : get_param, log_version, param_file_type use MOM_hor_index, only : hor_index_type use MOM_string_functions, only : uppercase use MOM_variables, only : thermo_var_ptrs @@ -428,13 +427,13 @@ subroutine int_density_dz_generic_plm(k, tv, T_t, T_b, S_t, S_b, e, rho_ref, & real :: rho_scale ! A scaling factor for densities from kg m-3 to R [R m3 kg-1 ~> 1] real :: rho_ref_mks ! The reference density in MKS units, never rescaled from kg m-3 [kg m-3] real :: dz(HI%iscB:HI%iecB+1) ! Layer thicknesses at tracer points [Z ~> m] - real :: dz_x(5,HI%iscB:HI%iecB) ! Layer thicknesses along an x-line of subrid locations [Z ~> m] - real :: dz_y(5,HI%isc:HI%iec) ! Layer thicknesses along a y-line of subrid locations [Z ~> m] + real :: dz_x(5,HI%iscB:HI%iecB) ! Layer thicknesses along an x-line of subgrid locations [Z ~> m] + real :: dz_y(5,HI%isc:HI%iec) ! Layer thicknesses along a y-line of subgrid locations [Z ~> m] real :: massWeightToggle ! A non-dimensional toggle factor (0 or 1) [nondim] real :: Ttl, Tbl, Ttr, Tbr ! Temperatures at the velocity cell corners [degC] real :: Stl, Sbl, Str, Sbr ! Salinities at the velocity cell corners [ppt] real :: z0pres ! The height at which the pressure is zero [Z ~> m] - real :: hWght ! A topographically limited thicknes weight [Z ~> m] + real :: hWght ! A topographically limited thickness weight [Z ~> m] real :: hL, hR ! Thicknesses to the left and right [Z ~> m] real :: iDenom ! The denominator of the thickness weight expressions [Z-2 ~> m-2] logical :: use_stanley_eos ! True is SGS variance fields exist in tv. @@ -864,7 +863,7 @@ subroutine int_density_dz_generic_ppm(k, tv, T_t, T_b, S_t, S_b, e, & real :: T_top, T_mn, T_bot ! Left edge, cell mean and right edge values used in PPM reconstructions of T real :: S_top, S_mn, S_bot ! Left edge, cell mean and right edge values used in PPM reconstructions of S real :: z0pres ! The height at which the pressure is zero [Z ~> m] - real :: hWght ! A topographically limited thicknes weight [Z ~> m] + real :: hWght ! A topographically limited thickness weight [Z ~> m] real :: hL, hR ! Thicknesses to the left and right [Z ~> m] real :: iDenom ! The denominator of the thickness weight expressions [Z-2 ~> m-2] integer :: Isq, Ieq, Jsq, Jeq, i, j, m, n @@ -1455,7 +1454,10 @@ subroutine int_spec_vol_dp_generic_plm(T_t, T_b, S_t, S_b, p_t, p_b, alpha_ref, real :: p15(15) ! Pressures at fifteen quadrature points, scaled back to Pa as necessary [Pa] real :: a15(15) ! Specific volumes at fifteen quadrature points [R-1 ~> m3 kg-1] or [m3 kg-1] real :: wt_t(5), wt_b(5) ! Weights of top and bottom values at quadrature points [nondim] - real :: T_top, T_bot, S_top, S_bot, P_top, P_bot + real :: T_top, T_bot ! Horizontally interpolated temperature at the cell top and bottom [degC] + real :: S_top, S_bot ! Horizontally interpolated salinity at the cell top and bottom [ppt] + real :: P_top, P_bot ! Horizontally interpolated pressure at the cell top and bottom, + ! scaled back to Pa as necessary [Pa] real :: alpha_anom ! The depth averaged specific density anomaly [R-1 ~> m3 kg-1] or [m3 kg-1] real :: dp ! The pressure change through a layer [R L2 T-2 ~> Pa] diff --git a/src/core/MOM_dynamics_split_RK2.F90 b/src/core/MOM_dynamics_split_RK2.F90 index 1c045926de..d0a324f96f 100644 --- a/src/core/MOM_dynamics_split_RK2.F90 +++ b/src/core/MOM_dynamics_split_RK2.F90 @@ -148,15 +148,15 @@ module MOM_dynamics_split_RK2 !! dynamically. real :: be !< A nondimensional number from 0.5 to 1 that controls - !! the backward weighting of the time stepping scheme. + !! the backward weighting of the time stepping scheme [nondim] real :: begw !< A nondimensional number from 0 to 1 that controls !! the extent to which the treatment of gravity waves !! is forward-backward (0) or simulated backward - !! Euler (1). 0 is almost always used. + !! Euler (1) [nondim]. 0 is often used. logical :: debug !< If true, write verbose checksums for debugging purposes. logical :: debug_OBC !< If true, do debugging calls for open boundary conditions. - logical :: module_is_initialized = .false. !< Record whether this mouled has been initialzed. + logical :: module_is_initialized = .false. !< Record whether this module has been initialized. !>@{ Diagnostic IDs integer :: id_uh = -1, id_vh = -1 @@ -254,41 +254,41 @@ module MOM_dynamics_split_RK2 subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_surf_begin, p_surf_end, & uh, vh, uhtr, vhtr, eta_av, G, GV, US, CS, calc_dtbt, VarMix, & MEKE, thickness_diffuse_CSp, pbv, Waves) - type(ocean_grid_type), intent(inout) :: G !< ocean grid structure - type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure + type(ocean_grid_type), intent(inout) :: G !< Ocean grid structure + type(verticalGrid_type), intent(in) :: GV !< Ocean vertical grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), & - target, intent(inout) :: u !< zonal velocity [L T-1 ~> m s-1] + target, intent(inout) :: u !< Zonal velocity [L T-1 ~> m s-1] real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), & - target, intent(inout) :: v !< merid velocity [L T-1 ~> m s-1] + target, intent(inout) :: v !< Meridional velocity [L T-1 ~> m s-1] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), & - intent(inout) :: h !< layer thickness [H ~> m or kg m-2] - type(thermo_var_ptrs), intent(in) :: tv !< thermodynamic type - type(vertvisc_type), intent(inout) :: visc !< vertical visc, bottom drag, and related - type(time_type), intent(in) :: Time_local !< model time at end of time step - real, intent(in) :: dt !< time step [T ~> s] + intent(inout) :: h !< Layer thickness [H ~> m or kg m-2] + type(thermo_var_ptrs), intent(in) :: tv !< Thermodynamic type + type(vertvisc_type), intent(inout) :: visc !< Vertical visc, bottom drag, and related + type(time_type), intent(in) :: Time_local !< Model time at end of time step + real, intent(in) :: dt !< Baroclinic dynamics time step [T ~> s] type(mech_forcing), intent(in) :: forces !< A structure with the driving mechanical forces - real, dimension(:,:), pointer :: p_surf_begin !< surf pressure at the start of this dynamic + real, dimension(:,:), pointer :: p_surf_begin !< Surface pressure at the start of this dynamic !! time step [R L2 T-2 ~> Pa] - real, dimension(:,:), pointer :: p_surf_end !< surf pressure at the end of this dynamic + real, dimension(:,:), pointer :: p_surf_end !< Surface pressure at the end of this dynamic !! time step [R L2 T-2 ~> Pa] real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), & - target, intent(inout) :: uh !< zonal volume/mass transport + target, intent(inout) :: uh !< Zonal volume or mass transport !! [H L2 T-1 ~> m3 s-1 or kg s-1] real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), & - target, intent(inout) :: vh !< merid volume/mass transport + target, intent(inout) :: vh !< Meridional volume or mass transport !! [H L2 T-1 ~> m3 s-1 or kg s-1] real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)), & - intent(inout) :: uhtr !< accumulatated zonal volume/mass transport + intent(inout) :: uhtr !< Accumulated zonal volume or mass transport !! since last tracer advection [H L2 ~> m3 or kg] real, dimension(SZI_(G),SZJB_(G),SZK_(GV)), & - intent(inout) :: vhtr !< accumulatated merid volume/mass transport + intent(inout) :: vhtr !< Accumulated meridional volume or mass transport !! since last tracer advection [H L2 ~> m3 or kg] - real, dimension(SZI_(G),SZJ_(G)), intent(out) :: eta_av !< free surface height or column mass time + real, dimension(SZI_(G),SZJ_(G)), intent(out) :: eta_av !< Free surface height or column mass !! averaged over time step [H ~> m or kg m-2] - type(MOM_dyn_split_RK2_CS), pointer :: CS !< module control structure - logical, intent(in) :: calc_dtbt !< if true, recalculate barotropic time step - type(VarMix_CS), intent(inout) :: VarMix !< Variable mixing control struct + type(MOM_dyn_split_RK2_CS), pointer :: CS !< Module control structure + logical, intent(in) :: calc_dtbt !< If true, recalculate the barotropic time step + type(VarMix_CS), intent(inout) :: VarMix !< Variable mixing control structure type(MEKE_type), intent(inout) :: MEKE !< MEKE fields type(thickness_diffuse_CS), intent(inout) :: thickness_diffuse_CSp !< Pointer to a structure containing !! interface height diffusivities @@ -324,12 +324,22 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s real :: pres_to_eta ! A factor that converts pressures to the units of eta ! [H T2 R-1 L-2 ~> m Pa-1 or kg m-2 Pa-1] real, pointer, dimension(:,:) :: & - p_surf => NULL(), eta_PF_start => NULL(), & - taux_bot => NULL(), tauy_bot => NULL(), & - eta => NULL() + p_surf => NULL(), & ! A pointer to the surface pressure [R L2 T-2 ~> Pa] + eta_PF_start => NULL(), & ! The value of eta that corresponds to the starting pressure + ! for the barotropic solver [H ~> m or kg m-2] + taux_bot => NULL(), & ! A pointer to the zonal bottom stress in some cases [R L Z T-2 ~> Pa] + tauy_bot => NULL(), & ! A pointer to the meridional bottom stress in some cases [R L Z T-2 ~> Pa] + ! This pointer is just used as shorthand for CS%eta. + eta => NULL() ! A pointer to the instantaneous free surface height (in Boussinesq + ! mode) or column mass anomaly (in non-Boussinesq mode) [H ~> m or kg m-2] real, pointer, dimension(:,:,:) :: & - uh_ptr => NULL(), u_ptr => NULL(), vh_ptr => NULL(), v_ptr => NULL(), & + ! These pointers are used to alter which fields are passed to btstep with various options: + u_ptr => NULL(), & ! A pointer to a zonal velocity [L T-1] + v_ptr => NULL(), & ! A pointer to a meridional velocity [L T-1] + uh_ptr => NULL(), & ! A pointer to a zonal volume or mass transport [H L2 T-1 ~> m3 s-1 or kg s-1] + vh_ptr => NULL(), & ! A pointer to a meridional volume or mass transport [H L2 T-1 ~> m3 s-1 or kg s-1] + ! These pointers are just used as shorthand for CS%u_av, CS%v_av, and CS%h_av. u_av, & ! The zonal velocity time-averaged over a time step [L T-1 ~> m s-1]. v_av, & ! The meridional velocity time-averaged over a time step [L T-1 ~> m s-1]. h_av ! The layer thickness time-averaged over a time step [H ~> m or kg m-2]. @@ -339,12 +349,12 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s ! hf_CAu, hf_CAv, & ! Coriolis force accel. x fract. thickness [L T-2 ~> m s-2]. ! hf_u_BT_accel, hf_v_BT_accel ! barotropic correction accel. x fract. thickness [L T-2 ~> m s-2]. ! 3D diagnostics hf_PFu etc. are commented because there is no clarity on proper remapping grid option. - ! The code is retained for degugging purposes in the future. + ! The code is retained for debugging purposes in the future. real, allocatable, dimension(:,:) :: & - hf_PFu_2d, hf_PFv_2d, & ! Depth integeral of hf_PFu, hf_PFv [L T-2 ~> m s-2]. - hf_CAu_2d, hf_CAv_2d, & ! Depth integeral of hf_CAu, hf_CAv [L T-2 ~> m s-2]. - hf_u_BT_accel_2d, hf_v_BT_accel_2d ! Depth integeral of hf_u_BT_accel, hf_v_BT_accel + hf_PFu_2d, hf_PFv_2d, & ! Depth integral of hf_PFu, hf_PFv [L T-2 ~> m s-2]. + hf_CAu_2d, hf_CAv_2d, & ! Depth integral of hf_CAu, hf_CAv [L T-2 ~> m s-2]. + hf_u_BT_accel_2d, hf_v_BT_accel_2d ! Depth integral of hf_u_BT_accel, hf_v_BT_accel ! Diagnostics for thickness x momentum budget terms real, allocatable, dimension(:,:,:) :: & @@ -352,7 +362,7 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s h_CAu, h_CAv, & ! Coriolis force accel. x thickness [H L T-2 ~> m2 s-2]. h_u_BT_accel, h_v_BT_accel ! barotropic correction accel. x thickness [H L T-2 ~> m2 s-2]. - ! Dignostics for layer-sum of thickness x momentum budget terms + ! Diagnostics for layer-sum of thickness x momentum budget terms real, dimension(SZIB_(G),SZJ_(G)) :: & intz_PFu_2d, intz_CAu_2d, intz_u_BT_accel_2d ! [H L T-2 ~> m2 s-2]. real, dimension(SZI_(G),SZJB_(G)) :: & @@ -888,6 +898,9 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s ! The time-averaged free surface height has already been set by the last ! call to btstep. + ! Deallocate this memory to avoid a memory leak. ###We should also revisit how this array is declared. - RWH + !### if (dyn_p_surf .and. associated(eta_PF_start)) deallocate(eta_PF_start) + ! Here various terms used in to update the momentum equations are ! offered for time averaging. if (CS%id_PFu > 0) call post_data(CS%id_PFu, CS%PFu, CS%diag) @@ -922,7 +935,7 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s ! Diagnostics for terms multiplied by fractional thicknesses ! 3D diagnostics hf_PFu etc. are commented because there is no clarity on proper remapping grid option. - ! The code is retained for degugging purposes in the future. + ! The code is retained for debugging purposes in the future. !if (CS%id_hf_PFu > 0) then ! allocate(hf_PFu(G%IsdB:G%IedB,G%jsd:G%jed,GV%ke)) ! do k=1,nz ; do j=js,je ; do I=Isq,Ieq @@ -1173,18 +1186,18 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s end subroutine step_MOM_dyn_split_RK2 !> This subroutine sets up any auxiliary restart variables that are specific -!! to the unsplit time stepping scheme. All variables registered here should +!! to the split-explicit time stepping scheme. All variables registered here should !! have the ability to be recreated if they are not present in a restart file. subroutine register_restarts_dyn_split_RK2(HI, GV, param_file, CS, restart_CS, uh, vh) type(hor_index_type), intent(in) :: HI !< Horizontal index structure type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure type(param_file_type), intent(in) :: param_file !< parameter file type(MOM_dyn_split_RK2_CS), pointer :: CS !< module control structure - type(MOM_restart_CS), intent(inout) :: restart_CS !< MOM restart control struct + type(MOM_restart_CS), intent(inout) :: restart_CS !< MOM restart control structure real, dimension(SZIB_(HI),SZJ_(HI),SZK_(GV)), & - target, intent(inout) :: uh !< zonal volume/mass transport [H L2 T-1 ~> m3 s-1 or kg s-1] + target, intent(inout) :: uh !< zonal volume or mass transport [H L2 T-1 ~> m3 s-1 or kg s-1] real, dimension(SZI_(HI),SZJB_(HI),SZK_(GV)), & - target, intent(inout) :: vh !< merid volume/mass transport [H L2 T-1 ~> m3 s-1 or kg s-1] + target, intent(inout) :: vh !< merid volume or mass transport [H L2 T-1 ~> m3 s-1 or kg s-1] type(vardesc) :: vd(2) character(len=40) :: mdl = "MOM_dynamics_split_RK2" ! This module's name. @@ -1270,7 +1283,7 @@ subroutine initialize_dyn_split_RK2(u, v, h, uh, vh, eta, Time, G, GV, US, param type(param_file_type), intent(in) :: param_file !< parameter file for parsing type(diag_ctrl), target, intent(inout) :: diag !< to control diagnostics type(MOM_dyn_split_RK2_CS), pointer :: CS !< module control structure - type(MOM_restart_CS), intent(in) :: restart_CS !< MOM restart control struct + type(MOM_restart_CS), intent(in) :: restart_CS !< MOM restart control structure real, intent(in) :: dt !< time step [T ~> s] type(accel_diag_ptrs), target, intent(inout) :: Accel_diag !< points to momentum equation terms for !! budget analysis @@ -1280,7 +1293,7 @@ subroutine initialize_dyn_split_RK2(u, v, h, uh, vh, eta, Time, G, GV, US, param type(VarMix_CS), intent(inout) :: VarMix !< points to spatially variable viscosities type(MEKE_type), intent(inout) :: MEKE !< MEKE fields type(thickness_diffuse_CS), intent(inout) :: thickness_diffuse_CSp !< Pointer to the control structure - !! used for the isopycnal height diffusive transport. + !! used for the isopycnal height diffusive transport. type(ocean_OBC_type), pointer :: OBC !< points to OBC related fields type(update_OBC_CS), pointer :: update_OBC_CSp !< points to OBC update related fields type(ALE_CS), pointer :: ALE_CSp !< points to ALE control structure @@ -1296,19 +1309,19 @@ subroutine initialize_dyn_split_RK2(u, v, h, uh, vh, eta, Time, G, GV, US, param !! from the continuity solver. ! local variables - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_tmp + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_tmp ! A temporary copy of the layer thicknesses [H ~> m or kg m-2] character(len=40) :: mdl = "MOM_dynamics_split_RK2" ! This module's name. ! This include declares and sets the variable "version". # include "version_variable.h" character(len=48) :: thickness_units, flux_units, eta_rest_name - real :: H_rescale ! A rescaling factor for thicknesses from the representation in - ! a restart file to the internal representation in this run. - real :: vel_rescale ! A rescaling factor for velocities from the representation in - ! a restart file to the internal representation in this run. - real :: uH_rescale ! A rescaling factor for thickness transports from the representation in - ! a restart file to the internal representation in this run. - real :: accel_rescale ! A rescaling factor for accelerations from the representation in - ! a restart file to the internal representation in this run. + real :: H_rescale ! A rescaling factor for thicknesses from the representation in a + ! restart file to the internal representation in this run [various units ~> 1] + real :: vel_rescale ! A rescaling factor for velocities from the representation in a + ! restart file to the internal representation in this run [various units ~> 1] + real :: uH_rescale ! A rescaling factor for thickness transports from the representation in a + ! restart file to the internal representation in this run [various units ~> 1] + real :: accel_rescale ! A rescaling factor for accelerations from the representation in a + ! restart file to the internal representation in this run [various units ~> 1] type(group_pass_type) :: pass_av_h_uvh logical :: use_tides, debug_truncations diff --git a/src/core/MOM_dynamics_unsplit.F90 b/src/core/MOM_dynamics_unsplit.F90 index 6edf5aeacc..88a11e071c 100644 --- a/src/core/MOM_dynamics_unsplit.F90 +++ b/src/core/MOM_dynamics_unsplit.F90 @@ -12,7 +12,7 @@ module MOM_dynamics_unsplit !* Runge-Kutta time stepping scheme for the momentum and a forward- * !* backward coupling between the momentum and continuity equations. * !* This was the orignal unsplit time stepping scheme used in early * -!* versions of HIM and its precuror. While it is very simple and * +!* versions of HIM and its precursor. While it is very simple and * !* accurate, it is much less efficient that the split time stepping * !* scheme for realistic oceanographic applications. It has been * !* retained for all of these years primarily to verify that the split * @@ -121,7 +121,7 @@ module MOM_dynamics_unsplit !! for viscosity. The default should be true, but it is false. logical :: debug !< If true, write verbose checksums for debugging purposes. - logical :: module_is_initialized = .false. !< Record whether this mouled has been initialzed. + logical :: module_is_initialized = .false. !< Record whether this module has been initialized. !>@{ Diagnostic IDs integer :: id_uh = -1, id_vh = -1 @@ -215,19 +215,19 @@ subroutine step_MOM_dyn_unsplit(u, v, h, tv, visc, Time_local, dt, forces, & !! column mass [H ~> m or kg m-2]. type(MOM_dyn_unsplit_CS), pointer :: CS !< The control structure set up by !! initialize_dyn_unsplit. - type(VarMix_CS), intent(inout) :: VarMix !< Variable mixing control struct + type(VarMix_CS), intent(inout) :: VarMix !< Variable mixing control structure type(MEKE_type), intent(inout) :: MEKE !< MEKE fields type(porous_barrier_ptrs), intent(in) :: pbv !< porous barrier fractional cell metrics type(wave_parameters_CS), optional, pointer :: Waves !< A pointer to a structure containing !! fields related to the surface wave conditions ! Local variables - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_av, hp ! Prediced or averaged layer thicknesses [H ~> m or kg m-2] + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_av, hp ! Predicted or averaged layer thicknesses [H ~> m or kg m-2] real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)) :: up, upp ! Predicted zonal velocities [L T-1 ~> m s-1] real, dimension(SZI_(G),SZJB_(G),SZK_(GV)) :: vp, vpp ! Predicted meridional velocities [L T-1 ~> m s-1] real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)) :: ueffA ! Effective Area of U-Faces [H L ~> m2] real, dimension(SZI_(G),SZJB_(G),SZK_(GV)) :: veffA ! Effective Area of V-Faces [H L ~> m2] - real, dimension(:,:), pointer :: p_surf => NULL() + real, dimension(:,:), pointer :: p_surf => NULL() ! A pointer to the surface pressure [R L2 T-2 ~> Pa] real :: dt_pred ! The time step for the predictor part of the baroclinic time stepping [T ~> s]. real :: dt_visc ! The time step for a part of the update due to viscosity [T ~> s]. logical :: dyn_p_surf @@ -603,7 +603,7 @@ subroutine initialize_dyn_unsplit(u, v, h, Time, G, GV, US, param_file, diag, CS !! the appropriate control structure. type(ALE_CS), pointer :: ALE_CSp !< This points to the ALE control !! structure. - type(set_visc_CS), target, intent(in) :: set_visc !< set_visc control struct + type(set_visc_CS), target, intent(in) :: set_visc !< set_visc control structure type(vertvisc_type), intent(inout) :: visc !< A structure containing vertical !! viscosities, bottom drag !! viscosities, and related fields. diff --git a/src/core/MOM_dynamics_unsplit_RK2.F90 b/src/core/MOM_dynamics_unsplit_RK2.F90 index 862082d567..26bd00aaf5 100644 --- a/src/core/MOM_dynamics_unsplit_RK2.F90 +++ b/src/core/MOM_dynamics_unsplit_RK2.F90 @@ -113,17 +113,17 @@ module MOM_dynamics_unsplit_RK2 !! to the seafloor [R L Z T-2 ~> Pa] real :: be !< A nondimensional number from 0.5 to 1 that controls - !! the backward weighting of the time stepping scheme. + !! the backward weighting of the time stepping scheme [nondim]. real :: begw !< A nondimensional number from 0 to 1 that controls !! the extent to which the treatment of gravity waves !! is forward-backward (0) or simulated backward - !! Euler (1). 0 is almost always used. + !! Euler (1) [nondim]. 0 is often used. logical :: use_correct_dt_visc !< If true, use the correct timestep in the calculation of the !! turbulent mixed layer properties for viscosity. !! The default should be true, but it is false. logical :: debug !< If true, write verbose checksums for debugging purposes. - logical :: module_is_initialized = .false. !< Record whether this mouled has been initialzed. + logical :: module_is_initialized = .false. !< Record whether this module has been initialized. !>@{ Diagnostic IDs integer :: id_uh = -1, id_vh = -1 @@ -226,18 +226,19 @@ subroutine step_MOM_dyn_unsplit_RK2(u_in, v_in, h_in, tv, visc, Time_local, dt, !! or column mass [H ~> m or kg m-2]. type(MOM_dyn_unsplit_RK2_CS), pointer :: CS !< The control structure set up by !! initialize_dyn_unsplit_RK2. - type(VarMix_CS), intent(inout) :: VarMix !< Variable mixing control struct + type(VarMix_CS), intent(inout) :: VarMix !< Variable mixing control structure type(MEKE_type), intent(inout) :: MEKE !< MEKE fields !! fields related to the Mesoscale !! Eddy Kinetic Energy. type(porous_barrier_ptrs), intent(in) :: pbv !< porous barrier fractional cell metrics ! Local variables - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_av, hp + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_av ! Averaged layer thicknesses [H ~> m or kg m-2] + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: hp ! Predicted layer thicknesses [H ~> m or kg m-2] real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)) :: up ! Predicted zonal velocities [L T-1 ~> m s-1] real, dimension(SZI_(G),SZJB_(G),SZK_(GV)) :: vp ! Predicted meridional velocities [L T-1 ~> m s-1] real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)) :: ueffA ! Effective Area of U-Faces [H L ~> m2] real, dimension(SZI_(G),SZJB_(G),SZK_(GV)) :: veffA ! Effective Area of V-Faces [H L ~> m2] - real, dimension(:,:), pointer :: p_surf => NULL() + real, dimension(:,:), pointer :: p_surf => NULL() ! A pointer to the surface pressure [R L2 T-2 ~> Pa] real :: dt_pred ! The time step for the predictor part of the baroclinic time stepping [T ~> s] real :: dt_visc ! The time step for a part of the update due to viscosity [T ~> s] logical :: dyn_p_surf @@ -548,7 +549,7 @@ subroutine initialize_dyn_unsplit_RK2(u, v, h, Time, G, GV, US, param_file, diag !! to the appropriate control structure. type(ALE_CS), pointer :: ALE_CSp !< This points to the ALE !! control structure. - type(set_visc_CS), target, intent(in) :: set_visc !< set visc control struct + type(set_visc_CS), target, intent(in) :: set_visc !< set visc control structure type(vertvisc_type), intent(inout) :: visc !< A structure containing !! vertical viscosities, bottom drag !! viscosities, and related fields. diff --git a/src/core/MOM_forcing_type.F90 b/src/core/MOM_forcing_type.F90 index 9cec8587db..5d9d319a49 100644 --- a/src/core/MOM_forcing_type.F90 +++ b/src/core/MOM_forcing_type.F90 @@ -104,8 +104,11 @@ module MOM_forcing_type lrunoff => NULL(), & !< liquid river runoff entering ocean [R Z T-1 ~> kg m-2 s-1] frunoff => NULL(), & !< frozen river runoff (calving) entering ocean [R Z T-1 ~> kg m-2 s-1] seaice_melt => NULL(), & !< snow/seaice melt (positive) or formation (negative) [R Z T-1 ~> kg m-2 s-1] - netMassIn => NULL(), & !< Sum of water mass flux out of the ocean [kg m-2 s-1] - netMassOut => NULL(), & !< Net water mass flux into of the ocean [kg m-2 s-1] + netMassIn => NULL(), & !< Sum of water mass flux out of the ocean integrated over a + !! forcing timestep [H ~> m or kg m-2] + netMassOut => NULL(), & !< Net water mass flux into of the ocean integrated over a + !! forcing timestep [H ~> m or kg m-2] + !### Net salt is used with inconsistent units and only in one place and should be eliminated as unneeded. netSalt => NULL() !< Net salt entering the ocean [kgSalt m-2 s-1] ! heat associated with water crossing ocean surface @@ -152,14 +155,14 @@ module MOM_forcing_type ! iceberg related inputs real, pointer, dimension(:,:) :: & ustar_berg => NULL(), & !< iceberg contribution to top ustar [Z T-1 ~> m s-1]. - area_berg => NULL(), & !< area of ocean surface covered by icebergs [m2 m-2] + area_berg => NULL(), & !< fractional area of ocean surface covered by icebergs [nondim] mass_berg => NULL() !< mass of icebergs [R Z ~> kg m-2] ! land ice-shelf related inputs real, pointer, dimension(:,:) :: ustar_shelf => NULL() !< Friction velocity under ice-shelves [Z T-1 ~> m s-1]. !! as computed by the ocean at the previous time step. real, pointer, dimension(:,:) :: frac_shelf_h => NULL() !< Fractional ice shelf coverage of - !! h-cells, nondimensional from 0 to 1. This is only + !! h-cells, from 0 to 1 [nondim]. This is only !! associated if ice shelves are enabled, and are !! exactly 0 away from shelves or on land. real, pointer, dimension(:,:) :: iceshelf_melt => NULL() !< Ice shelf melt rate (positive) @@ -177,7 +180,7 @@ module MOM_forcing_type !! fluxes have been applied to the ocean. real :: dt_buoy_accum = -1.0 !< The amount of time over which the buoyancy fluxes !! should be applied [T ~> s]. If negative, this forcing - !! type variable has not yet been inialized. + !! type variable has not yet been initialized. logical :: gustless_accum_bug = .true. !< If true, use an incorrect expression in the time !! average of the gustless wind stress. real :: C_p !< heat capacity of seawater [Q degC-1 ~> J kg-1 degC-1]. @@ -231,7 +234,7 @@ module MOM_forcing_type ! iceberg related inputs real, pointer, dimension(:,:) :: & - area_berg => NULL(), & !< fractional area of ocean surface covered by icebergs [m2 m-2] + area_berg => NULL(), & !< fractional area of ocean surface covered by icebergs [nondim] mass_berg => NULL() !< mass of icebergs per unit ocean area [R Z ~> kg m-2] ! land ice-shelf related inputs @@ -257,15 +260,15 @@ module MOM_forcing_type !! ice needs to be accumulated, and the rigidity explicitly !! reset to zero at the driver level when appropriate. real, pointer, dimension(:,:) :: & - ustk0 => NULL(), & !< Surface Stokes drift, zonal [m/s] - vstk0 => NULL() !< Surface Stokes drift, meridional [m/s] + ustk0 => NULL(), & !< Surface Stokes drift, zonal [m s-1] + vstk0 => NULL() !< Surface Stokes drift, meridional [m s-1] real, pointer, dimension(:) :: & - stk_wavenumbers => NULL() !< The central wave number of Stokes bands [rad/m] + stk_wavenumbers => NULL() !< The central wave number of Stokes bands [rad m-1] real, pointer, dimension(:,:,:) :: & - ustkb => NULL(), & !< Stokes Drift spectrum, zonal [m/s] + ustkb => NULL(), & !< Stokes Drift spectrum, zonal [m s-1] !! Horizontal - u points !! 3rd dimension - wavenumber - vstkb => NULL() !< Stokes Drift spectrum, meridional [m/s] + vstkb => NULL() !< Stokes Drift spectrum, meridional [m s-1] !! Horizontal - v points !! 3rd dimension - wavenumber @@ -460,7 +463,7 @@ subroutine extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, dt, & ! [degC H T-1 ~> degC m s-1 or degC kg m-2 s-1] real :: Ih_limit ! inverse depth at which surface fluxes start to be limited ! or 0 for no limiting [H-1 ~> m-1 or m2 kg-1] - real :: scale ! scale scales away fluxes if depth < FluxRescaleDepth + real :: scale ! scale scales away fluxes if depth < FluxRescaleDepth [nondim] real :: I_Cp ! 1.0 / C_p [degC Q-1 ~> kg degC J-1] real :: I_Cp_Hconvert ! Unit conversion factors divided by the heat capacity ! [degC H R-1 Z-1 Q-1 ~> degC m3 J-1 or kg degC J-1] @@ -946,23 +949,21 @@ subroutine calculateBuoyancyFlux1d(G, GV, US, fluxes, optics, nsw, h, Temp, Salt logical :: useRiverHeatContent logical :: useCalvingHeatContent - real :: depthBeforeScalingFluxes ! A depth scale [H ~> m or kg m-2] real :: GoRho ! The gravitational acceleration divided by mean density times a ! unit conversion factor [L2 H-1 R-1 T-2 ~> m4 kg-1 s-2 or m7 kg-2 s-2] - real :: H_limit_fluxes ! Another depth scale [H ~> m or kg m-2] + real :: H_limit_fluxes ! A depth scale that specifies when the ocean is shallow that + ! it is necessary to eliminate fluxes [H ~> m or kg m-2] integer :: i, k ! smg: what do we do when have heat fluxes from calving and river? useRiverHeatContent = .False. useCalvingHeatContent = .False. - depthBeforeScalingFluxes = max( GV%Angstrom_H, 1.e-30*GV%m_to_H ) + H_limit_fluxes = max( GV%Angstrom_H, 1.e-30*GV%m_to_H ) pressure(:) = 0. if (associated(tv%p_surf)) then ; do i=G%isc,G%iec ; pressure(i) = tv%p_surf(i,j) ; enddo ; endif GoRho = (GV%g_Earth * GV%H_to_Z) / GV%Rho0 - H_limit_fluxes = depthBeforeScalingFluxes - ! The surface forcing is contained in the fluxes type. ! We aggregate the thermodynamic forcing for a time step into the following: ! netH = water added/removed via surface fluxes [H T-1 ~> m s-1 or kg m-2 s-1] @@ -971,7 +972,7 @@ subroutine calculateBuoyancyFlux1d(G, GV, US, fluxes, optics, nsw, h, Temp, Salt ! Note that unlike other calls to extractFLuxes1d() that return the time-integrated flux ! this call returns the rate because dt=1 (in arbitrary time units) call extractFluxes1d(G, GV, US, fluxes, optics, nsw, j, 1.0, & - depthBeforeScalingFluxes, useRiverHeatContent, useCalvingHeatContent, & + H_limit_fluxes, useRiverHeatContent, useCalvingHeatContent, & h(:,j,:), Temp(:,j,:), netH, netEvap, netHeatMinusSW, & netSalt, penSWbnd, tv, .false.) @@ -1421,11 +1422,10 @@ subroutine register_forcing_type_diags(Time, diag, US, use_temperature, handles, handles%id_massout_flux = register_diag_field('ocean_model', 'massout_flux', diag%axesT1, Time, & 'Net mass flux of freshwater out of the ocean (used in the boundary flux calculation)', & 'kg m-2', conversion=diag%GV%H_to_kg_m2) - ! This diagnostic is calculated in MKS units. handles%id_massin_flux = register_diag_field('ocean_model', 'massin_flux', diag%axesT1, Time, & - 'Net mass flux of freshwater into the ocean (used in boundary flux calculation)', 'kg m-2') - ! This diagnostic is calculated in MKS units. + 'Net mass flux of freshwater into the ocean (used in boundary flux calculation)', & + 'kg m-2', conversion=diag%GV%H_to_kg_m2) !========================================================================= ! area integrated surface mass transport, all are rescaled to MKS units before area integration. @@ -1981,12 +1981,12 @@ subroutine fluxes_accumulate(flux_tmp, fluxes, G, wt2, forces) real, intent(out) :: wt2 !< The relative weight of the new fluxes type(mech_forcing), optional, intent(in) :: forces !< A structure with the driving mechanical forces - ! This subroutine copies mechancal forcing from flux_tmp to fluxes and + ! This subroutine copies mechanical forcing from flux_tmp to fluxes and ! stores the time-weighted averages of the various buoyancy fluxes in fluxes, ! and increments the amount of time over which the buoyancy forcing in fluxes should be ! applied based on the time interval stored in flux_tmp. - real :: wt1 + real :: wt1 ! The relative weight of the previous fluxes [nondim] integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq, i0, j0 integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, isr, ier, jsr, jer is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec @@ -2342,17 +2342,18 @@ subroutine forcing_diagnostics(fluxes_in, sfc_state, G_in, US, time_end, diag, h type(diag_ctrl), intent(inout) :: diag !< diagnostic regulator type(forcing_diags), intent(inout) :: handles !< diagnostic ids - ! local + ! local variables type(ocean_grid_type), pointer :: G ! Grid metric on model index map type(forcing), pointer :: fluxes ! Fluxes on the model index map - real, dimension(SZI_(diag%G),SZJ_(diag%G)) :: res - real :: total_transport ! for diagnosing integrated boundary transport - real :: ave_flux ! for diagnosing averaged boundary flux + real, dimension(SZI_(diag%G),SZJ_(diag%G)) :: res ! A temporary array for rescaled combinations + ! of fluxes in MKS units, like [kg m-2 s-1] or [W m-2] + real :: total_transport ! for diagnosing integrated boundary transport, in MKS units like [kg s-1] or [W] + real :: ave_flux ! for diagnosing averaged boundary flux, in MKS units like [kg m-2 s-1] or [W m-2] real :: RZ_T_conversion ! A combination of scaling factors for mass fluxes [kg T m-2 s-1 R-1 Z-1 ~> 1] real :: I_dt ! inverse time step [T-1 ~> s-1] - real :: ppt2mks ! conversion between ppt and mks + real :: ppt2mks ! conversion between ppt and mks units [nondim] integer :: turns ! Number of index quarter turns - integer :: i,j,is,ie,js,je + integer :: i, j, is, ie, js, je call cpu_clock_begin(handles%id_clock_forcing) @@ -3306,8 +3307,8 @@ end subroutine deallocate_mech_forcing !< Rotate the fluxes by a set number of quarter turns subroutine rotate_forcing(fluxes_in, fluxes, turns) - type(forcing), intent(in) :: fluxes_in !< Input forcing struct - type(forcing), intent(inout) :: fluxes !< Rotated forcing struct + type(forcing), intent(in) :: fluxes_in !< Input forcing structure + type(forcing), intent(inout) :: fluxes !< Rotated forcing structure integer, intent(in) :: turns !< Number of quarter turns logical :: do_ustar, do_water, do_heat, do_salt, do_press, do_shelf, & @@ -3504,7 +3505,7 @@ end subroutine rotate_mech_forcing !! \subsection subsection_mass_fluxes Surface boundary mass fluxes !! !! The ocean gains or loses mass through evaporation, precipitation, -!! sea ice melt/form, and and river runoff. Positive mass fluxes +!! sea ice melt/form, and river runoff. Positive mass fluxes !! add mass to the liquid ocean. The boundary mass flux units are !! (kilogram per square meter per sec: kg/(m2/sec)). !! diff --git a/src/core/MOM_grid.F90 b/src/core/MOM_grid.F90 index 1e63d033f2..d9ed8ffee4 100644 --- a/src/core/MOM_grid.F90 +++ b/src/core/MOM_grid.F90 @@ -59,8 +59,8 @@ module MOM_grid integer :: JsgB !< The start j-index of cell vertices within the global domain integer :: JegB !< The end j-index of cell vertices within the global domain - integer :: isd_global !< The value of isd in the global index space (decompoistion invariant). - integer :: jsd_global !< The value of isd in the global index space (decompoistion invariant). + integer :: isd_global !< The value of isd in the global index space (decomposition invariant). + integer :: jsd_global !< The value of isd in the global index space (decomposition invariant). integer :: idg_offset !< The offset between the corresponding global and local i-indices. integer :: jdg_offset !< The offset between the corresponding global and local j-indices. integer :: ke !< The number of layers in the vertical. @@ -206,7 +206,7 @@ subroutine MOM_grid_init(G, param_file, US, HI, global_indexing, bathymetry_at_v !! are entirely determined from thickness points. ! Local variables - real :: mean_SeaLev_scale + real :: mean_SeaLev_scale ! A scaling factor for the reference height variable [1] or [Z m-1 ~> 1] integer :: isd, ied, jsd, jed, nk integer :: IsdB, IedB, JsdB, JedB integer :: ied_max, jed_max @@ -398,10 +398,10 @@ end subroutine MOM_grid_init subroutine rescale_grid_bathymetry(G, m_in_new_units) type(ocean_grid_type), intent(inout) :: G !< The horizontal grid structure real, intent(in) :: m_in_new_units !< The new internal representation of 1 m depth. - ! It appears that this routine is never called. + !### It appears that this routine is never called. ! Local variables - real :: rescale + real :: rescale ! A unit rescaling factor [various combinations of units ~> 1] integer :: i, j, isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed @@ -489,14 +489,16 @@ logical function isPointInCell(G, i, j, x, y) real, intent(in) :: x !< x coordinate of point real, intent(in) :: y !< y coordinate of point ! Local variables - real :: xNE, xNW, xSE, xSW, yNE, yNW, ySE, ySW - real :: p0, p1, p2, p3, l0, l1, l2, l3 + real :: xNE, xNW, xSE, xSW ! Longitudes of cell corners [degLon] + real :: yNE, yNW, ySE, ySW ! Latitudes of cell corners [degLat] + real :: l0, l1, l2, l3 ! Crossed products of differences in position [degLon degLat] + real :: p0, p1, p2, p3 ! Trinary unitary values reflecting the signs of the crossed products [nondim] isPointInCell = .false. xNE = G%geoLonBu(i ,j ) ; yNE = G%geoLatBu(i ,j ) xNW = G%geoLonBu(i-1,j ) ; yNW = G%geoLatBu(i-1,j ) xSE = G%geoLonBu(i ,j-1) ; ySE = G%geoLatBu(i ,j-1) xSW = G%geoLonBu(i-1,j-1) ; ySW = G%geoLatBu(i-1,j-1) - ! This is a crude calculation that assume a geographic coordinate system + ! This is a crude calculation that assumes a geographic coordinate system if (xmax(xNE,xNW,xSE,xSW) .or. & ymax(yNE,yNW,ySE,ySW) ) then return ! Avoid the more complicated calculation diff --git a/src/core/MOM_interface_heights.F90 b/src/core/MOM_interface_heights.F90 index 17729e586c..7047dd6421 100644 --- a/src/core/MOM_interface_heights.F90 +++ b/src/core/MOM_interface_heights.F90 @@ -17,7 +17,7 @@ module MOM_interface_heights public find_eta -!> Calculates the heights of sruface or all interfaces from layer thicknesses. +!> Calculates the heights of the free surface or all interfaces from layer thicknesses. interface find_eta module procedure find_eta_2d, find_eta_3d end interface find_eta diff --git a/src/core/MOM_isopycnal_slopes.F90 b/src/core/MOM_isopycnal_slopes.F90 index 40777c8227..80d94ec7fe 100644 --- a/src/core/MOM_isopycnal_slopes.F90 +++ b/src/core/MOM_isopycnal_slopes.F90 @@ -102,15 +102,15 @@ subroutine calc_isoneutral_slopes(G, GV, US, h, e, tv, dt_kappa_smooth, & real :: h_neglect ! A thickness that is so small it is usually lost ! in roundoff and can be neglected [H ~> m or kg m-2]. real :: h_neglect2 ! h_neglect^2 [H2 ~> m2 or kg2 m-4]. - real :: dz_neglect ! A change in interface heighs that is so small it is usually lost + real :: dz_neglect ! A change in interface heights that is so small it is usually lost ! in roundoff and can be neglected [Z ~> m]. logical :: use_EOS ! If true, density is calculated from T & S using an equation of state. real :: G_Rho0 ! The gravitational acceleration divided by density [L2 Z-1 T-2 R-1 ~> m4 s-2 kg-1] real :: Z_to_L ! A conversion factor between from units for e to the - ! units for lateral distances. + ! units for lateral distances [L Z-1 ~> 1] real :: L_to_Z ! A conversion factor between from units for lateral distances - ! to the units for e. - real :: H_to_Z ! A conversion factor from thickness units to the units of e. + ! to the units for e [Z L-1 ~> 1] + real :: H_to_Z ! A conversion factor from thickness units to the units of e [Z H-1 ~> 1 or m3 kg-1] logical :: present_N2_u, present_N2_v integer, dimension(2) :: EOSdom_u, EOSdom_v ! Domains for the equation of state calculations at u and v points @@ -457,7 +457,7 @@ subroutine vert_fill_TS(h, T_in, S_in, kappa_dt, T_f, S_f, G, GV, halo_here, lar real :: kap_dt_x2 ! The 2*kappa_dt converted to H units [H2 ~> m2 or kg2 m-4]. real :: h_neglect ! A negligible thickness [H ~> m or kg m-2], to allow for zero thicknesses. real :: h0 ! A negligible thickness to allow for zero thickness layers without - ! completely decouping groups of layers [H ~> m or kg m-2]. + ! completely decoupling groups of layers [H ~> m or kg m-2]. ! Often 0 < h_neglect << h0. real :: h_tr ! h_tr is h at tracer points with a tiny thickness ! added to ensure positive definiteness [H ~> m or kg m-2]. diff --git a/src/core/MOM_open_boundary.F90 b/src/core/MOM_open_boundary.F90 index cb9422a412..6d8696216a 100644 --- a/src/core/MOM_open_boundary.F90 +++ b/src/core/MOM_open_boundary.F90 @@ -203,8 +203,8 @@ module MOM_open_boundary type(segment_tracer_registry_type), pointer :: tr_Reg=> NULL()!< A pointer to the tracer registry for the segment. type(hor_index_type) :: HI !< Horizontal index ranges real :: Tr_InvLscale_out !< An effective inverse length scale for restoring - !! the tracer concentration in a ficticious - !! reservior towards interior values when flow + !! the tracer concentration in a fictitious + !! reservoir towards interior values when flow !! is exiting the domain [L-1 ~> m-1] real :: Tr_InvLscale_in !< An effective inverse length scale for restoring !! the tracer concentration towards an externally @@ -283,7 +283,7 @@ module MOM_open_boundary ! The following parameters are used in the baroclinic radiation code: real :: gamma_uv !< The relative weighting for the baroclinic radiation !! velocities (or speed of characteristics) at the - !! new time level (1) or the running mean (0) for velocities. + !! new time level (1) or the running mean (0) for velocities [nondim]. !! Valid values range from 0 to 1, with a default of 0.3. real :: rx_max !< The maximum magnitude of the baroclinic radiation velocity (or speed of !! characteristics) in units of grid points per timestep [nondim]. @@ -303,13 +303,12 @@ module MOM_open_boundary !! the independence of the OBCs to this external data [H ~> m or kg m-2]. real :: silly_u !< A silly value of velocity outside of the domain that can be used to test !! the independence of the OBCs to this external data [L T-1 ~> m s-1]. - logical :: ramp = .false. !< If True, ramp from zero to the external values - !! for SSH. + logical :: ramp = .false. !< If True, ramp from zero to the external values for SSH. logical :: ramping_is_activated = .false. !< True if the ramping has been initialized - real :: ramp_timescale !< If ramp is True, use this timescale for ramping. - real :: trunc_ramp_time !< If ramp is True, time after which ramp is done. + real :: ramp_timescale !< If ramp is True, use this timescale for ramping [s]. + real :: trunc_ramp_time !< If ramp is True, time after which ramp is done [s]. real :: ramp_value !< If ramp is True, where we are on the ramp from - !! zero to one. + !! zero to one [nondim]. type(time_type) :: ramp_start_time !< Time when model was started. end type ocean_OBC_type @@ -653,7 +652,7 @@ subroutine initialize_segment_data(G, OBC, PF) character(len=256) :: filename character(len=20) :: segnam, suffix character(len=32) :: varnam, fieldname - real :: value + real :: value ! A value that is parsed from the segment data string [various units] character(len=32), dimension(MAX_OBC_FIELDS) :: fields ! segment field names character(len=128) :: inputdir type(OBC_segment_type), pointer :: segment => NULL() ! pointer to segment type list @@ -896,7 +895,7 @@ subroutine initialize_segment_data(G, OBC, PF) segment%field(m)%value = value segment%field(m)%name = trim(fields(m)) ! Check if this is a tidal field. If so, the number - ! of expected constiuents must be 1. + ! of expected constituents must be 1. if ((index(segment%field(m)%name, 'phase') > 0) .or. (index(segment%field(m)%name, 'amp') > 0)) then if (OBC%n_tide_constituents .gt. 1 .and. OBC%add_tide_constituents) then call MOM_error(FATAL, 'Only one constituent is supported when specifying '//& @@ -1025,7 +1024,7 @@ subroutine initialize_obc_tides(OBC, param_file) " are true, or if OBC_TIDE_N_CONSTITUENTS > 0 and "//trim(OBC%tide_names(c))//& " is in OBC_TIDE_CONSTITUENTS.", units="s-1", default=tidal_frequency(trim(OBC%tide_names(c)))) - ! Find equilibrum phase if needed + ! Find equilibrium phase if needed if (OBC%add_eq_phase) then OBC%tide_eq_phases(c) = eq_phase(trim(OBC%tide_names(c)), OBC%tidal_longitudes) else @@ -1179,7 +1178,7 @@ subroutine setup_u_point_obc(OBC, G, US, segment_str, l_seg, PF, reentrant_y) integer :: j, a_loop character(len=32) :: action_str(8) character(len=128) :: segment_param_str - real, allocatable, dimension(:) :: tnudge + real, allocatable, dimension(:) :: tnudge ! Nudging timescales [T ~> s] ! This returns the global indices for the segment call parse_segment_str(G%ieg, G%jeg, segment_str, I_obc, Js_obc, Je_obc, action_str, reentrant_y) @@ -1319,7 +1318,7 @@ subroutine setup_v_point_obc(OBC, G, US, segment_str, l_seg, PF, reentrant_x) integer :: i, a_loop character(len=32) :: action_str(8) character(len=128) :: segment_param_str - real, allocatable, dimension(:) :: tnudge + real, allocatable, dimension(:) :: tnudge ! Nudging timescales [T ~> s] ! This returns the global indices for the segment call parse_segment_str(G%ieg, G%jeg, segment_str, J_obc, Is_obc, Ie_obc, action_str, reentrant_x) @@ -1470,7 +1469,7 @@ subroutine parse_segment_str(ni_global, nj_global, segment_str, l, m, n, action_ mn_max = nj_global if (.not. (word2(1:2)=='J=')) call MOM_error(FATAL, "MOM_open_boundary.F90, parse_segment_str: "//& "Second word of string '"//trim(segment_str)//"' must start with 'J='.") - elseif (word1(1:2)=='J=') then ! Note that the file_parser uniformaly expands "=" to " = " + elseif (word1(1:2)=='J=') then ! Note that the file_parser uniformly expands "=" to " = " l_max = nj_global mn_max = ni_global if (.not. (word2(1:2)=='I=')) call MOM_error(FATAL, "MOM_open_boundary.F90, parse_segment_str: "//& @@ -1640,7 +1639,7 @@ subroutine parse_for_tracer_reservoirs(OBC, PF, use_temperature) character(len=256) :: filename character(len=20) :: segnam, suffix character(len=32) :: varnam, fieldname - real :: value + real :: value ! A value that is parsed from the segment data string [various units] character(len=32), dimension(MAX_OBC_FIELDS) :: fields ! segment field names type(OBC_segment_type), pointer :: segment => NULL() ! pointer to segment type list character(len=256) :: mesg ! Message for error messages. @@ -1910,7 +1909,7 @@ subroutine open_boundary_end(OBC) call open_boundary_dealloc(OBC) end subroutine open_boundary_end -!> Sets the slope of bathymetry normal to an open bounndary to zero. +!> Sets the slope of bathymetry normal to an open boundary to zero. subroutine open_boundary_impose_normal_slope(OBC, G, depth) type(ocean_OBC_type), pointer :: OBC !< Open boundary control structure type(dyn_horgrid_type), intent(in) :: G !< Ocean grid structure @@ -2236,7 +2235,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, do k=1,nz ; do j=segment%HI%jsd,segment%HI%jed if (segment%radiation) then dhdt = (u_old(I-1,j,k) - u_new(I-1,j,k)) !old-new - dhdx = (u_new(I-1,j,k) - u_new(I-2,j,k)) !in new time backward sasha for I-1 + dhdx = (u_new(I-1,j,k) - u_new(I-2,j,k)) !in new time backward sashay for I-1 rx_new = 0.0 if (dhdt*dhdx > 0.0) rx_new = min( (dhdt/dhdx), rx_max) ! outward phase speed if (gamma_u < 1.0) then @@ -2256,7 +2255,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, endif elseif (segment%oblique) then dhdt = (u_old(I-1,j,k) - u_new(I-1,j,k)) !old-new - dhdx = (u_new(I-1,j,k) - u_new(I-2,j,k)) !in new time backward sasha for I-1 + dhdx = (u_new(I-1,j,k) - u_new(I-2,j,k)) !in new time backward sashay for I-1 if (dhdt*(segment%grad_normal(J,1,k) + segment%grad_normal(J-1,1,k)) > 0.0) then dhdy = segment%grad_normal(J-1,1,k) elseif (dhdt*(segment%grad_normal(J,1,k) + segment%grad_normal(J-1,1,k)) == 0.0) then @@ -2319,7 +2318,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, else do J=segment%HI%JsdB,segment%HI%JedB dhdt = v_old(i,J,k)-v_new(i,J,k) !old-new - dhdx = v_new(i,J,k)-v_new(i-1,J,k) !in new time backward sasha for I-1 + dhdx = v_new(i,J,k)-v_new(i-1,J,k) !in new time backward sashay for I-1 rx_tang_rad(I,J,k) = 0.0 if (dhdt*dhdx > 0.0) rx_tang_rad(I,J,k) = min( (dhdt/dhdx), rx_max) ! outward phase speed enddo @@ -2398,7 +2397,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, else do J=segment%HI%JsdB,segment%HI%JedB dhdt = v_old(i,J,k)-v_new(i,J,k) !old-new - dhdx = v_new(i,J,k)-v_new(i-1,J,k) !in new time backward sasha for I-1 + dhdx = v_new(i,J,k)-v_new(i-1,J,k) !in new time backward sashay for I-1 if (dhdt*(segment%grad_tan(j,1,k) + segment%grad_tan(j+1,1,k)) > 0.0) then dhdy = segment%grad_tan(j,1,k) elseif (dhdt*(segment%grad_tan(j,1,k) + segment%grad_tan(j+1,1,k)) == 0.0) then @@ -2480,7 +2479,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, do k=1,nz ; do j=segment%HI%jsd,segment%HI%jed if (segment%radiation) then dhdt = (u_old(I+1,j,k) - u_new(I+1,j,k)) !old-new - dhdx = (u_new(I+1,j,k) - u_new(I+2,j,k)) !in new time forward sasha for I+1 + dhdx = (u_new(I+1,j,k) - u_new(I+2,j,k)) !in new time forward sashay for I+1 rx_new = 0.0 if (dhdt*dhdx > 0.0) rx_new = min( (dhdt/dhdx), rx_max) if (gamma_u < 1.0) then @@ -2500,7 +2499,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, endif elseif (segment%oblique) then dhdt = (u_old(I+1,j,k) - u_new(I+1,j,k)) !old-new - dhdx = (u_new(I+1,j,k) - u_new(I+2,j,k)) !in new time forward sasha for I+1 + dhdx = (u_new(I+1,j,k) - u_new(I+2,j,k)) !in new time forward sashay for I+1 if (dhdt*(segment%grad_normal(J,1,k) + segment%grad_normal(J-1,1,k)) > 0.0) then dhdy = segment%grad_normal(J-1,1,k) elseif (dhdt*(segment%grad_normal(J,1,k) + segment%grad_normal(J-1,1,k)) == 0.0) then @@ -2564,7 +2563,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, else do J=segment%HI%JsdB,segment%HI%JedB dhdt = v_old(i+1,J,k)-v_new(i+1,J,k) !old-new - dhdx = v_new(i+1,J,k)-v_new(i+2,J,k) !in new time backward sasha for I-1 + dhdx = v_new(i+1,J,k)-v_new(i+2,J,k) !in new time backward sashay for I-1 rx_tang_rad(I,J,k) = 0.0 if (dhdt*dhdx > 0.0) rx_tang_rad(I,J,k) = min( (dhdt/dhdx), rx_max) ! outward phase speed enddo @@ -2643,7 +2642,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, else do J=segment%HI%JsdB,segment%HI%JedB dhdt = v_old(i+1,J,k)-v_new(i+1,J,k) !old-new - dhdx = v_new(i+1,J,k)-v_new(i+2,J,k) !in new time backward sasha for I-1 + dhdx = v_new(i+1,J,k)-v_new(i+2,J,k) !in new time backward sashay for I-1 if (dhdt*(segment%grad_tan(j,1,k) + segment%grad_tan(j+1,1,k)) > 0.0) then dhdy = segment%grad_tan(j,1,k) elseif (dhdt*(segment%grad_tan(j,1,k) + segment%grad_tan(j+1,1,k)) == 0.0) then @@ -2725,7 +2724,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, do k=1,nz ; do i=segment%HI%isd,segment%HI%ied if (segment%radiation) then dhdt = (v_old(i,J-1,k) - v_new(i,J-1,k)) !old-new - dhdy = (v_new(i,J-1,k) - v_new(i,J-2,k)) !in new time backward sasha for J-1 + dhdy = (v_new(i,J-1,k) - v_new(i,J-2,k)) !in new time backward sashay for J-1 ry_new = 0.0 if (dhdt*dhdy > 0.0) ry_new = min( (dhdt/dhdy), ry_max) if (gamma_u < 1.0) then @@ -2745,7 +2744,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, endif elseif (segment%oblique) then dhdt = (v_old(i,J-1,k) - v_new(i,J-1,k)) !old-new - dhdy = (v_new(i,J-1,k) - v_new(i,J-2,k)) !in new time backward sasha for J-1 + dhdy = (v_new(i,J-1,k) - v_new(i,J-2,k)) !in new time backward sashay for J-1 if (dhdt*(segment%grad_normal(I,1,k) + segment%grad_normal(I-1,1,k)) > 0.0) then dhdx = segment%grad_normal(I-1,1,k) elseif (dhdt*(segment%grad_normal(I,1,k) + segment%grad_normal(I-1,1,k)) == 0.0) then @@ -2808,7 +2807,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, else do I=segment%HI%IsdB,segment%HI%IedB dhdt = u_old(I,j-1,k)-u_new(I,j-1,k) !old-new - dhdy = u_new(I,j-1,k)-u_new(I,j-2,k) !in new time backward sasha for I-1 + dhdy = u_new(I,j-1,k)-u_new(I,j-2,k) !in new time backward sashay for I-1 ry_tang_rad(I,J,k) = 0.0 if (dhdt*dhdy > 0.0) ry_tang_rad(I,J,k) = min( (dhdt/dhdy), rx_max) ! outward phase speed enddo @@ -2887,7 +2886,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, else do I=segment%HI%IsdB,segment%HI%IedB dhdt = u_old(I,j,k)-u_new(I,j,k) !old-new - dhdy = u_new(I,j,k)-u_new(I,j-1,k) !in new time backward sasha for I-1 + dhdy = u_new(I,j,k)-u_new(I,j-1,k) !in new time backward sashay for I-1 if (dhdt*(segment%grad_tan(i,1,k) + segment%grad_tan(i+1,1,k)) > 0.0) then dhdx = segment%grad_tan(i,1,k) elseif (dhdt*(segment%grad_tan(i,1,k) + segment%grad_tan(i+1,1,k)) == 0.0) then @@ -2969,7 +2968,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, do k=1,nz ; do i=segment%HI%isd,segment%HI%ied if (segment%radiation) then dhdt = (v_old(i,J+1,k) - v_new(i,J+1,k)) !old-new - dhdy = (v_new(i,J+1,k) - v_new(i,J+2,k)) !in new time backward sasha for J-1 + dhdy = (v_new(i,J+1,k) - v_new(i,J+2,k)) !in new time backward sashay for J-1 ry_new = 0.0 if (dhdt*dhdy > 0.0) ry_new = min( (dhdt/dhdy), ry_max) if (gamma_u < 1.0) then @@ -2989,7 +2988,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, endif elseif (segment%oblique) then dhdt = (v_old(i,J+1,k) - v_new(i,J+1,k)) !old-new - dhdy = (v_new(i,J+1,k) - v_new(i,J+2,k)) !in new time backward sasha for J-1 + dhdy = (v_new(i,J+1,k) - v_new(i,J+2,k)) !in new time backward sashay for J-1 if (dhdt*(segment%grad_normal(I,1,k) + segment%grad_normal(I-1,1,k)) > 0.0) then dhdx = segment%grad_normal(I-1,1,k) elseif (dhdt*(segment%grad_normal(I,1,k) + segment%grad_normal(I-1,1,k)) == 0.0) then @@ -3053,7 +3052,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, else do I=segment%HI%IsdB,segment%HI%IedB dhdt = u_old(I,j+1,k)-u_new(I,j+1,k) !old-new - dhdy = u_new(I,j+1,k)-u_new(I,j+2,k) !in new time backward sasha for I-1 + dhdy = u_new(I,j+1,k)-u_new(I,j+2,k) !in new time backward sashay for I-1 ry_tang_rad(I,J,k) = 0.0 if (dhdt*dhdy > 0.0) ry_tang_rad(I,J,k) = min( (dhdt/dhdy), rx_max) ! outward phase speed enddo @@ -3132,7 +3131,7 @@ subroutine radiation_open_bdry_conds(OBC, u_new, u_old, v_new, v_old, G, GV, US, else do I=segment%HI%IsdB,segment%HI%IedB dhdt = u_old(I,j+1,k)-u_new(I,j+1,k) !old-new - dhdy = u_new(I,j+1,k)-u_new(I,j+2,k) !in new time backward sasha for I-1 + dhdy = u_new(I,j+1,k)-u_new(I,j+2,k) !in new time backward sashay for I-1 if (dhdt*(segment%grad_tan(i,1,k) + segment%grad_tan(i+1,1,k)) > 0.0) then dhdx = segment%grad_tan(i,1,k) elseif (dhdt*(segment%grad_tan(i,1,k) + segment%grad_tan(i+1,1,k)) == 0.0) then @@ -3420,19 +3419,8 @@ subroutine set_tracer_data(OBC, tv, h, G, GV, PF) real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: h !< Thickness type(param_file_type), intent(in) :: PF !< Parameter file handle - integer :: i, j, k, itt, is, ie, js, je, isd, ied, jsd, jed, nz, n - integer :: isd_off, jsd_off - integer :: IsdB, IedB, JsdB, JedB type(OBC_segment_type), pointer :: segment => NULL() ! pointer to segment type list - character(len=40) :: mdl = "set_tracer_data" ! This subroutine's name. - character(len=200) :: filename, OBC_file, inputdir ! Strings for file/path - - real :: temp_u(G%domain%niglobal+1,G%domain%njglobal) - real :: temp_v(G%domain%niglobal,G%domain%njglobal+1) - - is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke - isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed - IsdB = G%IsdB ; IedB = G%IedB ; JsdB = G%JsdB ; JedB = G%JedB + integer :: i, j, k, n ! For now, there are no radiation conditions applied to the thicknesses, since ! the thicknesses might not be physically motivated. Instead, sponges should be @@ -3723,23 +3711,23 @@ subroutine update_OBC_segment_data(G, GV, US, OBC, tv, h, Time) character(len=200) :: filename, OBC_file, inputdir ! Strings for file/path type(OBC_segment_type), pointer :: segment => NULL() integer, dimension(4) :: siz - real, dimension(:,:,:), pointer :: tmp_buffer_in => NULL() ! Unrotated input + real, dimension(:,:,:), pointer :: tmp_buffer_in => NULL() ! Unrotated input [various units] integer :: ni_seg, nj_seg ! number of src gridpoints along the segments integer :: ni_buf, nj_buf ! Number of filled values in tmp_buffer integer :: i2, j2 ! indices for referencing local domain array integer :: is_obc, ie_obc, js_obc, je_obc ! segment indices within local domain integer :: ishift, jshift ! offsets for staggered locations - real, dimension(:,:,:), allocatable, target :: tmp_buffer + real, dimension(:,:,:), allocatable, target :: tmp_buffer ! A buffer for input data [various units] real, dimension(:), allocatable :: h_stack ! Thicknesses at corner points [H ~> m or kg m-2] integer :: is_obc2, js_obc2 real :: net_H_src ! Total thickness of the incoming flow in the source field [H ~> m or kg m-2] real :: net_H_int ! Total thickness of the incoming flow in the model [H ~> m or kg m-2] - real :: scl_fac ! A nondimensional scaling factor [nondim] - real :: tidal_vel ! Tangential tidal velocity [m s-1] - real :: tidal_elev ! Tidal elevation at an OBC point [m] + real :: scl_fac ! A scaling factor to compensate for differences in total thicknesses [nondim] + real :: tidal_vel ! Interpolated tidal velocity at the OBC points [m s-1] + real :: tidal_elev ! Interpolated tidal elevation at the OBC points [m] real, allocatable :: normal_trans_bt(:,:) ! barotropic transport [H L2 T-1 ~> m3 s-1] - integer :: turns ! Number of index quarter turns - real :: time_delta ! Time since tidal reference date [s] + integer :: turns ! Number of index quarter turns + real :: time_delta ! Time since tidal reference date [s] is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed @@ -3857,15 +3845,15 @@ subroutine update_OBC_segment_data(G, GV, US, OBC, tv, h, Time) ! NOTE: buffer is sized for vertex points, but may be used for faces if (siz(1)==1) then if (OBC%brushcutter_mode) then - allocate(tmp_buffer(1,nj_seg*2-1,segment%field(m)%nk_src)) ! segment data is currrently on supergrid + allocate(tmp_buffer(1,nj_seg*2-1,segment%field(m)%nk_src)) ! segment data is currently on supergrid else - allocate(tmp_buffer(1,nj_seg,segment%field(m)%nk_src)) ! segment data is currrently on supergrid + allocate(tmp_buffer(1,nj_seg,segment%field(m)%nk_src)) ! segment data is currently on supergrid endif else if (OBC%brushcutter_mode) then - allocate(tmp_buffer(ni_seg*2-1,1,segment%field(m)%nk_src)) ! segment data is currrently on supergrid + allocate(tmp_buffer(ni_seg*2-1,1,segment%field(m)%nk_src)) ! segment data is currently on supergrid else - allocate(tmp_buffer(ni_seg,1,segment%field(m)%nk_src)) ! segment data is currrently on supergrid + allocate(tmp_buffer(ni_seg,1,segment%field(m)%nk_src)) ! segment data is currently on supergrid endif endif @@ -4351,11 +4339,12 @@ end subroutine update_OBC_segment_data subroutine update_OBC_ramp(Time, OBC, activate) type(time_type), target, intent(in) :: Time !< Current model time type(ocean_OBC_type), intent(inout) :: OBC !< Open boundary structure - logical, optional, intent(in) :: activate !< Specifiy whether to record the value of + logical, optional, intent(in) :: activate !< Specify whether to record the value of !! Time as the beginning of the ramp period ! Local variables - real :: deltaTime, wghtA + real :: deltaTime ! The time since start of ramping [s] + real :: wghtA ! A temporary variable used to set OBC%ramp_value [nondim] character(len=12) :: msg if (.not. OBC%ramp) return ! This indicates the ramping is turned off @@ -4502,7 +4491,7 @@ subroutine segment_tracer_registry_init(param_file, segment) end subroutine segment_tracer_registry_init -!> Register a tracer array that is active on an OBC segment, potentially also specifing how the +!> Register a tracer array that is active on an OBC segment, potentially also specifying how the !! tracer inflow values are specified. subroutine register_segment_tracer(tr_ptr, param_file, GV, segment, & OBC_scalar, OBC_array) @@ -4974,7 +4963,7 @@ subroutine open_boundary_register_restarts(HI, GV, OBC, Reg, param_file, restart else ! This would be coming from user code such as DOME. if (OBC%ntr /= Reg%ntr) then -! call MOM_error(FATAL, "open_boundary_regiser_restarts: Inconsistent value for ntr") +! call MOM_error(FATAL, "open_boundary_register_restarts: Inconsistent value for ntr") write(mesg,'("Inconsistent values for ntr ", I8," and ",I8,".")') OBC%ntr, Reg%ntr call MOM_error(WARNING, 'open_boundary_register_restarts: '//mesg) endif @@ -5115,7 +5104,7 @@ subroutine adjustSegmentEtaToFitBathymetry(G, GV, US, segment,fld) integer :: n real, allocatable, dimension(:,:,:) :: eta ! Segment source data interface heights [Z ~> m] real :: hTolerance = 0.1 !< Tolerance to exceed adjustment criteria [Z ~> m] - real :: hTmp, eTmp, dilate + ! real :: dilate ! A factor by which to dilate the water column [nondim] character(len=100) :: mesg hTolerance = 0.1*US%m_to_Z diff --git a/src/core/MOM_verticalGrid.F90 b/src/core/MOM_verticalGrid.F90 index 2b597e355f..b856cff3dc 100644 --- a/src/core/MOM_verticalGrid.F90 +++ b/src/core/MOM_verticalGrid.F90 @@ -90,7 +90,7 @@ subroutine verticalGridInit( param_file, GV, US ) ! Local variables integer :: nk, H_power - real :: H_rescale_factor + real :: H_rescale_factor ! The integer power of 2 by which thicknesses are rescaled [nondim] ! This include declares and sets the variable "version". # include "version_variable.h" character(len=16) :: mdl = 'MOM_verticalGrid' diff --git a/src/tracer/MOM_CFC_cap.F90 b/src/tracer/MOM_CFC_cap.F90 index 216b74c735..267d16a174 100644 --- a/src/tracer/MOM_CFC_cap.F90 +++ b/src/tracer/MOM_CFC_cap.F90 @@ -366,7 +366,7 @@ function CFC_cap_stock(h, stocks, G, GV, CS, names, units, stock_index) integer :: CFC_cap_stock !< The number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or nondim] + real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] real :: mass ! The cell volume or mass [H L2 ~> m3 or kg] integer :: i, j, k, is, ie, js, je, nz is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke diff --git a/src/tracer/MOM_OCMIP2_CFC.F90 b/src/tracer/MOM_OCMIP2_CFC.F90 index df86300351..008f6386c9 100644 --- a/src/tracer/MOM_OCMIP2_CFC.F90 +++ b/src/tracer/MOM_OCMIP2_CFC.F90 @@ -506,7 +506,7 @@ function OCMIP2_CFC_stock(h, stocks, G, GV, CS, names, units, stock_index) integer :: OCMIP2_CFC_stock !< The number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or nondim] + real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] real :: mass ! The cell volume or mass [H L2 ~> m3 or kg] integer :: i, j, k, is, ie, js, je, nz is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke diff --git a/src/tracer/MOM_generic_tracer.F90 b/src/tracer/MOM_generic_tracer.F90 index ca5fc65f37..4627d0ec80 100644 --- a/src/tracer/MOM_generic_tracer.F90 +++ b/src/tracer/MOM_generic_tracer.F90 @@ -583,7 +583,7 @@ function MOM_generic_tracer_stock(h, stocks, G, GV, CS, names, units, stock_inde !! number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or nondim] + real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] type(g_tracer_type), pointer :: g_tracer, g_tracer_next real, dimension(:,:,:,:), pointer :: tr_field real, dimension(:,:,:), pointer :: tr_ptr diff --git a/src/tracer/MOM_lateral_boundary_diffusion.F90 b/src/tracer/MOM_lateral_boundary_diffusion.F90 index ef27146a18..19d40f2db1 100644 --- a/src/tracer/MOM_lateral_boundary_diffusion.F90 +++ b/src/tracer/MOM_lateral_boundary_diffusion.F90 @@ -174,8 +174,8 @@ subroutine lateral_boundary_diffusion(G, GV, US, h, Coef_x, Coef_y, dt, Reg, CS) real, dimension(SZI_(G),SZJ_(G)) :: tracer_end !< integrated tracer after LBD is applied. !! [conc H L2 ~> conc m3 or conc kg] integer :: i, j, k, m !< indices to loop over - real :: Idt !< inverse of the time step [s-1] - real :: tmp1, tmp2 !< temporary variables + real :: Idt !< inverse of the time step [T-1 ~> s-1] + real :: tmp1, tmp2 !< temporary variables [conc H L2 ~> conc m3 or conc kg] call cpu_clock_begin(id_clock_lbd) Idt = 1./dt @@ -569,51 +569,50 @@ end subroutine boundary_k_range !> Calculate the lateral boundary diffusive fluxes using the layer by layer method. !! See \ref section_method subroutine fluxes_layer_method(boundary, ke, hbl_L, hbl_R, h_L, h_R, phi_L, phi_R, & - khtr_u, F_layer, area_L, area_R, CS) - - integer, intent(in ) :: boundary !< Which boundary layer SURFACE or BOTTOM [nondim] - integer, intent(in ) :: ke !< Number of layers in the native grid [nondim] - real, intent(in ) :: hbl_L !< Thickness of the boundary boundary - !! layer (left) [H ~> m or kg m-2] - real, intent(in ) :: hbl_R !< Thickness of the boundary boundary - !! layer (right) [H ~> m or kg m-2] - real, dimension(ke), intent(in ) :: h_L !< Thicknesses in the native grid (left) [H ~> m or kg m-2] - real, dimension(ke), intent(in ) :: h_R !< Thicknesses in the native grid (right) [H ~> m or kg m-2] - real, dimension(ke), intent(in ) :: phi_L !< Tracer values in the native grid (left) [conc] - real, dimension(ke), intent(in ) :: phi_R !< Tracer values in the native grid (right) [conc] - real, intent(in ) :: khtr_u !< Horizontal diffusivities times delta t - !! at a velocity point [L2 ~> m2] - real, dimension(ke), intent( out) :: F_layer !< Layerwise diffusive flux at U- or V-point in the native - !! grid [H L2 conc ~> m3 conc] - real, intent(in ) :: area_L !< Area of the horizontal grid (left) [L2 ~> m2] - real, intent(in ) :: area_R !< Area of the horizontal grid (right) [L2 ~> m2] - type(lbd_CS), pointer :: CS !< Lateral diffusion control structure - !! the boundary layer + khtr_u, F_layer, area_L, area_R, CS) + + integer, intent(in ) :: boundary !< Which boundary layer SURFACE or BOTTOM [nondim] + integer, intent(in ) :: ke !< Number of layers in the native grid [nondim] + real, intent(in ) :: hbl_L !< Thickness of the boundary boundary + !! layer (left) [H ~> m or kg m-2] + real, intent(in ) :: hbl_R !< Thickness of the boundary boundary + !! layer (right) [H ~> m or kg m-2] + real, dimension(ke), intent(in ) :: h_L !< Thicknesses in the native grid (left) [H ~> m or kg m-2] + real, dimension(ke), intent(in ) :: h_R !< Thicknesses in the native grid (right) [H ~> m or kg m-2] + real, dimension(ke), intent(in ) :: phi_L !< Tracer values in the native grid (left) [conc] + real, dimension(ke), intent(in ) :: phi_R !< Tracer values in the native grid (right) [conc] + real, intent(in ) :: khtr_u !< Horizontal diffusivities times the time step + !! at a velocity point [L2 ~> m2] + real, dimension(ke), intent( out) :: F_layer !< Layerwise diffusive flux at U- or V-point + !! in the native grid [H L2 conc ~> m3 conc] + real, intent(in ) :: area_L !< Area of the horizontal grid (left) [L2 ~> m2] + real, intent(in ) :: area_R !< Area of the horizontal grid (right) [L2 ~> m2] + type(lbd_CS), pointer :: CS !< Lateral diffusion control structure + ! Local variables - real, dimension(:), allocatable :: dz_top !< The LBD z grid to be created [L ~ m] - real, dimension(:), allocatable :: phi_L_z !< Tracer values in the ztop grid (left) [conc] - real, dimension(:), allocatable :: phi_R_z !< Tracer values in the ztop grid (right) [conc] - real, dimension(:), allocatable :: F_layer_z !< Diffusive flux at U/V-point in the ztop grid [H L2 conc ~> m3 conc] - real, dimension(ke) :: h_vel !< Thicknesses at u- and v-points in the native grid - !! The harmonic mean is used to avoid zero values [H ~> m or kg m-2] - real :: khtr_avg !< Thickness-weighted diffusivity at the u-point [m^2 s^-1] + real, allocatable :: dz_top(:) !< The LBD z grid to be created [H ~> m or kg m-2] + real, allocatable :: phi_L_z(:) !< Tracer values in the ztop grid (left) [conc] + real, allocatable :: phi_R_z(:) !< Tracer values in the ztop grid (right) [conc] + real, allocatable :: F_layer_z(:) !< Diffusive flux at U/V-point in the ztop grid [H L2 conc ~> m3 conc] + real :: h_vel(ke) !< Thicknesses at u- and v-points in the native grid + !! The harmonic mean is used to avoid zero values [H ~> m or kg m-2] + real :: khtr_avg !< Thickness-weighted diffusivity at the velocity-point [L2 T-1 ~> m s-1] !! This is just to remind developers that khtr_avg should be !! computed once khtr is 3D. - real :: htot !< Total column thickness [H ~> m or kg m-2] + real :: htot !< Total column thickness [H ~> m or kg m-2] integer :: k, k_bot_min, k_top_max !< k-indices, min and max for bottom and top, respectively integer :: k_bot_max, k_top_min !< k-indices, max and min for bottom and top, respectively integer :: k_bot_diff, k_top_diff !< different between left and right k-indices for bottom and top, respectively integer :: k_top_L, k_bot_L !< k-indices left native grid integer :: k_top_R, k_bot_R !< k-indices right native grid real :: zeta_top_L, zeta_top_R !< distance from the top of a layer to the boundary - !! layer depth in the native grid [nondim] + !! layer depth in the native grid [nondim] real :: zeta_bot_L, zeta_bot_R !< distance from the bottom of a layer to the boundary - !!layer depth in the native grid [nondim] - real :: hbl_min !< minimum BLD (left and right) [m] - real :: wgt !< weight to be used in the linear transition to the interior [nondim] + !! layer depth in the native grid [nondim] + real :: wgt !< weight to be used in the linear transition to the interior [nondim] real :: a !< coefficient to be used in the linear transition to the interior [nondim] - real :: tmp1, tmp2 !< dummy variables - real :: htot_max !< depth below which no fluxes should be applied + real :: tmp1, tmp2 !< dummy variables [H ~> m or kg m-2] + real :: htot_max !< depth below which no fluxes should be applied [H ~> m or kg m-2] integer :: nk !< number of layers in the LBD grid F_layer(:) = 0.0 @@ -743,16 +742,16 @@ logical function near_boundary_unit_tests( verbose ) integer, parameter :: nk = 2 ! Number of layers real, dimension(nk+1) :: eta1 ! Updated interfaces with one extra value [m] real, dimension(:), allocatable :: h1 ! Upates layer thicknesses [m] - real, dimension(nk) :: phi_L, phi_R ! Tracer values (left and right column) [ nondim m^-3 ] + real, dimension(nk) :: phi_L, phi_R ! Tracer values (left and right column) [conc] real, dimension(nk) :: h_L, h_R ! Layer thickness (left and right) [m] - real :: khtr_u ! Horizontal diffusivities at U-point [m^2 s^-1] + real :: khtr_u ! Horizontal diffusivities at U-point [m2 s-1] real :: hbl_L, hbl_R ! Depth of the boundary layer (left and right) [m] - real, dimension(nk) :: F_layer ! Diffusive flux within each layer at U-point [nondim s^-1] + real, dimension(nk) :: F_layer ! Diffusive flux within each layer at U-point [conc m3 s-1] character(len=120) :: test_name ! Title of the unit test integer :: k_top ! Index of cell containing top of boundary - real :: zeta_top ! Nondimension position + real :: zeta_top ! Nondimension position [nondim] integer :: k_bot ! Index of cell containing bottom of boundary - real :: zeta_bot ! Nondimension position + real :: zeta_bot ! Nondimension position [nondim] type(lbd_CS), pointer :: CS allocate(CS) diff --git a/src/tracer/MOM_neutral_diffusion.F90 b/src/tracer/MOM_neutral_diffusion.F90 index a06af6cd57..fd479eeaf3 100644 --- a/src/tracer/MOM_neutral_diffusion.F90 +++ b/src/tracer/MOM_neutral_diffusion.F90 @@ -969,9 +969,9 @@ subroutine find_neutral_surface_positions_continuous(nk, Pl, Tl, Sl, dRdTl, dRdS real, dimension(nk+1), intent(in) :: dRdTr !< Left-column dRho/dT [R degC-1 ~> kg m-3 degC-1] real, dimension(nk+1), intent(in) :: dRdSr !< Left-column dRho/dS [R ppt-1 ~> kg m-3 ppt-1] real, dimension(2*nk+2), intent(inout) :: PoL !< Fractional position of neutral surface within - !! layer KoL of left column + !! layer KoL of left column [nondim] real, dimension(2*nk+2), intent(inout) :: PoR !< Fractional position of neutral surface within - !! layer KoR of right column + !! layer KoR of right column [nondim] integer, dimension(2*nk+2), intent(inout) :: KoL !< Index of first left interface above neutral surface integer, dimension(2*nk+2), intent(inout) :: KoR !< Index of first right interface above neutral surface real, dimension(2*nk+1), intent(inout) :: hEff !< Effective thickness between two neutral surfaces @@ -986,7 +986,6 @@ subroutine find_neutral_surface_positions_continuous(nk, Pl, Tl, Sl, dRdTl, dRdS integer :: k_surface ! Index of neutral surface integer :: kl ! Index of left interface integer :: kr ! Index of right interface - real :: dRdT, dRdS ! dRho/dT [kg m-3 degC-1] and dRho/dS [kg m-3 ppt-1] for the neutral surface logical :: searching_left_column ! True if searching for the position of a right interface in the left column logical :: searching_right_column ! True if searching for the position of a left interface in the right column logical :: reached_bottom ! True if one of the bottom-most interfaces has been used as the target @@ -1246,7 +1245,7 @@ subroutine find_neutral_surface_positions_discontinuous(CS, nk, & integer, optional, intent(in) :: k_bot_L !< k-index for the boundary layer (left) [nondim] integer, optional, intent(in) :: k_bot_R !< k-index for the boundary layer (right) [nondim] logical, optional, intent(in) :: hard_fail_heff !< If true (default) bring down the model if the - !! neutral surfaces ever cross [logical] + !! neutral surfaces ever cross ! Local variables integer :: ns ! Number of neutral surfaces integer :: k_surface ! Index of neutral surface diff --git a/src/tracer/MOM_offline_aux.F90 b/src/tracer/MOM_offline_aux.F90 index 9486e87369..af8b422238 100644 --- a/src/tracer/MOM_offline_aux.F90 +++ b/src/tracer/MOM_offline_aux.F90 @@ -704,9 +704,9 @@ subroutine update_offline_from_files(G, GV, nk_input, mean_file, sum_file, snap_ fluxes%netMassOut(:,:) = 0.0 fluxes%netMassIn(:,:) = 0.0 call MOM_read_data(surf_file,'massout_flux_sum',fluxes%netMassOut, G%Domain, & - timelevel=ridx_sum) + timelevel=ridx_sum, scale=GV%kg_m2_to_H) call MOM_read_data(surf_file,'massin_flux_sum', fluxes%netMassIn, G%Domain, & - timelevel=ridx_sum) + timelevel=ridx_sum, scale=GV%kg_m2_to_H) do j=js,je ; do i=is,ie if (G%mask2dT(i,j)<1.0) then diff --git a/src/tracer/MOM_tracer_diabatic.F90 b/src/tracer/MOM_tracer_diabatic.F90 index c1e39598cc..a701c98553 100644 --- a/src/tracer/MOM_tracer_diabatic.F90 +++ b/src/tracer/MOM_tracer_diabatic.F90 @@ -425,7 +425,7 @@ subroutine applyTracerBoundaryFluxesInOut(G, GV, Tr, dt, fluxes, h, evap_CFL_lim type(ocean_grid_type), intent(in ) :: G !< Grid structure type(verticalGrid_type), intent(in ) :: GV !< ocean vertical grid structure - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: Tr !< Tracer concentration on T-cell + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: Tr !< Tracer concentration on T-cell [conc] real, intent(in ) :: dt !< Time-step over which forcing is applied [T ~> s] type(forcing), intent(in ) :: fluxes !< Surface fluxes container real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(inout) :: h !< Layer thickness [H ~> m or kg m-2] @@ -436,30 +436,36 @@ subroutine applyTracerBoundaryFluxesInOut(G, GV, Tr, dt, fluxes, h, evap_CFL_lim !! which fluxes can be applied [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G)), optional, intent(in ) :: in_flux_optional !< The total time-integrated !! amount of tracer that enters with freshwater + !! [conc H ~> conc m or conc kg m-2] real, dimension(SZI_(G),SZJ_(G)), optional, intent(in) :: out_flux_optional !< The total time-integrated !! amount of tracer that leaves with freshwater + !! [conc H ~> conc m or conc kg m-2] logical, optional, intent(in) :: update_h_opt !< Optional flag to determine whether !! h should be updated integer, parameter :: maxGroundings = 5 integer :: numberOfGroundings, iGround(maxGroundings), jGround(maxGroundings) - real :: H_limit_fluxes, IforcingDepthScale - real :: dThickness, dTracer - real :: fractionOfForcing, hOld, Ithickness - real :: RivermixConst ! A constant used in implementing river mixing [Pa s]. - real, dimension(SZI_(G)) :: & - netMassInOut, & ! surface water fluxes [H ~> m or kg m-2] over time step - netMassIn, & ! mass entering ocean surface [H ~> m or kg m-2] over a time step - netMassOut ! mass leaving ocean surface [H ~> m or kg m-2] over a time step - - real, dimension(SZI_(G),SZK_(GV)) :: h2d, Tr2d - real, dimension(SZI_(G),SZJ_(G)) :: in_flux ! The total time-integrated amount of tracer! - ! that enters with freshwater - real, dimension(SZI_(G),SZJ_(G)) :: out_flux ! The total time-integrated amount of tracer! - ! that leaves with freshwater - real, dimension(SZI_(G)) :: in_flux_1d, out_flux_1d - real :: hGrounding(maxGroundings) - real :: Tr_in + real :: IforcingDepthScale ! The inverse of the scale over which to apply forcing [H-1 ~> m-1 or m2 kg-1] + real :: dThickness ! The change in a layer's thickness [H ~> m or kg m-2] + real :: dTracer ! The change in the integrated tracer content of a layer [conc H ~> conc m or conc kg m-2] + real :: fractionOfForcing ! The fraction of the forcing to apply to a layer [nondim] + real :: hOld ! The layer thickness before surface forcing is applied [H ~> m or kg m-2] + real :: Ithickness ! The inverse of the new layer thickness [H-1 ~> m-1 or m2 kg-1] + + real :: h2d(SZI_(G),SZK_(GV)) ! A 2-d work copy of layer thicknesses [H ~> m or kg m-2] + real :: Tr2d(SZI_(G),SZK_(GV)) ! A 2-d work copy of tracer concentrations [conc] + real :: in_flux(SZI_(G),SZJ_(G)) ! The total time-integrated amount of tracer that + ! enters with freshwater [conc H ~> conc m or conc kg m-2] + real :: out_flux(SZI_(G),SZJ_(G)) ! The total time-integrated amount of tracer that + ! leaves with freshwater [conc H ~> conc m or conc kg m-2] + real :: netMassIn(SZI_(G)) ! The remaining mass entering ocean surface [H ~> m or kg m-2] + real :: netMassOut(SZI_(G)) ! The remaining mass leaving ocean surface [H ~> m or kg m-2] + real :: in_flux_1d(SZI_(G)) ! The remaining amount of tracer that enters with + ! the freshwater [conc H ~> conc m or conc kg m-2] + real :: out_flux_1d(SZI_(G)) ! The remaining amount of tracer that leaves with + ! the freshwater [conc H ~> conc m or conc kg m-2] + real :: hGrounding(maxGroundings) ! The remaining fresh water flux that was not able to be + ! supplied from a column that grounded out [H ~> m or kg m-2] logical :: update_h integer :: i, j, is, ie, js, je, k, nz, n, nsw character(len=45) :: mesg @@ -493,10 +499,9 @@ subroutine applyTracerBoundaryFluxesInOut(G, GV, Tr, dt, fluxes, h, evap_CFL_lim !$OMP IforcingDepthScale,minimum_forcing_depth, & !$OMP numberOfGroundings,iGround,jGround,update_h, & !$OMP in_flux,out_flux,hGrounding,evap_CFL_limit) & -!$OMP private(h2d,Tr2d,netMassInOut,netMassOut, & +!$OMP private(h2d,Tr2d,netMassIn,netMassOut, & !$OMP in_flux_1d,out_flux_1d,fractionOfForcing, & -!$OMP dThickness,dTracer,hOld,Ithickness, & -!$OMP netMassIn, Tr_in) +!$OMP dThickness,dTracer,hOld,Ithickness) ! Work in vertical slices for efficiency do j=js,je @@ -521,8 +526,8 @@ subroutine applyTracerBoundaryFluxesInOut(G, GV, Tr, dt, fluxes, h, evap_CFL_lim ! Note here that the aggregateFW flag has already been taken care of in the call to ! applyBoundaryFluxesInOut do i=is,ie - netMassOut(i) = fluxes%netMassOut(i,j) - netMassIn(i) = fluxes%netMassIn(i,j) + netMassOut(i) = fluxes%netMassOut(i,j) + netMassIn(i) = fluxes%netMassIn(i,j) enddo ! Apply the surface boundary fluxes in three steps: diff --git a/src/tracer/MOM_tracer_registry.F90 b/src/tracer/MOM_tracer_registry.F90 index dedfe4e30a..c729231927 100644 --- a/src/tracer/MOM_tracer_registry.F90 +++ b/src/tracer/MOM_tracer_registry.F90 @@ -826,7 +826,7 @@ subroutine MOM_tracer_chkinv(mesg, G, GV, h, Tr, ntr) integer, intent(in) :: ntr !< number of registered tracers ! Local variables - real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> nondim or m3 kg-1] + real :: vol_scale ! The dimensional scaling factor to convert volumes to m3 [m3 H-1 L-2 ~> 1 or m3 kg-1] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: tr_inv ! Volumetric tracer inventory in each cell [conc m3] real :: total_inv ! The total amount of tracer [conc m3] integer :: is, ie, js, je, nz diff --git a/src/tracer/advection_test_tracer.F90 b/src/tracer/advection_test_tracer.F90 index d6d1ac25fe..b713803182 100644 --- a/src/tracer/advection_test_tracer.F90 +++ b/src/tracer/advection_test_tracer.F90 @@ -358,7 +358,7 @@ function advection_test_stock(h, stocks, G, GV, CS, names, units, stock_index) integer :: advection_test_stock !< the number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or nondim] + real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] integer :: i, j, k, is, ie, js, je, nz, m is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke diff --git a/src/tracer/boundary_impulse_tracer.F90 b/src/tracer/boundary_impulse_tracer.F90 index bc5d19b4fb..18e9b8dc8e 100644 --- a/src/tracer/boundary_impulse_tracer.F90 +++ b/src/tracer/boundary_impulse_tracer.F90 @@ -302,7 +302,7 @@ function boundary_impulse_stock(h, stocks, G, GV, CS, names, units, stock_index) ! is present, only the stock corresponding to that coded index is returned. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or nondim] + real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] integer :: i, j, k, is, ie, js, je, nz, m is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke diff --git a/src/tracer/dye_example.F90 b/src/tracer/dye_example.F90 index 9a3ca019bd..91806bb94e 100644 --- a/src/tracer/dye_example.F90 +++ b/src/tracer/dye_example.F90 @@ -341,7 +341,7 @@ function dye_stock(h, stocks, G, GV, CS, names, units, stock_index) !! calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or nondim] + real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] integer :: i, j, k, is, ie, js, je, nz, m is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke diff --git a/src/tracer/ideal_age_example.F90 b/src/tracer/ideal_age_example.F90 index 60d9c02aa0..ca47a8ca1d 100644 --- a/src/tracer/ideal_age_example.F90 +++ b/src/tracer/ideal_age_example.F90 @@ -385,7 +385,7 @@ function ideal_age_stock(h, stocks, G, GV, CS, names, units, stock_index) integer :: ideal_age_stock !< The number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or nondim] + real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] integer :: i, j, k, is, ie, js, je, nz, m is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke diff --git a/src/tracer/oil_tracer.F90 b/src/tracer/oil_tracer.F90 index 0ebf9dcfc9..862209a688 100644 --- a/src/tracer/oil_tracer.F90 +++ b/src/tracer/oil_tracer.F90 @@ -327,7 +327,7 @@ subroutine oil_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US ! Local variables real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_work ! Used so that h can be modified [H ~> m or kg m-2] real :: Isecs_per_year = 1.0 / (365.0*86400.0) - real :: vol_scale ! A conversion factor for volumes into m3 [m3 H-1 L-2 ~> nondim or m3 kg-1] + real :: vol_scale ! A conversion factor for volumes into m3 [m3 H-1 L-2 ~> 1 or m3 kg-1] real :: year, h_total, ldecay integer :: i, j, k, is, ie, js, je, nz, m, k_max is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke @@ -414,7 +414,7 @@ function oil_stock(h, stocks, G, GV, CS, names, units, stock_index) integer :: oil_stock !< The number of stocks calculated here. ! Local variables - real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or nondim] + real :: stock_scale ! The dimensional scaling factor to convert stocks to kg [kg H-1 L-2 ~> kg m-3 or 1] integer :: i, j, k, is, ie, js, je, nz, m is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke