MOM6
Functions/Subroutines
mom_isopycnal_slopes Module Reference

Detailed Description

Calculations of isoneutral slopes and stratification.

Functions/Subroutines

subroutine, public calc_isoneutral_slopes (G, GV, h, e, tv, dt_kappa_smooth, slope_x, slope_y, N2_u, N2_v, halo)
 
subroutine vert_fill_ts (h, T_in, S_in, kappa, dt, T_f, S_f, G, GV, halo_here)
 

Function/Subroutine Documentation

◆ calc_isoneutral_slopes()

subroutine, public mom_isopycnal_slopes::calc_isoneutral_slopes ( type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke+1), intent(in)  e,
type(thermo_var_ptrs), intent(in)  tv,
real, intent(in)  dt_kappa_smooth,
real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke+1), intent(inout)  slope_x,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke+1), intent(inout)  slope_y,
real, dimension( g %isdb: g %iedb, g %jsd: g %jed, g %ke+1), intent(inout), optional  N2_u,
real, dimension( g %isd: g %ied, g %jsdb: g %jedb, g %ke+1), intent(inout), optional  N2_v,
integer, intent(in), optional  halo 
)
Parameters
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[in]hLayer thicknesses, in H (usually m or kg m-2)
[in]tvA structure pointing to various thermodynamic variables

Definition at line 21 of file MOM_isopycnal_slopes.F90.

References mom_eos::calculate_density_derivs(), and vert_fill_ts().

Referenced by mom_lateral_mixing_coeffs::calc_slope_functions().

21  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
22  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
23  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses, in H (usually m or kg m-2)
24  real, dimension(SZI_(G),SZJ_(G),SZK_(G)+1), intent(in) :: e
25  type(thermo_var_ptrs), intent(in) :: tv !< A structure pointing to various thermodynamic variables
26  real, intent(in) :: dt_kappa_smooth
27  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)+1), intent(inout) :: slope_x
28  real, dimension(SZI_(G),SZJB_(G),SZK_(G)+1), intent(inout) :: slope_y
29  real, dimension(SZIB_(G),SZJ_(G),SZK_(G)+1), intent(inout) :: n2_u
30  real, dimension(SZI_(G),SZJB_(G),SZK_(G)+1), intent(inout) :: n2_v
31  optional :: n2_u, n2_v
32  integer, optional, intent(in) :: halo
33  ! Local variables
34  real, dimension(SZI_(G), SZJ_(G), SZK_(G)) :: &
35  t, & ! The temperature (or density) in C, with the values in
36  ! in massless layers filled vertically by diffusion.
37  s, & ! The filled salinity, in PSU, with the values in
38  ! in massless layers filled vertically by diffusion.
39  rho ! Density itself, when a nonlinear equation of state is
40  ! not in use.
41  real, dimension(SZI_(G), SZJ_(G), SZK_(G)+1) :: &
42  pres ! The pressure at an interface, in Pa.
43  real, dimension(SZIB_(G)) :: &
44  drho_dt_u, & ! The derivatives of density with temperature and
45  drho_ds_u ! salinity at u points, in kg m-3 K-1 and kg m-3 psu-1.
46  real, dimension(SZI_(G)) :: &
47  drho_dt_v, & ! The derivatives of density with temperature and
48  drho_ds_v ! salinity at v points, in kg m-3 K-1 and kg m-3 psu-1.
49  real, dimension(SZIB_(G)) :: &
50  t_u, s_u, & ! Temperature, salinity, and pressure on the interface at
51  pres_u ! the u-point in the horizontal.
52  real, dimension(SZI_(G)) :: &
53  t_v, s_v, & ! Temperature, salinity, and pressure on the interface at
54  pres_v ! the v-point in the horizontal.
55  real :: drdia, drdib ! Along layer zonal- and meridional- potential density
56  real :: drdja, drdjb ! gradients in the layers above (A) and below(B) the
57  ! interface times the grid spacing, in kg m-3.
58  real :: drdkl, drdkr ! Vertical density differences across an interface,
59  ! in kg m-3.
60  real :: hg2a, hg2b, hg2l, hg2r
61  real :: haa, hab, hal, har
62  real :: dzal, dzar
63  real :: wta, wtb, wtl, wtr
64  real :: drdx, drdy, drdz ! Zonal, meridional, and vertical density gradients,
65  ! in units of kg m-4.
66  real :: slope ! The slope of density surfaces, calculated in a way
67  ! that is always between -1 and 1.
68  real :: mag_grad2 ! The squared magnitude of the 3-d density gradient, in kg2 m-8.
69  real :: slope2_ratio ! The ratio of the slope squared to slope_max squared.
70  real :: h_neglect ! A thickness that is so small it is usually lost
71  ! in roundoff and can be neglected, in H.
72  real :: h_neglect2 ! h_neglect^2, in H2.
73  real :: dz_neglect ! A thickness in m that is so small it is usually lost
74  ! in roundoff and can be neglected, in m.
75  logical :: use_eos ! If true, density is calculated from T & S using an
76  ! equation of state.
77  real :: g_rho0, n2, dzn2, h_x(szib_(g)), h_y(szi_(g))
78 
79  logical :: present_n2_u, present_n2_v
80  integer :: is, ie, js, je, nz, isdb
81  integer :: i, j, k
82 
83  if (present(halo)) then
84  is = g%isc-halo ; ie = g%iec+halo ; js = g%jsc-halo ; je = g%jec+halo
85  else
86  is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
87  endif
88  nz = g%ke ; isdb = g%IsdB
89 
90  h_neglect = gv%H_subroundoff ; h_neglect2 = h_neglect**2
91  dz_neglect = gv%H_subroundoff*gv%H_to_m
92 
93  use_eos = associated(tv%eqn_of_state)
94 
95  present_n2_u = PRESENT(n2_u)
96  present_n2_v = PRESENT(n2_v)
97  g_rho0 = gv%g_Earth / gv%Rho0
98  if (present_n2_u) then
99  do j=js,je ; do i=is-1,ie
100  n2_u(i,j,1) = 0.
101  n2_u(i,j,nz+1) = 0.
102  enddo ; enddo
103  endif
104  if (present_n2_v) then
105  do j=js-1,je ; do i=is,ie
106  n2_v(i,j,1) = 0.
107  n2_v(i,j,nz+1) = 0.
108  enddo ; enddo
109  endif
110 
111  if (use_eos) then
112  if (present(halo)) then
113  call vert_fill_ts(h, tv%T, tv%S, dt_kappa_smooth, 1.0, t, s, g, gv, halo+1)
114  else
115  call vert_fill_ts(h, tv%T, tv%S, dt_kappa_smooth, 1.0, t, s, g, gv, 1)
116  endif
117  endif
118 
119  ! Find the maximum and minimum permitted streamfunction.
120 !$OMP parallel default(none) shared(is,ie,js,je,pres,GV,h,nz)
121 !$OMP do
122  do j=js-1,je+1 ; do i=is-1,ie+1
123  pres(i,j,1) = 0.0 ! ### This should be atmospheric pressure.
124  pres(i,j,2) = pres(i,j,1) + gv%H_to_Pa*h(i,j,1)
125  enddo ; enddo
126 !$OMP do
127  do j=js-1,je+1
128  do k=2,nz ; do i=is-1,ie+1
129  pres(i,j,k+1) = pres(i,j,k) + gv%H_to_Pa*h(i,j,k)
130  enddo ; enddo
131  enddo
132 !$OMP end parallel
133 
134 !$OMP parallel do default(none) shared(nz,is,ie,js,je,use_EOS,G,GV,pres,T,S, &
135 !$OMP IsdB,tv,h,h_neglect,e,dz_neglect, &
136 !$OMP h_neglect2,present_N2_u,G_Rho0,N2_u,slope_x) &
137 !$OMP private(drdiA,drdiB,drdkL,drdkR,pres_u,T_u,S_u, &
138 !$OMP drho_dT_u,drho_dS_u,hg2A,hg2B,hg2L,hg2R,haA, &
139 !$OMP haB,haL,haR,dzaL,dzaR,wtA,wtB,wtL,wtR,drdz, &
140 !$OMP drdx,mag_grad2,Slope,slope2_Ratio)
141  do j=js,je ; do k=nz,2,-1
142  if (.not.(use_eos)) then
143  drdia = 0.0 ; drdib = 0.0
144  drdkl = gv%Rlay(k)-gv%Rlay(k-1) ; drdkr = gv%Rlay(k)-gv%Rlay(k-1)
145  endif
146 
147  ! Calculate the zonal isopycnal slope.
148  if (use_eos) then
149  do i=is-1,ie
150  pres_u(i) = 0.5*(pres(i,j,k) + pres(i+1,j,k))
151  t_u(i) = 0.25*((t(i,j,k) + t(i+1,j,k)) + (t(i,j,k-1) + t(i+1,j,k-1)))
152  s_u(i) = 0.25*((s(i,j,k) + s(i+1,j,k)) + (s(i,j,k-1) + s(i+1,j,k-1)))
153  enddo
154  call calculate_density_derivs(t_u, s_u, pres_u, drho_dt_u, &
155  drho_ds_u, (is-isdb+1)-1, ie-is+2, tv%eqn_of_state)
156  endif
157 
158  do i=is-1,ie
159  if (use_eos) then
160  ! Estimate the horizontal density gradients along layers.
161  drdia = drho_dt_u(i) * (t(i+1,j,k-1)-t(i,j,k-1)) + &
162  drho_ds_u(i) * (s(i+1,j,k-1)-s(i,j,k-1))
163  drdib = drho_dt_u(i) * (t(i+1,j,k)-t(i,j,k)) + &
164  drho_ds_u(i) * (s(i+1,j,k)-s(i,j,k))
165 
166  ! Estimate the vertical density gradients times the grid spacing.
167  drdkl = (drho_dt_u(i) * (t(i,j,k)-t(i,j,k-1)) + &
168  drho_ds_u(i) * (s(i,j,k)-s(i,j,k-1)))
169  drdkr = (drho_dt_u(i) * (t(i+1,j,k)-t(i+1,j,k-1)) + &
170  drho_ds_u(i) * (s(i+1,j,k)-s(i+1,j,k-1)))
171  endif
172 
173 
174  if (use_eos) then
175  hg2a = h(i,j,k-1)*h(i+1,j,k-1) + h_neglect2
176  hg2b = h(i,j,k)*h(i+1,j,k) + h_neglect2
177  hg2l = h(i,j,k-1)*h(i,j,k) + h_neglect2
178  hg2r = h(i+1,j,k-1)*h(i+1,j,k) + h_neglect2
179  haa = 0.5*(h(i,j,k-1) + h(i+1,j,k-1))
180  hab = 0.5*(h(i,j,k) + h(i+1,j,k)) + h_neglect
181  hal = 0.5*(h(i,j,k-1) + h(i,j,k)) + h_neglect
182  har = 0.5*(h(i+1,j,k-1) + h(i+1,j,k)) + h_neglect
183  if (gv%Boussinesq) then
184  dzal = hal * gv%H_to_m ; dzar = har * gv%H_to_m
185  else
186  dzal = 0.5*(e(i,j,k-1) - e(i,j,k+1)) + dz_neglect
187  dzar = 0.5*(e(i+1,j,k-1) - e(i+1,j,k+1)) + dz_neglect
188  endif
189  ! Use the harmonic mean thicknesses to weight the horizontal gradients.
190  ! These unnormalized weights have been rearranged to minimize divisions.
191  wta = hg2a*hab ; wtb = hg2b*haa
192  wtl = hg2l*(har*dzar) ; wtr = hg2r*(hal*dzal)
193 
194  drdz = (wtl * drdkl + wtr * drdkr) / (dzal*wtl + dzar*wtr)
195  ! The expression for drdz above is mathematically equivalent to:
196  ! drdz = ((hg2L/haL) * drdkL/dzaL + (hg2R/haR) * drdkR/dzaR) / &
197  ! ((hg2L/haL) + (hg2R/haR))
198  ! This is the gradient of density along geopotentials.
199  drdx = ((wta * drdia + wtb * drdib) / (wta + wtb) - &
200  drdz * (e(i,j,k)-e(i+1,j,k))) * g%IdxCu(i,j)
201 
202  ! This estimate of slope is accurate for small slopes, but bounded
203  ! to be between -1 and 1.
204  mag_grad2 = drdx**2 + drdz**2
205  if (mag_grad2 > 0.0) then
206  slope_x(i,j,k) = drdx / sqrt(mag_grad2)
207  else ! Just in case mag_grad2 = 0 ever.
208  slope_x(i,j,k) = 0.0
209  endif
210 
211  if (present_n2_u) n2_u(i,j,k) = g_rho0 * drdz * g%mask2dCu(i,j) ! Square of Brunt-Vaisala frequency (s-2)
212 
213  else ! With .not.use_EOS, the layers are constant density.
214  slope_x(i,j,k) = ((e(i,j,k)-e(i+1,j,k))*g%IdxCu(i,j)) * gv%m_to_H
215  endif
216 
217  enddo ! I
218  enddo ; enddo ! end of j-loop
219 
220  ! Calculate the meridional isopycnal slope.
221 !$OMP parallel do default(none) shared(nz,is,ie,js,je,use_EOS,G,GV,pres,T,S, &
222 !$OMP IsdB,tv,h,h_neglect,e,dz_neglect, &
223 !$OMP h_neglect2,present_N2_v,G_Rho0,N2_v,slope_y) &
224 !$OMP private(drdjA,drdjB,drdkL,drdkR,pres_v,T_v,S_v, &
225 !$OMP drho_dT_v,drho_dS_v,hg2A,hg2B,hg2L,hg2R,haA, &
226 !$OMP haB,haL,haR,dzaL,dzaR,wtA,wtB,wtL,wtR,drdz, &
227 !$OMP drdy,mag_grad2,Slope,slope2_Ratio)
228  do j=js-1,je ; do k=nz,2,-1
229  if (.not.(use_eos)) then
230  drdja = 0.0 ; drdjb = 0.0
231  drdkl = gv%Rlay(k)-gv%Rlay(k-1) ; drdkr = gv%Rlay(k)-gv%Rlay(k-1)
232  endif
233 
234  if (use_eos) then
235  do i=is,ie
236  pres_v(i) = 0.5*(pres(i,j,k) + pres(i,j+1,k))
237  t_v(i) = 0.25*((t(i,j,k) + t(i,j+1,k)) + (t(i,j,k-1) + t(i,j+1,k-1)))
238  s_v(i) = 0.25*((s(i,j,k) + s(i,j+1,k)) + (s(i,j,k-1) + s(i,j+1,k-1)))
239  enddo
240  call calculate_density_derivs(t_v, s_v, pres_v, drho_dt_v, &
241  drho_ds_v, is, ie-is+1, tv%eqn_of_state)
242  endif
243  do i=is,ie
244  if (use_eos) then
245  ! Estimate the horizontal density gradients along layers.
246  drdja = drho_dt_v(i) * (t(i,j+1,k-1)-t(i,j,k-1)) + &
247  drho_ds_v(i) * (s(i,j+1,k-1)-s(i,j,k-1))
248  drdjb = drho_dt_v(i) * (t(i,j+1,k)-t(i,j,k)) + &
249  drho_ds_v(i) * (s(i,j+1,k)-s(i,j,k))
250 
251  ! Estimate the vertical density gradients times the grid spacing.
252  drdkl = (drho_dt_v(i) * (t(i,j,k)-t(i,j,k-1)) + &
253  drho_ds_v(i) * (s(i,j,k)-s(i,j,k-1)))
254  drdkr = (drho_dt_v(i) * (t(i,j+1,k)-t(i,j+1,k-1)) + &
255  drho_ds_v(i) * (s(i,j+1,k)-s(i,j+1,k-1)))
256  endif
257 
258  if (use_eos) then
259  hg2a = h(i,j,k-1)*h(i,j+1,k-1) + h_neglect2
260  hg2b = h(i,j,k)*h(i,j+1,k) + h_neglect2
261  hg2l = h(i,j,k-1)*h(i,j,k) + h_neglect2
262  hg2r = h(i,j+1,k-1)*h(i,j+1,k) + h_neglect2
263  haa = 0.5*(h(i,j,k-1) + h(i,j+1,k-1)) + h_neglect
264  hab = 0.5*(h(i,j,k) + h(i,j+1,k)) + h_neglect
265  hal = 0.5*(h(i,j,k-1) + h(i,j,k)) + h_neglect
266  har = 0.5*(h(i,j+1,k-1) + h(i,j+1,k)) + h_neglect
267  if (gv%Boussinesq) then
268  dzal = hal * gv%H_to_m ; dzar = har * gv%H_to_m
269  else
270  dzal = 0.5*(e(i,j,k-1) - e(i,j,k+1)) + dz_neglect
271  dzar = 0.5*(e(i,j+1,k-1) - e(i,j+1,k+1)) + dz_neglect
272  endif
273  ! Use the harmonic mean thicknesses to weight the horizontal gradients.
274  ! These unnormalized weights have been rearranged to minimize divisions.
275  wta = hg2a*hab ; wtb = hg2b*haa
276  wtl = hg2l*(har*dzar) ; wtr = hg2r*(hal*dzal)
277 
278  drdz = (wtl * drdkl + wtr * drdkr) / (dzal*wtl + dzar*wtr)
279  ! The expression for drdz above is mathematically equivalent to:
280  ! drdz = ((hg2L/haL) * drdkL/dzaL + (hg2R/haR) * drdkR/dzaR) / &
281  ! ((hg2L/haL) + (hg2R/haR))
282  ! This is the gradient of density along geopotentials.
283  drdy = ((wta * drdja + wtb * drdjb) / (wta + wtb) - &
284  drdz * (e(i,j,k)-e(i,j+1,k))) * g%IdyCv(i,j)
285 
286  ! This estimate of slope is accurate for small slopes, but bounded
287  ! to be between -1 and 1.
288  mag_grad2 = drdy**2 + drdz**2
289  if (mag_grad2 > 0.0) then
290  slope_y(i,j,k) = drdy / sqrt(mag_grad2)
291  else ! Just in case mag_grad2 = 0 ever.
292  slope_y(i,j,k) = 0.0
293  endif
294 
295  if (present_n2_v) n2_v(i,j,k) = g_rho0 * drdz * g%mask2dCv(i,j) ! Square of Brunt-Vaisala frequency (s-2)
296 
297  else ! With .not.use_EOS, the layers are constant density.
298  slope_y(i,j,k) = ((e(i,j,k)-e(i,j+1,k))*g%IdyCv(i,j)) * gv%m_to_H
299  endif
300 
301  enddo ! i
302  enddo ; enddo ! end of j-loop
303 
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:19
Here is the call graph for this function:
Here is the caller graph for this function:

◆ vert_fill_ts()

subroutine mom_isopycnal_slopes::vert_fill_ts ( real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  h,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  T_in,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(in)  S_in,
real, intent(in)  kappa,
real, intent(in)  dt,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(out)  T_f,
real, dimension( g %isd: g %ied, g %jsd: g %jed, g %ke), intent(out)  S_f,
type(ocean_grid_type), intent(in)  G,
type(verticalgrid_type), intent(in)  GV,
integer, intent(in), optional  halo_here 
)
private
Parameters
[in]gThe ocean's grid structure
[in]gvThe ocean's vertical grid structure
[in]hLayer thicknesses, in H (usually m or kg m-2)
[in]dtThe time increment, in s.

Definition at line 307 of file MOM_isopycnal_slopes.F90.

Referenced by calc_isoneutral_slopes().

307  type(ocean_grid_type), intent(in) :: g !< The ocean's grid structure
308  type(verticalgrid_type), intent(in) :: gv !< The ocean's vertical grid structure
309  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: h !< Layer thicknesses, in H (usually m or kg m-2)
310  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: t_in
311  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(in) :: s_in
312  real, intent(in) :: kappa
313  real, intent(in) :: dt !< The time increment, in s.
314  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(out) :: t_f
315  real, dimension(SZI_(G),SZJ_(G),SZK_(G)), intent(out) :: s_f
316  integer, optional, intent(in) :: halo_here
317 ! This subroutine fills massless layers with sensible values of two
318 !* tracer arrays (nominally temperature and salinity) by diffusing
319 !* vertically with a (small?) constant diffusivity.
320 
321 ! Arguments: h - Layer thickness, in m or kg m-2.
322 ! (in) T_in - The input temperature, in K.
323 ! (in) S_in - The input salinity, in psu.
324 ! (in) kappa - The diapycnal diffusivity, in m2 s-1.
325 ! (in) dt - Time increment in s.
326 ! (out) T_f - The filled temperature, in K.
327 ! (out) S_f - The filled salinity, in psu.
328 ! (in) G - The ocean's grid structure.
329 ! (in) GV - The ocean's vertical grid structure.
330 ! (in,opt) halo_here - the number of halo points to work on, 0 by default.
331 
332  real :: ent(szi_(g),szk_(g)+1) ! The diffusive entrainment (kappa*dt)/dz
333  ! between layers in a timestep in m or kg m-2.
334  real :: b1(szi_(g)), d1(szi_(g)) ! b1, c1, and d1 are variables used by the
335  real :: c1(szi_(g),szk_(g)) ! tridiagonal solver.
336  real :: kap_dt_x2 ! The product of 2*kappa*dt, converted to
337  ! the same units as h, in m2 or kg2 m-4.
338  real :: h_neglect ! A negligible thickness, in m or kg m-2, to
339  ! allow for zero thicknesses.
340  integer :: i, j, k, is, ie, js, je, nz, halo
341 
342  halo=0 ; if (present(halo_here)) halo = max(halo_here,0)
343 
344  is = g%isc-halo ; ie = g%iec+halo ; js = g%jsc-halo ; je = g%jec+halo
345  nz = g%ke
346 
347  kap_dt_x2 = (2.0*kappa*dt)*gv%m_to_H**2
348  h_neglect = gv%H_subroundoff
349 
350  if (kap_dt_x2 <= 0.0) then
351 !$OMP parallel do default(none) shared(is,ie,js,je,nz,T_f,T_in,S_f,S_in)
352  do k=1,nz ; do j=js,je ; do i=is,ie
353  t_f(i,j,k) = t_in(i,j,k) ; s_f(i,j,k) = s_in(i,j,k)
354  enddo ; enddo ; enddo
355  else
356 !$OMP parallel do default(none) private(ent,b1,d1,c1) &
357 !$OMP shared(is,ie,js,je,nz,kap_dt_x2,h,h_neglect,T_f,S_f,T_in,S_in)
358  do j=js,je
359  do i=is,ie
360  ent(i,2) = kap_dt_x2 / ((h(i,j,1)+h(i,j,2)) + h_neglect)
361  b1(i) = 1.0 / (h(i,j,1)+ent(i,2))
362  d1(i) = b1(i) * h(i,j,1)
363  t_f(i,j,1) = (b1(i)*h(i,j,1))*t_in(i,j,1)
364  s_f(i,j,1) = (b1(i)*h(i,j,1))*s_in(i,j,1)
365  enddo
366  do k=2,nz-1 ; do i=is,ie
367  ent(i,k+1) = kap_dt_x2 / ((h(i,j,k)+h(i,j,k+1)) + h_neglect)
368  c1(i,k) = ent(i,k) * b1(i)
369  b1(i) = 1.0 / ((h(i,j,k) + d1(i)*ent(i,k)) + ent(i,k+1))
370  d1(i) = b1(i) * (h(i,j,k) + d1(i)*ent(i,k))
371  t_f(i,j,k) = b1(i) * (h(i,j,k)*t_in(i,j,k) + ent(i,k)*t_f(i,j,k-1))
372  s_f(i,j,k) = b1(i) * (h(i,j,k)*s_in(i,j,k) + ent(i,k)*s_f(i,j,k-1))
373  enddo ; enddo
374  do i=is,ie
375  c1(i,nz) = ent(i,nz) * b1(i)
376  b1(i) = 1.0 / (h(i,j,nz) + d1(i)*ent(i,nz) + h_neglect)
377  t_f(i,j,nz) = b1(i) * (h(i,j,nz)*t_in(i,j,nz) + ent(i,nz)*t_f(i,j,nz-1))
378  s_f(i,j,nz) = b1(i) * (h(i,j,nz)*s_in(i,j,nz) + ent(i,nz)*s_f(i,j,nz-1))
379  enddo
380  do k=nz-1,1,-1 ; do i=is,ie
381  t_f(i,j,k) = t_f(i,j,k) + c1(i,k+1)*t_f(i,j,k+1)
382  s_f(i,j,k) = s_f(i,j,k) + c1(i,k+1)*s_f(i,j,k+1)
383  enddo ; enddo
384  enddo
385  endif
386 
mom_grid::ocean_grid_type
Ocean grid type. See mom_grid for details.
Definition: MOM_grid.F90:19
Here is the caller graph for this function: