benchmark_initialization Module Reference

Detailed Description

Initialization for the "bench mark" configuration.

Functions/Subroutines
subroutine, public	benchmark_initialize_topography (D, G, param_file, max_depth, US)
	This subroutine sets up the benchmark test case topography. More...
subroutine, public	benchmark_initialize_thickness (h, G, GV, US, param_file, eqn_of_state, P_ref, just_read_params)
	Initializes layer thicknesses for the benchmark test case, by finding the depths of interfaces in a specified latitude-dependent temperature profile with an exponentially decaying thermocline on top of a linear stratification. More...
subroutine, public	benchmark_init_temperature_salinity (T, S, G, GV, US, param_file, eqn_of_state, P_Ref, just_read_params)
	Initializes layer temperatures and salinities for benchmark. More...

Function/Subroutine Documentation

benchmark_init_temperature_salinity()

subroutine, public benchmark_initialization::benchmark_init_temperature_salinity	(	real, dimension(szi_(g),szj_(g), szk_(g)), intent(out)	T,
		real, dimension(szi_(g),szj_(g), szk_(g)), intent(out)	S,
		type(ocean_grid_type), intent(in)	G,
		type(verticalgrid_type), intent(in)	GV,
		type(unit_scale_type), intent(in)	US,
		type(param_file_type), intent(in)	param_file,
		type(eos_type), pointer	eqn_of_state,
		real, intent(in)	P_Ref,
		logical, intent(in), optional	just_read_params
	)

Initializes layer temperatures and salinities for benchmark.

Parameters

[in]	g	The ocean's grid structure.
[in]	gv	The ocean's vertical grid structure.
[out]	t	The potential temperature that is being initialized.
[out]	s	The salinity that is being initialized.
[in]	us	A dimensional unit scaling type
[in]	param_file	A structure indicating the open file to parse for model parameter values.
	eqn_of_state	integer that selects the equation of state.
[in]	p_ref	The coordinate-density reference pressure [Pa].
[in]	just_read_params	If present and true, this call will only read parameters without changing h.

Definition at line 217 of file benchmark_initialization.F90.

  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_(G)), intent(out) :: T      !< The potential temperature
                                                                   !! that is being initialized.
  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: S      !< The salinity that is being
                                                                   !! initialized.
  type(unit_scale_type),                     intent(in)  :: US     !< A dimensional unit scaling type
  type(param_file_type),               intent(in)  :: param_file   !< A structure indicating the
                                                                   !! open file to parse for
                                                                   !! model parameter values.
  type(EOS_type),                      pointer     :: eqn_of_state !< integer that selects the
                                                                   !! equation of state.
  real,                                intent(in)  :: P_Ref        !< The coordinate-density
                                                                   !! reference pressure [Pa].
  logical,       optional, intent(in)  :: just_read_params !< If present and true, this call will
                                                      !! only read parameters without changing h.
  ! Local variables
  real :: T0(SZK_(G)), S0(SZK_(G))
  real :: pres(SZK_(G))      ! Reference pressure [kg m-3].
  real :: drho_dT(SZK_(G))   ! Derivative of density with temperature [R degC-1 ~> kg m-3 degC-1].
  real :: drho_dS(SZK_(G))   ! Derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1].
  real :: rho_guess(SZK_(G)) ! Potential density at T0 & S0 [R ~> kg m-3].
  real :: PI        ! 3.1415926... calculated as 4*atan(1)
  real :: SST       !  The initial sea surface temperature [degC].
  real :: lat
  logical :: just_read    ! If true, just read parameters but set nothing.
  character(len=40)  :: mdl = "benchmark_init_temperature_salinity" ! This subroutine's name.
  integer :: i, j, k, k1, is, ie, js, je, nz, itt
 
  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
 
  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
 
  if (just_read) return ! All run-time parameters have been read, so return.
 
  k1 = gv%nk_rho_varies + 1
 
  do k=1,nz
    pres(k) = p_ref ; s0(k) = 35.0
  enddo
 
  t0(k1) = 29.0
  call calculate_density(t0(k1),s0(k1),pres(k1),rho_guess(k1),eqn_of_state, scale=us%kg_m3_to_R)
  call calculate_density_derivs(t0,s0,pres,drho_dt,drho_ds,k1,1,eqn_of_state, scale=us%kg_m3_to_R)
 
! A first guess of the layers' temperatures.                         !
  do k=1,nz
    t0(k) = t0(k1) + (gv%Rlay(k) - rho_guess(k1)) / drho_dt(k1)
  enddo
 
! Refine the guesses for each layer.                                 !
  do itt = 1,6
    call calculate_density(t0,s0,pres,rho_guess,1,nz,eqn_of_state, scale=us%kg_m3_to_R)
    call calculate_density_derivs(t0,s0,pres,drho_dt,drho_ds,1,nz,eqn_of_state, scale=us%kg_m3_to_R)
    do k=1,nz
      t0(k) = t0(k) + (gv%Rlay(k) - rho_guess(k)) / drho_dt(k)
    enddo
  enddo
 
  do k=1,nz ; do i=is,ie ; do j=js,je
    t(i,j,k) = t0(k)
    s(i,j,k) = s0(k)
  enddo ; enddo ; enddo
  pi = 4.0*atan(1.0)
  do i=is,ie ; do j=js,je
    sst = 0.5*(t0(k1)+t0(nz)) - 0.9*0.5*(t0(k1)-t0(nz)) * &
                               cos(pi*(g%geoLatT(i,j)-g%south_lat)/(g%len_lat))
    do k=1,k1-1
      t(i,j,k) = sst
    enddo
  enddo ; enddo
 

Referenced by mom_state_initialization::mom_initialize_state().

Here is the caller graph for this function:

benchmark_initialize_thickness()

subroutine, public benchmark_initialization::benchmark_initialize_thickness	(	real, dimension(szi_(g),szj_(g),szk_(gv)), intent(out)	h,
		type(ocean_grid_type), intent(in)	G,
		type(verticalgrid_type), intent(in)	GV,
		type(unit_scale_type), intent(in)	US,
		type(param_file_type), intent(in)	param_file,
		type(eos_type), pointer	eqn_of_state,
		real, intent(in)	P_ref,
		logical, intent(in), optional	just_read_params
	)

Initializes layer thicknesses for the benchmark test case, by finding the depths of interfaces in a specified latitude-dependent temperature profile with an exponentially decaying thermocline on top of a linear stratification.

Parameters

[in]	g	The ocean's grid structure.
[in]	gv	The ocean's vertical grid structure.
[in]	us	A dimensional unit scaling type
[out]	h	The thickness that is being initialized [H ~> m or kg m-2].
[in]	param_file	A structure indicating the open file to parse for model parameter values.
	eqn_of_state	integer that selects the equation of state.
[in]	p_ref	The coordinate-density reference pressure [Pa].
[in]	just_read_params	If present and true, this call will only read parameters without changing h.

Definition at line 87 of file benchmark_initialization.F90.

  type(ocean_grid_type),   intent(in)  :: G           !< The ocean's grid structure.
  type(verticalGrid_type), intent(in)  :: GV          !< The ocean's vertical grid structure.
  type(unit_scale_type),   intent(in)  :: US !< A dimensional unit scaling type
  real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), &
                           intent(out) :: h           !< The thickness that is being initialized [H ~> m or kg m-2].
  type(param_file_type),   intent(in)  :: param_file  !< A structure indicating the open file
                                                      !! to parse for model parameter values.
  type(EOS_type),          pointer     :: eqn_of_state !< integer that selects the
                                                      !! equation of state.
  real,                    intent(in)  :: P_Ref       !< The coordinate-density
                                                      !! reference pressure [Pa].
  logical,       optional, intent(in)  :: just_read_params !< If present and true, this call will
                                                      !! only read parameters without changing h.
  ! Local variables
  real :: e0(SZK_(GV)+1)     ! The resting interface heights, in depth units [Z ~> m],
                             ! usually negative because it is positive upward.
  real :: e_pert(SZK_(GV)+1) ! Interface height perturbations, positive upward,
                             ! in depth units [Z ~> m].
  real :: eta1D(SZK_(GV)+1)  ! Interface height relative to the sea surface
                             ! positive upward, in depth units [Z ~> m].
  real :: SST       !  The initial sea surface temperature [degC].
  real :: T_int     !  The initial temperature of an interface [degC].
  real :: ML_depth  !  The specified initial mixed layer depth, in depth units [Z ~> m].
  real :: thermocline_scale ! The e-folding scale of the thermocline, in depth units [Z ~> m].
  real, dimension(SZK_(GV)) :: &
    T0, pres, S0, & ! drho
    rho_guess, &    ! Potential density at T0 & S0 [R ~> kg m-3].
    drho_dT, &      ! Derivative of density with temperature [R degC-1 ~> kg m-3 degC-1].
    drho_dS         ! Derivative of density with salinity [R ppt-1 ~> kg m-3 ppt-1].
  real :: a_exp     ! The fraction of the overall stratification that is exponential.
  real :: I_ts, I_md ! Inverse lengthscales [Z-1 ~> m-1].
  real :: T_frac    ! A ratio of the interface temperature to the range
                    ! between SST and the bottom temperature.
  real :: err, derr_dz  ! The error between the profile's temperature and the
                    ! interface temperature for a given z and its derivative.
  real :: pi, z
  logical :: just_read
  ! This include declares and sets the variable "version".
# include "version_variable.h"
  character(len=40)  :: mdl = "benchmark_initialize_thickness" ! This subroutine's name.
  integer :: i, j, k, k1, is, ie, js, je, nz, itt
 
  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
 
  just_read = .false. ; if (present(just_read_params)) just_read = just_read_params
  if (.not.just_read) call log_version(param_file, mdl, version, "")
  call get_param(param_file, mdl, "BENCHMARK_ML_DEPTH_IC", ml_depth, &
                 "Initial mixed layer depth in the benchmark test case.", &
                 units='m', default=50.0, scale=us%m_to_Z, do_not_log=just_read)
  call get_param(param_file, mdl, "BENCHMARK_THERMOCLINE_SCALE", thermocline_scale, &
                 "Initial thermocline depth scale in the benchmark test case.", &
                 default=500.0, units="m", scale=us%m_to_Z, do_not_log=just_read)
 
  if (just_read) return ! This subroutine has no run-time parameters.
 
  call mom_mesg("  benchmark_initialization.F90, benchmark_initialize_thickness: setting thickness", 5)
 
  k1 = gv%nk_rho_varies + 1
 
  a_exp = 0.9
 
! This block calculates T0(k) for the purpose of diagnosing where the
! interfaces will be found.
  do k=1,nz
    pres(k) = p_ref ; s0(k) = 35.0
  enddo
  t0(k1) = 29.0
  call calculate_density(t0(k1), s0(k1), pres(k1), rho_guess(k1), eqn_of_state, scale=us%kg_m3_to_R)
  call calculate_density_derivs(t0, s0, pres, drho_dt, drho_ds, k1, 1, eqn_of_state, scale=us%kg_m3_to_R)
 
! A first guess of the layers' temperatures.
  do k=1,nz
    t0(k) = t0(k1) + (gv%Rlay(k) - rho_guess(k1)) / drho_dt(k1)
  enddo
 
! Refine the guesses for each layer.
  do itt=1,6
    call calculate_density(t0, s0, pres, rho_guess, 1, nz, eqn_of_state, scale=us%kg_m3_to_R)
    call calculate_density_derivs(t0, s0, pres, drho_dt, drho_ds, 1, nz, eqn_of_state, scale=us%kg_m3_to_R)
    do k=1,nz
      t0(k) = t0(k) + (gv%Rlay(k) - rho_guess(k)) / drho_dt(k)
    enddo
  enddo
 
  pi = 4.0*atan(1.0)
  i_ts = 1.0 / thermocline_scale
  i_md = 1.0 / g%max_depth
  do j=js,je ; do i=is,ie
    sst = 0.5*(t0(k1)+t0(nz)) - 0.9*0.5*(t0(k1)-t0(nz)) * &
                               cos(pi*(g%geoLatT(i,j)-g%south_lat)/(g%len_lat))
 
    do k=1,nz ; e_pert(k) = 0.0 ; enddo
 
    !   This sets the initial thickness (in [H ~> m or kg m-2]) of the layers.  The thicknesses
    ! are set to insure that: 1. each layer is at least  Gv%Angstrom_m thick, and
    ! 2. the interfaces are where they should be based on the resting depths and interface
    ! height perturbations, as long at this doesn't interfere with 1.
    eta1d(nz+1) = -g%bathyT(i,j)
 
    do k=nz,2,-1
      t_int = 0.5*(t0(k) + t0(k-1))
      t_frac = (t_int - t0(nz)) / (sst - t0(nz))
      ! Find the z such that T_frac = a exp(z/thermocline_scale) + (1-a) (z+D)/D
      z = 0.0
      do itt=1,6
        err = a_exp * exp(z*i_ts) + (1.0 - a_exp) * (z*i_md + 1.0) - t_frac
        derr_dz = a_exp * i_ts * exp(z*i_ts) + (1.0 - a_exp) * i_md
        z = z - err / derr_dz
      enddo
      e0(k) = z
!       e0(K) = -ML_depth + thermocline_scale * log((T_int - T0(nz)) / (SST - T0(nz)))
 
      eta1d(k) = e0(k) + e_pert(k)
 
      if (eta1d(k) > -ml_depth) eta1d(k) = -ml_depth
 
      if (eta1d(k) < eta1d(k+1) + gv%Angstrom_Z) &
        eta1d(k) = eta1d(k+1) + gv%Angstrom_Z
 
      h(i,j,k) = max(gv%Z_to_H * (eta1d(k) - eta1d(k+1)), gv%Angstrom_H)
    enddo
    h(i,j,1) = max(gv%Z_to_H * (0.0 - eta1d(2)), gv%Angstrom_H)
 
  enddo ; enddo
 

References mom_error_handler::mom_mesg().

Here is the call graph for this function:

benchmark_initialize_topography()

subroutine, public benchmark_initialization::benchmark_initialize_topography	(	real, dimension(g%isd:g%ied,g%jsd:g%jed), intent(out)	D,
		type(dyn_horgrid_type), intent(in)	G,
		type(param_file_type), intent(in)	param_file,
		real, intent(in)	max_depth,
		type(unit_scale_type), intent(in), optional	US
	)

This subroutine sets up the benchmark test case topography.

Parameters

[in]	g	The dynamic horizontal grid type
[out]	d	Ocean bottom depth in m or [Z ~> m] if US is present
[in]	param_file	Parameter file structure
[in]	max_depth	Maximum model depth in the units of D
[in]	us	A dimensional unit scaling type

Definition at line 35 of file benchmark_initialization.F90.

  type(dyn_horgrid_type),          intent(in)  :: G !< The dynamic horizontal grid type
  real, dimension(G%isd:G%ied,G%jsd:G%jed), &
                                   intent(out) :: D !< Ocean bottom depth in m or [Z ~> m] if US is present
  type(param_file_type),           intent(in)  :: param_file !< Parameter file structure
  real,                            intent(in)  :: max_depth !< Maximum model depth in the units of D
  type(unit_scale_type), optional, intent(in)  :: US !< A dimensional unit scaling type
 
  ! Local variables
  real :: min_depth            ! The minimum and maximum depths [Z ~> m].
  real :: PI                   ! 3.1415926... calculated as 4*atan(1)
  real :: D0                   ! A constant to make the maximum     !
                               ! basin depth MAXIMUM_DEPTH.         !
  real :: m_to_Z  ! A dimensional rescaling factor.
  real :: x, y
  ! This include declares and sets the variable "version".
# include "version_variable.h"
  character(len=40)  :: mdl = "benchmark_initialize_topography" ! This subroutine's name.
  integer :: i, j, is, ie, js, je, isd, ied, jsd, jed
  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
 
  call mom_mesg("  benchmark_initialization.F90, benchmark_initialize_topography: setting topography", 5)
 
  m_to_z = 1.0 ; if (present(us)) m_to_z = us%m_to_Z
 
  call log_version(param_file, mdl, version, "")
  call get_param(param_file, mdl, "MINIMUM_DEPTH", min_depth, &
                 "The minimum depth of the ocean.", units="m", default=0.0, scale=m_to_z)
 
  pi = 4.0*atan(1.0)
  d0 = max_depth / 0.5
 
!  Calculate the depth of the bottom.
  do j=js,je ; do i=is,ie
    x = (g%geoLonT(i,j)-g%west_lon) / g%len_lon
    y = (g%geoLatT(i,j)-g%south_lat) / g%len_lat
!  This sets topography that has a reentrant channel to the south.
    d(i,j) = -d0 * ( y*(1.0 + 0.6*cos(4.0*pi*x)) &
                   + 0.75*exp(-6.0*y) &
                   + 0.05*cos(10.0*pi*x) - 0.7 )
    if (d(i,j) > max_depth) d(i,j) = max_depth
    if (d(i,j) < min_depth) d(i,j) = 0.
  enddo ; enddo
 

References mom_error_handler::mom_mesg().

Referenced by mom_fixed_initialization::mom_initialize_topography().

Here is the call graph for this function:

Here is the caller graph for this function: