Detailed Description

By Robert Hallberg, May 2012.

This file contains a subroutine that increments the diapycnal diffusivity in a specified band of latitudes and densities.

A small fragment of the grid is shown below:

j+1 x ^ x ^ x At x: q j+1 > o > o > At ^: v j x ^ x ^ x At >: u j > o > o > At o: h, T, S, Kd, etc. j-1 x ^ x ^ x i-1 i i+1 At x & ^: i i+1 At > & o:

The boundaries always run through q grid points (x).

Data Types
type	user_change_diff_cs

Functions/Subroutines
subroutine, public	user_change_diff (h, tv, G, CS, Kd, Kd_int, T_f, S_f, Kd_int_add)
	This subroutine provides an interface for a user to use to modify the main code to alter the diffusivities as needed. The specific example implemented here augments the diffusivity for a specified range of latitude and coordinate potential density. More...

logical function	range_ok (range)
	This subroutine checks whether the 4 values of range are in ascending order. More...

real function	val_weights (val, range)
	This subroutine returns a value that goes smoothly from 0 to 1, stays at 1, and then goes smoothly back to 0 at the four values of range. The transitions are cubic, and have zero first derivatives where the curves hit 0 and 1. The values in range must be in ascending order, as can be checked by calling range_OK. More...

subroutine, public	user_change_diff_init (Time, G, param_file, diag, CS)
	Set up the module control structure. More...

subroutine, public	user_change_diff_end (CS)
	Clean up the module control structure. More...

Function/Subroutine Documentation

◆ range_ok()

logical function user_change_diffusivity::range_ok ( real, dimension(4), intent(in) range )

private

This subroutine checks whether the 4 values of range are in ascending order.

Parameters

[in] range Four values to check.

Returns: Return value.

Definition at line 160 of file user_change_diffusivity.F90.

Referenced by user_change_diff(), and user_change_diff_init().

   real, dimension(4), intent(in) :: range  !< Four values to check.
   logical                        :: ok     !< Return value.
 
   ok = ((range(1) <= range(2)) .and. (range(2) <= range(3)) .and. &
         (range(3) <= range(4)))

Here is the caller graph for this function:

◆ user_change_diff()

subroutine, public user_change_diffusivity::user_change_diff	(	real, dimension(szi_(g),szj_(g),szk_(g)), intent(in)	h,
		type(thermo_var_ptrs), intent(in)	tv,
		type(ocean_grid_type), intent(in)	G,
		type(user_change_diff_cs), pointer	CS,
		real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(inout), optional	Kd,
		real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke+1), intent(inout), optional	Kd_int,
		real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in), optional	T_f,
		real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in), optional	S_f,
		real, dimension(:,:,:), optional, pointer	Kd_int_add
	)

This subroutine provides an interface for a user to use to modify the main code to alter the diffusivities as needed. The specific example implemented here augments the diffusivity for a specified range of latitude and coordinate potential density.

Parameters

[in]	g	The ocean's grid structure.
[in]	h	Layer thickness, in m or kg m-2.
[in]	tv	A structure containing pointers to any available thermodynamic fields. Absent fields have NULL ptrs.
	cs	This module's control structure.
[in,out]	kd	The diapycnal diffusivity of each layer in m2 s-1.
[in,out]	kd_int	The diapycnal diffusivity at each interface in m2 s-1.
[in]	t_f	Temperature with massless layers filled in vertically.
[in]	s_f	Salinity with massless layers filled in vertically.
	kd_int_add	The diapycnal diffusivity that is being added at each interface in m2 s-1.

Definition at line 57 of file user_change_diffusivity.F90.

References mom_error_handler::mom_error(), range_ok(), and val_weights().

   type(ocean_grid_type),                    intent(in)    :: g   !< The ocean's grid structure.
   real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in)    :: h   !< Layer thickness, in m or kg m-2.
   type(thermo_var_ptrs),                    intent(in)    :: tv  !< A structure containing pointers
                                                                  !! to any available thermodynamic
                                                                  !! fields. Absent fields have NULL ptrs.
   type(user_change_diff_cs),                pointer       :: cs  !< This module's control structure.
   real, dimension(SZI_(G),SZJ_(G),SZK_(G)),   optional, intent(inout) :: kd !< The diapycnal diffusivity of
                                                                   !! each layer in m2 s-1.
   real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1), optional, intent(inout) :: kd_int !< The diapycnal diffusivity
                                                                   !! at each interface in m2 s-1.
   real, dimension(SZI_(G),SZJ_(G),SZK_(G)),   optional, intent(in)    :: t_f !< Temperature with massless
                                                                   !! layers filled in vertically.
   real, dimension(SZI_(G),SZJ_(G),SZK_(G)),   optional, intent(in)    :: s_f !< Salinity with massless
                                                                   !! layers filled in vertically.
   real, dimension(:,:,:),                     optional, pointer       :: kd_int_add !< The diapycnal
                                                                   !! diffusivity that is being added at
                                                                   !! each interface in m2 s-1.
 
   real :: rcv(szi_(g),szk_(g)) ! The coordinate density in layers in kg m-3.
   real :: p_ref(szi_(g))       ! An array of tv%P_Ref pressures.
   real :: rho_fn      ! The density dependence of the input function, 0-1, ND.
   real :: lat_fn      ! The latitude dependence of the input function, 0-1, ND.
   logical :: use_eos  ! If true, density is calculated from T & S using an
                       ! equation of state.
   logical :: store_kd_add  ! Save the added diffusivity as a diagnostic if true.
   integer :: i, j, k, is, ie, js, je, nz
   integer :: isd, ied, jsd, jed
 
   real :: kappa_fill  ! diffusivity used to fill massless layers
   real :: dt_fill     ! timestep used to fill massless layers
   character(len=200) :: mesg
 
   is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
   isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
 
   if (.not.associated(cs)) call mom_error(fatal,"user_set_diffusivity: "//&
          "Module must be initialized before it is used.")
 
   use_eos = associated(tv%eqn_of_state)
   if (.not.use_eos) return
   store_kd_add = .false.
   if (present(kd_int_add)) store_kd_add = associated(kd_int_add)
 
   if (.not.range_ok(cs%lat_range)) then
     write(mesg, '(4(1pe15.6))') cs%lat_range(1:4)
     call mom_error(fatal, "user_set_diffusivity: bad latitude range: \n  "//&
                     trim(mesg))
   endif
   if (.not.range_ok(cs%rho_range)) then
     write(mesg, '(4(1pe15.6))') cs%rho_range(1:4)
     call mom_error(fatal, "user_set_diffusivity: bad density range: \n  "//&
                     trim(mesg))
   endif
 
   if (store_kd_add) kd_int_add(:,:,:) = 0.0
 
   do i=is,ie ; p_ref(i) = tv%P_Ref ; enddo
   do j=js,je
     if (present(t_f) .and. present(s_f)) then
       do k=1,nz
         call calculate_density(t_f(:,j,k),s_f(:,j,k),p_ref,rcv(:,k),&
                                is,ie-is+1,tv%eqn_of_state)
       enddo
     else
       do k=1,nz
         call calculate_density(tv%T(:,j,k),tv%S(:,j,k),p_ref,rcv(:,k),&
                                is,ie-is+1,tv%eqn_of_state)
       enddo
     endif
 
     if (present(kd)) then
       do k=1,nz ; do i=is,ie
         if (cs%use_abs_lat) then
           lat_fn = val_weights(abs(g%geoLatT(i,j)), cs%lat_range)
         else
           lat_fn = val_weights(g%geoLatT(i,j), cs%lat_range)
         endif
         rho_fn = val_weights(rcv(i,k), cs%rho_range)
         if (rho_fn * lat_fn > 0.0) &
           kd(i,j,k) = kd(i,j,k) + cs%Kd_add * rho_fn * lat_fn
       enddo ; enddo
     endif
     if (present(kd_int)) then
       do k=2,nz ; do i=is,ie
         if (cs%use_abs_lat) then
           lat_fn = val_weights(abs(g%geoLatT(i,j)), cs%lat_range)
         else
           lat_fn = val_weights(g%geoLatT(i,j), cs%lat_range)
         endif
         !   rho_int = 0.5*(Rcv(i,k-1) + Rcv(i,k))
         rho_fn = val_weights( 0.5*(rcv(i,k-1) + rcv(i,k)), cs%rho_range)
         if (rho_fn * lat_fn > 0.0) then
           kd_int(i,j,k) = kd_int(i,j,k) + cs%Kd_add * rho_fn * lat_fn
           if (store_kd_add) kd_int_add(i,j,k) = cs%Kd_add * rho_fn * lat_fn
         endif
       enddo ; enddo
     endif
   enddo
 

Here is the call graph for this function:

◆ user_change_diff_end()

subroutine, public user_change_diffusivity::user_change_diff_end ( type(user_change_diff_cs), pointer CS )

Clean up the module control structure.

Parameters

cs	A pointer that is set to point to the control structure for this module.

Definition at line 269 of file user_change_diffusivity.F90.

   type(user_change_diff_cs), pointer :: cs         !< A pointer that is set to
                                                    !! point to the control
                                                    !! structure for this module.
 
   if (associated(cs)) deallocate(cs)
 

◆ user_change_diff_init()

subroutine, public user_change_diffusivity::user_change_diff_init	(	type(time_type), intent(in)	Time,
		type(ocean_grid_type), intent(in)	G,
		type(param_file_type), intent(in)	param_file,
		type(diag_ctrl), intent(inout), target	diag,
		type(user_change_diff_cs), pointer	CS
	)

Set up the module control structure.

Parameters

[in]	time	The current model time.
[in]	g	The ocean's grid structure.
[in]	param_file	A structure indicating the open file to parse for model parameter values.
[in,out]	diag	A structure that is used to regulate diagnostic output.
	cs	A pointer that is set to point to the control structure for this module.

Definition at line 200 of file user_change_diffusivity.F90.

References mom_error_handler::mom_error(), and range_ok().

   type(time_type),           intent(in)    :: time       !< The current model time.
   type(ocean_grid_type),     intent(in)    :: g          !< The ocean's grid structure.
   type(param_file_type),     intent(in)    :: param_file !< A structure indicating the
                                                          !! open file to parse for
                                                          !! model parameter values.
   type(diag_ctrl), target,   intent(inout) :: diag       !< A structure that is used to
                                                          !! regulate diagnostic output.
   type(user_change_diff_cs), pointer       :: cs         !< A pointer that is set to
                                                          !! point to the control
                                                          !! structure for this module.
 
 ! This include declares and sets the variable "version".
 #include "version_variable.h"
   character(len=40)  :: mdl = "user_set_diffusivity"  ! This module's name.
   character(len=200) :: mesg
   integer :: i, j, is, ie, js, je
 
   if (associated(cs)) then
     call mom_error(warning, "diabatic_entrain_init called with an associated "// &
                             "control structure.")
     return
   endif
   allocate(cs)
 
   is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
 
   cs%diag => diag
 
   ! Read all relevant parameters and write them to the model log.
   call log_version(param_file, mdl, version, "")
   call get_param(param_file, mdl, "USER_KD_ADD", cs%Kd_add, &
                  "A user-specified additional diffusivity over a range of \n"//&
                  "latitude and density.", units="m2 s-1", default=0.0)
   if (cs%Kd_add /= 0.0) then
     call get_param(param_file, mdl, "USER_KD_ADD_LAT_RANGE", cs%lat_range(:), &
                  "Four successive values that define a range of latitudes \n"//&
                  "over which the user-specified extra diffusivity is \n"//&
                  "applied.  The four values specify the latitudes at \n"//&
                  "which the extra diffusivity starts to increase from 0, \n"//&
                  "hits its full value, starts to decrease again, and is \n"//&
                  "back to 0.", units="degree", default=-1.0e9)
     call get_param(param_file, mdl, "USER_KD_ADD_RHO_RANGE", cs%rho_range(:), &
                  "Four successive values that define a range of potential \n"//&
                  "densities over which the user-given extra diffusivity \n"//&
                  "is applied.  The four values specify the density at \n"//&
                  "which the extra diffusivity starts to increase from 0, \n"//&
                  "hits its full value, starts to decrease again, and is \n"//&
                  "back to 0.", units="kg m-3", default=-1.0e9)
     call get_param(param_file, mdl, "USER_KD_ADD_USE_ABS_LAT", cs%use_abs_lat, &
                  "If true, use the absolute value of latitude when \n"//&
                  "checking whether a point fits into range of latitudes.", &
                  default=.false.)
   endif
 
  if (.not.range_ok(cs%lat_range)) then
     write(mesg, '(4(1pe15.6))') cs%lat_range(1:4)
     call mom_error(fatal, "user_set_diffusivity: bad latitude range: \n  "//&
                     trim(mesg))
   endif
   if (.not.range_ok(cs%rho_range)) then
     write(mesg, '(4(1pe15.6))') cs%rho_range(1:4)
     call mom_error(fatal, "user_set_diffusivity: bad density range: \n  "//&
                     trim(mesg))
   endif
 

Here is the call graph for this function:

◆ val_weights()

real function user_change_diffusivity::val_weights	(	real, intent(in)	val,
		real, dimension(4), intent(in)	range
	)

private

This subroutine returns a value that goes smoothly from 0 to 1, stays at 1, and then goes smoothly back to 0 at the four values of range. The transitions are cubic, and have zero first derivatives where the curves hit 0 and 1. The values in range must be in ascending order, as can be checked by calling range_OK.

Parameters

[in]	val	Value for which we need an answer.
[in]	range	Range over which the answer is non-zero.

Returns: Return value.

Definition at line 175 of file user_change_diffusivity.F90.

Referenced by user_change_diff().

   real,               intent(in) :: val    !< Value for which we need an answer.
   real, dimension(4), intent(in) :: range  !< Range over which the answer is non-zero.
   real                           :: ans    !< Return value.
 
   real :: x   ! A nondimensional number between 0 and 1.
 
   ans = 0.0
   if ((val > range(1)) .and. (val < range(4))) then
     if (val < range(2)) then
       ! x goes from 0 to 1; ans goes from 0 to 1, with 0 derivatives at the ends.
       x = (val - range(1)) / (range(2) - range(1))
       ans = x**2 * (3.0 - 2.0 * x)
     elseif (val > range(3)) then
       ! x goes from 0 to 1; ans goes from 0 to 1, with 0 derivatives at the ends.
       x = (range(4) - val) / (range(4) - range(3))
       ans = x**2 * (3.0 - 2.0 * x)
     else
       ans = 1.0
     endif
   endif
 

Here is the caller graph for this function:

Detailed Description

Data Types

Functions/Subroutines

Function/Subroutine Documentation

◆ range_ok()

◆ user_change_diff()

◆ user_change_diff_end()

◆ user_change_diff_init()

◆ val_weights()