seamount_initialization Module Reference

Detailed Description

Configures the model for the idealized seamount test case.

Functions/Subroutines
subroutine, public	seamount_initialize_topography (D, G, param_file, max_depth)
	Initialization of topography. More...
subroutine, public	seamount_initialize_thickness (h, G, GV, US, param_file, just_read_params)
	Initialization of thicknesses. This subroutine initializes the layer thicknesses to be uniform. More...
subroutine, public	seamount_initialize_temperature_salinity (T, S, h, G, GV, param_file, eqn_of_state, just_read_params)
	Initial values for temperature and salinity. More...

Variables
character(len=40)	mdl = "seamount_initialization"
	This module's name. More...

Function/Subroutine Documentation

seamount_initialize_temperature_salinity()

subroutine, public seamount_initialization::seamount_initialize_temperature_salinity	(	real, dimension(szi_(g),szj_(g), szk_(g)), intent(out)	T,
		real, dimension(szi_(g),szj_(g), szk_(g)), intent(out)	S,
		real, dimension(szi_(g),szj_(g), szk_(g)), intent(in)	h,
		type(ocean_grid_type), intent(in)	G,
		type(verticalgrid_type), intent(in)	GV,
		type(param_file_type), intent(in)	param_file,
		type(eos_type), pointer	eqn_of_state,
		logical, intent(in), optional	just_read_params
	)

Initial values for temperature and salinity.

Parameters

[in]	g	Ocean grid structure
[in]	gv	Vertical grid structure
[out]	t	Potential temperature [degC]
[out]	s	Salinity [ppt]
[in]	h	Layer thickness [H ~> m or kg m-2]
[in]	param_file	Parameter file structure
	eqn_of_state	Equation of state structure
[in]	just_read_params	If present and true, this call will only read parameters without changing h.

Definition at line 196 of file seamount_initialization.F90.

  type(ocean_grid_type),                     intent(in)  :: G !< Ocean grid structure
  type(verticalGrid_type),                   intent(in)  :: GV !< Vertical grid structure
  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: T !< Potential temperature [degC]
  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(out) :: S !< Salinity [ppt]
  real, dimension(SZI_(G),SZJ_(G), SZK_(G)), intent(in)  :: h !< Layer thickness [H ~> m or kg m-2]
  type(param_file_type),                     intent(in)  :: param_file !< Parameter file structure
  type(EOS_type),                            pointer     :: eqn_of_state !< Equation of state structure
  logical,       optional, intent(in)  :: just_read_params !< If present and true, this call will
                                                      !! only read parameters without changing h.
 
  ! Local variables
  integer :: i, j, k, is, ie, js, je, nz, k_light
  real    :: xi0, xi1, dxi, r, S_surf, T_surf, S_range, T_range
  real    :: T_ref, T_Light, T_Dense, S_ref, S_Light, S_Dense, a1, frac_dense, k_frac, res_rat
  logical :: just_read    ! If true, just read parameters but set nothing.
  character(len=20) :: verticalCoordinate, density_profile
 
  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
 
  call get_param(param_file, mdl, "REGRIDDING_COORDINATE_MODE", verticalcoordinate, &
                 default=default_coordinate_mode, do_not_log=just_read)
  call get_param(param_file, mdl,"INITIAL_DENSITY_PROFILE", density_profile, &
                 'Initial profile shape. Valid values are "linear", "parabolic" '//&
                 'and "exponential".', default='linear', do_not_log=just_read)
  call get_param(param_file, mdl,"INITIAL_SSS", s_surf, &
                 'Initial surface salinity', units='1e-3', default=34., do_not_log=just_read)
  call get_param(param_file, mdl,"INITIAL_SST", t_surf, &
                 'Initial surface temperature', units='C', default=0., do_not_log=just_read)
  call get_param(param_file, mdl,"INITIAL_S_RANGE", s_range, &
                 'Initial salinity range (bottom - surface)', units='1e-3', &
                 default=2., do_not_log=just_read)
  call get_param(param_file, mdl,"INITIAL_T_RANGE", t_range, &
                 'Initial temperature range (bottom - surface)', units='C', &
                 default=0., do_not_log=just_read)
 
  select case ( coordinatemode(verticalcoordinate) )
    case ( regridding_layer ) ! Initial thicknesses for layer isopycnal coordinates
      ! These parameters are used in MOM_fixed_initialization.F90 when CONFIG_COORD="ts_range"
      call get_param(param_file, mdl, "T_REF", t_ref, default=10.0, do_not_log=.true.)
      call get_param(param_file, mdl, "TS_RANGE_T_LIGHT", t_light, default=t_ref, do_not_log=.true.)
      call get_param(param_file, mdl, "TS_RANGE_T_DENSE", t_dense, default=t_ref, do_not_log=.true.)
      call get_param(param_file, mdl, "S_REF", s_ref, default=35.0, do_not_log=.true.)
      call get_param(param_file, mdl, "TS_RANGE_S_LIGHT", s_light, default = s_ref, do_not_log=.true.)
      call get_param(param_file, mdl, "TS_RANGE_S_DENSE", s_dense, default = s_ref, do_not_log=.true.)
      call get_param(param_file, mdl, "TS_RANGE_RESOLN_RATIO", res_rat, default=1.0, do_not_log=.true.)
      if (just_read) return ! All run-time parameters have been read, so return.
 
      ! Emulate the T,S used in the "ts_range" coordinate configuration code
      k_light = gv%nk_rho_varies + 1
      do j=js,je ; do i=is,ie
        t(i,j,k_light) = t_light ; s(i,j,k_light) = s_light
      enddo ; enddo
      a1 = 2.0 * res_rat / (1.0 + res_rat)
      do k=k_light+1,nz
        k_frac = real(k-k_light)/real(nz-k_light)
        frac_dense = a1 * k_frac + (1.0 - a1) * k_frac**2
        do j=js,je ; do i=is,ie
          t(i,j,k) = frac_dense * (t_dense - t_light) + t_light
          s(i,j,k) = frac_dense * (s_dense - s_light) + s_light
        enddo ; enddo
      enddo
    case ( regridding_sigma, regridding_zstar, regridding_rho ) ! All other coordinate use FV initialization
      if (just_read) return ! All run-time parameters have been read, so return.
      do j=js,je ; do i=is,ie
        xi0 = 0.0
        do k = 1,nz
          xi1 = xi0 + gv%H_to_Z * h(i,j,k) / g%max_depth
          select case ( trim(density_profile) )
            case ('linear')
             !S(i,j,k) = S_surf + S_range * 0.5 * (xi0 + xi1)
              s(i,j,k) = s_surf + ( 0.5 * s_range ) * (xi0 + xi1) ! Coded this way to reproduce old hard-coded answers
              t(i,j,k) = t_surf + t_range * 0.5 * (xi0 + xi1)
            case ('parabolic')
              s(i,j,k) = s_surf + s_range * (2.0 / 3.0) * (xi1**3 - xi0**3) / (xi1 - xi0)
              t(i,j,k) = t_surf + t_range * (2.0 / 3.0) * (xi1**3 - xi0**3) / (xi1 - xi0)
            case ('exponential')
              r = 0.8 ! small values give sharp profiles
              s(i,j,k) = s_surf + s_range * (exp(xi1/r)-exp(xi0/r)) / (xi1 - xi0)
              t(i,j,k) = t_surf + t_range * (exp(xi1/r)-exp(xi0/r)) / (xi1 - xi0)
            case default
              call mom_error(fatal, 'Unknown value for "INITIAL_DENSITY_PROFILE"')
          end select
          xi0 = xi1
        enddo
      enddo ; enddo
  end select
 

References mdl, mom_error_handler::mom_error(), regrid_consts::regridding_rho, and regrid_consts::regridding_sigma.

Referenced by mom_state_initialization::mom_initialize_state().

Here is the call graph for this function:

Here is the caller graph for this function:

seamount_initialize_thickness()

subroutine, public seamount_initialization::seamount_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,
		logical, intent(in), optional	just_read_params
	)

Initialization of thicknesses. This subroutine initializes the layer thicknesses to be uniform.

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.
[in]	just_read_params	If present and true, this call will only read parameters without changing h.

Definition at line 81 of file seamount_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.
  logical,       optional, intent(in)  :: just_read_params !< If present and true, this call will
                                                      !! only read parameters without changing h.
 
  real :: e0(SZK_(G)+1)   ! The resting interface heights [Z ~> m], usually
                          ! negative because it is positive upward.
  real :: eta1D(SZK_(G)+1)! Interface height relative to the sea surface
                          ! positive upward [Z ~> m].
  real :: min_thickness   ! The minimum layer thicknesses [Z ~> m].
  real :: S_surf, S_range, S_ref, S_light, S_dense ! Various salinities [ppt].
  real :: eta_IC_quanta   ! The granularity of quantization of intial interface heights [Z-1 ~> m-1].
  character(len=20) :: verticalCoordinate
  logical :: just_read    ! If true, just read parameters but set nothing.
  integer :: i, j, k, is, ie, js, je, nz
 
  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 mom_mesg("MOM_initialization.F90, initialize_thickness_uniform: setting thickness")
 
  call get_param(param_file, mdl,"MIN_THICKNESS",min_thickness, &
                'Minimum thickness for layer', &
                 units='m', default=1.0e-3, do_not_log=just_read, scale=us%m_to_Z)
  call get_param(param_file, mdl,"REGRIDDING_COORDINATE_MODE",verticalcoordinate, &
                 default=default_coordinate_mode, do_not_log=just_read)
 
  ! WARNING: this routine specifies the interface heights so that the last layer
  !          is vanished, even at maximum depth. In order to have a uniform
  !          layer distribution, use this line of code within the loop:
  !          e0(k) = -G%max_depth * real(k-1) / real(nz)
  !          To obtain a thickness distribution where the last layer is
  !          vanished and the other thicknesses uniformly distributed, use:
  !          e0(k) = -G%max_depth * real(k-1) / real(nz-1)
  !do k=1,nz+1
  !  e0(k) = -G%max_depth * real(k-1) / real(nz)
  !enddo
 
  select case ( coordinatemode(verticalcoordinate) )
 
  case ( regridding_layer, regridding_rho ) ! Initial thicknesses for isopycnal coordinates
    call get_param(param_file, mdl,"INITIAL_SSS", s_surf, default=34., do_not_log=.true.)
    call get_param(param_file, mdl,"INITIAL_S_RANGE", s_range, default=2., do_not_log=.true.)
    call get_param(param_file, mdl, "S_REF", s_ref, default=35.0, do_not_log=.true.)
    call get_param(param_file, mdl, "TS_RANGE_S_LIGHT", s_light, default = s_ref, do_not_log=.true.)
    call get_param(param_file, mdl, "TS_RANGE_S_DENSE", s_dense, default = s_ref, do_not_log=.true.)
    call get_param(param_file, mdl, "INTERFACE_IC_QUANTA", eta_ic_quanta, &
                   "The granularity of initial interface height values "//&
                   "per meter, to avoid sensivity to order-of-arithmetic changes.", &
                   default=2048.0, units="m-1", scale=us%Z_to_m, do_not_log=just_read)
    if (just_read) return ! All run-time parameters have been read, so return.
 
    do k=1,nz+1
      ! Salinity of layer k is S_light + (k-1)/(nz-1) * (S_dense - S_light)
      ! Salinity of interface K is S_light + (K-3/2)/(nz-1) * (S_dense - S_light)
      ! Salinity at depth z should be S(z) = S_surf - S_range * z/max_depth
      ! Equating: S_surf - S_range * z/max_depth = S_light + (K-3/2)/(nz-1) * (S_dense - S_light)
      ! Equating: - S_range * z/max_depth = S_light - S_surf + (K-3/2)/(nz-1) * (S_dense - S_light)
      ! Equating: z/max_depth = - ( S_light - S_surf + (K-3/2)/(nz-1) * (S_dense - S_light) ) / S_range
      e0(k) = - g%max_depth * ( ( s_light  - s_surf ) + ( s_dense - s_light ) * &
                              ( (real(k)-1.5) / real(nz-1) ) ) / s_range
      ! Force round numbers ... the above expression has irrational factors ...
      if (eta_ic_quanta > 0.0) &
        e0(k) = nint(eta_ic_quanta*e0(k)) / eta_ic_quanta
      e0(k) = min(real(1-k)*gv%Angstrom_Z, e0(k)) ! Bound by surface
      e0(k) = max(-g%max_depth, e0(k)) ! Bound by bottom
    enddo
    do j=js,je ; do i=is,ie
      eta1d(nz+1) = -g%bathyT(i,j)
      do k=nz,1,-1
        eta1d(k) = e0(k)
        if (eta1d(k) < (eta1d(k+1) + gv%Angstrom_Z)) then
          eta1d(k) = eta1d(k+1) + gv%Angstrom_Z
          h(i,j,k) = gv%Angstrom_H
        else
          h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
        endif
      enddo
    enddo ; enddo
 
  case ( regridding_zstar )                       ! Initial thicknesses for z coordinates
    if (just_read) return ! All run-time parameters have been read, so return.
    do j=js,je ; do i=is,ie
      eta1d(nz+1) = -g%bathyT(i,j)
      do k=nz,1,-1
        eta1d(k) =  -g%max_depth * real(k-1) / real(nz)
        if (eta1d(k) < (eta1d(k+1) + min_thickness)) then
          eta1d(k) = eta1d(k+1) + min_thickness
          h(i,j,k) = gv%Z_to_H * min_thickness
        else
          h(i,j,k) = gv%Z_to_H * (eta1d(k) - eta1d(k+1))
        endif
      enddo
    enddo ; enddo
 
  case ( regridding_sigma )             ! Initial thicknesses for sigma coordinates
    if (just_read) return ! All run-time parameters have been read, so return.
    do j=js,je ; do i=is,ie
      h(i,j,:) = gv%Z_to_H * g%bathyT(i,j) / dfloat(nz)
    enddo ; enddo
 
end select
 

References mdl, mom_error_handler::mom_mesg(), regrid_consts::regridding_rho, and regrid_consts::regridding_sigma.

Here is the call graph for this function:

seamount_initialize_topography()

subroutine, public seamount_initialization::seamount_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
	)

Initialization of topography.

Parameters

[in]	g	The dynamic horizontal grid type
[out]	d	Ocean bottom depth in the units of depth_max
[in]	param_file	Parameter file structure
[in]	max_depth	Maximum ocean depth in arbitrary units

Definition at line 42 of file seamount_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 the units of depth_max
  type(param_file_type),   intent(in)  :: param_file !< Parameter file structure
  real,                    intent(in)  :: max_depth !< Maximum ocean depth in arbitrary units
 
  ! Local variables
  integer   :: i, j
  real      :: x, y, delta, Lx, rLx, Ly, rLy
 
  call get_param(param_file, mdl,"SEAMOUNT_DELTA",delta, &
                 "Non-dimensional height of seamount.", &
                 units="non-dim", default=0.5)
  call get_param(param_file, mdl,"SEAMOUNT_X_LENGTH_SCALE",lx, &
                 "Length scale of seamount in x-direction. "//&
                 "Set to zero make topography uniform in the x-direction.", &
                 units="Same as x,y", default=20.)
  call get_param(param_file, mdl,"SEAMOUNT_Y_LENGTH_SCALE",ly, &
                 "Length scale of seamount in y-direction. "//&
                 "Set to zero make topography uniform in the y-direction.", &
                 units="Same as x,y", default=0.)
 
  lx = lx / g%len_lon
  ly = ly / g%len_lat
  rlx = 0. ; if (lx>0.) rlx = 1. / lx
  rly = 0. ; if (ly>0.) rly = 1. / ly
 
  do j=g%jsc,g%jec ; do i=g%isc,g%iec
    ! Compute normalized zonal coordinates (x,y=0 at center of domain)
    x = ( g%geoLonT(i,j) - g%west_lon ) / g%len_lon - 0.5
    y = ( g%geoLatT(i,j) - g%south_lat ) / g%len_lat - 0.5
    d(i,j) = g%max_depth * ( 1.0 - delta * exp(-(rlx*x)**2 -(rly*y)**2) )
  enddo ; enddo
 

References mdl.

Referenced by mom_fixed_initialization::mom_initialize_topography().

Here is the caller graph for this function:

Variable Documentation

mdl

character(len=40) seamount_initialization::mdl = "seamount_initialization"

private

This module's name.

Definition at line 26 of file seamount_initialization.F90.

character(len=40) :: mdl = "seamount_initialization" !< This module's name.

Referenced by seamount_initialize_temperature_salinity(), seamount_initialize_thickness(), and seamount_initialize_topography().