118 type(adapt_CS),
intent(in) :: CS
119 type(ocean_grid_type),
intent(in) :: G
120 type(verticalGrid_type),
intent(in) :: GV
121 type(thermo_var_ptrs),
intent(in) :: tv
123 integer,
intent(in) :: i
124 integer,
intent(in) :: j
125 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)+1),
intent(in) :: zInt
126 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)+1),
intent(in) :: tInt
127 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)+1),
intent(in) :: sInt
128 real,
dimension(SZI_(G),SZJ_(G),SZK_(GV)),
intent(in) :: h
129 real,
dimension(SZK_(GV)+1),
intent(inout) :: zNext
133 real :: h_up, b1, b_denom_1, d1, depth, drdz, nominal_z, stretching
134 real,
dimension(SZK_(GV)+1) :: alpha, beta, del2sigma
135 real,
dimension(SZK_(GV)) :: kGrid, c1
141 znext(nz+1) = zint(i,j,nz+1)
144 depth = g%bathyT(i,j) * gv%Z_to_H
154 if (g%mask2dT(i,j-1) > 0.)
then
155 call calculate_density_derivs( &
156 0.5 * (tint(i,j,2:nz) + tint(i,j-1,2:nz)), &
157 0.5 * (sint(i,j,2:nz) + sint(i,j-1,2:nz)), &
158 0.5 * (zint(i,j,2:nz) + zint(i,j-1,2:nz)) * gv%H_to_Pa, &
159 alpha, beta, 2, nz - 1, tv%eqn_of_state)
161 del2sigma(2:nz) = del2sigma(2:nz) + &
162 (alpha(2:nz) * (tint(i,j-1,2:nz) - tint(i,j,2:nz)) + &
163 beta(2:nz) * (sint(i,j-1,2:nz) - sint(i,j,2:nz)))
166 if (g%mask2dT(i,j+1) > 0.)
then
167 call calculate_density_derivs( &
168 0.5 * (tint(i,j,2:nz) + tint(i,j+1,2:nz)), &
169 0.5 * (sint(i,j,2:nz) + sint(i,j+1,2:nz)), &
170 0.5 * (zint(i,j,2:nz) + zint(i,j+1,2:nz)) * gv%H_to_Pa, &
171 alpha, beta, 2, nz - 1, tv%eqn_of_state)
173 del2sigma(2:nz) = del2sigma(2:nz) + &
174 (alpha(2:nz) * (tint(i,j+1,2:nz) - tint(i,j,2:nz)) + &
175 beta(2:nz) * (sint(i,j+1,2:nz) - sint(i,j,2:nz)))
178 if (g%mask2dT(i-1,j) > 0.)
then
179 call calculate_density_derivs( &
180 0.5 * (tint(i,j,2:nz) + tint(i-1,j,2:nz)), &
181 0.5 * (sint(i,j,2:nz) + sint(i-1,j,2:nz)), &
182 0.5 * (zint(i,j,2:nz) + zint(i-1,j,2:nz)) * gv%H_to_Pa, &
183 alpha, beta, 2, nz - 1, tv%eqn_of_state)
185 del2sigma(2:nz) = del2sigma(2:nz) + &
186 (alpha(2:nz) * (tint(i-1,j,2:nz) - tint(i,j,2:nz)) + &
187 beta(2:nz) * (sint(i-1,j,2:nz) - sint(i,j,2:nz)))
190 if (g%mask2dT(i+1,j) > 0.)
then
191 call calculate_density_derivs( &
192 0.5 * (tint(i,j,2:nz) + tint(i+1,j,2:nz)), &
193 0.5 * (sint(i,j,2:nz) + sint(i+1,j,2:nz)), &
194 0.5 * (zint(i,j,2:nz) + zint(i+1,j,2:nz)) * gv%H_to_Pa, &
195 alpha, beta, 2, nz - 1, tv%eqn_of_state)
197 del2sigma(2:nz) = del2sigma(2:nz) + &
198 (alpha(2:nz) * (tint(i+1,j,2:nz) - tint(i,j,2:nz)) + &
199 beta(2:nz) * (sint(i+1,j,2:nz) - sint(i,j,2:nz)))
208 call calculate_density_derivs(tint(i,j,:), sint(i,j,:), zint(i,j,:) * gv%H_to_Pa, &
209 alpha, beta, 1, nz + 1, tv%eqn_of_state)
212 del2sigma(k) = del2sigma(k) * (0.5 * (h(i,j,k-1) + h(i,j,k))) / &
213 max(alpha(k) * (tv%T(i,j,k) - tv%T(i,j,k-1)) + &
214 beta(k) * (tv%S(i,j,k) - tv%S(i,j,k-1)), 1e-20)
219 h_up = merge(h(i,j,k), h(i,j,k-1), del2sigma(k) > 0.)
220 del2sigma(k) = 0.5 * cs%adaptAlpha * &
221 sign(min(abs(del2sigma(k)), 0.5 * h_up), del2sigma(k))
224 znext(k) = zint(i,j,k) + del2sigma(k)
235 drdz = 0.5 * (alpha(k) + alpha(k+1)) * (tint(i,j,k+1) - tint(i,j,k)) + &
236 0.5 * (beta(k) + beta(k+1)) * (sint(i,j,k+1) - sint(i,j,k))
238 drdz = drdz / (znext(k) - znext(k+1) + gv%H_subroundoff)
243 kgrid(k) = (cs%adaptTimeRatio * nz**2 * depth) * &
244 (cs%adaptZoomCoeff / (cs%adaptZoom + 0.5*(znext(k) + znext(k+1))) + &
245 (cs%adaptBuoyCoeff * drdz / cs%adaptDrho0) + &
246 max(1.0 - cs%adaptZoomCoeff - cs%adaptBuoyCoeff, 0.0) / depth)
256 b_denom_1 = 1. + d1 * kgrid(k-1)
258 b1 = 1.0 / (b_denom_1 + kgrid(k))
261 c1(k) = kgrid(k) * b1
266 znext(k) = b1 * (znext(k) + kgrid(k-1)*znext(k-1))
270 znext(k) = znext(k) + c1(k)*znext(k+1)
273 if (cs%adaptDoMin)
then
275 stretching = zint(i,j,nz+1) / depth
278 nominal_z = nominal_z + cs%coordinateResolution(k-1) * stretching
280 znext(k) = max(znext(k), nominal_z)
282 znext(k) = min(znext(k), zint(i,j,nz+1))