ROMS
Loading...
Searching...
No Matches
rp_rho_eos.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined TL_IOMS && defined SOLVE3D
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine computes "in situ" density and other associated !
13! quantitites as a function of potential temperature, salinity, !
14! and pressure from a polynomial expression (Jackett and McDougall, !
15! 1992). The polynomial expression was found from fitting to 248 !
16! values in the oceanographic ranges of salinity, potential !
17! temperature, and pressure. It assumes no pressure variation !
18! along geopotential surfaces, that is, depth (meters; negative) !
19! and pressure (dbar; assumed negative here) are interchangeable. !
20! !
21! Check Values: (T=3 C, S=35.5 PSU, Z=-5000 m) !
22! !
23! alpha = 2.1014611551470d-04 (1/Celsius) !
24! beta = 7.2575037309946d-04 (1/PSU) !
25! gamma = 3.9684764511766d-06 (1/Pa) !
26! den = 1050.3639165364 (kg/m3) !
27! den1 = 1028.2845117925 (kg/m3) !
28! sound = 1548.8815240223 (m/s) !
29! bulk = 23786.056026320 (Pa) !
30! !
31! Reference: !
32! !
33! Jackett, D. R. and T. J. McDougall, 1995, Minimal Adjustment of !
34! Hydrostatic Profiles to Achieve Static Stability, J. of Atmos. !
35! and Oceanic Techn., vol. 12, pp. 381-389. !
36! !
37!=======================================================================
38!
39 implicit none
40!
41 PRIVATE
42 PUBLIC :: rp_rho_eos
43!
44 CONTAINS
45!
46!***********************************************************************
47 SUBROUTINE rp_rho_eos (ng, tile, model)
48!***********************************************************************
49!
50 USE mod_param
51 USE mod_coupling
52 USE mod_grid
53 USE mod_mixing
54 USE mod_ocean
55 USE mod_stepping
56!
57! Imported variable declarations.
58!
59 integer, intent(in) :: ng, tile, model
60!
61! Local variable declarations.
62!
63 character (len=*), parameter :: myfile = &
64 & __FILE__
65!
66# include "tile.h"
67!
68# ifdef PROFILE
69 CALL wclock_on (ng, model, 14, __line__, myfile)
70# endif
71 CALL rp_rho_eos_tile (ng, tile, model, &
72 & lbi, ubi, lbj, ubj, &
73 & imins, imaxs, jmins, jmaxs, &
74 & nrhs(ng), &
75# ifdef MASKING
76 & grid(ng) % rmask, &
77# endif
78# ifdef VAR_RHO_2D_NOT_YET
79 & grid(ng) % Hz, &
80 & grid(ng) % tl_Hz, &
81# endif
82 & grid(ng) % z_r, &
83 & grid(ng) % tl_z_r, &
84 & grid(ng) % z_w, &
85 & grid(ng) % tl_z_w, &
86 & ocean(ng) % t, &
87 & ocean(ng) % tl_t, &
88# ifdef VAR_RHO_2D_NOT_YET
89 & coupling(ng) % rhoA, &
90 & coupling(ng) % tl_rhoA, &
91 & coupling(ng) % rhoS, &
92 & coupling(ng) % tl_rhoS, &
93# endif
94# ifdef BV_FREQUENCY_NOT_YET
95 & mixing(ng) % tl_bvf, &
96# endif
97# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
98 defined bulk_fluxes
99 & mixing(ng) % alpha, &
100 & mixing(ng) % tl_alpha, &
101 & mixing(ng) % beta, &
102 & mixing(ng) % tl_beta, &
103# ifdef LMD_DDMIX_NOT_YET
104 & mixing(ng) % tl_alfaobeta, &
105# endif
106# endif
107# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
108 & ocean(ng) % tl_pden, &
109# endif
110 & ocean(ng) % rho, &
111 & ocean(ng) % tl_rho)
112# ifdef PROFILE
113 CALL wclock_off (ng, model, 14, __line__, myfile)
114# endif
115!
116 RETURN
117 END SUBROUTINE rp_rho_eos
118
119# ifdef NONLIN_EOS
120!
121!***********************************************************************
122 SUBROUTINE rp_rho_eos_tile (ng, tile, model, &
123 & LBi, UBi, LBj, UBj, &
124 & IminS, ImaxS, JminS, JmaxS, &
125 & nrhs, &
126# ifdef MASKING
127 & rmask, &
128# endif
129# ifdef VAR_RHO_2D_NOT_YET
130 & Hz, tl_Hz, &
131# endif
132 & z_r, tl_z_r, &
133 & z_w, tl_z_w, &
134 & t, tl_t, &
135# ifdef VAR_RHO_2D_NOT_YET
136 & rhoA, tl_rhoA, &
137 & rhoS, tl_rhoS, &
138# endif
139# ifdef BV_FREQUENCY_NOT_YET
140 & tl_bvf, &
141# endif
142# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
143 defined BULK_FLUXES
144 & alpha, tl_alpha, &
145 & beta, tl_beta, &
146# ifdef LMD_DDMIX_NOT_YET
147 & tl_alfaobeta, &
148# endif
149# endif
150# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
151 & tl_pden, &
152# endif
153 & rho, tl_rho)
154!***********************************************************************
155!
156 USE mod_param
157 USE mod_eoscoef
158 USE mod_scalars
159# ifdef SEDIMENT_NOT_YET
160 USE mod_sediment
161# endif
162!
165# ifdef DISTRIBUTE
167# endif
168!
169! Imported variable declarations.
170!
171 integer, intent(in) :: ng, tile, model
172 integer, intent(in) :: LBi, UBi, LBj, UBj
173 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
174 integer, intent(in) :: nrhs
175!
176# ifdef ASSUMED_SHAPE
177# ifdef MASKING
178 real(r8), intent(in) :: rmask(LBi:,LBj:)
179# endif
180# ifdef VAR_RHO_2D_NOT_YET
181 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
182# endif
183 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
184 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
185 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
186# ifdef VAR_RHO_2D_NOT_YET
187 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
188# endif
189 real(r8), intent(in) :: tl_z_r(LBi:,LBj:,:)
190 real(r8), intent(in) :: tl_z_w(LBi:,LBj:,0:)
191 real(r8), intent(in) :: tl_t(LBi:,LBj:,:,:,:)
192# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
193 defined bulk_fluxes
194 real(r8), intent(inout) :: alpha(LBi:,LBj:)
195 real(r8), intent(inout) :: beta(LBi:,LBj:)
196# endif
197# ifdef VAR_RHO_2D_NOT_YET
198 real(r8), intent(out) :: rhoA(LBi:,LBj:)
199 real(r8), intent(out) :: rhoS(LBi:,LBj:)
200 real(r8), intent(out) :: tl_rhoA(LBi:,LBj:)
201 real(r8), intent(out) :: tl_rhoS(LBi:,LBj:)
202# endif
203# ifdef BV_FREQUENCY_NOT_YET
204 real(r8), intent(out) :: tl_bvf(LBi:,LBj:,0:)
205# endif
206# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
207 defined bulk_fluxes
208 real(r8), intent(out) :: tl_alpha(LBi:,LBj:)
209 real(r8), intent(out) :: tl_beta(LBi:,LBj:)
210# ifdef LMD_DDMIX_NOT_YET
211 real(r8), intent(out) :: tl_alfaobeta(LBi:,LBj:,0:)
212# endif
213# endif
214# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
215 real(r8), intent(out) :: tl_pden(LBi:,LBj:,:)
216# endif
217 real(r8), intent(out) :: rho(LBi:,LBj:,:)
218 real(r8), intent(out) :: tl_rho(LBi:,LBj:,:)
219# else
220# ifdef MASKING
221 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
222# endif
223# ifdef VAR_RHO_2D_NOT_YET
224 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
225# endif
226 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
227 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
228 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
229# ifdef VAR_RHO_2D_NOT_YET
230 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
231# endif
232 real(r8), intent(in) :: tl_z_r(LBi:UBi,LBj:UBj,N(ng))
233 real(r8), intent(in) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
234 real(r8), intent(in) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
235# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
236 defined bulk_fluxes
237 real(r8), intent(inout) :: alpha(LBi:UBi,LBj:UBj)
238 real(r8), intent(inout) :: beta(LBi:UBi,LBj:UBj)
239# endif
240# ifdef VAR_RHO_2D_NOT_YET
241 real(r8), intent(out) :: rhoA(LBi:UBi,LBj:UBj)
242 real(r8), intent(out) :: rhoS(LBi:UBi,LBj:UBj)
243 real(r8), intent(out) :: tl_rhoA(LBi:UBi,LBj:UBj)
244 real(r8), intent(out) :: tl_rhoS(LBi:UBi,LBj:UBj)
245# endif
246# ifdef BV_FREQUENCY_NOT_YET
247 real(r8), intent(out) :: tl_bvf(LBi:UBi,LBj:UBj,0:N(ng))
248# endif
249# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
250 defined bulk_fluxes
251 real(r8), intent(out) :: tl_alpha(LBi:UBi,LBj:UBj)
252 real(r8), intent(out) :: tl_beta(LBi:UBi,LBj:UBj)
253# ifdef LMD_DDMIX_NOT_YET
254 real(r8), intent(out) :: tl_alfaobeta(LBi:UBi,LBj:UBj,0:N(ng))
255# endif
256# endif
257# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
258 real(r8), intent(out) :: tl_pden(LBi:UBi,LBj:UBj,N(ng))
259# endif
260 real(r8), intent(out) :: rho(LBi:UBi,LBj:UBj,N(ng))
261 real(r8), intent(out) :: tl_rho(LBi:UBi,LBj:UBj,N(ng))
262# endif
263!
264! Local variable declarations.
265!
266 integer :: i, ised, itrc, j, k
267
268 real(r8) :: SedDen, Tp, Tpr10, Ts, Tt, sqrtTs
269 real(r8) :: tl_SedDen, tl_Tp, tl_Tpr10, tl_Ts, tl_Tt, tl_sqrtTs
270# ifdef BV_FREQUENCY_NOT_YET
271 real(r8) :: bulk_dn, bulk_up, den_dn, den_up
272 real(r8) :: tl_bulk_dn, tl_bulk_up, tl_den_dn, tl_den_up
273# endif
274 real(r8) :: cff, cff1, cff2, cff3
275 real(r8) :: tl_cff, tl_cff1, tl_cff2, tl_cff3
276
277 real(r8), dimension(0:9) :: C
278 real(r8), dimension(0:9) :: tl_C
279# ifdef EOS_TDERIVATIVE
280 real(r8), dimension(0:9) :: dCdT(0:9)
281 real(r8), dimension(0:9) :: tl_dCdT(0:9)
282 real(r8), dimension(0:9) :: d2Cd2T(0:9)
283
284 real(r8), dimension(IminS:ImaxS,N(ng)) :: DbulkDS
285 real(r8), dimension(IminS:ImaxS,N(ng)) :: DbulkDT
286 real(r8), dimension(IminS:ImaxS,N(ng)) :: Dden1DS
287 real(r8), dimension(IminS:ImaxS,N(ng)) :: Dden1DT
288 real(r8), dimension(IminS:ImaxS,N(ng)) :: Scof
289 real(r8), dimension(IminS:ImaxS,N(ng)) :: Tcof
290 real(r8), dimension(IminS:ImaxS,N(ng)) :: wrk
291
292 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_DbulkDS
293 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_DbulkDT
294 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_Dden1DS
295 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_Dden1DT
296 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_Scof
297 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_Tcof
298 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_wrk
299# endif
300 real(r8), dimension(IminS:ImaxS,N(ng)) :: bulk
301 real(r8), dimension(IminS:ImaxS,N(ng)) :: bulk0
302 real(r8), dimension(IminS:ImaxS,N(ng)) :: bulk1
303 real(r8), dimension(IminS:ImaxS,N(ng)) :: bulk2
304 real(r8), dimension(IminS:ImaxS,N(ng)) :: den
305 real(r8), dimension(IminS:ImaxS,N(ng)) :: den1
306
307 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_bulk
308 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_bulk0
309 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_bulk1
310 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_bulk2
311 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_den
312 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_den1
313
314# include "set_bounds.h"
315!
316!=======================================================================
317! Nonlinear equation of state. Notice that this equation of state
318! is only valid for potential temperature range of -2C to 40C and
319! a salinity range of 0 PSU to 42 PSU.
320!=======================================================================
321!
322 DO j=jstrt,jendt
323 DO k=1,n(ng)
324 DO i=istrt,iendt
325!
326! Check temperature and salinity lower values. Assign depth to the
327! pressure.
328!
329 tt=max(-2.0_r8,t(i,j,k,nrhs,itemp))
330 tl_tt=(0.5_r8-sign(0.5_r8,-2.0_r8-t(i,j,k,nrhs,itemp)))* &
331 & tl_t(i,j,k,nrhs,itemp)- &
332# ifdef TL_IOMS
333 & 2.0_r8*(0.5_r8+sign(0.5_r8, &
334 & -2.0_r8-t(i,j,k,nrhs,itemp)))
335# endif
336# ifdef SALINITY
337 ts=max(0.0_r8,t(i,j,k,nrhs,isalt))
338 tl_ts=(0.5_r8-sign(0.5_r8,-t(i,j,k,nrhs,isalt)))* &
339 & tl_t(i,j,k,nrhs,isalt)
340 sqrtts=sqrt(ts)
341 IF (ts.ne.0.0_r8) THEN
342 tl_sqrtts=0.5_r8*(tl_ts/sqrtts)+ &
343# ifdef TL_IOMS
344 & 0.5_r8*sqrtts
345# endif
346 ELSE
347 tl_sqrtts=0.0_r8
348 END IF
349# else
350 ts=0.0_r8
351 tl_ts=0.0_r8
352 sqrtts=0.0_r8
353 tl_sqrtts=0.0_r8
354# endif
355 tp=z_r(i,j,k)
356 tl_tp=tl_z_r(i,j,k)
357 tpr10=0.1_r8*tp
358 tl_tpr10=0.1_r8*tl_tp
359!
360!-----------------------------------------------------------------------
361! Compute BASIC STATE and tangent linear density (kg/m3) at standard
362! one atmosphere pressure.
363!-----------------------------------------------------------------------
364!
365 c(0)=q00+tt*(q01+tt*(q02+tt*(q03+tt*(q04+tt*q05))))
366 c(1)=u00+tt*(u01+tt*(u02+tt*(u03+tt*u04)))
367 c(2)=v00+tt*(v01+tt*v02)
368# ifdef EOS_TDERIVATIVE
369!
370 dcdt(0)=q01+tt*(2.0_r8*q02+tt*(3.0_r8*q03+tt*(4.0_r8*q04+ &
371 & tt*5.0_r8*q05)))
372 dcdt(1)=u01+tt*(2.0_r8*u02+tt*(3.0_r8*u03+tt*4.0_r8*u04))
373 dcdt(2)=v01+tt*2.0_r8*v02
374# endif
375 tl_c(0)=tl_tt*dcdt(0)+ &
376# ifdef TL_IOMS
377 & q00-tt*tt*(q02+tt*(2.0_r8*q03+tt*(3.0_r8*q04+ &
378 & tt*4.0_r8*q05)))
379# endif
380 tl_c(1)=tl_tt*dcdt(1)+ &
381# ifdef TL_IOMS
382 & u00-tt*tt*(u02+tt*(2.0_r8*u03+tt*3.0_r8*u04))
383# endif
384 tl_c(2)=tl_tt*dcdt(2)+ &
385# ifdef TL_IOMS
386 & v00-v02*tt*tt
387# endif
388!
389 den1(i,k)=c(0)+ts*(c(1)+sqrtts*c(2)+ts*w00)
390 tl_den1(i,k)=tl_c(0)+ &
391 & tl_ts*(c(1)+sqrtts*c(2)+ts*w00)+ &
392 & ts*(tl_c(1)+tl_sqrtts*c(2)+ &
393 & sqrtts*tl_c(2)+tl_ts*w00)- &
394# ifdef TL_IOMS
395 & ts*(c(1)+2.0_r8*sqrtts*c(2)+ts*w00)
396# endif
397# ifdef EOS_TDERIVATIVE
398!
399! Compute d(den1)/d(S) and d(den1)/d(T) derivatives used in the
400! computation of thermal expansion and saline contraction
401! coefficients.
402!
403 d2cd2t(0)=2.0_r8*q02+tt*(6.0_r8*q03+tt*(12.0_r8*q04+ &
404 & tt*20.0_r8*q05))
405 d2cd2t(1)=2.0_r8*u02+tt*(6.0_r8*u03+tt*12.0_r8*u04)
406 d2cd2t(2)=2.0_r8*v02
407!
408 tl_dcdt(0)=tl_tt*d2cd2t(0)+ &
409# ifdef TL_IOMS
410 & q01-tt*tt*(3.0_r8*q03+tt*(8.0_r8*q04+ &
411 & tt*15.0_r8*q05*tt))
412# endif
413 tl_dcdt(1)=tl_tt*d2cd2t(1)+ &
414# ifdef TL_IOMS
415 & u01-tt*tt*(3.0_r8*u03+tt*8.0_r8*u04)
416# endif
417 tl_dcdt(2)=tl_tt*d2cd2t(2)+ &
418# ifdef TL_IOMS
419 & v01
420# endif
421!
422 dden1ds(i,k)=c(1)+1.5_r8*c(2)*sqrtts+2.0_r8*w00*ts
423 dden1dt(i,k)=dcdt(0)+ts*(dcdt(1)+sqrtts*dcdt(2))
424!
425 tl_dden1ds(i,k)=tl_c(1)+ &
426 & 1.5_r8*(tl_c(2)*sqrtts+ &
427 & c(2)*tl_sqrtts)+ &
428 & 2.0_r8*w00*tl_ts- &
429# ifdef TL_IOMS
430 & 1.5_r8*c(2)*sqrtts
431# endif
432 tl_dden1dt(i,k)=tl_dcdt(0)+ &
433 & tl_ts*(dcdt(1)+sqrtts*dcdt(2))+ &
434 & ts*(tl_dcdt(1)+tl_sqrtts*dcdt(2)+ &
435 & sqrtts*tl_dcdt(2))- &
436# ifdef TL_IOMS
437 & ts*(dcdt(1)+2.0_r8*sqrtts*dcdt(2))
438# endif
439# endif
440!
441!-----------------------------------------------------------------------
442! Compute BASIC STATE and tangent linear secant bulk modulus.
443!-----------------------------------------------------------------------
444!
445 c(3)=a00+tt*(a01+tt*(a02+tt*(a03+tt*a04)))
446 c(4)=b00+tt*(b01+tt*(b02+tt*b03))
447 c(5)=d00+tt*(d01+tt*d02)
448 c(6)=e00+tt*(e01+tt*(e02+tt*e03))
449 c(7)=f00+tt*(f01+tt*f02)
450 c(8)=g01+tt*(g02+tt*g03)
451 c(9)=h00+tt*(h01+tt*h02)
452# ifdef EOS_TDERIVATIVE
453!
454 dcdt(3)=a01+tt*(2.0_r8*a02+tt*(3.0_r8*a03+tt*4.0_r8*a04))
455 dcdt(4)=b01+tt*(2.0_r8*b02+tt*3.0_r8*b03)
456 dcdt(5)=d01+tt*2.0_r8*d02
457 dcdt(6)=e01+tt*(2.0_r8*e02+tt*3.0_r8*e03)
458 dcdt(7)=f01+tt*2.0_r8*f02
459 dcdt(8)=g02+tt*2.0_r8*g03
460 dcdt(9)=h01+tt*2.0_r8*h02
461# endif
462!
463 tl_c(3)=tl_tt*dcdt(3)+ &
464# ifdef TL_IOMS
465 & a00-tt*tt*(a02+tt*(2.0_r8*a03+tt*3.0_r8*a04))
466# endif
467 tl_c(4)=tl_tt*dcdt(4)+ &
468# ifdef TL_IOMS
469 & b00-tt*tt*(b02+tt*2.0_r8*b03)
470# endif
471 tl_c(5)=tl_tt*dcdt(5)+ &
472# ifdef TL_IOMS
473 & d00-tt*tt*d02
474# endif
475 tl_c(6)=tl_tt*dcdt(6)+ &
476# ifdef TL_IOMS
477 & e00-tt*tt*(e02+tt*2.0_r8*e03)
478# endif
479 tl_c(7)=tl_tt*dcdt(7)+ &
480# ifdef TL_IOMS
481 & f00-tt*tt*f02
482# endif
483 tl_c(8)=tl_tt*dcdt(8)+ &
484# ifdef TL_IOMS
485 & g01-tt*tt*g03
486# endif
487 tl_c(9)=tl_tt*dcdt(9)+ &
488# ifdef TL_IOMS
489 & h00-tt*tt*h02
490# endif
491!
492 bulk0(i,k)=c(3)+ts*(c(4)+sqrtts*c(5))
493 bulk1(i,k)=c(6)+ts*(c(7)+sqrtts*g00)
494 bulk2(i,k)=c(8)+ts*c(9)
495 bulk(i,k)=bulk0(i,k)-tp*(bulk1(i,k)-tp*bulk2(i,k))
496!
497 tl_bulk0(i,k)=tl_c(3)+ &
498 & tl_ts*(c(4)+sqrtts*c(5))+ &
499 & ts*(tl_c(4)+tl_sqrtts*c(5)+ &
500 & sqrtts*tl_c(5))- &
501# ifdef TL_IOMS
502 & ts*(c(4)+2.0_r8*sqrtts*c(5))
503# endif
504 tl_bulk1(i,k)=tl_c(6)+ &
505 & tl_ts*(c(7)+sqrtts*g00)+ &
506 & ts*(tl_c(7)+tl_sqrtts*g00)- &
507# ifdef TL_IOMS
508 & ts*(c(7)+sqrtts*g00)
509# endif
510 tl_bulk2(i,k)=tl_c(8)+tl_ts*c(9)+ts*tl_c(9)- &
511# ifdef TL_IOMS
512 & ts*c(9)
513# endif
514 tl_bulk(i,k)=tl_bulk0(i,k)- &
515 & tl_tp*(bulk1(i,k)-tp*bulk2(i,k))- &
516 & tp*(tl_bulk1(i,k)- &
517 & tl_tp*bulk2(i,k)- &
518 & tp*tl_bulk2(i,k))+ &
519# ifdef TL_IOMS
520 & tp*(bulk1(i,k)-2.0_r8*tp*bulk2(i,k))
521# endif
522
523# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
524 defined bulk_fluxes
525!
526! Compute d(bulk)/d(S) and d(bulk)/d(T) derivatives used
527! in the computation of thermal expansion and saline contraction
528! coefficients.
529!
530 d2cd2t(3)=2.0_r8*a02+tt*(6.0_r8*a03+tt*12.0_r8*a04)
531 d2cd2t(4)=2.0_r8*b02+tt*6.0_r8*b03
532 d2cd2t(5)=2.0_r8*d02
533 d2cd2t(6)=2.0_r8*e02+tt*6.0_r8*e03
534 d2cd2t(7)=2.0_r8*f02
535 d2cd2t(8)=2.0_r8*g03
536 d2cd2t(9)=2.0_r8*h02
537!
538 tl_dcdt(3)=tl_tt*d2cd2t(3)+ &
539# ifdef TL_IOMS
540 & a01-tt*tt*(3.0_r8*a03+tt*8.0_r8*a04)
541# endif
542 tl_dcdt(4)=tl_tt*d2cd2t(4)+ &
543# ifdef TL_IOMS
544 & b01-tt*tt*3.0_r8*b03
545# endif
546 tl_dcdt(5)=tl_tt*d2cd2t(5)+ &
547# ifdef TL_IOMS
548 & d01
549# endif
550 tl_dcdt(6)=tl_tt*d2cd2t(6)+ &
551# ifdef TL_IOMS
552 & e01-tt*tt*3.0_r8*e03
553# endif
554 tl_dcdt(7)=tl_tt*d2cd2t(7)+ &
555# ifdef TL_IOMS
556 & f01
557# endif
558 tl_dcdt(8)=tl_tt*d2cd2t(8)+ &
559# ifdef TL_IOMS
560 & g02
561# endif
562 tl_dcdt(9)=tl_tt*d2cd2t(9)+ &
563# ifdef TL_IOMS
564 & h01
565# endif
566!
567 dbulkds(i,k)=c(4)+sqrtts*1.5_r8*c(5)- &
568 & tp*(c(7)+sqrtts*1.5_r8*g00-tp*c(9))
569 dbulkdt(i,k)=dcdt(3)+ts*(dcdt(4)+sqrtts*dcdt(5))- &
570 & tp*(dcdt(6)+ts*dcdt(7)- &
571 & tp*(dcdt(8)+ts*dcdt(9)))
572!
573 tl_dbulkds(i,k)=tl_c(4)+ &
574 & 1.5_r8*(tl_sqrtts*c(5)+ &
575 & sqrtts*tl_c(5))- &
576 & tl_tp*(c(7)+sqrtts*1.5_r8*g00- &
577 & tp*c(9))- &
578 & tp*(tl_c(7)+tl_sqrtts*1.5_r8*g00- &
579 & tl_tp*c(9)-tp*tl_c(9))- &
580# ifdef TL_IOMS
581 & sqrtts*1.5_r8*c(5)+ &
582 & tp*(c(7)+sqrtts*1.5_r8*g00-2.0_r8*tp*c(9))
583# endif
584 tl_dbulkdt(i,k)=tl_dcdt(3)+ &
585 & tl_ts*(dcdt(4)+sqrtts*dcdt(5))+ &
586 & ts*(tl_dcdt(4)+tl_sqrtts*dcdt(5)+ &
587 & sqrtts*tl_dcdt(5))- &
588 & tl_tp*(dcdt(6)+ts*dcdt(7)- &
589 & tp*(dcdt(8)+ts*dcdt(9)))- &
590 & tp*(tl_dcdt(6)+tl_ts*dcdt(7)+ts*tl_dcdt(7)- &
591 & tl_tp*(dcdt(8)+ts*dcdt(9))- &
592 & tp*(tl_dcdt(8)+tl_ts*dcdt(9)+ &
593 & ts*tl_dcdt(9)))- &
594# ifdef TL_IOMS
595 & ts*(dcdt(4)+2.0_r8*sqrtts*dcdt(5))+ &
596 & tp*(dcdt(6)+2.0_r8*ts*dcdt(7)- &
597 & tp*(2.0_r8*dcdt(8)+ &
598 & 3.0_r8*ts*dcdt(9)))
599# endif
600# endif
601!
602!-----------------------------------------------------------------------
603! Compute local "in situ" density anomaly (kg/m3 - 1000).
604!-----------------------------------------------------------------------
605!
606 cff=1.0_r8/(bulk(i,k)+tpr10)
607 tl_cff=-cff*cff*(tl_bulk(i,k)+tl_tpr10)+ &
608# ifdef TL_IOMS
609 & 2.0_r8*cff
610# endif
611 den(i,k)=den1(i,k)*bulk(i,k)*cff
612 tl_den(i,k)=tl_den1(i,k)*bulk(i,k)*cff+ &
613 & den1(i,k)*(tl_bulk(i,k)*cff+ &
614 & bulk(i,k)*tl_cff)- &
615# ifdef TL_IOMS
616 & 2.0_r8*den(i,k)
617# endif
618# if defined SEDIMENT_NOT_YET && defined SED_DENS_NOT_YET
619 sedden=0.0_r8
620 tl_sedden=0.0_r8
621 DO ised=1,nst
622 itrc=idsed(ised)
623 cff1=1.0_r8/srho(ised,ng)
624 sedden=sedden+ &
625 & t(i,j,k,nrhs,itrc)* &
626 & (srho(ised,ng)-den(i,k))*cff1
627 tl_sedden=tl_sedden+ &
628 & (tl_t(i,j,k,nrhs,itrc)* &
629 & (srho(ised,ng)-den(i,k))- &
630 & t(i,j,k,nrhs,itrc)* &
631 & tl_den(i,k))*cff1+ &
632# ifdef TL_IOMS
633 & t(i,j,k,nrhs,itrc)*den(i,k)*cff1
634# endif
635 END DO
636 den(i,k)=den(i,k)+sedden
637 tl_den(i,k)=tl_den(i,k)+tl_sedden
638# endif
639 den(i,k)=den(i,k)-1000.0_r8
640# ifdef TL_IOMS
641 tl_den(i,k)=tl_den(i,k)-1000.0_r8
642# endif
643# ifdef MASKING
644 den(i,k)=den(i,k)*rmask(i,j)
645 tl_den(i,k)=tl_den(i,k)*rmask(i,j)
646# endif
647 END DO
648 END DO
649
650# ifdef VAR_RHO_2D_NOT_YET
651!
652!-----------------------------------------------------------------------
653! Compute vertical averaged density (rhoA) and density perturbation
654! (rhoS) used in barotropic pressure gradient.
655!-----------------------------------------------------------------------
656!
657 DO i=istrt,iendt
658 cff1=den(i,n(ng))*hz(i,j,n(ng))
659 tl_cff1=tl_den(i,n(ng))*hz(i,j,n(ng))+ &
660 & den(i,n(ng))*tl_hz(i,j,n(ng))- &
661# ifdef TL_IOMS
662 & cff1
663# endif
664 rhos(i,j)=0.5_r8*cff1*hz(i,j,n(ng))
665 tl_rhos(i,j)=0.5_r8*(tl_cff1*hz(i,j,n(ng))+ &
666 & cff1*tl_hz(i,j,n(ng)))- &
667# ifdef TL_IOMS
668 & rhos(i,j)
669# endif
670 rhoa(i,j)=cff1
671 tl_rhoa(i,j)=tl_cff1
672 END DO
673 DO k=n(ng)-1,1,-1
674 DO i=istrt,iendt
675 cff1=den(i,k)*hz(i,j,k)
676 tl_cff1=tl_den(i,k)*hz(i,j,k)+ &
677 & den(i,k)*tl_hz(i,j,k)- &
678# ifdef TL_IOMS
679 & cff1
680# endif
681 rhos(i,j)=rhos(i,j)+hz(i,j,k)*(rhoa(i,j)+0.5_r8*cff1)
682 tl_rhos(i,j)=tl_rhos(i,j)+ &
683 & tl_hz(i,j,k)*(rhoa(i,j)+0.5_r8*cff1)+ &
684 & hz(i,j,k)*(tl_rhoa(i,j)+0.5_r8*tl_cff1)- &
685# ifdef TL_IOMS
686 & hz(i,j,k)*(rhoa(i,j)+0.5_r8*cff1)
687# endif
688 rhoa(i,j)=rhoa(i,j)+cff1
689 tl_rhoa(i,j)=tl_rhoa(i,j)+tl_cff1
690 END DO
691 END DO
692 cff2=1.0_r8/rho0
693 DO i=istrt,iendt
694 cff1=1.0_r8/(z_w(i,j,n(ng))-z_w(i,j,0))
695 tl_cff1=-cff1*cff1*(tl_z_w(i,j,n(ng))-tl_z_w(i,j,0))+ &
696# ifdef TL_IOMS
697 & 2.0_r8*cff1
698# endif
699!
700! Here we reverse the order of the NL and TL operations since an
701! intermeridiate value of rhoA and rhoS is needed because they are
702! recursive.
703!
704 tl_rhoa(i,j)=cff2*(tl_cff1*rhoa(i,j)+cff1*tl_rhoa(i,j))
705 rhoa(i,j)=cff2*cff1*rhoa(i,j)
706# ifdef TL_IOMS
707 tl_rhoa(i,j)=tl_rhoa(i,j)-rhoa(i,j)
708# endif
709 tl_rhos(i,j)=2.0_r8*cff2* &
710 & cff1*(2.0_r8*tl_cff1*rhos(i,j)+ &
711 & cff1*tl_rhos(i,j))
712 rhos(i,j)=2.0_r8*cff1*cff1*cff2*rhos(i,j)
713# ifdef TL_IOMS
714 tl_rhos(i,j)=tl_rhos(i,j)-2.0_r8*rhos(i,j)
715# endif
716 END DO
717# endif
718
719# if defined BV_FREQUENCY_NOT_YET
720!
721!-----------------------------------------------------------------------
722! Compute Brunt-Vaisala frequency (1/s2) at horizontal RHO-points
723! and vertical W-points:
724!
725! bvf = - g/rho d(rho)/d(z).
726!
727! The density anomaly difference is computed by lowering/rising the
728! water parcel above/below adiabatically at W-point depth "z_w".
729!-----------------------------------------------------------------------
730!
731 DO k=1,n(ng)-1
732 DO i=istrt,iendt
733 bulk_up=bulk0(i,k+1)- &
734 & z_w(i,j,k)*(bulk1(i,k+1)- &
735 & bulk2(i,k+1)*z_w(i,j,k))
736 tl_bulk_up=tl_bulk0(i,k+1)- &
737 & tl_z_w(i,j,k)*(bulk1(i,k+1)- &
738 & bulk2(i,k+1)*z_w(i,j,k))- &
739 & z_w(i,j,k)*(tl_bulk1(i,k+1)- &
740 & tl_bulk2(i,k+1)*z_w(i,j,k)- &
741 & bulk2(i,k+1)*tl_z_w(i,j,k))+ &
742# ifdef TL_IOMS
743 & z_w(i,j,k)*(bulk1(i,k+1)- &
744 & 2.0_r8*bulk2(i,k+1)*z_w(i,j,k))
745# endif
746 bulk_dn=bulk0(i,k )- &
747 & z_w(i,j,k)*(bulk1(i,k )- &
748 & bulk2(i,k )*z_w(i,j,k))
749 tl_bulk_dn=tl_bulk0(i,k )- &
750 & tl_z_w(i,j,k)*(bulk1(i,k )- &
751 & bulk2(i,k )*z_w(i,j,k))- &
752 & z_w(i,j,k)*(tl_bulk1(i,k )- &
753 & tl_bulk2(i,k )*z_w(i,j,k)- &
754 & bulk2(i,k )*tl_z_w(i,j,k))+ &
755# ifdef TL_IOMS
756 & z_w(i,j,k)*(bulk1(i,k )- &
757 & 2.0_r8*bulk2(i,k )*z_w(i,j,k))
758# endif
759 cff1=1.0_r8/(bulk_up+0.1_r8*z_w(i,j,k))
760 cff2=1.0_r8/(bulk_dn+0.1_r8*z_w(i,j,k))
761 tl_cff1=-cff1*cff1*(tl_bulk_up+0.1_r8*tl_z_w(i,j,k))+ &
762# ifdef TL_IOMS
763 & 2.0_r8*cff1
764# endif
765 tl_cff2=-cff2*cff2*(tl_bulk_dn+0.1_r8*tl_z_w(i,j,k))+ &
766# ifdef TL_IOMS
767 & 2.0_r8*cff2
768# endif
769 den_up=cff1*(den1(i,k+1)*bulk_up)
770 den_dn=cff2*(den1(i,k )*bulk_dn)
771 tl_den_up=tl_cff1*(den1(i,k+1)*bulk_up)+ &
772 & cff1*(tl_den1(i,k+1)*bulk_up+ &
773 & den1(i,k+1)*tl_bulk_up)- &
774# ifdef TL_IOMS
775 & 2.0_r8*den_up
776# endif
777 tl_den_dn=tl_cff2*(den1(i,k )*bulk_dn)+ &
778 & cff2*(tl_den1(i,k )*bulk_dn+ &
779 & den1(i,k )*tl_bulk_dn)- &
780# ifdef TL_IOMS
781 & 2.0_r8*den_dn
782# endif
783!^ bvf(i,j,k)=-g*(den_up-den_dn)/ &
784!^ & (0.5_r8*(den_up+den_dn)* &
785!^ & (z_r(i,j,k+1)-z_r(i,j,k)))
786!^
787 cff3=1.0_r8/(0.5_r8*(den_up+den_dn)* &
788 & (z_r(i,j,k+1)-z_r(i,j,k)))
789 tl_cff3=-cff3*cff3* &
790 & 0.5_r8*((tl_den_up+tl_den_dn)* &
791 & (z_r(i,j,k+1)-z_r(i,j,k))+ &
792 & (den_up+den_dn)* &
793 & (tl_z_r(i,j,k+1)-tl_z_r(i,j,k)))+ &
794# ifdef TL_IOMS
795 & 3.0_r8*cff3
796# endif
797 tl_bvf(i,j,k)=-g*((tl_den_up-tl_den_dn)*cff3+ &
798 & (den_up-den_dn)*tl_cff3)+ &
799# ifdef TL_IOMS
800 & 2.0_r8*g*(den_up-den_dn)*cff3
801# endif
802 END DO
803 END DO
804 DO i=istrt,iendt
805!^ bvf(i,j,0)=0.0_r8
806!^
807 tl_bvf(i,j,0)=0.0_r8
808!^ bvf(i,j,N(ng))=0.0_r8
809!^
810 tl_bvf(i,j,n(ng))=0.0_r8
811 END DO
812# endif
813
814# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
815 defined bulk_fluxes
816!
817!-----------------------------------------------------------------------
818! Compute thermal expansion (1/Celsius) and saline contraction
819! (1/PSU) coefficients.
820!-----------------------------------------------------------------------
821!
822# ifdef LMD_DDMIX_NOT_YET
823 DO k=1,n(ng)
824# else
825 DO k=n(ng),n(ng)
826# endif
827 DO i=istrt,iendt
828 tpr10=0.1_r8*z_r(i,j,k)
829 tl_tpr10=0.1_r8*tl_z_r(i,j,k)
830!
831! Compute thermal expansion and saline contraction coefficients.
832!
833 cff=bulk(i,k)+tpr10
834 tl_cff=tl_bulk(i,k)+tl_tpr10
835 cff1=tpr10*den1(i,k)
836 tl_cff1=tl_tpr10*den1(i,k)+tpr10*tl_den1(i,k)- &
837# ifdef TL_IOMS
838 & cff1
839# endif
840 cff2=bulk(i,k)*cff
841 tl_cff2=tl_bulk(i,k)*cff+bulk(i,k)*tl_cff- &
842# ifdef TL_IOMS
843 & cff2
844# endif
845 wrk(i,k)=(den(i,k)+1000.0_r8)*cff*cff
846 tl_wrk(i,k)=cff*(cff*tl_den(i,k)+ &
847 & 2.0_r8*tl_cff*(den(i,k)+1000.0_r8))- &
848# ifdef TL_IOMS
849 & cff*cff*(2.0_r8*den(i,k)+1000.0_r8)
850# endif
851 tcof(i,k)=-(dbulkdt(i,k)*cff1+ &
852 & dden1dt(i,k)*cff2)
853 tl_tcof(i,k)=-(tl_dbulkdt(i,k)*cff1+ &
854 & dbulkdt(i,k)*tl_cff1+ &
855 & tl_dden1dt(i,k)*cff2+ &
856 & dden1dt(i,k)*tl_cff2)- &
857# ifdef TL_IOMS
858 & tcof(i,k)
859# endif
860 scof(i,k)= (dbulkds(i,k)*cff1+ &
861 & dden1ds(i,k)*cff2)
862 tl_scof(i,k)= (tl_dbulkds(i,k)*cff1+ &
863 & dbulkds(i,k)*tl_cff1+ &
864 & tl_dden1ds(i,k)*cff2+ &
865 & dden1ds(i,k)*tl_cff2)- &
866# ifdef TL_IOMS
867 & scof(i,k)
868# endif
869# ifdef LMD_DDMIX_NOT_YET
870!^ alfaobeta(i,j,k)=Tcof(i,k)/Scof(i,k)
871!^
872 tl_alfaobeta(i,j,k)=(tl_tcof(i,k)*scof(i,k)- &
873 & tcof(i,k)*tl_scof(i,k))/ &
874 & (scof(i,k)*scof(i,k))+ &
875# ifdef TL_IOMS
876 & tcof(i,k)/scof(i,k)
877# endif
878# endif
879 END DO
880 IF (k.eq.n(ng)) THEN
881 DO i=istrt,iendt
882 cff=1.0_r8/wrk(i,n(ng))
883 tl_cff=-cff*cff*tl_wrk(i,n(ng))+ &
884# ifdef TL_IOMS
885 & 2.0_r8*cff
886# endif
887 alpha(i,j)=cff*tcof(i,n(ng))
888 tl_alpha(i,j)=tl_cff*tcof(i,n(ng))+cff*tl_tcof(i,n(ng))- &
889# ifdef TL_IOMS
890 & alpha(i,j)
891# endif
892 beta(i,j)=cff*scof(i,n(ng))
893 tl_beta(i,j)=tl_cff*scof(i,n(ng))+cff*tl_scof(i,n(ng))- &
894# ifdef TL_IOMS
895 & beta(i,j)
896# endif
897 END DO
898 END IF
899 END DO
900# endif
901!
902!-----------------------------------------------------------------------
903! Load "in situ" density anomaly (kg/m3 - 1000) and potential
904! density anomaly (kg/m3 - 1000) referenced to the surface into global
905! arrays. Notice that this is done in a separate (i,k) DO-loops to
906! facilitate the adjoint.
907!-----------------------------------------------------------------------
908!
909 DO k=1,n(ng)
910 DO i=istrt,iendt
911 rho(i,j,k)=den(i,k)
912 tl_rho(i,j,k)=tl_den(i,k)
913# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
914!^ pden(i,j,k)=(den1(i,k)-1000.0_r8) ! This gives a fatal
915!^ ! result in 4D-Var
916 tl_pden(i,j,k)=tl_den1(i,k)- &
917# ifdef TL_IOMS
918 & 1000.0_r8 ! posterior error...
919# endif
920# ifdef MASKING
921!^ pden(i,j,k)=pden(i,k)*rmask(i,j)
922!^
923 tl_pden(i,j,k)=tl_pden(i,k)*rmask(i,j)
924# endif
925# endif
926 END DO
927 END DO
928 END DO
929!
930!-----------------------------------------------------------------------
931! Exchange boundary data.
932!-----------------------------------------------------------------------
933!
934 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
935 CALL exchange_r3d_tile (ng, tile, &
936 & lbi, ubi, lbj, ubj, 1, n(ng), &
937 & rho)
938 CALL exchange_r3d_tile (ng, tile, &
939 & lbi, ubi, lbj, ubj, 1, n(ng), &
940 & tl_rho)
941
942# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
943!^ CALL exchange_r3d_tile (ng, tile, &
944!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
945!^ & pden)
946!^
947 CALL exchange_r3d_tile (ng, tile, &
948 & lbi, ubi, lbj, ubj, 1, n(ng), &
949 & tl_pden)
950# endif
951
952# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
953 defined bulk_fluxes_not_yet
954# ifdef LMD_DDMIX_NOT_YET
955!^ CALL exchange_w3d_tile (ng, tile, &
956!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
957!^ & alfaobeta)
958!^
959 CALL exchange_w3d_tile (ng, tile, &
960 & lbi, ubi, lbj, ubj, 0, n(ng), &
961 & tl_alfaobeta)
962# endif
963 CALL exchange_r2d_tile (ng, tile, &
964 & lbi, ubi, lbj, ubj, &
965 & alpha)
966 CALL exchange_r2d_tile (ng, tile, &
967 & lbi, ubi, lbj, ubj, &
968 & tl_alpha)
969 CALL exchange_r2d_tile (ng, tile, &
970 & lbi, ubi, lbj, ubj, &
971 & beta)
972 CALL exchange_r2d_tile (ng, tile, &
973 & lbi, ubi, lbj, ubj, &
974 & tl_beta)
975# endif
976
977# ifdef VAR_RHO_2D_NOT_YET
978 CALL exchange_r2d_tile (ng, tile, &
979 & lbi, ubi, lbj, ubj, &
980 & rhoa)
981 CALL exchange_r2d_tile (ng, tile, &
982 & lbi, ubi, lbj, ubj, &
983 & tl_rhoa)
984 CALL exchange_r2d_tile (ng, tile, &
985 & lbi, ubi, lbj, ubj, &
986 & rhos)
987 CALL exchange_r2d_tile (ng, tile, &
988 & lbi, ubi, lbj, ubj, &
989 & tl_rhos)
990# endif
991
992# ifdef BV_FREQUENCY_NOT_YET
993!^ CALL exchange_w3d_tile (ng, tile, &
994!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
995!^ & bvf)
996!^
997 CALL exchange_w3d_tile (ng, tile, &
998 & lbi, ubi, lbj, ubj, 0, n(ng), &
999 & tl_bvf)
1000# endif
1001 END IF
1002
1003# ifdef DISTRIBUTE
1004!
1005 CALL mp_exchange3d (ng, tile, model, 1, &
1006 & lbi, ubi, lbj, ubj, 1, n(ng), &
1007 & nghostpoints, &
1008 & ewperiodic(ng), nsperiodic(ng), &
1009 & rho)
1010 CALL mp_exchange3d (ng, tile, model, 1, &
1011 & lbi, ubi, lbj, ubj, 1, n(ng), &
1012 & nghostpoints, &
1013 & ewperiodic(ng), nsperiodic(ng), &
1014 & tl_rho)
1015
1016# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1017!^ CALL mp_exchange3d (ng, tile, model, 1, &
1018!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
1019!^ & NghostPoints, &
1020!^ & EWperiodic(ng), NSperiodic(ng), &
1021!^ & pden)
1022!^
1023 CALL mp_exchange3d (ng, tile, model, 1, &
1024 & lbi, ubi, lbj, ubj, 1, n(ng), &
1025 & nghostpoints, &
1026 & ewperiodic(ng), nsperiodic(ng), &
1027 & tl_pden)
1028# endif
1029
1030# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1031 defined bulk_fluxes
1032# ifdef LMD_DDMIX_NOT_YET
1033!^ CALL mp_exchange3d (ng, tile, model, 1, &
1034!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1035!^ & NghostPoints, &
1036!^ & EWperiodic(ng), NSperiodic(ng), &
1037!^ & alfaobeta)
1038!^
1039 CALL mp_exchange3d (ng, tile, model, 1, &
1040 & lbi, ubi, lbj, ubj, 0, n(ng), &
1041 & nghostpoints, &
1042 & ewperiodic(ng), nsperiodic(ng), &
1043 & tl_alfaobeta)
1044# endif
1045 CALL mp_exchange2d (ng, tile, model, 2, &
1046 & lbi, ubi, lbj, ubj, &
1047 & nghostpoints, &
1048 & ewperiodic(ng), nsperiodic(ng), &
1049 & alpha, beta)
1050 CALL mp_exchange2d (ng, tile, model, 2, &
1051 & lbi, ubi, lbj, ubj, &
1052 & nghostpoints, &
1053 & ewperiodic(ng), nsperiodic(ng), &
1054 & tl_alpha, tl_beta)
1055# endif
1056
1057# ifdef VAR_RHO_2D_NOT_YET
1058 CALL mp_exchange2d (ng, tile, model, 2, &
1059 & lbi, ubi, lbj, ubj, &
1060 & nghostpoints, &
1061 & ewperiodic(ng), nsperiodic(ng), &
1062 & rhoa, rhos)
1063 CALL mp_exchange2d (ng, tile, model, 2, &
1064 & lbi, ubi, lbj, ubj, &
1065 & nghostpoints, &
1066 & ewperiodic(ng), nsperiodic(ng), &
1067 & tl_rhoa, tl_rhos)
1068# endif
1069
1070# ifdef BV_FREQUENCY_NOT_YET
1071!^ CALL mp_exchange3d (ng, tile, model, 1, &
1072!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1073!^ & NghostPoints, &
1074!^ & EWperiodic(ng), NSperiodic(ng), &
1075!^ & bvf)
1076!^
1077 CALL mp_exchange3d (ng, tile, model, 1, &
1078 & lbi, ubi, lbj, ubj, 0, n(ng), &
1079 & nghostpoints, &
1080 & ewperiodic(ng), nsperiodic(ng), &
1081 & tl_bvf)
1082# endif
1083# endif
1084!
1085 RETURN
1086 END SUBROUTINE rp_rho_eos_tile
1087# endif
1088
1089# ifndef NONLIN_EOS
1090!
1091!***********************************************************************
1092 SUBROUTINE rp_rho_eos_tile (ng, tile, model, &
1093 & LBi, UBi, LBj, UBj, &
1094 & IminS, ImaxS, JminS, JmaxS, &
1095 & nrhs, &
1096# ifdef MASKING
1097 & rmask, &
1098# endif
1099# ifdef VAR_RHO_2D_NOT_YET
1100 & Hz, tl_Hz, &
1101# endif
1102 & z_r, tl_z_r, &
1103 & z_w, tl_z_w, &
1104 & t, tl_t, &
1105# ifdef VAR_RHO_2D_NOT_YET
1106 & rhoA, tl_rhoA, &
1107 & rhoS, tl_rhoS, &
1108# endif
1109# ifdef BV_FREQUENCY_NOT_YET
1110 & tl_bvf, &
1111# endif
1112# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1113 defined BULK_FLUXES
1114 & alpha, tl_alpha, &
1115 & beta, tl_beta, &
1116# ifdef LMD_DDMIX_NOT_YET
1117 & tl_alfaobeta, &
1118# endif
1119# endif
1120# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1121 & tl_pden, &
1122# endif
1123 & rho, tl_rho)
1124!***********************************************************************
1125!
1126 USE mod_param
1127 USE mod_scalars
1128# ifdef SEDIMENT_NOT_YET
1129 USE mod_sediment
1130# endif
1131!
1132 USE exchange_2d_mod
1133 USE exchange_3d_mod
1134# ifdef DISTRIBUTE
1136# endif
1137!
1138! Imported variable declarations.
1139!
1140 integer, intent(in) :: ng, tile, model
1141 integer, intent(in) :: LBi, UBi, LBj, UBj
1142 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
1143 integer, intent(in) :: nrhs
1144!
1145# ifdef ASSUMED_SHAPE
1146# ifdef MASKING
1147 real(r8), intent(in) :: rmask(LBi:,LBj:)
1148# endif
1149# ifdef VAR_RHO_2D_NOT_YET
1150 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
1151# endif
1152 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
1153 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
1154 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
1155# ifdef VAR_RHO_2D_NOT_YET
1156 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
1157# endif
1158 real(r8), intent(in) :: tl_z_r(LBi:,LBj:,:)
1159 real(r8), intent(in) :: tl_z_w(LBi:,LBj:,0:)
1160 real(r8), intent(in) :: tl_t(LBi:,LBj:,:,:,:)
1161# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1162 defined bulk_fluxes
1163 real(r8), intent(inout) :: alpha(LBi:,LBj:)
1164 real(r8), intent(inout) :: beta(LBi:,LBj:)
1165# endif
1166# ifdef VAR_RHO_2D_NOT_YET
1167 real(r8), intent(out) :: rhoA(LBi:,LBj:)
1168 real(r8), intent(out) :: rhoS(LBi:,LBj:)
1169 real(r8), intent(out) :: tl_rhoA(LBi:,LBj:)
1170 real(r8), intent(out) :: tl_rhoS(LBi:,LBj:)
1171# endif
1172# ifdef BV_FREQUENCY_NOT_YET
1173 real(r8), intent(out) :: tl_bvf(LBi:,LBj:,0:)
1174# endif
1175# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1176 defined bulk_fluxes
1177 real(r8), intent(out) :: tl_alpha(LBi:,LBj:)
1178 real(r8), intent(out) :: tl_beta(LBi:,LBj:)
1179# ifdef LMD_DDMIX_NOT_YET
1180 real(r8), intent(out) :: tl_alfaobeta(LBi:,LBj:,0:)
1181# endif
1182# endif
1183# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1184 real(r8), intent(out) :: tl_pden(LBi:,LBj:,:)
1185# endif
1186 real(r8), intent(out) :: rho(LBi:,LBj:,:)
1187 real(r8), intent(out) :: tl_rho(LBi:,LBj:,:)
1188# else
1189# ifdef MASKING
1190 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
1191# endif
1192# ifdef VAR_RHO_2D_NOT_YET
1193 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
1194# endif
1195 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
1196 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
1197 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
1198# ifdef VAR_RHO_2D_NOT_YET
1199 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
1200# endif
1201 real(r8), intent(in) :: tl_z_r(LBi:UBi,LBj:UBj,N(ng))
1202 real(r8), intent(in) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
1203 real(r8), intent(in) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
1204# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1205 defined bulk_fluxes
1206 real(r8), intent(inout) :: alpha(LBi:UBi,LBj:UBj)
1207 real(r8), intent(inout) :: beta(LBi:UBi,LBj:UBj)
1208# endif
1209# ifdef VAR_RHO_2D_NOT_YET
1210 real(r8), intent(out) :: rhoA(LBi:UBi,LBj:UBj)
1211 real(r8), intent(out) :: rhoS(LBi:UBi,LBj:UBj)
1212 real(r8), intent(out) :: tl_rhoA(LBi:UBi,LBj:UBj)
1213 real(r8), intent(out) :: tl_rhoS(LBi:UBi,LBj:UBj)
1214# endif
1215# ifdef BV_FREQUENCY_NOT_YET
1216 real(r8), intent(out) :: tl_bvf(LBi:UBi,LBj:UBj,0:N(ng))
1217# endif
1218# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1219 defined bulk_fluxes
1220 real(r8), intent(out) :: tl_alpha(LBi:UBi,LBj:UBj)
1221 real(r8), intent(out) :: tl_beta(LBi:UBi,LBj:UBj)
1222# ifdef LMD_DDMIX_NOT_YET
1223 real(r8), intent(out) :: tl_alfaobeta(LBi:UBi,LBj:UBj,0:N(ng))
1224# endif
1225# endif
1226# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1227 real(r8), intent(out) :: tl_pden(LBi:UBi,LBj:UBj,N(ng))
1228# endif
1229 real(r8), intent(out) :: rho(LBi:UBi,LBj:UBj,N(ng))
1230 real(r8), intent(out) :: tl_rho(LBi:UBi,LBj:UBj,N(ng))
1231# endif
1232!
1233! Local variable declarations.
1234!
1235 integer :: i, ised, itrc, j, k
1236
1237 real(r8) :: SedDen, cff, cff1, cff2
1238 real(r8) :: tl_SedDen, tl_cff, tl_cff1
1239
1240# include "set_bounds.h"
1241!
1242!=======================================================================
1243! Tangent linear, linear equation of state.
1244!=======================================================================
1245!
1246!-----------------------------------------------------------------------
1247! Compute "in situ" density anomaly (kg/m3 - 1000) using the linear
1248! equation of state.
1249!-----------------------------------------------------------------------
1250!
1251 DO j=jstrt,jendt
1252 DO k=1,n(ng)
1253 DO i=istrt,iendt
1254 rho(i,j,k)=r0(ng)- &
1255 & r0(ng)*tcoef(ng)*(t(i,j,k,nrhs,itemp)-t0(ng))
1256 tl_rho(i,j,k)=-r0(ng)*tcoef(ng)*tl_t(i,j,k,nrhs,itemp)+ &
1257# ifdef TL_IOMS
1258 & r0(ng)+r0(ng)*tcoef(ng)*t0(ng)
1259# endif
1260# ifdef SALINITY
1261 rho(i,j,k)=rho(i,j,k)+ &
1262 & r0(ng)*scoef(ng)*(t(i,j,k,nrhs,isalt)-s0(ng))
1263 tl_rho(i,j,k)=tl_rho(i,j,k)+ &
1264 & r0(ng)*scoef(ng)*tl_t(i,j,k,nrhs,isalt)- &
1265# ifdef TL_IOMS
1266 & r0(ng)*scoef(ng)*s0(ng)
1267# endif
1268# endif
1269# if defined SEDIMENT_NOT_YET && defined SED_DENS_NOT_YET
1270 sedden=0.0_r8
1271 tl_sedden=0.0_r8
1272 DO ised=1,nst
1273 itrc=idsed(ised)
1274 cff1=1.0_r8/srho(ised,ng)
1275 sedden=sedden+ &
1276 & t(i,j,k,nrhs,itrc)* &
1277 & (srho(ised,ng)-rho(i,j,k))*cff1
1278 tl_sedden=tl_sedden+ &
1279 & (tl_t(i,j,k,nrhs,itrc)* &
1280 & (srho(ised,ng)-rho(i,j,k))- &
1281 & t(i,j,k,nrhs,itrc)* &
1282 & tl_rho(i,j,k))*cff1+ &
1283# ifdef TL_IOMS
1284 & t(i,j,k,nrhs,itrc)*rho(i,j,k)*cff1
1285# endif
1286 END DO
1287 rho(i,j,k)=rho(i,j,k)+sedden
1288 tl_rho(i,j,k)=tl_rho(i,j,k)+tl_sedden
1289# endif
1290 rho(i,j,k)=rho(i,j,k)-1000.0_r8
1291# ifdef TL_IOMS
1292 tl_rho(i,j,k)=tl_rho(i,j,k)-1000.0_r8
1293# endif
1294# ifdef MASKING
1295 rho(i,j,k)=rho(i,j,k)*rmask(i,j)
1296 tl_rho(i,j,k)=tl_rho(i,j,k)*rmask(i,j)
1297# endif
1298# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1299!^ pden(i,j,k)=rho(i,j,k)
1300!^
1301 tl_pden(i,j,k)=tl_rho(i,j,k)
1302# endif
1303 END DO
1304 END DO
1305
1306# ifdef VAR_RHO_2D_NOT_YET
1307!
1308!-----------------------------------------------------------------------
1309! Compute vertical averaged density (rhoA) and density perturbation
1310! used (rhoS) in barotropic pressure gradient.
1311!-----------------------------------------------------------------------
1312!
1313 DO i=istrt,iendt
1314 cff1=rho(i,j,n(ng))*hz(i,j,n(ng))
1315 tl_cff1=tl_rho(i,j,n(ng))*hz(i,j,n(ng))+ &
1316 & rho(i,j,n(ng))*tl_hz(i,j,n(ng))- &
1317# ifdef TL_IOMS
1318 & cff1
1319# endif
1320 rhos(i,j)=0.5_r8*cff1*hz(i,j,n(ng))
1321 tl_rhos(i,j)=0.5_r8*(tl_cff1*hz(i,j,n(ng))+ &
1322 & cff1*tl_hz(i,j,n(ng)))- &
1323# ifdef TL_IOMS
1324 & rhos(i,j)
1325# endif
1326 rhoa(i,j)=cff1
1327 tl_rhoa(i,j)=tl_cff1
1328 END DO
1329 DO k=n(ng)-1,1,-1
1330 DO i=istrt,iendt
1331 cff1=rho(i,j,k)*hz(i,j,k)
1332 tl_cff1=tl_rho(i,j,k)*hz(i,j,k)+ &
1333 & rho(i,j,k)*tl_hz(i,j,k)- &
1334# ifdef TL_IOMS
1335 & cff1
1336# endif
1337 rhos(i,j)=rhos(i,j)+hz(i,j,k)*(rhoa(i,j)+0.5_r8*cff1)
1338 tl_rhos(i,j)=tl_rhos(i,j)+ &
1339 & tl_hz(i,j,k)*(rhoa(i,j)+0.5_r8*cff1)+ &
1340 & hz(i,j,k)*(tl_rhoa(i,j)+0.5_r8*tl_cff1)- &
1341# ifdef TL_IOMS
1342 & hz(i,j,k)*(rhoa(i,j)+0.5_r8*cff1)
1343# endif
1344 rhoa(i,j)=rhoa(i,j)+cff1
1345 tl_rhoa(i,j)=tl_rhoa(i,j)+tl_cff1
1346 END DO
1347 END DO
1348 cff2=1.0_r8/rho0
1349 DO i=istrt,iendt
1350 cff1=1.0_r8/(z_w(i,j,n(ng))-z_w(i,j,0))
1351 tl_cff1=-cff1*cff1*(tl_z_w(i,j,n(ng))-tl_z_w(i,j,0))+ &
1352# ifdef TL_IOMS
1353 & 2.0_r8*cff1
1354# endif
1355!
1356! Here we reverse the order of the NL and TL operations since an
1357! intermeridiate value of rhoA and rhoS is needed because they are
1358! recursive.
1359!
1360 tl_rhoa(i,j)=cff2*(tl_cff1*rhoa(i,j)+cff1*tl_rhoa(i,j))
1361 rhoa(i,j)=cff2*cff1*rhoa(i,j)
1362# ifdef TL_IOMS
1363 tl_rhoa(i,j)=tl_rhoa(i,j)-rhoa(i,j)
1364# endif
1365 tl_rhos(i,j)=2.0_r8*cff2* &
1366 & cff1*(2.0_r8*tl_cff1*rhos(i,j)+ &
1367 & cff1*tl_rhos(i,j))
1368 rhos(i,j)=2.0_r8*cff1*cff1*cff2*rhos(i,j)
1369# ifdef TL_IOMS
1370 tl_rhos(i,j)=tl_rhos(i,j)-2.0_r8*rhos(i,j)
1371# endif
1372 END DO
1373# endif
1374
1375# ifdef BV_FREQUENCY_NOT_YET
1376!
1377!-----------------------------------------------------------------------
1378! Compute Brunt-Vaisala frequency (1/s2) at horizontal RHO-points
1379! and vertical W-points.
1380!-----------------------------------------------------------------------
1381!
1382 DO k=1,n(ng)-1
1383 DO i=istrt,iendt
1384!^ bvf(i,j,k)=-gorho0*(rho(i,j,k+1)-rho(i,j,k))/ &
1385!^ & (z_r(i,j,k+1)-z_r(i,j,k))
1386!^
1387 cff=1.0_r8/(z_r(i,j,k+1)-z_r(i,j,k))
1388 tl_cff=-cff*cff*(tl_z_r(i,j,k+1)-tl_z_r(i,j,k))+ &
1389# ifdef TL_IOMS
1390 & 2.0*cff
1391# endif
1392 tl_bvf(i,j,k)=-gorho0* &
1393 & (cff*(tl_rho(i,j,k+1)-tl_rho(i,j,k))+ &
1394 & tl_cff*(rho(i,j,k+1)-rho(i,j,k)))+ &
1395# ifdef TL_IOMS
1396 & gorho0*(rho(i,j,k+1)-rho(i,j,k))*cff
1397# endif
1398 END DO
1399 END DO
1400# endif
1401
1402# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1403 defined bulk_fluxes
1404!
1405!-----------------------------------------------------------------------
1406! Compute thermal expansion (1/Celsius) and saline contraction
1407! (1/PSU) coefficients.
1408!-----------------------------------------------------------------------
1409!
1410 DO i=istrt,iendt
1411 alpha(i,j)=abs(tcoef(ng))
1412# ifdef TL_IOMS
1413 tl_alpha(i,j)=abs(tcoef(ng))
1414# else
1415 tl_alpha(i,j)=0.0_r8
1416# endif
1417# ifdef SALINITY
1418 beta(i,j)=abs(scoef(ng))
1419
1420# ifdef TL_IOMS
1421 tl_beta(i,j)=abs(scoef(ng))
1422# else
1423 tl_beta(i,j)=0.0_r8
1424# endif
1425# else
1426 beta(i,j)=0.0_r8
1427 tl_beta(i,j)=0.0_r8
1428# endif
1429 END DO
1430# ifdef LMD_DDMIX_NOT_YET
1431!
1432! Compute ratio of thermal expansion and saline contraction
1433! coefficients.
1434!
1435 IF (scoef(ng).eq.0.0_r8) THEN
1436 cff=1.0_r8
1437 ELSE
1438 cff=1.0_r8/scoef(ng)
1439 END IF
1440 DO k=1,n(ng)
1441 DO i=istrt,iendt
1442!^ alfaobeta(i,j,k)=cff*Tcoef(ng)
1443!^
1444# ifdef TL_IOMS
1445 tl_alfaobeta(i,j,k)=cff*tcoef(ng)
1446# else
1447 tl_alfaobeta(i,j,k)=0.0_r8
1448# endif
1449 END DO
1450 END DO
1451# endif
1452# endif
1453 END DO
1454!
1455!-----------------------------------------------------------------------
1456! Exchange boundary data.
1457!-----------------------------------------------------------------------
1458!
1459 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1460 CALL exchange_r3d_tile (ng, tile, &
1461 & lbi, ubi, lbj, ubj, 1, n(ng), &
1462 & rho)
1463 CALL exchange_r3d_tile (ng, tile, &
1464 & lbi, ubi, lbj, ubj, 1, n(ng), &
1465 & tl_rho)
1466
1467# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1468!^ CALL exchange_r3d_tile (ng, tile, &
1469!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
1470!^ & pden)
1471!^
1472 CALL exchange_r3d_tile (ng, tile, &
1473 & lbi, ubi, lbj, ubj, 1, n(ng), &
1474 & tl_pden)
1475# endif
1476
1477# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1478 defined bulk_fluxes_not_yet
1479# ifdef LMD_DDMIX_NOT_YET
1480!^ CALL exchange_w3d_tile (ng, tile, &
1481!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1482!^ & alfaobeta)
1483!^
1484 CALL exchange_w3d_tile (ng, tile, &
1485 & lbi, ubi, lbj, ubj, 0, n(ng), &
1486 & tl_alfaobeta)
1487# endif
1488 CALL exchange_r2d_tile (ng, tile, &
1489 & lbi, ubi, lbj, ubj, &
1490 & alpha)
1491 CALL exchange_r2d_tile (ng, tile, &
1492 & lbi, ubi, lbj, ubj, &
1493 & tl_alpha)
1494 CALL exchange_r2d_tile (ng, tile, &
1495 & lbi, ubi, lbj, ubj, &
1496 & beta)
1497 CALL exchange_r2d_tile (ng, tile, &
1498 & lbi, ubi, lbj, ubj, &
1499 & tl_beta)
1500# endif
1501
1502# ifdef VAR_RHO_2D_NOT_YET
1503 CALL exchange_r2d_tile (ng, tile, &
1504 & lbi, ubi, lbj, ubj, &
1505 & rhoa)
1506 CALL exchange_r2d_tile (ng, tile, &
1507 & lbi, ubi, lbj, ubj, &
1508 & tl_rhoa)
1509 CALL exchange_r2d_tile (ng, tile, &
1510 & lbi, ubi, lbj, ubj, &
1511 & rhos)
1512 CALL exchange_r2d_tile (ng, tile, &
1513 & lbi, ubi, lbj, ubj, &
1514 & tl_rhos)
1515# endif
1516
1517# ifdef BV_FREQUENCY_NOT_YET
1518!^ CALL exchange_w3d_tile (ng, tile, &
1519!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1520!^ & bvf)
1521!^
1522 CALL exchange_w3d_tile (ng, tile, &
1523 & lbi, ubi, lbj, ubj, 0, n(ng), &
1524 & tl_bvf)
1525# endif
1526 END IF
1527
1528# ifdef DISTRIBUTE
1529!
1530 CALL mp_exchange3d (ng, tile, model, 1, &
1531 & lbi, ubi, lbj, ubj, 1, n(ng), &
1532 & nghostpoints, &
1533 & ewperiodic(ng), nsperiodic(ng), &
1534 & rho)
1535 CALL mp_exchange3d (ng, tile, model, 1, &
1536 & lbi, ubi, lbj, ubj, 1, n(ng), &
1537 & nghostpoints, &
1538 & ewperiodic(ng), nsperiodic(ng), &
1539 & tl_rho)
1540
1541# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1542!^ CALL mp_exchange3d (ng, tile, model, 1, &
1543!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
1544!^ & NghostPoints, &
1545!^ & EWperiodic(ng), NSperiodic(ng), &
1546!^ & pden)
1547!^
1548 CALL mp_exchange3d (ng, tile, model, 1, &
1549 & lbi, ubi, lbj, ubj, 1, n(ng), &
1550 & nghostpoints, &
1551 & ewperiodic(ng), nsperiodic(ng), &
1552 & tl_pden)
1553# endif
1554
1555# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1556 defined bulk_fluxes
1557# ifdef LMD_DDMIX_NOT_YET
1558!^ CALL mp_exchange3d (ng, tile, model, 1, &
1559!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1560!^ & NghostPoints, &
1561!^ & EWperiodic(ng), NSperiodic(ng), &
1562!^ & alfaobeta)
1563!^
1564 CALL mp_exchange3d (ng, tile, model, 1, &
1565 & lbi, ubi, lbj, ubj, 0, n(ng), &
1566 & nghostpoints, &
1567 & ewperiodic(ng), nsperiodic(ng), &
1568 & tl_alfaobeta)
1569# endif
1570 CALL mp_exchange2d (ng, tile, model, 2, &
1571 & lbi, ubi, lbj, ubj, &
1572 & nghostpoints, &
1573 & ewperiodic(ng), nsperiodic(ng), &
1574 & alpha, beta)
1575 CALL mp_exchange2d (ng, tile, model, 2, &
1576 & lbi, ubi, lbj, ubj, &
1577 & nghostpoints, &
1578 & ewperiodic(ng), nsperiodic(ng), &
1579 & tl_alpha, tl_beta)
1580# endif
1581
1582# ifdef VAR_RHO_2D_NOT_YET
1583 CALL mp_exchange2d (ng, tile, model, 2, &
1584 & lbi, ubi, lbj, ubj, &
1585 & nghostpoints, &
1586 & ewperiodic(ng), nsperiodic(ng), &
1587 & rhoa, rhos)
1588 CALL mp_exchange2d (ng, tile, model, 2, &
1589 & lbi, ubi, lbj, ubj, &
1590 & nghostpoints, &
1591 & ewperiodic(ng), nsperiodic(ng), &
1592 & tl_rhoa, tl_rhos)
1593# endif
1594
1595# ifdef BV_FREQUENCY_NOT_YET
1596!^ CALL mp_exchange3d (ng, tile, model, 1, &
1597!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1598!^ & NghostPoints, &
1599!^ & EWperiodic(ng), NSperiodic(ng), &
1600!^ & bvf)
1601!^
1602 CALL mp_exchange3d (ng, tile, model, 1, &
1603 & lbi, ubi, lbj, ubj, 0, n(ng), &
1604 & nghostpoints, &
1605 & ewperiodic(ng), nsperiodic(ng), &
1606 & tl_bvf)
1607# endif
1608# endif
1609!
1610 RETURN
1611 END SUBROUTINE rp_rho_eos_tile
1612# endif
1613#endif
1614 END MODULE rp_rho_eos_mod
subroutine exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_w3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_coupling), dimension(:), allocatable coupling
real(r8), parameter b03
Definition mod_eoscoef.F:32
real(r8), parameter e02
Definition mod_eoscoef.F:38
real(r8), parameter h02
Definition mod_eoscoef.F:49
real(r8), parameter q03
Definition mod_eoscoef.F:53
real(r8), parameter a03
Definition mod_eoscoef.F:27
real(r8), parameter q04
Definition mod_eoscoef.F:54
real(r8), parameter q05
Definition mod_eoscoef.F:55
real(r8), parameter u01
Definition mod_eoscoef.F:57
real(r8), parameter f00
Definition mod_eoscoef.F:40
real(r8), parameter v01
Definition mod_eoscoef.F:62
real(r8), parameter g01
Definition mod_eoscoef.F:44
real(r8), parameter a04
Definition mod_eoscoef.F:28
real(r8), parameter f01
Definition mod_eoscoef.F:41
real(r8), parameter g00
Definition mod_eoscoef.F:43
real(r8), parameter u02
Definition mod_eoscoef.F:58
real(r8), parameter d01
Definition mod_eoscoef.F:34
real(r8), parameter w00
Definition mod_eoscoef.F:64
real(r8), parameter b00
Definition mod_eoscoef.F:29
real(r8), parameter u03
Definition mod_eoscoef.F:59
real(r8), parameter d00
Definition mod_eoscoef.F:33
real(r8), parameter d02
Definition mod_eoscoef.F:35
real(r8), parameter g03
Definition mod_eoscoef.F:46
real(r8), parameter b01
Definition mod_eoscoef.F:30
real(r8), parameter u04
Definition mod_eoscoef.F:60
real(r8), parameter u00
Definition mod_eoscoef.F:56
real(r8), parameter v02
Definition mod_eoscoef.F:63
real(r8), parameter q01
Definition mod_eoscoef.F:51
real(r8), parameter e00
Definition mod_eoscoef.F:36
real(r8), parameter h00
Definition mod_eoscoef.F:47
real(r8), parameter g02
Definition mod_eoscoef.F:45
real(r8), parameter e03
Definition mod_eoscoef.F:39
real(r8), parameter a00
Definition mod_eoscoef.F:24
real(r8), parameter a01
Definition mod_eoscoef.F:25
real(r8), parameter a02
Definition mod_eoscoef.F:26
real(r8), parameter b02
Definition mod_eoscoef.F:31
real(r8), parameter h01
Definition mod_eoscoef.F:48
real(r8), parameter e01
Definition mod_eoscoef.F:37
real(r8), parameter f02
Definition mod_eoscoef.F:42
real(r8), parameter v00
Definition mod_eoscoef.F:61
real(r8), parameter q02
Definition mod_eoscoef.F:52
real(r8), parameter q00
Definition mod_eoscoef.F:50
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer nghostpoints
Definition mod_param.F:710
integer nst
Definition mod_param.F:521
real(r8), dimension(:), allocatable t0
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(r8), dimension(:), allocatable s0
real(r8), dimension(:), allocatable tcoef
real(dp) gorho0
integer isalt
integer itemp
real(dp) g
real(dp) rho0
real(r8), dimension(:), allocatable r0
real(r8), dimension(:), allocatable scoef
real(r8), dimension(:,:), allocatable srho
integer, dimension(:), allocatable idsed
integer, dimension(:), allocatable nrhs
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public rp_rho_eos(ng, tile, model)
Definition rp_rho_eos.F:48
subroutine rp_rho_eos_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, rmask, z_r, tl_z_r, z_w, tl_z_w, t, tl_t, alpha, tl_alpha, beta, tl_beta, rho, tl_rho)
Definition rp_rho_eos.F:154
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3