place stochastic array in fluxes container and make SPPT specific arr…

…ays allocatable

config_src/nuopc_driver/mom_cap.F90

@@ -1379,6 +1379,8 @@ subroutine ModelAdvance(gcomp, rc)
   character(len=*),parameter             :: subname='(MOM_cap:ModelAdvance)'
   character(len=8)                       :: suffix
   integer                                :: num_rest_files
+  logical                                :: do_sppt = .false.
+  logical                                :: pert_epbl = .false.
 
   rc = ESMF_SUCCESS
   if(profile_memory) call ESMF_VMLogMemInfo("Entering MOM Model_ADVANCE: ")
@@ -1587,15 +1589,17 @@ subroutine ModelAdvance(gcomp, rc)
         ! write restart file(s)
         call ocean_model_restart(ocean_state, restartname=restartname)
 
-        ! write stochastic physics restart file if active
-        if (ESMF_AlarmIsRinging(stop_alarm, rc=rc)) then
-           write(restartname,'(A)')"ocn_stoch.res.nc"
-        else
-           write(restartname,'(A,I4.4,"-",I2.2,"-",I2.2,"-",I2.2,"-",I2.2,"-",I2.2,A)') &
-                "ocn_stoch.res.", year, month, day, hour, minute, seconds,".nc"
+        if (ocean_state%do_sppt .OR. ocean_state%pert_epbl) then
+          if (ESMF_AlarmIsRinging(stop_alarm, rc=rc)) then
+             write(restartname,'(A)')"ocn_stoch.res.nc"
+          else
+             write(restartname,'(A,I4.4,"-",I2.2,"-",I2.2,"-",I2.2,"-",I2.2,"-",I2.2,A)') &
+                  "ocn_stoch.res.", year, month, day, hour, minute, seconds,".nc"
+          endif
+          call ESMF_LogWrite("MOM_cap: Writing stoch restart :  "//trim(restartname), &
+                             ESMF_LOGMSG_INFO)
+          call write_stoch_restart_ocn('RESTART/'//trim(restartname))
         endif
-        call ESMF_LogWrite("MOM_cap: Writing restart :  "//trim(restartname), ESMF_LOGMSG_INFO)
-        call write_stoch_restart_ocn('RESTART/'//trim(restartname))
      endif
 
      if (is_root_pe()) then

config_src/nuopc_driver/mom_ocean_model_nuopc.F90

@@ -63,6 +63,7 @@ module MOM_ocean_model_nuopc
 use MOM_surface_forcing_nuopc, only : ice_ocean_boundary_type, surface_forcing_CS
 use MOM_surface_forcing_nuopc, only : forcing_save_restart
 use MOM_domains,               only : root_PE,PE_here,num_PEs
+use MOM_coms,                  only : Get_PElist
 use stochastic_physics,        only : init_stochastic_physics_ocn, run_stochastic_physics_ocn
 
 #include <MOM_memory.h>
@@ -176,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, allocate array for SPPT
+  logical,public :: pert_epbl !< If true, allocate arrays for energetic PBL perturbations
 
   real :: eps_omesh           !< Max allowable difference between ESMF mesh and MOM6
                               !! domain coordinates
@@ -253,6 +254,7 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i
                       !! If HFrz <= 0 (default), melt potential will not be computed.
   logical :: use_melt_pot!< If true, allocate melt_potential array
 ! stochastic physics
+  integer,allocatable :: pelist(:) ! list of pes for this instance of the ocean
   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
@@ -439,7 +441,9 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i
                  default=.false.)
   if (OS%do_sppt .OR. OS%pert_epbl) then
      num_procs=num_PEs()
-     call mpp_get_pelist(Ocean_sfc%domain, mom_comm)
+     allocate(pelist(num_procs))
+     !call mpp_get_pelist(pelist, commID=mom_comm)
+     call Get_PElist(pelist,commID = mom_comm)
      me=PE_here()
      master=root_PE()
 

src/core/MOM.F90

@@ -123,7 +123,7 @@ module MOM
 use MOM_variables,             only : surface, allocate_surface_state, deallocate_surface_state
 use MOM_variables,             only : thermo_var_ptrs, vertvisc_type
 use MOM_variables,             only : accel_diag_ptrs, cont_diag_ptrs, ocean_internal_state
-use MOM_variables,             only : rotate_surface_state,stochastic_pattern
+use MOM_variables,             only : rotate_surface_state
 use MOM_verticalGrid,          only : verticalGrid_type, verticalGridInit, verticalGridEnd
 use MOM_verticalGrid,          only : fix_restart_scaling
 use MOM_verticalGrid,          only : get_thickness_units, get_flux_units, get_tr_flux_units

src/parameterizations/vertical/MOM_diabatic_driver.F90

@@ -289,9 +289,9 @@ subroutine diabatic(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Time_end, &
   integer :: i, j, k, m, is, ie, js, je, nz
   logical :: showCallTree ! If true, show the call tree
 
-  real, allocatable(:,:,:)    :: h_in  ! thickness before thermodynamics
-  real, allocatable(:,:,:)    :: t_in  ! temperature before thermodynamics
-  real, allocatable(:,:,:)    :: s_in  ! salinity before thermodynamics
+  real, allocatable, dimension(:,:,:)    :: h_in  ! thickness before thermodynamics
+  real, allocatable, dimension(:,:,:)    :: t_in  ! temperature before thermodynamics
+  real, allocatable, dimension(:,:,:)    :: s_in  ! salinity before thermodynamics
   real :: t_tend,s_tend,h_tend                  ! holder for tendencey needed for SPPT
   real :: t_pert,s_pert,h_pert                  ! holder for perturbations needed for SPPT
 
@@ -302,9 +302,9 @@ subroutine diabatic(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Time_end, &
     allocate(h_in(G%isd:G%ied, G%jsd:G%jed,G%ke))
     allocate(t_in(G%isd:G%ied, G%jsd:G%jed,G%ke))
     allocate(s_in(G%isd:G%ied, G%jsd:G%jed,G%ke))
-    h_in(:,:) = h(:,:)
-    t_in(:,:) = tv%T(:,:)
-    s_in(:,:) = tv%S(:,:)
+    h_in(:,:,:)=h(:,:,:)
+    t_in(:,:,:)=tv%T(:,:,:)
+    s_in(:,:,:)=tv%S(:,:,:)
 
     if (CS%id_sppt_wts > 0) then
       call post_data(CS%id_sppt_wts, fluxes%sppt_wts, CS%diag)
@@ -458,6 +458,7 @@ subroutine diabatic(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Time_end, &
   if (CS%debugConservation) call MOM_state_stats('leaving diabatic', u, v, h, tv%T, tv%S, G, GV, US)
 
   if (CS%do_sppt) then
+  ! perturb diabatic tendecies
     do k=1,nz
       do j=js,je
         do i=is,ie
@@ -468,17 +469,20 @@ subroutine diabatic(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Time_end, &
           t_pert=t_tend+t_in(i,j,k)
           s_pert=s_tend+s_in(i,j,k)
           if (h_pert > GV%Angstrom_H) then
-             h(i,j,k)=h_pert
+            h(i,j,k) = h_pert
           else
-             h(i,j,k)=GV%Angstrom_H
+            h(i,j,k) = GV%Angstrom_H
           endif
-          tv%T(i,j,k)=t_pert
+          tv%T(i,j,k) = t_pert
           if (s_pert > 0.0) then
-             tv%S(i,j,k)=s_pert
+            tv%S(i,j,k) = s_pert
           endif
         enddo
       enddo
     enddo
+    deallocate(h_in)
+    deallocate(t_in)
+    deallocate(s_in)
   endif
 
 end subroutine diabatic

src/parameterizations/vertical/MOM_energetic_PBL.F90

@@ -461,7 +461,8 @@ subroutine energetic_PBL(h_3d, u_3d, v_3d, tv, fluxes, dt, Kd_int, G, GV, US, CS
         call ePBL_column(h, u, v, T0, S0, dSV_dT_1d, dSV_dS_1d, TKE_forcing, B_flux, absf, &
                          u_star, u_star_mean, dt, MLD_io, Kd, mixvel, mixlen, GV, &
                          US, CS, eCD, dt_diag=dt_diag, Waves=Waves, G=G, &
-                         epbl1_wt=epbl1_wts(i,j),epbl2_wt=epbl2_wts(i,j),i=i, j=j)
+                         epbl1_wt=fluxes%epbl1_wts(i,j),epbl2_wt=fluxes%epbl2_wts(i,j), &
+                         i=i, j=j)
       else
         call ePBL_column(h, u, v, T0, S0, dSV_dT_1d, dSV_dS_1d, TKE_forcing, B_flux, absf, &
                          u_star, u_star_mean, dt, MLD_io, Kd, mixvel, mixlen, GV, &
@@ -540,8 +541,8 @@ subroutine energetic_PBL(h_3d, u_3d, v_3d, tv, fluxes, dt, Kd_int, G, GV, US, CS
     ! only write random patterns if running with stochastic physics, otherwise the
     ! array is unallocated and will give an error
     if  (CS%pert_epbl) then
-      if (CS%id_epbl1_wts > 0)    call post_data(CS%id_epbl1_wts, stochastics%epbl1_wts, CS%diag)
-      if (CS%id_epbl2_wts > 0)    call post_data(CS%id_epbl2_wts, stochastics%epbl2_wts, CS%diag)
+      if (CS%id_epbl1_wts > 0)    call post_data(CS%id_epbl1_wts, fluxes%epbl1_wts, CS%diag)
+      if (CS%id_epbl2_wts > 0)    call post_data(CS%id_epbl2_wts, fluxes%epbl2_wts, CS%diag)
     endif
   endif