ROMS
Loading...
Searching...
No Matches
ad_rho_eos.F
Go to the documentation of this file.
1#include "cppdefs.h"
2
4#if defined ADJOINT && defined SOLVE3D
5!
6!git $Id$
7!================================================== Hernan G. Arango ===
8! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
9! Licensed under a MIT/X style license !
10! See License_ROMS.md !
11!=======================================================================
12! !
13! This routine computes "in situ" density and other associated !
14! quantitites as a function of potential temperature, salinity, !
15! and pressure from a polynomial expression (Jackett and McDougall, !
16! 1992). The polynomial expression was found from fitting to 248 !
17! values in the oceanographic ranges of salinity, potential !
18! temperature, and pressure. It assumes no pressure variation !
19! along geopotential surfaces, that is, depth (meters; negative) !
20! and pressure (dbar; assumed negative here) are interchangeable. !
21! !
22! Check Values: (T=3 C, S=35.5 PSU, Z=-5000 m) !
23! !
24! alpha = 2.1014611551470d-04 (1/Celsius) !
25! beta = 7.2575037309946d-04 (1/PSU) !
26! gamma = 3.9684764511766d-06 (1/Pa) !
27! den = 1050.3639165364 (kg/m3) !
28! den1 = 1028.2845117925 (kg/m3) !
29! sound = 1548.8815240223 (m/s) !
30! bulk = 23786.056026320 (Pa) !
31! !
32! Reference: !
33! !
34! Jackett, D. R. and T. J. McDougall, 1995, Minimal Adjustment of !
35! Hydrostatic Profiles to Achieve Static Stability, J. of Atmos. !
36! and Oceanic Techn., vol. 12, pp. 381-389. !
37! !
38!=======================================================================
39!
40 implicit none
41!
42 PRIVATE
43 PUBLIC :: ad_rho_eos
44!
45 CONTAINS
46!
47!***********************************************************************
48 SUBROUTINE ad_rho_eos (ng, tile, model)
49!***********************************************************************
50!
51 USE mod_param
52 USE mod_parallel
53 USE mod_coupling
54 USE mod_grid
55 USE mod_mixing
56 USE mod_ocean
57 USE mod_stepping
58!
59! Imported variable declarations.
60!
61 integer, intent(in) :: ng, tile, model
62!
63! Local variable declarations.
64!
65 character (len=*), parameter :: myfile = &
66 & __FILE__
67!
68# include "tile.h"
69!
70# ifdef PROFILE
71 CALL wclock_on (ng, model, 14, __line__, myfile)
72# endif
73 CALL ad_rho_eos_tile (ng, tile, model, &
74 & lbi, ubi, lbj, ubj, &
75 & imins, imaxs, jmins, jmaxs, &
76 & nrhs(ng), &
77# ifdef MASKING
78 & grid(ng) % rmask, &
79# endif
80# ifdef VAR_RHO_2D_NOT_YET
81 & grid(ng) % Hz, &
82 & grid(ng) % ad_Hz, &
83# endif
84 & grid(ng) % z_r, &
85 & grid(ng) % ad_z_r, &
86# if defined BV_FREQUENCY_NOT_YET || defined VAR_RHO_2D_NOT_YET
87 & grid(ng) % z_w, &
88 & grid(ng) % ad_z_w, &
89# endif
90 & ocean(ng) % t, &
91 & ocean(ng) % ad_t, &
92# ifdef VAR_RHO_2D_NOT_YET
93 & coupling(ng) % ad_rhoA, &
94 & coupling(ng) % ad_rhoS, &
95# endif
96# ifdef BV_FREQUENCY_NOT_YET
97 & mixing(ng) % ad_bvf, &
98# endif
99# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
100 defined bulk_fluxes_not_yet
101 & mixing(ng) % ad_alpha, &
102 & mixing(ng) % ad_beta, &
103# ifdef LMD_DDMIX_NOT_YET
104 & mixing(ng) % ad_alfaobeta, &
105# endif
106# endif
107# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
108 & ocean(ng) % ad_pden, &
109# endif
110 & ocean(ng) % rho, &
111 & ocean(ng) % ad_rho)
112# ifdef PROFILE
113 CALL wclock_off (ng, model, 14, __line__, myfile)
114# endif
115!
116 RETURN
117 END SUBROUTINE ad_rho_eos
118
119# ifdef NONLIN_EOS
120!
121!***********************************************************************
122 SUBROUTINE ad_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, ad_Hz, &
131# endif
132 & z_r, ad_z_r, &
133# if defined BV_FREQUENCY_NOT_YET || defined VAR_RHO_2D_NOT_YET
134 & z_w, ad_z_w, &
135# endif
136 & t, ad_t, &
137# ifdef VAR_RHO_2D_NOT_YET
138 & ad_rhoA, ad_rhoS, &
139# endif
140# ifdef BV_FREQUENCY_NOT_YET
141 & ad_bvf, &
142# endif
143# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
144 defined BULK_FLUXES_NOT_YET
145 & ad_alpha, ad_beta, &
146# ifdef LMD_DDMIX_NOT_YET
147 & ad_alfaobeta, &
148# endif
149# endif
150# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
151 & ad_pden, &
152# endif
153 & rho, ad_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# if defined BV_FREQUENCY_NOT_YET || defined VAR_RHO_2D_NOT_YET
185 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
186# endif
187 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
188 real(r8), intent(in) :: rho(LBi:,LBj:,:)
189
190# ifdef VAR_RHO_2D_NOT_YET
191 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
192# endif
193 real(r8), intent(inout) :: ad_z_r(LBi:,LBj:,:)
194# if defined BV_FREQUENCY_NOT_YET || defined VAR_RHO_2D_NOT_YET
195 real(r8), intent(inout) :: ad_z_w(LBi:,LBj:,0:)
196# endif
197 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
198# ifdef VAR_RHO_2D_NOT_YET
199 real(r8), intent(inout) :: ad_rhoA(LBi:,LBj:)
200 real(r8), intent(inout) :: ad_rhoS(LBi:,LBj:)
201# endif
202# ifdef BV_FREQUENCY_NOT_YET
203 real(r8), intent(inout) :: ad_bvf(LBi:,LBj:,0:)
204# endif
205# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
206 defined bulk_fluxes_not_yet
207 real(r8), intent(inout) :: ad_alpha(LBi:,LBj:)
208 real(r8), intent(inout) :: ad_beta(LBi:,LBj:)
209# ifdef LMD_DDMIX_NOT_YET
210 real(r8), intent(inout) :: ad_alfaobeta(LBi:,LBj:,0:)
211# endif
212# endif
213# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
214 real(r8), intent(inout) :: ad_pden(LBi:,LBj:,:)
215# endif
216 real(r8), intent(inout) :: ad_rho(LBi:,LBj:,:)
217# else
218# ifdef MASKING
219 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
220# endif
221# ifdef VAR_RHO_2D_NOT_YET
222 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
223# endif
224 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
225# if defined BV_FREQUENCY_NOT_YET || defined VAR_RHO_2D_NOT_YET
226 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
227# endif
228 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
229 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
230
231# ifdef VAR_RHO_2D_NOT_YET
232 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
233# endif
234 real(r8), intent(inout) :: ad_z_r(LBi:UBi,LBj:UBj,N(ng))
235# if defined BV_FREQUENCY_NOT_YET || defined VAR_RHO_2D_NOT_YET
236 real(r8), intent(inout) :: ad_z_w(LBi:UBi,LBj:UBj,0:N(ng))
237# endif
238 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
239
240# ifdef VAR_RHO_2D_NOT_YET
241 real(r8), intent(inout) :: ad_rhoA(LBi:UBi,LBj:UBj)
242 real(r8), intent(inout) :: ad_rhoS(LBi:UBi,LBj:UBj)
243# endif
244# ifdef BV_FREQUENCY_NOT_YET
245 real(r8), intent(inout) :: ad_bvf(LBi:UBi,LBj:UBj,0:N(ng))
246# endif
247# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
248 defined bulk_fluxes_not_yet
249 real(r8), intent(inout) :: ad_alpha(LBi:UBi,LBj:UBj)
250 real(r8), intent(inout) :: ad_beta(LBi:UBi,LBj:UBj)
251# ifdef LMD_DDMIX_NOT_YET
252 real(r8), intent(inout) :: ad_alfaobeta(LBi:UBi,LBj:UBj,0:N(ng))
253# endif
254# endif
255# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
256 real(r8), intent(inout) :: ad_pden(LBi:UBi,LBj:UBj,N(ng))
257# endif
258 real(r8), intent(inout) :: ad_rho(LBi:UBi,LBj:UBj,N(ng))
259# endif
260!
261! Local variable declarations.
262!
263 integer :: i, ised, itrc, j, k
264
265 real(r8) :: SedDen, Tp, Tpr10, Ts, Tt, sqrtTs
266 real(r8) :: ad_SedDen, ad_Tp, ad_Tpr10, ad_Ts, ad_Tt, ad_sqrtTs
267# ifdef BV_FREQUENCY_NOT_YET
268 real(r8) :: bulk_dn, bulk_up, den_dn, den_up
269 real(r8) :: ad_bulk_dn, ad_bulk_up, ad_den_dn, ad_den_up
270# endif
271 real(r8) :: cff, cff1, cff2, cff3
272 real(r8) :: ad_cff, ad_cff1, ad_cff2, ad_cff3
273 real(r8) :: adfac, adfac1, adfac2, adfac3
274
275 real(r8), dimension(0:9) :: C
276 real(r8), dimension(0:9) :: ad_C
277# ifdef EOS_TDERIVATIVE
278 real(r8), dimension(0:9) :: dCdT(0:9)
279 real(r8), dimension(0:9) :: ad_dCdT(0:9)
280 real(r8), dimension(0:9) :: d2Cd2T(0:9)
281
282 real(r8), dimension(IminS:ImaxS,N(ng)) :: DbulkDS
283 real(r8), dimension(IminS:ImaxS,N(ng)) :: DbulkDT
284 real(r8), dimension(IminS:ImaxS,N(ng)) :: Dden1DS
285 real(r8), dimension(IminS:ImaxS,N(ng)) :: Dden1DT
286 real(r8), dimension(IminS:ImaxS,N(ng)) :: Scof
287 real(r8), dimension(IminS:ImaxS,N(ng)) :: Tcof
288 real(r8), dimension(IminS:ImaxS,N(ng)) :: wrk
289
290 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_DbulkDS
291 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_DbulkDT
292 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_Dden1DS
293 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_Dden1DT
294 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_Scof
295 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_Tcof
296 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_wrk
297# endif
298 real(r8), dimension(IminS:ImaxS,N(ng)) :: bulk
299 real(r8), dimension(IminS:ImaxS,N(ng)) :: bulk0
300 real(r8), dimension(IminS:ImaxS,N(ng)) :: bulk1
301 real(r8), dimension(IminS:ImaxS,N(ng)) :: bulk2
302 real(r8), dimension(IminS:ImaxS,N(ng)) :: den
303 real(r8), dimension(IminS:ImaxS,N(ng)) :: den1
304
305 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_bulk
306 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_bulk0
307 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_bulk1
308 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_bulk2
309 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_den
310 real(r8), dimension(IminS:ImaxS,N(ng)) :: ad_den1
311
312# ifdef VAR_RHO_2D_NOT_YET
313 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: rhoA1
314 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: rhoS1
315# endif
316
317# include "set_bounds.h"
318!
319!-------------------------------------------------------------------------
320! Initialize adjoint private variables.
321!-------------------------------------------------------------------------
322!
323 ad_tt=0.0_r8
324 ad_ts=0.0_r8
325 ad_tp=0.0_r8
326 ad_tpr10=0.0_r8
327# ifdef BV_FREQUENCY_NOT_YET
328 ad_bulk_dn=0.0_r8
329 ad_bulk_up=0.0_r8
330 ad_den_dn=0.0_r8
331 ad_den_up=0.0_r8
332# endif
333 ad_sqrtts=0.0_r8
334 ad_cff=0.0_r8
335 ad_cff1=0.0_r8
336 ad_cff2=0.0_r8
337 ad_cff3=0.0_r8
338
339 ad_c=0.0_r8
340 ad_dcdt=0.0_r8
341
342 DO k=1,n(ng)
343 DO i=imins,imaxs
344# ifdef EOS_TDERIVATIVE
345 ad_dbulkds(i,k)=0.0_r8
346 ad_dbulkdt(i,k)=0.0_r8
347 ad_dden1ds(i,k)=0.0_r8
348 ad_dden1dt(i,k)=0.0_r8
349 ad_scof(i,k)=0.0_r8
350 ad_tcof(i,k)=0.0_r8
351 ad_wrk(i,k)=0.0_r8
352# endif
353 ad_bulk(i,k)=0.0_r8
354 ad_bulk0(i,k)=0.0_r8
355 ad_bulk1(i,k)=0.0_r8
356 ad_bulk2(i,k)=0.0_r8
357 ad_den(i,k)=0.0_r8
358 ad_den1(i,k)=0.0_r8
359 END DO
360 END DO
361!
362!=======================================================================
363! Adjoint nonlinear equation of state. Notice that this equation
364! of state is only valid for potential temperature range of -2C to 40C
365! and a salinity range of 0 PSU to 42 PSU.
366!=======================================================================
367!
368!-----------------------------------------------------------------------
369! Exchange boundary data.
370!-----------------------------------------------------------------------
371!
372# ifdef DISTRIBUTE
373# ifdef BV_FREQUENCY_NOT_YET
374!^ CALL mp_exchange3d (ng, tile, model, 1, &
375!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
376!^ & NghostPoints, &
377!^ & EWperiodic(ng), NSperiodic(ng), &
378!^ & tl_bvf)
379!^
380 CALL ad_mp_exchange3d (ng, tile, model, 1, &
381 & lbi, ubi, lbj, ubj, 0, n(ng), &
382 & nghostpoints, &
383 & ewperiodic(ng), nsperiodic(ng), &
384 & ad_bvf)
385# endif
386# ifdef VAR_RHO_2D_NOT_YET
387!^ CALL mp_exchange2d (ng, tile, model, 2, &
388!^ & LBi, UBi, LBj, UBj, &
389!^ & NghostPoints, &
390!^ & EWperiodic(ng), NSperiodic(ng), &
391!^ & tl_rhoA, tl_rhoS)
392!^
393 CALL ad_mp_exchange2d (ng, tile, model, 2, &
394 & lbi, ubi, lbj, ubj, &
395 & nghostpoints, &
396 & ewperiodic(ng), nsperiodic(ng), &
397 & ad_rhoa, ad_rhos)
398# endif
399# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
400 defined bulk_fluxes_not_yet
401!^ CALL mp_exchange2d (ng, tile, model, 2, &
402!^ & LBi, UBi, LBj, UBj, &
403!^ & NghostPoints, &
404!^ & EWperiodic(ng), NSperiodic(ng), &
405!^ & tl_alpha, tl_beta)
406!^
407 CALL ad_mp_exchange2d (ng, tile, model, 2, &
408 & lbi, ubi, lbj, ubj, &
409 & nghostpoints, &
410 & ewperiodic(ng), nsperiodic(ng), &
411 & ad_alpha, ad_beta)
412# ifdef LMD_DDMIX_NOT_YET
413!^ CALL mp_exchange3d (ng, tile, model, 1, &
414!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
415!^ & NghostPoints, &
416!^ & EWperiodic(ng), NSperiodic(ng), &
417!^ & tl_alfaobeta)
418!^
419 CALL ad_mp_exchange3d (ng, tile, model, 1, &
420 & lbi, ubi, lbj, ubj, 0, n(ng), &
421 & nghostpoints, &
422 & ewperiodic(ng), nsperiodic(ng), &
423 & ad_alfaobeta)
424# endif
425# endif
426# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
427!^ CALL mp_exchange3d (ng, tile, model, 1, &
428!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
429!^ & NghostPoints, &
430!^ & EWperiodic(ng), NSperiodic(ng), &
431!^ & tl_pden)
432!^
433 CALL ad_mp_exchange3d (ng, tile, model, 1, &
434 & lbi, ubi, lbj, ubj, 1, n(ng), &
435 & nghostpoints, &
436 & ewperiodic(ng), nsperiodic(ng), &
437 & ad_pden)
438# endif
439!^ CALL mp_exchange3d (ng, tile, model, 1, &
440!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
441!^ & NghostPoints, &
442!^ & EWperiodic(ng), NSperiodic(ng), &
443!^ & tl_rho)
444!^
445 CALL ad_mp_exchange3d (ng, tile, model, 1, &
446 & lbi, ubi, lbj, ubj, 1, n(ng), &
447 & nghostpoints, &
448 & ewperiodic(ng), nsperiodic(ng), &
449 & ad_rho)
450!
451# endif
452
453 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
454# ifdef BV_FREQUENCY_NOT_YET
455!^ CALL exchange_w3d_tile (ng, tile, &
456!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
457!^ & tl_bvf)
458!^
459 CALL ad_exchange_w3d_tile (ng, tile, &
460 & lbi, ubi, lbj, ubj, 0, n(ng), &
461 & ad_bvf)
462# endif
463# ifdef VAR_RHO_2D_NOT_YET
464!^ CALL exchange_r2d_tile (ng, tile, &
465!^ & LBi, UBi, LBj, UBj, &
466!^ & tl_rhoS)
467!^
468 CALL ad_exchange_r2d_tile (ng, tile, &
469 & lbi, ubi, lbj, ubj, &
470 & ad_rhos)
471!^ CALL exchange_r2d_tile (ng, tile, &
472!^ & LBi, UBi, LBj, UBj, &
473!^ & tl_rhoA)
474!^
475 CALL ad_exchange_r2d_tile (ng, tile, &
476 & lbi, ubi, lbj, ubj, &
477 & ad_rhoa)
478# endif
479# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
480 defined bulk_fluxes_not_yet
481!^ CALL exchange_r2d_tile (ng, tile, &
482!^ & LBi, UBi, LBj, UBj, &
483!^ & tl_beta)
484!^
485 CALL ad_exchange_r2d_tile (ng, tile, &
486 & lbi, ubi, lbj, ubj, &
487 & ad_beta)
488!^ CALL exchange_r2d_tile (ng, tile, &
489!^ & LBi, UBi, LBj, UBj, &
490!^ & tl_alpha)
491!^
492 CALL ad_exchange_r2d_tile (ng, tile, &
493 & lbi, ubi, lbj, ubj, &
494 & ad_alpha)
495# ifdef LMD_DDMIX_NOT_YET
496!^ CALL exchange_w3d_tile (ng, tile, &
497!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
498!^ & tl_alfaobeta)
499!^
500 CALL ad_exchange_w3d_tile (ng, tile, &
501 & lbi, ubi, lbj, ubj, 0, n(ng), &
502 & ad_alfaobeta)
503# endif
504# endif
505# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
506!^ CALL exchange_r3d_tile (ng, tile, &
507!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
508!^ & tl_pden)
509!^
510 CALL ad_exchange_r3d_tile (ng, tile, &
511 & lbi, ubi, lbj, ubj, 1, n(ng), &
512 & ad_pden)
513# endif
514!^ CALL exchange_r3d_tile (ng, tile, &
515!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
516!^ & tl_rho)
517!^
518 CALL ad_exchange_r3d_tile (ng, tile, &
519 & lbi, ubi, lbj, ubj, 1, n(ng), &
520 & ad_rho)
521 END IF
522!
523!-----------------------------------------------------------------------
524! Compute BASIC STATE related variables.
525!-----------------------------------------------------------------------
526!
527 DO j=jstrt,jendt
528 DO k=1,n(ng)
529 DO i=istrt,iendt
530 tt=max(-2.0_r8,t(i,j,k,nrhs,itemp))
531# ifdef SALINITY
532 ts=max(0.0_r8,t(i,j,k,nrhs,isalt))
533 sqrtts=sqrt(ts)
534# else
535 ts=0.0_r8
536 sqrtts=0.0_r8
537# endif
538 tp=z_r(i,j,k)
539 tpr10=0.1_r8*tp
540!
541! Compute local nonlinear equation of state coefficients and their
542! derivatives when appropriate.
543!
544 c(0)=q00+tt*(q01+tt*(q02+tt*(q03+tt*(q04+tt*q05))))
545 c(1)=u00+tt*(u01+tt*(u02+tt*(u03+tt*u04)))
546 c(2)=v00+tt*(v01+tt*v02)
547 c(3)=a00+tt*(a01+tt*(a02+tt*(a03+tt*a04)))
548 c(4)=b00+tt*(b01+tt*(b02+tt*b03))
549 c(5)=d00+tt*(d01+tt*d02)
550 c(6)=e00+tt*(e01+tt*(e02+tt*e03))
551 c(7)=f00+tt*(f01+tt*f02)
552 c(8)=g01+tt*(g02+tt*g03)
553 c(9)=h00+tt*(h01+tt*h02)
554# ifdef EOS_TDERIVATIVE
555!
556 dcdt(0)=q01+tt*(2.0_r8*q02+tt*(3.0_r8*q03+tt*(4.0_r8*q04+ &
557 & tt*5.0_r8*q05)))
558 dcdt(1)=u01+tt*(2.0_r8*u02+tt*(3.0_r8*u03+tt*4.0_r8*u04))
559 dcdt(2)=v01+tt*2.0_r8*v02
560 dcdt(3)=a01+tt*(2.0_r8*a02+tt*(3.0_r8*a03+tt*4.0_r8*a04))
561 dcdt(4)=b01+tt*(2.0_r8*b02+tt*3.0_r8*b03)
562 dcdt(5)=d01+tt*2.0_r8*d02
563 dcdt(6)=e01+tt*(2.0_r8*e02+tt*3.0_r8*e03)
564 dcdt(7)=f01+tt*2.0_r8*f02
565 dcdt(8)=g02+tt*2.0_r8*g03
566 dcdt(9)=h01+tt*2.0_r8*h02
567!
568 d2cd2t(0)=2.0_r8*q02+tt*(6.0_r8*q03+tt*(12.0_r8*q04+ &
569 & tt*20.0_r8*q05))
570 d2cd2t(1)=2.0_r8*u02+tt*(6.0_r8*u03+tt*12.0_r8*u04)
571 d2cd2t(2)=2.0_r8*v02
572 d2cd2t(3)=2.0_r8*a02+tt*(6.0_r8*a03+tt*12.0_r8*a04)
573 d2cd2t(4)=2.0_r8*b02+tt*6.0_r8*b03
574 d2cd2t(5)=2.0_r8*d02
575 d2cd2t(6)=2.0_r8*e02+tt*6.0_r8*e03
576 d2cd2t(7)=2.0_r8*f02
577 d2cd2t(8)=2.0_r8*g03
578 d2cd2t(9)=2.0_r8*h02
579# endif
580!
581! Compute BASIC STATE density (kg/m3) at standard one atmosphere
582! pressure.
583!
584 den1(i,k)=c(0)+ts*(c(1)+sqrtts*c(2)+ts*w00)
585
586# ifdef EOS_TDERIVATIVE
587!
588! Compute BASIC STATE d(den1)/d(S) and d(den1)/d(T) derivatives used
589! in the computation of thermal expansion and saline contraction
590! coefficients.
591!
592 dden1ds(i,k)=c(1)+1.5_r8*c(2)*sqrtts+2.0_r8*w00*ts
593 dden1dt(i,k)=dcdt(0)+ts*(dcdt(1)+sqrtts*dcdt(2))
594# endif
595!
596! Compute BASIC STATE secant bulk modulus.
597!
598 bulk0(i,k)=c(3)+ts*(c(4)+sqrtts*c(5))
599 bulk1(i,k)=c(6)+ts*(c(7)+sqrtts*g00)
600 bulk2(i,k)=c(8)+ts*c(9)
601 bulk(i,k)=bulk0(i,k)-tp*(bulk1(i,k)-tp*bulk2(i,k))
602
603# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
604 defined bulk_fluxes_not_yet
605!
606! Compute BASIC STATE d(bulk)/d(S) and d(bulk)/d(T) derivatives used
607! in the computation of thermal expansion and saline contraction
608! coefficients.
609!
610 dbulkds(i,k)=c(4)+sqrtts*1.5_r8*c(5)- &
611 & tp*(c(7)+sqrtts*1.5_r8*g00-tp*c(9))
612 dbulkdt(i,k)=dcdt(3)+ts*(dcdt(4)+sqrtts*dcdt(5))- &
613 & tp*(dcdt(6)+ts*dcdt(7)- &
614 & tp*(dcdt(8)+ts*dcdt(9)))
615# endif
616!
617! Compute local "in situ" density anomaly (kg/m3 - 1000). The (i,k)
618! DO-loop is closed here because of the adjoint to facilitate vertical
619! integrals of the BASIC STATE.
620!
621 cff=1.0_r8/(bulk(i,k)+tpr10)
622 den(i,k)=den1(i,k)*bulk(i,k)*cff
623# if defined SEDIMENT_NOT_YET && defined SED_DENS_NOT_YET
624 sedden=0.0_r8
625 DO ised=1,nst
626 cff1=1.0_r8/srho(ised,ng)
627 sedden=sedden+ &
628 & t(i,j,k,nrhs,idsed(ised))* &
629 & (srho(ised,ng)-den(i,k))*cff1
630 END DO
631 den(i,k)=den(i,k)+sedden
632# endif
633 den(i,k)=den(i,k)-1000.0_r8
634# ifdef MASKING
635 den(i,k)=den(i,k)*rmask(i,j)
636# endif
637 END DO
638 END DO
639
640# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
641 defined bulk_fluxes_not_yet
642!
643!-----------------------------------------------------------------------
644! Compute BASIC STATE thermal expansion (1/Celsius) and saline
645! contraction (1/PSU) coefficients.
646!-----------------------------------------------------------------------
647!
648# ifdef LMD_DDMIX_NOT_YET
649 DO k=1,n(ng)
650# else
651 DO k=n(ng),n(ng)
652# endif
653 DO i=istrt,iendt
654 tpr10=0.1_r8*z_r(i,j,k)
655!
656! Compute thermal expansion and saline contraction coefficients.
657!
658 cff=bulk(i,k)+tpr10
659 cff1=tpr10*den1(i,k)
660 cff2=bulk(i,k)*cff
661 wrk(i,k)=(den(i,k)+1000.0_r8)*cff*cff
662 tcof(i,k)=-(dbulkdt(i,k)*cff1+ &
663 & dden1dt(i,k)*cff2)
664 scof(i,k)= (dbulkds(i,k)*cff1+ &
665 & dden1ds(i,k)*cff2)
666 END DO
667 END DO
668# endif
669!
670!-----------------------------------------------------------------------
671! Load adjoint "in situ" density anomaly (kg/m3 - 1000) and adjoint
672! potential density anomaly (kg/m3 - 1000) referenced to the surface
673! into global arrays.
674!-----------------------------------------------------------------------
675!
676 DO k=1,n(ng)
677 DO i=istrt,iendt
678# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
679# ifdef MASKING
680!^ tl_pden(i,j,k)=tl_pden(i,k)*rmask(i,j)
681!^
682 ad_pden(i,j,k)=ad_pden(i,k)*rmask(i,j)
683# endif
684!^ tl_pden(i,j,k)=tl_den1(i,k) ! This gives a fatal
685!^ ! result in 4D-Var
686 ad_den1(i,k)=ad_den1(i,k)+ad_pden(i,j,k)! posterior error...
687 ad_pden(i,j,k)=0.0_r8
688# endif
689!^ tl_rho(i,j,k)=tl_den(i,k)
690!^
691 ad_den(i,k)=ad_den(i,k)+ad_rho(i,j,k)
692 ad_rho(i,j,k)=0.0_r8
693 END DO
694 END DO
695
696# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
697 defined bulk_fluxes_not_yet
698!
699!-----------------------------------------------------------------------
700! Compute adjoint thermal expansion (1/Celsius) and adjoint saline
701! contraction (1/PSU) coefficients.
702!-----------------------------------------------------------------------
703!
704# ifdef LMD_DDMIX_NOT_YET
705!^ DO k=1,N(ng)
706!^
707 DO k=n(ng),1,-1
708# else
709 DO k=n(ng),n(ng)
710# endif
711 IF (k.eq.n(ng)) THEN
712 DO i=istrt,iendt
713 cff=1.0_r8/wrk(i,n(ng))
714!^ tl_beta (i,j)=tl_cff*Scof(i,N(ng))+cff*tl_Scof(i,N(ng))
715!^ tl_alpha(i,j)=tl_cff*Tcof(i,N(ng))+cff*tl_Tcof(i,N(ng))
716!^
717 ad_scof(i,n(ng))=ad_scof(i,n(ng))+cff*ad_beta(i,j)
718 ad_tcof(i,n(ng))=ad_tcof(i,n(ng))+cff*ad_alpha(i,j)
719 ad_cff=ad_beta(i,j)*scof(i,n(ng))+ &
720 & ad_alpha(i,j)*tcof(i,n(ng))
721 ad_beta(i,j)=0.0_r8
722 ad_alpha(i,j)=0.0_r8
723!^ tl_cff=-cff*cff*tl_wrk(i,N(ng))
724!^
725 ad_wrk(i,n(ng))=ad_wrk(i,n(ng))-cff*cff*ad_cff
726 ad_cff=0.0_r8
727 END DO
728 END IF
729 DO i=istrt,iendt
730 tpr10=0.1_r8*z_r(i,j,k)
731 cff=bulk(i,k)+tpr10
732 cff1=tpr10*den1(i,k)
733 cff2=bulk(i,k)*cff
734# ifdef LMD_DDMIX_NOT_YET
735!^ tl_alfaobeta(i,j,k)=(tl_Tcof(i,k)*Scof(i,k)- &
736!^ & Tcof(i,k)*tl_Scof(i,k))/ &
737!^ & (Scof(i,k)*Scof(i,k))
738!^
739 adfac=ad_alfaobeta(i,j,k)/(scof(i,k)*scof(i,k))
740 ad_tcof(i,k)=ad_tcof(i,k)+scof(i,k)*adfac
741 ad_scof(i,k)=ad_scof(i,k)-tcof(i,k)*adfac
742 ad_alfaobeta(i,j,k)=0.0_r8
743# endif
744!^ tl_Scof(i,k)= (tl_DbulkDS(i,k)*cff1+ &
745!^ & DbulkDS(i,k)*tl_cff1+ &
746!^ & tl_Dden1DS(i,k)*cff2+ &
747!^ & Dden1DS(i,k)*tl_cff2)
748!^
749 ad_dbulkds(i,k)=ad_dbulkds(i,k)+ad_scof(i,k)*cff1
750 ad_dden1ds(i,k)=ad_dden1ds(i,k)+ad_scof(i,k)*cff2
751 ad_cff1=dbulkds(i,k)*ad_scof(i,k)
752 ad_cff2=dden1ds(i,k)*ad_scof(i,k)
753 ad_scof(i,k)=0.0_r8
754!^ tl_Tcof(i,k)=-(tl_DbulkDT(i,k)*cff1+ &
755!^ & DbulkDT(i,k)*tl_cff1+ &
756!^ & tl_Dden1DT(i,k)*cff2+ &
757!^ & Dden1DT(i,k)*tl_cff2)
758!^
759 ad_dbulkdt(i,k)=ad_dbulkdt(i,k)-ad_tcof(i,k)*cff1
760 ad_dden1dt(i,k)=ad_dden1dt(i,k)-ad_tcof(i,k)*cff2
761 ad_cff1=ad_cff1-dbulkdt(i,k)*ad_tcof(i,k)
762 ad_cff2=ad_cff2-dden1dt(i,k)*ad_tcof(i,k)
763 ad_tcof(i,k)=0.0_r8
764!^ tl_wrk(i,k)=cff*(cff*tl_den(i,k)+ &
765!^ & 2.0_r8*tl_cff*(den(i,k)+1000.0_r8))
766!^
767 adfac=cff*ad_wrk(i,k)
768 ad_den(i,k)=ad_den(i,k)+cff*adfac
769 ad_cff=ad_cff+2.0_r8*(den(i,k)+1000.0_r8)*adfac
770 ad_wrk(i,k)=0.0_r8
771!^ tl_cff2=tl_bulk(i,k)*cff+bulk(i,k)*tl_cff
772!^
773 ad_bulk(i,k)=ad_bulk(i,k)+ad_cff2*cff
774 ad_cff=ad_cff+bulk(i,k)*ad_cff2
775 ad_cff2=0.0_r8
776!^ tl_cff1=tl_Tpr10*den1(i,k)+Tpr10*tl_den1(i,k)
777!^
778 ad_tpr10=ad_tpr10+ad_cff1*den1(i,k)
779 ad_den1(i,k)=ad_den1(i,k)+tpr10*ad_cff1
780 ad_cff1=0.0_r8
781!^ tl_cff=tl_bulk(i,k)+tl_Tpr10
782!^
783 ad_bulk(i,k)=ad_bulk(i,k)+ad_cff
784 ad_tpr10=ad_tpr10+ad_cff
785 ad_cff=0.0_r8
786!^ tl_Tpr10=0.1_r8*tl_z_r(i,j,k)
787!^
788 ad_z_r(i,j,k)=ad_z_r(i,j,k)+0.1_r8*ad_tpr10
789 ad_tpr10=0.0_r8
790 END DO
791 END DO
792# endif
793
794# if defined BV_FREQUENCY_NOT_YET
795!
796!-----------------------------------------------------------------------
797! Compute adjoint Brunt-Vaisala frequency (1/s2) at horizontal
798! RHO-points and vertical W-points.
799!-----------------------------------------------------------------------
800!
801 DO i=istrt,iendt
802!^ tl_bvf(i,j,N(ng))=0.0_r8
803!^
804 ad_bvf(i,j,n(ng))=0.0_r8
805!^ tl_bvf(i,j,0)=0.0_r8
806!^
807 ad_bvf(i,j,0)=0.0_r8
808 END DO
809 DO k=1,n(ng)-1
810 DO i=istrt,iendt
811 bulk_up=bulk0(i,k+1)- &
812 & z_w(i,j,k)*(bulk1(i,k+1)- &
813 & bulk2(i,k+1)*z_w(i,j,k))
814 bulk_dn=bulk0(i,k )- &
815 & z_w(i,j,k)*(bulk1(i,k )- &
816 & bulk2(i,k )*z_w(i,j,k))
817 cff1=1.0_r8/(bulk_up+0.1_r8*z_w(i,j,k))
818 cff2=1.0_r8/(bulk_dn+0.1_r8*z_w(i,j,k))
819 den_up=cff1*(den1(i,k+1)*bulk_up)
820 den_dn=cff2*(den1(i,k )*bulk_dn)
821 cff3=1.0_r8/(0.5_r8*(den_up+den_dn)* &
822 & (z_r(i,j,k+1)-z_r(i,j,k)))
823!^ tl_bvf(i,j,k)=-g*((tl_den_up-tl_den_dn)*cff3+ &
824!^ & (den_up-den_dn)*tl_cff3)
825!^
826 adfac=-g*ad_bvf(i,j,k)
827 adfac1=adfac*cff3
828 ad_cff3=ad_cff3+(den_up-den_dn)*adfac
829 ad_den_up=ad_den_up+adfac1
830 ad_den_dn=ad_den_dn-adfac1
831 ad_bvf(i,j,k)=0.0_r8
832!^ tl_cff3=-cff3*cff3* &
833!^ & 0.5_r8*((tl_den_up+tl_den_dn)* &
834!^ & (z_r(i,j,k+1)-z_r(i,j,k))+ &
835!^ & (den_up+den_dn)* &
836!^ & (tl_z_r(i,j,k+1)-tl_z_r(i,j,k)))
837!^
838 adfac=-cff3*cff3*0.5_r8*ad_cff3
839 adfac1=adfac*(z_r(i,j,k+1)-z_r(i,j,k))
840 adfac2=adfac*(den_up+den_dn)
841 ad_z_r(i,j,k )=ad_z_r(i,j,k )-adfac2
842 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)+adfac2
843 ad_den_up=ad_den_up+adfac1
844 ad_den_dn=ad_den_dn+adfac1
845 ad_cff3=0.0_r8
846!^ tl_den_dn=tl_cff2*(den1(i,k )*bulk_dn)+ &
847!^ & cff2*(tl_den1(i,k )*bulk_dn+ &
848!^ & den1(i,k )*tl_bulk_dn)
849!^ tl_den_up=tl_cff1*(den1(i,k+1)*bulk_up)+ &
850!^ & cff1*(tl_den1(i,k+1)*bulk_up+ &
851!^ & den1(i,k+1)*tl_bulk_up)
852!^
853 adfac1=cff2*ad_den_dn
854 adfac2=cff1*ad_den_up
855 ad_cff2=ad_cff2+(den1(i,k )*bulk_dn)*ad_den_dn
856 ad_cff1=ad_cff1+(den1(i,k+1)*bulk_up)*ad_den_up
857 ad_den1(i,k )=ad_den1(i,k )+bulk_dn*adfac1
858 ad_den1(i,k+1)=ad_den1(i,k+1)+bulk_up*adfac2
859 ad_bulk_dn=ad_bulk_dn+den1(i,k )*adfac1
860 ad_bulk_up=ad_bulk_up+den1(i,k+1)*adfac2
861 ad_den_dn=0.0_r8
862 ad_den_up=0.0_r8
863!^ tl_cff2=-cff2*cff2*(tl_bulk_dn+0.1_r8*tl_z_w(i,j,k))
864!^ tl_cff1=-cff1*cff1*(tl_bulk_up+0.1_r8*tl_z_w(i,j,k))
865!^
866 adfac1=-cff2*cff2*ad_cff2
867 adfac2=-cff1*cff1*ad_cff1
868 ad_bulk_dn=ad_bulk_dn+adfac1
869 ad_bulk_up=ad_bulk_up+adfac2
870 ad_z_w(i,j,k)=ad_z_w(i,j,k)+ &
871 & 0.1_r8*(adfac1+adfac2)
872 ad_cff2=0.0_r8
873 ad_cff1=0.0_r8
874!^ tl_bulk_dn=tl_bulk0(i,k )- &
875!^ & tl_z_w(i,j,k)*(bulk1(i,k )- &
876!^ & bulk2(i,k )*z_w(i,j,k))- &
877!^ & z_w(i,j,k)*(tl_bulk1(i,k )- &
878!^ & tl_bulk2(i,k )*z_w(i,j,k)- &
879!^ & bulk2(i,k )*tl_z_w(i,j,k))
880!^ tl_bulk_up=tl_bulk0(i,k+1)- &
881!^ & tl_z_w(i,j,k)*(bulk1(i,k+1)- &
882!^ & bulk2(i,k+1)*z_w(i,j,k))- &
883!^ & z_w(i,j,k)*(tl_bulk1(i,k+1)- &
884!^ & tl_bulk2(i,k+1)*z_w(i,j,k)- &
885!^ & bulk2(i,k+1)*tl_z_w(i,j,k))
886!^
887 adfac1=z_w(i,j,k)*ad_bulk_dn
888 adfac2=z_w(i,j,k)*ad_bulk_up
889 ad_bulk0(i,k )=ad_bulk0(i,k )+ad_bulk_dn
890 ad_bulk0(i,k+1)=ad_bulk0(i,k+1)+ad_bulk_up
891 ad_z_w(i,j,k)=ad_z_w(i,j,k)- &
892 & (bulk1(i,k )- &
893 & bulk2(i,k )*z_w(i,j,k)- &
894 & bulk2(i,k ))*ad_bulk_dn- &
895 & (bulk1(i,k+1)- &
896 & bulk2(i,k+1)*z_w(i,j,k)- &
897 & bulk2(i,k+1))*ad_bulk_up
898 ad_bulk1(i,k )=ad_bulk1(i,k )-adfac1
899 ad_bulk1(i,k+1)=ad_bulk1(i,k+1)-adfac2
900 ad_bulk2(i,k )=ad_bulk2(i,k )+z_w(i,j,k)*adfac1
901 ad_bulk2(i,k+1)=ad_bulk2(i,k+1)+z_w(i,j,k)*adfac2
902 ad_bulk_dn=0.0_r8
903 ad_bulk_up=0.0_r8
904 END DO
905 END DO
906# endif
907
908# ifdef VAR_RHO_2D_NOT_YET
909!
910!-----------------------------------------------------------------------
911! Compute adjoint vertical averaged density (ad_rhoA) and adjoint
912! density perturbation (ad_rhoS) used in adjoint barotropic pressure
913! gradient.
914!-----------------------------------------------------------------------
915!
916! Compute temporary intermediate BASIC STATE "rhoS1" and "rhoA1".
917!
918 DO i=istrt,iendt
919 cff1=den(i,n(ng))*hz(i,j,n(ng))
920 rhos1(i,j)=0.5_r8*cff1*hz(i,j,n(ng))
921 rhoa1(i,j)=cff1
922 END DO
923 DO k=n(ng)-1,1,-1
924 DO i=istrt,iendt
925 cff1=den(i,k)*hz(i,j,k)
926 rhos1(i,j)=rhos1(i,j)+hz(i,j,k)*(rhoa1(i,j)+0.5_r8*cff1)
927 rhoa1(i,j)=rhoa1(i,j)+cff1
928 END DO
929 END DO
930!
931 cff2=1.0_r8/rho0
932 DO i=istrt,iendt
933 cff1=1.0_r8/(z_w(i,j,n(ng))-z_w(i,j,0))
934!^ tl_rhoS(i,j)=2.0_r8*cff2* &
935!^ & cff1*(2.0_r8*tl_cff1*rhoS1(i,j)+ &
936!^ & cff1*tl_rhoS(i,j))
937!^
938 adfac=2.0_r8*cff2*cff1*ad_rhos(i,j)
939 ad_cff1=2.0_r8*rhos1(i,j)*adfac
940 ad_rhos(i,j)=cff1*adfac
941!^ tl_rhoA(i,j)=cff2*(tl_cff1*rhoA1(i,j)+cff1*tl_rhoA(i,j))
942!^
943 adfac=cff2*ad_rhoa(i,j)
944 ad_cff1=ad_cff1+rhoa1(i,j)*adfac
945 ad_rhoa(i,j)=cff1*adfac
946!^ tl_cff1=-cff1*cff1*(tl_z_w(i,j,N(ng))-tl_z_w(i,j,0))
947!^
948 adfac=-cff1*cff1*ad_cff1
949 ad_z_w(i,j,n(ng))=ad_z_w(i,j,n(ng))+adfac
950 ad_z_w(i,j,0 )=ad_z_w(i,j,0 )-adfac
951 ad_cff1=0.0_r8
952 END DO
953!
954! Compute appropriate intermediate BASIC STATE "rhoA1".
955!
956 DO i=istrt,iendt
957 cff1=den(i,n(ng))*hz(i,j,n(ng))
958 rhoa1(i,j)=cff1
959 END DO
960 DO k=n(ng)-1,1,-1
961 DO i=istrt,iendt
962 cff1=den(i,k)*hz(i,j,k)
963!^ tl_rhoA(i,j)=tl_rhoA(i,j)+tl_cff1
964!^
965 ad_cff1=ad_rhoa(i,j)
966!^ tl_rhoS(i,j)=tl_rhoS(i,j)+ &
967!^ & tl_Hz(i,j,k)*(rhoA1(i,j)+0.5_r8*cff1)+ &
968!^ & Hz(i,j,k)*(tl_rhoA(i,j)+0.5_r8*tl_cff1)
969!^
970 adfac=hz(i,j,k)*ad_rhos(i,j)
971 ad_rhoa(i,j)=ad_rhoa(i,j)+adfac
972 ad_cff1=ad_cff1+0.5_r8*adfac
973 ad_hz(i,j,k)=ad_hz(i,j,k)+ &
974 & (rhoa1(i,j)+0.5_r8*cff1)*ad_rhos(i,j)
975!^ tl_cff1=tl_den(i,k)*Hz(i,j,k)+ &
976!^ & den(i,k)*tl_Hz(i,j,k)
977!^
978 ad_den(i,k)=ad_den(i,k)+hz(i,j,k)*ad_cff1
979 ad_hz(i,j,k)=ad_hz(i,j,k)+den(i,k)*ad_cff1
980 ad_cff1=0.0_r8
981 rhoa1(i,j)=rhoa1(i,j)+cff1
982 END DO
983 END DO
984 DO i=istrt,iendt
985 cff1=den(i,n(ng))*hz(i,j,n(ng))
986!^ tl_rhoA(i,j)=tl_cff1
987!^
988 ad_cff1=ad_rhoa(i,j)
989 ad_rhoa(i,j)=0.0_r8
990!^ tl_rhoS(i,j)=0.5_r8*(tl_cff1*Hz(i,j,N(ng))+ &
991!^ & cff1*tl_Hz(i,j,N(ng)))
992!^
993 adfac=0.5_r8*ad_rhos(i,j)
994 ad_cff1=ad_cff1+hz(i,j,n(ng))*adfac
995 ad_hz(i,j,n(ng))=ad_hz(i,j,n(ng))+cff1*adfac
996 ad_rhos(i,j)=0.0_r8
997!^ tl_cff1=tl_den(i,N(ng))*Hz(i,j,N(ng))+ &
998!^ & den(i,N(ng))*tl_Hz(i,j,N(ng))
999!^
1000 ad_den(i,n(ng))=ad_den(i,n(ng))+hz(i,j,n(ng))*ad_cff1
1001 ad_hz(i,j,n(ng))=ad_hz(i,j,n(ng))+den(i,n(ng))*ad_cff1
1002 ad_cff1=0.0_r8
1003 END DO
1004# endif
1005!
1006!-----------------------------------------------------------------------
1007! Adjoint nonlinear equation of state.
1008!-----------------------------------------------------------------------
1009!
1010 DO k=1,n(ng)
1011 DO i=istrt,iendt
1012!
1013! Check temperature and salinity minimum valid values. Assign depth
1014! to the pressure.
1015!
1016 tt=max(-2.0_r8,t(i,j,k,nrhs,itemp))
1017# ifdef SALINITY
1018 ts=max(0.0_r8,t(i,j,k,nrhs,isalt))
1019 sqrtts=sqrt(ts)
1020# else
1021 ts=0.0_r8
1022 sqrtts=0.0_r8
1023# endif
1024 tp=z_r(i,j,k)
1025 tpr10=0.1_r8*tp
1026!
1027! Compute local nonlinear equation of state coefficients and their
1028! derivatives when appropriate. These coefficients can be stored
1029! in slab (i,k) arrays to avoid recompute them twice. However, the
1030! equivalent of 50 slabs arrays are required.
1031!
1032 c(0)=q00+tt*(q01+tt*(q02+tt*(q03+tt*(q04+tt*q05))))
1033 c(1)=u00+tt*(u01+tt*(u02+tt*(u03+tt*u04)))
1034 c(2)=v00+tt*(v01+tt*v02)
1035 c(3)=a00+tt*(a01+tt*(a02+tt*(a03+tt*a04)))
1036 c(4)=b00+tt*(b01+tt*(b02+tt*b03))
1037 c(5)=d00+tt*(d01+tt*d02)
1038 c(6)=e00+tt*(e01+tt*(e02+tt*e03))
1039 c(7)=f00+tt*(f01+tt*f02)
1040 c(8)=g01+tt*(g02+tt*g03)
1041 c(9)=h00+tt*(h01+tt*h02)
1042# ifdef EOS_TDERIVATIVE
1043!
1044 dcdt(0)=q01+tt*(2.0_r8*q02+tt*(3.0_r8*q03+tt*(4.0_r8*q04+ &
1045 & tt*5.0_r8*q05)))
1046 dcdt(1)=u01+tt*(2.0_r8*u02+tt*(3.0_r8*u03+tt*4.0_r8*u04))
1047 dcdt(2)=v01+tt*2.0_r8*v02
1048 dcdt(3)=a01+tt*(2.0_r8*a02+tt*(3.0_r8*a03+tt*4.0_r8*a04))
1049 dcdt(4)=b01+tt*(2.0_r8*b02+tt*3.0_r8*b03)
1050 dcdt(5)=d01+tt*2.0_r8*d02
1051 dcdt(6)=e01+tt*(2.0_r8*e02+tt*3.0_r8*e03)
1052 dcdt(7)=f01+tt*2.0_r8*f02
1053 dcdt(8)=g02+tt*2.0_r8*g03
1054 dcdt(9)=h01+tt*2.0_r8*h02
1055!
1056 d2cd2t(0)=2.0_r8*q02+tt*(6.0_r8*q03+tt*(12.0_r8*q04+ &
1057 & tt*20.0_r8*q05))
1058 d2cd2t(1)=2.0_r8*u02+tt*(6.0_r8*u03+tt*12.0_r8*u04)
1059 d2cd2t(2)=2.0_r8*v02
1060 d2cd2t(3)=2.0_r8*a02+tt*(6.0_r8*a03+tt*12.0_r8*a04)
1061 d2cd2t(4)=2.0_r8*b02+tt*6.0_r8*b03
1062 d2cd2t(5)=2.0_r8*d02
1063 d2cd2t(6)=2.0_r8*e02+tt*6.0_r8*e03
1064 d2cd2t(7)=2.0_r8*f02
1065 d2cd2t(8)=2.0_r8*g03
1066 d2cd2t(9)=2.0_r8*h02
1067# endif
1068!
1069!-----------------------------------------------------------------------
1070! Compute local adjoint "in situ" density anomaly (kg/m3 - 1000).
1071!-----------------------------------------------------------------------
1072!
1073 cff=1.0_r8/(bulk(i,k)+tpr10)
1074# ifdef MASKING
1075!^ tl_den(i,k)=tl_den(i,k)*rmask(i,j)
1076!^
1077 ad_den(i,k)=ad_den(i,k)*rmask(i,j)
1078# endif
1079# if defined SEDIMENT_NOT_YET && defined SED_DENS_NOT_YET
1080!^ tl_den(i,k)=tl_den(i,k)+tl_SedDen
1081!^
1082 ad_sedden=ad_sedden+ad_den(i,k)
1083 DO ised=1,nst
1084 itrc=idsed(ised)
1085 cff1=1.0_r8/srho(ised,ng)
1086!^ tl_SedDen=tl_SedDen+ &
1087!^ & (tl_t(i,j,k,nrhs,itrc)* &
1088!^ & (Srho(ised,ng)-den(i,k))- &
1089!^ & t(i,j,k,nrhs,itrc)* &
1090!^ & tl_den(i,k))*cff1
1091!^
1092 adfac=cff1*ad_sedden
1093 ad_den(i,k)=ad_den(i,k)- &
1094 & t(i,j,k,nrhs,idsed(ised))*adfac
1095 ad_t(i,j,k,nrhs,itrc)=ad_t(i,j,k,nrhs,itrc)+ &
1096 & (srho(ised,ng)-den(i,k))*adfac
1097 END DO
1098!^ tl_SedDen=0.0_r8
1099!^
1100 ad_sedden=0.0_r8
1101# endif
1102!^ tl_den(i,k)=tl_den1(i,k)*bulk(i,k)*cff+ &
1103!^ & den1(i,k)*(tl_bulk(i,k)*cff+ &
1104!^ & bulk(i,k)*tl_cff)
1105!^
1106 adfac1=den1(i,k)*ad_den(i,k)
1107 ad_den1(i,k)=ad_den1(i,k)+bulk(i,k)*cff*ad_den(i,k)
1108 ad_bulk(i,k)=ad_bulk(i,k)+cff*adfac1
1109 ad_cff=ad_cff+bulk(i,k)*adfac1
1110 ad_den(i,k)=0.0_r8
1111!^ tl_cff=-cff*cff*(tl_bulk(i,k)+tl_Tpr10)
1112!^
1113 adfac=-cff*cff*ad_cff
1114 ad_bulk(i,k)=ad_bulk(i,k)+adfac
1115 ad_tpr10=ad_tpr10+adfac
1116 ad_cff=0.0_r8
1117
1118# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1119 defined bulk_fluxes_not_yet
1120!
1121! Compute d(bulk)/d(S) and d(bulk)/d(T) derivatives used
1122! in the computation of thermal expansion and saline contraction
1123! coefficients.
1124!
1125!^ tl_DbulkDT(i,k)=tl_dCdT(3)+ &
1126!^ & tl_Ts*(dCdT(4)+sqrtTs*dCdT(5))+ &
1127!^ & Ts*(tl_dCdT(4)+tl_sqrtTs*dCdT(5)+ &
1128!^ & sqrtTs*tl_dCdT(5))- &
1129!^ & tl_Tp*(dCdT(6)+Ts*dCdT(7)- &
1130!^ & Tp*(dCdT(8)+Ts*dCdT(9)))- &
1131!^ & Tp*(tl_dCdT(6)+tl_Ts*dCdT(7)+Ts*tl_dCdT(7)- &
1132!^ & tl_Tp*(dCdT(8)+Ts*dCdT(9))- &
1133!^ & Tp*(tl_dCdT(8)+tl_Ts*dCdT(9)+ &
1134!^ & Ts*tl_dCdT(9)))
1135!^
1136 adfac1=ts*ad_dbulkdt(i,k)
1137 adfac2=tp*ad_dbulkdt(i,k)
1138 adfac3=adfac2*tp
1139 ad_dcdt(3)=ad_dcdt(3)+ad_dbulkdt(i,k)
1140 ad_dcdt(4)=ad_dcdt(4)+adfac1
1141 ad_dcdt(5)=ad_dcdt(5)+sqrtts*adfac1
1142 ad_dcdt(6)=ad_dcdt(6)-adfac2
1143 ad_dcdt(7)=ad_dcdt(7)-ts*adfac2
1144 ad_dcdt(8)=ad_dcdt(8)+adfac3
1145 ad_dcdt(9)=ad_dcdt(9)+ts*adfac3
1146 ad_sqrtts=ad_sqrtts+dcdt(5)*adfac1
1147 ad_ts=ad_ts+ &
1148 & ad_dbulkdt(i,k)* &
1149 & (dcdt(4)+sqrtts*dcdt(5)- &
1150 & tp*(dcdt(7)-tp*dcdt(9)))
1151 ad_tp=ad_tp- &
1152 & ad_dbulkdt(i,k)* &
1153 & (dcdt(6)+ts*dcdt(7)- &
1154 & 2.0_r8*tp*(dcdt(8)+ts*dcdt(9)))
1155 ad_dbulkdt(i,k)=0.0_r8
1156!^ tl_DbulkDS(i,k)=tl_C(4)+ &
1157!^ & 1.5_r8*(tl_sqrtTs*C(5)+ &
1158!^ & sqrtTs*tl_C(5))- &
1159!^ & tl_Tp*(C(7)+sqrtTs*1.5_r8*G00- &
1160!^ & Tp*C(9))- &
1161!^ & Tp*(tl_C(7)+tl_sqrtTs*1.5_r8*G00- &
1162!^ & tl_Tp*C(9)-Tp*tl_C(9))
1163!^
1164 adfac1=1.5_r8*ad_dbulkds(i,k)
1165 adfac2=tp*ad_dbulkds(i,k)
1166 ad_c(4)=ad_c(4)+ad_dbulkds(i,k)
1167 ad_c(5)=ad_c(5)+sqrtts*adfac1
1168 ad_c(7)=ad_c(7)-adfac2
1169 ad_c(9)=ad_c(9)+tp*adfac2
1170 ad_sqrtts=ad_sqrtts+ &
1171 & (c(5)-tp*g00)*adfac1
1172 ad_tp=ad_tp- &
1173 & ad_dbulkds(i,k)* &
1174 & (c(7)+sqrtts*1.5_r8*g00-tp*c(9)- &
1175 & tp*c(9))
1176 ad_dbulkds(i,k)=0.0_r8
1177!^ tl_dCdT(9)=tl_Tt*d2Cd2T(9)
1178!^ tl_dCdT(8)=tl_Tt*d2Cd2T(8)
1179!^ tl_dCdT(7)=tl_Tt*d2Cd2T(7)
1180!^ tl_dCdT(6)=tl_Tt*d2Cd2T(6)
1181!^ tl_dCdT(5)=tl_Tt*d2Cd2T(5)
1182!^ tl_dCdT(4)=tl_Tt*d2Cd2T(4)
1183!^ tl_dCdT(3)=tl_Tt*d2Cd2T(3)
1184!^
1185 ad_tt=ad_tt+d2cd2t(9)*ad_dcdt(9)+ &
1186 & d2cd2t(8)*ad_dcdt(8)+ &
1187 & d2cd2t(7)*ad_dcdt(7)+ &
1188 & d2cd2t(6)*ad_dcdt(6)+ &
1189 & d2cd2t(5)*ad_dcdt(5)+ &
1190 & d2cd2t(4)*ad_dcdt(4)+ &
1191 & d2cd2t(3)*ad_dcdt(3)
1192 ad_dcdt(9)=0.0_r8
1193 ad_dcdt(8)=0.0_r8
1194 ad_dcdt(7)=0.0_r8
1195 ad_dcdt(6)=0.0_r8
1196 ad_dcdt(5)=0.0_r8
1197 ad_dcdt(4)=0.0_r8
1198 ad_dcdt(3)=0.0_r8
1199# endif
1200!
1201! Compute adjoint secant bulk modulus.
1202!
1203!^ tl_bulk (i,k)=tl_bulk0(i,k)- &
1204!^ & tl_Tp*(bulk1(i,k)-Tp*bulk2(i,k))- &
1205!^ & Tp*(tl_bulk1(i,k)- &
1206!^ & tl_Tp*bulk2(i,k)- &
1207!^ & Tp*tl_bulk2(i,k))
1208!^
1209 adfac=tp*ad_bulk(i,k)
1210 ad_bulk0(i,k)=ad_bulk0(i,k)+ad_bulk(i,k)
1211 ad_bulk1(i,k)=ad_bulk1(i,k)-adfac
1212 ad_bulk2(i,k)=ad_bulk2(i,k)+adfac*tp
1213 ad_tp=ad_tp- &
1214 & ad_bulk(i,k)*(bulk1(i,k)-tp*bulk2(i,k))+ &
1215 & adfac*bulk2(i,k)
1216 ad_bulk(i,k)=0.0_r8
1217!^ tl_bulk2(i,k)=tl_C(8)+tl_Ts*C(9)+Ts*tl_C(9)
1218!^
1219 ad_c(8)=ad_c(8)+ad_bulk2(i,k)
1220 ad_c(9)=ad_c(9)+ts*ad_bulk2(i,k)
1221 ad_ts=ad_ts+ad_bulk2(i,k)*c(9)
1222 ad_bulk2(i,k)=0.0_r8
1223!^ tl_bulk1(i,k)=tl_C(6)+ &
1224!^ & tl_Ts*(C(7)+sqrtTs*G00)+ &
1225!^ & Ts*(tl_C(7)+tl_sqrtTs*G00)
1226!^
1227 adfac=ts*ad_bulk1(i,k)
1228 ad_c(6)=ad_c(6)+ad_bulk1(i,k)
1229 ad_c(7)=ad_c(7)+adfac
1230 ad_ts=ad_ts+ad_bulk1(i,k)*(c(7)+sqrtts*g00)
1231 ad_sqrtts=ad_sqrtts+adfac*g00
1232 ad_bulk1(i,k)=0.0_r8
1233!^ tl_bulk0(i,k)=tl_C(3)+ &
1234!^ & tl_Ts*(C(4)+sqrtTs*C(5))+ &
1235!^ & Ts*(tl_C(4)+tl_sqrtTs*C(5)+ &
1236!^ & sqrtTs*tl_C(5))
1237!^
1238 adfac=ts*ad_bulk0(i,k)
1239 ad_c(3)=ad_c(3)+ad_bulk0(i,k)
1240 ad_c(4)=ad_c(4)+adfac
1241 ad_c(5)=ad_c(5)+sqrtts*adfac
1242 ad_ts=ad_ts+ad_bulk0(i,k)*(c(4)+sqrtts*c(5))
1243 ad_sqrtts=ad_sqrtts+c(5)*adfac
1244 ad_bulk0(i,k)=0.0_r8
1245!^ tl_C(9)=tl_Tt*dCdT(9)
1246!^ tl_C(8)=tl_Tt*dCdT(8)
1247!^ tl_C(7)=tl_Tt*dCdT(7)
1248!^ tl_C(6)=tl_Tt*dCdT(6)
1249!^ tl_C(5)=tl_Tt*dCdT(5)
1250!^ tl_C(4)=tl_Tt*dCdT(4)
1251!^ tl_C(3)=tl_Tt*dCdT(3)
1252!^
1253 ad_tt=ad_tt+ad_c(9)*dcdt(9)+ &
1254 & ad_c(8)*dcdt(8)+ &
1255 & ad_c(7)*dcdt(7)+ &
1256 & ad_c(6)*dcdt(6)+ &
1257 & ad_c(5)*dcdt(5)+ &
1258 & ad_c(4)*dcdt(4)+ &
1259 & ad_c(3)*dcdt(3)
1260 ad_c(9)=0.0_r8
1261 ad_c(8)=0.0_r8
1262 ad_c(7)=0.0_r8
1263 ad_c(6)=0.0_r8
1264 ad_c(5)=0.0_r8
1265 ad_c(4)=0.0_r8
1266 ad_c(3)=0.0_r8
1267
1268# ifdef EOS_TDERIVATIVE
1269!
1270! Compute d(den1)/d(S) and d(den1)/d(T) derivatives used in the
1271! computation of thermal expansion and saline contraction
1272! coefficients.
1273!
1274!^ tl_Dden1DT(i,k)=tl_dCdT(0)+ &
1275!^ & tl_Ts*(dCdT(1)+sqrtTs*dCdT(2))+ &
1276!^ & Ts*(tl_dCdT(1)+tl_sqrtTs*dCdT(2)+ &
1277!^ & sqrtTs*tl_dCdT(2))
1278!^
1279 adfac1=ts*ad_dden1dt(i,k)
1280 ad_dcdt(0)=ad_dcdt(0)+ad_dden1dt(i,k)
1281 ad_dcdt(1)=ad_dcdt(1)+adfac1
1282 ad_dcdt(2)=ad_dcdt(2)+sqrtts*adfac1
1283 ad_ts=ad_ts+ &
1284 & (dcdt(1)+sqrtts*dcdt(2))*ad_dden1dt(i,k)
1285 ad_sqrtts=ad_sqrtts+dcdt(2)*adfac1
1286 ad_dden1dt(i,k)=0.0_r8
1287!^ tl_Dden1DS(i,k)=tl_C(1)+ &
1288!^ & 1.5_r8*(tl_C(2)*sqrtTs+ &
1289!^ & C(2)*tl_sqrtTs)+ &
1290!^ & 2.0_r8*W00*tl_Ts
1291!^
1292 adfac1=1.5_r8*ad_dden1ds(i,k)
1293 ad_c(1)=ad_c(1)+ad_dden1ds(i,k)
1294 ad_c(2)=ad_c(2)+sqrtts*adfac1
1295 ad_ts=ad_ts+2.0_r8*w00*ad_dden1ds(i,k)
1296 ad_sqrtts=ad_sqrtts+c(2)*adfac1
1297 ad_dden1ds(i,k)=0.0_r8
1298!^ tl_dCdT(2)=tl_Tt*d2Cd2T(2)
1299!^ tl_dCdT(1)=tl_Tt*d2Cd2T(1)
1300!^ tl_dCdT(0)=tl_Tt*d2Cd2T(0)
1301!^
1302 ad_tt=ad_tt+d2cd2t(2)*ad_dcdt(2)+ &
1303 & d2cd2t(1)*ad_dcdt(1)+ &
1304 & d2cd2t(0)*ad_dcdt(0)
1305 ad_dcdt(2)=0.0_r8
1306 ad_dcdt(1)=0.0_r8
1307 ad_dcdt(0)=0.0_r8
1308# endif
1309!
1310! Compute basic state and tangent linear density (kg/m3) at standard
1311! one atmosphere pressure.
1312!
1313!^ tl_den1(i,k)=tl_C(0)+ &
1314!^ & tl_Ts*(C(1)+sqrtTs*C(2)+Ts*W00)+ &
1315!^ & Ts*(tl_C(1)+tl_sqrtTs*C(2)+ &
1316!^ & sqrtTs*tl_C(2)+tl_Ts*W00)
1317!^
1318 adfac=ts*ad_den1(i,k)
1319 ad_c(0)=ad_c(0)+ad_den1(i,k)
1320 ad_c(1)=ad_c(1)+adfac
1321 ad_c(2)=ad_c(2)+adfac*sqrtts
1322 ad_ts=ad_ts+ &
1323 & ad_den1(i,k)*(c(1)+sqrtts*c(2)+ts*w00)+ &
1324 & adfac*w00
1325 ad_sqrtts=ad_sqrtts+adfac*c(2)
1326 ad_den1(i,k)=0.0_r8
1327!^ tl_C(2)=tl_Tt*dCdT(2)
1328!^ tl_C(1)=tl_Tt*dCdT(1)
1329!^ tl_C(0)=tl_Tt*dCdT(0)
1330!^
1331 ad_tt=ad_tt+ad_c(2)*dcdt(2)+ &
1332 & ad_c(1)*dcdt(1)+ &
1333 & ad_c(0)*dcdt(0)
1334 ad_c(2)=0.0_r8
1335 ad_c(1)=0.0_r8
1336 ad_c(0)=0.0_r8
1337!
1338! Check temperature and salinity minimum valid values. Assign depth
1339! to the pressure.
1340!
1341!^ tl_Tpr10=0.1_r8*tl_Tp
1342!^
1343 ad_tp=ad_tp+0.1_r8*ad_tpr10
1344 ad_tpr10=0.0_r8
1345!^ tl_Tp=tl_z_r(i,j,k)
1346!^
1347 ad_z_r(i,j,k)=ad_z_r(i,j,k)+ad_tp
1348 ad_tp=0.0_r8
1349
1350# ifdef SALINITY
1351 IF (ts.ne.0.0_r8) THEN
1352!^ tl_sqrtTs=0.5_r8*tl_Ts/SQRT(Ts)
1353!^
1354 ad_ts=ad_ts+0.5_r8*ad_sqrtts/sqrt(ts)
1355 ad_sqrtts=0.0_r8
1356 ELSE
1357!^ tl_sqrtTs=0.0_r8
1358!^
1359 ad_sqrtts=0.0_r8
1360 END IF
1361!^ tl_Ts=(0.5_r8-SIGN(0.5_r8,-t(i,j,k,nrhs,isalt)))*
1362!^ & tl_t(i,j,k,nrhs,isalt)
1363!^
1364 ad_t(i,j,k,nrhs,isalt)=ad_t(i,j,k,nrhs,isalt)+ &
1365 & (0.5_r8-sign(0.5_r8, &
1366 & -t(i,j,k,nrhs,isalt)))* &
1367 & ad_ts
1368 ad_ts=0.0_r8
1369# else
1370!^ tl_sqrtTs=0.0_r8
1371!^
1372 ad_sqrtts=0.0_r8
1373!^ tl_Ts=0.0_r8
1374!^
1375 ad_ts=0.0_r8
1376# endif
1377!^ tl_Tt=(0.5_r8-SIGN(0.5_r8,-2.0_r8-t(i,j,k,nrhs,itemp)))*
1378!^ & tl_t(i,j,k,nrhs,itemp)
1379!^
1380 ad_t(i,j,k,nrhs,itemp)=ad_t(i,j,k,nrhs,itemp)+ &
1381 & (0.5_r8-sign(0.5_r8,-2.0_r8- &
1382 & t(i,j,k,nrhs,itemp)))* &
1383 & ad_tt
1384 ad_tt=0.0_r8
1385 END DO
1386 END DO
1387 END DO
1388!
1389 RETURN
1390 END SUBROUTINE ad_rho_eos_tile
1391# endif
1392
1393# ifndef NONLIN_EOS
1394!
1395!***********************************************************************
1396 SUBROUTINE ad_rho_eos_tile (ng, tile, model, &
1397 & LBi, UBi, LBj, UBj, &
1398 & IminS, ImaxS, JminS, JmaxS, &
1399 & nrhs, &
1400# ifdef MASKING
1401 & rmask, &
1402# endif
1403# ifdef VAR_RHO_2D_NOT_YET
1404 & Hz, ad_Hz, &
1405# endif
1406 & z_r, ad_z_r, &
1407# ifdef VAR_RHO_2D_NOT_YET
1408 & z_w, ad_z_w, &
1409# endif
1410 & t, ad_t, &
1411# ifdef VAR_RHO_2D_NOT_YET
1412 & ad_rhoA, ad_rhoS, &
1413# endif
1414# ifdef BV_FREQUENCY_NOT_YET
1415 & ad_bvf, &
1416# endif
1417# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1418 defined BULK_FLUXES_NOT_YET
1419 & ad_alpha, ad_beta, &
1420# ifdef LMD_DDMIX_NOT_YET
1421 & ad_alfaobeta, &
1422# endif
1423# endif
1424# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1425 & ad_pden, &
1426# endif
1427 & rho, ad_rho)
1428!***********************************************************************
1429!
1430 USE mod_param
1431 USE mod_parallel
1432 USE mod_scalars
1433# ifdef SEDIMENT_NOT_YET
1434 USE mod_sediment
1435# endif
1436!
1439# ifdef DISTRIBUTE
1441# endif
1442!
1443! Imported variable declarations.
1444!
1445 integer, intent(in) :: ng, tile, model
1446 integer, intent(in) :: LBi, UBi, LBj, UBj
1447 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
1448 integer, intent(in) :: nrhs
1449!
1450# ifdef ASSUMED_SHAPE
1451# ifdef MASKING
1452 real(r8), intent(in) :: rmask(LBi:,LBj:)
1453# endif
1454# ifdef VAR_RHO_2D_NOT_YET
1455 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
1456# endif
1457 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
1458# ifdef VAR_RHO_2D_NOT_YET
1459 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
1460# endif
1461 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
1462
1463 real(r8), intent(in) :: rho(LBi:,LBj:,:)
1464# ifdef VAR_RHO_2D_NOT_YET
1465 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
1466# endif
1467 real(r8), intent(inout) :: ad_z_r(LBi:,LBj:,:)
1468# ifdef VAR_RHO_2D_NOT_YET
1469 real(r8), intent(inout) :: ad_z_w(LBi:,LBj:,0:)
1470# endif
1471 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
1472# ifdef VAR_RHO_2D_NOT_YET
1473 real(r8), intent(inout) :: ad_rhoA(LBi:,LBj:)
1474 real(r8), intent(inout) :: ad_rhoS(LBi:,LBj:)
1475# endif
1476# ifdef BV_FREQUENCY_NOT_YET
1477 real(r8), intent(inout) :: tl_bvf(LBi:,LBj:,0:)
1478# endif
1479# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1480 defined bulk_fluxes_not_yet
1481 real(r8), intent(inout) :: ad_alpha(LBi:,LBj:)
1482 real(r8), intent(inout) :: ad_beta(LBi:,LBj:)
1483# ifdef LMD_DDMIX_NOT_YET
1484 real(r8), intent(inout) :: ad_alfaobeta(LBi:,LBj:,0:)
1485# endif
1486# endif
1487# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1488 real(r8), intent(inout) :: ad_pden(LBi:,LBj:,:)
1489# endif
1490 real(r8), intent(inout) :: ad_rho(LBi:,LBj:,:)
1491# else
1492# ifdef MASKING
1493 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
1494# endif
1495# ifdef VAR_RHO_2D_NOT_YET
1496 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
1497# endif
1498 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
1499# ifdef VAR_RHO_2D_NOT_YET
1500 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
1501# endif
1502 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
1503
1504 real(r8), intent(inout) :: rho(LBi:UBi,LBj:UBj,N(ng))
1505# ifdef VAR_RHO_2D_NOT_YET
1506 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
1507# endif
1508 real(r8), intent(inout) :: ad_z_r(LBi:UBi,LBj:UBj,N(ng))
1509# ifdef VAR_RHO_2D_NOT_YET
1510 real(r8), intent(inout) :: ad_z_w(LBi:UBi,LBj:UBj,0:N(ng))
1511# endif
1512 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
1513# ifdef VAR_RHO_2D_NOT_YET
1514 real(r8), intent(inout) :: ad_rhoA(LBi:UBi,LBj:UBj)
1515 real(r8), intent(inout) :: ad_rhoS(LBi:UBi,LBj:UBj)
1516# endif
1517# ifdef BV_FREQUENCY_NOT_YET
1518 real(r8), intent(inout) :: ad_bvf(LBi:UBi,LBj:UBj,0:N(ng))
1519# endif
1520# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1521 defined bulk_fluxes_not_yet
1522 real(r8), intent(inout) :: ad_alpha(LBi:UBi,LBj:UBj)
1523 real(r8), intent(inout) :: ad_beta(LBi:UBi,LBj:UBj)
1524# ifdef LMD_DDMIX_NOT_YET
1525 real(r8), intent(inout) :: ad_alfaobeta(LBi:UBi,LBj:UBj,0:N(ng))
1526# endif
1527# endif
1528# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1529 real(r8), intent(inout) :: ad_pden(LBi:UBi,LBj:UBj,N(ng))
1530# endif
1531 real(r8), intent(inout) :: ad_rho(LBi:UBi,LBj:UBj,N(ng))
1532# endif
1533!
1534! Local variable declarations.
1535!
1536 integer :: i, ised, itrc, j, k
1537
1538 real(r8) :: SedDen, cff, cff1, cff2
1539 real(r8) :: ad_SedDen, ad_cff, ad_cff1
1540 real(r8) :: adfac, adfac1
1541
1542# ifdef VAR_RHO_2D_NOT_YET
1543 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: rhoA1
1544 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: rhoS1
1545# endif
1546
1547# include "set_bounds.h"
1548!
1549!=======================================================================
1550! Adjoint linear equation of state.
1551!=======================================================================
1552!
1553!-----------------------------------------------------------------------
1554! Exchange boundary data.
1555!-----------------------------------------------------------------------
1556!
1557# ifdef DISTRIBUTE
1558# ifdef BV_FREQUENCY_NOT_YET
1559!^ CALL mp_exchange3d (ng, tile, model, 1, &
1560!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1561!^ & NghostPoints, &
1562!^ & EWperiodic(ng), NSperiodic(ng), &
1563!^ & tl_bvf)
1564!^
1565 CALL ad_mp_exchange3d (ng, tile, model, 1, &
1566 & lbi, ubi, lbj, ubj, 0, n(ng), &
1567 & nghostpoints, &
1568 & ewperiodic(ng), nsperiodic(ng), &
1569 & ad_bvf)
1570# endif
1571# ifdef VAR_RHO_2D_NOT_YET
1572!^ CALL mp_exchange2d (ng, tile, model, 2, &
1573!^ & LBi, UBi, LBj, UBj, &
1574!^ & NghostPoints, &
1575!^ & EWperiodic(ng), NSperiodic(ng), &
1576!^ & tl_rhoA, tl_rhoS)
1577!^
1578 CALL ad_mp_exchange2d (ng, tile, model, 2, &
1579 & lbi, ubi, lbj, ubj, &
1580 & nghostpoints, &
1581 & ewperiodic(ng), nsperiodic(ng), &
1582 & ad_rhoa, ad_rhos)
1583# endif
1584# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1585 defined bulk_fluxes_not_yet
1586!^ CALL mp_exchange2d (ng, tile, model, 2, &
1587!^ & LBi, UBi, LBj, UBj, &
1588!^ & NghostPoints, &
1589!^ & EWperiodic(ng), NSperiodic(ng), &
1590!^ & tl_alpha, tl_beta)
1591!^
1592 CALL ad_mp_exchange2d (ng, tile, model, 2, &
1593 & lbi, ubi, lbj, ubj, &
1594 & nghostpoints, &
1595 & ewperiodic(ng), nsperiodic(ng), &
1596 & ad_alpha, ad_beta)
1597# ifdef LMD_DDMIX_NOT_YET
1598!^ CALL mp_exchange3d (ng, tile, model, 1, &
1599!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1600!^ & NghostPoints, &
1601!^ & EWperiodic(ng), NSperiodic(ng), &
1602!^ & tl_alfaobeta)
1603!^
1604 CALL ad_mp_exchange3d (ng, tile, model, 1, &
1605 & lbi, ubi, lbj, ubj, 0, n(ng), &
1606 & nghostpoints, &
1607 & ewperiodic(ng), nsperiodic(ng), &
1608 & ad_alfaobeta)
1609# endif
1610# endif
1611# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1612!^ CALL mp_exchange3d (ng, tile, model, 1, &
1613!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
1614!^ & NghostPoints, &
1615!^ & EWperiodic(ng), NSperiodic(ng), &
1616!^ & tl_pden)
1617!^
1618 CALL ad_mp_exchange3d (ng, tile, model, 1, &
1619 & lbi, ubi, lbj, ubj, 1, n(ng), &
1620 & nghostpoints, &
1621 & ewperiodic(ng), nsperiodic(ng), &
1622 & ad_pden)
1623# endif
1624!^ CALL mp_exchange3d (ng, tile, model, 1, &
1625!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
1626!^ & NghostPoints, &
1627!^ & EWperiodic(ng), NSperiodic(ng), &
1628!^ & tl_rho)
1629!^
1630 CALL ad_mp_exchange3d (ng, tile, model, 1, &
1631 & lbi, ubi, lbj, ubj, 1, n(ng), &
1632 & nghostpoints, &
1633 & ewperiodic(ng), nsperiodic(ng), &
1634 & ad_rho)
1635!
1636# endif
1637
1638 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1639# ifdef BV_FREQUENCY_NOT_YET
1640!^ CALL exchange_w3d_tile (ng, tile, &
1641!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1642!^ & tl_bvf)
1643!^
1644 CALL ad_exchange_w3d_tile (ng, tile, &
1645 & lbi, ubi, lbj, ubj, 0, n(ng), &
1646 & ad_bvf)
1647# endif
1648# ifdef VAR_RHO_2D_NOT_YET
1649!^ CALL exchange_r2d_tile (ng, tile, &
1650!^ & LBi, UBi, LBj, UBj, &
1651!^ & tl_rhoS)
1652!^
1653 CALL ad_exchange_r2d_tile (ng, tile, &
1654 & lbi, ubi, lbj, ubj, &
1655 & ad_rhos)
1656!^ CALL exchange_r2d_tile (ng, tile, &
1657!^ & LBi, UBi, LBj, UBj, &
1658!^ & tl_rhoA)
1659!^
1660 CALL ad_exchange_r2d_tile (ng, tile, &
1661 & lbi, ubi, lbj, ubj, &
1662 & ad_rhoa)
1663# endif
1664# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1665 defined bulk_fluxes_not_yet
1666!^ CALL exchange_r2d_tile (ng, tile, &
1667!^ & LBi, UBi, LBj, UBj, &
1668!^ & tl_beta)
1669!^
1670 CALL ad_exchange_r2d_tile (ng, tile, &
1671 & lbi, ubi, lbj, ubj, &
1672 & ad_beta)
1673!^ CALL exchange_r2d_tile (ng, tile, &
1674!^ & LBi, UBi, LBj, UBj, &
1675!^ & tl_alpha)
1676!^
1677 CALL ad_exchange_r2d_tile (ng, tile, &
1678 & lbi, ubi, lbj, ubj, &
1679 & ad_alpha)
1680# ifdef LMD_DDMIX_NOT_YET
1681!^ CALL exchange_w3d_tile (ng, tile, &
1682!^ & LBi, UBi, LBj, UBj, 0, N(ng), &
1683!^ & tl_alfaobeta)
1684!^
1685 CALL ad_exchange_w3d_tile (ng, tile, &
1686 & lbi, ubi, lbj, ubj, 0, n(ng), &
1687 & ad_alfaobeta)
1688# endif
1689# endif
1690# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1691!^ CALL exchange_r3d_tile (ng, tile, &
1692!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
1693!^ & tl_pden)
1694!^
1695 CALL ad_exchange_r3d_tile (ng, tile, &
1696 & lbi, ubi, lbj, ubj, 1, n(ng), &
1697 & ad_pden)
1698# endif
1699!^ CALL exchange_r3d_tile (ng, tile, &
1700!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
1701!^ & tl_rho)
1702!^
1703 CALL ad_exchange_r3d_tile (ng, tile, &
1704 & lbi, ubi, lbj, ubj, 1, n(ng), &
1705 & ad_rho)
1706 END IF
1707!
1708 DO j=jstrt,jendt
1709
1710# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET || \
1711 defined bulk_fluxes_not_yet
1712!
1713!-----------------------------------------------------------------------
1714! Compute adjoint thermal expansion (1/Celsius) and saline contraction
1715! (1/PSU) coefficients.
1716!-----------------------------------------------------------------------
1717
1718# ifdef LMD_DDMIX_NOT_YET
1719!
1720! Compute ratio of thermal expansion and saline contraction
1721! coefficients.
1722!
1723 IF (scoef(ng).eq.0.0_r8) THEN
1724 cff=1.0_r8
1725 ELSE
1726 cff=1.0_r8/scoef(ng)
1727 END IF
1728 DO k=1,n(ng)
1729 DO i=istrt,iendt
1730!^ tl_alfaobeta(i,j,k)=0.0_r8
1731!^
1732 ad_alfaobeta(i,j,k)=0.0_r8
1733 END DO
1734 END DO
1735# endif
1736!
1737! Set thermal expansion and saline contraction coefficients.
1738!
1739 DO i=istrt,iendt
1740# ifdef SALINITY
1741!^ tl_beta(i,j)=0.0_r8
1742!^
1743 ad_beta(i,j)=0.0_r8
1744# else
1745!^ beta(i,j)=0.0_r8
1746!^
1747 ad_beta(i,j)=0.0_r8
1748# endif
1749!^ tl_alpha(i,j)=0.0_r8
1750!^
1751 ad_alpha(i,j)=0.0_r8
1752 END DO
1753# endif
1754
1755# ifdef BV_FREQUENCY_NOT_YET
1756!
1757!-----------------------------------------------------------------------
1758! Compute Brunt-Vaisala frequency (1/s2) at horizontal RHO-points
1759! and vertical W-points.
1760!-----------------------------------------------------------------------
1761!
1762 ad_cff=0.0_r8
1763 DO k=1,n(ng)
1764 DO i=istrt,iendt
1765 cff=1.0_r8/(z_r(i,j,k+1)-z_r(i,j,k))
1766!^ tl_bvf(i,j,k)=-gorho0* &
1767!^ & (cff*(tl_rho(i,j,k+1)-tl_rho(i,j,k))+ &
1768!^ & tl_cff*(rho(i,j,k+1)-rho(i,j,k)))
1769!^
1770 adfac=-gorho0*ad_bvf(i,j,k)
1771 adfac1=adfac*cff
1772 ad_rho(i,j,k+1)=ad_rho(i,j,k+1)+adfac1
1773 ad_rho(i,j,k )=ad_rho(i,j,k )-adfac1
1774 ad_cff=ad_cff+(rho(i,j,k+1)-rho(i,j,k))*adfac
1775 ad_bvf(i,j,k)=0.0_r8
1776!^ tl_cff=-cff*cff*(tl_z_r(i,j,k+1)-tl_z_r(i,j,k))
1777!^
1778 adfac=-cff*cff*ad_cff
1779 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)+adfac
1780 ad_z_r(i,j,k )=ad_z_r(i,j,k )-adfac
1781 ad_cff=0.0_r8
1782 END DO
1783 END DO
1784# endif
1785
1786# ifdef VAR_RHO_2D_NOT_YET
1787!
1788!---------------------------------------------------------------------
1789! Compute adjoint vertical averaged density (ad_rhoA) and adjoint
1790! density perturbation (ad_rhoS) used in adjoint barotropic pressure
1791! gradient.
1792!---------------------------------------------------------------------
1793!
1794! Compute intermediate BASIC STATE rhoS1 and rhoA1.
1795!
1796 DO i=istrt,iendt
1797 cff1=rho(i,j,n(ng))*hz(i,j,n(ng))
1798 rhos1(i,j)=0.5_r8*cff1*hz(i,j,n(ng))
1799 rhoa1(i,j)=cff1
1800 END DO
1801 DO k=n(ng)-1,1,-1
1802 DO i=istrt,iendt
1803 cff1=rho(i,j,k)*hz(i,j,k)
1804 rhos1(i,j)=rhos1(i,j)+hz(i,j,k)*(rhoa1(i,j)+0.5_r8*cff1)
1805 rhoa1(i,j)=rhoa1(i,j)+cff1
1806 END DO
1807 END DO
1808!
1809 cff2=1.0_r8/rho0
1810 DO i=istrt,iendt
1811 cff1=1.0_r8/(z_w(i,j,n(ng))-z_w(i,j,0))
1812!^ tl_rhoS(i,j)=2.0_r8*cff2* &
1813!^ & cff1*(2.0_r8*tl_cff1*rhoS1(i,j)+ &
1814!^ & cff1*tl_rhoS(i,j))
1815!^
1816 adfac=2.0_r8*cff2*cff1*ad_rhos(i,j)
1817 ad_cff1=2.0_r8*rhos1(i,j)*adfac
1818 ad_rhos(i,j)=cff1*adfac
1819!^ tl_rhoA(i,j)=cff2*(tl_cff1*rhoA1(i,j)+cff1*tl_rhoA(i,j))
1820!^
1821 adfac=cff2*ad_rhoa(i,j)
1822 ad_cff1=ad_cff1+rhoa1(i,j)*adfac
1823 ad_rhoa(i,j)=cff1*adfac
1824!^ tl_cff1=-cff1*cff1*(tl_z_w(i,j,N(ng))-tl_z_w(i,j,0))
1825!^
1826 adfac=-cff1*cff1*ad_cff1
1827 ad_z_w(i,j,n(ng))=ad_z_w(i,j,n(ng))+adfac
1828 ad_z_w(i,j,0 )=ad_z_w(i,j,0 )-adfac
1829 ad_cff1=0.0_r8
1830 END DO
1831!
1832! Compute appropriate intermediate BASIC STATE "rhoA1".
1833!
1834 DO i=istrt,iendt
1835 cff1=rho(i,j,n(ng))*hz(i,j,n(ng))
1836 rhoa1(i,j)=cff1
1837 END DO
1838 DO k=n(ng)-1,1,-1
1839 DO i=istrt,iendt
1840 cff1=rho(i,j,k)*hz(i,j,k)
1841!^ tl_rhoA(i,j)=tl_rhoA(i,j)+tl_cff1
1842!^
1843 ad_cff1=ad_rhoa(i,j)
1844!^ tl_rhoS(i,j)=tl_rhoS(i,j)+ &
1845!^ & tl_Hz(i,j,k)*(rhoA1(i,j)+0.5_r8*cff1)+ &
1846!^ & Hz(i,j,k)*(tl_rhoA(i,j)+0.5_r8*tl_cff1)
1847!^
1848 adfac=hz(i,j,k)*ad_rhos(i,j)
1849 ad_rhoa(i,j)=ad_rhoa(i,j)+adfac
1850 ad_cff1=ad_cff1+0.5_r8*adfac
1851 ad_hz(i,j,k)=ad_hz(i,j,k)+ &
1852 & (rhoa1(i,j)+0.5_r8*cff1)*ad_rhos(i,j)
1853!^ tl_cff1=tl_rho(i,j,k)*Hz(i,j,k)+ &
1854!^ & rho(i,j,k)*tl_Hz(i,j,k)
1855!^
1856 ad_rho(i,j,k)=ad_rho(i,j,k)+hz(i,j,k)*ad_cff1
1857 ad_hz(i,j,k)=ad_hz(i,j,k)+rho(i,j,k)*ad_cff1
1858 ad_cff1=0.0_r8
1859 rhoa1(i,j)=rhoa1(i,j)+cff1
1860 END DO
1861 END DO
1862 DO i=istrt,iendt
1863 cff1=rho(i,j,n(ng))*hz(i,j,n(ng))
1864!^ tl_rhoA(i,j)=tl_cff1
1865!^
1866 ad_cff1=ad_rhoa(i,j)
1867 ad_rhoa(i,j)=0.0_r8
1868!^ tl_rhoS(i,j)=0.5_r8*(tl_cff1*Hz(i,j,N(ng))+
1869!^ & cff1*tl_Hz(i,j,N(ng)))
1870!^
1871 adfac=0.5_r8*ad_rhos(i,j)
1872 ad_cff1=ad_cff1+hz(i,j,n(ng))*adfac
1873 ad_hz(i,j,n(ng))=ad_hz(i,j,n(ng))+cff1*adfac
1874 ad_rhos(i,j)=0.0_r8
1875!^ tl_cff1=tl_rho(i,j,N(ng))*Hz(i,j,N(ng))+
1876!^ & rho(i,j,N(ng))*tl_Hz(i,j,N(ng))
1877!^
1878 ad_rho(i,j,n(ng))=ad_rho(i,j,n(ng))+hz(i,j,n(ng))*ad_cff1
1879 ad_hz(i,j,n(ng))=ad_hz(i,j,n(ng))+rho(i,j,n(ng))*ad_cff1
1880 ad_cff1=0.0_r8
1881 END DO
1882# endif
1883!
1884!-----------------------------------------------------------------------
1885! Compute adjoint "in situ" density anomaly (kg/m3 - 1000) using linear
1886! equation of state.
1887!-----------------------------------------------------------------------
1888!
1889 DO k=1,n(ng)
1890 DO i=istrt,iendt
1891# if defined LMD_SKPP_NOT_YET || defined LMD_BKPP_NOT_YET
1892!^ tl_pden(i,j,k)=tl_rho(i,j,k)
1893!^
1894 ad_rho(i,j,k)=ad_rho(i,j,k)+ad_pden(i,j,k)
1895 ad_pden(i,j,k)=0.0
1896# endif
1897# ifdef MASKING
1898!^ tl_rho(i,j,k)=tl_rho(i,j,k)*rmask(i,j)
1899!^
1900 ad_rho(i,j,k)=ad_rho(i,j,k)*rmask(i,j)
1901# endif
1902# if defined SEDIMENT_NOT_YET && defined SED_DENS_NOT_YET
1903!^ tl_rho(i,j,k)=tl_rho(i,j,k)+tl_SedDen
1904!^
1905 ad_sedden=ad_sedden+tl_rho(i,j,k)
1906 DO ised=1,nst
1907 itrc=idsed(ised)
1908 cff1=1.0_r8/srho(ised,ng)
1909!^ tl_SedDen=tl_SedDen+ &
1910!^ & (tl_t(i,j,k,nrhs,itrc)* &
1911!^ & (Srho(ised,ng)-rho(i,j,k))- &
1912!^ & t(i,j,k,nrhs,itrc)* &
1913!^ & tl_rho(i,j,k))*cff1
1914!^
1915 adfac=cff1*ad_sedden
1916 ad_rho(i,j,k)=ad_rho(i,j,k)- &
1917 & t(i,j,k,nrhs,itrc)*adfac
1918 tl_t(i,j,k,nrhs,itrc)=tl_t(i,j,k,nrhs,itrc)+ &
1919 & (srho(ised,ng)-rho(i,j,k))*adfac
1920 END DO
1921!^ tl_SedDen=0.0_r8
1922!^
1923 ad_sedden=0.0_r8
1924# endif
1925# ifdef SALINITY
1926!^ tl_rho(i,j,k)=tl_rho(i,j,k)+ &
1927!^ & R0(ng)*Scoef(ng)*t(i,j,k,nrhs,isalt)
1928!^
1929 ad_t(i,j,k,nrhs,isalt)=ad_t(i,j,k,nrhs,isalt)+ &
1930 & r0(ng)*scoef(ng)*ad_rho(i,j,k)
1931# endif
1932!^ tl_rho(i,j,k)=-R0(ng)*Tcoef*tl_t(i,j,k,nrhs,itemp)
1933!^
1934 ad_t(i,j,k,nrhs,itemp)=ad_t(i,j,k,nrhs,itemp)- &
1935 & r0(ng)*tcoef(ng)*ad_rho(i,j,k)
1936 ad_rho(i,j,k)=0.0_r8
1937 END DO
1938 END DO
1939 END DO
1940!
1941 RETURN
1942 END SUBROUTINE ad_rho_eos_tile
1943# endif
1944#endif
1945 END MODULE ad_rho_eos_mod
subroutine ad_exchange_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, ad_a)
subroutine ad_exchange_w3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, ad_a)
subroutine ad_exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, ad_a)
subroutine, public ad_rho_eos(ng, tile, model)
Definition ad_rho_eos.F:49
subroutine ad_rho_eos_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, rmask, z_r, ad_z_r, t, ad_t, ad_alpha, ad_beta, rho, ad_rho)
Definition ad_rho_eos.F:154
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
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
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 ad_mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, ad_a, ad_b, ad_c, ad_d)
subroutine ad_mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, ad_a, ad_b, ad_c, ad_d)
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