Skip to content

Commit

Permalink
A new tracer that keeps track of "mixed layer age" has been added to …
Browse files Browse the repository at this point in the history 
…the ideal age

module. This PR also adds the ability to use the actual BL depth that is diagnosed by the
active BL scheme inside the ideal age module (for all ideal age tracers).
  • Loading branch information
@sdbachman
sdbachman committed Jul 29, 2022
1 parent 9e27b52 commit bb02a51
Show file tree
Hide file tree
Showing 2 changed files with 199 additions and 33 deletions.
5 changes: 3 additions & 2 deletions src/tracer/MOM_tracer_flow_control.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -468,7 +468,8 @@ subroutine call_tracer_column_fns(h_old, h_new, ea, eb, fluxes, Hml, dt, G, GV,
call ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, &
G, GV, US, CS%ideal_age_tracer_CSp, &
evap_CFL_limit=evap_CFL_limit, &
minimum_forcing_depth=minimum_forcing_depth)
minimum_forcing_depth=minimum_forcing_depth, &
Hml=Hml)
if (CS%use_regional_dyes) &
call dye_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, &
G, GV, US, CS%dye_tracer_CSp, &
Expand Down Expand Up @@ -544,7 +545,7 @@ subroutine call_tracer_column_fns(h_old, h_new, ea, eb, fluxes, Hml, dt, G, GV,
G, GV, US, CS%RGC_tracer_CSp)
if (CS%use_ideal_age) &
call ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, &
G, GV, US, CS%ideal_age_tracer_CSp)
G, GV, US, CS%ideal_age_tracer_CSp, Hml=Hml)
if (CS%use_regional_dyes) &
call dye_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, &
G, GV, US, CS%dye_tracer_CSp)
Expand Down
227 changes: 196 additions & 31 deletions src/tracer/ideal_age_example.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ module ideal_age_example
use MOM_coms, only : EFP_type
use MOM_coupler_types, only : set_coupler_type_data, atmos_ocn_coupler_flux
use MOM_diag_mediator, only : diag_ctrl
use MOM_error_handler, only : MOM_error, FATAL, WARNING
use MOM_error_handler, only : MOM_error, FATAL, WARNING, NOTE
use MOM_file_parser, only : get_param, log_param, log_version, param_file_type
use MOM_forcing_type, only : forcing
use MOM_grid, only : ocean_grid_type
Expand All @@ -31,8 +31,9 @@ module ideal_age_example
public register_ideal_age_tracer, initialize_ideal_age_tracer
public ideal_age_tracer_column_physics, ideal_age_tracer_surface_state
public ideal_age_stock, ideal_age_example_end
public count_ML_layers

integer, parameter :: NTR_MAX = 3 !< the maximum number of tracers in this module.
integer, parameter :: NTR_MAX = 4 !< the maximum number of tracers in this module.

!> The control structure for the ideal_age_tracer package
type, public :: ideal_age_tracer_CS ; private
Expand All @@ -49,9 +50,12 @@ module ideal_age_example
real, dimension(NTR_MAX) :: IC_val = 0.0 !< The (uniform) initial condition value.
real, dimension(NTR_MAX) :: young_val = 0.0 !< The value assigned to tr at the surface.
real, dimension(NTR_MAX) :: land_val = -1.0 !< The value of tr used where land is masked out.
real, dimension(NTR_MAX) :: sfc_growth_rate !< The exponential growth rate for the surface value [year-1].
real, dimension(NTR_MAX) :: growth_rate !< The exponential growth rate for the young value [year-1].
real, dimension(NTR_MAX) :: tracer_start_year !< The year in which tracers start aging, or at which the
!! surface value equals young_val, in years.
logical :: use_real_BL_depth !< If true, uses the BL scheme to determine the number of
!! layers above the BL depth instead of the fixed nkml value.
integer :: ML_residence_num !! The tracer number assigned to the ML residence tracer in this module
logical :: tracers_may_reinit !< If true, these tracers be set up via the initialization code if
!! they are not found in the restart files.
logical :: tracer_ages(NTR_MAX) !< Indicates whether each tracer ages.
Expand All @@ -64,6 +68,7 @@ module ideal_age_example
type(MOM_restart_CS), pointer :: restart_CSp => NULL() !< A pointer to the restart controls structure

type(vardesc) :: tr_desc(NTR_MAX) !< Descriptions and metadata for the tracers

end type ideal_age_tracer_CS

contains
Expand All @@ -87,7 +92,7 @@ function register_ideal_age_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
character(len=48) :: var_name ! The variable's name.
real, pointer :: tr_ptr(:,:,:) => NULL()
logical :: register_ideal_age_tracer
logical :: do_ideal_age, do_vintage, do_ideal_age_dated
logical :: do_ideal_age, do_vintage, do_ideal_age_dated, do_ML_residence
integer :: isd, ied, jsd, jed, nz, m
isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed ; nz = GV%ke

Expand All @@ -114,8 +119,14 @@ function register_ideal_age_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
"the standard ideal age tracer - i.e. is set to 0 age in "//&
"the mixed layer and ages at unit rate in the interior.", &
default=.false.)


call get_param(param_file, mdl, "DO_ML_RESIDENCE", do_ML_residence, &
"If true, use a residence tracer that is set to 0 age "//&
"in the interior and ages at unit rate in the mixed layer.", &
default=.false.)
call get_param(param_file, mdl, "USE_REAL_BL_DEPTH", CS%use_real_BL_depth, &
"If true, the ideal age tracers will use the boundary layer "//&
"depth diagnosed from the BL or bulkmixedlayer scheme.", &
default=.false.)
call get_param(param_file, mdl, "AGE_IC_FILE", CS%IC_file, &
"The file in which the age-tracer initial values can be "//&
"found, or an empty string for internal initialization.", &
Expand All @@ -139,15 +150,15 @@ function register_ideal_age_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
if (do_ideal_age) then
CS%ntr = CS%ntr + 1 ; m = CS%ntr
CS%tr_desc(m) = var_desc("age", "yr", "Ideal Age Tracer", cmor_field_name="agessc", caller=mdl)
CS%tracer_ages(m) = .true. ; CS%sfc_growth_rate(m) = 0.0
CS%tracer_ages(m) = .true. ; CS%growth_rate(m) = 0.0
CS%IC_val(m) = 0.0 ; CS%young_val(m) = 0.0 ; CS%tracer_start_year(m) = 0.0
endif

if (do_vintage) then
CS%ntr = CS%ntr + 1 ; m = CS%ntr
CS%tr_desc(m) = var_desc("vintage", "yr", "Exponential Vintage Tracer", &
caller=mdl)
CS%tracer_ages(m) = .false. ; CS%sfc_growth_rate(m) = 1.0/30.0
CS%tracer_ages(m) = .false. ; CS%growth_rate(m) = 1.0/30.0
CS%IC_val(m) = 0.0 ; CS%young_val(m) = 1e-20 ; CS%tracer_start_year(m) = 0.0
call get_param(param_file, mdl, "IDEAL_VINTAGE_START_YEAR", CS%tracer_start_year(m), &
"The date at which the ideal vintage tracer starts.", &
Expand All @@ -158,13 +169,21 @@ function register_ideal_age_tracer(HI, GV, param_file, CS, tr_Reg, restart_CS)
CS%ntr = CS%ntr + 1 ; m = CS%ntr
CS%tr_desc(m) = var_desc("age_dated","yr","Ideal Age Tracer with a Start Date",&
caller=mdl)
CS%tracer_ages(m) = .true. ; CS%sfc_growth_rate(m) = 0.0
CS%tracer_ages(m) = .true. ; CS%growth_rate(m) = 0.0
CS%IC_val(m) = 0.0 ; CS%young_val(m) = 0.0 ; CS%tracer_start_year(m) = 0.0
call get_param(param_file, mdl, "IDEAL_AGE_DATED_START_YEAR", CS%tracer_start_year(m), &
"The date at which the dated ideal age tracer starts.", &
units="years", default=0.0)
endif

CS%ML_residence_num = 0
if (do_ML_residence) then
CS%ntr = CS%ntr + 1 ; m = CS%ntr; CS%ML_residence_num = CS%ntr
CS%tr_desc(m) = var_desc("ML_age", "yr", "ML Residence Time Tracer", caller=mdl)
CS%tracer_ages(m) = .true. ; CS%growth_rate(m) = 0.0
CS%IC_val(m) = 0.0 ; CS%young_val(m) = 0.0 ; CS%tracer_start_year(m) = 0.0
endif

allocate(CS%tr(isd:ied,jsd:jed,nz,CS%ntr), source=0.0)

do m=1,CS%ntr
Expand Down Expand Up @@ -220,6 +239,7 @@ subroutine initialize_ideal_age_tracer(restart, day, G, GV, US, h, diag, OBC, CS
logical :: OK
integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz, m
integer :: IsdB, IedB, JsdB, JedB
logical :: use_real_BL_depth

if (.not.associated(CS)) return
if (CS%ntr < 1) return
Expand Down Expand Up @@ -277,7 +297,7 @@ end subroutine initialize_ideal_age_tracer

!> Applies diapycnal diffusion, aging and regeneration at the surface to the ideal age tracers
subroutine ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G, GV, US, CS, &
evap_CFL_limit, minimum_forcing_depth)
evap_CFL_limit, minimum_forcing_depth, Hml)
type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
Expand All @@ -302,20 +322,33 @@ subroutine ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G,
!! be fluxed out of the top layer in a timestep [nondim]
real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over which
!! fluxes can be applied [H ~> m or kg m-2]
real, dimension(SZI_(G),SZJ_(G)), optional, intent(in) :: Hml !< Mixed layer depth [Z ~> m]

! This subroutine applies diapycnal diffusion and any other column
! tracer physics or chemistry to the tracers from this file.
! This is a simple example of a set of advected passive tracers.

! The arguments to this subroutine are redundant in that
! h_new(k) = h_old(k) + ea(k) - eb(k-1) + eb(k) - ea(k+1)
! Local variables
real, dimension(SZI_(G),SZJ_(G)) :: ML_layers ! Stores number of layers in mixed layer
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: h_work ! Used so that h can be modified
real :: sfc_val ! The surface value for the tracers.
real :: young_val ! The "young" value for the tracers.
real :: Isecs_per_year ! The inverse of the amount of time in a year [T-1 ~> s-1]
real :: year ! The time in years.
integer :: i, j, k, is, ie, js, je, nz, m
real :: layer_frac
integer :: i, j, k, is, ie, js, je, nz, m, nk
character(len=255) :: msg
is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke

if (CS%use_real_BL_depth .and. .not. present(Hml)) then
call MOM_error(FATAL,"Attempting to use real boundary layer depth for ideal age tracers, but no valid boundary layer scheme was found")
endif

if (CS%use_real_BL_depth .and. present(Hml)) then
call count_ML_layers(G, GV, h_old, Hml, ML_layers)
endif

if (.not.associated(CS)) return
if (CS%ntr < 1) return

Expand All @@ -340,27 +373,123 @@ subroutine ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, G,
year = US%s_to_T*time_type_to_real(CS%Time) * Isecs_per_year

do m=1,CS%ntr
if (CS%sfc_growth_rate(m) == 0.0) then
sfc_val = CS%young_val(m)

if (CS%growth_rate(m) == 0.0) then
young_val = CS%young_val(m)
else
sfc_val = CS%young_val(m) * &
exp((year-CS%tracer_start_year(m)) * CS%sfc_growth_rate(m))
young_val = CS%young_val(m) * &
exp((year-CS%tracer_start_year(m)) * CS%growth_rate(m))
endif
do k=1,CS%nkml ; do j=js,je ; do i=is,ie
if (G%mask2dT(i,j) > 0.0) then
CS%tr(i,j,k,m) = sfc_val
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif
enddo ; enddo ; enddo
enddo
do m=1,CS%ntr ; if (CS%tracer_ages(m) .and. &
(year>=CS%tracer_start_year(m))) then
!$OMP parallel do default(none) shared(is,ie,js,je,CS,nz,G,dt,Isecs_per_year,m)
do k=CS%nkml+1,nz ; do j=js,je ; do i=is,ie
CS%tr(i,j,k,m) = CS%tr(i,j,k,m) + G%mask2dT(i,j)*dt*Isecs_per_year
enddo ; enddo ; enddo
endif ; enddo

if (m == CS%ML_residence_num) then

if (CS%use_real_BL_depth) then
do j=js,je ; do i=is,ie
nk = floor(ML_layers(i,j))

do k=1,nk
if (G%mask2dT(i,j) > 0.0) then
CS%tr(i,j,k,m) = CS%tr(i,j,k,m) + G%mask2dT(i,j)*dt*Isecs_per_year
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif
enddo

k = MIN(nk+1,nz)

write(msg,*) TRIM("ML_layers= "),ML_layers(i,j), TRIM(", k= "),(k)
call MOM_error(NOTE,msg)

if (G%mask2dT(i,j) > 0.0) then
layer_frac = ML_layers(i,j)-nk
layer_frac = 0.9
CS%tr(i,j,k,m) = layer_frac * (CS%tr(i,j,k,m) + G%mask2dT(i,j)*dt*Isecs_per_year) + (1.-layer_frac) * young_val
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif


do k=nk+2,nz
if (G%mask2dT(i,j) > 0.0) then
CS%tr(i,j,k,m) = young_val
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif
enddo
enddo ; enddo

else ! use real BL depth
do j=js,je ; do i=is,ie
do k=1,CS%nkml
if (G%mask2dT(i,j) > 0.0) then
CS%tr(i,j,k,m) = CS%tr(i,j,k,m) + G%mask2dT(i,j)*dt*Isecs_per_year
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif
enddo

do k=CS%nkml+1,nz
if (G%mask2dT(i,j) > 0.0) then
CS%tr(i,j,k,m) = young_val
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif
enddo
enddo ; enddo

endif ! use real BL depth

else ! if ML residence tracer

if (CS%use_real_BL_depth) then
do j=js,je ; do i=is,ie
nk = floor(ML_layers(i,j))
do k=1,nk
if (G%mask2dT(i,j) > 0.0) then
CS%tr(i,j,k,m) = young_val
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif
enddo

k = MIN(nk+1,nz)
if (G%mask2dT(i,j) > 0.0) then
layer_frac = ML_layers(i,j)-nk
CS%tr(i,j,k,m) = (1.-layer_frac) * (CS%tr(i,j,k,m) + G%mask2dT(i,j)*dt*Isecs_per_year) + layer_frac * young_val
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif

do k=nk+2,nz
if (G%mask2dT(i,j) > 0.0) then
CS%tr(i,j,k,m) = CS%tr(i,j,k,m) + G%mask2dT(i,j)*dt*Isecs_per_year
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif
enddo
enddo ; enddo

else ! use real BL depth
do k=1,CS%nkml ; do j=js,je ; do i=is,ie
if (G%mask2dT(i,j) > 0.0) then
CS%tr(i,j,k,m) = young_val
else
CS%tr(i,j,k,m) = CS%land_val(m)
endif
enddo ; enddo ; enddo

if (CS%tracer_ages(m) .and. (year>=CS%tracer_start_year(m))) then
!$OMP parallel do default(none) shared(is,ie,js,je,CS,nz,G,dt,Isecs_per_year,m)
do k=CS%nkml+1,nz ; do j=js,je ; do i=is,ie
CS%tr(i,j,k,m) = CS%tr(i,j,k,m) + G%mask2dT(i,j)*dt*Isecs_per_year
enddo ; enddo ; enddo
endif


endif ! if use real BL depth
endif ! if ML residence tracer

enddo ! loop over all tracers

end subroutine ideal_age_tracer_column_physics

Expand Down Expand Up @@ -448,6 +577,42 @@ subroutine ideal_age_example_end(CS)
endif
end subroutine ideal_age_example_end

subroutine count_ML_layers(G, GV, h, Hml, ML_layers)
type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
intent(in) :: h !< Layer thicknesses [H ~> m or kg m-2].
real, dimension(SZI_(G),SZJ_(G)), intent(in) :: Hml !< Mixed layer depth [Z ~> m]
real, dimension(SZI_(G),SZJ_(G)), intent(out) :: ML_layers !< Number of model layers in the mixed layer

real :: current_depth
integer :: i, j, k, is, ie, js, je, nz, m, nk
character(len=255) :: msg
is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke

ML_layers(:,:) = 0.
do j=js,je ; do i=is,ie

! write(msg,*) TRIM("Hml= "),Hml(i,j)
! call MOM_error(NOTE,msg)
current_depth = 0.
do k=1,nz
current_depth = current_depth + h(i,j,k)*GV%H_to_Z
if (Hml(i,j) <= current_depth) then
ML_layers(i,j) = ML_layers(i,j) + (1.0 - (current_depth - Hml(i,j)) / (h(i,j,k)*GV%H_to_Z))
! write(msg,*) TRIM("ML_layers(i,j) found = "),ML_layers(i,j)
! call MOM_error(NOTE,msg)
exit
else
ML_layers(i,j) = ML_layers(i,j) + 1.0
! write(msg,*) TRIM("ML_layers(i,j) adding = "),ML_layers(i,j)
! call MOM_error(NOTE,msg)
endif
enddo
enddo ; enddo

end subroutine count_ML_layers

!> \namespace ideal_age_example
!!
!! Originally by Robert Hallberg, 2002
Expand Down

0 comments on commit bb02a51

Please sign in to comment.