MOM6
mom_ocean_model_mct Module Reference

Detailed Description

Top-level module for the MOM6 ocean model in coupled mode.

Data Types

type  ocean_public_type
 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". More...
 
type  ocean_state_type
 The ocean_state_type contains all information about the state of the ocean, with a format that is private so it can be readily changed without disrupting other coupled components. More...
 

Functions/Subroutines

subroutine, public ocean_model_init (Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, input_restart_file)
 ocean_model_init initializes the ocean model, including registering fields for restarts and reading restart files if appropriate. More...
 
subroutine, public update_ocean_model (Ice_ocean_boundary, OS, Ocean_sfc, time_start_update, Ocean_coupling_time_step, update_dyn, update_thermo, Ocn_fluxes_used)
 update_ocean_model uses the forcing in Ice_ocean_boundary to advance the ocean model's state from the input value of Ocean_state (which must be for time time_start_update) for a time interval of Ocean_coupling_time_step, returning the publicly visible ocean surface properties in Ocean_sfc and storing the new ocean properties in Ocean_state. More...
 
subroutine, public ocean_model_restart (OS, timestamp, restartname)
 This subroutine writes out the ocean model restart file. More...
 
subroutine, public ocean_model_end (Ocean_sfc, Ocean_state, Time)
 ocean_model_end terminates the model run, saving the ocean state in a restart and deallocating any data associated with the ocean. More...
 
subroutine, public ocean_model_save_restart (OS, Time, directory, filename_suffix)
 ocean_model_save_restart causes restart files associated with the ocean to be written out. More...
 
subroutine initialize_ocean_public_type (input_domain, Ocean_sfc, diag, maskmap, gas_fields_ocn)
 Initialize the public ocean type. More...
 
subroutine convert_state_to_ocean_type (sfc_state, Ocean_sfc, G, US, patm, press_to_z)
 This subroutine translates the coupler's ocean_data_type into MOM's surface state variable. This may eventually be folded into the MOM code that calculates the surface state in the first place. Note the offset in the arrays because the ocean_data_type has no halo points in its arrays and always uses absolute indicies. More...
 
subroutine, public ocean_model_init_sfc (OS, Ocean_sfc)
 This subroutine extracts the surface properties from the ocean's internal state and stores them in the ocean type returned to the calling ice model. It has to be separate from the ocean_initialization call because the coupler module allocates the space for some of these variables. More...
 
subroutine, public ocean_model_flux_init (OS, verbosity)
 ocean_model_flux_init is used to initialize properties of the air-sea fluxes as determined by various run-time parameters. It can be called from non-ocean PEs, or PEs that have not yet been initialzed, and it can safely be called multiple times. More...
 
subroutine, public ocean_stock_pe (OS, index, value, time_index)
 Ocean_stock_pe - returns the integrated stocks of heat, water, etc. for conservation checks. Because of the way FMS is coded, only the root PE has the integrated amount, while all other PEs get 0. More...
 
subroutine, public ocean_public_type_chksum (id, timestep, ocn)
 Write out FMS-format checsums on fields from the ocean surface state. More...
 
subroutine, public get_ocean_grid (OS, Gridp)
 

Function/Subroutine Documentation

◆ convert_state_to_ocean_type()

subroutine mom_ocean_model_mct::convert_state_to_ocean_type ( type(surface), intent(inout)  sfc_state,
type(ocean_public_type), intent(inout), target  Ocean_sfc,
type(ocean_grid_type), intent(inout)  G,
type(unit_scale_type), intent(in)  US,
real, dimension(:,:), intent(in), optional  patm,
real, intent(in), optional  press_to_z 
)

This subroutine translates the coupler's ocean_data_type into MOM's surface state variable. This may eventually be folded into the MOM code that calculates the surface state in the first place. Note the offset in the arrays because the ocean_data_type has no halo points in its arrays and always uses absolute indicies.

Parameters
[in,out]sfc_stateA structure containing fields that describe the surface state of the ocean.
[in,out]ocean_sfcA structure containing various publicly
[in,out]gThe ocean's grid structure
[in]usA dimensional unit scaling type
[in]patmThe pressure at the ocean surface, in Pa.
[in]press_to_zA conversion factor between pressure and ocean depth in m, usually 1/(rho_0*g), in m Pa-1.

Definition at line 844 of file mom_ocean_model_mct.F90.

844  type(surface), intent(inout) :: sfc_state !< A structure containing fields that
845  !! describe the surface state of the ocean.
846  type(ocean_public_type), &
847  target, intent(inout) :: Ocean_sfc !< A structure containing various publicly
848  !! visible ocean surface fields, whose elements
849  !! have their data set here.
850  type(ocean_grid_type), intent(inout) :: G !< The ocean's grid structure
851  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
852  real, optional, intent(in) :: patm(:,:) !< The pressure at the ocean surface, in Pa.
853  real, optional, intent(in) :: press_to_z !< A conversion factor between pressure and
854  !! ocean depth in m, usually 1/(rho_0*g), in m Pa-1.
855 
856  ! Local variables
857  real :: IgR0
858  character(len=48) :: val_str
859  integer :: isc_bnd, iec_bnd, jsc_bnd, jec_bnd
860  integer :: i, j, i0, j0, is, ie, js, je
861 
862  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
863  call pass_vector(sfc_state%u,sfc_state%v,g%Domain)
864 
865  call mpp_get_compute_domain(ocean_sfc%Domain, isc_bnd, iec_bnd, &
866  jsc_bnd, jec_bnd)
867  if (present(patm)) then
868  ! Check that the inidicies in patm are (isc_bnd:iec_bnd,jsc_bnd:jec_bnd).
869  if (.not.present(press_to_z)) call mom_error(fatal, &
870  'convert_state_to_ocean_type: press_to_z must be present if patm is.')
871  endif
872 
873  i0 = is - isc_bnd ; j0 = js - jsc_bnd
874  if (sfc_state%T_is_conT) then
875  ! Convert the surface T from conservative T to potential T.
876  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
877  ocean_sfc%t_surf(i,j) = gsw_pt_from_ct(sfc_state%SSS(i+i0,j+j0), &
878  sfc_state%SST(i+i0,j+j0)) + celsius_kelvin_offset
879  enddo ; enddo
880  else
881  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
882  ocean_sfc%t_surf(i,j) = sfc_state%SST(i+i0,j+j0) + celsius_kelvin_offset
883  enddo ; enddo
884  endif
885  if (sfc_state%S_is_absS) then
886  ! Convert the surface S from absolute salinity to practical salinity.
887  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
888  ocean_sfc%s_surf(i,j) = gsw_sp_from_sr(sfc_state%SSS(i+i0,j+j0))
889  enddo ; enddo
890  else
891  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
892  ocean_sfc%s_surf(i,j) = sfc_state%SSS(i+i0,j+j0)
893  enddo ; enddo
894  endif
895 
896  if (present(patm)) then
897  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
898  ocean_sfc%sea_lev(i,j) = sfc_state%sea_lev(i+i0,j+j0) + patm(i,j) * press_to_z
899  ocean_sfc%area(i,j) = us%L_to_m**2*g%areaT(i+i0,j+j0)
900  enddo ; enddo
901  else
902  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
903  ocean_sfc%sea_lev(i,j) = sfc_state%sea_lev(i+i0,j+j0)
904  ocean_sfc%area(i,j) = us%L_to_m**2*g%areaT(i+i0,j+j0)
905  enddo ; enddo
906  endif
907 
908  if (associated(sfc_state%frazil)) then
909  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
910  ocean_sfc%frazil(i,j) = sfc_state%frazil(i+i0,j+j0)
911  enddo ; enddo
912  endif
913 
914  if (allocated(sfc_state%melt_potential)) then
915  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
916  ocean_sfc%melt_potential(i,j) = sfc_state%melt_potential(i+i0,j+j0)
917  enddo ; enddo
918  endif
919 
920  if (allocated(sfc_state%Hml)) then
921  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
922  ocean_sfc%OBLD(i,j) = sfc_state%Hml(i+i0,j+j0)
923  enddo ; enddo
924  endif
925 
926  if (ocean_sfc%stagger == agrid) then
927  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
928  ocean_sfc%u_surf(i,j) = g%mask2dT(i+i0,j+j0) * &
929  0.5*(sfc_state%u(i+i0,j+j0)+sfc_state%u(i-1+i0,j+j0))
930  ocean_sfc%v_surf(i,j) = g%mask2dT(i+i0,j+j0) * &
931  0.5*(sfc_state%v(i+i0,j+j0)+sfc_state%v(i+i0,j-1+j0))
932  enddo ; enddo
933  elseif (ocean_sfc%stagger == bgrid_ne) then
934  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
935  ocean_sfc%u_surf(i,j) = g%mask2dBu(i+i0,j+j0) * &
936  0.5*(sfc_state%u(i+i0,j+j0)+sfc_state%u(i+i0,j+j0+1))
937  ocean_sfc%v_surf(i,j) = g%mask2dBu(i+i0,j+j0) * &
938  0.5*(sfc_state%v(i+i0,j+j0)+sfc_state%v(i+i0+1,j+j0))
939  enddo ; enddo
940  elseif (ocean_sfc%stagger == cgrid_ne) then
941  do j=jsc_bnd,jec_bnd ; do i=isc_bnd,iec_bnd
942  ocean_sfc%u_surf(i,j) = g%mask2dCu(i+i0,j+j0)*sfc_state%u(i+i0,j+j0)
943  ocean_sfc%v_surf(i,j) = g%mask2dCv(i+i0,j+j0)*sfc_state%v(i+i0,j+j0)
944  enddo ; enddo
945  else
946  write(val_str, '(I8)') ocean_sfc%stagger
947  call mom_error(fatal, "convert_state_to_ocean_type: "//&
948  "Ocean_sfc%stagger has the unrecognized value of "//trim(val_str))
949  endif
950 
951  if (coupler_type_initialized(sfc_state%tr_fields)) then
952  if (.not.coupler_type_initialized(ocean_sfc%fields)) then
953  call mom_error(fatal, "convert_state_to_ocean_type: "//&
954  "Ocean_sfc%fields has not been initialized.")
955  endif
956  call coupler_type_copy_data(sfc_state%tr_fields, ocean_sfc%fields)
957  endif
958 

References mom_constants::celsius_kelvin_offset.

Referenced by ocean_model_init(), ocn_comp_mct::ocean_model_init_sfc(), ocean_model_init_sfc(), and update_ocean_model().

Here is the caller graph for this function:

◆ get_ocean_grid()

subroutine, public mom_ocean_model_mct::get_ocean_grid ( type(ocean_state_type OS,
type(ocean_grid_type), pointer  Gridp 
)

Definition at line 1073 of file mom_ocean_model_mct.F90.

1073  ! Obtain the ocean grid.
1074  type(ocean_state_type) :: OS
1075  type(ocean_grid_type) , pointer :: Gridp
1076 
1077  gridp => os%grid
1078  return

◆ initialize_ocean_public_type()

subroutine mom_ocean_model_mct::initialize_ocean_public_type ( type(domain2d), intent(in)  input_domain,
type(ocean_public_type), intent(inout)  Ocean_sfc,
type(diag_ctrl), intent(in)  diag,
logical, dimension(:,:), intent(in), optional  maskmap,
type(coupler_1d_bc_type), intent(in), optional  gas_fields_ocn 
)

Initialize the public ocean type.

Parameters
[in]input_domainThe ocean model domain description
[in,out]ocean_sfcA structure containing various publicly visible ocean surface properties after initialization, whose elements are allocated here.
[in]diagA structure that regulates diagnsotic output
[in]maskmapA mask indicating which virtual processors
[in]gas_fields_ocnIf present, this type describes the

Definition at line 783 of file mom_ocean_model_mct.F90.

783  type(domain2D), intent(in) :: input_domain !< The ocean model domain description
784  type(ocean_public_type), intent(inout) :: Ocean_sfc !< A structure containing various publicly
785  !! visible ocean surface properties after initialization, whose
786  !! elements are allocated here.
787  type(diag_ctrl), intent(in) :: diag !< A structure that regulates diagnsotic output
788  logical, dimension(:,:), &
789  optional, intent(in) :: maskmap !< A mask indicating which virtual processors
790  !! are actually in use. If missing, all are used.
791  type(coupler_1d_bc_type), &
792  optional, intent(in) :: gas_fields_ocn !< If present, this type describes the
793  !! ocean and surface-ice fields that will participate
794  !! in the calculation of additional gas or other
795  !! tracer fluxes.
796  integer :: xsz, ysz, layout(2)
797  ! ice-ocean-boundary fields are always allocated using absolute indicies
798  ! and have no halos.
799  integer :: isc, iec, jsc, jec
800 
801  call mpp_get_layout(input_domain,layout)
802  call mpp_get_global_domain(input_domain, xsize=xsz, ysize=ysz)
803  if(PRESENT(maskmap)) then
804  call mpp_define_domains((/1,xsz,1,ysz/),layout,ocean_sfc%Domain, maskmap=maskmap)
805  else
806  call mpp_define_domains((/1,xsz,1,ysz/),layout,ocean_sfc%Domain)
807  endif
808  call mpp_get_compute_domain(ocean_sfc%Domain, isc, iec, jsc, jec)
809 
810  allocate (ocean_sfc%t_surf (isc:iec,jsc:jec), &
811  ocean_sfc%s_surf (isc:iec,jsc:jec), &
812  ocean_sfc%u_surf (isc:iec,jsc:jec), &
813  ocean_sfc%v_surf (isc:iec,jsc:jec), &
814  ocean_sfc%sea_lev(isc:iec,jsc:jec), &
815  ocean_sfc%area (isc:iec,jsc:jec), &
816  ocean_sfc%OBLD (isc:iec,jsc:jec), &
817  ocean_sfc%melt_potential(isc:iec,jsc:jec), &
818  ocean_sfc%frazil (isc:iec,jsc:jec))
819 
820  ocean_sfc%t_surf = 0.0 ! time averaged sst (Kelvin) passed to atmosphere/ice model
821  ocean_sfc%s_surf = 0.0 ! time averaged sss (psu) passed to atmosphere/ice models
822  ocean_sfc%u_surf = 0.0 ! time averaged u-current (m/sec) passed to atmosphere/ice models
823  ocean_sfc%v_surf = 0.0 ! time averaged v-current (m/sec) passed to atmosphere/ice models
824  ocean_sfc%sea_lev = 0.0 ! time averaged thickness of top model grid cell (m) plus patm/rho0/grav
825  ocean_sfc%frazil = 0.0 ! time accumulated frazil (J/m^2) passed to ice model
826  ocean_sfc%melt_potential = 0.0 ! time accumulated melt potential (J/m^2) passed to ice model
827  ocean_sfc%OBLD = 0.0 ! ocean boundary layer depth, in m
828  ocean_sfc%area = 0.0
829  ocean_sfc%axes = diag%axesT1%handles !diag axes to be used by coupler tracer flux diagnostics
830 
831  if (present(gas_fields_ocn)) then
832  call coupler_type_spawn(gas_fields_ocn, ocean_sfc%fields, (/isc,isc,iec,iec/), &
833  (/jsc,jsc,jec,jec/), suffix = '_ocn', as_needed=.true.)
834  endif
835 

Referenced by ocean_model_init().

Here is the caller graph for this function:

◆ ocean_model_end()

subroutine, public mom_ocean_model_mct::ocean_model_end ( type(ocean_public_type), intent(inout)  Ocean_sfc,
type(ocean_state_type), pointer  Ocean_state,
type(time_type), intent(in)  Time 
)

ocean_model_end terminates the model run, saving the ocean state in a restart and deallocating any data associated with the ocean.

Parameters
[in,out]ocean_sfcAn ocean_public_type structure that is to be deallocated upon termination.
ocean_stateA pointer to the structure containing the internal ocean state to be deallocated upon termination.
[in]timeThe model time, used for writing restarts.

Definition at line 729 of file mom_ocean_model_mct.F90.

729  type(ocean_public_type), intent(inout) :: Ocean_sfc !< An ocean_public_type structure that is
730  !! to be deallocated upon termination.
731  type(ocean_state_type), pointer :: Ocean_state !< A pointer to the structure containing
732  !! the internal ocean state to be deallocated
733  !! upon termination.
734  type(time_type), intent(in) :: Time !< The model time, used for writing restarts.
735 
736  call ocean_model_save_restart(ocean_state, time)
737  call diag_mediator_end(time, ocean_state%diag, end_diag_manager=.true.)
738  ! print time stats
739  call mom_infra_end
740  call mom_end(ocean_state%MOM_CSp)
741  if (ocean_state%use_ice_shelf) call ice_shelf_end(ocean_state%Ice_shelf_CSp)

References mom_diag_mediator::diag_mediator_end(), mom_ice_shelf::ice_shelf_end(), and ocean_model_save_restart().

Here is the call graph for this function:

◆ ocean_model_flux_init()

subroutine, public mom_ocean_model_mct::ocean_model_flux_init ( type(ocean_state_type), optional, pointer  OS,
integer, intent(in), optional  verbosity 
)

ocean_model_flux_init is used to initialize properties of the air-sea fluxes as determined by various run-time parameters. It can be called from non-ocean PEs, or PEs that have not yet been initialzed, and it can safely be called multiple times.

Parameters
osAn optional pointer to the ocean state, used to figure out if this is an ocean PE that has already been initialized.
[in]verbosityA 0-9 integer indicating a level of verbosity.

Definition at line 987 of file mom_ocean_model_mct.F90.

987  type(ocean_state_type), optional, pointer :: OS !< An optional pointer to the ocean state,
988  !! used to figure out if this is an ocean PE that
989  !! has already been initialized.
990  integer, optional, intent(in) :: verbosity !< A 0-9 integer indicating a level of verbosity.
991 
992  logical :: OS_is_set
993  integer :: verbose
994 
995  os_is_set = .false. ; if (present(os)) os_is_set = associated(os)
996 
997  ! Use this to control the verbosity of output; consider rethinking this logic later.
998  verbose = 5 ; if (os_is_set) verbose = 3
999  if (present(verbosity)) verbose = verbosity
1000 
1001  call call_tracer_flux_init(verbosity=verbose)
1002 

References mom_tracer_flow_control::call_tracer_flux_init().

Here is the call graph for this function:

◆ ocean_model_init()

subroutine, public mom_ocean_model_mct::ocean_model_init ( type(ocean_public_type), intent(inout), target  Ocean_sfc,
type(ocean_state_type), pointer  OS,
type(time_type), intent(in)  Time_init,
type(time_type), intent(in)  Time_in,
type(coupler_1d_bc_type), intent(in), optional  gas_fields_ocn,
character(len=*), intent(in), optional  input_restart_file 
)

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 initialize_ocean_type.

Parameters
[in,out]ocean_sfcA structure containing various publicly
osA structure whose internal contents are private to ocean_model_mod that may be used to contain all information about the ocean's interior state.
[in]time_initThe start time for the coupled model's calendar
[in]time_inThe time at which to initialize the ocean model.
[in]gas_fields_ocnIf present, this type describes the
[in]input_restart_fileIf present, name of restart file to read

Definition at line 224 of file mom_ocean_model_mct.F90.

224  type(ocean_public_type), target, &
225  intent(inout) :: Ocean_sfc !< A structure containing various publicly
226  !! visible ocean surface properties after initialization,
227  !! the data in this type is intent out.
228  type(ocean_state_type), pointer :: OS !< A structure whose internal
229  !! contents are private to ocean_model_mod that may be used to
230  !! contain all information about the ocean's interior state.
231  type(time_type), intent(in) :: Time_init !< The start time for the coupled model's calendar
232  type(time_type), intent(in) :: Time_in !< The time at which to initialize the ocean model.
233  type(coupler_1d_bc_type), &
234  optional, intent(in) :: gas_fields_ocn !< If present, this type describes the
235  !! ocean and surface-ice fields that will participate
236  !! in the calculation of additional gas or other
237  !! tracer fluxes, and can be used to spawn related
238  !! internal variables in the ice model.
239  character(len=*), optional, intent(in) :: input_restart_file !< If present, name of restart file to read
240  ! Local variables
241  real :: Rho0 ! The Boussinesq ocean density, in kg m-3.
242  real :: G_Earth ! The gravitational acceleration in m s-2.
243  real :: HFrz !< If HFrz > 0 (m), melt potential will be computed.
244  !! The actual depth over which melt potential is computed will
245  !! min(HFrz, OBLD), where OBLD is the boundary layer depth.
246  !! If HFrz <= 0 (default), melt potential will not be computed.
247  logical :: use_melt_pot!< If true, allocate melt_potential array
248 
249 ! This include declares and sets the variable "version".
250 #include "version_variable.h"
251  character(len=40) :: mdl = "ocean_model_init" ! This module's name.
252  character(len=48) :: stagger
253  integer :: secs, days
254  type(param_file_type) :: param_file !< A structure to parse for run-time parameters
255  logical :: use_temperature
256 
257  call calltree_enter("ocean_model_init(), MOM_ocean_model_mct.F90")
258  if (associated(os)) then
259  call mom_error(warning, "ocean_model_init called with an associated "// &
260  "ocean_state_type structure. Model is already initialized.")
261  return
262  endif
263  allocate(os)
264 
265  os%is_ocean_pe = ocean_sfc%is_ocean_pe
266  if (.not.os%is_ocean_pe) return
267 
268  os%Time = time_in
269  call initialize_mom(os%Time, time_init, param_file, os%dirs, os%MOM_CSp, &
270  os%restart_CSp, time_in, offline_tracer_mode=os%offline_tracer_mode, &
271  input_restart_file=input_restart_file, &
272  diag_ptr=os%diag, count_calls=.true.)
273  call get_mom_state_elements(os%MOM_CSp, g=os%grid, gv=os%GV, us=os%US, c_p=os%C_p, &
274  use_temp=use_temperature)
275  os%fluxes%C_p = os%C_p
276 
277  ! Read all relevant parameters and write them to the model log.
278  call log_version(param_file, mdl, version, "")
279 
280  call get_param(param_file, mdl, "SINGLE_STEPPING_CALL", os%single_step_call, &
281  "If true, advance the state of MOM with a single step "//&
282  "including both dynamics and thermodynamics. If false, "//&
283  "the two phases are advanced with separate calls.", default=.true.)
284  call get_param(param_file, mdl, "DT", os%dt, &
285  "The (baroclinic) dynamics time step. The time-step that "//&
286  "is actually used will be an integer fraction of the "//&
287  "forcing time-step.", units="s", fail_if_missing=.true.)
288  call get_param(param_file, mdl, "DT_THERM", os%dt_therm, &
289  "The thermodynamic and tracer advection time step. "//&
290  "Ideally DT_THERM should be an integer multiple of DT "//&
291  "and less than the forcing or coupling time-step, unless "//&
292  "THERMO_SPANS_COUPLING is true, in which case DT_THERM "//&
293  "can be an integer multiple of the coupling timestep. By "//&
294  "default DT_THERM is set to DT.", units="s", default=os%dt)
295  call get_param(param_file, "MOM", "THERMO_SPANS_COUPLING", os%thermo_spans_coupling, &
296  "If true, the MOM will take thermodynamic and tracer "//&
297  "timesteps that can be longer than the coupling timestep. "//&
298  "The actual thermodynamic timestep that is used in this "//&
299  "case is the largest integer multiple of the coupling "//&
300  "timestep that is less than or equal to DT_THERM.", default=.false.)
301  call get_param(param_file, mdl, "DIABATIC_FIRST", os%diabatic_first, &
302  "If true, apply diabatic and thermodynamic processes, "//&
303  "including buoyancy forcing and mass gain or loss, "//&
304  "before stepping the dynamics forward.", default=.false.)
305 
306  call get_param(param_file, mdl, "RESTART_CONTROL", os%Restart_control, &
307  "An integer whose bits encode which restart files are "//&
308  "written. Add 2 (bit 1) for a time-stamped file, and odd "//&
309  "(bit 0) for a non-time-stamped file. A restart file "//&
310  "will be saved at the end of the run segment for any "//&
311  "non-negative value.", default=1)
312  call get_param(param_file, mdl, "OCEAN_SURFACE_STAGGER", stagger, &
313  "A case-insensitive character string to indicate the "//&
314  "staggering of the surface velocity field that is "//&
315  "returned to the coupler. Valid values include "//&
316  "'A', 'B', or 'C'.", default="C")
317  if (uppercase(stagger(1:1)) == 'A') then ; ocean_sfc%stagger = agrid
318  elseif (uppercase(stagger(1:1)) == 'B') then ; ocean_sfc%stagger = bgrid_ne
319  elseif (uppercase(stagger(1:1)) == 'C') then ; ocean_sfc%stagger = cgrid_ne
320  else ; call mom_error(fatal,"ocean_model_init: OCEAN_SURFACE_STAGGER = "// &
321  trim(stagger)//" is invalid.") ; endif
322 
323  call get_param(param_file, mdl, "RESTORE_SALINITY",os%restore_salinity, &
324  "If true, the coupled driver will add a globally-balanced "//&
325  "fresh-water flux that drives sea-surface salinity "//&
326  "toward specified values.", default=.false.)
327  call get_param(param_file, mdl, "RESTORE_TEMPERATURE",os%restore_temp, &
328  "If true, the coupled driver will add a "//&
329  "heat flux that drives sea-surface temperature "//&
330  "toward specified values.", default=.false.)
331  call get_param(param_file, mdl, "RHO_0", rho0, &
332  "The mean ocean density used with BOUSSINESQ true to "//&
333  "calculate accelerations and the mass for conservation "//&
334  "properties, or with BOUSSINSEQ false to convert some "//&
335  "parameters from vertical units of m to kg m-2.", &
336  units="kg m-3", default=1035.0)
337  call get_param(param_file, mdl, "G_EARTH", g_earth, &
338  "The gravitational acceleration of the Earth.", &
339  units="m s-2", default = 9.80)
340 
341  call get_param(param_file, mdl, "ICE_SHELF", os%use_ice_shelf, &
342  "If true, enables the ice shelf model.", default=.false.)
343 
344  call get_param(param_file, mdl, "ICEBERGS_APPLY_RIGID_BOUNDARY", os%icebergs_alter_ocean, &
345  "If true, allows icebergs to change boundary condition felt by ocean", default=.false.)
346 
347  os%press_to_z = 1.0/(rho0*g_earth)
348 
349  call get_param(param_file, mdl, "HFREEZE", hfrz, &
350  "If HFREEZE > 0, melt potential will be computed. The actual depth "//&
351  "over which melt potential is computed will be min(HFREEZE, OBLD), "//&
352  "where OBLD is the boundary layer depth. If HFREEZE <= 0 (default), "//&
353  "melt potential will not be computed.", units="m", default=-1.0, do_not_log=.true.)
354 
355  if (hfrz .gt. 0.0) then
356  use_melt_pot=.true.
357  else
358  use_melt_pot=.false.
359  endif
360 
361  ! Consider using a run-time flag to determine whether to do the diagnostic
362  ! vertical integrals, since the related 3-d sums are not negligible in cost.
363  call allocate_surface_state(os%sfc_state, os%grid, use_temperature, &
364  do_integrals=.true., gas_fields_ocn=gas_fields_ocn, use_meltpot=use_melt_pot)
365 
366  call surface_forcing_init(time_in, os%grid, os%US, param_file, os%diag, &
367  os%forcing_CSp, os%restore_salinity, os%restore_temp)
368 
369  if (os%use_ice_shelf) then
370  call initialize_ice_shelf(param_file, os%grid, os%Time, os%ice_shelf_CSp, &
371  os%diag, os%forces, os%fluxes)
372  endif
373 
374  if (os%icebergs_alter_ocean) then
375  call marine_ice_init(os%Time, os%grid, param_file, os%diag, os%marine_ice_CSp)
376  if (.not. os%use_ice_shelf) &
377  call allocate_forcing_type(os%grid, os%fluxes, shelf=.true.)
378  endif
379 
380  call get_param(param_file, mdl, "USE_WAVES", os%Use_Waves, &
381  "If true, enables surface wave modules.", default=.false.)
382  if (os%use_waves) then
383  call mom_wave_interface_init(os%Time, os%grid, os%GV, os%US, param_file, os%Waves, os%diag)
384  else
385  call mom_wave_interface_init_lite(param_file)
386  endif
387 
388  if (associated(os%grid%Domain%maskmap)) then
389  call initialize_ocean_public_type(os%grid%Domain%mpp_domain, ocean_sfc, &
390  os%diag, maskmap=os%grid%Domain%maskmap, &
391  gas_fields_ocn=gas_fields_ocn)
392  else
393  call initialize_ocean_public_type(os%grid%Domain%mpp_domain, ocean_sfc, &
394  os%diag, gas_fields_ocn=gas_fields_ocn)
395  endif
396 
397  ! This call can only occur here if the coupler_bc_type variables have been
398  ! initialized already using the information from gas_fields_ocn.
399  if (present(gas_fields_ocn)) then
400  call coupler_type_set_diags(ocean_sfc%fields, "ocean_sfc", &
401  ocean_sfc%axes(1:2), time_in)
402 
403  call extract_surface_state(os%MOM_CSp, os%sfc_state)
404 
405  call convert_state_to_ocean_type(os%sfc_state, ocean_sfc, os%grid, os%US)
406  endif
407 
408  call close_param_file(param_file)
409  call diag_mediator_close_registration(os%diag)
410 
411  if (is_root_pe()) &
412  write(*,'(/12x,a/)') '======== COMPLETED MOM INITIALIZATION ========'
413 
414  call calltree_leave("ocean_model_init(")

References mom_error_handler::calltree_enter(), mom_error_handler::calltree_leave(), mom_file_parser::close_param_file(), convert_state_to_ocean_type(), mom_diag_mediator::diag_mediator_close_registration(), initialize_ocean_public_type(), mom_marine_ice::marine_ice_init(), mom_wave_interface::mom_wave_interface_init_lite(), and mom_string_functions::uppercase().

Here is the call graph for this function:

◆ ocean_model_init_sfc()

subroutine, public mom_ocean_model_mct::ocean_model_init_sfc ( type(ocean_state_type), pointer  OS,
type(ocean_public_type), intent(inout)  Ocean_sfc 
)

This subroutine extracts the surface properties from the ocean's internal state and stores them in the ocean type returned to the calling ice model. It has to be separate from the ocean_initialization call because the coupler module allocates the space for some of these variables.

Parameters
osThe structure with the complete ocean state
[in,out]ocean_sfcA structure containing various publicly visible ocean surface properties after initialization, whose elements have their data set here.

Definition at line 966 of file mom_ocean_model_mct.F90.

966  type(ocean_state_type), pointer :: OS !< The structure with the complete ocean state
967  type(ocean_public_type), intent(inout) :: Ocean_sfc !< A structure containing various publicly
968  !! visible ocean surface properties after initialization, whose
969  !! elements have their data set here.
970  integer :: is, ie, js, je
971 
972  is = os%grid%isc ; ie = os%grid%iec ; js = os%grid%jsc ; je = os%grid%jec
973  call coupler_type_spawn(ocean_sfc%fields, os%sfc_state%tr_fields, &
974  (/is,is,ie,ie/), (/js,js,je,je/), as_needed=.true.)
975 
976  call extract_surface_state(os%MOM_CSp, os%sfc_state)
977 
978  call convert_state_to_ocean_type(os%sfc_state, ocean_sfc, os%grid, os%US)
979 

References convert_state_to_ocean_type().

Here is the call graph for this function:

◆ ocean_model_restart()

subroutine, public mom_ocean_model_mct::ocean_model_restart ( type(ocean_state_type), pointer  OS,
character(len=*), intent(in), optional  timestamp,
character(len=*), intent(in), optional  restartname 
)

This subroutine writes out the ocean model restart file.

Parameters
osA pointer to the structure containing the internal ocean state being saved to a restart file
[in]timestampAn optional timestamp string that should be prepended to the file name. (Currently this is unused.)
[in]restartnameName of restart file to use This option distinguishes the cesm interface from the non-cesm interface

Definition at line 677 of file mom_ocean_model_mct.F90.

677  type(ocean_state_type), pointer :: OS !< A pointer to the structure containing the
678  !! internal ocean state being saved to a restart file
679  character(len=*), optional, intent(in) :: timestamp !< An optional timestamp string that should be
680  !! prepended to the file name. (Currently this is unused.)
681  character(len=*), optional, intent(in) :: restartname !< Name of restart file to use
682  !! This option distinguishes the cesm interface from the
683  !! non-cesm interface
684 
685  if (.not.mom_state_is_synchronized(os%MOM_CSp)) &
686  call mom_error(warning, "End of MOM_main reached with inconsistent "//&
687  "dynamics and advective times. Additional restart fields "//&
688  "that have not been coded yet would be required for reproducibility.")
689  if (.not.os%fluxes%fluxes_used) call mom_error(fatal, "ocean_model_restart "//&
690  "was called with unused buoyancy fluxes. For conservation, the ocean "//&
691  "restart files can only be created after the buoyancy forcing is applied.")
692 
693  if (present(restartname)) then
694  call save_restart(os%dirs%restart_output_dir, os%Time, os%grid, &
695  os%restart_CSp, gv=os%GV, filename=restartname)
696  call forcing_save_restart(os%forcing_CSp, os%grid, os%Time, &
697  os%dirs%restart_output_dir) ! Is this needed?
698  if (os%use_ice_shelf) then
699  call ice_shelf_save_restart(os%Ice_shelf_CSp, os%Time, &
700  os%dirs%restart_output_dir)
701  endif
702  else
703  if (btest(os%Restart_control,1)) then
704  call save_restart(os%dirs%restart_output_dir, os%Time, os%grid, &
705  os%restart_CSp, .true., gv=os%GV)
706  call forcing_save_restart(os%forcing_CSp, os%grid, os%Time, &
707  os%dirs%restart_output_dir, .true.)
708  if (os%use_ice_shelf) then
709  call ice_shelf_save_restart(os%Ice_shelf_CSp, os%Time, os%dirs%restart_output_dir, .true.)
710  endif
711  endif
712  if (btest(os%Restart_control,0)) then
713  call save_restart(os%dirs%restart_output_dir, os%Time, os%grid, &
714  os%restart_CSp, gv=os%GV)
715  call forcing_save_restart(os%forcing_CSp, os%grid, os%Time, &
716  os%dirs%restart_output_dir)
717  if (os%use_ice_shelf) then
718  call ice_shelf_save_restart(os%Ice_shelf_CSp, os%Time, os%dirs%restart_output_dir)
719  endif
720  endif
721  end if
722 

References mom_surface_forcing_mct::forcing_save_restart(), mom_ice_shelf::ice_shelf_save_restart(), and mom_restart::save_restart().

Here is the call graph for this function:

◆ ocean_model_save_restart()

subroutine, public mom_ocean_model_mct::ocean_model_save_restart ( type(ocean_state_type), pointer  OS,
type(time_type), intent(in)  Time,
character(len=*), intent(in), optional  directory,
character(len=*), intent(in), optional  filename_suffix 
)

ocean_model_save_restart causes restart files associated with the ocean to be written out.

Parameters
osA pointer to the structure containing the internal ocean state (in).
[in]timeThe model time at this call, needed for mpp_write calls.
[in]directoryAn optional directory into which to write these restart files.
[in]filename_suffixAn optional suffix (e.g., a time-stamp) to append to the restart file names.

Definition at line 747 of file mom_ocean_model_mct.F90.

747  type(ocean_state_type), pointer :: OS !< A pointer to the structure containing the
748  !! internal ocean state (in).
749  type(time_type), intent(in) :: Time !< The model time at this call, needed for mpp_write calls.
750  character(len=*), optional, intent(in) :: directory !< An optional directory into which to
751  !! write these restart files.
752  character(len=*), optional, intent(in) :: filename_suffix !< An optional suffix (e.g., a time-stamp)
753  !! to append to the restart file names.
754 ! Note: This is a new routine - it will need to exist for the new incremental
755 ! checkpointing. It will also be called by ocean_model_end, giving the same
756 ! restart behavior as now in FMS.
757  character(len=200) :: restart_dir
758 
759  if (.not.mom_state_is_synchronized(os%MOM_CSp)) &
760  call mom_error(warning, "ocean_model_save_restart called with inconsistent "//&
761  "dynamics and advective times. Additional restart fields "//&
762  "that have not been coded yet would be required for reproducibility.")
763  if (.not.os%fluxes%fluxes_used) call mom_error(fatal, "ocean_model_save_restart "//&
764  "was called with unused buoyancy fluxes. For conservation, the ocean "//&
765  "restart files can only be created after the buoyancy forcing is applied.")
766 
767  if (present(directory)) then ; restart_dir = directory
768  else ; restart_dir = os%dirs%restart_output_dir ; endif
769 
770  call save_restart(restart_dir, time, os%grid, os%restart_CSp, gv=os%GV)
771 
772  call forcing_save_restart(os%forcing_CSp, os%grid, time, restart_dir)
773 
774  if (os%use_ice_shelf) then
775  call ice_shelf_save_restart(os%Ice_shelf_CSp, os%Time, os%dirs%restart_output_dir)
776  endif
777 

References mom_surface_forcing_mct::forcing_save_restart(), mom_ice_shelf::ice_shelf_save_restart(), and mom_restart::save_restart().

Referenced by ocean_model_end().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ ocean_public_type_chksum()

subroutine, public mom_ocean_model_mct::ocean_public_type_chksum ( character(len=*), intent(in)  id,
integer, intent(in)  timestep,
type(ocean_public_type), intent(in)  ocn 
)

Write out FMS-format checsums on fields from the ocean surface state.

Parameters
[in]idAn identifying string for this call
[in]timestepThe number of elapsed timesteps
[in]ocnA structure containing various publicly visible ocean surface fields.

Definition at line 1049 of file mom_ocean_model_mct.F90.

1049 
1050  character(len=*), intent(in) :: id !< An identifying string for this call
1051  integer, intent(in) :: timestep !< The number of elapsed timesteps
1052  type(ocean_public_type), intent(in) :: ocn !< A structure containing various publicly
1053  !! visible ocean surface fields.
1054  integer :: n, m, outunit
1055 
1056  outunit = stdout()
1057 
1058  write(outunit,*) "BEGIN CHECKSUM(ocean_type):: ", id, timestep
1059  write(outunit,100) 'ocean%t_surf ',mpp_chksum(ocn%t_surf )
1060  write(outunit,100) 'ocean%s_surf ',mpp_chksum(ocn%s_surf )
1061  write(outunit,100) 'ocean%u_surf ',mpp_chksum(ocn%u_surf )
1062  write(outunit,100) 'ocean%v_surf ',mpp_chksum(ocn%v_surf )
1063  write(outunit,100) 'ocean%sea_lev ',mpp_chksum(ocn%sea_lev)
1064  write(outunit,100) 'ocean%frazil ',mpp_chksum(ocn%frazil )
1065  write(outunit,100) 'ocean%melt_potential ',mpp_chksum(ocn%melt_potential)
1066 
1067  call coupler_type_write_chksums(ocn%fields, outunit, 'ocean%')
1068 100 FORMAT(" CHECKSUM::",a20," = ",z20)
1069 

◆ ocean_stock_pe()

subroutine, public mom_ocean_model_mct::ocean_stock_pe ( type(ocean_state_type), pointer  OS,
integer, intent(in)  index,
real, intent(out)  value,
integer, intent(in), optional  time_index 
)

Ocean_stock_pe - returns the integrated stocks of heat, water, etc. for conservation checks. Because of the way FMS is coded, only the root PE has the integrated amount, while all other PEs get 0.

Parameters
osA structure containing the internal ocean state. The data in OS is intent in.
[in]indexThe stock index for the quantity of interest.
[out]valueSum returned for the conservation quantity of interest.
[in]time_indexAn unused optional argument, present only for interfacial compatibility with other models.

Definition at line 1009 of file mom_ocean_model_mct.F90.

1009  use stock_constants_mod, only : istock_water, istock_heat,istock_salt
1010  type(ocean_state_type), pointer :: OS !< A structure containing the internal ocean state.
1011  !! The data in OS is intent in.
1012  integer, intent(in) :: index !< The stock index for the quantity of interest.
1013  real, intent(out) :: value !< Sum returned for the conservation quantity of interest.
1014  integer, optional, intent(in) :: time_index !< An unused optional argument, present only for
1015  !! interfacial compatibility with other models.
1016 ! Arguments: OS - A structure containing the internal ocean state.
1017 ! (in) index - Index of conservation quantity of interest.
1018 ! (in) value - Sum returned for the conservation quantity of interest.
1019 ! (in,opt) time_index - Index for time level to use if this is necessary.
1020 
1021  real :: salt
1022 
1023  value = 0.0
1024  if (.not.associated(os)) return
1025  if (.not.os%is_ocean_pe) return
1026 
1027  select case (index)
1028  case (istock_water) ! Return the mass of fresh water in the ocean in kg.
1029  if (os%GV%Boussinesq) then
1030  call get_ocean_stocks(os%MOM_CSp, mass=value, on_pe_only=.true.)
1031  else ! In non-Boussinesq mode, the mass of salt needs to be subtracted.
1032  call get_ocean_stocks(os%MOM_CSp, mass=value, salt=salt, on_pe_only=.true.)
1033  value = value - salt
1034  endif
1035  case (istock_heat) ! Return the heat content of the ocean in J.
1036  call get_ocean_stocks(os%MOM_CSp, heat=value, on_pe_only=.true.)
1037  case (istock_salt) ! Return the mass of the salt in the ocean in kg.
1038  call get_ocean_stocks(os%MOM_CSp, salt=value, on_pe_only=.true.)
1039  case default ; value = 0.0
1040  end select
1041  ! If the FMS coupler is changed so that Ocean_stock_PE is only called on
1042  ! ocean PEs, uncomment the following and eliminate the on_PE_only flags above.
1043  ! if (.not.is_root_pe()) value = 0.0
1044 

References mom::get_ocean_stocks().

Here is the call graph for this function:

◆ update_ocean_model()

subroutine, public mom_ocean_model_mct::update_ocean_model ( type(ice_ocean_boundary_type), intent(in)  Ice_ocean_boundary,
type(ocean_state_type), pointer  OS,
type(ocean_public_type), intent(inout)  Ocean_sfc,
type(time_type), intent(in)  time_start_update,
type(time_type), intent(in)  Ocean_coupling_time_step,
logical, intent(in), optional  update_dyn,
logical, intent(in), optional  update_thermo,
logical, intent(in), optional  Ocn_fluxes_used 
)

update_ocean_model uses the forcing in Ice_ocean_boundary to advance the ocean model's state from the input value of Ocean_state (which must be for time time_start_update) for a time interval of Ocean_coupling_time_step, returning the publicly visible ocean surface properties in Ocean_sfc and storing the new ocean properties in Ocean_state.

Parameters
[in]ice_ocean_boundaryA structure containing the
osA pointer to a private structure containing
[in,out]ocean_sfcA structure containing all the
[in]time_start_updateThe time at the beginning of the update step.
[in]ocean_coupling_time_stepThe amount of time over which to advance the ocean.
[in]update_dynIf present and false, do not do updates due to the ocean dynamics.
[in]update_thermoIf present and false, do not do updates due to the ocean thermodynamics or remapping.
[in]ocn_fluxes_usedIf present, this indicates whether the cumulative thermodynamic fluxes from the ocean, like frazil, have been used and should be reset.

Definition at line 425 of file mom_ocean_model_mct.F90.

425  type(ice_ocean_boundary_type), &
426  intent(in) :: Ice_ocean_boundary !< A structure containing the
427  !! various forcing fields coming from the ice.
428  type(ocean_state_type), &
429  pointer :: OS !< A pointer to a private structure containing
430  !! the internal ocean state.
431  type(ocean_public_type), &
432  intent(inout) :: Ocean_sfc !< A structure containing all the
433  !! publicly visible ocean surface fields after
434  !! a coupling time step. The data in this type is
435  !! intent out.
436  type(time_type), intent(in) :: time_start_update !< The time at the beginning of the update step.
437  type(time_type), intent(in) :: Ocean_coupling_time_step !< The amount of time over
438  !! which to advance the ocean.
439  logical, optional, intent(in) :: update_dyn !< If present and false, do not do updates
440  !! due to the ocean dynamics.
441  logical, optional, intent(in) :: update_thermo !< If present and false, do not do updates
442  !! due to the ocean thermodynamics or remapping.
443  logical, optional, intent(in) :: Ocn_fluxes_used !< If present, this indicates whether the
444  !! cumulative thermodynamic fluxes from the ocean,
445  !! like frazil, have been used and should be reset.
446  ! Local variables
447  type(time_type) :: Master_time ! This allows step_MOM to temporarily change
448  ! the time that is seen by internal modules.
449  type(time_type) :: Time1 ! The value of the ocean model's time at the
450  ! start of a call to step_MOM.
451  integer :: index_bnds(4) ! The computational domain index bounds in the
452  ! ice-ocean boundary type.
453  real :: weight ! Flux accumulation weight
454  real :: dt_coupling ! The coupling time step in seconds.
455  integer :: nts ! The number of baroclinic dynamics time steps
456  ! within dt_coupling.
457  real :: dt_therm ! A limited and quantized version of OS%dt_therm (sec)
458  real :: dt_dyn ! The dynamics time step in sec.
459  real :: dtdia ! The diabatic time step in sec.
460  real :: t_elapsed_seg ! The elapsed time in this update segment, in s.
461  integer :: n, n_max, n_last_thermo
462  type(time_type) :: Time2 ! A temporary time.
463  logical :: thermo_does_span_coupling ! If true, thermodynamic forcing spans
464  ! multiple dynamic timesteps.
465  logical :: do_dyn ! If true, step the ocean dynamics and transport.
466  logical :: do_thermo ! If true, step the ocean thermodynamics.
467  logical :: step_thermo ! If true, take a thermodynamic step.
468  integer :: secs, days
469  integer :: is, ie, js, je
470 
471  call calltree_enter("update_ocean_model(), MOM_ocean_model_mct.F90")
472  call get_time(ocean_coupling_time_step, secs, days)
473  dt_coupling = 86400.0*real(days) + real(secs)
474 
475  if (time_start_update /= os%Time) then
476  call mom_error(warning, "update_ocean_model: internal clock does not "//&
477  "agree with time_start_update argument.")
478  endif
479 
480  if (.not.associated(os)) then
481  call mom_error(fatal, "update_ocean_model called with an unassociated "// &
482  "ocean_state_type structure. ocean_model_init must be "// &
483  "called first to allocate this structure.")
484  return
485  endif
486 
487  do_dyn = .true. ; if (present(update_dyn)) do_dyn = update_dyn
488  do_thermo = .true. ; if (present(update_thermo)) do_thermo = update_thermo
489 
490  ! This is benign but not necessary if ocean_model_init_sfc was called or if
491  ! OS%sfc_state%tr_fields was spawned in ocean_model_init. Consider removing it.
492  is = os%grid%isc ; ie = os%grid%iec ; js = os%grid%jsc ; je = os%grid%jec
493  call coupler_type_spawn(ocean_sfc%fields, os%sfc_state%tr_fields, &
494  (/is,is,ie,ie/), (/js,js,je,je/), as_needed=.true.)
495 
496  ! Translate Ice_ocean_boundary into fluxes.
497  call mpp_get_compute_domain(ocean_sfc%Domain, index_bnds(1), index_bnds(2), &
498  index_bnds(3), index_bnds(4))
499  weight = 1.0
500 
501  call convert_iob_to_forces(ice_ocean_boundary, os%forces, index_bnds, os%Time, &
502  os%grid, os%US, os%forcing_CSp)
503 
504  if (os%fluxes%fluxes_used) then
505 
506  ! GMM, is enable_averaging needed now?
507  call enable_averaging(dt_coupling, os%Time + ocean_coupling_time_step, os%diag)
508 
509  if (do_thermo) &
510  call convert_iob_to_fluxes(ice_ocean_boundary, os%fluxes, index_bnds, os%Time, dt_coupling, &
511  os%grid, os%US, os%forcing_CSp, os%sfc_state, &
512  os%restore_salinity, os%restore_temp)
513 
514  ! Add ice shelf fluxes
515  if (os%use_ice_shelf) then
516  if (do_thermo) &
517  call shelf_calc_flux(os%sfc_state, os%fluxes, os%Time, dt_coupling, os%Ice_shelf_CSp)
518  if (do_dyn) &
519  call add_shelf_forces(os%grid, os%US, os%Ice_shelf_CSp, os%forces)
520  endif
521  if (os%icebergs_alter_ocean) then
522  if (do_dyn) &
523  call iceberg_forces(os%grid, os%forces, os%use_ice_shelf, &
524  os%sfc_state, dt_coupling, os%marine_ice_CSp)
525  if (do_thermo) &
526  call iceberg_fluxes(os%grid, os%US, os%fluxes, os%use_ice_shelf, &
527  os%sfc_state, dt_coupling, os%marine_ice_CSp)
528  endif
529 
530  ! Fields that exist in both the forcing and mech_forcing types must be copied.
531  call copy_common_forcing_fields(os%forces, os%fluxes, os%grid)
532 
533 #ifdef _USE_GENERIC_TRACER
534  call enable_averaging(dt_coupling, os%Time + ocean_coupling_time_step, os%diag) !Is this needed?
535  call mom_generic_tracer_fluxes_accumulate(os%fluxes, weight) !here weight=1, just saving the current fluxes
536 #endif
537 
538  else
539 
540  os%flux_tmp%C_p = os%fluxes%C_p
541 
542  if (do_thermo) &
543  call convert_iob_to_fluxes(ice_ocean_boundary, os%flux_tmp, index_bnds, os%Time, dt_coupling, &
544  os%grid, os%US, os%forcing_CSp, os%sfc_state, os%restore_salinity,os%restore_temp)
545 
546  if (os%use_ice_shelf) then
547  if (do_thermo) &
548  call shelf_calc_flux(os%sfc_state, os%flux_tmp, os%Time, dt_coupling, os%Ice_shelf_CSp)
549  if (do_dyn) &
550  call add_shelf_forces(os%grid, os%US, os%Ice_shelf_CSp, os%forces)
551  endif
552  if (os%icebergs_alter_ocean) then
553  if (do_dyn) &
554  call iceberg_forces(os%grid, os%forces, os%use_ice_shelf, &
555  os%sfc_state, dt_coupling, os%marine_ice_CSp)
556  if (do_thermo) &
557  call iceberg_fluxes(os%grid, os%US, os%flux_tmp, os%use_ice_shelf, &
558  os%sfc_state, dt_coupling, os%marine_ice_CSp)
559  endif
560 
561  call forcing_accumulate(os%flux_tmp, os%forces, os%fluxes, os%grid, weight)
562 
563  ! Some of the fields that exist in both the forcing and mech_forcing types
564  ! (e.g., ustar) are time-averages must be copied back to the forces type.
565  call copy_back_forcing_fields(os%fluxes, os%forces, os%grid)
566 
567 #ifdef _USE_GENERIC_TRACER
568  call mom_generic_tracer_fluxes_accumulate(os%flux_tmp, weight) !weight of the current flux in the running average
569 #endif
570  endif
571 
572  call set_derived_forcing_fields(os%forces, os%fluxes, os%grid, os%US, os%GV%Rho0)
573  call set_net_mass_forcing(os%fluxes, os%forces, os%grid, os%US)
574 
575  if (os%use_waves) then
576  call update_surface_waves(os%grid, os%GV, os%US, os%time, ocean_coupling_time_step, os%waves)
577  endif
578 
579  if (os%nstep==0) then
580  call finish_mom_initialization(os%Time, os%dirs, os%MOM_CSp, os%restart_CSp)
581  endif
582 
583  call disable_averaging(os%diag)
584  master_time = os%Time ; time1 = os%Time
585 
586  if(os%offline_tracer_mode) then
587  call step_offline(os%forces, os%fluxes, os%sfc_state, time1, dt_coupling, os%MOM_CSp)
588 
589  elseif ((.not.do_thermo) .or. (.not.do_dyn)) then
590  ! The call sequence is being orchestrated from outside of update_ocean_model.
591  call step_mom(os%forces, os%fluxes, os%sfc_state, time1, dt_coupling, os%MOM_CSp, &
592  waves=os%Waves, do_dynamics=do_thermo, do_thermodynamics=do_dyn, &
593  reset_therm=ocn_fluxes_used)
594 
595  elseif (os%single_step_call) then
596  call step_mom(os%forces, os%fluxes, os%sfc_state, time1, dt_coupling, os%MOM_CSp, waves=os%Waves)
597 
598  else
599  n_max = 1 ; if (dt_coupling > os%dt) n_max = ceiling(dt_coupling/os%dt - 0.001)
600  dt_dyn = dt_coupling / real(n_max)
601  thermo_does_span_coupling = (os%thermo_spans_coupling .and. &
602  (os%dt_therm > 1.5*dt_coupling))
603 
604  if (thermo_does_span_coupling) then
605  dt_therm = dt_coupling * floor(os%dt_therm / dt_coupling + 0.001)
606  nts = floor(dt_therm/dt_dyn + 0.001)
607  else
608  nts = max(1,min(n_max,floor(os%dt_therm/dt_dyn + 0.001)))
609  n_last_thermo = 0
610  endif
611 
612  time2 = time1 ; t_elapsed_seg = 0.0
613  do n=1,n_max
614  if (os%diabatic_first) then
615  if (thermo_does_span_coupling) call mom_error(fatal, &
616  "MOM is not yet set up to have restarts that work with "//&
617  "THERMO_SPANS_COUPLING and DIABATIC_FIRST.")
618  if (modulo(n-1,nts)==0) then
619  dtdia = dt_dyn*min(nts,n_max-(n-1))
620  call step_mom(os%forces, os%fluxes, os%sfc_state, time2, dtdia, os%MOM_CSp, &
621  waves=os%Waves, do_dynamics=.false., do_thermodynamics=.true., &
622  start_cycle=(n==1), end_cycle=.false., cycle_length=dt_coupling)
623  endif
624 
625  call step_mom(os%forces, os%fluxes, os%sfc_state, time2, dt_dyn, os%MOM_CSp, &
626  waves=os%Waves, do_dynamics=.true., do_thermodynamics=.false., &
627  start_cycle=.false., end_cycle=(n==n_max), cycle_length=dt_coupling)
628  else
629  call step_mom(os%forces, os%fluxes, os%sfc_state, time2, dt_dyn, os%MOM_CSp, &
630  waves=os%Waves, do_dynamics=.true., do_thermodynamics=.false., &
631  start_cycle=(n==1), end_cycle=.false., cycle_length=dt_coupling)
632 
633  step_thermo = .false.
634  if (thermo_does_span_coupling) then
635  dtdia = dt_therm
636  step_thermo = mom_state_is_synchronized(os%MOM_CSp, adv_dyn=.true.)
637  elseif ((modulo(n,nts)==0) .or. (n==n_max)) then
638  dtdia = dt_dyn*(n - n_last_thermo)
639  n_last_thermo = n
640  step_thermo = .true.
641  endif
642 
643  if (step_thermo) then
644  ! Back up Time2 to the start of the thermodynamic segment.
645  time2 = time2 - set_time(int(floor((dtdia - dt_dyn) + 0.5)))
646  call step_mom(os%forces, os%fluxes, os%sfc_state, time2, dtdia, os%MOM_CSp, &
647  waves=os%Waves, do_dynamics=.false., do_thermodynamics=.true., &
648  start_cycle=.false., end_cycle=(n==n_max), cycle_length=dt_coupling)
649  endif
650  endif
651 
652  t_elapsed_seg = t_elapsed_seg + dt_dyn
653  time2 = time1 + set_time(int(floor(t_elapsed_seg + 0.5)))
654  enddo
655  endif
656 
657  os%Time = master_time + ocean_coupling_time_step
658  os%nstep = os%nstep + 1
659 
660  call mech_forcing_diags(os%forces, dt_coupling, os%grid, os%Time, os%diag, os%forcing_CSp%handles)
661 
662  if (os%fluxes%fluxes_used) then
663  call forcing_diagnostics(os%fluxes, os%sfc_state, os%grid, os%US, os%Time, os%diag, os%forcing_CSp%handles)
664  endif
665 
666 ! Translate state into Ocean.
667 ! call convert_state_to_ocean_type(OS%sfc_state, Ocean_sfc, OS%grid, &
668 ! Ice_ocean_boundary%p, OS%press_to_z)
669  call convert_state_to_ocean_type(os%sfc_state, ocean_sfc, os%grid, os%US)
670  call coupler_type_send_data(ocean_sfc%fields, os%Time)
671 
672  call calltree_leave("update_ocean_model()")

References mom_ice_shelf::add_shelf_forces(), mom_error_handler::calltree_enter(), mom_error_handler::calltree_leave(), convert_state_to_ocean_type(), mom_forcing_type::forcing_diagnostics(), mom_marine_ice::iceberg_fluxes(), mom_marine_ice::iceberg_forces(), mom_forcing_type::mech_forcing_diags(), mom::step_offline(), and mom_wave_interface::update_surface_waves().

Here is the call graph for this function: