Skip to content

Commit

Permalink
stochastic physics re-write
Browse files Browse the repository at this point in the history
  • Loading branch information
@pjpegion @kshedstrom
pjpegion authored and kshedstrom committed Apr 6, 2022
1 parent 0932b9e commit 199126e
Show file tree
Hide file tree
Showing 7 changed files with 204 additions and 164 deletions.
38 changes: 3 additions & 35 deletions config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -64,7 +64,6 @@ module MOM_ocean_model_nuopc
use MOM_surface_forcing_nuopc, only : forcing_save_restart
use MOM_domains, only : root_PE,num_PEs
use MOM_coms, only : Get_PElist
use stochastic_physics, only : init_stochastic_physics_ocn, run_stochastic_physics_ocn

#include <MOM_memory.h>

Expand Down Expand Up @@ -178,8 +177,8 @@ module MOM_ocean_model_nuopc
!! steps can span multiple coupled time steps.
logical :: diabatic_first !< If true, apply diabatic and thermodynamic
!! processes before time stepping the dynamics.
logical :: do_sppt !< If true, allocate array for SPPT
logical :: pert_epbl !< If true, allocate arrays for energetic PBL perturbations
logical,public :: do_sppt !< If true, write stochastic physics restarts
logical,public :: pert_epbl !< If true, write stochastic physics restarts

real :: eps_omesh !< Max allowable difference between ESMF mesh and MOM6
!! domain coordinates
Expand Down Expand Up @@ -255,12 +254,6 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i
!! min(HFrz, OBLD), where OBLD is the boundary layer depth.
!! If HFrz <= 0 (default), melt potential will not be computed.
logical :: use_melt_pot!< If true, allocate melt_potential array
! stochastic physics
integer :: mom_comm ! list of pes for this instance of the ocean
integer :: num_procs ! number of processors to pass to stochastic physics
integer :: iret ! return code from stochastic physics
integer :: me ! my pe
integer :: master ! root pe

! This include declares and sets the variable "version".
#include "version_variable.h"
Expand Down Expand Up @@ -432,7 +425,7 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i

endif

! get number of processors and PE list for stocasthci physics initialization
! check to see if stochastic physics is active
call get_param(param_file, mdl, "DO_SPPT", OS%do_sppt, &
"If true, then stochastically perturb the thermodynamic "//&
"tendemcies of T,S, amd h. Amplitude and correlations are "//&
Expand All @@ -443,27 +436,7 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i
"production and dissipation terms. Amplitude and correlations are "//&
"controlled by the nam_stoch namelist in the UFS model only.", &
default=.false.)
if (OS%do_sppt .OR. OS%pert_epbl) then
num_procs=num_PEs()
allocate(pelist(num_procs))
call Get_PElist(pelist,commID = mom_comm)
master=root_PE()

call init_stochastic_physics_ocn(OS%dt_therm,OS%grid%geoLonT,OS%grid%geoLatT,OS%grid%ied-OS%grid%isd+1,OS%grid%jed-OS%grid%jsd+1,OS%grid%ke,&
OS%pert_epbl,OS%do_sppt,master,mom_comm,iret)
if (iret/=0) then
write(6,*) 'call to init_stochastic_physics_ocn failed'
call MOM_error(FATAL, "stochastic physics in enambled in MOM6 but "// &
"not activated in stochastic_physics namelist ")
return
endif

if (OS%do_sppt) allocate(OS%fluxes%sppt_wts(OS%grid%isd:OS%grid%ied,OS%grid%jsd:OS%grid%jed))
if (OS%pert_epbl) then
allocate(OS%fluxes%epbl1_wts(OS%grid%isd:OS%grid%ied,OS%grid%jsd:OS%grid%jed))
allocate(OS%fluxes%epbl2_wts(OS%grid%isd:OS%grid%ied,OS%grid%jsd:OS%grid%jed))
endif
endif
call close_param_file(param_file)
call diag_mediator_close_registration(OS%diag)

Expand Down Expand Up @@ -635,11 +608,6 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, &
call disable_averaging(OS%diag)
Master_time = OS%Time ; Time1 = OS%Time

! update stochastic physics patterns before running next time-step
if (OS%do_sppt .OR. OS%pert_epbl ) then
call run_stochastic_physics_ocn(OS%fluxes%sppt_wts,OS%fluxes%epbl1_wts,OS%fluxes%epbl2_wts)
endif

if (OS%offline_tracer_mode) then
call step_offline(OS%forces, OS%fluxes, OS%sfc_state, Time1, dt_coupling, OS%MOM_CSp)
elseif ((.not.do_thermo) .or. (.not.do_dyn)) then
Expand Down
25 changes: 13 additions & 12 deletions src/core/MOM.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -403,16 +403,6 @@ module MOM
type(ODA_CS), pointer :: odaCS => NULL() !< a pointer to the control structure for handling
!! ensemble model state vectors and data assimilation
!! increments and priors
type(porous_barrier_ptrs) :: pbv !< porous barrier fractional cell metrics
real ALLOCABLE_, dimension(NIMEMB_PTR_,NJMEM_,NKMEM_) &
:: por_face_areaU !< fractional open area of U-faces [nondim]
real ALLOCABLE_, dimension(NIMEM_,NJMEMB_PTR_,NKMEM_) &
:: por_face_areaV !< fractional open area of V-faces [nondim]
real ALLOCABLE_, dimension(NIMEMB_PTR_,NJMEM_,NK_INTERFACE_) &
:: por_layer_widthU !< fractional open width of U-faces [nondim]
real ALLOCABLE_, dimension(NIMEM_,NJMEMB_PTR_,NK_INTERFACE_) &
:: por_layer_widthV !< fractional open width of V-faces [nondim]
type(particles), pointer :: particles => NULL() !<Lagrangian particles
type(stochastic_CS), pointer :: stoch_CS => NULL() !< a pointer to the stochastics control structure
end type MOM_control_struct

Expand Down Expand Up @@ -673,7 +663,7 @@ subroutine step_MOM(forces_in, fluxes_in, sfc_state, Time_start, time_int_in, CS
call disable_averaging(CS%diag)
endif
endif
! advance the random pattern if stochastic physics is active

if (CS%stoch_CS%do_sppt .OR. CS%stoch_CS%pert_epbl) call update_stochastics(CS%stoch_CS)

if (do_dyn) then
Expand Down Expand Up @@ -825,6 +815,7 @@ subroutine step_MOM(forces_in, fluxes_in, sfc_state, Time_start, time_int_in, CS
enddo ; enddo
endif


call step_MOM_dynamics(forces, CS%p_surf_begin, CS%p_surf_end, dt, &
dt_therm_here, bbl_time_int, CS, &
Time_local, Waves=Waves)
Expand Down Expand Up @@ -1391,7 +1382,7 @@ subroutine step_MOM_thermo(CS, G, GV, US, u, v, h, tv, fluxes, dtdia, &
call cpu_clock_begin(id_clock_diabatic)

call diabatic(u, v, h, tv, CS%Hml, fluxes, CS%visc, CS%ADp, CS%CDp, dtdia, &
Time_end_thermo, G, GV, US, CS%diabatic_CSp, OBC=CS%OBC, Waves=Waves)
Time_end_thermo, G, GV, US, CS%diabatic_CSp, CS%stoch_CS,OBC=CS%OBC, Waves=Waves)
fluxes%fluxes_used = .true.

if (showCallTree) call callTree_waypoint("finished diabatic (step_MOM_thermo)")
Expand Down Expand Up @@ -2500,6 +2491,16 @@ subroutine initialize_MOM(Time, Time_init, param_file, dirs, CS, restart_CSp, &
else
call set_first_direction(G, modulo(first_direction, 2))
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)
call destroy_dyn_horgrid(dG_in)

if (.not. CS%rotate_index) &
G => G_in
! initialize stochastic physics
!call stochastics_init(CS%dt_therm, CS%G, CS%GV, CS%stoch_CS, param_file, diag, Time)
! Set a few remaining fields that are specific to the ocean grid type.
call set_first_direction(G, first_direction)
! Allocate the auxiliary non-symmetric domain for debugging or I/O purposes.
if (CS%debug .or. G%symmetric) then
call clone_MOM_domain(G%Domain, G%Domain_aux, symmetric=.false.)
Expand Down
2 changes: 0 additions & 2 deletions src/core/MOM_forcing_type.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -170,8 +170,6 @@ module MOM_forcing_type
!! exactly 0 away from shelves or on land.
real, pointer, dimension(:,:) :: iceshelf_melt => NULL() !< Ice shelf melt rate (positive)
!! or freezing (negative) [R Z T-1 ~> kg m-2 s-1]
real, pointer, dimension(:,:) :: shelf_sfc_mass_flux => NULL() !< Ice shelf surface mass flux
!! deposition from the atmosphere. [R Z T-1 ~> kg m-2 s-1]

! Scalars set by surface forcing modules
real :: vPrecGlobalAdj = 0. !< adjustment to restoring vprec to zero out global net [kg m-2 s-1]
Expand Down
50 changes: 30 additions & 20 deletions src/parameterizations/stochastic/MOM_stochastics.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,13 @@ module MOM_stochastics
! This file is part of MOM6. See LICENSE.md for the license.

! This is the top level module for the MOM6 ocean model. It contains routines
! for initialization, update, and writing restart of stochastic physics. This
! for initialization, termination and update of ocean model state. This
! particular version wraps all of the calls for MOM6 in the calls that had
! been used for MOM4.
!
! This code is a stop-gap wrapper of the MOM6 code to enable it to be called
! in the same way as MOM4.

use MOM_diag_mediator, only : register_diag_field, diag_ctrl, time_type
use MOM_grid, only : ocean_grid_type
use MOM_verticalGrid, only : verticalGrid_type
Expand All @@ -26,13 +29,11 @@ module MOM_stochastics

public stochastics_init, update_stochastics

!> This control structure holds parameters for the MOM_stochastics module
type, public:: stochastic_CS
logical :: do_sppt !< If true, stochastically perturb the diabatic
logical :: pert_epbl !< If true, then randomly perturb the KE dissipation and genration terms
integer :: id_sppt_wts = -1 !< Diagnostic id for SPPT
integer :: id_epbl1_wts=-1 !< Diagnostic id for epbl generation perturbation
integer :: id_epbl2_wts=-1 !< Diagnostic id for epbl dissipation perturbation
logical :: do_sppt !< If true, stochastically perturb the diabatic
logical :: pert_epbl !! If true, then randomly perturb the KE dissipation and genration terms
integer :: id_sppt_wts = -1
integer :: id_epbl1_wts=-1,id_epbl2_wts=-1
! stochastic patterns
real, allocatable :: sppt_wts(:,:) !< Random pattern for ocean SPPT
!! tendencies with a number between 0 and 2
Expand All @@ -42,14 +43,22 @@ module MOM_stochastics
type(time_type), pointer :: Time !< Pointer to model time (needed for sponges)
end type stochastic_CS

!> This type is used for communication with other components via the FMS coupler.
!! The element names and types can be changed only with great deliberation, hence
!! the persistnce of things like the cutsy element name "avg_kount".
contains

!! This subroutine initializes the stochastics physics control structure.
!> ocean_model_init initializes the ocean model, including registering fields
!! for restarts and reading restart files if appropriate.
!!
!! This subroutine initializes both the ocean state and the ocean surface type.
!! Because of the way that indicies and domains are handled, Ocean_sfc must have
!! been used in a previous call to cean_type.
subroutine stochastics_init(dt, grid, GV, CS, param_file, diag, Time)
real, intent(in) :: dt !< time step [T ~> s]
type(ocean_grid_type), intent(in) :: grid !< horizontal grid information
type(verticalGrid_type), intent(in) :: GV !< vertical grid structure
type(stochastic_CS), pointer, intent(inout):: CS !< stochastic control structure
type(ocean_grid_type), intent(in) :: grid ! horizontal grid information
type(verticalGrid_type), intent(in) :: GV ! vertical grid structure
type(stochastic_CS), pointer, intent(inout):: CS
type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
type(diag_ctrl), target, intent(inout) :: diag !< structure to regulate diagnostic output
type(time_type), target :: Time !< model time
Expand All @@ -59,7 +68,7 @@ subroutine stochastics_init(dt, grid, GV, CS, param_file, diag, Time)
integer :: num_procs ! number of processors to pass to stochastic physics
integer :: iret ! return code from stochastic physics
integer :: me ! my pe
integer :: pe_zero ! root pe
integer :: master ! root pe
integer :: nx ! number of x-points including halo
integer :: ny ! number of x-points including halo

Expand Down Expand Up @@ -95,13 +104,15 @@ subroutine stochastics_init(dt, grid, GV, CS, param_file, diag, Time)
num_procs=num_PEs()
allocate(pelist(num_procs))
call Get_PElist(pelist,commID = mom_comm)
pe_zero=root_PE()
nx = grid%ied - grid%isd + 1
ny = grid%jed - grid%jsd + 1
master=root_PE()
nx=grid%ied-grid%isd+1
ny=grid%jed-grid%jsd+1
call init_stochastic_physics_ocn(dt,grid%geoLonT,grid%geoLatT,nx,ny,GV%ke, &
CS%pert_epbl,CS%do_sppt,pe_zero,mom_comm,iret)
CS%pert_epbl,CS%do_sppt,master,mom_comm,iret)
if (iret/=0) then
call MOM_error(FATAL, "call to init_stochastic_physics_ocn failed")
write(6,*) 'call to init_stochastic_physics_ocn failed'
call MOM_error(FATAL, "stochastic physics in enambled in MOM6 but "// &
"not activated in stochastic_physics namelist ")
return
endif

Expand All @@ -122,7 +133,6 @@ subroutine stochastics_init(dt, grid, GV, CS, param_file, diag, Time)
write(*,'(/12x,a/)') '=== COMPLETED MOM STOCHASTIC INITIALIZATION ====='

call callTree_leave("ocean_model_init(")
return
end subroutine stochastics_init

!> update_ocean_model uses the forcing in Ice_ocean_boundary to advance the
Expand All @@ -135,9 +145,9 @@ subroutine update_stochastics(CS)
call callTree_enter("update_stochastics(), MOM_stochastics.F90")

! update stochastic physics patterns before running next time-step
call run_stochastic_physics_ocn(CS%sppt_wts,CS%epbl1_wts,CS%epbl2_wts)
call run_stochastic_physics_ocn(CS%sppt_wts,CS%epbl1_wts,CS%epbl2_wts)
print*,'in update_stoch',minval(CS%sppt_wts),maxval(CS%sppt_wts),minval(CS%epbl1_wts),maxval(CS%epbl1_wts)

return
end subroutine update_stochastics

end module MOM_stochastics
Expand Down
64 changes: 64 additions & 0 deletions src/parameterizations/stochastic/MOM_stochastics_stub.F90
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,64 @@
!> Top-level module for the MOM6 ocean model in coupled mode.
module MOM_stochastics

! This file is part of MOM6. See LICENSE.md for the license.

! This is the top level module for the MOM6 ocean model. It contains routines
! for initialization, termination and update of ocean model state. This
! particular version wraps all of the calls for MOM6 in the calls that had
! been used for MOM4.
!
! This code is a stop-gap wrapper of the MOM6 code to enable it to be called
! in the same way as MOM4.

use MOM_diag_mediator, only : register_diag_field, diag_ctrl, time_type
use MOM_grid, only : ocean_grid_type
use MOM_verticalGrid, only : verticalGrid_type
use MOM_error_handler, only : MOM_error, FATAL, WARNING, is_root_pe
use MOM_error_handler, only : callTree_enter, callTree_leave
use MOM_file_parser, only : get_param, log_version, close_param_file, param_file_type
use mpp_domains_mod, only : domain2d, mpp_get_layout, mpp_get_global_domain
use mpp_domains_mod, only : mpp_define_domains, mpp_get_compute_domain, mpp_get_data_domain
use MOM_domains, only : root_PE,num_PEs
use MOM_coms, only : Get_PElist

#include <MOM_memory.h>

implicit none ; private

public stochastics_init, update_stochastics

type, public:: stochastic_CS
logical :: do_sppt !< If true, stochastically perturb the diabatic
logical :: pert_epbl !! If true, then randomly perturb the KE dissipation and genration terms
integer :: id_sppt_wts = -1
integer :: id_epbl1_wts=-1,id_epbl2_wts=-1
! stochastic patterns
real, allocatable :: sppt_wts(:,:) !< Random pattern for ocean SPPT
!! tendencies with a number between 0 and 2
real, allocatable :: epbl1_wts(:,:) !< Random pattern for K.E. generation
real, allocatable :: epbl2_wts(:,:) !< Random pattern for K.E. dissipation
type(diag_ctrl), pointer :: diag !< structure used to regulate timing of diagnostic output
type(time_type), pointer :: Time !< Pointer to model time (needed for sponges)
end type stochastic_CS

contains

subroutine stochastics_init(dt, grid, GV, CS, param_file, diag, Time)
real, intent(in) :: dt !< time step [T ~> s]
type(ocean_grid_type), intent(in) :: grid ! horizontal grid information
type(verticalGrid_type), intent(in) :: GV ! vertical grid structure
type(stochastic_CS), pointer, intent(inout):: CS
type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
type(diag_ctrl), target, intent(inout) :: diag !< structure to regulate diagnostic output
type(time_type), target :: Time !< model time
return
end subroutine stochastics_init

subroutine update_stochastics(CS)
type(stochastic_CS), intent(inout) :: CS !< diabatic control structure
return
end subroutine update_stochastics

end module MOM_stochastics

Loading

0 comments on commit 199126e

Please sign in to comment.