MOM6
user_change_diffusivity Module Reference

Detailed Description

Increments the diapycnal diffusivity in a specified band of latitudes and densities.

Data Types

type  user_change_diff_cs
 Control structure for user_change_diffusivity. More...
 

Functions/Subroutines

subroutine, public user_change_diff (h, tv, G, GV, CS, Kd_lay, 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, GV, US, 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]rangeFour values to check.
Returns
Return value.

Definition at line 152 of file user_change_diffusivity.F90.

152  real, dimension(4), intent(in) :: range !< Four values to check.
153  logical :: OK !< Return value.
154 
155  ok = ((range(1) <= range(2)) .and. (range(2) <= range(3)) .and. &
156  (range(3) <= range(4)))
157 

Referenced by user_change_diff(), and user_change_diff_init().

Here is the caller graph for this function:

◆ user_change_diff()

subroutine, public user_change_diffusivity::user_change_diff ( real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
type(thermo_var_ptrs), intent(in)  tv,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
type(user_change_diff_cs), pointer  CS,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(inout), optional  Kd_lay,
real, dimension(szi_(g),szj_(g),szk_(g)+1), intent(inout), optional  Kd_int,
real, dimension(szi_(g),szj_(g),szk_(g)), intent(in), optional  T_f,
real, dimension(szi_(g),szj_(g),szk_(g)), 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]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure
[in]hLayer thickness [H ~> m or kg m-2].
[in]tvA structure containing pointers to any available thermodynamic fields. Absent fields have NULL ptrs.
csThis module's control structure.
[in,out]kd_layThe diapycnal diffusivity of each layer [Z2 T-1 ~> m2 s-1].
[in,out]kd_intThe diapycnal diffusivity at each interface [Z2 T-1 ~> m2 s-1].
[in]t_fTemperature with massless layers filled in vertically [degC].
[in]s_fSalinity with massless layers filled in vertically [ppt].
kd_int_addThe diapycnal diffusivity that is being added at each interface [Z2 T-1 ~> m2 s-1].

Definition at line 48 of file user_change_diffusivity.F90.

48  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
49  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
50  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thickness [H ~> m or kg m-2].
51  type(thermo_var_ptrs), intent(in) :: tv !< A structure containing pointers
52  !! to any available thermodynamic
53  !! fields. Absent fields have NULL ptrs.
54  type(user_change_diff_CS), pointer :: CS !< This module's control structure.
55  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), optional, intent(inout) :: Kd_lay !< The diapycnal diffusivity of
56  !! each layer [Z2 T-1 ~> m2 s-1].
57  real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1), optional, intent(inout) :: Kd_int !< The diapycnal diffusivity
58  !! at each interface [Z2 T-1 ~> m2 s-1].
59  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), optional, intent(in) :: T_f !< Temperature with massless
60  !! layers filled in vertically [degC].
61  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), optional, intent(in) :: S_f !< Salinity with massless
62  !! layers filled in vertically [ppt].
63  real, dimension(:,:,:), optional, pointer :: Kd_int_add !< The diapycnal
64  !! diffusivity that is being added at
65  !! each interface [Z2 T-1 ~> m2 s-1].
66  ! Local variables
67  real :: Rcv(SZI_(G),SZK_(G)) ! The coordinate density in layers [kg m-3].
68  real :: p_ref(SZI_(G)) ! An array of tv%P_Ref pressures.
69  real :: rho_fn ! The density dependence of the input function, 0-1 [nondim].
70  real :: lat_fn ! The latitude dependence of the input function, 0-1 [nondim].
71  logical :: use_EOS ! If true, density is calculated from T & S using an
72  ! equation of state.
73  logical :: store_Kd_add ! Save the added diffusivity as a diagnostic if true.
74  integer :: i, j, k, is, ie, js, je, nz
75  integer :: isd, ied, jsd, jed
76 
77  real :: kappa_fill ! diffusivity used to fill massless layers
78  real :: dt_fill ! timestep used to fill massless layers
79  character(len=200) :: mesg
80 
81  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
82  isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
83 
84  if (.not.associated(cs)) call mom_error(fatal,"user_set_diffusivity: "//&
85  "Module must be initialized before it is used.")
86 
87  use_eos = associated(tv%eqn_of_state)
88  if (.not.use_eos) return
89  store_kd_add = .false.
90  if (present(kd_int_add)) store_kd_add = associated(kd_int_add)
91 
92  if (.not.range_ok(cs%lat_range)) then
93  write(mesg, '(4(1pe15.6))') cs%lat_range(1:4)
94  call mom_error(fatal, "user_set_diffusivity: bad latitude range: \n "//&
95  trim(mesg))
96  endif
97  if (.not.range_ok(cs%rho_range)) then
98  write(mesg, '(4(1pe15.6))') cs%rho_range(1:4)
99  call mom_error(fatal, "user_set_diffusivity: bad density range: \n "//&
100  trim(mesg))
101  endif
102 
103  if (store_kd_add) kd_int_add(:,:,:) = 0.0
104 
105  do i=is,ie ; p_ref(i) = tv%P_Ref ; enddo
106  do j=js,je
107  if (present(t_f) .and. present(s_f)) then
108  do k=1,nz
109  call calculate_density(t_f(:,j,k),s_f(:,j,k),p_ref,rcv(:,k),&
110  is,ie-is+1,tv%eqn_of_state)
111  enddo
112  else
113  do k=1,nz
114  call calculate_density(tv%T(:,j,k),tv%S(:,j,k),p_ref,rcv(:,k),&
115  is,ie-is+1,tv%eqn_of_state)
116  enddo
117  endif
118 
119  if (present(kd_lay)) then
120  do k=1,nz ; do i=is,ie
121  if (cs%use_abs_lat) then
122  lat_fn = val_weights(abs(g%geoLatT(i,j)), cs%lat_range)
123  else
124  lat_fn = val_weights(g%geoLatT(i,j), cs%lat_range)
125  endif
126  rho_fn = val_weights(rcv(i,k), cs%rho_range)
127  if (rho_fn * lat_fn > 0.0) &
128  kd_lay(i,j,k) = kd_lay(i,j,k) + cs%Kd_add * rho_fn * lat_fn
129  enddo ; enddo
130  endif
131  if (present(kd_int)) then
132  do k=2,nz ; do i=is,ie
133  if (cs%use_abs_lat) then
134  lat_fn = val_weights(abs(g%geoLatT(i,j)), cs%lat_range)
135  else
136  lat_fn = val_weights(g%geoLatT(i,j), cs%lat_range)
137  endif
138  ! rho_int = 0.5*(Rcv(i,k-1) + Rcv(i,k))
139  rho_fn = val_weights( 0.5*(rcv(i,k-1) + rcv(i,k)), cs%rho_range)
140  if (rho_fn * lat_fn > 0.0) then
141  kd_int(i,j,k) = kd_int(i,j,k) + cs%Kd_add * rho_fn * lat_fn
142  if (store_kd_add) kd_int_add(i,j,k) = cs%Kd_add * rho_fn * lat_fn
143  endif
144  enddo ; enddo
145  endif
146  enddo
147 

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

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
csA pointer that is set to point to the control structure for this module.

Definition at line 263 of file user_change_diffusivity.F90.

263  type(user_change_diff_CS), pointer :: CS !< A pointer that is set to point to the control
264  !! structure for this module.
265 
266  if (associated(cs)) deallocate(cs)
267 

Referenced by mom_set_diffusivity::set_diffusivity_end().

Here is the caller graph for this function:

◆ 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(verticalgrid_type), intent(in)  GV,
type(unit_scale_type), intent(in)  US,
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]timeThe current model time.
[in]gThe ocean's grid structure.
[in]gvThe ocean's vertical grid structure
[in]usA dimensional unit scaling type
[in]param_fileA structure indicating the open file to parse for model parameter values.
[in,out]diagA structure that is used to regulate diagnostic output.
csA pointer that is set to point to the control structure for this module.

Definition at line 191 of file user_change_diffusivity.F90.

191  type(time_type), intent(in) :: Time !< The current model time.
192  type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure.
193  type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure
194  type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
195  type(param_file_type), intent(in) :: param_file !< A structure indicating the
196  !! open file to parse for
197  !! model parameter values.
198  type(diag_ctrl), target, intent(inout) :: diag !< A structure that is used to
199  !! regulate diagnostic output.
200  type(user_change_diff_CS), pointer :: CS !< A pointer that is set to
201  !! point to the control
202  !! structure for this module.
203 
204 ! This include declares and sets the variable "version".
205 #include "version_variable.h"
206  character(len=40) :: mdl = "user_set_diffusivity" ! This module's name.
207  character(len=200) :: mesg
208  integer :: i, j, is, ie, js, je
209 
210  if (associated(cs)) then
211  call mom_error(warning, "diabatic_entrain_init called with an associated "// &
212  "control structure.")
213  return
214  endif
215  allocate(cs)
216 
217  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
218 
219  cs%diag => diag
220 
221  ! Read all relevant parameters and write them to the model log.
222  call log_version(param_file, mdl, version, "")
223  call get_param(param_file, mdl, "USER_KD_ADD", cs%Kd_add, &
224  "A user-specified additional diffusivity over a range of "//&
225  "latitude and density.", default=0.0, units="m2 s-1", &
226  scale=us%m2_s_to_Z2_T)
227  if (cs%Kd_add /= 0.0) then
228  call get_param(param_file, mdl, "USER_KD_ADD_LAT_RANGE", cs%lat_range(:), &
229  "Four successive values that define a range of latitudes "//&
230  "over which the user-specified extra diffusivity is "//&
231  "applied. The four values specify the latitudes at "//&
232  "which the extra diffusivity starts to increase from 0, "//&
233  "hits its full value, starts to decrease again, and is "//&
234  "back to 0.", units="degree", default=-1.0e9)
235  call get_param(param_file, mdl, "USER_KD_ADD_RHO_RANGE", cs%rho_range(:), &
236  "Four successive values that define a range of potential "//&
237  "densities over which the user-given extra diffusivity "//&
238  "is applied. The four values specify the density at "//&
239  "which the extra diffusivity starts to increase from 0, "//&
240  "hits its full value, starts to decrease again, and is "//&
241  "back to 0.", units="kg m-3", default=-1.0e9)
242  call get_param(param_file, mdl, "USER_KD_ADD_USE_ABS_LAT", cs%use_abs_lat, &
243  "If true, use the absolute value of latitude when "//&
244  "checking whether a point fits into range of latitudes.", &
245  default=.false.)
246  endif
247 
248  if (.not.range_ok(cs%lat_range)) then
249  write(mesg, '(4(1pe15.6))') cs%lat_range(1:4)
250  call mom_error(fatal, "user_set_diffusivity: bad latitude range: \n "//&
251  trim(mesg))
252  endif
253  if (.not.range_ok(cs%rho_range)) then
254  write(mesg, '(4(1pe15.6))') cs%rho_range(1:4)
255  call mom_error(fatal, "user_set_diffusivity: bad density range: \n "//&
256  trim(mesg))
257  endif
258 

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

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]valValue for which we need an answer.
[in]rangeRange over which the answer is non-zero.
Returns
Return value.

Definition at line 166 of file user_change_diffusivity.F90.

166  real, intent(in) :: val !< Value for which we need an answer.
167  real, dimension(4), intent(in) :: range !< Range over which the answer is non-zero.
168  real :: ans !< Return value.
169  ! Local variables
170  real :: x ! A nondimensional number between 0 and 1.
171 
172  ans = 0.0
173  if ((val > range(1)) .and. (val < range(4))) then
174  if (val < range(2)) then
175  ! x goes from 0 to 1; ans goes from 0 to 1, with 0 derivatives at the ends.
176  x = (val - range(1)) / (range(2) - range(1))
177  ans = x**2 * (3.0 - 2.0 * x)
178  elseif (val > range(3)) then
179  ! x goes from 0 to 1; ans goes from 0 to 1, with 0 derivatives at the ends.
180  x = (range(4) - val) / (range(4) - range(3))
181  ans = x**2 * (3.0 - 2.0 * x)
182  else
183  ans = 1.0
184  endif
185  endif
186 

Referenced by user_change_diff().

Here is the caller graph for this function: