ROMS
Loading...
Searching...
No Matches
mod_ocean.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE mod_ocean
3!
4!git $Id$
5!================================================== Hernan G. Arango ===
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md !
9!=======================================================================
10! !
11! 2D Primitive Variables. !
12! !
13! rubar Right-hand-side of 2D U-momentum equation (m4/s2). !
14! rvbar Right-hand-side of 2D V-momentum equation (m4/s2). !
15! rzeta Right-hand-side of free surface equation (m3/s). !
16! ubar Vertically integrated U-momentum component (m/s). !
17! vbar Vertically integrated V-momentum component (m/s). !
18! zeta Free surface (m). !
19#if defined TIDE_GENERATING_FORCES
20! eq_tide Equilibrium tides (m) used in tide generating forces. !
21#endif
22! !
23! 3D Primitive Variables. !
24! !
25! pden Potential Density anomaly (kg/m3). !
26! rho Density anomaly (kg/m3). !
27! ru Right-hand-side of 3D U-momentum equation (m4/s2). !
28! rv Right hand side of 3D V-momentum equation (m4/s2). !
29! t Tracer type variables (active and passive). !
30! u 3D U-momentum component (m/s). !
31! ua 3D U-momentun component (m/s) at RHO-points, A-grid. !
32! v 3D V-momentum component (m/s). !
33! va 3D U-momentun component (m/s) at RHO-points, A-grid. !
34! W S-coordinate (omega*Hz/mn) vertical velocity (m3/s). !
35#ifdef BIOLOGY
36! !
37! Biology Variables. !
38! !
39# if defined BIO_FENNEL && defined CARBON
40! pH Surface concentration of hydrogen ions. !
41# endif
42# ifdef RED_TIDE
43! CystIni Red tide dinoflagellate bottom cyst initial !
44! concentration (cysts/m2). !
45! DIN_obs Observed Dissolved Inorganic Nutrient (micromoles). !
46# endif
47#endif
48# ifdef HYPOXIA_SRM
49! repiration Total biological respiration rate (1/day). !
50# endif
51! !
52#ifdef WEC
53! Waves Effect on Currents variables, !
54! !
55! zetat zeta total = zetaw + qsp + bh !
56! zetaw Quasi-static sea level adjustment !
57! qsp Quasi-static pressure !
58! bh Bernoulli head !
59! !
60! Stokes velocities: !
61! !
62! ubar_stokes 2D U-Stokes velocity (m/s). !
63! vbar_stokes 2D V-Stokes velocity (m/s). !
64! u_stokes 3D U-Stokes velocity (m/s). !
65! v_stokes 3D V-Stokes velocity (m/s). !
66! W_stokes 3D W-Stokes S-coordindate drift velocity (m3/s). !
67! wstvel 3D W-Stokes Z-coordinate drift velocity (m/s). !
68! !
69#endif
70!=======================================================================
71!
72 USE mod_kinds
73!
74 implicit none
75!
76 PUBLIC :: allocate_ocean
77 PUBLIC :: deallocate_ocean
78 PUBLIC :: initialize_ocean
79!
80!-----------------------------------------------------------------------
81! Define T_OCEAN structure.
82!-----------------------------------------------------------------------
83!
85!
86! Nonlinear model state.
87!
88 real(r8), pointer :: rubar(:,:,:)
89 real(r8), pointer :: rvbar(:,:,:)
90 real(r8), pointer :: rzeta(:,:,:)
91 real(r8), pointer :: ubar(:,:,:)
92 real(r8), pointer :: vbar(:,:,:)
93 real(r8), pointer :: zeta(:,:,:)
94#if defined TIDE_GENERATING_FORCES
95 real(r8), pointer :: eq_tide(:,:)
96#endif
97#if defined WEC_VF
98 real(r8), pointer :: zetat(:,:)
99 real(r8), pointer :: zetaw(:,:)
100 real(r8), pointer :: qsp(:,:)
101 real(r8), pointer :: bh(:,:)
102#endif
103#if defined WEC
104 real(r8), pointer :: ubar_stokes(:,:)
105 real(r8), pointer :: vbar_stokes(:,:)
106#endif
107#ifdef SOLVE3D
108 real(r8), pointer :: pden(:,:,:)
109 real(r8), pointer :: rho(:,:,:)
110 real(r8), pointer :: ru(:,:,:,:)
111 real(r8), pointer :: rv(:,:,:,:)
112 real(r8), pointer :: t(:,:,:,:,:)
113 real(r8), pointer :: u(:,:,:,:)
114 real(r8), pointer :: ua(:,:,:)
115 real(r8), pointer :: v(:,:,:,:)
116 real(r8), pointer :: va(:,:,:)
117 real(r8), pointer :: w(:,:,:)
118# if defined OMEGA_IMPLICIT
119 real(r8), pointer :: wi(:,:,:)
120# endif
121 real(r8), pointer :: wvel(:,:,:)
122# if defined WEC
123 real(r8), pointer :: u_stokes(:,:,:)
124 real(r8), pointer :: v_stokes(:,:,:)
125 real(r8), pointer :: w_stokes(:,:,:)
126 real(r8), pointer :: wstvel(:,:,:)
127# endif
128# ifdef UV_KIRBY
129 real(r8), pointer :: uwave(:,:)
130 real(r8), pointer :: vwave(:,:)
131# endif
132# if defined BIO_FENNEL && defined CARBON
133 real(r8), pointer :: ph(:,:)
134# endif
135# ifdef RED_TIDE
136 real(r8), pointer :: cystini(:,:)
137 real(r8), pointer :: din_obs(:,:,:)
138 real(r8), pointer :: din_obsg(:,:,:,:)
139# endif
140# ifdef HYPOXIA_SRM
141 real(r8), pointer :: respiration(:,:,:)
142# ifndef ANA_RESPIRATION
143 real(r8), pointer :: respirationg(:,:,:,:)
144# endif
145# endif
146#endif
147
148#if defined TANGENT || defined TL_IOMS
149!
150! Tangent linear model state.
151!
152 real(r8), pointer :: tl_rubar(:,:,:)
153 real(r8), pointer :: tl_rvbar(:,:,:)
154 real(r8), pointer :: tl_rzeta(:,:,:)
155 real(r8), pointer :: tl_ubar(:,:,:)
156 real(r8), pointer :: tl_vbar(:,:,:)
157 real(r8), pointer :: tl_zeta(:,:,:)
158# if defined TIDE_GENERATING_FORCES
159 real(r8), pointer :: tl_eq_tide(:,:)
160# endif
161# if defined WEC_VF
162 real(r8), pointer :: tl_zetat(:,:)
163 real(r8), pointer :: tl_zetaw(:,:)
164 real(r8), pointer :: tl_qsp(:,:)
165 real(r8), pointer :: tl_bh(:,:)
166# endif
167# if defined WEC
168 real(r8), pointer :: tl_ubar_stokes(:,:)
169 real(r8), pointer :: tl_vbar_stokes(:,:)
170# endif
171# ifdef SOLVE3D
172 real(r8), pointer :: tl_pden(:,:,:)
173 real(r8), pointer :: tl_rho(:,:,:)
174 real(r8), pointer :: tl_ru(:,:,:,:)
175 real(r8), pointer :: tl_rv(:,:,:,:)
176 real(r8), pointer :: tl_t(:,:,:,:,:)
177 real(r8), pointer :: tl_u(:,:,:,:)
178 real(r8), pointer :: tl_ua(:,:,:)
179 real(r8), pointer :: tl_v(:,:,:,:)
180 real(r8), pointer :: tl_va(:,:,:)
181 real(r8), pointer :: tl_w(:,:,:)
182# if defined OMEGA_IMPLICIT
183 real(r8), pointer :: tl_wi(:,:,:)
184# endif
185# ifdef NOT_YET
186 real(r8), pointer :: tl_wvel(:,:,:)
187# endif
188# if defined WEC
189 real(r8), pointer :: tl_u_stokes(:,:,:)
190 real(r8), pointer :: tl_v_stokes(:,:,:)
191 real(r8), pointer :: tl_w_stokes(:,:,:)
192# endif
193# endif
194#endif
195
196#ifdef ADJOINT
197!
198! Adjoint model state.
199!
200 real(r8), pointer :: ad_rubar(:,:,:)
201 real(r8), pointer :: ad_rvbar(:,:,:)
202 real(r8), pointer :: ad_rzeta(:,:,:)
203 real(r8), pointer :: ad_ubar(:,:,:)
204 real(r8), pointer :: ad_vbar(:,:,:)
205 real(r8), pointer :: ad_zeta(:,:,:)
206 real(r8), pointer :: ad_ubar_sol(:,:)
207 real(r8), pointer :: ad_vbar_sol(:,:)
208 real(r8), pointer :: ad_zeta_sol(:,:)
209# if defined TIDE_GENERATING_FORCES
210 real(r8), pointer :: ad_eq_tide(:,:)
211# endif
212# if defined WEC_VF
213 real(r8), pointer :: ad_zetat(:,:)
214 real(r8), pointer :: ad_zetaw(:,:)
215 real(r8), pointer :: ad_qsp(:,:)
216 real(r8), pointer :: ad_bh(:,:)
217# endif
218# if defined WEC
219 real(r8), pointer :: ad_ubar_stokes(:,:)
220 real(r8), pointer :: ad_vbar_stokes(:,:)
221# endif
222# ifdef SOLVE3D
223 real(r8), pointer :: ad_pden(:,:,:)
224 real(r8), pointer :: ad_rho(:,:,:)
225 real(r8), pointer :: ad_ru(:,:,:,:)
226 real(r8), pointer :: ad_rv(:,:,:,:)
227 real(r8), pointer :: ad_t(:,:,:,:,:)
228 real(r8), pointer :: ad_u(:,:,:,:)
229 real(r8), pointer :: ad_ua(:,:,:)
230 real(r8), pointer :: ad_v(:,:,:,:)
231 real(r8), pointer :: ad_va(:,:,:)
232 real(r8), pointer :: ad_w(:,:,:)
233 real(r8), pointer :: ad_wvel(:,:,:)
234
235 real(r8), pointer :: ad_t_sol(:,:,:,:)
236 real(r8), pointer :: ad_u_sol(:,:,:)
237 real(r8), pointer :: ad_v_sol(:,:,:)
238 real(r8), pointer :: ad_w_sol(:,:,:)
239# if defined OMEGA_IMPLICIT
240 real(r8), pointer :: ad_wi(:,:,:)
241# endif
242# if defined WEC
243 real(r8), pointer :: ad_u_stokes(:,:,:)
244 real(r8), pointer :: ad_v_stokes(:,:,:)
245 real(r8), pointer :: ad_w_stokes(:,:,:)
246# endif
247# endif
248#endif
249
250#if defined FOUR_DVAR || defined IMPULSE || \
251 (defined hessian_sv && defined bnorm)
252!
253! Working arrays to store adjoint impulse forcing, error covariance,
254! standard deviations, or descent conjugate vectors (directions).
255!
256 real(r8), pointer :: b_ubar(:,:,:)
257 real(r8), pointer :: b_vbar(:,:,:)
258 real(r8), pointer :: b_zeta(:,:,:)
259# ifdef SOLVE3D
260 real(r8), pointer :: b_t(:,:,:,:,:)
261 real(r8), pointer :: b_u(:,:,:,:)
262 real(r8), pointer :: b_v(:,:,:,:)
263# endif
264# if defined FOUR_DVAR || (defined HESSIAN_SV && defined BNORM)
265 real(r8), pointer :: d_ubar(:,:)
266 real(r8), pointer :: d_vbar(:,:)
267 real(r8), pointer :: d_zeta(:,:)
268# ifdef SOLVE3D
269 real(r8), pointer :: d_t(:,:,:,:)
270 real(r8), pointer :: d_u(:,:,:)
271 real(r8), pointer :: d_v(:,:,:)
272# endif
273 real(r8), pointer :: e_ubar(:,:,:)
274 real(r8), pointer :: e_vbar(:,:,:)
275 real(r8), pointer :: e_zeta(:,:,:)
276# ifdef SOLVE3D
277 real(r8), pointer :: e_t(:,:,:,:,:)
278 real(r8), pointer :: e_u(:,:,:,:)
279 real(r8), pointer :: e_v(:,:,:,:)
280# endif
281# endif
282# ifdef WEAK_CONSTRAINT
283 real(r8), pointer :: f_ubar(:,:)
284 real(r8), pointer :: f_vbar(:,:)
285 real(r8), pointer :: f_zeta(:,:)
286# ifdef SOLVE3D
287 real(r8), pointer :: f_t(:,:,:,:)
288 real(r8), pointer :: f_u(:,:,:)
289 real(r8), pointer :: f_v(:,:,:)
290# endif
291# ifdef TIME_CONV
292 real(r8), pointer :: f_ubars(:,:,:)
293 real(r8), pointer :: f_vbars(:,:,:)
294 real(r8), pointer :: f_zetas(:,:,:)
295# ifdef SOLVE3D
296 real(r8), pointer :: f_ts(:,:,:,:,:)
297 real(r8), pointer :: f_us(:,:,:,:)
298 real(r8), pointer :: f_vs(:,:,:,:)
299# endif
300# endif
301# endif
302#endif
303#if defined FORCING_SV || defined HESSIAN_FSV
304 real(r8), pointer :: f_ubar(:,:)
305 real(r8), pointer :: f_vbar(:,:)
306 real(r8), pointer :: f_zeta(:,:)
307# ifdef SOLVE3D
308 real(r8), pointer :: f_t(:,:,:,:)
309 real(r8), pointer :: f_u(:,:,:)
310 real(r8), pointer :: f_v(:,:,:)
311# endif
312#endif
313
314#if defined FORWARD_READ && \
315 (defined tangent || defined tl_ioms || defined adjoint)
316!
317! Latest two records of the nonlinear trajectory used to interpolate
318! the background state in the tangent linear and adjoint models.
319!
320# ifdef FORWARD_RHS
321 real(r8), pointer :: rubarg(:,:,:)
322 real(r8), pointer :: rvbarg(:,:,:)
323 real(r8), pointer :: rzetag(:,:,:)
324# endif
325 real(r8), pointer :: ubarg(:,:,:)
326 real(r8), pointer :: vbarg(:,:,:)
327 real(r8), pointer :: zetag(:,:,:)
328# ifdef SOLVE3D
329# ifdef FORWARD_RHS
330 real(r8), pointer :: rug(:,:,:,:)
331 real(r8), pointer :: rvg(:,:,:,:)
332# endif
333 real(r8), pointer :: tg(:,:,:,:,:)
334 real(r8), pointer :: ug(:,:,:,:)
335 real(r8), pointer :: vg(:,:,:,:)
336# endif
337# ifdef WEAK_CONSTRAINT
338 real(r8), pointer :: f_zetag(:,:,:)
339# ifdef SOLVE3D
340 real(r8), pointer :: f_tg(:,:,:,:,:)
341 real(r8), pointer :: f_ug(:,:,:,:)
342 real(r8), pointer :: f_vg(:,:,:,:)
343# endif
344 real(r8), pointer :: f_ubarg(:,:,:)
345 real(r8), pointer :: f_vbarg(:,:,:)
346# endif
347#endif
348
349 END TYPE t_ocean
350!
351 TYPE (t_ocean), allocatable :: ocean(:)
352!
353 CONTAINS
354!
355 SUBROUTINE allocate_ocean (ng, LBi, UBi, LBj, UBj)
356!
357!=======================================================================
358! !
359! This routine allocates all variables in the module for all nested !
360! grids. !
361! !
362!=======================================================================
363!
364 USE mod_param
365!
366! Imported variable declarations.
367!
368 integer, intent(in) :: ng, lbi, ubi, lbj, ubj
369!
370! Local variable declarations.
371!
372 real(r8) :: size2d
373!
374!-----------------------------------------------------------------------
375! Allocate and initialize module variables.
376!-----------------------------------------------------------------------
377!
378 IF (ng.eq.1) allocate ( ocean(ngrids) )
379!
380! Set horizontal array size.
381!
382 size2d=real((ubi-lbi+1)*(ubj-lbj+1),r8)
383!
384! Nonlinear model state.
385!
386 allocate ( ocean(ng) % rubar(lbi:ubi,lbj:ubj,2) )
387 dmem(ng)=dmem(ng)+2.0_r8*size2d
388
389 allocate ( ocean(ng) % rvbar(lbi:ubi,lbj:ubj,2) )
390 dmem(ng)=dmem(ng)+2.0_r8*size2d
391
392 allocate ( ocean(ng) % rzeta(lbi:ubi,lbj:ubj,2) )
393 dmem(ng)=dmem(ng)+2.0_r8*size2d
394
395 allocate ( ocean(ng) % ubar(lbi:ubi,lbj:ubj,3) )
396 dmem(ng)=dmem(ng)+3.0_r8*size2d
397
398 allocate ( ocean(ng) % vbar(lbi:ubi,lbj:ubj,3) )
399 dmem(ng)=dmem(ng)+3.0_r8*size2d
400
401 allocate ( ocean(ng) % zeta(lbi:ubi,lbj:ubj,3) )
402 dmem(ng)=dmem(ng)+3.0_r8*size2d
403
404#if defined TIDE_GENERATING_FORCES
405 allocate ( ocean(ng) % eq_tide(lbi:ubi,lbj:ubj) )
406 dmem(ng)=dmem(ng)+size2d
407#endif
408#if defined WEC_VF
409 allocate ( ocean(ng) % zetat(lbi:ubi,lbj:ubj) )
410 dmem(ng)=dmem(ng)+size2d
411 allocate ( ocean(ng) % zetaw(lbi:ubi,lbj:ubj) )
412 dmem(ng)=dmem(ng)+size2d
413 allocate ( ocean(ng) % qsp(lbi:ubi,lbj:ubj) )
414 dmem(ng)=dmem(ng)+size2d
415 allocate ( ocean(ng) % bh(lbi:ubi,lbj:ubj) )
416 dmem(ng)=dmem(ng)+size2d
417#endif
418#if defined WEC
419 allocate ( ocean(ng) % ubar_stokes(lbi:ubi,lbj:ubj) )
420 dmem(ng)=dmem(ng)+size2d
421
422 allocate ( ocean(ng) % vbar_stokes(lbi:ubi,lbj:ubj) )
423 dmem(ng)=dmem(ng)+size2d
424#endif
425
426#ifdef SOLVE3D
427 allocate ( ocean(ng) % pden(lbi:ubi,lbj:ubj,n(ng)) )
428 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
429
430 allocate ( ocean(ng) % rho(lbi:ubi,lbj:ubj,n(ng)) )
431 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
432
433 allocate ( ocean(ng) % ru(lbi:ubi,lbj:ubj,0:n(ng),2) )
434 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)+1,r8)*size2d
435
436 allocate ( ocean(ng) % rv(lbi:ubi,lbj:ubj,0:n(ng),2) )
437 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)+1,r8)*size2d
438
439 allocate ( ocean(ng) % t(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)) )
440 dmem(ng)=dmem(ng)+3.0_r8*real(n(ng)*nt(ng),r8)*size2d
441
442 allocate ( ocean(ng) % u(lbi:ubi,lbj:ubj,n(ng),2) )
443 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
444
445 allocate ( ocean(ng) % ua(lbi:ubi,lbj:ubj,n(ng)) )
446 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
447
448 allocate ( ocean(ng) % v(lbi:ubi,lbj:ubj,n(ng),2) )
449 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
450
451 allocate ( ocean(ng) % va(lbi:ubi,lbj:ubj,n(ng)) )
452 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
453
454 allocate ( ocean(ng) % W(lbi:ubi,lbj:ubj,0:n(ng)) )
455 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
456
457# if defined OMEGA_IMPLICIT
458 allocate ( ocean(ng) % Wi(lbi:ubi,lbj:ubj,0:n(ng)) )
459 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
460# endif
461
462 allocate ( ocean(ng) % wvel(lbi:ubi,lbj:ubj,0:n(ng)) )
463 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
464
465# if defined WEC
466 allocate ( ocean(ng) % u_stokes(lbi:ubi,lbj:ubj,n(ng)) )
467 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
468
469 allocate ( ocean(ng) % v_stokes(lbi:ubi,lbj:ubj,n(ng)) )
470 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
471
472 allocate ( ocean(ng) % W_stokes(lbi:ubi,lbj:ubj,0:n(ng)) )
473 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
474
475 allocate ( ocean(ng) % wstvel(lbi:ubi,lbj:ubj,0:n(ng)) )
476 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
477# endif
478
479# ifdef UV_KIRBY
480 allocate ( ocean(ng) % uwave(lbi:ubi,lbj:ubj) )
481 dmem(ng)=dmem(ng)+size2d
482
483 allocate ( ocean(ng) % vwave(lbi:ubi,lbj:ubj) )
484 dmem(ng)=dmem(ng)+size2d
485# endif
486
487# if defined BIO_FENNEL && defined CARBON
488 allocate ( ocean(ng) % pH(lbi:ubi,lbj:ubj) )
489 dmem(ng)=dmem(ng)+size2d
490# endif
491
492# ifdef RED_TIDE
493 allocate ( ocean(ng) % CystIni(lbi:ubi,lbj:ubj) )
494 dmem(ng)=dmem(ng)+size2d
495
496 allocate ( ocean(ng) % DIN_obs(lbi:ubi,lbj:ubj,n(ng)) )
497 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
498
499 allocate ( ocean(ng) % DIN_obsG(lbi:ubi,lbj:ubj,n(ng),2) )
500 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
501# endif
502
503# ifdef HYPOXIA_SRM
504 allocate ( ocean(ng) % respiration(lbi:ubi,lbj:ubj,n(ng)) )
505 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
506
507# ifndef ANA_RESPIRATION
508 allocate ( ocean(ng) % respirationG(lbi:ubi,lbj:ubj,n(ng),2) )
509 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
510# endif
511# endif
512#endif
513
514#if defined TANGENT || defined TL_IOMS
515!
516! Tangent linear model state.
517!
518 allocate ( ocean(ng) % tl_rubar(lbi:ubi,lbj:ubj,2) )
519 dmem(ng)=dmem(ng)+2.0_r8*size2d
520
521 allocate ( ocean(ng) % tl_rvbar(lbi:ubi,lbj:ubj,2) )
522 dmem(ng)=dmem(ng)+2.0_r8*size2d
523
524 allocate ( ocean(ng) % tl_rzeta(lbi:ubi,lbj:ubj,2) )
525 dmem(ng)=dmem(ng)+2.0_r8*size2d
526
527 allocate ( ocean(ng) % tl_ubar(lbi:ubi,lbj:ubj,3) )
528 dmem(ng)=dmem(ng)+3.0_r8*size2d
529
530 allocate ( ocean(ng) % tl_vbar(lbi:ubi,lbj:ubj,3) )
531 dmem(ng)=dmem(ng)+3.0_r8*size2d
532
533 allocate ( ocean(ng) % tl_zeta(lbi:ubi,lbj:ubj,3) )
534 dmem(ng)=dmem(ng)+3.0_r8*size2d
535
536# if defined TIDE_GENERATING_FORCES
537 allocate ( ocean(ng) % tl_eq_tide(lbi:ubi,lbj:ubj) )
538 dmem(ng)=dmem(ng)+3.0_r8*size2d
539# endif
540
541# if defined WEC_VF
542 allocate ( ocean(ng) % tl_zetat(lbi:ubi,lbj:ubj,3) )
543 dmem(ng)=dmem(ng)+3.0_r8*size2d
544
545 allocate ( ocean(ng) % tl_zetaw(lbi:ubi,lbj:ubj,3) )
546 dmem(ng)=dmem(ng)+3.0_r8*size2d
547
548 allocate ( ocean(ng) % tl_qsp(lbi:ubi,lbj:ubj,3) )
549 dmem(ng)=dmem(ng)+3.0_r8*size2d
550
551 allocate ( ocean(ng) % tl_bh(lbi:ubi,lbj:ubj,3) )
552 dmem(ng)=dmem(ng)+3.0_r8*size2d
553# endif
554
555# if defined WEC
556 allocate ( ocean(ng) % tl_ubar_stokes(lbi:ubi,lbj:ubj) )
557 dmem(ng)=dmem(ng)+size2d
558
559 allocate ( ocean(ng) % tl_vbar_stokes(lbi:ubi,lbj:ubj) )
560 dmem(ng)=dmem(ng)+size2d
561# endif
562
563# ifdef SOLVE3D
564 allocate ( ocean(ng) % tl_pden(lbi:ubi,lbj:ubj,n(ng)) )
565 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
566
567 allocate ( ocean(ng) % tl_rho(lbi:ubi,lbj:ubj,n(ng)) )
568 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
569
570 allocate ( ocean(ng) % tl_ru(lbi:ubi,lbj:ubj,0:n(ng),2) )
571 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)+1,r8)*size2d
572
573 allocate ( ocean(ng) % tl_rv(lbi:ubi,lbj:ubj,0:n(ng),2) )
574 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)+1,r8)*size2d
575
576 allocate ( ocean(ng) % tl_t(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)) )
577 dmem(ng)=dmem(ng)+3.0_r8*real(n(ng)*nt(ng),r8)*size2d
578
579 allocate ( ocean(ng) % tl_u(lbi:ubi,lbj:ubj,n(ng),2) )
580 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
581
582 allocate ( ocean(ng) % tl_ua(lbi:ubi,lbj:ubj,n(ng)) )
583 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
584
585 allocate ( ocean(ng) % tl_v(lbi:ubi,lbj:ubj,n(ng),2) )
586 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
587
588 allocate ( ocean(ng) % tl_va(lbi:ubi,lbj:ubj,n(ng)) )
589 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
590
591 allocate ( ocean(ng) % tl_W(lbi:ubi,lbj:ubj,0:n(ng)) )
592 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
593
594# if defined OMEGA_IMPLICIT
595 allocate ( ocean(ng) % tl_Wi(lbi:ubi,lbj:ubj,0:n(ng)) )
596 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
597# endif
598
599# ifdef NOT_YET
600 allocate ( ocean(ng) % tl_wvel(lbi:ubi,lbj:ubj,0:n(ng)) )
601 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
602# endif
603
604# if defined WEC
605 allocate ( ocean(ng) % tl_u_stokes(lbi:ubi,lbj:ubj,n(ng)) )
606 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
607
608 allocate ( ocean(ng) % tl_v_stokes(lbi:ubi,lbj:ubj,n(ng)) )
609 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
610
611 allocate ( ocean(ng) % tl_W_stokes(lbi:ubi,lbj:ubj,0:n(ng)) )
612 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
613# endif
614# endif
615#endif
616
617#ifdef ADJOINT
618!
619! Adjoint model state.
620!
621 allocate ( ocean(ng) % ad_rubar(lbi:ubi,lbj:ubj,2) )
622 dmem(ng)=dmem(ng)+2.0_r8*size2d
623
624 allocate ( ocean(ng) % ad_rvbar(lbi:ubi,lbj:ubj,2) )
625 dmem(ng)=dmem(ng)+2.0_r8*size2d
626
627 allocate ( ocean(ng) % ad_rzeta(lbi:ubi,lbj:ubj,2) )
628 dmem(ng)=dmem(ng)+2.0_r8*size2d
629
630 allocate ( ocean(ng) % ad_ubar(lbi:ubi,lbj:ubj,3) )
631 dmem(ng)=dmem(ng)+3.0_r8*size2d
632
633 allocate ( ocean(ng) % ad_vbar(lbi:ubi,lbj:ubj,3) )
634 dmem(ng)=dmem(ng)+3.0_r8*size2d
635
636 allocate ( ocean(ng) % ad_zeta(lbi:ubi,lbj:ubj,3) )
637 dmem(ng)=dmem(ng)+3.0_r8*size2d
638
639 allocate ( ocean(ng) % ad_ubar_sol(lbi:ubi,lbj:ubj) )
640 dmem(ng)=dmem(ng)+size2d
641
642 allocate ( ocean(ng) % ad_vbar_sol(lbi:ubi,lbj:ubj) )
643 dmem(ng)=dmem(ng)+size2d
644
645 allocate ( ocean(ng) % ad_zeta_sol(lbi:ubi,lbj:ubj) )
646 dmem(ng)=dmem(ng)+size2d
647
648# if defined TIDE_GENERATING_FORCES
649 allocate ( ocean(ng) % ad_eq_tide(lbi:ubi,lbj:ubj) )
650 dmem(ng)=dmem(ng)+3.0_r8*size2d
651# endif
652
653# if defined WEC_VF
654 allocate ( ocean(ng) % ad_zetat(lbi:ubi,lbj:ubj,3) )
655 dmem(ng)=dmem(ng)+3.0_r8*size2d
656
657 allocate ( ocean(ng) % ad_zetaw(lbi:ubi,lbj:ubj,3) )
658 dmem(ng)=dmem(ng)+3.0_r8*size2d
659
660 allocate ( ocean(ng) % ad_qsp(lbi:ubi,lbj:ubj,3) )
661 dmem(ng)=dmem(ng)+3.0_r8*size2d
662
663 allocate ( ocean(ng) % ad_bh(lbi:ubi,lbj:ubj,3) )
664 dmem(ng)=dmem(ng)+3.0_r8*size2d
665# endif
666
667# if defined WEC
668 allocate ( ocean(ng) % ad_ubar_stokes(lbi:ubi,lbj:ubj) )
669 dmem(ng)=dmem(ng)+size2d
670
671 allocate ( ocean(ng) % ad_vbar_stokes(lbi:ubi,lbj:ubj) )
672 dmem(ng)=dmem(ng)+size2d
673# endif
674
675# ifdef SOLVE3D
676 allocate ( ocean(ng) % ad_pden(lbi:ubi,lbj:ubj,n(ng)) )
677 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
678
679 allocate ( ocean(ng) % ad_rho(lbi:ubi,lbj:ubj,n(ng)) )
680 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
681
682 allocate ( ocean(ng) % ad_ru(lbi:ubi,lbj:ubj,0:n(ng),2) )
683 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)+1,r8)*size2d
684
685 allocate ( ocean(ng) % ad_rv(lbi:ubi,lbj:ubj,0:n(ng),2) )
686 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)+1,r8)*size2d
687
688 allocate ( ocean(ng) % ad_t(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)) )
689 dmem(ng)=dmem(ng)+3.0_r8*real(n(ng)*nt(ng),r8)*size2d
690
691 allocate ( ocean(ng) % ad_u(lbi:ubi,lbj:ubj,n(ng),2) )
692 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
693
694 allocate ( ocean(ng) % ad_ua(lbi:ubi,lbj:ubj,n(ng)) )
695 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
696
697 allocate ( ocean(ng) % ad_v(lbi:ubi,lbj:ubj,n(ng),2) )
698 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
699
700 allocate ( ocean(ng) % ad_va(lbi:ubi,lbj:ubj,n(ng)) )
701 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
702
703 allocate ( ocean(ng) % ad_W(lbi:ubi,lbj:ubj,0:n(ng)) )
704 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
705
706# if defined OMEGA_IMPLICIT
707 allocate ( ocean(ng) % ad_Wi(lbi:ubi,lbj:ubj,0:n(ng)) )
708 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
709# endif
710
711 allocate ( ocean(ng) % ad_wvel(lbi:ubi,lbj:ubj,0:n(ng)) )
712 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
713
714 allocate ( ocean(ng) % ad_t_sol(lbi:ubi,lbj:ubj,n(ng),nt(ng)) )
715 dmem(ng)=dmem(ng)+real(n(ng)*nt(ng),r8)*size2d
716
717 allocate ( ocean(ng) % ad_u_sol(lbi:ubi,lbj:ubj,n(ng)) )
718 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
719
720 allocate ( ocean(ng) % ad_v_sol(lbi:ubi,lbj:ubj,n(ng)) )
721 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
722
723 allocate ( ocean(ng) % ad_W_sol(lbi:ubi,lbj:ubj,0:n(ng)) )
724 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
725
726# if defined WEC
727 allocate ( ocean(ng) % ad_u_stokes(lbi:ubi,lbj:ubj,n(ng)) )
728 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
729
730 allocate ( ocean(ng) % ad_v_stokes(lbi:ubi,lbj:ubj,n(ng)) )
731 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
732
733 allocate ( ocean(ng) % ad_W_stokes(lbi:ubi,lbj:ubj,0:n(ng)) )
734 dmem(ng)=dmem(ng)+real(n(ng)+1,r8)*size2d
735# endif
736# endif
737#endif
738
739#if defined FOUR_DVAR || defined IMPULSE || \
740 (defined hessian_sv && defined bnorm)
741!
742! Working arrays to store adjoint impulse forcing, background error
743! covariance, background-error standard deviations, or descent
744! conjugate vectors (directions).
745!
746 allocate ( ocean(ng) % b_ubar(lbi:ubi,lbj:ubj,nsa) )
747 dmem(ng)=dmem(ng)+real(nsa,r8)*size2d
748
749 allocate ( ocean(ng) % b_vbar(lbi:ubi,lbj:ubj,nsa) )
750 dmem(ng)=dmem(ng)+real(nsa,r8)*size2d
751
752 allocate ( ocean(ng) % b_zeta(lbi:ubi,lbj:ubj,nsa) )
753 dmem(ng)=dmem(ng)+real(nsa,r8)*size2d
754
755# ifdef SOLVE3D
756 allocate ( ocean(ng) % b_t(lbi:ubi,lbj:ubj,n(ng),nsa,nt(ng)) )
757 dmem(ng)=dmem(ng)+real(n(ng)*nt(ng)*nsa,r8)*size2d
758
759 allocate ( ocean(ng) % b_u(lbi:ubi,lbj:ubj,n(ng),nsa) )
760 dmem(ng)=dmem(ng)+real(n(ng)*nsa,r8)*size2d
761
762 allocate ( ocean(ng) % b_v(lbi:ubi,lbj:ubj,n(ng),nsa) )
763 dmem(ng)=dmem(ng)+real(n(ng)*nsa,r8)*size2d
764# endif
765
766# if defined FOUR_DVAR || (defined HESSIAN_SV && defined BNORM)
767 allocate ( ocean(ng) % d_ubar(lbi:ubi,lbj:ubj) )
768 dmem(ng)=dmem(ng)+size2d
769
770 allocate ( ocean(ng) % d_vbar(lbi:ubi,lbj:ubj) )
771 dmem(ng)=dmem(ng)+size2d
772
773 allocate ( ocean(ng) % d_zeta(lbi:ubi,lbj:ubj) )
774 dmem(ng)=dmem(ng)+size2d
775
776# ifdef SOLVE3D
777 allocate ( ocean(ng) % d_t(lbi:ubi,lbj:ubj,n(ng),nt(ng)) )
778 dmem(ng)=dmem(ng)+real(n(ng)*nt(ng),r8)*size2d
779
780 allocate ( ocean(ng) % d_u(lbi:ubi,lbj:ubj,n(ng)) )
781 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
782
783 allocate ( ocean(ng) % d_v(lbi:ubi,lbj:ubj,n(ng)) )
784 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
785# endif
786
787 allocate ( ocean(ng) % e_ubar(lbi:ubi,lbj:ubj,nsa) )
788 dmem(ng)=dmem(ng)+real(nsa,r8)*size2d
789
790 allocate ( ocean(ng) % e_vbar(lbi:ubi,lbj:ubj,nsa) )
791 dmem(ng)=dmem(ng)+real(nsa,r8)*size2d
792
793 allocate ( ocean(ng) % e_zeta(lbi:ubi,lbj:ubj,nsa) )
794 dmem(ng)=dmem(ng)+real(nsa,r8)*size2d
795
796# ifdef SOLVE3D
797 allocate ( ocean(ng) % e_t(lbi:ubi,lbj:ubj,n(ng),nsa,nt(ng)) )
798 dmem(ng)=dmem(ng)+real(n(ng)*nt(ng)*nsa,r8)*size2d
799
800 allocate ( ocean(ng) % e_u(lbi:ubi,lbj:ubj,n(ng),nsa) )
801 dmem(ng)=dmem(ng)+real(n(ng)*nsa,r8)*size2d
802
803 allocate ( ocean(ng) % e_v(lbi:ubi,lbj:ubj,n(ng),nsa) )
804 dmem(ng)=dmem(ng)+real(n(ng)*nsa,r8)*size2d
805# endif
806
807# ifdef WEAK_CONSTRAINT
808# ifdef SOLVE3D
809 allocate ( ocean(ng) % f_t(lbi:ubi,lbj:ubj,n(ng),nt(ng)) )
810 dmem(ng)=dmem(ng)+real(n(ng)*nt(ng),r8)*size2d
811
812 allocate ( ocean(ng) % f_u(lbi:ubi,lbj:ubj,n(ng)) )
813 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
814
815 allocate ( ocean(ng) % f_v(lbi:ubi,lbj:ubj,n(ng)) )
816 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
817# endif
818 allocate ( ocean(ng) % f_ubar(lbi:ubi,lbj:ubj) )
819 dmem(ng)=dmem(ng)+size2d
820
821 allocate ( ocean(ng) % f_vbar(lbi:ubi,lbj:ubj) )
822 dmem(ng)=dmem(ng)+size2d
823
824 allocate ( ocean(ng) % f_zeta(lbi:ubi,lbj:ubj) )
825 dmem(ng)=dmem(ng)+size2d
826# ifdef TIME_CONV
827# ifdef SOLVE3D
828 allocate ( ocean(ng) % f_tS(lbi:ubi,lbj:ubj,n(ng),nrectc(ng), &
829 & nt(ng)) )
830 dmem(ng)=dmem(ng)+real(n(ng)*nrectc(ng)*nt(ng),r8)*size2d
831
832 allocate ( ocean(ng) % f_uS(lbi:ubi,lbj:ubj,n(ng),nrectc(ng)) )
833 dmem(ng)=dmem(ng)+real(n(ng)*nrectc(ng),r8)*size2d
834
835 allocate ( ocean(ng) % f_vS(lbi:ubi,lbj:ubj,n(ng),nrectc(ng)) )
836 dmem(ng)=dmem(ng)+real(n(ng)*nrectc(ng),r8)*size2d
837# endif
838 allocate ( ocean(ng) % f_ubarS(lbi:ubi,lbj:ubj,nrectc(ng)) )
839 dmem(ng)=dmem(ng)+real(nrectc(ng),r8)*size2d
840
841 allocate ( ocean(ng) % f_vbarS(lbi:ubi,lbj:ubj,nrectc(ng)) )
842 dmem(ng)=dmem(ng)+real(nrectc(ng),r8)*size2d
843
844 allocate ( ocean(ng) % f_zetaS(lbi:ubi,lbj:ubj,nrectc(ng)) )
845 dmem(ng)=dmem(ng)+real(nrectc(ng),r8)*size2d
846# endif
847# endif
848# endif
849#endif
850#if defined FORCING_SV || defined HESSIAN_FSV
851 allocate ( ocean(ng) % f_ubar(lbi:ubi,lbj:ubj) )
852 dmem(ng)=dmem(ng)+size2d
853
854 allocate ( ocean(ng) % f_vbar(lbi:ubi,lbj:ubj) )
855 dmem(ng)=dmem(ng)+size2d
856
857 allocate ( ocean(ng) % f_zeta(lbi:ubi,lbj:ubj) )
858 dmem(ng)=dmem(ng)+size2d
859
860# ifdef SOLVE3D
861 allocate ( ocean(ng) % f_t(lbi:ubi,lbj:ubj,n(ng),nt(ng)) )
862 dmem(ng)=dmem(ng)+real(n(ng)*nt(ng),r8)*size2d
863
864 allocate ( ocean(ng) % f_u(lbi:ubi,lbj:ubj,n(ng)) )
865 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
866
867 allocate ( ocean(ng) % f_v(lbi:ubi,lbj:ubj,n(ng)) )
868 dmem(ng)=dmem(ng)+real(n(ng),r8)*size2d
869# endif
870#endif
871
872#if defined FORWARD_READ && \
873 (defined tangent || defined tl_ioms || defined adjoint)
874!
875! Latest two records of the nonlinear trajectory used to interpolate
876! the background state in the tangent linear and adjoint models.
877!
878# ifdef FORWARD_RHS
879 allocate ( ocean(ng) % rubarG(lbi:ubi,lbj:ubj,2) )
880 dmem(ng)=dmem(ng)+2.0_r8*size2d
881
882 allocate ( ocean(ng) % rvbarG(lbi:ubi,lbj:ubj,2) )
883 dmem(ng)=dmem(ng)+2.0_r8*size2d
884
885 allocate ( ocean(ng) % rzetaG(lbi:ubi,lbj:ubj,2) )
886 dmem(ng)=dmem(ng)+2.0_r8*size2d
887# endif
888 allocate ( ocean(ng) % ubarG(lbi:ubi,lbj:ubj,2) )
889 dmem(ng)=dmem(ng)+2.0_r8*size2d
890
891 allocate ( ocean(ng) % vbarG(lbi:ubi,lbj:ubj,2) )
892 dmem(ng)=dmem(ng)+2.0_r8*size2d
893
894 allocate ( ocean(ng) % zetaG(lbi:ubi,lbj:ubj,2) )
895 dmem(ng)=dmem(ng)+2.0_r8*size2d
896
897# ifdef SOLVE3D
898# ifdef FORWARD_RHS
899 allocate ( ocean(ng) % ruG(lbi:ubi,lbj:ubj,0:n(ng),2) )
900 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)+1,r8)*size2d
901
902 allocate ( ocean(ng) % rvG(lbi:ubi,lbj:ubj,0:n(ng),2) )
903 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)+1,r8)*size2d
904# endif
905 allocate ( ocean(ng) % tG(lbi:ubi,lbj:ubj,n(ng),2,nt(ng)) )
906 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)*nt(ng),r8)*size2d
907
908 allocate ( ocean(ng) % uG(lbi:ubi,lbj:ubj,n(ng),2) )
909 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
910
911 allocate ( ocean(ng) % vG(lbi:ubi,lbj:ubj,n(ng),2) )
912 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
913
914# endif
915# ifdef WEAK_CONSTRAINT
916# ifdef SOLVE3D
917 allocate ( ocean(ng) % f_tG(lbi:ubi,lbj:ubj,n(ng),2,nt(ng)) )
918 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng)*nt(ng),r8)*size2d
919
920 allocate ( ocean(ng) % f_uG(lbi:ubi,lbj:ubj,n(ng),2) )
921 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
922
923 allocate ( ocean(ng) % f_vG(lbi:ubi,lbj:ubj,n(ng),2) )
924 dmem(ng)=dmem(ng)+2.0_r8*real(n(ng),r8)*size2d
925# endif
926 allocate ( ocean(ng) % f_ubarG(lbi:ubi,lbj:ubj,2) )
927 dmem(ng)=dmem(ng)+2.0_r8*size2d
928
929 allocate ( ocean(ng) % f_vbarG(lbi:ubi,lbj:ubj,2) )
930 dmem(ng)=dmem(ng)+2.0_r8*size2d
931
932 allocate ( ocean(ng) % f_zetaG(lbi:ubi,lbj:ubj,2) )
933 dmem(ng)=dmem(ng)+2.0_r8*size2d
934# endif
935#endif
936!
937 RETURN
938 END SUBROUTINE allocate_ocean
939!
940 SUBROUTINE deallocate_ocean (ng)
941!
942!=======================================================================
943! !
944! This routine deallocates all variables in the module for all nested !
945! grids. !
946! !
947!=======================================================================
948!
949 USE mod_param, ONLY : ngrids
950#ifdef SUBOBJECT_DEALLOCATION
951 USE destroy_mod, ONLY : destroy
952#endif
953!
954! Imported variable declarations.
955!
956 integer, intent(in) :: ng
957!
958! Local variable declarations.
959!
960 character (len=*), parameter :: myfile = &
961 & __FILE__//", deallocate_ocean"
962
963#ifdef SUBOBJECT_DEALLOCATION
964!
965!-----------------------------------------------------------------------
966! Deallocate each variable in the derived-type T_OCEAN structure
967! separately.
968!-----------------------------------------------------------------------
969!
970! Nonlinear model state.
971!
972 IF (.not.destroy(ng, ocean(ng)%rubar, myfile, &
973 & __line__, 'OCEAN(ng)%rubar')) RETURN
974
975 IF (.not.destroy(ng, ocean(ng)%rvbar, myfile, &
976 & __line__, 'OCEAN(ng)%rvbar')) RETURN
977
978 IF (.not.destroy(ng, ocean(ng)%rzeta, myfile, &
979 & __line__, 'OCEAN(ng)%rzeta')) RETURN
980
981 IF (.not.destroy(ng, ocean(ng)%ubar, myfile, &
982 & __line__, 'OCEAN(ng)%ubar')) RETURN
983
984 IF (.not.destroy(ng, ocean(ng)%vbar, myfile, &
985 & __line__, 'OCEAN(ng)%vbar')) RETURN
986
987 IF (.not.destroy(ng, ocean(ng)%zeta, myfile, &
988 & __line__, 'OCEAN(ng)%zeta')) RETURN
989
990# if defined TIDE_GENERATING_FORCES
991 IF (.not.destroy(ng, ocean(ng)%eq_tide, myfile, &
992 & __line__, 'OCEAN(ng)%eq_tide')) RETURN
993# endif
994
995# if defined WEC_VF
996 IF (.not.destroy(ng, ocean(ng)%zetat, myfile, &
997 & __line__, 'OCEAN(ng)%zetat')) RETURN
998
999 IF (.not.destroy(ng, ocean(ng)%zetaw, myfile, &
1000 & __line__, 'OCEAN(ng)%zetaw')) RETURN
1001
1002 IF (.not.destroy(ng, ocean(ng)%qsp, myfile, &
1003 & __line__, 'OCEAN(ng)%qsp')) RETURN
1004
1005 IF (.not.destroy(ng, ocean(ng)%bh, myfile, &
1006 & __line__, 'OCEAN(ng)%bh')) RETURN
1007# endif
1008
1009# if defined WEC
1010 IF (.not.destroy(ng, ocean(ng)%ubar_stokes, myfile, &
1011 & __line__, 'OCEAN(ng)%ubar_stokes')) RETURN
1012
1013 IF (.not.destroy(ng, ocean(ng)%vbar_stokes, myfile, &
1014 & __line__, 'OCEAN(ng)%vbar_stokes')) RETURN
1015# endif
1016
1017# ifdef SOLVE3D
1018 IF (.not.destroy(ng, ocean(ng)%pden, myfile, &
1019 & __line__, 'OCEAN(ng)%pden')) RETURN
1020
1021 IF (.not.destroy(ng, ocean(ng)%rho, myfile, &
1022 & __line__, 'OCEAN(ng)%rho')) RETURN
1023
1024 IF (.not.destroy(ng, ocean(ng)%ru, myfile, &
1025 & __line__, 'OCEAN(ng)%ru')) RETURN
1026
1027 IF (.not.destroy(ng, ocean(ng)%rv, myfile, &
1028 & __line__, 'OCEAN(ng)%rv')) RETURN
1029
1030 IF (.not.destroy(ng, ocean(ng)%t, myfile, &
1031 & __line__, 'OCEAN(ng)%t')) RETURN
1032
1033 IF (.not.destroy(ng, ocean(ng)%ua, myfile, &
1034 & __line__, 'OCEAN(ng)%u')) RETURN
1035
1036 IF (.not.destroy(ng, ocean(ng)%u, myfile, &
1037 & __line__, 'OCEAN(ng)%u')) RETURN
1038
1039 IF (.not.destroy(ng, ocean(ng)%v, myfile, &
1040 & __line__, 'OCEAN(ng)%v')) RETURN
1041
1042 IF (.not.destroy(ng, ocean(ng)%W, myfile, &
1043 & __line__, 'OCEAN(ng)%W')) RETURN
1044
1045# if defined OMEGA_IMPLICIT
1046 IF (.not.destroy(ng, ocean(ng)%Wi, myfile, &
1047 & __line__, 'OCEAN(ng)%Wi')) RETURN
1048# endif
1049
1050 IF (.not.destroy(ng, ocean(ng)%wvel, myfile, &
1051 & __line__, 'OCEAN(ng)%wvel')) RETURN
1052
1053# if defined WEC
1054 IF (.not.destroy(ng, ocean(ng)%u_stokes, myfile, &
1055 & __line__, 'OCEAN(ng)%u_stokes')) RETURN
1056
1057 IF (.not.destroy(ng, ocean(ng)%v_stokes, myfile, &
1058 & __line__, 'OCEAN(ng)%v_stokes')) RETURN
1059
1060
1061 IF (.not.destroy(ng, ocean(ng)%W_stokes, myfile, &
1062 & __line__, 'OCEAN(ng)%W_stokes')) RETURN
1063
1064 IF (.not.destroy(ng, ocean(ng)%wstvel, myfile, &
1065 & __line__, 'OCEAN(ng)%wstvel')) RETURN
1066# endif
1067
1068# if defined BIO_FENNEL && defined CARBON
1069 IF (.not.destroy(ng, ocean(ng)%pH, myfile, &
1070 & __line__, 'OCEAN(ng)%pH')) RETURN
1071# endif
1072
1073# ifdef RED_TIDE
1074 IF (.not.destroy(ng, ocean(ng)%CystIni, myfile, &
1075 & __line__, 'OCEAN(ng)%CystIni')) RETURN
1076
1077 IF (.not.destroy(ng, ocean(ng)%DIN_obs, myfile, &
1078 & __line__, 'OCEAN(ng)%DIN_obs')) RETURN
1079
1080 IF (.not.destroy(ng, ocean(ng)%DIN_obsG, myfile, &
1081 & __line__, 'OCEAN(ng)%DIN_obsG')) RETURN
1082# endif
1083
1084# ifdef HYPOXIA_SRM
1085 IF (.not.destroy(ng, ocean(ng)%respiration, myfile, &
1086 & __line__, 'OCEAN(ng)%respiration')) RETURN
1087
1088# ifndef ANA_RESPIRATION
1089 IF (.not.destroy(ng, ocean(ng)%respirationG, myfile, &
1090 & __line__, 'OCEAN(ng)%respirationG')) RETURN
1091# endif
1092# endif
1093# endif
1094
1095# if defined TANGENT || defined TL_IOMS
1096!
1097! Tangent linear model state.
1098!
1099 IF (.not.destroy(ng, ocean(ng)%tl_rubar, myfile, &
1100 & __line__, 'OCEAN(ng)%tl_rubar')) RETURN
1101
1102 IF (.not.destroy(ng, ocean(ng)%tl_rvbar, myfile, &
1103 & __line__, 'OCEAN(ng)%tl_rvbar')) RETURN
1104
1105 IF (.not.destroy(ng, ocean(ng)%tl_rzeta, myfile, &
1106 & __line__, 'OCEAN(ng)%tl_rzeta')) RETURN
1107
1108 IF (.not.destroy(ng, ocean(ng)%tl_ubar, myfile, &
1109 & __line__, 'OCEAN(ng)%tl_ubar')) RETURN
1110
1111 IF (.not.destroy(ng, ocean(ng)%tl_vbar, myfile, &
1112 & __line__, 'OCEAN(ng)%tl_vbar')) RETURN
1113
1114 IF (.not.destroy(ng, ocean(ng)%tl_zeta, myfile, &
1115 & __line__, 'OCEAN(ng)%tl_zeta')) RETURN
1116
1117# if defined TIDE_GENERATING_FORCES
1118 IF (.not.destroy(ng, ocean(ng)%tl_eq_tide, myfile, &
1119 & __line__, 'OCEAN(ng)%tl_eq_tide')) RETURN
1120# endif
1121
1122# if defined WEC_VF
1123 IF (.not.destroy(ng, ocean(ng)%tl_zetat, myfile, &
1124 & __line__, 'OCEAN(ng)%tl_zetat')) RETURN
1125
1126 IF (.not.destroy(ng, ocean(ng)%tl_zetaw, myfile, &
1127 & __line__, 'OCEAN(ng)%tl_zetaw')) RETURN
1128
1129 IF (.not.destroy(ng, ocean(ng)%tl_qsp, myfile, &
1130 & __line__, 'OCEAN(ng)%tl_qsp')) RETURN
1131
1132 IF (.not.destroy(ng, ocean(ng)%tl_bh, myfile, &
1133 & __line__, 'OCEAN(ng)%tl_bh')) RETURN
1134# endif
1135
1136# if defined WEC
1137 IF (.not.destroy(ng, ocean(ng)%tl_ubar_stokes, myfile, &
1138 & __line__, 'OCEAN(ng)%tl_ubar_stokes')) RETURN
1139
1140 IF (.not.destroy(ng, ocean(ng)%tl_vbar_stokes, myfile, &
1141 & __line__, 'OCEAN(ng)%tl_vbar_stokes')) RETURN
1142# endif
1143
1144# ifdef SOLVE3D
1145 IF (.not.destroy(ng, ocean(ng)%tl_pden, myfile, &
1146 & __line__, 'OCEAN(ng)%tl_pden')) RETURN
1147
1148 IF (.not.destroy(ng, ocean(ng)%tl_rho, myfile, &
1149 & __line__, 'OCEAN(ng)%tl_rho')) RETURN
1150
1151 IF (.not.destroy(ng, ocean(ng)%tl_ru, myfile, &
1152 & __line__, 'OCEAN(ng)%tl_ru')) RETURN
1153
1154 IF (.not.destroy(ng, ocean(ng)%tl_rv, myfile, &
1155 & __line__, 'OCEAN(ng)%tl_rv')) RETURN
1156
1157 IF (.not.destroy(ng, ocean(ng)%tl_t, myfile, &
1158 & __line__, 'OCEAN(ng)%tl_t')) RETURN
1159
1160 IF (.not.destroy(ng, ocean(ng)%tl_u, myfile, &
1161 & __line__, 'OCEAN(ng)%tl_u')) RETURN
1162
1163 IF (.not.destroy(ng, ocean(ng)%tl_v, myfile, &
1164 & __line__, 'OCEAN(ng)%tl_v')) RETURN
1165
1166 IF (.not.destroy(ng, ocean(ng)%tl_W, myfile, &
1167 & __line__, 'OCEAN(ng)%tl_W')) RETURN
1168
1169# if defined OMEGA_IMPLICIT
1170 IF (.not.destroy(ng, ocean(ng)%tl_Wi, myfile, &
1171 & __line__, 'OCEAN(ng)%tl_Wi')) RETURN
1172# endif
1173
1174# ifdef NOT_YET
1175 IF (.not.destroy(ng, ocean(ng)%tl_wvel, myfile, &
1176 & __line__, 'OCEAN(ng)%tl_wvel')) RETURN
1177# endif
1178
1179# if defined WEC
1180 IF (.not.destroy(ng, ocean(ng)%tl_u_stokes, myfile, &
1181 & __line__, 'OCEAN(ng)%tl_u_stokes')) RETURN
1182
1183 IF (.not.destroy(ng, ocean(ng)%tl_v_stokes, myfile, &
1184 & __line__, 'OCEAN(ng)%tl_v_stokes')) RETURN
1185
1186 IF (.not.destroy(ng, ocean(ng)%tl_W_stokes, myfile, &
1187 & __line__, 'OCEAN(ng)%tl_W_stokes')) RETURN
1188# endif
1189
1190# endif
1191# endif
1192
1193# ifdef ADJOINT
1194!
1195! Adjoint model state.
1196!
1197 IF (.not.destroy(ng, ocean(ng)%ad_rubar, myfile, &
1198 & __line__, 'OCEAN(ng)%ad_rubar')) RETURN
1199
1200 IF (.not.destroy(ng, ocean(ng)%ad_rvbar, myfile, &
1201 & __line__, 'OCEAN(ng)%ad_rvbar')) RETURN
1202
1203 IF (.not.destroy(ng, ocean(ng)%ad_rzeta, myfile, &
1204 & __line__, 'OCEAN(ng)%ad_rzeta')) RETURN
1205
1206 IF (.not.destroy(ng, ocean(ng)%ad_ubar, myfile, &
1207 & __line__, 'OCEAN(ng)%ad_ubar')) RETURN
1208
1209 IF (.not.destroy(ng, ocean(ng)%ad_vbar, myfile, &
1210 & __line__, 'OCEAN(ng)%ad_vbar')) RETURN
1211
1212 IF (.not.destroy(ng, ocean(ng)%ad_zeta, myfile, &
1213 & __line__, 'OCEAN(ng)%ad_zeta')) RETURN
1214
1215 IF (.not.destroy(ng, ocean(ng)%ad_ubar_sol, myfile, &
1216 & __line__, 'OCEAN(ng)%ad_ubar_sol')) RETURN
1217
1218 IF (.not.destroy(ng, ocean(ng)%ad_vbar_sol, myfile, &
1219 & __line__, 'OCEAN(ng)%ad_vbar_sol')) RETURN
1220
1221 IF (.not.destroy(ng, ocean(ng)%ad_zeta_sol, myfile, &
1222 & __line__, 'OCEAN(ng)%ad_zeta_sol')) RETURN
1223
1224# if defined TIDE_GENERATING_FORCES
1225 IF (.not.destroy(ng, ocean(ng)%ad_eq_tide, myfile, &
1226 & __line__, 'OCEAN(ng)%ad_eq_tide')) RETURN
1227# endif
1228
1229# if defined WEC_VF
1230 IF (.not.destroy(ng, ocean(ng)%ad_zetat, myfile, &
1231 & __line__, 'OCEAN(ng)%ad_zetat')) RETURN
1232
1233 IF (.not.destroy(ng, ocean(ng)%ad_zetaw, myfile, &
1234 & __line__, 'OCEAN(ng)%ad_zetaw')) RETURN
1235
1236 IF (.not.destroy(ng, ocean(ng)%ad_qsp, myfile, &
1237 & __line__, 'OCEAN(ng)%ad_qsp')) RETURN
1238
1239 IF (.not.destroy(ng, ocean(ng)%ad_bh, myfile, &
1240 & __line__, 'OCEAN(ng)%ad_bh')) RETURN
1241# endif
1242
1243# if defined WEC
1244 IF (.not.destroy(ng, ocean(ng)%ad_ubar_stokes, myfile, &
1245 & __line__, 'OCEAN(ng)%ad_ubar_stokes')) RETURN
1246
1247 IF (.not.destroy(ng, ocean(ng)%ad_vbar_stokes, myfile, &
1248 & __line__, 'OCEAN(ng)%ad_vbar_stokes')) RETURN
1249# endif
1250
1251# ifdef SOLVE3D
1252 IF (.not.destroy(ng, ocean(ng)%ad_pden, myfile, &
1253 & __line__, 'OCEAN(ng)%ad_pden')) RETURN
1254
1255 IF (.not.destroy(ng, ocean(ng)%ad_rho, myfile, &
1256 & __line__, 'OCEAN(ng)%ad_rho')) RETURN
1257
1258 IF (.not.destroy(ng, ocean(ng)%ad_ru, myfile, &
1259 & __line__, 'OCEAN(ng)%ad_ru')) RETURN
1260
1261 IF (.not.destroy(ng, ocean(ng)%ad_rv, myfile, &
1262 & __line__, 'OCEAN(ng)%ad_rv')) RETURN
1263
1264 IF (.not.destroy(ng, ocean(ng)%ad_t, myfile, &
1265 & __line__, 'OCEAN(ng)%ad_t')) RETURN
1266
1267 IF (.not.destroy(ng, ocean(ng)%ad_u, myfile, &
1268 & __line__, 'OCEAN(ng)%ad_u')) RETURN
1269
1270 IF (.not.destroy(ng, ocean(ng)%ad_v, myfile, &
1271 & __line__, 'OCEAN(ng)%ad_v')) RETURN
1272
1273 IF (.not.destroy(ng, ocean(ng)%ad_W, myfile, &
1274 & __line__, 'OCEAN(ng)%ad_W')) RETURN
1275
1276# if defined OMEGA_IMPLICIT
1277 IF (.not.destroy(ng, ocean(ng)%ad_Wi, myfile, &
1278 & __line__, 'OCEAN(ng)%ad_Wi')) RETURN
1279# endif
1280
1281 IF (.not.destroy(ng, ocean(ng)%ad_wvel, myfile, &
1282 & __line__, 'OCEAN(ng)%ad_wvel')) RETURN
1283
1284# if defined WEC
1285 IF (.not.destroy(ng, ocean(ng)%ad_u_stokes, myfile, &
1286 & __line__, 'OCEAN(ng)%ad_u_stokes')) RETURN
1287
1288 IF (.not.destroy(ng, ocean(ng)%tl_v_stokes, myfile, &
1289 & __line__, 'OCEAN(ng)%ad_v_stokes')) RETURN
1290
1291 IF (.not.destroy(ng, ocean(ng)%tl_W_stokes, myfile, &
1292 & __line__, 'OCEAN(ng)%ad_W_stokes')) RETURN
1293# endif
1294
1295 IF (.not.destroy(ng, ocean(ng)%ad_t_sol, myfile, &
1296 & __line__, 'OCEAN(ng)%ad_t_sol')) RETURN
1297
1298 IF (.not.destroy(ng, ocean(ng)%ad_u_sol, myfile, &
1299 & __line__, 'OCEAN(ng)%ad_u_sol')) RETURN
1300
1301 IF (.not.destroy(ng, ocean(ng)%ad_v_sol, myfile, &
1302 & __line__, 'OCEAN(ng)%ad_v_sol')) RETURN
1303
1304 IF (.not.destroy(ng, ocean(ng)%ad_W_sol, myfile, &
1305 & __line__, 'OCEAN(ng)%ad_W_sol')) RETURN
1306
1307# endif
1308# endif
1309
1310# if defined FOUR_DVAR || defined IMPULSE || \
1311 (defined hessian_sv && defined bnorm)
1312!
1313! Working arrays to store adjoint impulse forcing, background error
1314! covariance, background-error standard deviations, or descent
1315! conjugate vectors (directions).
1316!
1317 IF (.not.destroy(ng, ocean(ng)%b_ubar, myfile, &
1318 & __line__, 'OCEAN(ng)%b_ubar')) RETURN
1319
1320 IF (.not.destroy(ng, ocean(ng)%b_vbar, myfile, &
1321 & __line__, 'OCEAN(ng)%b_vbar')) RETURN
1322
1323 IF (.not.destroy(ng, ocean(ng)%b_zeta, myfile, &
1324 & __line__, 'OCEAN(ng)%b_zeta')) RETURN
1325
1326# ifdef SOLVE3D
1327 IF (.not.destroy(ng, ocean(ng)%b_t, myfile, &
1328 & __line__, 'OCEAN(ng)%b_t')) RETURN
1329
1330 IF (.not.destroy(ng, ocean(ng)%b_u, myfile, &
1331 & __line__, 'OCEAN(ng)%b_u')) RETURN
1332
1333 IF (.not.destroy(ng, ocean(ng)%b_v, myfile, &
1334 & __line__, 'OCEAN(ng)%b_v')) RETURN
1335# endif
1336
1337# if defined FOUR_DVAR || (defined HESSIAN_SV && defined BNORM)
1338 IF (.not.destroy(ng, ocean(ng)%d_ubar, myfile, &
1339 & __line__, 'OCEAN(ng)%d_ubar')) RETURN
1340
1341 IF (.not.destroy(ng, ocean(ng)%d_vbar, myfile, &
1342 & __line__, 'OCEAN(ng)%d_vbar')) RETURN
1343
1344 IF (.not.destroy(ng, ocean(ng)%d_zeta, myfile, &
1345 & __line__, 'OCEAN(ng)%d_zeta')) RETURN
1346
1347# ifdef SOLVE3D
1348 IF (.not.destroy(ng, ocean(ng)%d_t, myfile, &
1349 & __line__, 'OCEAN(ng)%d_t')) RETURN
1350
1351 IF (.not.destroy(ng, ocean(ng)%d_u, myfile, &
1352 & __line__, 'OCEAN(ng)%d_u')) RETURN
1353
1354 IF (.not.destroy(ng, ocean(ng)%d_v, myfile, &
1355 & __line__, 'OCEAN(ng)%d_v')) RETURN
1356# endif
1357
1358 IF (.not.destroy(ng, ocean(ng)%e_ubar, myfile, &
1359 & __line__, 'OCEAN(ng)%e_ubar')) RETURN
1360
1361 IF (.not.destroy(ng, ocean(ng)%e_vbar, myfile, &
1362 & __line__, 'OCEAN(ng)%e_vbar')) RETURN
1363
1364 IF (.not.destroy(ng, ocean(ng)%e_zeta, myfile, &
1365 & __line__, 'OCEAN(ng)%e_zeta')) RETURN
1366
1367# ifdef SOLVE3D
1368 IF (.not.destroy(ng, ocean(ng)%e_t, myfile, &
1369 & __line__, 'OCEAN(ng)%e_t')) RETURN
1370
1371 IF (.not.destroy(ng, ocean(ng)%e_u, myfile, &
1372 & __line__, 'OCEAN(ng)%e_u')) RETURN
1373
1374 IF (.not.destroy(ng, ocean(ng)%e_v, myfile, &
1375 & __line__, 'OCEAN(ng)%e_v')) RETURN
1376# endif
1377
1378# ifdef WEAK_CONSTRAINT
1379# ifdef SOLVE3D
1380 IF (.not.destroy(ng, ocean(ng)%f_t, myfile, &
1381 & __line__, 'OCEAN(ng)%f_t')) RETURN
1382
1383 IF (.not.destroy(ng, ocean(ng)%f_u, myfile, &
1384 & __line__, 'OCEAN(ng)%f_u')) RETURN
1385
1386 IF (.not.destroy(ng, ocean(ng)%f_v, myfile, &
1387 & __line__, 'OCEAN(ng)%f_v')) RETURN
1388# endif
1389 IF (.not.destroy(ng, ocean(ng)%f_ubar, myfile, &
1390 & __line__, 'OCEAN(ng)%f_ubar')) RETURN
1391
1392 IF (.not.destroy(ng, ocean(ng)%f_vbar, myfile, &
1393 & __line__, 'OCEAN(ng)%f_vbar')) RETURN
1394
1395 IF (.not.destroy(ng, ocean(ng)%f_zeta, myfile, &
1396 & __line__, 'OCEAN(ng)%f_zeta')) RETURN
1397
1398# ifdef TIME_CONV
1399# ifdef SOLVE3D
1400 IF (.not.destroy(ng, ocean(ng)%f_tS, myfile, &
1401 & __line__, 'OCEAN(ng)%f_tS')) RETURN
1402
1403 IF (.not.destroy(ng, ocean(ng)%f_uS, myfile, &
1404 & __line__, 'OCEAN(ng)%f_uS')) RETURN
1405
1406 IF (.not.destroy(ng, ocean(ng)%f_vS, myfile, &
1407 & __line__, 'OCEAN(ng)%f_vS')) RETURN
1408# endif
1409
1410 IF (.not.destroy(ng, ocean(ng)%f_ubarS, myfile, &
1411 & __line__, 'OCEAN(ng)%f_ubarS')) RETURN
1412
1413 IF (.not.destroy(ng, ocean(ng)%f_vbarS, myfile, &
1414 & __line__, 'OCEAN(ng)%f_vbarS')) RETURN
1415
1416 IF (.not.destroy(ng, ocean(ng)%f_zetaS, myfile, &
1417 & __line__, 'OCEAN(ng)%f_zetaS')) RETURN
1418# endif
1419# endif
1420# endif
1421# endif
1422
1423# if defined FORCING_SV || defined HESSIAN_FSV
1424 IF (.not.destroy(ng, ocean(ng)%f_ubar, myfile, &
1425 & __line__, 'OCEAN(ng)%f_ubar')) RETURN
1426
1427 IF (.not.destroy(ng, ocean(ng)%f_vbar, myfile, &
1428 & __line__, 'OCEAN(ng)%f_vbar')) RETURN
1429
1430 IF (.not.destroy(ng, ocean(ng)%f_zeta, myfile, &
1431 & __line__, 'OCEAN(ng)%f_zeta')) RETURN
1432
1433# ifdef SOLVE3D
1434 IF (.not.destroy(ng, ocean(ng)%f_t, myfile, &
1435 & __line__, 'OCEAN(ng)%f_t')) RETURN
1436
1437 IF (.not.destroy(ng, ocean(ng)%f_u, myfile, &
1438 & __line__, 'OCEAN(ng)%f_u')) RETURN
1439
1440 IF (.not.destroy(ng, ocean(ng)%f_v, myfile, &
1441 & __line__, 'OCEAN(ng)%f_v')) RETURN
1442# endif
1443# endif
1444
1445# if defined FORWARD_READ && \
1446 (defined tangent || defined tl_ioms || defined adjoint)
1447!
1448! Latest two records of the nonlinear trajectory used to interpolate
1449! the background state in the tangent linear and adjoint models.
1450!
1451# ifdef FORWARD_RHS
1452 IF (.not.destroy(ng, ocean(ng)%rubarG, myfile, &
1453 & __line__, 'OCEAN(ng)%rubarG')) RETURN
1454
1455 IF (.not.destroy(ng, ocean(ng)%rvbarG, myfile, &
1456 & __line__, 'OCEAN(ng)%rvbarG')) RETURN
1457
1458 IF (.not.destroy(ng, ocean(ng)%rzetaG, myfile, &
1459 & __line__, 'OCEAN(ng)%rzetaG')) RETURN
1460# endif
1461
1462 IF (.not.destroy(ng, ocean(ng)%ubarG, myfile, &
1463 & __line__, 'OCEAN(ng)%ubarG')) RETURN
1464
1465 IF (.not.destroy(ng, ocean(ng)%vbarG, myfile, &
1466 & __line__, 'OCEAN(ng)%vbarG')) RETURN
1467
1468 IF (.not.destroy(ng, ocean(ng)%zetaG, myfile, &
1469 & __line__, 'OCEAN(ng)%zetaG')) RETURN
1470
1471# ifdef SOLVE3D
1472# ifdef FORWARD_RHS
1473 IF (.not.destroy(ng, ocean(ng)%ruG, myfile, &
1474 & __line__, 'OCEAN(ng)%ruG')) RETURN
1475
1476 IF (.not.destroy(ng, ocean(ng)%rvG, myfile, &
1477 & __line__, 'OCEAN(ng)%rvG')) RETURN
1478# endif
1479
1480 IF (.not.destroy(ng, ocean(ng)%tG, myfile, &
1481 & __line__, 'OCEAN(ng)%tG')) RETURN
1482
1483 IF (.not.destroy(ng, ocean(ng)%uG, myfile, &
1484 & __line__, 'OCEAN(ng)%uG')) RETURN
1485
1486 IF (.not.destroy(ng, ocean(ng)%vG, myfile, &
1487 & __line__, 'OCEAN(ng)%vG')) RETURN
1488# endif
1489
1490# ifdef WEAK_CONSTRAINT
1491# ifdef SOLVE3D
1492 IF (.not.destroy(ng, ocean(ng)%f_tG, myfile, &
1493 & __line__, 'OCEAN(ng)%f_tG')) RETURN
1494
1495 IF (.not.destroy(ng, ocean(ng)%f_uG, myfile, &
1496 & __line__, 'OCEAN(ng)%f_uG')) RETURN
1497
1498 IF (.not.destroy(ng, ocean(ng)%f_vG, myfile, &
1499 & __line__, 'OCEAN(ng)%f_vG')) RETURN
1500# endif
1501
1502 IF (.not.destroy(ng, ocean(ng)%f_ubarG, myfile, &
1503 & __line__, 'OCEAN(ng)%f_ubarG')) RETURN
1504
1505 IF (.not.destroy(ng, ocean(ng)%f_vbarG, myfile, &
1506 & __line__, 'OCEAN(ng)%f_vbarG')) RETURN
1507
1508 IF (.not.destroy(ng, ocean(ng)%f_zetaG, myfile, &
1509 & __line__, 'OCEAN(ng)%f_zetaG')) RETURN
1510# endif
1511# endif
1512#endif
1513!
1514!-----------------------------------------------------------------------
1515! Deallocate derived-type OCEAN structure.
1516!-----------------------------------------------------------------------
1517!
1518 IF (ng.eq.ngrids) THEN
1519 IF (allocated(ocean)) deallocate ( ocean )
1520 END IF
1521!
1522 RETURN
1523 END SUBROUTINE deallocate_ocean
1524!
1525 SUBROUTINE initialize_ocean (ng, tile, model)
1526!
1527!=======================================================================
1528! !
1529! This routine initialize all variables in the module using first !
1530! touch distribution policy. In shared-memory configuration, this !
1531! operation actually performs propagation of the "shared arrays" !
1532! across the cluster, unless another policy is specified to !
1533! override the default. !
1534! !
1535!=======================================================================
1536!
1537 USE mod_param
1538!
1539! Imported variable declarations.
1540!
1541 integer, intent(in) :: ng, tile, model
1542!
1543! Local variable declarations.
1544!
1545 integer :: imin, imax, jmin, jmax
1546 integer :: i, j, rec
1547#ifdef SOLVE3D
1548 integer :: itrc, k
1549#endif
1550
1551 real(r8), parameter :: inival = 0.0_r8
1552
1553#include "set_bounds.h"
1554!
1555! Set array initialization range.
1556!
1557#ifdef DISTRIBUTE
1558 imin=bounds(ng)%LBi(tile)
1559 imax=bounds(ng)%UBi(tile)
1560 jmin=bounds(ng)%LBj(tile)
1561 jmax=bounds(ng)%UBj(tile)
1562#else
1563 IF (domain(ng)%Western_Edge(tile)) THEN
1564 imin=bounds(ng)%LBi(tile)
1565 ELSE
1566 imin=istr
1567 END IF
1568 IF (domain(ng)%Eastern_Edge(tile)) THEN
1569 imax=bounds(ng)%UBi(tile)
1570 ELSE
1571 imax=iend
1572 END IF
1573 IF (domain(ng)%Southern_Edge(tile)) THEN
1574 jmin=bounds(ng)%LBj(tile)
1575 ELSE
1576 jmin=jstr
1577 END IF
1578 IF (domain(ng)%Northern_Edge(tile)) THEN
1579 jmax=bounds(ng)%UBj(tile)
1580 ELSE
1581 jmax=jend
1582 END IF
1583#endif
1584!
1585!-----------------------------------------------------------------------
1586! Initialize module variables.
1587!-----------------------------------------------------------------------
1588!
1589! Nonlinear model state.
1590!
1591 IF ((model.eq.0).or.(model.eq.inlm)) THEN
1592 DO j=jmin,jmax
1593 DO i=imin,imax
1594 ocean(ng) % rubar(i,j,1) = inival
1595 ocean(ng) % rubar(i,j,2) = inival
1596 ocean(ng) % rvbar(i,j,1) = inival
1597 ocean(ng) % rvbar(i,j,2) = inival
1598 ocean(ng) % rzeta(i,j,1) = inival
1599 ocean(ng) % rzeta(i,j,2) = inival
1600
1601 ocean(ng) % ubar(i,j,1) = inival
1602 ocean(ng) % ubar(i,j,2) = inival
1603 ocean(ng) % ubar(i,j,3) = inival
1604 ocean(ng) % vbar(i,j,1) = inival
1605 ocean(ng) % vbar(i,j,2) = inival
1606 ocean(ng) % vbar(i,j,3) = inival
1607 ocean(ng) % zeta(i,j,1) = inival
1608 ocean(ng) % zeta(i,j,2) = inival
1609 ocean(ng) % zeta(i,j,3) = inival
1610#if defined TIDE_GENERATING_FORCES
1611 ocean(ng) % eq_tide(i,j) = inival
1612#endif
1613#if defined WEC_VF
1614 ocean(ng) % zetat(i,j) = inival
1615 ocean(ng) % zetaw(i,j) = inival
1616 ocean(ng) % qsp(i,j) = inival
1617 ocean(ng) % bh(i,j) = inival
1618#endif
1619#if defined WEC
1620 ocean(ng) % ubar_stokes(i,j) = inival
1621 ocean(ng) % vbar_stokes(i,j) = inival
1622#endif
1623#ifdef UV_KIRBY
1624 ocean(ng) % uwave(i,j) = inival
1625 ocean(ng) % vwave(i,j) = inival
1626#endif
1627#if defined BIO_FENNEL && defined CARBON && defined SOLVE3D
1628 ocean(ng) % pH(i,j) = 8.0_r8
1629#endif
1630#ifdef RED_TIDE
1631 ocean(ng) % CystIni(i,j) = inival
1632#endif
1633 END DO
1634#ifdef SOLVE3D
1635 DO k=1,n(ng)
1636 DO i=imin,imax
1637 ocean(ng) % pden(i,j,k) = inival
1638 ocean(ng) % rho(i,j,k) = inival
1639
1640 ocean(ng) % u(i,j,k,1) = inival
1641 ocean(ng) % u(i,j,k,2) = inival
1642 ocean(ng) % ua(i,j,k) = inival
1643 ocean(ng) % v(i,j,k,1) = inival
1644 ocean(ng) % v(i,j,k,2) = inival
1645 ocean(ng) % va(i,j,k) = inival
1646# if defined WEC
1647 ocean(ng) % u_stokes(i,j,k) = inival
1648 ocean(ng) % v_stokes(i,j,k) = inival
1649# endif
1650# ifdef RED_TIDE
1651 ocean(ng) % DIN_obs(i,j,k) = inival
1652 ocean(ng) % DIN_obsG(i,j,k,1) = inival
1653 ocean(ng) % DIN_obsG(i,j,k,2) = inival
1654# endif
1655# ifdef HYPOXIA_SRM
1656 ocean(ng) % respiration(i,j,k) = inival
1657# ifndef ANA_RESPIRATION
1658 ocean(ng) % respirationG(i,j,k,1) = inival
1659 ocean(ng) % respirationG(i,j,k,2) = inival
1660# endif
1661# endif
1662 END DO
1663 END DO
1664 DO k=0,n(ng)
1665 DO i=imin,imax
1666 ocean(ng) % ru(i,j,k,1) = inival
1667 ocean(ng) % ru(i,j,k,2) = inival
1668 ocean(ng) % rv(i,j,k,1) = inival
1669 ocean(ng) % rv(i,j,k,2) = inival
1670
1671 ocean(ng) % W(i,j,k) = inival
1672# if defined OMEGA_IMPLICIT
1673 ocean(ng) % Wi(i,j,k) = inival
1674# endif
1675 ocean(ng) % wvel(i,j,k) = inival
1676# if defined WEC
1677 ocean(ng) % W_stokes(i,j,k) = inival
1678 ocean(ng) % wstvel(i,j,k) = inival
1679# endif
1680 END DO
1681 END DO
1682 DO itrc=1,nt(ng)
1683 DO k=1,n(ng)
1684 DO i=imin,imax
1685 ocean(ng) % t(i,j,k,1,itrc) = inival
1686 ocean(ng) % t(i,j,k,2,itrc) = inival
1687 ocean(ng) % t(i,j,k,3,itrc) = inival
1688 END DO
1689 END DO
1690 END DO
1691#endif
1692 END DO
1693 END IF
1694
1695#if defined TANGENT || defined TL_IOMS
1696!
1697! Tangent linear model state.
1698!
1699 IF ((model.eq.0).or.(model.eq.itlm).or.(model.eq.irpm)) THEN
1700 DO j=jmin,jmax
1701 DO i=imin,imax
1702 ocean(ng) % tl_rubar(i,j,1) = inival
1703 ocean(ng) % tl_rubar(i,j,2) = inival
1704 ocean(ng) % tl_rvbar(i,j,1) = inival
1705 ocean(ng) % tl_rvbar(i,j,2) = inival
1706 ocean(ng) % tl_rzeta(i,j,1) = inival
1707 ocean(ng) % tl_rzeta(i,j,2) = inival
1708
1709 ocean(ng) % tl_ubar(i,j,1) = inival
1710 ocean(ng) % tl_ubar(i,j,2) = inival
1711 ocean(ng) % tl_ubar(i,j,3) = inival
1712 ocean(ng) % tl_vbar(i,j,1) = inival
1713 ocean(ng) % tl_vbar(i,j,2) = inival
1714 ocean(ng) % tl_vbar(i,j,3) = inival
1715 ocean(ng) % tl_zeta(i,j,1) = inival
1716 ocean(ng) % tl_zeta(i,j,2) = inival
1717 ocean(ng) % tl_zeta(i,j,3) = inival
1718
1719# if defined TIDE_GENERATING_FORCES
1720 ocean(ng) % tl_eq_tide(i,j) = inival
1721# endif
1722
1723# if defined FORCING_SV || defined HESSIAN_FSV
1724 ocean(ng) % f_ubar(i,j) = inival
1725 ocean(ng) % f_vbar(i,j) = inival
1726 ocean(ng) % f_zeta(i,j) = inival
1727# endif
1728
1729# if defined WEC_VF
1730 ocean(ng) % tl_zetat(i,j) = inival
1731 ocean(ng) % tl_zetaw(i,j) = inival
1732 ocean(ng) % tl_qsp(i,j) = inival
1733 ocean(ng) % tl_bh(i,j) = inival
1734# endif
1735
1736# if defined WEC
1737 ocean(ng) % tl_ubar_stokes(i,j) = inival
1738 ocean(ng) % tl_vbar_stokes(i,j) = inival
1739# endif
1740 END DO
1741# ifdef SOLVE3D
1742 DO k=1,n(ng)
1743 DO i=imin,imax
1744 ocean(ng) % tl_pden(i,j,k) = inival
1745 ocean(ng) % tl_rho(i,j,k) = inival
1746
1747 ocean(ng) % tl_u(i,j,k,1) = inival
1748 ocean(ng) % tl_u(i,j,k,2) = inival
1749 ocean(ng) % tl_ua(i,j,k) = inival
1750 ocean(ng) % tl_v(i,j,k,1) = inival
1751 ocean(ng) % tl_v(i,j,k,2) = inival
1752 ocean(ng) % tl_va(i,j,k) = inival
1753# if defined FORCING_SV || defined HESSIAN_FSV
1754 ocean(ng) % f_u(i,j,k) = inival
1755 ocean(ng) % f_v(i,j,k) = inival
1756# endif
1757# if defined WEC
1758 ocean(ng) % tl_u_stokes(i,j,k) = inival
1759 ocean(ng) % tl_v_stokes(i,j,k) = inival
1760# endif
1761 END DO
1762 END DO
1763 DO k=0,n(ng)
1764 DO i=imin,imax
1765 ocean(ng) % tl_ru(i,j,k,1) = inival
1766 ocean(ng) % tl_ru(i,j,k,2) = inival
1767 ocean(ng) % tl_rv(i,j,k,1) = inival
1768 ocean(ng) % tl_rv(i,j,k,2) = inival
1769
1770 ocean(ng) % tl_W(i,j,k) = inival
1771# if defined OMEGA_IMPLICIT
1772 ocean(ng) % tl_Wi(i,j,k) = inival
1773# endif
1774# ifdef NOT_YET
1775 ocean(ng) % tl_wvel(i,j,k) = inival
1776# endif
1777# if defined WEC
1778 ocean(ng) % tl_W_stokes(i,j,k) = inival
1779# endif
1780 END DO
1781 END DO
1782 DO itrc=1,nt(ng)
1783 DO k=1,n(ng)
1784 DO i=imin,imax
1785 ocean(ng) % tl_t(i,j,k,1,itrc) = inival
1786 ocean(ng) % tl_t(i,j,k,2,itrc) = inival
1787 ocean(ng) % tl_t(i,j,k,3,itrc) = inival
1788# if defined FORCING_SV || defined HESSIAN_FSV
1789 ocean(ng) % f_t(i,j,k,itrc) = inival
1790# endif
1791 END DO
1792 END DO
1793 END DO
1794# endif
1795 END DO
1796 END IF
1797#endif
1798
1799#ifdef ADJOINT
1800!
1801! Adjoint model state.
1802!
1803 IF ((model.eq.0).or.(model.eq.iadm)) THEN
1804 DO j=jmin,jmax
1805 DO i=imin,imax
1806 ocean(ng) % ad_rubar(i,j,1) = inival
1807 ocean(ng) % ad_rubar(i,j,2) = inival
1808 ocean(ng) % ad_rvbar(i,j,1) = inival
1809 ocean(ng) % ad_rvbar(i,j,2) = inival
1810 ocean(ng) % ad_rzeta(i,j,1) = inival
1811 ocean(ng) % ad_rzeta(i,j,2) = inival
1812
1813 ocean(ng) % ad_ubar(i,j,1) = inival
1814 ocean(ng) % ad_ubar(i,j,2) = inival
1815 ocean(ng) % ad_ubar(i,j,3) = inival
1816 ocean(ng) % ad_vbar(i,j,1) = inival
1817 ocean(ng) % ad_vbar(i,j,2) = inival
1818 ocean(ng) % ad_vbar(i,j,3) = inival
1819 ocean(ng) % ad_zeta(i,j,1) = inival
1820 ocean(ng) % ad_zeta(i,j,2) = inival
1821 ocean(ng) % ad_zeta(i,j,3) = inival
1822
1823# if defined TIDE_GENERATING_FORCES
1824 ocean(ng) % ad_eq_tide(i,j) = inival
1825# endif
1826
1827# if defined FORCING_SV | defined HESSIAN_FSV
1828 ocean(ng) % f_ubar(i,j) = inival
1829 ocean(ng) % f_vbar(i,j) = inival
1830 ocean(ng) % f_zeta(i,j) = inival
1831# endif
1832
1833 ocean(ng) % ad_ubar_sol(i,j) = inival
1834 ocean(ng) % ad_vbar_sol(i,j) = inival
1835 ocean(ng) % ad_zeta_sol(i,j) = inival
1836
1837# if defined WEC_VF
1838 ocean(ng) % ad_zetat(i,j) = inival
1839 ocean(ng) % ad_zetaw(i,j) = inival
1840 ocean(ng) % ad_qsp(i,j) = inival
1841 ocean(ng) % ad_bh(i,j) = inival
1842# endif
1843
1844# if defined WEC
1845 ocean(ng) % ad_ubar_stokes(i,j) = inival
1846 ocean(ng) % ad_vbar_stokes(i,j) = inival
1847# endif
1848 END DO
1849# ifdef SOLVE3D
1850 DO k=1,n(ng)
1851 DO i=imin,imax
1852 ocean(ng) % ad_pden(i,j,k) = inival
1853 ocean(ng) % ad_rho(i,j,k) = inival
1854
1855 ocean(ng) % ad_u(i,j,k,1) = inival
1856 ocean(ng) % ad_u(i,j,k,2) = inival
1857 ocean(ng) % ad_ua(i,j,k) = inival
1858 ocean(ng) % ad_v(i,j,k,1) = inival
1859 ocean(ng) % ad_v(i,j,k,2) = inival
1860 ocean(ng) % ad_va(i,j,k) = inival
1861 ocean(ng) % ad_u_sol(i,j,k) = inival
1862 ocean(ng) % ad_v_sol(i,j,k) = inival
1863# if defined FORCING_SV || defined HESSIAN_FSV
1864 ocean(ng) % f_u(i,j,k) = inival
1865 ocean(ng) % f_v(i,j,k) = inival
1866# endif
1867# if defined WEC
1868 ocean(ng) % ad_u_stokes(i,j,k) = inival
1869 ocean(ng) % ad_v_stokes(i,j,k) = inival
1870# endif
1871 END DO
1872 END DO
1873 DO k=0,n(ng)
1874 DO i=imin,imax
1875 ocean(ng) % ad_ru(i,j,k,1) = inival
1876 ocean(ng) % ad_ru(i,j,k,2) = inival
1877 ocean(ng) % ad_rv(i,j,k,1) = inival
1878 ocean(ng) % ad_rv(i,j,k,2) = inival
1879
1880 ocean(ng) % ad_W(i,j,k) = inival
1881 ocean(ng) % ad_W_sol(i,j,k) = inival
1882# if defined OMEGA_IMPLICIT
1883 ocean(ng) % ad_Wi(i,j,k) = inival
1884# endif
1885 ocean(ng) % ad_wvel(i,j,k) = inival
1886# if defined WEC
1887 ocean(ng) % ad_W_stokes(i,j,k) = inival
1888# endif
1889 END DO
1890 END DO
1891 DO itrc=1,nt(ng)
1892 DO k=1,n(ng)
1893 DO i=imin,imax
1894 ocean(ng) % ad_t(i,j,k,1,itrc) = inival
1895 ocean(ng) % ad_t(i,j,k,2,itrc) = inival
1896 ocean(ng) % ad_t(i,j,k,3,itrc) = inival
1897 ocean(ng) % ad_t_sol(i,j,k,itrc) = inival
1898# if defined FORCING_SV || defined HESSIAN_FSV
1899 ocean(ng) % f_t(i,j,k,itrc) = inival
1900# endif
1901 END DO
1902 END DO
1903 END DO
1904# endif
1905 END DO
1906 END IF
1907#endif
1908
1909#if defined FOUR_DVAR || defined IMPULSE || \
1910 (defined hessian_sv && defined bnorm)
1911!
1912! Working arrays to store adjoint impulse forcing, background error
1913! covariance, background-error standard deviations, or descent
1914! conjugate vectors (directions).
1915!
1916 IF (model.eq.0) THEN
1917 DO j=jmin,jmax
1918 DO rec=1,nsa
1919 DO i=imin,imax
1920 ocean(ng) % b_ubar(i,j,rec) = inival
1921 ocean(ng) % b_vbar(i,j,rec) = inival
1922 ocean(ng) % b_zeta(i,j,rec) = inival
1923# if defined FOUR_DVAR || (defined HESSIAN_SV && defined BNORM)
1924 ocean(ng) % e_ubar(i,j,rec) = inival
1925 ocean(ng) % e_vbar(i,j,rec) = inival
1926 ocean(ng) % e_zeta(i,j,rec) = inival
1927# endif
1928 END DO
1929 END DO
1930# ifdef FOUR_DVAR
1931 DO i=imin,imax
1932 ocean(ng) % d_ubar(i,j) = inival
1933 ocean(ng) % d_vbar(i,j) = inival
1934 ocean(ng) % d_zeta(i,j) = inival
1935
1936# ifdef WEAK_CONSTRAINT
1937 ocean(ng) % f_ubar(i,j) = inival
1938 ocean(ng) % f_vbar(i,j) = inival
1939 ocean(ng) % f_zeta(i,j) = inival
1940# endif
1941 END DO
1942# endif
1943# ifdef SOLVE3D
1944 DO rec=1,nsa
1945 DO k=1,n(ng)
1946 DO i=imin,imax
1947 ocean(ng) % b_u(i,j,k,rec) = inival
1948 ocean(ng) % b_v(i,j,k,rec) = inival
1949# ifdef FOUR_DVAR
1950 ocean(ng) % e_u(i,j,k,rec) = inival
1951 ocean(ng) % e_v(i,j,k,rec) = inival
1952# endif
1953 END DO
1954 END DO
1955 END DO
1956# ifdef FOUR_DVAR
1957 DO k=1,n(ng)
1958 DO i=imin,imax
1959 ocean(ng) % d_u(i,j,k) = inival
1960 ocean(ng) % d_v(i,j,k) = inival
1961
1962# ifdef WEAK_CONSTRAINT
1963 ocean(ng) % f_u(i,j,k) = inival
1964 ocean(ng) % f_v(i,j,k) = inival
1965# endif
1966 END DO
1967 END DO
1968# endif
1969 DO itrc=1,nt(ng)
1970 DO rec=1,nsa
1971 DO k=1,n(ng)
1972 DO i=imin,imax
1973 ocean(ng) % b_t(i,j,k,rec,itrc) = inival
1974# ifdef FOUR_DVAR
1975 ocean(ng) % e_t(i,j,k,rec,itrc) = inival
1976# endif
1977 END DO
1978 END DO
1979 END DO
1980# ifdef FOUR_DVAR
1981 DO k=1,n(ng)
1982 DO i=imin,imax
1983 ocean(ng) % d_t(i,j,k,itrc) = inival
1984
1985# ifdef WEAK_CONSTRAINT
1986 ocean(ng) % f_t(i,j,k,itrc) = inival
1987# endif
1988 END DO
1989 END DO
1990# endif
1991 END DO
1992# endif
1993# if defined TIME_CONV && defined WEAK_CONSTRAINT
1994 DO rec=1,nrectc(ng)
1995 DO i=imin,imax
1996 ocean(ng) % f_ubarS(i,j,rec) = inival
1997 ocean(ng) % f_zetaS(i,j,rec) = inival
1998 ocean(ng) % f_vbarS(i,j,rec) = inival
1999 END DO
2000 END DO
2001# ifdef SOLVE3D
2002 DO rec=1,nrectc(ng)
2003 DO k=1,n(ng)
2004 DO i=imin,imax
2005 ocean(ng) % f_uS(i,j,k,rec) = inival
2006 ocean(ng) % f_vS(i,j,k,rec) = inival
2007 END DO
2008 END DO
2009 END DO
2010 DO itrc=1,nt(ng)
2011 DO rec=1,nrectc(ng)
2012 DO k=1,n(ng)
2013 DO i=imin,imax
2014 ocean(ng) % f_tS(i,j,k,rec,itrc) = inival
2015 END DO
2016 END DO
2017 END DO
2018 END DO
2019# endif
2020# endif
2021 END DO
2022 END IF
2023#endif
2024
2025#if defined FORWARD_READ && \
2026 (defined tangent || defined tl_ioms || defined adjoint)
2027!
2028! Latest two records of the nonlinear trajectory used to interpolate
2029! the background state in the tangent linear and adjoint models.
2030!
2031 IF (model.eq.0) THEN
2032 DO j=jmin,jmax
2033 DO i=imin,imax
2034# ifdef FORWARD_RHS
2035 ocean(ng) % rubarG(i,j,1) = inival
2036 ocean(ng) % rubarG(i,j,2) = inival
2037 ocean(ng) % rvbarG(i,j,1) = inival
2038 ocean(ng) % rvbarG(i,j,2) = inival
2039 ocean(ng) % rzetaG(i,j,1) = inival
2040 ocean(ng) % rzetaG(i,j,2) = inival
2041# endif
2042 ocean(ng) % ubarG(i,j,1) = inival
2043 ocean(ng) % ubarG(i,j,2) = inival
2044 ocean(ng) % vbarG(i,j,1) = inival
2045 ocean(ng) % vbarG(i,j,2) = inival
2046 ocean(ng) % zetaG(i,j,1) = inival
2047 ocean(ng) % zetaG(i,j,2) = inival
2048# ifdef WEAK_CONSTRAINT
2049 ocean(ng) % f_zetaG(i,j,1) = inival
2050 ocean(ng) % f_zetaG(i,j,2) = inival
2051 ocean(ng) % f_ubarG(i,j,1) = inival
2052 ocean(ng) % f_ubarG(i,j,2) = inival
2053 ocean(ng) % f_vbarG(i,j,1) = inival
2054 ocean(ng) % f_vbarG(i,j,2) = inival
2055# endif
2056 END DO
2057# ifdef SOLVE3D
2058 DO k=1,n(ng)
2059 DO i=imin,imax
2060 ocean(ng) % uG(i,j,k,1) = inival
2061 ocean(ng) % uG(i,j,k,2) = inival
2062 ocean(ng) % vG(i,j,k,1) = inival
2063 ocean(ng) % vG(i,j,k,2) = inival
2064# ifdef WEAK_CONSTRAINT
2065 ocean(ng) % f_uG(i,j,k,1) = inival
2066 ocean(ng) % f_uG(i,j,k,2) = inival
2067 ocean(ng) % f_vG(i,j,k,1) = inival
2068 ocean(ng) % f_vG(i,j,k,2) = inival
2069# endif
2070 END DO
2071 END DO
2072# ifdef FORWARD_RHS
2073 DO k=0,n(ng)
2074 DO i=imin,imax
2075 ocean(ng) % ruG(i,j,k,1) = inival
2076 ocean(ng) % ruG(i,j,k,2) = inival
2077 ocean(ng) % rvG(i,j,k,1) = inival
2078 ocean(ng) % rvG(i,j,k,2) = inival
2079 END DO
2080 END DO
2081# endif
2082 DO itrc=1,nt(ng)
2083 DO k=1,n(ng)
2084 DO i=imin,imax
2085 ocean(ng) % tG(i,j,k,1,itrc) = inival
2086 ocean(ng) % tG(i,j,k,2,itrc) = inival
2087# ifdef WEAK_CONSTRAINT
2088 ocean(ng) % f_tG(i,j,k,1,itrc) = inival
2089 ocean(ng) % f_tG(i,j,k,2,itrc) = inival
2090# endif
2091 END DO
2092 END DO
2093 END DO
2094# endif
2095 END DO
2096 END IF
2097#endif
2098 RETURN
2099 END SUBROUTINE initialize_ocean
2100
2101 END MODULE mod_ocean
integer, parameter r8
Definition mod_kinds.F:28
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
subroutine, public deallocate_ocean(ng)
Definition mod_ocean.F:941
subroutine, public initialize_ocean(ng, tile, model)
Definition mod_ocean.F:1526
subroutine, public allocate_ocean(ng, lbi, ubi, lbj, ubj)
Definition mod_ocean.F:356
integer, dimension(:), allocatable nrectc
Definition mod_param.F:627
integer, parameter inlm
Definition mod_param.F:662
integer, parameter irpm
Definition mod_param.F:664
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
real(r8), dimension(:), allocatable dmem
Definition mod_param.F:137
integer, parameter iadm
Definition mod_param.F:665
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer ngrids
Definition mod_param.F:113
integer, parameter itlm
Definition mod_param.F:663
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer nsa
Definition mod_param.F:612