ROMS
Loading...
Searching...
No Matches
mod_forces.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE mod_forces
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! Surface momentum stresses. !
12! !
13! sustr Surface momentum flux (wind stress) in the !
14! XI-direction (m2/s2) at horizontal U-points. !
15! sustrG Latest two-time snapshots of input "sustr" grided !
16! data used for interpolation. !
17! svstr Surface momentum flux (wind stress) in the !
18! ETA-direction (m2/s2) at horizontal V-points. !
19! svstrG Latest two-time snapshots of input "svstr" grided !
20! data used for interpolation. !
21! !
22! Bottom momentum stresses. !
23! !
24! bustr Bottom momentum flux (bottom stress) in the !
25! XI-direction (m2/s2) at horizontal U-points. !
26! bvstr Bottom momentum flux (bottom stress) in the !
27! ETA-direction (m2/s2) at horizontal V-points. !
28! !
29! Surface wind induced waves. !
30! !
31! Hwave Surface wind induced wave height (m). !
32! HwaveG Latest two-time snapshots of input "Hwave" grided !
33! data used for interpolation. !
34! Dwave Surface wind induced mean wave direction (radians). !
35! DwaveG Latest two-time snapshots of input "Dwave" grided !
36! data used for interpolation. !
37! Dwavep Surface wind induced peak wave direction (radians). !
38! DwavepG Latest two-time snapshots of input "Dwavep" grided !
39! data used for interpolation. !
40! Lwave Mean surface wavelength read in from swan output !
41! LwaveG Latest two-time snapshots of input "Lwave" grided !
42! data used for interpolation. !
43! Lwavep Peak surface wavelength read in from swan output !
44! LwavepG Latest two-time snapshots of input "Lwavep" grided !
45! data used for interpolation. !
46! Pwave_top Wind induced surface wave period (s). !
47! Pwave_topG Latest two-time snapshots of input "Pwave_top" grided !
48! data used for interpolation. !
49! Pwave_bot Wind induced bottom wave period (s). !
50! Pwave_botG Latest two-time snapshots of input "Pwave_bot" grided !
51! data used for interpolation. !
52! Uwave_rms Bottom orbital velocity read in from swan output !
53! Uwave_rmsG Latest two-time snapshots of input "Uwave_rms" grided !
54! data used for interpolation. !
55! Wave_break Percent of wave breaking for use with roller model. !
56! Wave_breakG Latest two-time snapshots of input "wave_break" !
57! gridded data used for interpolation. !
58! Wave_ds Wave directional spreading. !
59! Wave_qp Wave spectrum peakedness. !
60! !
61! Solar shortwave radiation flux. !
62! !
63! srflx Surface shortwave solar radiation flux (degC m/s) !
64! at horizontal RHO-points !
65! srflxG Latest two-time snapshots of input "srflx" grided !
66! data used for interpolation. !
67! !
68! Cloud fraction. !
69! !
70! cloud Cloud fraction (percentage/100). !
71! cloudG Latest two-time snapshots of input "cloud" grided !
72! data used for interpolation. !
73! !
74! Surface heat fluxes, Atmosphere-Ocean bulk parameterization. !
75! !
76! lhflx Latent heat flux (degC m/s). !
77! lrflx Longwave radiation (degC m/s). !
78! shflx Sensible heat flux (degC m/s). !
79! !
80! Surface air humidity. !
81! !
82! Hair Surface air specific (g/kg) or relative humidity !
83! (percentage). !
84! HairG Latest two-time snapshots of input "Hair" grided !
85! data used for interpolation. !
86! !
87! Surface air pressure. !
88! !
89! Pair Surface air pressure (mb). !
90! PairG Latest two-time snapshots of input "Pair" grided !
91! data used for interpolation. !
92! !
93! Surface air temperature. !
94! !
95! Tair Surface air temperature (Celsius) !
96! TairG Latest two-time snapshots of input "Tair" grided !
97! data used for interpolation. !
98! Surface Winds. !
99! !
100! Uwind Surface wind in the XI-direction (m/s) at !
101! horizontal RHO-points. !
102! UwindG Latest two-time snapshots of input "Uwind" grided !
103! data used for interpolation. !
104! Vwind Surface wind in the ETA-direction (m/s) at !
105! horizontal RHO-points. !
106! VwindG Latest two-time snapshots of input "Vwind" grided !
107! data used for interpolation. !
108! !
109! Rain fall rate. !
110! !
111! evap Evaporation rate (kg/m2/s). !
112! rain Rain fall rate (kg/m2/s). !
113! rainG Latest two-time snapshots of input "rain" grided !
114! data used for interpolation. !
115! !
116! Surface tracer fluxes. !
117! !
118! stflux Forcing surface flux of tracer type variables from !
119! data, coupling, bulk flux parameterization, or !
120! analytical formulas. !
121! !
122! stflux(:,:,itemp) surface net heat flux !
123! stflux(:,:,isalt) surface net freshwater flux (E-P) !
124! !
125! stfluxG Latest two-time snapshots of input "stflux" grided !
126! data used for interpolation. !
127! !
128! stflx ROMS state surface flux of tracer type variables !
129! (TracerUnits m/s) at horizontal RHO-points, as used !
130! in the governing equations. !
131
132#if defined ADJUST_STFLUX && defined FOUR_DVAR
133! !
134! tflux Surface tracer flux adjustments (:,:,Nfrec,2,itrc) !
135#endif
136#if defined ICE_MODEL && defined ICE_BULK_FLUXES && defined BULK_FLUXES
137! !
138! Seaice fluxes. !
139! !
140! Qnet_ai Surface net heat flux at the Air/Ice inteface !
141! Qnet_ao Surface net heat flux at the Air/Ocean inteface !
142! snow Snow fall rate on ice !
143! srflx_ice Shortwave radiation flux on ice !
144! sustr_ai Surface U-momentum stress at the Air/Ice interface !
145! svstr_ai Surface V-momentum stress at the Air/Ice interface !
146! sustr_ao Surface U-momentum stress at the Air/Ocean interface !
147! svstr_ao Surface V-momentum stress at the Air/Ocean interface !
148#endif
149! !
150! Bottom tracer fluxes. !
151! !
152! btflux Forcing bottom flux of tracer type variables from !
153! data or analytical formulas. Usually, the bottom !
154! flux of tracer is zero. !
155! !
156! btflux(:,:,itemp) bottom heat flux !
157! btflux(:,:,isalt) bottom freshwater flux !
158! !
159! btfluxG Latest two-time snapshots of input "vtflux" grided !
160! data used for interpolation. !
161! !
162! btflx ROMS state bottom flux of tracer type variables !
163! (TracerUnits m/s) at horizontal RHO-points, as used !
164! in the governing equations. !
165! !
166! Surface heat flux correction. !
167! !
168! dqdt Surface net heat flux sensitivity to SST, !
169! d(Q)/d(SST), (m/s). !
170! dqdtG Latest two-time snapshots of input "dqdt" grided !
171! data used for interpolation. !
172! sst Sea surface temperature (Celsius). !
173! sstG Latest two-time snapshots of input "sst" grided !
174! data used for interpolation. !
175! !
176! Surface freshwater flux correction. !
177! !
178! sss Sea surface salinity (PSU). !
179! sssG Latest two-time snapshots of input "sss" grided !
180! data used for interpolation. !
181! !
182! Surface spectral downwelling irradiance. !
183! !
184! SpecIr Spectral irradiance (NBands) from 400-700 nm at !
185! 5 nm bandwidth. !
186! avcos Cosine of average zenith angle of downwelling !
187! spectral photons. !
188! !
189!=======================================================================
190!
191 USE mod_kinds
192!
193 implicit none
194!
195 PUBLIC :: allocate_forces
196 PUBLIC :: deallocate_forces
197 PUBLIC :: initialize_forces
198!
199!-----------------------------------------------------------------------
200! Define T_FORCES structure.
201!-----------------------------------------------------------------------
202!
204!
205! Nonlinear model state.
206!
207 real(r8), pointer :: sustr(:,:)
208 real(r8), pointer :: svstr(:,:)
209#if !defined ANA_SMFLUX && !defined BULK_FLUXES || \
210 defined forward_fluxes
211 real(r8), pointer :: sustrg(:,:,:)
212 real(r8), pointer :: svstrg(:,:,:)
213#endif
214#ifdef ADJUST_WSTRESS
215 real(r8), pointer :: ustr(:,:,:,:)
216 real(r8), pointer :: vstr(:,:,:,:)
217#endif
218 real(r8), pointer :: bustr(:,:)
219 real(r8), pointer :: bvstr(:,:)
220#if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS || \
221 defined spectral_light
222 real(r8), pointer :: pair(:,:)
223# ifndef ANA_PAIR
224 real(r8), pointer :: pairg(:,:,:)
225# endif
226#endif
227
228#ifdef WAVES_DIR
229 real(r8), pointer :: dwave(:,:)
230# ifndef ANA_WWAVE
231 real(r8), pointer :: dwaveg(:,:,:)
232# endif
233#endif
234
235#ifdef WAVES_DIRP
236 real(r8), pointer :: dwavep(:,:)
237# ifndef ANA_WWAVE
238 real(r8), pointer :: dwavepg(:,:,:)
239# endif
240#endif
241
242#ifdef WAVES_HEIGHT
243 real(r8), pointer :: hwave(:,:)
244# ifndef ANA_WWAVE
245 real(r8), pointer :: hwaveg(:,:,:)
246# endif
247#endif
248
249#ifdef WAVES_LENGTH
250 real(r8), pointer :: lwave(:,:)
251# ifndef ANA_WWAVE
252 real(r8), pointer :: lwaveg(:,:,:)
253# endif
254#endif
255
256#ifdef WAVES_LENGTHP
257 real(r8), pointer :: lwavep(:,:)
258# ifndef ANA_WWAVE
259 real(r8), pointer :: lwavepg(:,:,:)
260# endif
261#endif
262
263#ifdef WAVES_TOP_PERIOD
264 real(r8), pointer :: pwave_top(:,:)
265# ifndef ANA_WWAVE
266 real(r8), pointer :: pwave_topg(:,:,:)
267# endif
268#endif
269
270#ifdef WAVES_BOT_PERIOD
271 real(r8), pointer :: pwave_bot(:,:)
272# ifndef ANA_WWAVE
273 real(r8), pointer :: pwave_botg(:,:,:)
274# endif
275#endif
276
277#if defined BBL_MODEL || defined SED_BEDLOAD_VANDERA || \
278 defined wav_coupling
279 real(r8), pointer :: uwave_rms(:,:)
280# ifndef ANA_WWAVE
281 real(r8), pointer :: uwave_rmsg(:,:,:)
282# endif
283#endif
284
285# if defined TKE_WAVEDISS || defined WAV_COUPLING || \
286 defined wdiss_thorguza || defined wdiss_churthor || \
287 defined waves_diss || defined wdiss_inwave
288 real(r8), pointer :: dissip_break(:,:)
289 real(r8), pointer :: dissip_wcap(:,:)
290# ifndef ANA_WWAVE
291 real(r8), pointer :: dissip_breakg(:,:,:)
292 real(r8), pointer :: dissip_wcapg(:,:,:)
293# endif
294#endif
295#if defined WAV_COUPLING || (defined WEC_VF && defined BOTTOM_STREAMING)
296 real(r8), pointer :: dissip_fric(:,:)
297# ifndef ANA_WWAVE
298 real(r8), pointer :: dissip_fricg(:,:,:)
299# endif
300#endif
301
302#if defined WEC_ROLLER
303 real(r8), pointer :: dissip_roller(:,:)
304# ifndef ANA_WWAVE
305 real(r8), pointer :: dissip_rollerg(:,:,:)
306# endif
307#endif
308
309#if defined ROLLER_SVENDSEN
310 real(r8), pointer :: wave_break(:,:)
311# ifndef ANA_WWAVE
312 real(r8), pointer :: wave_breakg(:,:,:)
313# endif
314#endif
315
316#if defined WAVES_DSPR
317 real(r8), pointer :: wave_ds(:,:)
318 real(r8), pointer :: wave_qp(:,:)
319# ifndef ANA_WWAVE
320 real(r8), pointer :: wave_dsg(:,:,:)
321 real(r8), pointer :: wave_qpg(:,:,:)
322# endif
323#endif
324
325#if defined WEC_ROLLER
326 real(r8), pointer :: rolla(:,:)
327#endif
328
329#ifdef SPECTRUM_STOKES
330 real(r8), pointer :: spec_wn(:,:,:)
331 real(r8), pointer :: spec_us(:,:,:)
332 real(r8), pointer :: spec_vs(:,:,:)
333#endif
334
335#ifdef SOLVE3D
336
337# ifdef SHORTWAVE
338 real(r8), pointer :: srflx(:,:)
339# ifndef ANA_SRFLUX
340 real(r8), pointer :: srflxg(:,:,:)
341# endif
342# ifdef ICE_ALBEDO
343 real(r8), pointer :: albedo(:,:)
344# ifdef ICE_MODEL
345 real(r8), pointer :: albedo_ice(:,:)
346# endif
347# endif
348# endif
349
350# if defined RED_TIDE && defined DAILY_SHORTWAVE
351 real(r8), pointer :: srflx_avg(:,:)
352 real(r8), pointer :: srflxg_avg(:,:,:)
353# endif
354
355# ifdef CLOUDS
356 real(r8), pointer :: cloud(:,:)
357# ifndef ANA_CLOUD
358 real(r8), pointer :: cloudg(:,:,:)
359# endif
360# endif
361# if defined BULK_FLUXES || defined FRC_COUPLING
362 real(r8), pointer :: lhflx(:,:)
363 real(r8), pointer :: lrflx(:,:)
364 real(r8), pointer :: shflx(:,:)
365# endif
366
367# if !defined LONGWAVE && defined BULK_FLUXES
368 real(r8), pointer :: lrflxg(:,:,:)
369# endif
370
371# if defined BULK_FLUXES || defined ECOSIM || \
372 (defined shortwave && defined ana_srflux) || \
373 defined spectral_light
374 real(r8), pointer :: hair(:,:)
375# ifndef ANA_HUMIDITY
376 real(r8), pointer :: hairg(:,:,:)
377# endif
378 real(r8), pointer :: tair(:,:)
379# ifndef ANA_TAIR
380 real(r8), pointer :: tairg(:,:,:)
381# endif
382# endif
383
384# if defined BULK_FLUXES || defined ECOSIM || \
385 defined spectral_light
386 real(r8), pointer :: uwind(:,:)
387 real(r8), pointer :: vwind(:,:)
388# ifndef ANA_WINDS
389 real(r8), pointer :: uwindg(:,:,:)
390 real(r8), pointer :: vwindg(:,:,:)
391# endif
392# endif
393
394# ifdef BULK_FLUXES
395 real(r8), pointer :: rain(:,:)
396# ifndef ANA_RAIN
397 real(r8), pointer :: raing(:,:,:)
398# endif
399# ifdef EMINUSP
400 real(r8), pointer :: evap(:,:)
401# endif
402# endif
403
404# if defined ICE_MODEL && \
405 defined ice_bulk_fluxes && defined bulk_fluxes
406 real(r8), pointer :: qnet_ai(:,:)
407 real(r8), pointer :: qnet_ao(:,:)
408 real(r8), pointer :: srflx_ice(:,:)
409 real(r8), pointer :: sustr_ao(:,:)
410 real(r8), pointer :: svstr_ao(:,:)
411 real(r8), pointer :: sustr_ai(:,:)
412 real(r8), pointer :: svstr_ai(:,:)
413 real(r8), pointer :: snow(:,:)
414# endif
415
416 real(r8), pointer :: stflux(:,:,:)
417# if !defined ANA_STFLUX || !defined ANA_SSFLUX || \
418 !defined ANA_SPFLUX
419 real(r8), pointer :: stfluxg(:,:,:,:)
420# endif
421 real(r8), pointer :: stflx(:,:,:)
422# ifdef ADJUST_STFLUX
423 real(r8), pointer :: tflux(:,:,:,:,:)
424# endif
425
426 real(r8), pointer :: btflux(:,:,:)
427# if !defined ANA_BTFLUX || !defined ANA_BSFLUX || \
428 !defined ANA_BPFLUX
429 real(r8), pointer :: btfluxg(:,:,:,:)
430# endif
431 real(r8), pointer :: btflx(:,:,:)
432
433# ifdef QCORRECTION
434 real(r8), pointer :: dqdt(:,:)
435 real(r8), pointer :: sst(:,:)
436# ifndef ANA_SST
437 real(r8), pointer :: dqdtg(:,:,:)
438 real(r8), pointer :: sstg(:,:,:)
439# endif
440# endif
441
442# if defined SALINITY && (defined SCORRECTION || defined SRELAXATION)
443 real(r8), pointer :: sss(:,:)
444# ifndef ANA_SSS
445 real(r8), pointer :: sssg(:,:,:)
446# endif
447# endif
448
449# if defined ECOSIM || defined SPECTRAL_LIGHT
450 real(r8), pointer :: specir(:,:,:)
451 real(r8), pointer :: avcos(:,:,:)
452# endif
453#endif
454
455#if defined TANGENT || defined TL_IOMS
456!
457! Tangent linear model state.
458!
459 real(r8), pointer :: tl_sustr(:,:)
460 real(r8), pointer :: tl_svstr(:,:)
461# ifdef ADJUST_WSTRESS
462 real(r8), pointer :: tl_ustr(:,:,:,:)
463 real(r8), pointer :: tl_vstr(:,:,:,:)
464# endif
465 real(r8), pointer :: tl_bustr(:,:)
466 real(r8), pointer :: tl_bvstr(:,:)
467# ifdef SOLVE3D
468 real(r8), pointer :: tl_stflux(:,:,:)
469 real(r8), pointer :: tl_stflx (:,:,:)
470 real(r8), pointer :: tl_btflux(:,:,:)
471 real(r8), pointer :: tl_btflx (:,:,:)
472# ifdef ADJUST_STFLUX
473 real(r8), pointer :: tl_tflux(:,:,:,:,:)
474# endif
475# ifdef SHORTWAVE
476 real(r8), pointer :: tl_srflx(:,:)
477# endif
478# ifdef BULK_FLUXES
479 real(r8), pointer :: tl_lhflx(:,:)
480 real(r8), pointer :: tl_lrflx(:,:)
481 real(r8), pointer :: tl_shflx(:,:)
482# ifdef EMINUSP
483 real(r8), pointer :: tl_evap(:,:)
484# endif
485# endif
486# endif
487#endif
488
489#ifdef ADJOINT
490!
491! Adjoint model state.
492!
493 real(r8), pointer :: ad_sustr(:,:)
494 real(r8), pointer :: ad_svstr(:,:)
495# ifdef ADJUST_WSTRESS
496 real(r8), pointer :: ad_ustr(:,:,:,:)
497 real(r8), pointer :: ad_vstr(:,:,:,:)
498# endif
499 real(r8), pointer :: ad_bustr(:,:)
500 real(r8), pointer :: ad_bvstr(:,:)
501 real(r8), pointer :: ad_bustr_sol(:,:)
502 real(r8), pointer :: ad_bvstr_sol(:,:)
503# ifdef SOLVE3D
504 real(r8), pointer :: ad_stflx(:,:,:)
505 real(r8), pointer :: ad_btflx(:,:,:)
506# ifdef ADJUST_STFLUX
507 real(r8), pointer :: ad_tflux(:,:,:,:,:)
508# endif
509# ifdef SHORTWAVE
510 real(r8), pointer :: ad_srflx(:,:)
511# endif
512# ifdef BULK_FLUXES
513 real(r8), pointer :: ad_lhflx(:,:)
514 real(r8), pointer :: ad_lrflx(:,:)
515 real(r8), pointer :: ad_shflx(:,:)
516# ifdef EMINUSP
517 real(r8), pointer :: ad_evap(:,:)
518# endif
519# endif
520# endif
521#endif
522
523#if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
524!
525! Working arrays to store adjoint impulse forcing, background error
526! covariance, background-error standard deviations, or descent
527! conjugate vectors (directions).
528!
529# if defined FOUR_DVAR || defined IMPULSE
530# ifdef ADJUST_WSTRESS
531 real(r8), pointer :: b_sustr(:,:)
532 real(r8), pointer :: b_svstr(:,:)
533# endif
534# if defined ADJUST_STFLUX && defined SOLVE3D
535 real(r8), pointer :: b_stflx(:,:,:)
536# endif
537# ifdef FOUR_DVAR
538# ifdef ADJUST_WSTRESS
539 real(r8), pointer :: d_sustr(:,:,:)
540 real(r8), pointer :: d_svstr(:,:,:)
541 real(r8), pointer :: e_sustr(:,:)
542 real(r8), pointer :: e_svstr(:,:)
543# endif
544# if defined ADJUST_STFLUX && defined SOLVE3D
545 real(r8), pointer :: d_stflx(:,:,:,:)
546 real(r8), pointer :: e_stflx(:,:,:)
547# endif
548# endif
549# endif
550#endif
551
552 END TYPE t_forces
553!
554 TYPE (t_forces), allocatable :: forces(:)
555!
556 CONTAINS
557!
558 SUBROUTINE allocate_forces (ng, LBi, UBi, LBj, UBj)
559!
560!=======================================================================
561! !
562! This routine allocates all variables in the module for all nested !
563! grids. !
564! !
565!=======================================================================
566!
567 USE mod_param
568#ifdef BIOLOGY
569 USE mod_biology
570#endif
571#if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
572 USE mod_scalars
573#endif
574!
575! Local variable declarations.
576!
577 integer, intent(in) :: ng, lbi, ubi, lbj, ubj
578!
579! Local variable declarations.
580!
581 real(r8) :: size2d
582!
583!-----------------------------------------------------------------------
584! Allocate module variables.
585!-----------------------------------------------------------------------
586!
587 IF (ng.eq.1) allocate ( forces(ngrids) )
588!
589! Set horizontal array size.
590!
591 size2d=real((ubi-lbi+1)*(ubj-lbj+1),r8)
592!
593! Nonlinear model state
594!
595 allocate ( forces(ng) % sustr(lbi:ubi,lbj:ubj) )
596 dmem(ng)=dmem(ng)+size2d
597
598 allocate ( forces(ng) % svstr(lbi:ubi,lbj:ubj) )
599 dmem(ng)=dmem(ng)+size2d
600
601#if !defined ANA_SMFLUX && !defined BULK_FLUXES || \
602 defined forward_fluxes
603 allocate ( forces(ng) % sustrG(lbi:ubi,lbj:ubj,2) )
604 dmem(ng)=dmem(ng)+2.0_r8*size2d
605
606 allocate ( forces(ng) % svstrG(lbi:ubi,lbj:ubj,2) )
607 dmem(ng)=dmem(ng)+2.0_r8*size2d
608#endif
609#ifdef ADJUST_WSTRESS
610 allocate ( forces(ng) % ustr(lbi:ubi,lbj:ubj,nfrec(ng),2) )
611 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng),r8)*size2d
612
613 allocate ( forces(ng) % vstr(lbi:ubi,lbj:ubj,nfrec(ng),2) )
614 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng),r8)*size2d
615#endif
616 allocate ( forces(ng) % bustr(lbi:ubi,lbj:ubj) )
617 dmem(ng)=dmem(ng)+size2d
618
619 allocate ( forces(ng) % bvstr(lbi:ubi,lbj:ubj) )
620 dmem(ng)=dmem(ng)+size2d
621
622#if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS || \
623 defined spectral_light
624 allocate ( forces(ng) % Pair(lbi:ubi,lbj:ubj) )
625 dmem(ng)=dmem(ng)+size2d
626# ifndef ANA_PAIR
627 allocate ( forces(ng) % PairG(lbi:ubi,lbj:ubj,2) )
628 dmem(ng)=dmem(ng)+2.0_r8*size2d
629# endif
630#endif
631
632#ifdef WAVES_DIR
633 allocate ( forces(ng) % Dwave(lbi:ubi,lbj:ubj) )
634 dmem(ng)=dmem(ng)+size2d
635# ifndef ANA_WWAVE
636 allocate ( forces(ng) % DwaveG(lbi:ubi,lbj:ubj,2) )
637 dmem(ng)=dmem(ng)+2.0_r8*size2d
638# endif
639#endif
640
641#ifdef WAVES_DIRP
642 allocate ( forces(ng) % Dwavep(lbi:ubi,lbj:ubj) )
643 dmem(ng)=dmem(ng)+size2d
644# ifndef ANA_WWAVE
645 allocate ( forces(ng) % DwavepG(lbi:ubi,lbj:ubj,2) )
646 dmem(ng)=dmem(ng)+2.0_r8*size2d
647# endif
648#endif
649
650#ifdef WAVES_HEIGHT
651 allocate ( forces(ng) % Hwave(lbi:ubi,lbj:ubj) )
652 dmem(ng)=dmem(ng)+size2d
653# ifndef ANA_WWAVE
654 allocate ( forces(ng) % HwaveG(lbi:ubi,lbj:ubj,2) )
655 dmem(ng)=dmem(ng)+2.0_r8*size2d
656# endif
657#endif
658
659#ifdef WAVES_LENGTH
660 allocate ( forces(ng) % Lwave(lbi:ubi,lbj:ubj) )
661 dmem(ng)=dmem(ng)+size2d
662# ifndef ANA_WWAVE
663 allocate ( forces(ng) % LwaveG(lbi:ubi,lbj:ubj,2) )
664 dmem(ng)=dmem(ng)+2.0_r8*size2d
665# endif
666#endif
667
668#ifdef WAVES_LENGTHP
669 allocate ( forces(ng) % Lwavep(lbi:ubi,lbj:ubj) )
670 dmem(ng)=dmem(ng)+size2d
671# ifndef ANA_WWAVE
672 allocate ( forces(ng) % LwavepG(lbi:ubi,lbj:ubj,2) )
673 dmem(ng)=dmem(ng)+2.0_r8*size2d
674# endif
675#endif
676
677#ifdef WAVES_TOP_PERIOD
678 allocate ( forces(ng) % Pwave_top(lbi:ubi,lbj:ubj) )
679 dmem(ng)=dmem(ng)+size2d
680# ifndef ANA_WWAVE
681 allocate ( forces(ng) % Pwave_topG(lbi:ubi,lbj:ubj,2) )
682 dmem(ng)=dmem(ng)+2.0_r8*size2d
683# endif
684#endif
685
686#ifdef WAVES_BOT_PERIOD
687 allocate ( forces(ng) % Pwave_bot(lbi:ubi,lbj:ubj) )
688 dmem(ng)=dmem(ng)+size2d
689# ifndef ANA_WWAVE
690 allocate ( forces(ng) % Pwave_botG(lbi:ubi,lbj:ubj,2) )
691 dmem(ng)=dmem(ng)+2.0_r8*size2d
692# endif
693#endif
694
695#if defined BBL_MODEL || defined WAV_COUPLING
696 allocate ( forces(ng) % Uwave_rms(lbi:ubi,lbj:ubj) )
697 dmem(ng)=dmem(ng)+size2d
698# ifndef ANA_WWAVE
699 allocate ( forces(ng) % Uwave_rmsG(lbi:ubi,lbj:ubj,2) )
700 dmem(ng)=dmem(ng)+2.0_r8*size2d
701# endif
702#endif
703
704#if defined TKE_WAVEDISS || defined WAV_COUPLING || \
705 defined wdiss_thorguza || defined wdiss_churthor || \
706 defined waves_diss || defined wdiss_inwave
707 allocate ( forces(ng) % Dissip_break(lbi:ubi,lbj:ubj) )
708 dmem(ng)=dmem(ng)+size2d
709 allocate ( forces(ng) % Dissip_wcap(lbi:ubi,lbj:ubj) )
710 dmem(ng)=dmem(ng)+size2d
711# ifndef ANA_WWAVE
712 allocate ( forces(ng) % Dissip_breakG(lbi:ubi,lbj:ubj,2) )
713 dmem(ng)=dmem(ng)+2.0_r8*size2d
714 allocate ( forces(ng) % Dissip_wcapG(lbi:ubi,lbj:ubj,2) )
715 dmem(ng)=dmem(ng)+2.0_r8*size2d
716# endif
717#endif
718
719#if defined WAV_COUPLING || (defined WEC_VF && defined BOTTOM_STREAMING)
720 allocate ( forces(ng) % Dissip_fric(lbi:ubi,lbj:ubj) )
721 dmem(ng)=dmem(ng)+size2d
722# ifndef ANA_WWAVE
723 allocate ( forces(ng) % Dissip_fricG(lbi:ubi,lbj:ubj,2) )
724 dmem(ng)=dmem(ng)+2.0_r8*size2d
725# endif
726#endif
727
728#if defined WEC_ROLLER
729 allocate ( forces(ng) % Dissip_roller(lbi:ubi,lbj:ubj) )
730 dmem(ng)=dmem(ng)+size2d
731# ifndef ANA_WWAVE
732 allocate ( forces(ng) % Dissip_rollerG(lbi:ubi,lbj:ubj,2) )
733 dmem(ng)=dmem(ng)+2.0_r8*size2d
734# endif
735#endif
736
737#if defined ROLLER_SVENDSEN
738 allocate ( forces(ng) % wave_break(lbi:ubi,lbj:ubj) )
739 dmem(ng)=dmem(ng)+size2d
740# ifndef ANA_WWAVE
741 allocate ( forces(ng) % Wave_breakG(lbi:ubi,lbj:ubj,2) )
742 dmem(ng)=dmem(ng)+2.0_r8*size2d
743# endif
744#endif
745
746#if defined WAVES_DSPR
747 allocate ( forces(ng) % Wave_ds(lbi:ubi,lbj:ubj) )
748 dmem(ng)=dmem(ng)+size2d
749 allocate ( forces(ng) % Wave_qp(lbi:ubi,lbj:ubj) )
750 dmem(ng)=dmem(ng)+size2d
751# ifndef ANA_WWAVE
752 allocate ( forces(ng) % Wave_dsG(lbi:ubi,lbj:ubj,2) )
753 dmem(ng)=dmem(ng)+2.0_r8*size2d
754 allocate ( forces(ng) % Wave_qpG(lbi:ubi,lbj:ubj,2) )
755 dmem(ng)=dmem(ng)+2.0_r8*size2d
756# endif
757#endif
758
759#if defined WEC_ROLLER
760 allocate ( forces(ng) % rollA(lbi:ubi,lbj:ubj) )
761 dmem(ng)=dmem(ng)+size2d
762#endif
763
764#ifdef SPECTRUM_STOKES
765 allocate ( forces(ng) % spec_wn(lbi:ubi,lbj:ubj,mscs) )
766 dmem(ng)=dmem(ng)+size2d*real(mscs,r8)
767 allocate ( forces(ng) % spec_us(lbi:ubi,lbj:ubj,mscs) )
768 dmem(ng)=dmem(ng)+size2d*real(mscs,r8)
769 allocate ( forces(ng) % spec_vs(lbi:ubi,lbj:ubj,mscs) )
770 dmem(ng)=dmem(ng)+size2d*real(mscs,r8)
771#endif
772
773#ifdef SOLVE3D
774
775# ifdef SHORTWAVE
776 allocate ( forces(ng) % srflx(lbi:ubi,lbj:ubj) )
777 dmem(ng)=dmem(ng)+size2d
778
779# ifndef ANA_SRFLUX
780 allocate ( forces(ng) % srflxG(lbi:ubi,lbj:ubj,2) )
781 dmem(ng)=dmem(ng)+size2d
782# endif
783# ifdef ICE_ALBEDO
784 allocate ( forces(ng) % albedo(lbi:ubi,lbj:ubj) )
785 dmem(ng)=dmem(ng)+size2d
786# ifdef ICE_MODEL
787 allocate ( forces(ng) % albedo_ice(lbi:ubi,lbj:ubj) )
788 dmem(ng)=dmem(ng)+size2d
789# endif
790# endif
791# endif
792
793# if defined RED_TIDE && defined DAILY_SHORTWAVE
794 allocate ( forces(ng) % srflx_avg(lbi:ubi,lbj:ubj) )
795 dmem(ng)=dmem(ng)+size2d
796
797 allocate ( forces(ng) % srflxG_avg(lbi:ubi,lbj:ubj,2) )
798 dmem(ng)=dmem(ng)+2.0_r8*size2d
799# endif
800
801# ifdef CLOUDS
802 allocate ( forces(ng) % cloud(lbi:ubi,lbj:ubj) )
803 dmem(ng)=dmem(ng)+size2d
804
805# ifndef ANA_CLOUD
806 allocate ( forces(ng) % cloudG(lbi:ubi,lbj:ubj,2) )
807 dmem(ng)=dmem(ng)+2.0_r8*size2d
808# endif
809# endif
810
811# if defined BULK_FLUXES || defined FRC_COUPLING
812 allocate ( forces(ng) % lhflx(lbi:ubi,lbj:ubj) )
813 dmem(ng)=dmem(ng)+size2d
814
815 allocate ( forces(ng) % lrflx(lbi:ubi,lbj:ubj) )
816 dmem(ng)=dmem(ng)+size2d
817
818 allocate ( forces(ng) % shflx(lbi:ubi,lbj:ubj) )
819 dmem(ng)=dmem(ng)+size2d
820# endif
821
822# if !defined LONGWAVE && defined BULK_FLUXES
823 allocate ( forces(ng) % lrflxG(lbi:ubi,lbj:ubj,2) )
824 dmem(ng)=dmem(ng)+2.0_r8*size2d
825# endif
826
827# if defined BULK_FLUXES || defined ECOSIM || \
828 (defined shortwave && defined ana_srflux) || \
829 defined spectral_light
830 allocate ( forces(ng) % Hair(lbi:ubi,lbj:ubj) )
831 dmem(ng)=dmem(ng)+size2d
832
833# ifndef ANA_HUMIDITY
834 allocate ( forces(ng) % HairG(lbi:ubi,lbj:ubj,2) )
835 dmem(ng)=dmem(ng)+2.0_r8*size2d
836# endif
837
838 allocate ( forces(ng) % Tair(lbi:ubi,lbj:ubj) )
839 dmem(ng)=dmem(ng)+size2d
840
841# ifndef ANA_TAIR
842 allocate ( forces(ng) % TairG(lbi:ubi,lbj:ubj,2) )
843 dmem(ng)=dmem(ng)+2.0_r8*size2d
844# endif
845# endif
846
847# if defined BULK_FLUXES || defined ECOSIM || \
848 defined spectral_light
849 allocate ( forces(ng) % Uwind(lbi:ubi,lbj:ubj) )
850 dmem(ng)=dmem(ng)+size2d
851
852 allocate ( forces(ng) % Vwind(lbi:ubi,lbj:ubj) )
853 dmem(ng)=dmem(ng)+size2d
854
855# ifndef ANA_WINDS
856 allocate ( forces(ng) % UwindG(lbi:ubi,lbj:ubj,2) )
857 dmem(ng)=dmem(ng)+2.0_r8*size2d
858
859 allocate ( forces(ng) % VwindG(lbi:ubi,lbj:ubj,2) )
860 dmem(ng)=dmem(ng)+2.0_r8*size2d
861# endif
862# endif
863
864# ifdef BULK_FLUXES
865 allocate ( forces(ng) % rain(lbi:ubi,lbj:ubj) )
866 dmem(ng)=dmem(ng)+size2d
867
868# ifndef ANA_RAIN
869 allocate ( forces(ng) % rainG(lbi:ubi,lbj:ubj,2) )
870 dmem(ng)=dmem(ng)+2.0_r8*size2d
871
872# endif
873# ifdef EMINUSP
874 allocate ( forces(ng) % evap(lbi:ubi,lbj:ubj) )
875 dmem(ng)=dmem(ng)+size2d
876# endif
877# endif
878
879# if defined ICE_MODEL && defined ICE_BULK_FLUXES && defined BULK_FLUXES
880 allocate ( forces(ng) % Qnet_ai(lbi:ubi,lbj:ubj) )
881 dmem(ng)=dmem(ng)+size2d
882
883 allocate ( forces(ng) % Qnet_ao(lbi:ubi,lbj:ubj) )
884 dmem(ng)=dmem(ng)+size2d
885
886 allocate ( forces(ng) % snow(lbi:ubi,lbj:ubj) )
887 dmem(ng)=dmem(ng)+size2d
888
889 allocate ( forces(ng) % srflx_ice(lbi:ubi,lbj:ubj) )
890 dmem(ng)=dmem(ng)+size2d
891
892 allocate ( forces(ng) % sustr_ao(lbi:ubi,lbj:ubj) )
893 dmem(ng)=dmem(ng)+size2d
894
895 allocate ( forces(ng) % svstr_ao(lbi:ubi,lbj:ubj) )
896 dmem(ng)=dmem(ng)+size2d
897
898 allocate ( forces(ng) % sustr_ai(lbi:ubi,lbj:ubj) )
899 dmem(ng)=dmem(ng)+size2d
900
901 allocate ( forces(ng) % svstr_ai(lbi:ubi,lbj:ubj) )
902 dmem(ng)=dmem(ng)+size2d
903# endif
904
905 allocate ( forces(ng) % stflux(lbi:ubi,lbj:ubj,nt(ng)) )
906 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
907
908# if !defined ANA_STFLUX || !defined ANA_SSFLUX || \
909 !defined ANA_SPFLUX
910 allocate ( forces(ng) % stfluxG(lbi:ubi,lbj:ubj,2,nt(ng)) )
911 dmem(ng)=dmem(ng)+2.0_r8*real(nt(ng),r8)*size2d
912
913# endif
914
915 allocate ( forces(ng) % stflx(lbi:ubi,lbj:ubj,nt(ng)) )
916 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
917
918# ifdef ADJUST_STFLUX
919 allocate ( forces(ng) % tflux(lbi:ubi,lbj:ubj,nfrec(ng), &
920 & 2,nt(ng)) )
921 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng)*nt(ng),r8)*size2d
922
923# endif
924
925 allocate ( forces(ng) % btflux(lbi:ubi,lbj:ubj,nt(ng)) )
926 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
927
928# if !defined ANA_BTFLUX || !defined ANA_BSFLUX || \
929 !defined ANA_BPFLUX
930 allocate ( forces(ng) % btfluxG(lbi:ubi,lbj:ubj,2,nt(ng)) )
931 dmem(ng)=dmem(ng)+2.0_r8*real(nt(ng),r8)*size2d
932# endif
933
934 allocate ( forces(ng) % btflx(lbi:ubi,lbj:ubj,nt(ng)) )
935 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
936
937# ifdef QCORRECTION
938 allocate ( forces(ng) % dqdt(lbi:ubi,lbj:ubj) )
939 dmem(ng)=dmem(ng)+size2d
940
941 allocate ( forces(ng) % sst(lbi:ubi,lbj:ubj) )
942 dmem(ng)=dmem(ng)+size2d
943
944# ifndef ANA_SST
945 allocate ( forces(ng) % dqdtG(lbi:ubi,lbj:ubj,2) )
946 dmem(ng)=dmem(ng)+2.0_r8*size2d
947
948 allocate ( forces(ng) % sstG(lbi:ubi,lbj:ubj,2) )
949 dmem(ng)=dmem(ng)+2.0_r8*size2d
950# endif
951# endif
952
953# if defined SALINITY && (defined SCORRECTION || defined SRELAXATION)
954 allocate ( forces(ng) % sss(lbi:ubi,lbj:ubj) )
955 dmem(ng)=dmem(ng)+size2d
956
957# ifndef ANA_SSS
958 allocate ( forces(ng) % sssG(lbi:ubi,lbj:ubj,2) )
959 dmem(ng)=dmem(ng)+2.0_r8*size2d
960# endif
961# endif
962
963# if defined ECOSIM || defined SPECTRAL_LIGHT
964 allocate ( forces(ng) % SpecIr(lbi:ubi,lbj:ubj,nbands) )
965 dmem(ng)=dmem(ng)+real(nbands,r8)*size2d
966
967 allocate ( forces(ng) % avcos(lbi:ubi,lbj:ubj,nbands) )
968 dmem(ng)=dmem(ng)+real(nbands,r8)*size2d
969# endif
970
971#endif
972
973#if defined TANGENT || defined TL_IOMS
974!
975! Tangent linear model state
976!
977 allocate ( forces(ng) % tl_sustr(lbi:ubi,lbj:ubj) )
978 dmem(ng)=dmem(ng)+size2d
979
980 allocate ( forces(ng) % tl_svstr(lbi:ubi,lbj:ubj) )
981 dmem(ng)=dmem(ng)+size2d
982# ifdef ADJUST_WSTRESS
983 allocate ( forces(ng) % tl_ustr(lbi:ubi,lbj:ubj,nfrec(ng),2) )
984 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng),r8)*size2d
985
986 allocate ( forces(ng) % tl_vstr(lbi:ubi,lbj:ubj,nfrec(ng),2) )
987 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng),r8)*size2d
988# endif
989 allocate ( forces(ng) % tl_bustr(lbi:ubi,lbj:ubj) )
990 dmem(ng)=dmem(ng)+size2d
991
992 allocate ( forces(ng) % tl_bvstr(lbi:ubi,lbj:ubj) )
993 dmem(ng)=dmem(ng)+size2d
994# ifdef SOLVE3D
995 allocate ( forces(ng) % tl_stflux(lbi:ubi,lbj:ubj,nt(ng)) )
996 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
997
998 allocate ( forces(ng) % tl_stflx(lbi:ubi,lbj:ubj,nt(ng)) )
999 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
1000
1001 allocate ( forces(ng) % tl_btflux(lbi:ubi,lbj:ubj,nt(ng)) )
1002 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
1003
1004 allocate ( forces(ng) % tl_btflx(lbi:ubi,lbj:ubj,nt(ng)) )
1005 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
1006# ifdef ADJUST_STFLUX
1007 allocate ( forces(ng) % tl_tflux(lbi:ubi,lbj:ubj,nfrec(ng), &
1008 & 2,nt(ng)) )
1009 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng)*nt(ng),r8)*size2d
1010# endif
1011# ifdef SHORTWAVE
1012 allocate ( forces(ng) % tl_srflx(lbi:ubi,lbj:ubj) )
1013 dmem(ng)=dmem(ng)+size2d
1014# endif
1015# ifdef BULK_FLUXES
1016 allocate ( forces(ng) % tl_lhflx(lbi:ubi,lbj:ubj) )
1017 dmem(ng)=dmem(ng)+size2d
1018
1019 allocate ( forces(ng) % tl_lrflx(lbi:ubi,lbj:ubj) )
1020 dmem(ng)=dmem(ng)+size2d
1021
1022 allocate ( forces(ng) % tl_shflx(lbi:ubi,lbj:ubj) )
1023 dmem(ng)=dmem(ng)+size2d
1024# ifdef EMINUSP
1025 allocate ( forces(ng) % tl_evap(lbi:ubi,lbj:ubj) )
1026 dmem(ng)=dmem(ng)+size2d
1027# endif
1028# endif
1029# endif
1030#endif
1031
1032#ifdef ADJOINT
1033!
1034! Adjoint model state
1035!
1036 allocate ( forces(ng) % ad_sustr(lbi:ubi,lbj:ubj) )
1037 dmem(ng)=dmem(ng)+size2d
1038
1039 allocate ( forces(ng) % ad_svstr(lbi:ubi,lbj:ubj) )
1040 dmem(ng)=dmem(ng)+size2d
1041
1042# ifdef ADJUST_WSTRESS
1043 allocate ( forces(ng) % ad_ustr(lbi:ubi,lbj:ubj,nfrec(ng),2) )
1044 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng),r8)*size2d
1045
1046 allocate ( forces(ng) % ad_vstr(lbi:ubi,lbj:ubj,nfrec(ng),2) )
1047 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng),r8)*size2d
1048# endif
1049 allocate ( forces(ng) % ad_bustr(lbi:ubi,lbj:ubj) )
1050 dmem(ng)=dmem(ng)+size2d
1051
1052 allocate ( forces(ng) % ad_bvstr(lbi:ubi,lbj:ubj) )
1053 dmem(ng)=dmem(ng)+size2d
1054
1055 allocate ( forces(ng) % ad_bustr_sol(lbi:ubi,lbj:ubj) )
1056 dmem(ng)=dmem(ng)+size2d
1057
1058 allocate ( forces(ng) % ad_bvstr_sol(lbi:ubi,lbj:ubj) )
1059 dmem(ng)=dmem(ng)+size2d
1060# ifdef SOLVE3D
1061 allocate ( forces(ng) % ad_stflx(lbi:ubi,lbj:ubj,nt(ng)) )
1062 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
1063
1064 allocate ( forces(ng) % ad_btflx(lbi:ubi,lbj:ubj,nt(ng)) )
1065 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
1066# ifdef ADJUST_STFLUX
1067 allocate ( forces(ng) % ad_tflux(lbi:ubi,lbj:ubj,nfrec(ng), &
1068 & 2,nt(ng)) )
1069 dmem(ng)=dmem(ng)+2.0_r8*real(nfrec(ng)*nt(ng),r8)*size2d
1070# endif
1071# ifdef SHORTWAVE
1072 allocate ( forces(ng) % ad_srflx(lbi:ubi,lbj:ubj) )
1073 dmem(ng)=dmem(ng)+size2d
1074# endif
1075# ifdef BULK_FLUXES
1076 allocate ( forces(ng) % ad_lhflx(lbi:ubi,lbj:ubj) )
1077 dmem(ng)=dmem(ng)+size2d
1078
1079 allocate ( forces(ng) % ad_lrflx(lbi:ubi,lbj:ubj) )
1080 dmem(ng)=dmem(ng)+size2d
1081
1082 allocate ( forces(ng) % ad_shflx(lbi:ubi,lbj:ubj) )
1083 dmem(ng)=dmem(ng)+size2d
1084# ifdef EMINUSP
1085 allocate ( forces(ng) % ad_evap(lbi:ubi,lbj:ubj) )
1086 dmem(ng)=dmem(ng)+size2d
1087# endif
1088# endif
1089# endif
1090#endif
1091
1092#if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
1093!
1094! Working arrays to store adjoint impulse forcing, background error
1095! covariance, background-error standard deviations, or descent
1096! conjugate vectors (directions).
1097!
1098# if defined FOUR_DVAR || defined IMPULSE
1099# ifdef ADJUST_WSTRESS
1100 allocate ( forces(ng) % b_sustr(lbi:ubi,lbj:ubj) )
1101 dmem(ng)=dmem(ng)+size2d
1102
1103 allocate ( forces(ng) % b_svstr(lbi:ubi,lbj:ubj) )
1104 dmem(ng)=dmem(ng)+size2d
1105# endif
1106# if defined ADJUST_STFLUX && defined SOLVE3D
1107 allocate ( forces(ng) % b_stflx(lbi:ubi,lbj:ubj,nt(ng)) )
1108 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
1109# endif
1110# endif
1111# ifdef FOUR_DVAR
1112# ifdef ADJUST_WSTRESS
1113 allocate ( forces(ng) % d_sustr(lbi:ubi,lbj:ubj,nfrec(ng)) )
1114 dmem(ng)=dmem(ng)+real(nfrec(ng),r8)*size2d
1115
1116 allocate ( forces(ng) % d_svstr(lbi:ubi,lbj:ubj,nfrec(ng)) )
1117 dmem(ng)=dmem(ng)+real(nfrec(ng),r8)*size2d
1118
1119 allocate ( forces(ng) % e_sustr(lbi:ubi,lbj:ubj) )
1120 dmem(ng)=dmem(ng)+size2d
1121
1122 allocate ( forces(ng) % e_svstr(lbi:ubi,lbj:ubj) )
1123 dmem(ng)=dmem(ng)+size2d
1124# endif
1125# if defined ADJUST_STFLUX && defined SOLVE3D
1126 allocate ( forces(ng) % d_stflx(lbi:ubi,lbj:ubj, &
1127 & nfrec(ng),nt(ng)) )
1128 dmem(ng)=dmem(ng)+real(nfrec(ng)*nt(ng),r8)*size2d
1129
1130 allocate ( forces(ng) % e_stflx(lbi:ubi,lbj:ubj,nt(ng)) )
1131 dmem(ng)=dmem(ng)+real(nt(ng),r8)*size2d
1132# endif
1133# endif
1134#endif
1135!
1136 RETURN
1137 END SUBROUTINE allocate_forces
1138!
1139 SUBROUTINE deallocate_forces (ng)
1140!
1141!=======================================================================
1142! !
1143! This routine deallocates all variables in the module for all nested !
1144! grids. !
1145! !
1146!=======================================================================
1147!
1148 USE mod_param, ONLY : ngrids
1149#ifdef SUBOBJECT_DEALLOCATION
1150 USE destroy_mod, ONLY : destroy
1151#endif
1152!
1153! Imported variable declarations.
1154!
1155 integer, intent(in) :: ng
1156!
1157! Local variable declarations.
1158!
1159 character (len=*), parameter :: myfile = &
1160 & __FILE__//", deallocate_forces"
1161
1162#ifdef SUBOBJECT_DEALLOCATION
1163!
1164!-----------------------------------------------------------------------
1165! Deallocate each variable in the derived-type T_FORCES structure
1166! separately.
1167!-----------------------------------------------------------------------
1168!
1169! Nonlinear model state
1170!
1171 IF (.not.destroy(ng, forces(ng)%sustr, myfile, &
1172 & __line__, 'FORCES(ng)%sustr')) RETURN
1173
1174 IF (.not.destroy(ng, forces(ng)%svstr, myfile, &
1175 & __line__, 'FORCES(ng)%svstr')) RETURN
1176
1177# if !defined ANA_SMFLUX && !defined BULK_FLUXES || \
1178 defined forward_fluxes
1179 IF (.not.destroy(ng, forces(ng)%sustrG, myfile, &
1180 & __line__, 'FORCES(ng)%sustrG')) RETURN
1181
1182 IF (.not.destroy(ng, forces(ng)%svstrG, myfile, &
1183 & __line__, 'FORCES(ng)%svstrG')) RETURN
1184# endif
1185
1186# ifdef ADJUST_WSTRESS
1187 IF (.not.destroy(ng, forces(ng)%ustr, myfile, &
1188 & __line__, 'FORCES(ng)%ustr')) RETURN
1189
1190 IF (.not.destroy(ng, forces(ng)%vstr, myfile, &
1191 & __line__, 'FORCES(ng)%vstr')) RETURN
1192# endif
1193
1194 IF (.not.destroy(ng, forces(ng)%bustr, myfile, &
1195 & __line__, 'FORCES(ng)%bustr')) RETURN
1196
1197 IF (.not.destroy(ng, forces(ng)%bvstr, myfile, &
1198 & __line__, 'FORCES(ng)%bvstr')) RETURN
1199
1200# if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS || \
1201 defined spectral_light
1202 IF (.not.destroy(ng, forces(ng)%Pair, myfile, &
1203 & __line__, 'FORCES(ng)%Pair')) RETURN
1204# ifndef ANA_PAIR
1205 IF (.not.destroy(ng, forces(ng)%PairG, myfile, &
1206 & __line__, 'FORCES(ng)%PairG')) RETURN
1207# endif
1208# endif
1209
1210# ifdef WAVES_DIR
1211 IF (.not.destroy(ng, forces(ng)%Dwave, myfile, &
1212 & __line__, 'FORCES(ng)%Dwave')) RETURN
1213# ifndef ANA_WWAVE
1214 IF (.not.destroy(ng, forces(ng)%DwaveG, myfile, &
1215 & __line__, 'FORCES(ng)%DwaveG')) RETURN
1216# endif
1217# endif
1218
1219# ifdef WAVES_DIRP
1220 IF (.not.destroy(ng, forces(ng)%Dwavep, myfile, &
1221 & __line__, 'FORCES(ng)%Dwavep')) RETURN
1222# ifndef ANA_WWAVE
1223 IF (.not.destroy(ng, forces(ng)%DwavepG, myfile, &
1224 & __line__, 'FORCES(ng)%DwavepG')) RETURN
1225# endif
1226# endif
1227
1228# ifdef WAVES_HEIGHT
1229 IF (.not.destroy(ng, forces(ng)%Hwave, myfile, &
1230 & __line__, 'FORCES(ng)%Hwave')) RETURN
1231# ifndef ANA_WWAVE
1232 IF (.not.destroy(ng, forces(ng)%HwaveG, myfile, &
1233 & __line__, 'FORCES(ng)%HwaveG')) RETURN
1234# endif
1235# endif
1236
1237# ifdef WAVES_LENGTH
1238 IF (.not.destroy(ng, forces(ng)%Lwave, myfile, &
1239 & __line__, 'FORCES(ng)%Lwave')) RETURN
1240# ifndef ANA_WWAVE
1241 IF (.not.destroy(ng, forces(ng)%LwaveG, myfile, &
1242 & __line__, 'FORCES(ng)%LwaveG')) RETURN
1243# endif
1244# endif
1245
1246# ifdef WAVES_LENGTHP
1247 IF (.not.destroy(ng, forces(ng)%Lwavep, myfile, &
1248 & __line__, 'FORCES(ng)%Lwavep')) RETURN
1249# ifndef ANA_WWAVE
1250 IF (.not.destroy(ng, forces(ng)%LwavepG, myfile, &
1251 & __line__, 'FORCES(ng)%LwavepG')) RETURN
1252# endif
1253# endif
1254
1255# ifdef WAVES_TOP_PERIOD
1256 IF (.not.destroy(ng, forces(ng)%Pwave_top, myfile, &
1257 & __line__, 'FORCES(ng)%Pwave_top')) RETURN
1258# ifndef ANA_WWAVE
1259 IF (.not.destroy(ng, forces(ng)%Pwave_topG, myfile, &
1260 & __line__, 'FORCES(ng)%Pwave_topG')) RETURN
1261# endif
1262# endif
1263
1264# ifdef WAVES_BOT_PERIOD
1265 IF (.not.destroy(ng, forces(ng)%Pwave_bot, myfile, &
1266 & __line__, 'FORCES(ng)%Pwave_bot')) RETURN
1267# ifndef ANA_WWAVE
1268 IF (.not.destroy(ng, forces(ng)%Pwave_botG, myfile, &
1269 & __line__, 'FORCES(ng)%Pwave_botG')) RETURN
1270# endif
1271# endif
1272
1273# if defined BBL_MODEL || defined SED_BEDLOAD_VANDERA || \
1274 defined wav_coupling
1275 IF (.not.destroy(ng, forces(ng)%Uwave_rms, myfile, &
1276 & __line__, 'FORCES(ng)%Uwave_rms')) RETURN
1277# ifndef ANA_WWAVE
1278 IF (.not.destroy(ng, forces(ng)%Uwave_rmsG, myfile, &
1279 & __line__, 'FORCES(ng)%Uwave_rmsG')) RETURN
1280# endif
1281# endif
1282
1283# if defined TKE_WAVEDISS || defined WAV_COUPLING || \
1284 defined wdiss_thorguza || defined wdiss_churthor || \
1285 defined waves_diss || defined wdiss_inwave
1286 IF (.not.destroy(ng, forces(ng)%Dissip_break, myfile, &
1287 & __line__, 'FORCES(ng)%Dissip_break')) RETURN
1288 IF (.not.destroy(ng, forces(ng)%Dissip_wcap, myfile, &
1289 & __line__, 'FORCES(ng)%Dissip_wcap')) RETURN
1290# ifndef ANA_WWAVE
1291 IF (.not.destroy(ng, forces(ng)%Dissip_breakG, myfile, &
1292 & __line__, 'FORCES(ng)%Dissip_breakG')) RETURN
1293 IF (.not.destroy(ng, forces(ng)%Dissip_wcapG, myfile, &
1294 & __line__, 'FORCES(ng)%Dissip_wcapG')) RETURN
1295# endif
1296# endif
1297# if defined WAV_COUPLING || (defined WEC_VF && defined BOTTOM_STREAMING)
1298 IF (.not.destroy(ng, forces(ng)%Dissip_fric, myfile, &
1299 & __line__, 'FORCES(ng)%Dissip_fric')) RETURN
1300# ifndef ANA_WWAVE
1301 IF (.not.destroy(ng, forces(ng)%Dissip_fricG, myfile, &
1302 & __line__, 'FORCES(ng)%Dissip_fricG')) RETURN
1303# endif
1304# endif
1305
1306# if defined WEC_ROLLER
1307 IF (.not.destroy(ng, forces(ng)%Dissip_roller, myfile, &
1308 & __line__, 'FORCES(ng)%Dissip_roller')) RETURN
1309 IF (.not.destroy(ng, forces(ng)%rollA, myfile, &
1310 & __line__, 'FORCES(ng)%rollA')) RETURN
1311# ifndef ANA_WWAVE
1312 IF (.not.destroy(ng, forces(ng)%Dissip_rollerG, myfile, &
1313 & __line__, 'FORCES(ng)%Dissip_rollerG')) RETURN
1314 IF (.not.destroy(ng, forces(ng)%rollAG, myfile, &
1315 & __line__, 'FORCES(ng)%rollAG')) RETURN
1316# endif
1317# endif
1318
1319# if defined ROLLER_SVENDSEN
1320 IF (.not.destroy(ng, forces(ng)%wave_break, myfile, &
1321 & __line__, 'FORCES(ng)%wave_break')) RETURN
1322# ifndef ANA_WWAVE
1323 IF (.not.destroy(ng, forces(ng)%Wave_breakG, myfile, &
1324 & __line__, 'FORCES(ng)%Wave_breakG')) RETURN
1325# endif
1326# endif
1327
1328# if defined WAVES_DSPR
1329 IF (.not.destroy(ng, forces(ng)%Wave_ds, myfile, &
1330 & __line__, 'FORCES(ng)%Wave_ds')) RETURN
1331 IF (.not.destroy(ng, forces(ng)%Wave_qp, myfile, &
1332 & __line__, 'FORCES(ng)%Wave_qp')) RETURN
1333# ifndef ANA_WWAVE
1334 IF (.not.destroy(ng, forces(ng)%Wave_dsG, myfile, &
1335 & __line__, 'FORCES(ng)%Wwave_dsG')) RETURN
1336 IF (.not.destroy(ng, forces(ng)%Wave_qpG, myfile, &
1337 & __line__, 'FORCES(ng)%Wave_qpG')) RETURN
1338# endif
1339# endif
1340
1341# ifdef SPECTRUM_STOKES
1342 IF (.not.destroy(ng, forces(ng)%spec_wn, myfile, &
1343 & __line__, 'FORCES(ng)%spec_wn')) RETURN
1344 IF (.not.destroy(ng, forces(ng)%spec_us, myfile, &
1345 & __line__, 'FORCES(ng)%spec_us')) RETURN
1346 IF (.not.destroy(ng, forces(ng)%spec_vs, myfile, &
1347 & __line__, 'FORCES(ng)%spec_vs')) RETURN
1348# endif
1349
1350# ifdef SOLVE3D
1351
1352# ifdef SHORTWAVE
1353 IF (.not.destroy(ng, forces(ng)%srflx, myfile, &
1354 & __line__, 'FORCES(ng)%srflx')) RETURN
1355# ifndef ANA_SRFLUX
1356 IF (.not.destroy(ng, forces(ng)%srflxG, myfile, &
1357 & __line__, 'FORCES(ng)%srflxG')) RETURN
1358# endif
1359
1360# ifdef ICE_ALBEDO
1361 IF (.not.destroy(ng, forces(ng)%albedo, myfile, &
1362 & __line__, 'FORCES(ng)%albedo')) RETURN
1363# ifdef ICE_MODEL
1364 IF (.not.destroy(ng, forces(ng)%albedo_ice, myfile, &
1365 & __line__, 'FORCES(ng)%albedo_ice')) RETURN
1366# endif
1367# endif
1368# endif
1369
1370# if defined RED_TIDE && defined DAILY_SHORTWAVE
1371 IF (.not.destroy(ng, forces(ng)%srflx_avg, myfile, &
1372 & __line__, 'FORCES(ng)%srflx_avg')) RETURN
1373
1374 IF (.not.destroy(ng, forces(ng)%srflxG_avg, myfile, &
1375 & __line__, 'FORCES(ng)%srflxG_avg')) RETURN
1376# endif
1377
1378# ifdef CLOUDS
1379 IF (.not.destroy(ng, forces(ng)%cloud, myfile, &
1380 & __line__, 'FORCES(ng)%cloud')) RETURN
1381# ifndef ANA_CLOUD
1382 IF (.not.destroy(ng, forces(ng)%cloudG, myfile, &
1383 & __line__, 'FORCES(ng)%cloudG')) RETURN
1384# endif
1385# endif
1386
1387# if defined BULK_FLUXES || defined FRC_COUPLING
1388 IF (.not.destroy(ng, forces(ng)%lhflx, myfile, &
1389 & __line__, 'FORCES(ng)%lhflx')) RETURN
1390
1391 IF (.not.destroy(ng, forces(ng)%lrflx, myfile, &
1392 & __line__, 'FORCES(ng)%lrflx')) RETURN
1393
1394 IF (.not.destroy(ng, forces(ng)%shflx, myfile, &
1395 & __line__, 'FORCES(ng)%shflx')) RETURN
1396# endif
1397
1398# if !defined LONGWAVE && defined BULK_FLUXES
1399 IF (.not.destroy(ng, forces(ng)%lrflxG, myfile, &
1400 & __line__, 'FORCES(ng)%lrflxG')) RETURN
1401# endif
1402
1403# if defined BULK_FLUXES || defined ECOSIM || \
1404 (defined shortwave && defined ana_srflux) || \
1405 defined spectral_light
1406 IF (.not.destroy(ng, forces(ng)%Hair, myfile, &
1407 & __line__, 'FORCES(ng)%Hair')) RETURN
1408
1409# ifndef ANA_HUMIDITY
1410 IF (.not.destroy(ng, forces(ng)%HairG, myfile, &
1411 & __line__, 'FORCES(ng)%HairG')) RETURN
1412# endif
1413
1414 IF (.not.destroy(ng, forces(ng)%Tair, myfile, &
1415 & __line__, 'FORCES(ng)%Tair')) RETURN
1416
1417# ifndef ANA_TAIR
1418 IF (.not.destroy(ng, forces(ng)%TairG, myfile, &
1419 & __line__, 'FORCES(ng)%TairG')) RETURN
1420# endif
1421# endif
1422
1423# if defined BULK_FLUXES || defined ECOSIM || \
1424 defined spectral_light
1425 IF (.not.destroy(ng, forces(ng)%Uwind, myfile, &
1426 & __line__, 'FORCES(ng)%Uwind')) RETURN
1427
1428 IF (.not.destroy(ng, forces(ng)%Vwind, myfile, &
1429 & __line__, 'FORCES(ng)%Vwind')) RETURN
1430
1431# ifndef ANA_WINDS
1432 IF (.not.destroy(ng, forces(ng)%UwindG, myfile, &
1433 & __line__, 'FORCES(ng)%UwindG')) RETURN
1434
1435 IF (.not.destroy(ng, forces(ng)%VwindG, myfile, &
1436 & __line__, 'FORCES(ng)%VwindG')) RETURN
1437# endif
1438# endif
1439
1440# ifdef BULK_FLUXES
1441 IF (.not.destroy(ng, forces(ng)%rain, myfile, &
1442 & __line__, 'FORCES(ng)%rain')) RETURN
1443
1444# ifndef ANA_RAIN
1445 IF (.not.destroy(ng, forces(ng)%rainG, myfile, &
1446 & __line__, 'FORCES(ng)%rainG')) RETURN
1447# endif
1448
1449# ifdef EMINUSP
1450 IF (.not.destroy(ng, forces(ng)%evap, myfile, &
1451 & __line__, 'FORCES(ng)%evap')) RETURN
1452# endif
1453# endif
1454
1455# if defined ICE_MODEL && \
1456 defined ice_bulk_fluxes && defined bulk_fluxes
1457 IF (.not.destroy(ng, forces(ng)%Qnet_ai, myfile, &
1458 & __line__, 'FORCES(ng)%Qnet_ai')) RETURN
1459
1460 IF (.not.destroy(ng, forces(ng)%Qnet_ao, myfile, &
1461 & __line__, 'FORCES(ng)%Qnet_ao')) RETURN
1462
1463 IF (.not.destroy(ng, forces(ng)%snow, myfile, &
1464 & __line__, 'FORCES(ng)%snow')) RETURN
1465
1466 IF (.not.destroy(ng, forces(ng)%srflx_ice, myfile, &
1467 & __line__, 'FORCES(ng)%srflx_ice')) RETURN
1468
1469 IF (.not.destroy(ng, forces(ng)%sustr_ai, myfile, &
1470 & __line__, 'FORCES(ng)%sustr_ai')) RETURN
1471
1472 IF (.not.destroy(ng, forces(ng)%svstr_ai, myfile, &
1473 & __line__, 'FORCES(ng)%svstr_ai')) RETURN
1474
1475 IF (.not.destroy(ng, forces(ng)%sustr_ao, myfile, &
1476 & __line__, 'FORCES(ng)%sustr_ao')) RETURN
1477
1478 IF (.not.destroy(ng, forces(ng)%svstr_ao, myfile, &
1479 & __line__, 'FORCES(ng)%svstr_ao')) RETURN
1480# endif
1481
1482 IF (.not.destroy(ng, forces(ng)%stflux, myfile, &
1483 & __line__, 'FORCES(ng)%stflux')) RETURN
1484
1485# if !defined ANA_STFLUX || !defined ANA_SSFLUX || \
1486 !defined ANA_SPFLUX
1487 IF (.not.destroy(ng, forces(ng)%stfluxG, myfile, &
1488 & __line__, 'FORCES(ng)%stfluxG')) RETURN
1489# endif
1490
1491 IF (.not.destroy(ng, forces(ng)%stflx, myfile, &
1492 & __line__, 'FORCES(ng)%stflx')) RETURN
1493
1494# ifdef ADJUST_STFLUX
1495 IF (.not.destroy(ng, forces(ng)%tflux, myfile, &
1496 & __line__, 'FORCES(ng)%tflux')) RETURN
1497# endif
1498
1499 IF (.not.destroy(ng, forces(ng)%btflux, myfile, &
1500 & __line__, 'FORCES(ng)%btflux')) RETURN
1501
1502# if !defined ANA_BTFLUX || !defined ANA_BSFLUX || \
1503 !defined ANA_BPFLUX
1504 IF (.not.destroy(ng, forces(ng)%btfluxG, myfile, &
1505 & __line__, 'FORCES(ng)%btfluxG')) RETURN
1506# endif
1507
1508 IF (.not.destroy(ng, forces(ng)%btflx, myfile, &
1509 & __line__, 'FORCES(ng)%btflx')) RETURN
1510
1511# ifdef QCORRECTION
1512 IF (.not.destroy(ng, forces(ng)%dqdt, myfile, &
1513 & __line__, 'FORCES(ng)%dqdt')) RETURN
1514
1515 IF (.not.destroy(ng, forces(ng)%sst, myfile, &
1516 & __line__, 'FORCES(ng)%sst')) RETURN
1517
1518# ifndef ANA_SST
1519 IF (.not.destroy(ng, forces(ng)%dqdtG, myfile, &
1520 & __line__, 'FORCES(ng)%dqdtG')) RETURN
1521
1522 IF (.not.destroy(ng, forces(ng)%sstG, myfile, &
1523 & __line__, 'FORCES(ng)%sstG')) RETURN
1524# endif
1525# endif
1526
1527# if defined SALINITY && (defined SCORRECTION || defined SRELAXATION)
1528 IF (.not.destroy(ng, forces(ng)%sss, myfile, &
1529 & __line__, 'FORCES(ng)%sss')) RETURN
1530
1531# ifndef ANA_SSS
1532 IF (.not.destroy(ng, forces(ng)%sssG, myfile, &
1533 & __line__, 'FORCES(ng)%sssG')) RETURN
1534# endif
1535# endif
1536
1537# if defined ECOSIM || defined SPECTRAL_LIGHT
1538 IF (.not.destroy(ng, forces(ng)%SpecIr, myfile, &
1539 & __line__, 'FORCES(ng)%SpecIr')) RETURN
1540
1541 IF (.not.destroy(ng, forces(ng)%avcos, myfile, &
1542 & __line__, 'FORCES(ng)%avcos')) RETURN
1543# endif
1544
1545# endif
1546
1547# if defined TANGENT || defined TL_IOMS
1548!
1549! Tangent linear model state
1550!
1551 IF (.not.destroy(ng, forces(ng)%tl_sustr, myfile, &
1552 & __line__, 'FORCES(ng)%tl_sustr')) RETURN
1553
1554 IF (.not.destroy(ng, forces(ng)%tl_svstr, myfile, &
1555 & __line__, 'FORCES(ng)%tl_svstr')) RETURN
1556
1557# ifdef ADJUST_WSTRESS
1558 IF (.not.destroy(ng, forces(ng)%tl_ustr, myfile, &
1559 & __line__, 'FORCES(ng)%tl_ustr')) RETURN
1560
1561 IF (.not.destroy(ng, forces(ng)%tl_vstr, myfile, &
1562 & __line__, 'FORCES(ng)%tl_vstr')) RETURN
1563# endif
1564
1565 IF (.not.destroy(ng, forces(ng)%tl_bustr, myfile, &
1566 & __line__, 'FORCES(ng)%tl_bustr')) RETURN
1567
1568 IF (.not.destroy(ng, forces(ng)%tl_bvstr, myfile, &
1569 & __line__, 'FORCES(ng)%tl_bvstr')) RETURN
1570
1571# ifdef SOLVE3D
1572 IF (.not.destroy(ng, forces(ng)%tl_stflux, myfile, &
1573 & __line__, 'FORCES(ng)%tl_stflux')) RETURN
1574
1575 IF (.not.destroy(ng, forces(ng)%tl_stflx, myfile, &
1576 & __line__, 'FORCES(ng)%tl_stflx')) RETURN
1577
1578 IF (.not.destroy(ng, forces(ng)%tl_btflux, myfile, &
1579 & __line__, 'FORCES(ng)%tl_btflux')) RETURN
1580
1581 IF (.not.destroy(ng, forces(ng)%tl_btflx, myfile, &
1582 & __line__, 'FORCES(ng)%tl_btflx')) RETURN
1583
1584# ifdef ADJUST_STFLUX
1585 IF (.not.destroy(ng, forces(ng)%tl_tflux, myfile, &
1586 & __line__, 'FORCES(ng)%tl_tflux')) RETURN
1587# endif
1588
1589# ifdef SHORTWAVE
1590 IF (.not.destroy(ng, forces(ng)%tl_srflx, myfile, &
1591 & __line__, 'FORCES(ng)%tl_srflx')) RETURN
1592# endif
1593
1594# ifdef BULK_FLUXES
1595 IF (.not.destroy(ng, forces(ng)%tl_lhflx, myfile, &
1596 & __line__, 'FORCES(ng)%tl_lhflx')) RETURN
1597
1598 IF (.not.destroy(ng, forces(ng)%tl_lrflx, myfile, &
1599 & __line__, 'FORCES(ng)%tl_lrflx')) RETURN
1600
1601 IF (.not.destroy(ng, forces(ng)%tl_shflx, myfile, &
1602 & __line__, 'FORCES(ng)%tl_shflx')) RETURN
1603
1604# ifdef EMINUSP
1605 IF (.not.destroy(ng, forces(ng)%tl_evap, myfile, &
1606 & __line__, 'FORCES(ng)%tl_evap')) RETURN
1607# endif
1608# endif
1609# endif
1610# endif
1611
1612# ifdef ADJOINT
1613!
1614! Adjoint model state
1615!
1616 IF (.not.destroy(ng, forces(ng)%ad_sustr, myfile, &
1617 & __line__, 'FORCES(ng)%ad_sustr')) RETURN
1618
1619 IF (.not.destroy(ng, forces(ng)%ad_svstr, myfile, &
1620 & __line__, 'FORCES(ng)%ad_svstr')) RETURN
1621
1622# ifdef ADJUST_WSTRESS
1623 IF (.not.destroy(ng, forces(ng)%ad_ustr, myfile, &
1624 & __line__, 'FORCES(ng)%ad_ustr')) RETURN
1625
1626 IF (.not.destroy(ng, forces(ng)%ad_vstr, myfile, &
1627 & __line__, 'FORCES(ng)%ad_vstr')) RETURN
1628# endif
1629
1630 IF (.not.destroy(ng, forces(ng)%ad_bustr, myfile, &
1631 & __line__, 'FORCES(ng)%ad_bustr')) RETURN
1632
1633 IF (.not.destroy(ng, forces(ng)%ad_bvstr, myfile, &
1634 & __line__, 'FORCES(ng)%ad_bvstr')) RETURN
1635
1636 IF (.not.destroy(ng, forces(ng)%ad_bustr_sol, myfile, &
1637 & __line__, 'FORCES(ng)%ad_bustr_sol')) RETURN
1638
1639 IF (.not.destroy(ng, forces(ng)%ad_bvstr_sol, myfile, &
1640 & __line__, 'FORCES(ng)%ad_bvstr_sol')) RETURN
1641
1642# ifdef SOLVE3D
1643 IF (.not.destroy(ng, forces(ng)%ad_stflx, myfile, &
1644 & __line__, 'FORCES(ng)%ad_stflx')) RETURN
1645
1646 IF (.not.destroy(ng, forces(ng)%ad_btflx, myfile, &
1647 & __line__, 'FORCES(ng)%ad_btflx')) RETURN
1648
1649# ifdef ADJUST_STFLUX
1650 IF (.not.destroy(ng, forces(ng)%ad_tflux, myfile, &
1651 & __line__, 'FORCES(ng)%ad_tflux')) RETURN
1652# endif
1653
1654# ifdef SHORTWAVE
1655 IF (.not.destroy(ng, forces(ng)%ad_srflx, myfile, &
1656 & __line__, 'FORCES(ng)%ad_srflx')) RETURN
1657# endif
1658
1659# ifdef BULK_FLUXES
1660 IF (.not.destroy(ng, forces(ng)%ad_lhflx, myfile, &
1661 & __line__, 'FORCES(ng)%ad_lhflx')) RETURN
1662
1663 IF (.not.destroy(ng, forces(ng)%ad_lrflx, myfile, &
1664 & __line__, 'FORCES(ng)%ad_lrflx')) RETURN
1665
1666 IF (.not.destroy(ng, forces(ng)%ad_shflx, myfile, &
1667 & __line__, 'FORCES(ng)%ad_shflx')) RETURN
1668
1669# ifdef EMINUSP
1670 IF (.not.destroy(ng, forces(ng)%ad_evap, myfile, &
1671 & __line__, 'FORCES(ng)%ad_evap')) RETURN
1672# endif
1673# endif
1674# endif
1675# endif
1676
1677# if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
1678!
1679! Working arrays to store adjoint impulse forcing, background error
1680! covariance, background-error standard deviations, or descent
1681! conjugate vectors (directions).
1682!
1683# if defined FOUR_DVAR || defined IMPULSE
1684# ifdef ADJUST_WSTRESS
1685 IF (.not.destroy(ng, forces(ng)%b_sustr, myfile, &
1686 & __line__, 'FORCES(ng)%b_sustr')) RETURN
1687
1688 IF (.not.destroy(ng, forces(ng)%b_svstr, myfile, &
1689 & __line__, 'FORCES(ng)%b_svstr')) RETURN
1690# endif
1691
1692# if defined ADJUST_STFLUX && defined SOLVE3D
1693 IF (.not.destroy(ng, forces(ng)%b_stflx, myfile, &
1694 & __line__, 'FORCES(ng)%b_stflx')) RETURN
1695# endif
1696# endif
1697
1698# ifdef FOUR_DVAR
1699# ifdef ADJUST_WSTRESS
1700 IF (.not.destroy(ng, forces(ng)%d_sustr, myfile, &
1701 & __line__, 'FORCES(ng)%d_sustr')) RETURN
1702
1703 IF (.not.destroy(ng, forces(ng)%d_svstr, myfile, &
1704 & __line__, 'FORCES(ng)%d_svstr')) RETURN
1705
1706 IF (.not.destroy(ng, forces(ng)%e_sustr, myfile, &
1707 & __line__, 'FORCES(ng)%e_sustr')) RETURN
1708
1709 IF (.not.destroy(ng, forces(ng)%e_svstr, myfile, &
1710 & __line__, 'FORCES(ng)%e_svstr')) RETURN
1711# endif
1712
1713# if defined ADJUST_STFLUX && defined SOLVE3D
1714 IF (.not.destroy(ng, forces(ng)%d_stflx, myfile, &
1715 & __line__, 'FORCES(ng)%d_stflx')) RETURN
1716
1717 IF (.not.destroy(ng, forces(ng)%e_stflx, myfile, &
1718 & __line__, 'FORCES(ng)%e_stflx')) RETURN
1719# endif
1720# endif
1721# endif
1722#endif
1723!
1724!-----------------------------------------------------------------------
1725! Deallocate derived-type FORCES structure.
1726!-----------------------------------------------------------------------
1727!
1728 IF (ng.eq.ngrids) THEN
1729 IF (allocated(forces)) deallocate ( forces )
1730 END IF
1731!
1732 RETURN
1733 END SUBROUTINE deallocate_forces
1734!
1735 SUBROUTINE initialize_forces (ng, tile, model)
1736!
1737!=======================================================================
1738! !
1739! This routine initialize all variables in the module using first !
1740! touch distribution policy. In shared-memory configuration, this !
1741! operation actually performs propagation of the "shared arrays" !
1742! across the cluster, unless another policy is specified to !
1743! override the default. !
1744! !
1745!=======================================================================
1746!
1747 USE mod_param
1748#ifdef BIOLOGY
1749 USE mod_biology
1750#endif
1751#if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
1752 USE mod_scalars
1753#endif
1754!
1755! Imported variable declarations.
1756!
1757 integer, intent(in) :: ng, tile, model
1758!
1759! Local variable declarations.
1760!
1761 integer :: imin, imax, jmin, jmax
1762 integer :: i, j, k
1763#ifdef SOLVE3D
1764 integer :: itrc
1765#endif
1766
1767 real(r8), parameter :: inival = 0.0_r8
1768
1769#include "set_bounds.h"
1770!
1771! Set array initialization range.
1772!
1773#ifdef DISTRIBUTE
1774 imin=bounds(ng)%LBi(tile)
1775 imax=bounds(ng)%UBi(tile)
1776 jmin=bounds(ng)%LBj(tile)
1777 jmax=bounds(ng)%UBj(tile)
1778#else
1779 IF (domain(ng)%Western_Edge(tile)) THEN
1780 imin=bounds(ng)%LBi(tile)
1781 ELSE
1782 imin=istr
1783 END IF
1784 IF (domain(ng)%Eastern_Edge(tile)) THEN
1785 imax=bounds(ng)%UBi(tile)
1786 ELSE
1787 imax=iend
1788 END IF
1789 IF (domain(ng)%Southern_Edge(tile)) THEN
1790 jmin=bounds(ng)%LBj(tile)
1791 ELSE
1792 jmin=jstr
1793 END IF
1794 IF (domain(ng)%Northern_Edge(tile)) THEN
1795 jmax=bounds(ng)%UBj(tile)
1796 ELSE
1797 jmax=jend
1798 END IF
1799#endif
1800!
1801!-----------------------------------------------------------------------
1802! Initialize module variables.
1803!-----------------------------------------------------------------------
1804!
1805! Nonlinear model state.
1806!
1807 IF ((model.eq.0).or.(model.eq.inlm)) THEN
1808 DO j=jmin,jmax
1809 DO i=imin,imax
1810#ifdef ADJUST_WSTRESS
1811 DO k=1,nfrec(ng)
1812 forces(ng) % ustr(i,j,k,1) = inival
1813 forces(ng) % ustr(i,j,k,2) = inival
1814 forces(ng) % vstr(i,j,k,1) = inival
1815 forces(ng) % vstr(i,j,k,2) = inival
1816 END DO
1817#endif
1818 forces(ng) % sustr(i,j) = inival
1819 forces(ng) % svstr(i,j) = inival
1820#if !defined ANA_SMFLUX && !defined BULK_FLUXES || \
1821 defined forward_fluxes
1822 forces(ng) % sustrG(i,j,1) = inival
1823 forces(ng) % sustrG(i,j,2) = inival
1824 forces(ng) % svstrG(i,j,1) = inival
1825 forces(ng) % svstrG(i,j,2) = inival
1826#endif
1827 forces(ng) % bustr(i,j) = inival
1828 forces(ng) % bvstr(i,j) = inival
1829#if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS || \
1830 defined spectral_light
1831 forces(ng) % Pair(i,j) = inival
1832# ifndef ANA_PAIR
1833 forces(ng) % PairG(i,j,1) = inival
1834 forces(ng) % PairG(i,j,2) = inival
1835# endif
1836#endif
1837#ifdef WAVES_DIR
1838 forces(ng) % Dwave(i,j) = inival
1839# ifndef ANA_WWAVE
1840 forces(ng) % DwaveG(i,j,1) = inival
1841 forces(ng) % DwaveG(i,j,2) = inival
1842# endif
1843#endif
1844#ifdef WAVES_DIRP
1845 forces(ng) % Dwavep(i,j) = inival
1846# ifndef ANA_WWAVE
1847 forces(ng) % DwavepG(i,j,1) = inival
1848 forces(ng) % DwavepG(i,j,2) = inival
1849# endif
1850#endif
1851#ifdef WAVES_HEIGHT
1852 forces(ng) % Hwave(i,j) = inival
1853# ifndef ANA_WWAVE
1854 forces(ng) % HwaveG(i,j,1) = inival
1855 forces(ng) % HwaveG(i,j,2) = inival
1856# endif
1857#endif
1858#ifdef WAVES_LENGTH
1859 forces(ng) % Lwave(i,j) = inival
1860# ifndef ANA_WWAVE
1861 forces(ng) % LwaveG(i,j,1) = inival
1862 forces(ng) % LwaveG(i,j,2) = inival
1863# endif
1864#endif
1865#ifdef WAVES_LENGTHP
1866 forces(ng) % Lwavep(i,j) = inival
1867# ifndef ANA_WWAVE
1868 forces(ng) % LwavepG(i,j,1) = inival
1869 forces(ng) % LwavepG(i,j,2) = inival
1870# endif
1871#endif
1872#ifdef WAVES_TOP_PERIOD
1873 forces(ng) % Pwave_top(i,j) = inival
1874# ifndef ANA_WWAVE
1875 forces(ng) % Pwave_topG(i,j,1) = inival
1876 forces(ng) % Pwave_topG(i,j,2) = inival
1877# endif
1878#endif
1879#ifdef WAVES_BOT_PERIOD
1880 forces(ng) % Pwave_bot(i,j) = inival
1881# ifndef ANA_WWAVE
1882 forces(ng) % Pwave_botG(i,j,1) = inival
1883 forces(ng) % Pwave_botG(i,j,2) = inival
1884# endif
1885#endif
1886#if defined BBL_MODEL || defined WAV_COUPLING
1887 forces(ng) % Uwave_rms(i,j) = inival
1888# ifndef ANA_WWAVE
1889 forces(ng) % Uwave_rmsG(i,j,1) = inival
1890 forces(ng) % Uwave_rmsG(i,j,2) = inival
1891# endif
1892#endif
1893
1894#if defined TKE_WAVEDISS || defined WAV_COUPLING || \
1895 defined wdiss_thorguza || defined wdiss_churthor || \
1896 defined waves_diss || defined wdiss_inwave
1897 forces(ng) % Dissip_break(i,j) = inival
1898 forces(ng) % Dissip_wcap(i,j) = inival
1899# ifndef ANA_WWAVE
1900 forces(ng) % Dissip_breakG(i,j,1) = inival
1901 forces(ng) % Dissip_breakG(i,j,2) = inival
1902 forces(ng) % Dissip_wcapG(i,j,1) = inival
1903 forces(ng) % Dissip_wcapG(i,j,2) = inival
1904# endif
1905#endif
1906#if defined WAV_COUPLING || (defined WEC_VF && defined BOTTOM_STREAMING)
1907 forces(ng) % Dissip_fric(i,j) = inival
1908# ifndef ANA_WWAVE
1909 forces(ng) % Dissip_fricG(i,j,1) = inival
1910 forces(ng) % Dissip_fricG(i,j,2) = inival
1911# endif
1912#endif
1913
1914#if defined WEC_ROLLER
1915 forces(ng) % Dissip_roller(i,j) = inival
1916# ifndef ANA_WWAVE
1917 forces(ng) % Dissip_rollerG(i,j,1) = inival
1918 forces(ng) % Dissip_rollerG(i,j,2) = inival
1919# endif
1920#endif
1921
1922#if defined ROLLER_SVENDSEN
1923 forces(ng) % Wave_break(i,j) = inival
1924# ifndef ANA_WWAVE
1925 forces(ng) % Wave_breakG(i,j,1) = inival
1926 forces(ng) % Wave_breakG(i,j,2) = inival
1927# endif
1928#endif
1929
1930#if defined WAVES_DSPR
1931 forces(ng) % Wave_ds(i,j) = inival
1932 forces(ng) % Wave_qp(i,j) = inival
1933# ifndef ANA_WWAVE
1934 forces(ng) % Wave_dsG(i,j,1) = inival
1935 forces(ng) % Wave_dsG(i,j,2) = inival
1936 forces(ng) % Wave_qpG(i,j,1) = inival
1937 forces(ng) % Wave_qpG(i,j,2) = inival
1938# endif
1939#endif
1940
1941#if defined WEC_ROLLER
1942 forces(ng) % rollA(i,j) = inival
1943#endif
1944
1945#ifdef SPECTRUM_STOKES
1946 DO k=1,mscs
1947 forces(ng) % spec_wn(i,j,k) = inival
1948 forces(ng) % spec_us(i,j,k) = inival
1949 forces(ng) % spec_vs(i,j,k) = inival
1950 END DO
1951#endif
1952
1953#ifdef SOLVE3D
1954# ifdef SHORTWAVE
1955 forces(ng) % srflx(i,j) = inival
1956# ifndef ANA_SRFLUX
1957 forces(ng) % srflxG(i,j,1) = inival
1958 forces(ng) % srflxG(i,j,2) = inival
1959# endif
1960# ifdef ICE_ALBEDO
1961 forces(ng) % albedo(i,j) = inival
1962# ifdef ICE_MODEL
1963 forces(ng) % albedo_ice(i,j) = inival
1964# endif
1965# endif
1966# endif
1967# if defined RED_TIDE && defined DAILY_SHORTWAVE
1968 forces(ng) % srflx_avg(i,j) = inival
1969 forces(ng) % srflxG_avg(i,j,1) = inival
1970 forces(ng) % srflxG_avg(i,j,2) = inival
1971# endif
1972# ifdef CLOUDS
1973 forces(ng) % cloud(i,j) = inival
1974# ifndef ANA_CLOUD
1975 forces(ng) % cloudG(i,j,1) = inival
1976 forces(ng) % cloudG(i,j,2) = inival
1977# endif
1978# endif
1979# if defined BULK_FLUXES || defined FRC_COUPLING
1980 forces(ng) % lhflx(i,j) = inival
1981 forces(ng) % lrflx(i,j) = inival
1982 forces(ng) % shflx(i,j) = inival
1983# endif
1984# if defined BULK_FLUXES || defined ECOSIM || \
1985 (defined shortwave && defined ana_srflux) || \
1986 defined spectral_light
1987 forces(ng) % Hair(i,j) = inival
1988 forces(ng) % Tair(i,j) = inival
1989# endif
1990# if defined BULK_FLUXES || defined ECOSIM || \
1991 defined spectral_light
1992 forces(ng) % Uwind(i,j) = inival
1993 forces(ng) % Vwind(i,j) = inival
1994# endif
1995# ifdef BULK_FLUXES
1996 forces(ng) % rain(i,j) = inival
1997# ifdef EMINUSP
1998 forces(ng) % evap(i,j) = inival
1999# endif
2000# endif
2001# if defined ICE_MODEL && \
2002 defined ice_bulk_fluxes && defined bulk_fluxes
2003 forces(ng) % Qnet_ai(i,j) = inival
2004 forces(ng) % Qnet_ao(i,j) = inival
2005 forces(ng) % snow(i,j) = inival
2006 forces(ng) % srflx_ice(i,j) = inival
2007 forces(ng) % sustr_ai(i,j) = inival
2008 forces(ng) % svstr_ai(i,j) = inival
2009 forces(ng) % sustr_ao(i,j) = inival
2010 forces(ng) % svstr_ao(i,j) = inival
2011# endif
2012# if !defined LONGWAVE && defined BULK_FLUXES
2013 forces(ng) % lrflxG(i,j,1) = inival
2014 forces(ng) % lrflxG(i,j,2) = inival
2015# endif
2016# if defined BULK_FLUXES || defined ECOSIM || \
2017 (defined shortwave && defined ana_srflux) || \
2018 defined spectral_light
2019# ifndef ANA_HUMIDITY
2020 forces(ng) % HairG(i,j,1) = inival
2021 forces(ng) % HairG(i,j,2) = inival
2022# endif
2023# endif
2024# if defined BULK_FLUXES || defined ECOSIM
2025# ifndef ANA_TAIR
2026 forces(ng) % TairG(i,j,1) = inival
2027 forces(ng) % TairG(i,j,2) = inival
2028# endif
2029# ifndef ANA_WINDS
2030 forces(ng) % UwindG(i,j,1) = inival
2031 forces(ng) % UwindG(i,j,2) = inival
2032 forces(ng) % VwindG(i,j,1) = inival
2033 forces(ng) % VwindG(i,j,2) = inival
2034# endif
2035# endif
2036# if !defined ANA_RAIN && defined BULK_FLUXES
2037 forces(ng) % rainG(i,j,1) = inival
2038 forces(ng) % rainG(i,j,2) = inival
2039# endif
2040# ifdef QCORRECTION
2041 forces(ng) % dqdt(i,j) = inival
2042 forces(ng) % sst(i,j) = inival
2043# ifndef ANA_SST
2044 forces(ng) % dqdtG(i,j,1) = inival
2045 forces(ng) % dqdtG(i,j,2) = inival
2046 forces(ng) % sstG(i,j,1) = inival
2047 forces(ng) % sstG(i,j,2) = inival
2048# endif
2049# endif
2050# if defined SALINITY && (defined SCORRECTION || defined SRELAXATION)
2051 forces(ng) % sss(i,j) = inival
2052# ifndef ANA_SSS
2053 forces(ng) % sssG(i,j,1) = inival
2054 forces(ng) % sssG(i,j,2) = inival
2055# endif
2056# endif
2057 DO itrc=1,nt(ng)
2058 forces(ng) % stflux(i,j,itrc) = inival
2059# if !defined ANA_STFLUX || !defined ANA_SSFLUX || \
2060 !defined ANA_SPFLUX
2061 forces(ng) % stfluxG(i,j,1,itrc) = inival
2062 forces(ng) % stfluxG(i,j,2,itrc) = inival
2063# endif
2064 forces(ng) % stflx(i,j,itrc) = inival
2065# ifdef ADJUST_STFLUX
2066 DO k=1,nfrec(ng)
2067 forces(ng) % tflux(i,j,k,1,itrc) = inival
2068 forces(ng) % tflux(i,j,k,2,itrc) = inival
2069 END DO
2070# endif
2071
2072 forces(ng) % btflux(i,j,itrc) = inival
2073# if !defined ANA_BTFLUX || !defined ANA_BSFLUX || \
2074 !defined ANA_BPFLUX
2075 forces(ng) % btfluxG(i,j,1,itrc) = inival
2076 forces(ng) % btfluxG(i,j,2,itrc) = inival
2077# endif
2078 forces(ng) % btflx(i,j,itrc) = inival
2079 END DO
2080# if defined ECOSIM || defined SPECTRAL_LIGHT
2081 DO itrc=1,nbands
2082 forces(ng) % SpecIr(i,j,itrc) = inival
2083 forces(ng) % avcos(i,j,itrc) = inival
2084 END DO
2085# endif
2086#endif
2087 END DO
2088 END DO
2089 END IF
2090
2091#if defined TANGENT || defined TL_IOMS
2092!
2093! Tangent linear model state.
2094!
2095 IF ((model.eq.0).or.(model.eq.itlm).or.(model.eq.irpm)) THEN
2096 DO j=jmin,jmax
2097 DO i=imin,imax
2098# ifdef ADJUST_WSTRESS
2099 DO k=1,nfrec(ng)
2100 forces(ng) % tl_ustr(i,j,k,1) = inival
2101 forces(ng) % tl_ustr(i,j,k,2) = inival
2102 forces(ng) % tl_vstr(i,j,k,1) = inival
2103 forces(ng) % tl_vstr(i,j,k,2) = inival
2104 END DO
2105# endif
2106 forces(ng) % tl_sustr(i,j) = inival
2107 forces(ng) % tl_svstr(i,j) = inival
2108 forces(ng) % tl_bustr(i,j) = inival
2109 forces(ng) % tl_bvstr(i,j) = inival
2110 END DO
2111# ifdef SOLVE3D
2112# ifdef SHORTWAVE
2113 DO i=imin,imax
2114 forces(ng) % tl_srflx(i,j) = inival
2115 END DO
2116# endif
2117# ifdef BULK_FLUXES
2118 DO i=imin,imax
2119 forces(ng) % tl_lhflx(i,j) = inival
2120 forces(ng) % tl_lrflx(i,j) = inival
2121 forces(ng) % tl_shflx(i,j) = inival
2122# ifdef EMINUSP
2123 forces(ng) % tl_evap(i,j) = inival
2124# endif
2125 END DO
2126# endif
2127 DO itrc=1,nt(ng)
2128 DO i=imin,imax
2129# ifdef ADJUST_STFLUX
2130 DO k=1,nfrec(ng)
2131 forces(ng) % tl_tflux(i,j,k,1,itrc) = inival
2132 forces(ng) % tl_tflux(i,j,k,2,itrc) = inival
2133 END DO
2134# endif
2135 forces(ng) % tl_stflux(i,j,itrc) = inival
2136 forces(ng) % tl_stflx (i,j,itrc) = inival
2137 forces(ng) % tl_btflux(i,j,itrc) = inival
2138 forces(ng) % tl_btflx (i,j,itrc) = inival
2139 END DO
2140 END DO
2141# endif
2142 END DO
2143 END IF
2144#endif
2145
2146#ifdef ADJOINT
2147!
2148! Adjoint model state.
2149!
2150 IF ((model.eq.0).or.(model.eq.iadm)) THEN
2151 DO j=jmin,jmax
2152 DO i=imin,imax
2153# ifdef ADJUST_WSTRESS
2154 DO k=1,nfrec(ng)
2155 forces(ng) % ad_ustr(i,j,k,1) = inival
2156 forces(ng) % ad_ustr(i,j,k,2) = inival
2157 forces(ng) % ad_vstr(i,j,k,1) = inival
2158 forces(ng) % ad_vstr(i,j,k,2) = inival
2159 END DO
2160# endif
2161 forces(ng) % ad_sustr(i,j) = inival
2162 forces(ng) % ad_svstr(i,j) = inival
2163 forces(ng) % ad_bustr(i,j) = inival
2164 forces(ng) % ad_bvstr(i,j) = inival
2165 forces(ng) % ad_bustr_sol(i,j) = inival
2166 forces(ng) % ad_bvstr_sol(i,j) = inival
2167 END DO
2168# ifdef SOLVE3D
2169# ifdef SHORTWAVE
2170 DO i=imin,imax
2171 forces(ng) % ad_srflx(i,j) = inival
2172 END DO
2173# endif
2174# ifdef BULK_FLUXES
2175 DO i=imin,imax
2176 forces(ng) % ad_lhflx(i,j) = inival
2177 forces(ng) % ad_lrflx(i,j) = inival
2178 forces(ng) % ad_shflx(i,j) = inival
2179# ifdef EMINUSP
2180 forces(ng) % ad_evap(i,j) = inival
2181# endif
2182 END DO
2183# endif
2184 DO itrc=1,nt(ng)
2185 DO i=imin,imax
2186# ifdef ADJUST_STFLUX
2187 DO k=1,nfrec(ng)
2188 forces(ng) % ad_tflux(i,j,k,1,itrc) = inival
2189 forces(ng) % ad_tflux(i,j,k,2,itrc) = inival
2190 END DO
2191# endif
2192 forces(ng) % ad_stflx(i,j,itrc) = inival
2193 forces(ng) % ad_btflx(i,j,itrc) = inival
2194 END DO
2195 END DO
2196# endif
2197 END DO
2198 END IF
2199#endif
2200
2201#if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
2202!
2203! Working arrays to store adjoint impulse forcing, background error
2204! covariance, background-error standard deviations, or descent
2205! conjugate vectors (directions).
2206!
2207# if defined FOUR_DVAR || defined IMPULSE
2208 IF (model.eq.0) THEN
2209# ifdef ADJUST_WSTRESS
2210 DO j=jmin,jmax
2211 DO i=imin,imax
2212 forces(ng) % b_sustr(i,j) = inival
2213 forces(ng) % b_svstr(i,j) = inival
2214# ifdef FOUR_DVAR
2215 forces(ng) % e_sustr(i,j) = inival
2216 forces(ng) % e_svstr(i,j) = inival
2217 DO k=1,nfrec(ng)
2218 forces(ng) % d_sustr(i,j,k) = inival
2219 forces(ng) % d_svstr(i,j,k) = inival
2220 END DO
2221# endif
2222 END DO
2223 END DO
2224# endif
2225# if defined ADJUST_STFLUX && defined SOLVE3D
2226 DO itrc=1,nt(ng)
2227 DO j=jmin,jmax
2228 DO i=imin,imax
2229 forces(ng) % b_stflx(i,j,itrc) = inival
2230# ifdef FOUR_DVAR
2231 forces(ng) % e_stflx(i,j,itrc) = inival
2232 DO k=1,nfrec(ng)
2233 forces(ng) % d_stflx(i,j,k,itrc) = inival
2234 END DO
2235# endif
2236 END DO
2237 END DO
2238 END DO
2239# endif
2240 END IF
2241# endif
2242#endif
2243!
2244 RETURN
2245 END SUBROUTINE initialize_forces
2246
2247 END MODULE mod_forces
integer, parameter nbands
Definition ecosim_mod.h:201
subroutine, public initialize_forces(ng, tile, model)
subroutine, public deallocate_forces(ng)
subroutine, public allocate_forces(ng, lbi, ubi, lbj, ubj)
Definition mod_forces.F:559
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
integer, parameter r8
Definition mod_kinds.F:28
integer, parameter inlm
Definition mod_param.F:662
integer, parameter irpm
Definition mod_param.F:664
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, dimension(:), allocatable nfrec