Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merge gfdl updates #2

Merged
merged 6 commits into from
Nov 17, 2021
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
39 changes: 24 additions & 15 deletions config_src/drivers/unit_drivers/MOM_sum_driver.F90
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
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