Skip to content

Commit

Permalink
Merge gfdl updates (#2)
Browse files Browse the repository at this point in the history 
* Add weighted d[uv]_dt_str diagnostics (mom-ocean#1539)

This adds four new diagnostics building on the wind stress acceleration diagnostics du_dt_str, dv_dt_str (from mom-ocean#1437)
- their thickness-weighted versions: h_du_dt_str, h_dv_dt_str (completing the set of diags from 3D thickness x momentum diagnostics mom-ocean#1398)
- their viscous remnant fraction: du_dt_str_visc_rem, dv_dt_str_visc_rem (completing the set of diags from Visc_rem_[uv] multiplied momentum budget diagnostics ocean-eddy-cpt#10)

Nora did some quick tests with the CPT NeverWorld2 setup, which confirm that online and offline multiplication by 1) h or 2) visc_rem_[uv] coincide (up to 1) interpolation error and 2) numerical noise, respectively), and illustrated this beautifully with some figures that accompanied the first of the three commits that were squashed into one.

* +Argument cleanup in vertical diffusivity code

  Cleaned up 26 falsely optional or unused arguments in the vertical
diffusivity code, and related changes.  Several descriptive comments were
also corrected, including the correction of the units of 10 variables related
to CVMix_KPP.  This commit includes:
 - Made the Kd_int arguments to set_diffusivity() and 3 subsidiary routines
   mandatory and reordered the arguments so that the non-optional arguments
   come before the grid types
 - Made the halo_TS and double_diffuse arguments to set_diffusivity_init()
   mandatory.
 - Made the Time argument to ALE_sponge() mandatory.
 - Made the Kd and Kv arguments to calculate_CVMIX_conv() mandatory.
 - Removed the unused halo argument to adjust_salt().
 - Removed the unused Kddt_convect argument to full_convection().
 - Made the halo arguments to full_convection()and smoothed_dRdT_dRdS()
   mandatory.
 - Made the useALEalgorithm argument to geothermal_init() mandatory.
 - Removed the unused initialize_all arguments to Calculate_kappa_shear() and
   Calc_kappa_shear_vertex().
 - Removed the unused I_Ld2_1d and dz_Int_1d arguments to kappa_shear_column().
 - Made 3 arguments to calculate_projected_state() mandatory and reordered the
   arguments accordingly.
 - Eliminating the unused skip_diags arguments to calculateBuoyancyFlux() and
   extractFluxes(), which are now effectively always false.
All answers are bitwise identical, and no output is changed.

* +Move rotate_dyn_horgrid to MOM_dyn_horgrid module

  Moved the routine rotate_dyngrid() from the MOM_transcribe_grid module to
rotate_dyn_horgrid() in the MOM_dyn_horgrid module so that this routine can also
be used at some point by SIS2 to implement rotational consistency testing, and
also to reflect that this routine only works with types from its new module.
The two routines are the same apart from some added comments, and the old name
of rotate_dyngrid() is still available from MOM_transcribe_grid via a module use
statement.  All answers are bitwise identical.

* +Reduce use of dyn_horgrid in initialize_MOM

  Minimized the dependence on dyn_horgrid in initialize_MOM by working directly
with the horizontal index type whereever possible and by moving the calls that
create the MOM_grid_type earlier in the routine, to limit the duration of the
dyn_horgrid_type, and to better co-locate grid-related parameters in the
parameter_doc files.  Also uses the new interface to rotate_dyn_horgrid from the
MOM_dyn_horgrid module in place of the rotate_dyngrid interface from the
MOM_transcribe_grid module.  All answers are bitwise identical, but the order of
some entries in the MOM_parameter_doc files has changed.

* (*)Fix compile-time issues with MOM_sum_driver.F90

  Modified drivers/unit_drivers/MOM_sum_driver.F90 to compile with the latest
version of the rest of the MOM6 code by using the proper types in the various
initialization calls, and verified that it runs as intended.

Co-authored-by: Nora Loose <NoraLoose@users.noreply.github.com>
Co-authored-by: Robert Hallberg <Robert.Hallberg@noaa.gov>
  • Loading branch information
@sditkovsky @NoraLoose @Hallberg-NOAA
3 people authored Nov 17, 2021
1 parent 557db30 commit 91c9caa
Show file tree
Hide file tree
Showing 15 changed files with 493 additions and 553 deletions.
39 changes: 24 additions & 15 deletions config_src/drivers/unit_drivers/MOM_sum_driver.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -18,13 +18,14 @@ program MOM_main
use MOM_coms, only : EFP_type, operator(+), operator(-), assignment(=), EFP_to_real, real_to_EFP
use MOM_cpu_clock, only : cpu_clock_id, cpu_clock_begin, cpu_clock_end
use MOM_cpu_clock, only : CLOCK_COMPONENT
use MOM_domains, only : MOM_domains_init, MOM_infra_init, MOM_infra_end
use MOM_domains, only : MOM_domain_type, MOM_domains_init, MOM_infra_init, MOM_infra_end
use MOM_dyn_horgrid, only : dyn_horgrid_type, create_dyn_horgrid, destroy_dyn_horgrid
use MOM_error_handler, only : MOM_error, MOM_mesg, WARNING, FATAL, is_root_pe
use MOM_error_handler, only : MOM_set_verbosity
use MOM_file_parser, only : read_param, get_param, log_param, log_version, param_file_type
use MOM_file_parser, only : open_param_file, close_param_file
use MOM_grid, only : MOM_grid_init, ocean_grid_type
use MOM_grid_initialize, only : set_grid_metrics
use MOM_hor_index, only : hor_index_type, hor_index_init
use MOM_io, only : MOM_io_init, file_exists, open_file, close_file
use MOM_io, only : check_nml_error, io_infra_init, io_infra_end
use MOM_io, only : APPEND_FILE, ASCII_FILE, READONLY_FILE, SINGLE_FILE
Expand All @@ -33,9 +34,10 @@ program MOM_main

#include <MOM_memory.h>

type(ocean_grid_type) :: grid ! A structure containing metrics and grid info.

type(param_file_type) :: param_file ! The structure indicating the file(s)
type(MOM_domain_type), pointer :: Domain => NULL() !< Ocean model domain
type(dyn_horgrid_type), pointer :: grid => NULL() ! A structure containing metrics and grid info
type(hor_index_type) :: HI ! A hor_index_type for array extents
type(param_file_type) :: param_file ! The structure indicating the file(s)
! containing all run-time parameters.
real :: max_depth ! The maximum ocean depth [m]
integer :: verbosity
Expand Down Expand Up @@ -76,14 +78,16 @@ program MOM_main
verbosity = 2 ; call read_param(param_file, "VERBOSITY", verbosity)
call MOM_set_verbosity(verbosity)

call MOM_domains_init(grid%domain, param_file)
call MOM_domains_init(Domain, param_file)

call MOM_io_init(param_file)
! call diag_mediator_init(param_file)
call MOM_grid_init(grid, param_file)
call hor_index_init(Domain, HI, param_file)
call create_dyn_horgrid(grid, HI)
grid%Domain => Domain

is = grid%isc ; ie = grid%iec ; js = grid%jsc ; je = grid%jec
isd = grid%isd ; ied = grid%ied ; jsd = grid%jsd ; jed = grid%jed
is = HI%isc ; ie = HI%iec ; js = HI%jsc ; je = HI%jec
isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed

! Read all relevant parameters and write them to the model log.
call log_version(param_file, "MOM", version, "")
Expand All @@ -99,7 +103,7 @@ program MOM_main
allocate(depth_tot_std(num_sums)) ; depth_tot_std(:) = 0.0
allocate(depth_tot_fastR(num_sums)) ; depth_tot_fastR(:) = 0.0

! Set up the parameters of the physical domain (i.e. the grid), G
! Set up the parameters of the physical grid
call set_grid_metrics(grid, param_file)

! Set up the bottom depth, grid%bathyT either analytically or from file
Expand Down Expand Up @@ -157,21 +161,24 @@ program MOM_main
endif
enddo

call destroy_dyn_horgrid(grid)
call io_infra_end ; call MOM_infra_end

contains

!> This subroutine sets up the benchmark test case topography for debugging
subroutine benchmark_init_topog_local(D, G, param_file, max_depth)
type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
real, dimension(SZI_(G),SZJ_(G)), intent(out) :: D !< The ocean bottom depth in m
type(dyn_horgrid_type), intent(in) :: G !< The dynamic horizontal grid type
real, dimension(G%isd:G%ied,G%jsd:G%jed), &
intent(out) :: D !< Ocean bottom depth in m or [Z ~> m] if US is present
type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
real, intent(in) :: max_depth !< The maximum ocean depth [m]

real :: min_depth ! The minimum ocean depth in m.
real :: PI ! 3.1415926... calculated as 4*atan(1)
real :: D0 ! A constant to make the maximum !
! basin depth MAXIMUM_DEPTH. !
real :: m_to_Z ! A dimensional rescaling factor.
real :: x, y
! This include declares and sets the variable "version".
# include "version_variable.h"
Expand All @@ -180,12 +187,14 @@ subroutine benchmark_init_topog_local(D, G, param_file, max_depth)
is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec
isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed

call log_version(param_file, mdl, version)
m_to_Z = 1.0 ! ; if (present(US)) m_to_Z = US%m_to_Z

call log_version(param_file, mdl, version, "")
call get_param(param_file, mdl, "MINIMUM_DEPTH", min_depth, &
"The minimum depth of the ocean.", units="m", default=0.0)
"The minimum depth of the ocean.", units="m", default=0.0, scale=m_to_Z)

PI = 4.0*atan(1.0)
D0 = max_depth / 0.5;
D0 = max_depth / 0.5

! Calculate the depth of the bottom.
do i=is,ie ; do j=js,je
Expand Down
94 changes: 44 additions & 50 deletions src/core/MOM.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -75,6 +75,7 @@ module MOM
use MOM_dynamics_unsplit_RK2, only : initialize_dyn_unsplit_RK2, end_dyn_unsplit_RK2
use MOM_dynamics_unsplit_RK2, only : MOM_dyn_unsplit_RK2_CS
use MOM_dyn_horgrid, only : dyn_horgrid_type, create_dyn_horgrid, destroy_dyn_horgrid
use MOM_dyn_horgrid, only : rotate_dyn_horgrid
use MOM_EOS, only : EOS_init, calculate_density, calculate_TFreeze, EOS_domain
use MOM_fixed_initialization, only : MOM_initialize_fixed
use MOM_forcing_type, only : allocate_forcing_type, allocate_mech_forcing
Expand Down Expand Up @@ -123,7 +124,6 @@ module MOM
use MOM_tracer_flow_control, only : tracer_flow_control_init, call_tracer_surface_state
use MOM_tracer_flow_control, only : tracer_flow_control_end
use MOM_transcribe_grid, only : copy_dyngrid_to_MOM_grid, copy_MOM_grid_to_dyngrid
use MOM_transcribe_grid, only : rotate_dyngrid
use MOM_unit_scaling, only : unit_scale_type, unit_scaling_init
use MOM_unit_scaling, only : unit_scaling_end, fix_restart_unit_scaling
use MOM_variables, only : surface, allocate_surface_state, deallocate_surface_state
Expand Down Expand Up @@ -1714,7 +1714,7 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
type(hor_index_type), pointer :: HI => NULL() ! A hor_index_type for array extents
type(hor_index_type), target :: HI_in ! HI on the input grid
type(verticalGrid_type), pointer :: GV => NULL()
type(dyn_horgrid_type), pointer :: dG => NULL()
type(dyn_horgrid_type), pointer :: dG => NULL(), test_dG => NULL()
type(dyn_horgrid_type), pointer :: dG_in => NULL()
type(diag_ctrl), pointer :: diag => NULL()
type(unit_scale_type), pointer :: US => NULL()
Expand Down Expand Up @@ -2174,8 +2174,7 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
! Swap axes for quarter and 3-quarter turns
if (CS%rotate_index) then
allocate(CS%G)
call clone_MOM_domain(G_in%Domain, CS%G%Domain, turns=turns, &
domain_name="MOM_rot")
call clone_MOM_domain(G_in%Domain, CS%G%Domain, turns=turns, domain_name="MOM_rot")
first_direction = modulo(first_direction + turns, 2)
else
CS%G => G_in
Expand All @@ -2200,38 +2199,50 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
local_indexing=.not.global_indexing)
call create_dyn_horgrid(dG_in, HI_in, bathymetry_at_vel=bathy_at_vel)
call clone_MOM_domain(G_in%Domain, dG_in%Domain)
! Also allocate the input ocean_grid_type type at this point based on the same information.
call MOM_grid_init(G_in, param_file, US, HI_in, bathymetry_at_vel=bathy_at_vel)

! Allocate initialize time-invariant MOM variables.
call MOM_initialize_fixed(dG_in, US, OBC_in, param_file, .false., dirs%output_directory)

! Copy the grid metrics and bathymetry to the ocean_grid_type
call copy_dyngrid_to_MOM_grid(dG_in, G_in, US)

call callTree_waypoint("returned from MOM_initialize_fixed() (initialize_MOM)")

! Determine HI and dG for the model index map.
call verticalGridInit( param_file, CS%GV, US )
GV => CS%GV

! Shift from using the temporary dynamic grid type to using the final (potentially static)
! and properly rotated ocean-specific grid type and horizontal index type.
if (CS%rotate_index) then
allocate(HI)
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.
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)
call rotate_dyngrid(dG_in, dG, US, turns)
call rotate_dyn_horgrid(dG_in, dG, US, turns)
call copy_dyngrid_to_MOM_grid(dG, G, US)

if (associated(OBC_in)) then
! TODO: General OBC index rotations is not yet supported.
if (modulo(turns, 4) /= 1) &
call MOM_error(FATAL, "OBC index rotation of 180 and 270 degrees is not yet supported.")
allocate(CS%OBC)
call rotate_OBC_config(OBC_in, dG_in, CS%OBC, dG, turns)
endif

call destroy_dyn_horgrid(dG)
else
! If not rotated, then G_in and G are the same grid.
HI => HI_in
dG => dG_in
G => G_in
CS%OBC => OBC_in
endif

call verticalGridInit( param_file, CS%GV, US )
GV => CS%GV

! Allocate the auxiliary non-symmetric domain for debugging or I/O purposes.
if (CS%debug .or. dG%symmetric) &
call clone_MOM_domain(dG%Domain, dG%Domain_aux, symmetric=.false.)
! dG_in is retained for now so that it can be used with write_ocean_geometry_file() below.

call callTree_waypoint("grids initialized (initialize_MOM)")

Expand All @@ -2240,9 +2251,9 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
call tracer_registry_init(param_file, CS%tracer_Reg)

! Allocate and initialize space for the primary time-varying MOM variables.
is = dG%isc ; ie = dG%iec ; js = dG%jsc ; je = dG%jec ; nz = GV%ke
isd = dG%isd ; ied = dG%ied ; jsd = dG%jsd ; jed = dG%jed
IsdB = dG%IsdB ; IedB = dG%IedB ; JsdB = dG%JsdB ; JedB = dG%JedB
is = HI%isc ; ie = HI%iec ; js = HI%jsc ; je = HI%jec ; nz = GV%ke
isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed
IsdB = HI%IsdB ; IedB = HI%IedB ; JsdB = HI%JsdB ; JedB = HI%JedB
ALLOC_(CS%u(IsdB:IedB,jsd:jed,nz)) ; CS%u(:,:,:) = 0.0
ALLOC_(CS%v(isd:ied,JsdB:JedB,nz)) ; CS%v(:,:,:) = 0.0
ALLOC_(CS%h(isd:ied,jsd:jed,nz)) ; CS%h(:,:,:) = GV%Angstrom_H
Expand Down Expand Up @@ -2279,12 +2290,12 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
else
conv2salt = GV%H_to_kg_m2
endif
call register_tracer(CS%tv%T, CS%tracer_Reg, param_file, dG%HI, GV, &
call register_tracer(CS%tv%T, CS%tracer_Reg, param_file, HI, GV, &
tr_desc=vd_T, registry_diags=.true., flux_nameroot='T', &
flux_units='W', flux_longname='Heat', &
flux_scale=conv2watt, convergence_units='W m-2', &
convergence_scale=conv2watt, CMOR_tendprefix="opottemp", diag_form=2)
call register_tracer(CS%tv%S, CS%tracer_Reg, param_file, dG%HI, GV, &
call register_tracer(CS%tv%S, CS%tracer_Reg, param_file, HI, GV, &
tr_desc=vd_S, registry_diags=.true., flux_nameroot='S', &
flux_units=S_flux_units, flux_longname='Salt', &
flux_scale=conv2salt, convergence_units='kg m-2 s-1', &
Expand Down Expand Up @@ -2364,24 +2375,24 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
call restart_init(param_file, restart_CSp)
call set_restart_fields(GV, US, param_file, CS, restart_CSp)
if (CS%split) then
call register_restarts_dyn_split_RK2(dG%HI, GV, param_file, &
call register_restarts_dyn_split_RK2(HI, GV, param_file, &
CS%dyn_split_RK2_CSp, restart_CSp, CS%uh, CS%vh)
elseif (CS%use_RK2) then
call register_restarts_dyn_unsplit_RK2(dG%HI, GV, param_file, &
call register_restarts_dyn_unsplit_RK2(HI, GV, param_file, &
CS%dyn_unsplit_RK2_CSp, restart_CSp)
else
call register_restarts_dyn_unsplit(dG%HI, GV, param_file, &
call register_restarts_dyn_unsplit(HI, GV, param_file, &
CS%dyn_unsplit_CSp, restart_CSp)
endif

! This subroutine calls user-specified tracer registration routines.
! Additional calls can be added to MOM_tracer_flow_control.F90.
call call_tracer_register(dG%HI, GV, US, param_file, CS%tracer_flow_CSp, &
call call_tracer_register(HI, GV, US, param_file, CS%tracer_flow_CSp, &
CS%tracer_Reg, restart_CSp)

call MEKE_alloc_register_restart(dG%HI, param_file, CS%MEKE, restart_CSp)
call set_visc_register_restarts(dG%HI, GV, param_file, CS%visc, restart_CSp)
call mixedlayer_restrat_register_restarts(dG%HI, param_file, &
call MEKE_alloc_register_restart(HI, param_file, CS%MEKE, restart_CSp)
call set_visc_register_restarts(HI, GV, param_file, CS%visc, restart_CSp)
call mixedlayer_restrat_register_restarts(HI, param_file, &
CS%mixedlayer_restrat_CSp, restart_CSp)

if (CS%rotate_index .and. associated(OBC_in) .and. use_temperature) then
Expand Down Expand Up @@ -2410,33 +2421,16 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &

! This needs the number of tracers and to have called any code that sets whether
! reservoirs are used.
call open_boundary_register_restarts(dg%HI, GV, CS%OBC, CS%tracer_Reg, &
call open_boundary_register_restarts(HI, GV, CS%OBC, CS%tracer_Reg, &
param_file, restart_CSp, use_temperature)
endif

call callTree_waypoint("restart registration complete (initialize_MOM)")
call restart_registry_lock(restart_CSp)

! Shift from using the temporary dynamic grid type to using the final
! (potentially static) ocean-specific grid type.
! The next line would be needed if G%Domain had not already been init'd above:
! call clone_MOM_domain(dG%Domain, G%Domain)

! NOTE: If indices are rotated, then G and G_in must both be initialized.
! If not rotated, then G_in and G are the same grid.
if (CS%rotate_index) then
call MOM_grid_init(G, param_file, US, HI, bathymetry_at_vel=bathy_at_vel)
call copy_dyngrid_to_MOM_grid(dG, G, US)
call destroy_dyn_horgrid(dG)
endif
call MOM_grid_init(G_in, param_file, US, HI_in, bathymetry_at_vel=bathy_at_vel)
call copy_dyngrid_to_MOM_grid(dG_in, G_in, US)
if (.not. CS%rotate_index) G => G_in

! Write out all of the grid data used by this run.
new_sim = determine_is_new_run(dirs%input_filename, dirs%restart_input_dir, G_in, restart_CSp)
write_geom_files = ((write_geom==2) .or. ((write_geom==1) .and. new_sim))

! Write out all of the grid data used by this run.
if (write_geom_files) call write_ocean_geometry_file(dG_in, param_file, dirs%output_directory, US=US)

call destroy_dyn_horgrid(dG_in)
Expand Down Expand Up @@ -2562,16 +2556,16 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &

if (test_grid_copy) then
! Copy the data from the temporary grid to the dyn_hor_grid to CS%G.
call create_dyn_horgrid(dG, G%HI)
call clone_MOM_domain(G%Domain, dG%Domain)
call create_dyn_horgrid(test_dG, G%HI)
call clone_MOM_domain(G%Domain, test_dG%Domain)

call clone_MOM_domain(G%Domain, CS%G%Domain)
call MOM_grid_init(CS%G, param_file, US)

call copy_MOM_grid_to_dyngrid(G, dg, US)
call copy_dyngrid_to_MOM_grid(dg, CS%G, US)
call copy_MOM_grid_to_dyngrid(G, test_dG, US)
call copy_dyngrid_to_MOM_grid(test_dG, CS%G, US)

call destroy_dyn_horgrid(dG)
call destroy_dyn_horgrid(test_dG)
call MOM_grid_end(G) ; deallocate(G)

G => CS%G
Expand Down
Loading

0 comments on commit 91c9caa

Please sign in to comment.