Skip to content

Commit

Permalink
Use dt_trans in post_transport_diagnotics
Browse files Browse the repository at this point in the history
  Actually use the time-increment argument that is passed into
post_transport_diagnotics, instead of pulling the time interval out of the
control structure; for now they are the same, but they might not always be. Also
added clarifying comments and improved variable names to MOM.F90.  All answers
are bitwise identical.
  • Loading branch information
Hallberg-NOAA committed Jun 8, 2017
1 parent b8881f6 commit 4d7e1d9
Showing 1 changed file with 40 additions and 19 deletions.
59 changes: 40 additions & 19 deletions src/core/MOM.F90
Original file line number Diff line number Diff line change
Expand Up @@ -220,9 +220,11 @@ module MOM
real :: dt_therm !< thermodynamics time step (seconds)
logical :: thermo_spans_coupling !< If true, thermodynamic and tracer time
!! steps can span multiple coupled time steps.
real :: t_dyn_rel_adv !< The elapsed time since updating the tracers and
!! applying diabatic processes (sec); may
!! span multiple timesteps.
real :: t_dyn_rel_adv !< The time of the dynamics relative to tracer
!! advection and lateral mixing (in seconds), or
!! equivalently the elapsed time since advectively
!! updating the tracers. t_dyn_rel_adv is invariably
!! positive and may span multiple coupling timesteps.
type(time_type) :: Z_diag_interval !< amount of time between calculating Z-space diagnostics
type(time_type) :: Z_diag_time !< next time to compute Z-space diagnostics
type(time_type), pointer :: Time !< pointer to ocean clock
Expand Down Expand Up @@ -457,7 +459,7 @@ subroutine step_MOM(fluxes, state, Time_start, time_interval, CS)
type(ocean_grid_type), pointer :: G ! pointer to a structure containing
! metrics and related information
type(verticalGrid_type), pointer :: GV => NULL()
integer, save :: nt = 1 ! running number of iterations
integer, save :: nt_debug = 1 ! running number of iterations, for debugging only.
integer :: ntstep ! time steps between tracer updates or diabatic forcing
integer :: n_max ! number of steps to take in this call

Expand All @@ -470,6 +472,13 @@ subroutine step_MOM(fluxes, state, Time_start, time_interval, CS)
real :: dt_therm ! a limited and quantized version of CS%dt_therm (sec)
real :: dtbt_reset_time ! value of CS%rel_time when DTBT was last calculated (sec)

real :: mass_src_time ! The amount of time for the surface mass source from
! precipitation-evaporation, rivers, etc., that should
! be applied to the start of the barotropic solver to
! avoid generating tsunamis, in s. This is negative
! if the precipation has already been applied to the
! layers, and positive if it will be applied later.

real :: wt_end, wt_beg

logical :: calc_dtbt ! Indicates whether the dynamically adjusted
Expand Down Expand Up @@ -649,7 +658,7 @@ subroutine step_MOM(fluxes, state, Time_start, time_interval, CS)

do n=1,n_max

nt = nt + 1
nt_debug = nt_debug + 1
call cpu_clock_begin(id_clock_other)

! Set the universally visible time to the middle of the time step
Expand All @@ -661,6 +670,8 @@ subroutine step_MOM(fluxes, state, Time_start, time_interval, CS)
if (showCallTree) call callTree_enter("DT cycles (step_MOM) n=",n)
call cpu_clock_end(id_clock_other)

!===========================================================================
! This is the first place where the diabatic processes and remapping could occur.
if (CS%diabatic_first .and. (CS%t_dyn_rel_adv==0.0)) then ! do thermodynamics.

if (thermo_does_span_coupling) then
Expand Down Expand Up @@ -692,6 +703,10 @@ subroutine step_MOM(fluxes, state, Time_start, time_interval, CS)
if (associated(CS%tv%S)) call hchksum(CS%tv%S, "Pre set_viscous_BBL S", G%HI, haloshift=1)
endif

! Calculate the BBL properties and store them inside visc (u,h).
! This is here so that CS%visc is updated before diabatic() when
! DIABATIC_FIRST=True. Otherwise diabatic() is called after the dynamics
! and set_viscous_BBL is called as a part of the dynamic stepping.
!call cpu_clock_begin(id_clock_vertvisc)
call set_viscous_BBL(u, v, h, CS%tv, CS%visc, G, GV, CS%set_visc_CSp)
!call cpu_clock_end(id_clock_vertvisc)
Expand All @@ -718,15 +733,14 @@ subroutine step_MOM(fluxes, state, Time_start, time_interval, CS)

endif ! end of block "(CS%diabatic_first .and. (CS%t_dyn_rel_adv==0.0))"

call cpu_clock_begin(id_clock_other)

call cpu_clock_begin(id_clock_pass)
call cpu_clock_begin(id_clock_other) ; call cpu_clock_begin(id_clock_pass)
if (do_pass_kd_kv_turb) call do_group_pass(CS%pass_kd_kv_turb, G%Domain)
call cpu_clock_end(id_clock_pass)

call cpu_clock_end(id_clock_other)
call cpu_clock_end(id_clock_pass) ; call cpu_clock_end(id_clock_other)


!===========================================================================
! This is the start of the dynamics stepping part of the algorithm.
call cpu_clock_begin(id_clock_dynamics)
call disable_averaging(CS%diag)

Expand Down Expand Up @@ -801,15 +815,18 @@ subroutine step_MOM(fluxes, state, Time_start, time_interval, CS)
dtbt_reset_time = CS%rel_time
endif

mass_src_time = CS%t_dyn_rel_adv
!### This should be mass_src_time = CS%t_dyn_rel_dia

if (CS%legacy_split) then
call step_MOM_dyn_legacy_split(u, v, h, CS%tv, CS%visc, &
Time_local, dt, fluxes, CS%p_surf_begin, CS%p_surf_end, &
CS%t_dyn_rel_adv, dt_therm, CS%uh, CS%vh, CS%uhtr, CS%vhtr, &
mass_src_time, dt_therm, CS%uh, CS%vh, CS%uhtr, CS%vhtr, &
eta_av, G, GV, CS%dyn_legacy_split_CSp, calc_dtbt, CS%VarMix, CS%MEKE)
else
call step_MOM_dyn_split_RK2(u, v, h, CS%tv, CS%visc, &
Time_local, dt, fluxes, CS%p_surf_begin, CS%p_surf_end, &
CS%t_dyn_rel_adv, dt_therm, CS%uh, CS%vh, CS%uhtr, CS%vhtr, &
mass_src_time, dt_therm, CS%uh, CS%vh, CS%uhtr, CS%vhtr, &
eta_av, G, GV, CS%dyn_split_RK2_CSp, calc_dtbt, CS%VarMix, CS%MEKE)
endif
if (showCallTree) call callTree_waypoint("finished step_MOM_dyn_split (step_MOM)")
Expand Down Expand Up @@ -884,17 +901,21 @@ subroutine step_MOM(fluxes, state, Time_start, time_interval, CS)
if (CS%useMEKE) call step_forward_MEKE(CS%MEKE, h, CS%VarMix%SN_u, CS%VarMix%SN_v, &
CS%visc, dt, G, GV, CS%MEKE_CSp, CS%uhtr, CS%vhtr)
call disable_averaging(CS%diag)
call cpu_clock_end(id_clock_dynamics)

! Advance the dynamics time by dt.
CS%t_dyn_rel_adv = CS%t_dyn_rel_adv + dt
call cpu_clock_end(id_clock_dynamics)

!===========================================================================
! This is the start of the tracer advection part of the algorithm.

if (thermo_does_span_coupling) then
do_advection = (CS%t_dyn_rel_adv + 0.5*dt > dt_therm)
else
do_advection = ((MOD(n,ntstep) == 0) .or. (n==n_max))
endif

! do advection and (perhaps) thermodynamics
if (do_advection) then
if (do_advection) then ! Do advective transport and lateral tracer mixing.

if (CS%debug) then
call cpu_clock_begin(id_clock_other)
Expand Down Expand Up @@ -2700,12 +2721,12 @@ end subroutine MOM_timing_init

!> This routine posts diagnostics of the transports, including the subgridscale
!! contributions.
subroutine post_transport_diagnostics(G, GV, CS, diag, dt, h, h_pre_dyn)
subroutine post_transport_diagnostics(G, GV, CS, diag, dt_trans, h, h_pre_dyn)
type(ocean_grid_type), intent(inout) :: G !< ocean grid structure
type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure
type(MOM_control_struct), intent(in) :: CS !< control structure
type(diag_ctrl), intent(inout) :: diag !< regulates diagnostic output
real , intent(in) :: dt !< total time step associated with the transports, in s.
real , intent(in) :: dt_trans !< total time step associated with the transports, in s.
real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
intent(in) :: h !< The updated layer thicknesses, in H
real, dimension(SZI_(G),SZJ_(G),SZK_(G)), &
Expand All @@ -2720,15 +2741,15 @@ subroutine post_transport_diagnostics(G, GV, CS, diag, dt, h, h_pre_dyn)
is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = G%ke

call cpu_clock_begin(id_clock_Z_diag)
call calculate_Z_transport(CS%uhtr, CS%vhtr, h, CS%t_dyn_rel_adv, G, GV, &
call calculate_Z_transport(CS%uhtr, CS%vhtr, h, dt_trans, G, GV, &
CS%diag_to_Z_CSp)
call cpu_clock_end(id_clock_Z_diag)

! Post mass transports, including SGS
! Build the remap grids using the layer thicknesses from before the dynamics
call diag_update_remap_grids(diag, alt_h = h_pre_dyn)

H_to_kg_m2_dt = GV%H_to_kg_m2 / CS%t_dyn_rel_adv
H_to_kg_m2_dt = GV%H_to_kg_m2 / dt_trans
if (CS%id_umo_2d > 0) then
umo2d(:,:) = 0.0
do k=1,nz ; do j=js,je ; do I=is-1,ie
Expand Down

0 comments on commit 4d7e1d9

Please sign in to comment.