WAVEWATCH III  beta 0.0.1
w3adatmd.F90
Go to the documentation of this file.
1 
8 
9 #include "w3macros.h"
10 !/ ------------------------------------------------------------------- /
26 MODULE w3adatmd
27  !/
28  !/ +-----------------------------------+
29  !/ | WAVEWATCH III NOAA/NCEP |
30  !/ | H. L. Tolman |
31  !/ | FORTRAN 90 |
32  !/ | Last update : 22-Mar-2021 |
33  !/ +-----------------------------------+
34  !/
35  !/ 28-Dec-2004 : Origination. ( version 3.06 )
36  !/ 04-May-2005 : Adding MPI_COMM_WAVE. ( version 3.07 )
37  !/ 20-Jul-2005 : Adding output fields. ( version 3.07 )
38  !/ 09-Nov-2005 : Removing soft boundary option. ( version 3.08 )
39  !/ 13-Jun-2006 : Splitting STORE in G/SSTORE. ( version 3.09 )
40  !/ 04-Oct-2006 : Add filter to array pointers. ( version 3.10 )
41  !/ 28_Mar-2007 : Add partitioned data arrays. ( version 3.11 )
42  !/ Add aditional undefined arrays.
43  !/ 22-Feb-2008 ; Modify MAPTH2 declaration. ( version 3.13 )
44  !/ 29-May-2009 : Preparing distribution version. ( version 3.14 )
45  !/ 29-Oct-2010 : Adding unstructured grid data. ( version 3.14 )
46  !/ (A. Roland and F. Ardhuin)
47  !/ 31-Oct-2010 : Adding output parameters ( version 3.14 )
48  !/ 12-Dec-2012 : Adding SMC grid. JG_Li ( version 4.08 )
49  !/ 26-Dec-2012 : Memory reduction for outputs. ( version 4.11 )
50  !/ Add W3XETA.
51  !/ 28-Jun-2013 : Bug fix initialization P2SMS. ( version 4.11 )
52  !/ 11-Nov-2013 : SMC and rotated grid incorporated in the main
53  !/ trunk ( version 4.13 )
54  !/ 14-Nov-2013 : Move orphaned arrays as scalar to W3SRCE.
55  !/ Here update of documentation only.
56  !/ (Z0S, CDS, EMN, FMN, WNM, AMX) ( version 4.13 )
57  !/ 30-Apr-2014 : Memory reduction for group3. ( version 5.00 )
58  !/ 10-Dec-2014 : Add checks for allocate status ( version 5.04 )
59  !/ 01-May-2017 : Adds directional MSS parameters ( version 6.02 )
60  !/ 30-Jul-2017 : Adds TWS parameter ( version 6.02 )
61  !/ 05-Jun-2018 : Adds PDLIB and MEMCHECK ( version 6.04 )
62  !/ 21-Aug-2018 : Add WBT parameter ( version 6.06 )
63  !/ 22-Mar-2021 : Adds TAUA, WNMEAN, TAUOC parameters ( version 7.13 )
64  !/ 06-May-2021 : SMC shares variables with PR2/3. ( version 7.13 )
65  !
66  !/
67  !/ Copyright 2009-2013 National Weather Service (NWS),
68  !/ National Oceanic and Atmospheric Administration. All rights
69  !/ reserved. WAVEWATCH III is a trademark of the NWS.
70  !/ No unauthorized use without permission.
71  !/
72  ! 1. Purpose :
73  !
74  ! Define data structures to set up wave model auxiliary data for
75  ! several models simultaneously.
76  !
77  ! 2. Variables and types :
78  !
79  ! Name Type Scope Description
80  ! ----------------------------------------------------------------
81  ! NADATA Int. Public Number of models in array dim.
82  ! IADATA Int. Public Selected model for output, init. at -1.
83  ! MPIBUF I.P. Public Number of buffer arrays for 'hidden'
84  ! MPI communications (no hiding for
85  ! MPIBUF = 1).
86  ! WADAT TYPE Public Basic data structure.
87  ! WADATS WADAT Public Array of data structures.
88  ! ----------------------------------------------------------------
89  !
90  ! All elements of WADAT are aliased to pointers with the same
91  ! name. These pointers are defined as :
92  !
93  ! Name Type Scope Description
94  ! ----------------------------------------------------------------
95  ! Internal model definition:
96  !
97  ! CG R.A. Public Group velocities for all wave model
98  ! sea points and frequencies.
99  ! WN R.A. Public Idem, wavenumbers.
100  !
101  ! Aux. arrays for model input:
102  !
103  ! CA0-I R.A. Public Absolute current velocity (initial
104  ! and inc.) in W3UCUR.
105  ! CD0-I R.A. Public Current direction (initial and
106  ! increment) in W3UCUR.
107  ! UA0-I R.A. Public Absolute wind speeds (initial and
108  ! incr.) in W3UWND (m/s)
109  ! UD0-I R.A. Public Wind direction (initial and incr.)
110  ! in W3UWND (rad)
111  ! AS0-I R.A. Public Stability par. (initial and incr.)
112  ! in W3UWND (degr)
113  ! MA0-I R.A. Public Absolute atmospheric momentum (initial
114  ! and inc.) in W3UTAU.
115  ! RA0-I R.A. Public Absolute air density (initial and inc.)
116  ! in W3URHO.
117  ! MD0-I R.A. Public Atmospheric momentum direction (initial and
118  ! increment) in W3UTAU.
119  ! ATRNX/Y R.A. Public Actual transparency info.
120  !
121  ! Fields of mean wave parameters:
122  !
123  ! DW R.A. Public Water depths.
124  ! UA R.A. Public Absolute wind speeds.
125  ! UD R.A. Public Absolute wind direction.
126  ! U10 R.A. Public Wind speed used.
127  ! U10D R.A. Public Wind direction used.
128  ! AS R.A. Public Stability parameter.
129  ! CX/Y R.A. Public Current components.
130  ! TAUA R.A. Public Absolute atmospheric momentum.
131  ! TAUADIR R.A. Public Absolute atmospheric momentum direction.
132  !
133  ! HS R.A. Public Wave Height.
134  ! WLM R.A. Public Mean wave length.
135  ! T02 R.A. Public Mean wave period (m0,2).
136  ! T0M1 R.A. Public Mean wave period (m0,-1).
137  ! T01 R.A. Public Mean wave period (m0,1).
138  ! FP0 R.A. Public Peak frequency.
139  ! THM R.A. Public Mean wave direction.
140  ! THS R.A. Public Mean directional spread.
141  ! THP0 R.A. Public Peak direction.
142  ! HSIG R.A. Public Height of infragravity waves
143  ! STMAXE R.A. Public Expected maximum surface elevation (crest)
144  ! STMAXD R.A. Public STD of maximum surface elevation
145  ! HMAXE R.A. Public Expected maximum wave height (from covariance)
146  ! HMAXD R.A. Public Std of HMAXE
147  ! HCMAXE R.A. Public Expected maximum wave height (from crest)
148  ! HCMAXD R.A. Public STD of HCMAXE
149  ! WBT R.A. Public Dominant wave breaking probability
150  ! (b_T in Babanin et al. (2001, JGR))
151  ! WNMEAN R.A. Public Mean wave number
152  !
153  ! CHARN R.A. Public Charnock parameter for air-sea friction.
154  ! TWS R.A. Public Wind sea period (used for flux parameterizations)
155  ! CGE R.A. Public Energy flux.
156  ! PHIAW R.A. Public Wind to wave energy flux.
157  ! TAUWIX/Y R.A. Public Wind to wave energy flux.
158  ! TAUWNX/Y R.A. Public Wind to wave energy flux.
159  ! WHITECAP R.A. Public 1 : Whitecap coverage
160  ! 2 : Whitecap thickness
161  ! 3 : Mean breaking height
162  ! 4 : Mean breaking height
163  !
164  ! Sxx R.A. Public Radiation stresses.
165  ! TAUOX/Y R.A. Public Wave-ocean momentum flux.
166  ! BHD R.A. Public Wave-induced pressure (J term, Smith JPO 2006)
167  ! PHIOC R.A. Public Waves to ocean energy flux.
168  ! TUSX/Y R.A. Public Volume transport associated to Stokes drift.
169  ! USSX/Y R.A. Public Surface Stokes drift.
170  ! TAUOCX/Y R.A. Public Total ocean momentum flux
171  ! TAUICE R.A. Public Wave-ice momentum flux.
172  ! PHICE R.A. Public Waves to ice energy flux.
173  !
174  ! US3D R.A. Public 3D Stokes drift.
175  ! USSP R.A. Public Partitioned Surface Stokes drift
176  !
177  ! ABA R.A. Public Near-bottom rms wave ex. amplitude.
178  ! ABD R.A. Public Corresponding direction.
179  ! UBA R.A. Public Near-bottom rms wave velocity.
180  ! UBD R.A. Public Corresponding direction.
181  ! BEDFORMS R.A. Public Bed for parameters
182  ! PHIBBL R.A. Public Energy loss in WBBL.
183  ! TAUBBL R.A. Public Momentum loss in WBBL.
184  !
185  ! MSSX/Y R.A. Public Surface mean square slopes in X and Y direction.
186  ! MSCX/Y R.A. Public Phillips constant.
187  ! MSSD R.A. Public Direction of MSSX
188  ! MSCD R.A. Public Direction of MSCX
189  ! QP R.A. Public Goda peakedness parameter.
190  ! QKK R.A. Public Spectral bandwidth (De Carlo et al. 2023)
191  ! SKEW R.A. Public skewness lambda_3,0,0 (Srokosz 1986)
192  !
193  ! DTDYN R.A. Public Mean dynamic time step (raw).
194  ! FCUT R.A. Public Cut-off frequency for tail.
195  ! CFLXYMAX R.A. Public Max. CFL number for spatial advection.
196  ! CFLTHMAX R.A. Public Max. CFL number for refraction.
197  ! CFLKMAX R.A. Public Max. CFL number for wavenumber shift.
198  !
199  ! Orphans, commented out here, now automatic arrays in W3WAVE, ....
200  !
201  ! DRAT R.A. Public Density ration air/water. Was
202  ! placeholder only. Now scalar in W3SRCE,
203  ! TAUWX/Y R.A. Public Stresses.
204  !
205  ! Derivatives in space ....
206  !
207  ! DDDx R.A. Public Spatial derivatives of the depth.
208  ! DCxDx R.A. Public Spatial dirivatives of the current.
209  !
210  ! Mean parameters from partitiones spectra, 2D array with el.
211  ! 0 holding wind sea data, and 1:NOSWLL holding swell fields.
212  ! Last two arrays are regular single-entry arrays.
213  !
214  ! PHS R.A. Public Wave height of partition.
215  ! PTP R.A. Public Peak period of partition.
216  ! PLP R.A. Public Peak wave leingth of partition.
217  ! PDIR R.A. Public Mean direction of partition.
218  ! PSI R.A. Public Mean spread of partition.
219  ! PWS R.A. Public Wind sea fraction of partition.
220  !
221  ! PWST R.A. Public Total wind sea fraction.
222  ! PNR R.A. Public Number of partitions found.
223  !
224  ! PTHP0 R.A. Public Peak wave direction of partition.
225  ! PQP R.A. Public Goda peakdedness parameter of partition.
226  ! PPE R.A. Public JONSWAP peak enhancement factor of partition.
227  ! PGW R.A. Public Gaussian frequency width of partition.
228  ! PSW R.A. Public Spectral width of partition.
229  ! PTM1 R.A. Public Mean wave period (m-1,0) of partition.
230  ! PT1 R.A. Public Mean wave period (m0,1) of partition.
231  ! PT2 R.A. Public Mean wave period (m0,2) of partition.
232  ! PEP R.A. Public Peak spectral density of partition.
233  !
234  ! Empty dummy fields (NOEXTR)
235  !
236  ! USERO R.A. Public Empty output arrays than can be
237  ! used by users as a simple means to
238  ! add output.
239  !
240  ! Map data for propagation schemes (1Up).
241  !
242  ! IS0/2 I.A. Public Spectral propagation maps.
243  ! FACVX/Y R.A. Public Spatial propagation factor map.
244  !
245  ! Map data for propagation schemes (UQ).
246  !
247  ! NMXn Int. Public Counters for MAPX2, see W3MAP3.
248  ! NMYn Int. Public
249  ! NMXY Int. Public Dimension of MAPXY.
250  ! NACTn Int. Public Dimension of MAPAXY.
251  ! NCENT Int. Public Dimension of MAPAXY.
252  ! MAPX2 I.A. Public Map for prop. in 'x' (longitude) dir.
253  ! MAPY2 I.A. Public Idem in y' (latitude) direction.
254  ! MAPXY I.A. Public
255  ! MAPAXY I.A. Public List of active points used in W3QCK1.
256  ! MAPCXY I.A. Public List of central points used in avg.
257  ! MAPTH2 I.A. Public Like MAPX2 for refraction (rotated
258  ! and shifted, see W3KTP3). Like MAPAXY.
259  ! MAPWN2 I.A. Public Like MAPX2 for wavenumber shift.
260  ! MAPTRN L.A. Public Map to block out GSE mitigation in
261  ! proper grid points.
262  !
263  ! Nonlinear interactions ( !/NL1 ) :
264  !
265  ! NFR Int. Public Nuber of frequencies ( NFR = NK )
266  ! NFRHGH Int. Public Auxiliary frequency counter.
267  ! NFRCHG Int. Public Id.
268  ! NSPECX-Y Int. Public Auxiliary spectral counter.
269  ! IPnn I.A. Public Spectral address for Snl.
270  ! IMnn I.A. Public Id.
271  ! ICnn I.A. Public Id.
272  ! DALn Real Public Lambda dependend weight factors.
273  ! AWGn Real Public Interpolation weights for Snl.
274  ! SWGn Real Public Interpolation weights for diag. term.
275  ! AF11 R.A. Public Scaling array (f**11)
276  ! NLINIT Log. Public Flag for initialization.
277  !
278  ! MPP / MPI variables :
279  !
280  ! IAPPRO I.A. Public Processor numbers for propagation calc.
281  ! for each spectral component.
282  ! MPI_COMM_WAVE
283  ! Int. Public Communicator used in the wave model.
284  ! MPI_COMM_WCMP
285  ! Int. Public Idem, computational proc. only.
286  ! WW3_FIELD_VEC, WW3_SPEC_VEC
287  ! Int. Public MPI derived vecor types.
288  ! NRQSG1 Int. Public Number of handles in IRQSG1.
289  ! NRQSG2 Int. Public Number of handles in IRQSG2.
290  ! IBFLOC Int. Public Present active buffer number.
291  ! ISPLOC Int. Public Corresponding local spectral bin number
292  ! (1,NSPLOC,1).
293  ! NSPLOC Int. Public Total number of spectral bins for which
294  ! prop. is performed on present CPU.
295  ! BSTAT I.A. Public Status of buffer (size MPIBUF):
296  ! 0: Inactive.
297  ! 1: A --> STORE (active or finished).
298  ! 2: STORE --> A (active or finished).
299  ! BISPL I.A. Public Local spectral bin number for buffer
300  ! (size MPIBUF).
301  ! IRQSG1 I.A. Public MPI request handles for scatters and
302  ! gathers to A() (persistent).
303  ! IRQSG2 I.A. Public MPI request handles for gathers and
304  ! scatters to STORE (persistent).
305  ! G/SSTORE R.A. Public Communication buffer (NSEA,MPIBUF).
306  ! SPPNT R.A. Public Point output buffer.
307  !
308  ! Other:
309  !
310  ! ITIME Int. Public Discrete time step counter.
311  ! IPASS Int. Public Pass counter for log file.
312  ! IDLAST Int. Public Last day ID for log file.
313  ! NSEALM Int. Public Maximum number of local sea points.
314  ! ALPHA R.A. Public Phillips' alpha.
315  ! FLCOLD Log. Public Flag for 'cold start' of model.
316  ! FLIWND Log. Public Flag for initialization of model
317  ! based on wind.
318  ! AINIT(2) Log. Public Flag for array initialization.
319  ! FL_ALL Log. Public Flag for all/partial initialization.
320  ! ----------------------------------------------------------------
321  !
322  ! 3. Subroutines and functions :
323  !
324  ! Name Type Scope Description
325  ! ----------------------------------------------------------------
326  ! W3NAUX Subr. Public Set number of grids/models.
327  ! W3DIMA Subr. Public Set dimensions of arrays.
328  ! W3DMNL Subr. Public Set dimensions of arrays. ( !/NL1 )
329  ! W3SETA Subr. Public Point to selected grid / model.
330  ! W3XETA Subr. Public Like W3SETA for expanded output arrays.
331  ! ----------------------------------------------------------------
332  !
333  ! 4. Subroutines and functions used :
334  !
335  ! Name Type Module Description
336  ! ----------------------------------------------------------------
337  ! W3SETG Subr. W3GDATMD Point to proper model grid.
338  ! STRACE Subr. W3SERVMD Subroutine tracing.
339  ! EXTCDE Subr. W3SERVMD Abort program with exit code.
340  ! ----------------------------------------------------------------
341  !
342  ! 5. Remarks :
343  !
344  ! - The number of grids is taken from W3GDATMD, and needs to be
345  ! set first with W3DIMG.
346  !
347  ! 6. Switches :
348  !
349  ! !/SHRD, !/DIST, !/MPI
350  ! Shared / distributed memory model
351  !
352  ! !/PRn Propagation scheme selection.
353  !
354  ! !/S Enable subroutine tracing.
355  ! !/T Enable test output
356  !
357  ! 7. Source code :
358  !
359  !/ ------------------------------------------------------------------- /
360 
361  use w3servmd, only : print_memcheck
362 
363  ! module default
364  implicit none
365 
366  PUBLIC
367  !/
368  !/ Module private variable for checking error returns
369  !/
370  INTEGER, PRIVATE :: ISTAT
371  !/
372  !/ Conventional declarations
373  !/
374  INTEGER :: nadata = -1, iadata = -1
375 #ifdef W3_MPI
376  INTEGER, PARAMETER :: mpibuf = 6
377 #endif
378  !/
379  !/ Data structure WADAT
380  !/
381  TYPE wadat
382  !
383  ! The grid
384  !
385  REAL, POINTER :: cg(:,:), wn(:,:)
386 #ifdef W3_IC3
387  REAL, POINTER :: ic3wn_r(:,:), ic3wn_i(:,:), ic3cg(:,:)
388 #endif
389  !
390  ! Arrays for processing model input
391  !
392  REAL, POINTER :: ca0(:), cai(:), cd0(:), cdi(:), &
393  ua0(:), uai(:), ud0(:), udi(:), &
394  ma0(:), mai(:), ra0(:), rai(:), &
395  md0(:), mdi(:), as0(:), asi(:), &
396  atrnx(:,:), atrny(:,:)
397  !
398  ! Output fields group 1)
399  !
400  REAL, POINTER :: dw(:), ua(:), ud(:), u10(:), u10d(:),&
401  as(:), cx(:), cy(:), taua(:), tauadir(:)
402  !
403  ! Output fields group 2)
404  !
405  REAL, POINTER :: hs(:), wlm(:), t02(:), t0m1(:), &
406  t01 (:), fp0(:), thm(:), &
407  ths(:), thp0(:), &
408  hsig(:), stmaxe(:), stmaxd(:), &
409  hmaxe(:), hcmaxe(:), hmaxd(:), &
410  hcmaxd(:), qp(:), wbt(:), wnmean(:)
411  REAL, POINTER :: xhs(:), xwlm(:), xt02(:), xt0m1(:), &
412  xt01 (:), xfp0(:), xthm(:), &
413  xths(:), xthp0(:), &
414  xhsig(:), xstmaxe(:), xstmaxd(:), &
415  xhmaxe(:), xhcmaxe(:), xhmaxd(:), &
416  xhcmaxd(:), xqp(:), xwbt(:), &
417  xwnmean(:)
418  !
419  ! Output fields group 3)
420  !
421  REAL, POINTER :: ef(:,:), th1m(:,:), sth1m(:,:), &
422  th2m(:,:), sth2m(:,:) !, WN(:,:)
423  REAL, POINTER :: xef(:,:), xth1m(:,:), xsth1m(:,:),&
424  xth2m(:,:), xsth2m(:,:) !, XWN(:,:)
425  !
426  ! Output fields group 4)
427  !
428  REAL, POINTER :: phs(:,:), ptp(:,:), plp(:,:), &
429  pdir(:,:), psi(:,:), pws(:,:), &
430  pwst(:), pnr(:), pgw(:,:), &
431  pthp0(:,:), pqp(:,:), ppe(:,:), &
432  psw(:,:), ptm1(:,:), pt1(:,:), &
433  pt2(:,:), pep(:,:)
434  REAL, POINTER :: xphs(:,:), xptp(:,:), xplp(:,:), &
435  xpdir(:,:), xpsi(:,:), xpws(:,:), &
436  xpwst(:), xpnr(:), xpgw(:,:), &
437  xpthp0(:,:), xpqp(:,:), xppe(:,:), &
438  xpsw(:,:), xptm1(:,:), xpt1(:,:), &
439  xpt2(:,:), xpep(:,:)
440  !
441  ! Output fields group 5)
442  !
443  REAL, POINTER :: charn(:), cge(:), phiaw(:), &
444  tauwix(:), tauwiy(:), tauwnx(:), &
445  tauwny(:), whitecap(:,:), tws(:)
446  REAL, POINTER :: xcharn(:), xcge(:), xphiaw(:), &
447  xtauwix(:), xtauwiy(:), xtauwnx(:), &
448  xtauwny(:), xwhitecap(:,:), xtws(:)
449  !
450  ! Output fields group 6)
451  !
452  REAL, POINTER :: sxx(:), syy(:), sxy(:), tauox(:),&
453  tauoy(:), bhd(:), phioc(:), &
454  tusx(:), tusy(:), ussx(:), &
455  ussy(:), tauocx(:), tauocy(:), &
456  prms(:), tpms(:), phice(:), &
457  tauice(:,:)
458  REAL, POINTER :: p2sms(:,:), us3d(:,:), ussp(:,:)
459  REAL, POINTER :: xsxx(:), xsyy(:), xsxy(:), xtauox(:),&
460  xtauoy(:), xbhd(:), xphioc(:), &
461  xtusx(:), xtusy(:), xussx(:), &
462  xussy(:), xtauocx(:), xtauocy(:), &
463  xprms(:), xtpms(:), xphice(:), &
464  xtauice(:,:)
465  REAL, POINTER :: xp2sms(:,:), xus3d(:,:), xussp(:,:)
466  !
467  ! Output fields group 7)
468  !
469  REAL, POINTER :: aba(:), abd(:), uba(:), ubd(:), &
470  bedforms(:,:), phibbl(:), &
471  taubbl(:,:)
472  REAL, POINTER :: xaba(:), xabd(:), xuba(:), xubd(:), &
473  xbedforms(:,:), xphibbl(:), &
474  xtaubbl(:,:)
475  !
476  ! Output fields group 8)
477  !
478  REAL, POINTER :: mssx(:), mssy(:), mssd(:), &
479  mscx(:), mscy(:), mscd(:), qkk(:), skew(:), embia1(:), embia2(:)
480  REAL, POINTER :: xmssx(:), xmssy(:), xmssd(:), &
481  xmscx(:), xmscy(:), xmscd(:), xqkk(:), &
482  xskew(:), xembia1(:), xembia2(:)
483  !
484  ! Output fields group 9)
485  !
486  REAL, POINTER :: dtdyn(:), fcut(:), cflxymax(:), &
487  cflthmax(:), cflkmax(:)
488  REAL, POINTER :: xdtdyn(:), xfcut(:), xcflxymax(:), &
489  xcflthmax(:), xcflkmax(:)
490  !
491  ! Output fields group 10)
492  !
493  REAL, POINTER :: usero(:,:)
494  REAL, POINTER :: xusero(:,:)
495  !
496  ! Spatial derivatives
497  !
498  REAL, POINTER :: dddx(:,:), dddy(:,:), dcxdx(:,:), &
499  dcydx(:,:), dcxdy(:,:), dcydy(:,:)
500  REAL, POINTER :: dcdx(:,:,:), dcdy(:,:,:)
501 #ifdef W3_SMC
502  REAL, POINTER :: dhdx(:), dhdy(:), dhlmt(:,:)
503 #endif
504  !
505 #ifdef W3_PR1
506  INTEGER, POINTER :: is0(:), is2(:)
507  REAL, POINTER :: facvx(:), facvy(:)
508 #endif
509  !
510 #ifdef W3_PR2
511  INTEGER :: nmx0, nmx1, nmx2, nmy0, nmy1, nmy2, &
512  nact, nmxy
513  INTEGER, POINTER :: mapx2(:), mapy2(:), mapaxy(:), &
514  mapxy(:), mapth2(:), mapwn2(:)
515 #endif
516  !
517 #ifdef W3_PR3
518  INTEGER :: nmx0, nmx1, nmx2, nmy0, nmy1, nmy2, &
519  nact, ncent
520  INTEGER, POINTER :: mapx2(:), mapy2(:), mapaxy(:), &
521  mapcxy(:), mapth2(:), mapwn2(:)
522  LOGICAL, POINTER :: maptrn(:)
523 #endif
524  !
525  ! Warning Defined but not set if UGTYPE .EQ. .T.
526  INTEGER, POINTER :: iter(:,:)
527  !
528 #ifdef W3_NL1
529  INTEGER :: nfr, nfrhgh, nfrchg, nspecx, nspecy
530  INTEGER, POINTER :: ip11(:), ip12(:), ip13(:), ip14(:), &
531  im11(:), im12(:), im13(:), im14(:), &
532  ip21(:), ip22(:), ip23(:), ip24(:), &
533  im21(:), im22(:), im23(:), im24(:), &
534  ic11(:), ic12(:), ic21(:), ic22(:), &
535  ic31(:), ic32(:), ic41(:), ic42(:), &
536  ic51(:), ic52(:), ic61(:), ic62(:), &
537  ic71(:), ic72(:), ic81(:), ic82(:)
538  REAL :: dal1, dal2, dal3, &
539  awg1, awg2, awg3, awg4, awg5, awg6, &
540  awg7, awg8, swg1, swg2, swg3, swg4, &
541  swg5, swg6, swg7, swg8
542  REAL, POINTER :: af11(:)
543  LOGICAL :: nlinit
544 #endif
545  !
546  INTEGER, POINTER :: iappro(:)
547 #ifdef W3_MPI
550  nrqsg1 = 0, nrqsg2, ibfloc, isploc, &
551  nsploc
552 #endif
553 #ifdef W3_PDLIB
554  INTEGER :: nbfield, pdlib_mpi_type
555 #endif
556 #ifdef W3_MPI
557  INTEGER :: bstat(mpibuf), bispl(mpibuf)
558  INTEGER, POINTER :: irqsg1(:,:), irqsg2(:,:)
559  REAL, POINTER :: gstore(:,:), sstore(:,:)
560 #endif
561  REAL, POINTER :: sppnt(:,:,:)
562  !
563  INTEGER :: itime, ipass, idlast, nsealm
564  REAL, POINTER :: alpha(:,:)
565  LOGICAL :: ainit, ainit2, fl_all, flcold, fliwnd
566  !
567  END TYPE wadat
568  !/
569  !/ Data storage
570  !/
571  TYPE(wadat), TARGET, ALLOCATABLE :: wadats(:)
572  !/
573  !/ Data aliases for structure WADAT(S)
574  !/
575  REAL, POINTER :: cg(:,:), wn(:,:)
576  REAL, POINTER :: ic3wn_r(:,:), ic3wn_i(:,:), ic3cg(:,:)
577  !
578  REAL, POINTER :: ca0(:), cai(:), cd0(:), cdi(:), &
579  ua0(:), uai(:), ud0(:), udi(:), &
580  ma0(:), mai(:), ra0(:), rai(:), &
581  md0(:), mdi(:), as0(:), asi(:), &
582  atrnx(:,:), atrny(:,:)
583  !
584  REAL, POINTER :: dw(:), ua(:), ud(:), u10(:), u10d(:),&
585  as(:), cx(:), cy(:), taua(:), tauadir(:)
586  !
587  REAL, POINTER :: hs(:), wlm(:), t02(:), t0m1(:), &
588  t01 (:), fp0(:), thm(:), ths(:), &
589  thp0(:), hsig(:), &
590  stmaxe(:), stmaxd(:), hmaxe(:), &
591  hcmaxe(:), hmaxd(:), hcmaxd(:), &
592  qp(:), wbt(:), wnmean(:)
593  !
594  REAL, POINTER :: ef(:,:), th1m(:,:), sth1m(:,:), &
595  th2m(:,:), sth2m(:,:)
596  !
597  REAL, POINTER :: phs(:,:), ptp(:,:), plp(:,:), &
598  pdir(:,:), psi(:,:), pws(:,:), &
599  pwst(:), pnr(:), pgw(:,:), psw(:,:), &
600  pthp0(:,:), pqp(:,:), ppe(:,:), &
601  ptm1(:,:), pt1(:,:), pt2(:,:),pep(:,:)
602  !
603  REAL, POINTER :: charn(:), cge(:), phiaw(:), &
604  tauwix(:), tauwiy(:), tauwnx(:), &
605  tauwny(:), whitecap(:,:), tws(:)
606  !
607  REAL, POINTER :: sxx(:), syy(:), sxy(:), tauox(:), &
608  tauoy(:), bhd(:), phioc(:), &
609  tusx(:), tusy(:), ussx(:), ussy(:), &
610  tauocx(:), tauocy(:), prms(:), &
611  tpms(:), phice(:), tauice(:,:)
612  REAL, POINTER :: p2sms(:,:), us3d(:,:), ussp(:,:)
613  !
614  REAL, POINTER :: aba(:), abd(:), uba(:), ubd(:), &
615  bedforms(:,:), phibbl(:), taubbl(:,:)
616  !
617  REAL, POINTER :: mssx(:), mssy(:), mssd(:), &
618  mscx(:), mscy(:), mscd(:), qkk(:), skew(:), embia1(:), embia2(:)
619  !
620  REAL, POINTER :: dtdyn(:), fcut(:), cflxymax(:), &
621  cflthmax(:), cflkmax(:)
622  !
623  REAL, POINTER :: usero(:,:)
624  !
625  ! REAL, POINTER :: TAUWX(:), TAUWY(:)
626  !
627  REAL, POINTER :: dddx(:,:), dddy(:,:), dcxdx(:,:), &
628  dcydx(:,:), dcxdy(:,:), dcydy(:,:)
629  REAL, POINTER :: dcdx(:,:,:), dcdy(:,:,:)
630 #ifdef W3_SMC
631  REAL, POINTER :: dhdx(:), dhdy(:), dhlmt(:,:)
632 #endif
633  !
634 #ifdef W3_PR1
635  INTEGER, POINTER :: is0(:), is2(:)
636  REAL, POINTER :: facvx(:), facvy(:)
637 #endif
638  !
639 #ifdef W3_PR2
640  INTEGER, POINTER :: nmx0, nmx1, nmx2, nmy0, nmy1, nmy2, &
641  nact, nmxy
642  INTEGER, POINTER :: mapx2(:), mapy2(:), mapaxy(:), &
643  mapxy(:), mapth2(:), mapwn2(:)
644 #endif
645  !
646 #ifdef W3_PR3
647  INTEGER, POINTER :: nmx0, nmx1, nmx2, nmy0, nmy1, nmy2, &
648  nact, ncent
649  INTEGER, POINTER :: mapx2(:), mapy2(:), mapaxy(:), &
650  mapcxy(:), mapth2(:), mapwn2(:)
651  LOGICAL, POINTER :: maptrn(:)
652 #endif
653  !
654  INTEGER, POINTER :: iter(:,:)
655  !
656 #ifdef W3_NL1
657  INTEGER, POINTER :: nfr, nfrhgh, nfrchg, nspecx, nspecy
658  INTEGER, POINTER :: ip11(:), ip12(:), ip13(:), ip14(:), &
659  im11(:), im12(:), im13(:), im14(:), &
660  ip21(:), ip22(:), ip23(:), ip24(:), &
661  im21(:), im22(:), im23(:), im24(:), &
662  ic11(:), ic12(:), ic21(:), ic22(:), &
663  ic31(:), ic32(:), ic41(:), ic42(:), &
664  ic51(:), ic52(:), ic61(:), ic62(:), &
665  ic71(:), ic72(:), ic81(:), ic82(:)
666  REAL, POINTER :: dal1, dal2, dal3, &
667  awg1, awg2, awg3, awg4, awg5, awg6, &
668  awg7, awg8, swg1, swg2, swg3, swg4, &
669  swg5, swg6, swg7, swg8
670  REAL, POINTER :: af11(:)
671  LOGICAL, POINTER :: nlinit
672 #endif
673  !
674  INTEGER, POINTER :: iappro(:)
675 #ifdef W3_MPI
676  INTEGER, POINTER :: mpi_comm_wave, mpi_comm_wcmp, &
678  nrqsg1, nrqsg2, ibfloc, isploc, &
679  nsploc
680  INTEGER, POINTER :: bstat(:), bispl(:)
681  INTEGER, POINTER :: irqsg1(:,:), irqsg2(:,:)
682  REAL, POINTER :: gstore(:,:), sstore(:,:)
683 #endif
684  REAL, POINTER :: sppnt(:,:,:)
685  !
686  INTEGER, POINTER :: itime, ipass, idlast, nsealm
687  REAL, POINTER :: alpha(:,:)
688  LOGICAL, POINTER :: ainit, ainit2, fl_all, flcold, fliwnd
689  !/
690 CONTAINS
691  !/ ------------------------------------------------------------------- /
703  SUBROUTINE w3naux ( NDSE, NDST )
704  !/
705  !/ +-----------------------------------+
706  !/ | WAVEWATCH III NOAA/NCEP |
707  !/ | H. L. Tolman |
708  !/ | FORTRAN 90 |
709  !/ | Last update : 10-Dec-2014 !
710  !/ +-----------------------------------+
711  !/
712  !/ 14-Dec-2004 : Origination. ( version 3.06 )
713  !/ 04-Oct-2006 : Add filter to array pointers. ( version 3.10 )
714  !/ 10-Dec-2014 : Add checks for allocate status ( version 5.04 )
715  !/
716  ! 1. Purpose :
717  !
718  ! Set up the number of grids to be used.
719  !
720  ! 2. Method :
721  !
722  ! Use data stored in NGRIDS in W3GDATMD.
723  !
724  ! 3. Parameters :
725  !
726  ! Parameter list
727  ! ----------------------------------------------------------------
728  ! NDSE Int. I Error output unit number.
729  ! NDST Int. I Test output unit number.
730  ! ----------------------------------------------------------------
731  !
732  ! 4. Subroutines used :
733  !
734  ! See module documentation.
735  !
736  ! 5. Called by :
737  !
738  ! Any program that uses this grid structure.
739  !
740  ! 6. Error messages :
741  !
742  ! - Error checks on previous setting of variable NGRIDS.
743  !
744  ! 7. Remarks :
745  !
746  ! 8. Structure :
747  !
748  ! 9. Switches :
749  !
750  ! !/S Enable subroutine tracing.
751  ! !/T Enable test output
752  !
753  ! 10. Source code :
754  !
755  !/ ------------------------------------------------------------------- /
756  USE w3gdatmd, ONLY: ngrids
757  USE w3servmd, ONLY: extcde
758  USE w3odatmd, ONLY: iaproc
759 #ifdef W3_S
760  USE w3servmd, ONLY: strace
761 #endif
762  !
763  !/
764  !/ ------------------------------------------------------------------- /
765  !/ Parameter list
766  !/
767  INTEGER, INTENT(IN) :: NDSE, NDST
768  !/
769  !/ ------------------------------------------------------------------- /
770  !/ Local parameters
771  !/
772  INTEGER :: I
773 #ifdef W3_S
774  INTEGER, SAVE :: IENT = 0
775  CALL strace (ient, 'W3NAUX')
776 #endif
777  !
778  ! -------------------------------------------------------------------- /
779  ! 1. Test input and module status
780  !
781  IF ( ngrids .EQ. -1 ) THEN
782  WRITE (ndse,1001) ngrids
783  CALL extcde (1)
784  END IF
785  !
786  ! -------------------------------------------------------------------- /
787  ! 2. Set variable and allocate arrays
788  !
789  ALLOCATE ( wadats(ngrids), stat=istat )
790  check_alloc_status( istat )
791  nadata = ngrids
792  !
793  ! -------------------------------------------------------------------- /
794  ! 3. Initialize parameters
795  !
796  DO i=1, ngrids
797  wadats(i)%ITIME = 0
798  wadats(i)%IPASS = 0
799  wadats(i)%IDLAST = 0
800  wadats(i)%NSEALM = 0
801  wadats(i)%FLCOLD = .false.
802  wadats(i)%FLIWND = .false.
803  wadats(i)%AINIT = .false.
804  wadats(i)%AINIT2 = .false.
805  wadats(i)%FL_ALL = .false.
806 #ifdef W3_NL1
807  wadats(i)%NLINIT = .false.
808 #endif
809  END DO
810  !
811 #ifdef W3_T
812  WRITE (ndst,9000) ngrids
813 #endif
814  !
815  RETURN
816  !
817  ! Formats
818  !
819 1001 FORMAT (/' *** ERROR W3NAUX : NGRIDS NOT YET SET *** '/ &
820  ' NGRIDS = ',i10/ &
821  ' RUN W3NMOD FIRST'/)
822  !
823 #ifdef W3_T
824 9000 FORMAT (' TEST W3NAUX : SETTING UP FOR ',i4,' GRIDS')
825 #endif
826  !/
827  !/ End of W3NAUX ----------------------------------------------------- /
828  !/
829  END SUBROUTINE w3naux
830  !/ ------------------------------------------------------------------- /
845  SUBROUTINE w3dima ( IMOD, NDSE, NDST, D_ONLY )
846  !/
847  !/ +-----------------------------------+
848  !/ | WAVEWATCH III NOAA/NCEP |
849  !/ | H. L. Tolman |
850  !/ | FORTRAN 90 |
851  !/ | Last update : 22-Mar-2021 !
852  !/ +-----------------------------------+
853  !/
854  !/ 28-Dec-2004 : Origination. ( version 3.06 )
855  !/ 20-Jul-2005 : Adding output fields. ( version 3.07 )
856  !/ 04-Oct-2006 : Add filter to array pointers. ( version 3.10 )
857  !/ 28-Mar-2007 : Add partitioned data arrays. ( version 3.11 )
858  !/ Add additional undefined arrays.
859  !/ 22-Feb-2008 ; Modify MAPTH2 declaration. ( version 3.14 )
860  !/ 31-Oct-2010 : Added initialization of CX,CY,DW ( version 3.14 )
861  !/ 25-Dec-2012 : Memory reduction for outputs. ( version 4.11 )
862  !/ 28-Jul-2013 : Bug fix initialization P2SMS. ( version 4.11 )
863  !/ 30-Apr-2014 : Memory reduction for group3. ( version 5.00 )
864  !/ 10-Dec-2014 : Add checks for allocate status ( version 5.04 )
865  !/ 22-Mar-2021 : Adds TAUA, WNMEAN, TAUOC parameters ( version 7.13 )
866  !/
867  ! 1. Purpose :
868  !
869  ! Initialize an individual data grid at the proper dimensions.
870  !
871  ! 2. Method :
872  !
873  ! Allocate directly into the structure array. Note that
874  ! this cannot be done through the pointer alias!
875  !
876  ! 3. Parameters :
877  !
878  ! Parameter list
879  ! ----------------------------------------------------------------
880  ! IMOD Int. I Model number to point to.
881  ! NDSE Int. I Error output unit number.
882  ! NDST Int. I Test output unit number.
883  ! D_ONLY L.O. I FLag for initializing data arrays only.
884  ! ----------------------------------------------------------------
885  !
886  ! 4. Subroutines used :
887  !
888  ! See module documentation.
889  !
890  ! 5. Called by :
891  !
892  ! Name Type Module Description
893  ! ----------------------------------------------------------------
894  ! W3IOGO Subr. W3IOGOMD Grid output IO routine.
895  ! WW3_SHEL Prog. N/A Wave model driver.
896  ! ----------------------------------------------------------------
897  !
898  ! 6. Error messages :
899  !
900  ! - Check on input parameters.
901  ! - Check on previous allocation.
902  !
903  ! 7. Remarks :
904  !
905  ! - W3SETA needs to be called after allocation to point to
906  ! proper allocated arrays.
907  !
908  ! 8. Structure :
909  !
910  ! See source code.
911  !
912  ! 9. Switches :
913  !
914  ! !/SHRD, !/DIST
915  ! Shared / distributed memory model
916  !
917  ! !/PRn Propagation scheme selection.
918  !
919  ! !/S Enable subroutine tracing.
920  ! !/T Enable test output
921  !
922  ! 10. Source code :
923  !
924  !/ ------------------------------------------------------------------- /
925  USE constants, ONLY : lpdlib
926  USE w3gdatmd, ONLY: ngrids, igrid, w3setg, nk, nx, ny, nsea, &
927  nseal, nspec, nth, e3df, p2msf, us3df, &
929  USE w3odatmd, ONLY: iaproc, naproc, ntproc, napfld, &
930  noswll, noextr, undef, flogrd, flogr2
931  USE w3idatmd, ONLY: flcur, flwind, fltaua, flrhoa
932  USE w3servmd, ONLY: extcde
933 #ifdef W3_S
934  USE w3servmd, ONLY: strace
935 #endif
936  !
937  !/
938  !/ ------------------------------------------------------------------- /
939  !/ Parameter list
940  !/
941  INTEGER, INTENT(IN) :: IMOD, NDSE, NDST
942  LOGICAL, INTENT(IN), OPTIONAL :: D_ONLY
943  !/
944  !/ ------------------------------------------------------------------- /
945  !/ Local parameters
946  !/
947  INTEGER :: JGRID, NXXX, NSEAL_tmp
948 #ifdef W3_S
949  INTEGER, SAVE :: IENT = 0
950  CALL strace (ient, 'W3DIMA')
951 #endif
952  integer :: memunit
953  !
954  ! -------------------------------------------------------------------- /
955  ! 1. Test input and module status
956  !
957  memunit = 30000+iaproc
958  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 0')
959 
960  IF ( PRESENT(d_only) ) THEN
961  fl_all = .NOT. d_only
962  ELSE
963  fl_all = .true.
964  END IF
965  !
966  IF ( ngrids .EQ. -1 ) THEN
967  WRITE (ndse,1001)
968  CALL extcde (1)
969  END IF
970  !
971  IF ( imod.LT.1 .OR. imod.GT.nadata ) THEN
972  WRITE (ndse,1002) imod, nadata
973  CALL extcde (2)
974  END IF
975  !
976  IF ( wadats(imod)%AINIT ) THEN
977  WRITE (ndse,1003)
978  CALL extcde (3)
979  END IF
980  !
981 #ifdef W3_T
982  WRITE (ndst,9000) imod
983 #endif
984  !
985  jgrid = igrid
986  IF ( jgrid .NE. imod ) CALL w3setg ( imod, ndse, ndst )
987 
988  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 1')
989  !
990  ! -------------------------------------------------------------------- /
991  ! 2. Allocate arrays
992  ! Call W3SETA to assure of pointes FLCUR, FLWND, and FLTAUA
993  !
994  CALL w3seta ( imod, ndse, ndst )
995 
996  !
997  !AR: Check this below more ...
998  nxxx = nsealm * naproc
999  !
1000  ! Output and input parameteres by output type
1001  !
1002  ! 1) Forcing fields (these arrays are always needed)
1003  !
1004  ALLOCATE ( wadats(imod)%DW(0:nsea) , stat=istat )
1005  check_alloc_status( istat )
1006  wadats(imod)%DW(:)=0.
1007  !
1008  ALLOCATE ( wadats(imod)%CX(0:nsea) , wadats(imod)%CY(0:nsea) , &
1009  stat=istat )
1010  check_alloc_status( istat )
1011  wadats(imod)%CX(:)=0.
1012  wadats(imod)%CY(:)=0.
1013  !
1014  ALLOCATE ( wadats(imod)%UA(0:nsea) , wadats(imod)%UD(0:nsea) , &
1015  wadats(imod)%U10(nsea) , wadats(imod)%U10D(nsea) , &
1016  wadats(imod)%AS(0:nsea) , stat=istat )
1017  check_alloc_status( istat )
1018  !
1019  ALLOCATE ( wadats(imod)%TAUA(0:nsea) , &
1020  wadats(imod)%TAUADIR(0:nsea), stat=istat )
1021  check_alloc_status( istat )
1022  wadats(imod)%TAUA(:) =0.
1023  wadats(imod)%TAUADIR(:)=0.
1024 
1025  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 2')
1026  !
1027  ! Water level WLV stored in W3WDATMD
1028  ! Ice concentration ICE stored in W3WDATMD
1029  ! Ice floe sizes ICEF and ICEDMAX stored in W3WDATMD
1030  ! Iceberg damping BERG stored in W3WDATMD
1031  !
1032  ! 2) Standard mean wave parameters
1033  ! Here, all short arrays are always allocated to reduce logical
1034  ! checks in all computations. The coresponding full size arrays
1035  ! are allocated in W3MPIO only as needed to keep the memory
1036  ! footprint down.
1037  !
1038  IF (nsealm .eq. 0) THEN
1039  nsealm=nsea
1040  END IF
1041  ALLOCATE ( wadats(imod)%HS (nsealm), wadats(imod)%WLM (nsealm), &
1042  wadats(imod)%T02 (nsealm), wadats(imod)%T0M1(nsealm), &
1043  wadats(imod)%T01 (nsealm), wadats(imod)%FP0 (nsealm), &
1044  wadats(imod)%THM (nsealm), wadats(imod)%THS (nsealm), &
1045  wadats(imod)%THP0 (nsealm), wadats(imod)%HSIG(nsealm), &
1046  wadats(imod)%STMAXE (nsealm), &
1047  wadats(imod)%STMAXD(nsealm), &
1048  wadats(imod)%HMAXE(nsealm), wadats(imod)%HMAXD(nsealm), &
1049  wadats(imod)%HCMAXE(nsealm), &
1050  wadats(imod)%HCMAXD(nsealm), wadats(imod)%QP(nsealm), &
1051  wadats(imod)%WBT(nsealm), &
1052  wadats(imod)%WNMEAN(nsealm), &
1053  stat=istat )
1054  check_alloc_status( istat )
1055  !
1056  wadats(imod)%HS = undef
1057  wadats(imod)%WLM = undef
1058  wadats(imod)%T02 = undef
1059  wadats(imod)%T0M1 = undef
1060  wadats(imod)%T01 = undef
1061  wadats(imod)%FP0 = undef
1062  wadats(imod)%THM = undef
1063  wadats(imod)%THS = undef
1064  wadats(imod)%THP0 = undef
1065  wadats(imod)%HSIG = undef
1066  wadats(imod)%STMAXE = undef
1067  wadats(imod)%STMAXD = undef
1068  wadats(imod)%HMAXE = undef
1069  wadats(imod)%HMAXD = undef
1070  wadats(imod)%HCMAXE = undef
1071  wadats(imod)%HCMAXD = undef
1072  wadats(imod)%QP = undef
1073  wadats(imod)%WBT = undef
1074  wadats(imod)%WNMEAN = undef
1075 
1076  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 3')
1077  !
1078  ! 3) Frequency-dependent standard parameters
1079  !
1080  ! For the 3D arrays: the allocation is performed only if these arrays are allowed
1081  ! by specific variables defined through the mod_def file
1082  ! and read by w3iogr, which is called before W3DIMA.
1083  IF ( e3df(1,1).GT.0 ) THEN
1084  ALLOCATE(wadats(imod)%EF(nsealm,e3df(2,1):e3df(3,1)), &
1085  stat=istat )
1086  check_alloc_status( istat )
1087  END IF
1088  IF ( e3df(1,2).GT.0 ) THEN
1089  ALLOCATE(wadats(imod)%TH1M(nsealm,e3df(2,2):e3df(3,2)), &
1090  stat=istat )
1091  check_alloc_status( istat )
1092  END IF
1093  IF ( e3df(1,3).GT.0 ) THEN
1094  ALLOCATE(wadats(imod)%STH1M(nsealm,e3df(2,3):e3df(3,3)), &
1095  stat=istat )
1096  check_alloc_status( istat )
1097  END IF
1098  IF ( e3df(1,4).GT.0 ) THEN
1099  ALLOCATE(wadats(imod)%TH2M(nsealm,e3df(2,4):e3df(3,4)), &
1100  stat=istat )
1101  check_alloc_status( istat )
1102  END IF
1103  IF ( e3df(1,5).GT.0 ) THEN
1104  ALLOCATE(wadats(imod)%STH2M(nsealm,e3df(2,5):e3df(3,5)), &
1105  stat=istat )
1106  check_alloc_status( istat )
1107  END IF
1108  !
1109  IF ( e3df(1,1).GT.0 ) wadats(imod)%EF = undef
1110  IF ( e3df(1,2).GT.0 ) wadats(imod)%TH1M = undef
1111  IF ( e3df(1,3).GT.0 ) wadats(imod)%STH1M = undef
1112  IF ( e3df(1,4).GT.0 ) wadats(imod)%TH2M = undef
1113  IF ( e3df(1,5).GT.0 ) wadats(imod)%STH2M = undef
1114 
1115  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 4')
1116  !
1117  ! 4) Spectral Partitions parameters
1118  !
1119  ALLOCATE ( wadats(imod)%PHS(nsealm,0:noswll), &
1120  wadats(imod)%PTP(nsealm,0:noswll), &
1121  wadats(imod)%PLP(nsealm,0:noswll), &
1122  wadats(imod)%PDIR(nsealm,0:noswll), &
1123  wadats(imod)%PSI(nsealm,0:noswll), &
1124  wadats(imod)%PWS(nsealm,0:noswll), &
1125  wadats(imod)%PWST(nsealm), &
1126  wadats(imod)%PNR(nsealm), &
1127  wadats(imod)%PTHP0(nsealm,0:noswll), &
1128  wadats(imod)%PQP(nsealm,0:noswll), &
1129  wadats(imod)%PPE(nsealm,0:noswll), &
1130  wadats(imod)%PGW(nsealm,0:noswll), &
1131  wadats(imod)%PSW(nsealm,0:noswll), &
1132  wadats(imod)%PTM1(nsealm,0:noswll), &
1133  wadats(imod)%PT1(nsealm,0:noswll), &
1134  wadats(imod)%PT2(nsealm,0:noswll), &
1135  wadats(imod)%PEP(nsealm,0:noswll), &
1136  stat=istat )
1137  check_alloc_status( istat )
1138  !
1139  wadats(imod)%PHS = undef
1140  wadats(imod)%PTP = undef
1141  wadats(imod)%PLP = undef
1142  wadats(imod)%PDIR = undef
1143  wadats(imod)%PSI = undef
1144  wadats(imod)%PWS = undef
1145  wadats(imod)%PWST = undef
1146  wadats(imod)%PNR = undef
1147  wadats(imod)%PTHP0 = undef
1148  wadats(imod)%PQP = undef
1149  wadats(imod)%PPE = undef
1150  wadats(imod)%PGW = undef
1151  wadats(imod)%PSW = undef
1152  wadats(imod)%PTM1 = undef
1153  wadats(imod)%PT1 = undef
1154  wadats(imod)%PT2 = undef
1155  wadats(imod)%PEP = undef
1156  !
1157  ! 5) Atmosphere-waves layer
1158  !
1159  ! Friction velocity UST and USTDIR in W3WDATMD
1160  !
1161  ALLOCATE ( wadats(imod)%CHARN (nsealm), &
1162  wadats(imod)%TWS (nsealm), &
1163  wadats(imod)%CGE (nsealm), &
1164  wadats(imod)%PHIAW (nsealm), &
1165  wadats(imod)%TAUWIX (nsealm), &
1166  wadats(imod)%TAUWIY (nsealm), &
1167  wadats(imod)%TAUWNX (nsealm), &
1168  wadats(imod)%TAUWNY (nsealm), &
1169  wadats(imod)%WHITECAP(nsealm,4), &
1170  stat=istat )
1171  check_alloc_status( istat )
1172  !
1173  wadats(imod)%CHARN = undef
1174  wadats(imod)%TWS = undef
1175  wadats(imod)%CGE = undef
1176  wadats(imod)%PHIAW = undef
1177  wadats(imod)%TAUWIX = undef
1178  wadats(imod)%TAUWIY = undef
1179  wadats(imod)%TAUWNX = undef
1180  wadats(imod)%TAUWNY = undef
1181  wadats(imod)%WHITECAP = undef
1182 
1183  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 5')
1184  !
1185  ! 6) Wave-ocean layer
1186  !
1187  ALLOCATE ( wadats(imod)%SXX (nsealm) , &
1188  wadats(imod)%SYY (nsealm) , &
1189  wadats(imod)%SXY (nsealm) , &
1190  wadats(imod)%TAUOX (nsealm) , &
1191  wadats(imod)%TAUOY (nsealm) , &
1192  wadats(imod)%BHD (nsealm) , &
1193  wadats(imod)%PHIOC (nsealm) , &
1194  wadats(imod)%TUSX (nsealm) , &
1195  wadats(imod)%TUSY (nsealm) , &
1196  wadats(imod)%USSX (nsealm) , &
1197  wadats(imod)%USSY (nsealm) , &
1198  wadats(imod)%TAUOCX(nsealm) , &
1199  wadats(imod)%TAUOCY(nsealm) , &
1200  wadats(imod)%PRMS (nsealm) , &
1201  wadats(imod)%TPMS (nsealm) , &
1202  wadats(imod)%PHICE (nsealm) , &
1203  wadats(imod)%TAUICE(nsealm,2), &
1204  stat=istat )
1205  check_alloc_status( istat )
1206  !
1207  ! For the 3D arrays: the allocation is performed only if these arrays are allowed
1208  ! by specific variables defined through the mod_def file
1209  ! and read by w3iogr, which is called before W3DIMA.
1210  IF ( p2msf(1).GT.0 ) THEN
1211  ALLOCATE(wadats(imod)%P2SMS(nsealm,p2msf(2):p2msf(3)), stat=istat )
1212  check_alloc_status( istat )
1213  END IF
1214  IF ( us3df(1).GT.0 ) THEN ! maybe use US3DF(2:3)
1215  ALLOCATE(wadats(imod)%US3D(nsealm,nk*2), stat=istat )
1216  check_alloc_status( istat )
1217  END IF
1218  IF ( usspf(1).GT.0 ) THEN
1219  ALLOCATE(wadats(imod)%USSP(nsealm,nk*2), stat=istat )
1220  check_alloc_status( istat )
1221  END IF
1222  !
1223  wadats(imod)%SXX = undef
1224  wadats(imod)%SYY = undef
1225  wadats(imod)%SXY = undef
1226  wadats(imod)%TAUOX = undef
1227  wadats(imod)%TAUOY = undef
1228  wadats(imod)%BHD = undef
1229  wadats(imod)%PHIOC = undef
1230  wadats(imod)%TUSX = undef
1231  wadats(imod)%TUSY = undef
1232  wadats(imod)%USSX = undef
1233  wadats(imod)%USSY = undef
1234  wadats(imod)%TAUOCX = undef
1235  wadats(imod)%TAUOCY = undef
1236  wadats(imod)%PRMS = undef
1237  wadats(imod)%TPMS = undef
1238  wadats(imod)%PHICE = undef
1239  wadats(imod)%TAUICE = undef
1240  IF ( p2msf(1).GT.0 ) wadats(imod)%P2SMS = undef
1241  IF ( us3df(1).GT.0 ) wadats(imod)%US3D = undef
1242  IF ( usspf(1).GT.0 ) wadats(imod)%USSP = undef
1243 
1244  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 6')
1245  !
1246  ! 7) Wave-bottom layer
1247  !
1248  ALLOCATE ( wadats(imod)%ABA(nsealm) , wadats(imod)%ABD(nsealm) , &
1249  wadats(imod)%UBA(nsealm) , wadats(imod)%UBD(nsealm) , &
1250  wadats(imod)%BEDFORMS(nsealm,3), &
1251  wadats(imod)%PHIBBL (nsealm) , &
1252  wadats(imod)%TAUBBL (nsealm,2), stat=istat )
1253  check_alloc_status( istat )
1254  !
1255  wadats(imod)%ABA = undef
1256  wadats(imod)%ABD = undef
1257  wadats(imod)%UBA = undef
1258  wadats(imod)%UBD = undef
1259  wadats(imod)%BEDFORMS = undef
1260  wadats(imod)%PHIBBL = undef
1261  wadats(imod)%TAUBBL = undef
1262 
1263  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 7')
1264  !
1265  ! 8) Spectrum parameters
1266  !
1267  ALLOCATE ( wadats(imod)%MSSX(nsealm), wadats(imod)%MSSY(nsealm), &
1268  wadats(imod)%MSCX(nsealm), wadats(imod)%MSCY(nsealm), &
1269  wadats(imod)%MSSD(nsealm), wadats(imod)%MSCD(nsealm), &
1270  wadats(imod)%QKK(nsealm), wadats(imod)%SKEW(nsealm), &
1271  wadats(imod)%EMBIA1(nsealm), wadats(imod)%EMBIA2(nsealm), &
1272  stat=istat )
1273  check_alloc_status( istat )
1274  !
1275  wadats(imod)%MSSX = undef
1276  wadats(imod)%MSSY = undef
1277  wadats(imod)%MSSD = undef
1278  wadats(imod)%MSCX = undef
1279  wadats(imod)%MSCY = undef
1280  wadats(imod)%MSCD = undef
1281  wadats(imod)%QKK = undef
1282  wadats(imod)%SKEW = undef
1283  wadats(imod)%EMBIA1 = undef
1284  wadats(imod)%EMBIA2 = undef
1285  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 8')
1286  !
1287  ! 9) Numerical diagnostics
1288  !
1289  !
1290  ALLOCATE ( wadats(imod)%DTDYN (nsealm) , &
1291  wadats(imod)%FCUT (nsealm) , &
1292  wadats(imod)%CFLXYMAX(nsealm) , &
1293  wadats(imod)%CFLTHMAX(nsealm) , &
1294  wadats(imod)%CFLKMAX (nsealm) , stat=istat )
1295  check_alloc_status( istat )
1296  !
1297  wadats(imod)%DTDYN = undef
1298  wadats(imod)%FCUT = undef
1299  wadats(imod)%CFLXYMAX = undef
1300  wadats(imod)%CFLTHMAX = undef
1301  wadats(imod)%CFLKMAX = undef
1302 
1303  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 9')
1304  !
1305  ! 10) User defined
1306  !
1307  ALLOCATE ( wadats(imod)%USERO(nsealm,noextr), stat=istat )
1308  check_alloc_status( istat )
1309  !
1310  wadats(imod)%USERO = undef
1311  !
1312  ALLOCATE ( wadats(imod)%WN(0:nk+1,0:nsea), stat=istat )
1313  check_alloc_status( istat )
1314 
1315 #ifdef W3_IC3
1316  ALLOCATE (wadats(imod)%IC3WN_R(0:nk+1,0:300), stat=istat )
1317  check_alloc_status( istat )
1318  ALLOCATE (wadats(imod)%IC3WN_I(0:nk+1,0:300), stat=istat )
1319  check_alloc_status( istat )
1320 #endif
1321  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 10')
1322  !
1323  IF ( fl_all ) THEN
1324  !
1325  ALLOCATE ( wadats(imod)%CG(0:nk+1,0:nsea), stat=istat )
1326  check_alloc_status( istat )
1327 
1328 #ifdef W3_IC3
1329  ALLOCATE (wadats(imod)%IC3CG(0:nk+1,0:300), stat=istat )
1330  check_alloc_status( istat )
1331 #endif
1332 
1333  !
1334  IF ( flcur ) THEN
1335  ALLOCATE ( wadats(imod)%CA0(nsea) , &
1336  wadats(imod)%CAI(nsea) , &
1337  wadats(imod)%CD0(nsea) , &
1338  wadats(imod)%CDI(nsea) , &
1339  stat=istat )
1340  check_alloc_status( istat )
1341  END IF
1342  !
1343  IF ( flwind ) THEN
1344  ALLOCATE ( wadats(imod)%UA0(nsea) , &
1345  wadats(imod)%UAI(nsea) , &
1346  wadats(imod)%UD0(nsea) , &
1347  wadats(imod)%UDI(nsea) , &
1348  wadats(imod)%AS0(nsea) , &
1349  wadats(imod)%ASI(nsea) , &
1350  stat=istat )
1351  check_alloc_status( istat )
1352  END IF
1353  !
1354  IF ( fltaua ) THEN
1355  ALLOCATE ( wadats(imod)%MA0(nsea) , &
1356  wadats(imod)%MAI(nsea) , &
1357  wadats(imod)%MD0(nsea) , &
1358  wadats(imod)%MDI(nsea) , &
1359  stat=istat )
1360  check_alloc_status( istat )
1361  END IF
1362  !
1363  IF ( flrhoa ) THEN
1364  ALLOCATE ( wadats(imod)%RA0(nsea) , &
1365  wadats(imod)%RAI(nsea) , &
1366  stat=istat )
1367  check_alloc_status( istat )
1368  END IF
1369  !
1370  ALLOCATE ( wadats(imod)%ATRNX(ny*nx,-1:1) , &
1371  wadats(imod)%ATRNY(ny*nx,-1:1) , stat=istat )
1372  check_alloc_status( istat )
1373  !
1374  IF (.NOT. lpdlib) THEN
1375  ALLOCATE ( wadats(imod)%DDDX(ny,nx) , &
1376  wadats(imod)%DDDY(ny,nx) , &
1377  wadats(imod)%DCDX(0:nk+1,ny,nx) , &
1378  wadats(imod)%DCDY(0:nk+1,ny,nx) , &
1379  wadats(imod)%DCXDX(ny,nx) , &
1380  wadats(imod)%DCYDX(ny,nx) , &
1381  wadats(imod)%DCXDY(ny,nx) , &
1382  wadats(imod)%DCYDY(ny,nx) , stat=istat )
1383  ELSE
1384  ALLOCATE ( wadats(imod)%DDDX(1,nseal) , &
1385  wadats(imod)%DDDY(1,nseal) , &
1386  wadats(imod)%DCDX(0:nk+1,1,nseal) , &
1387  wadats(imod)%DCDY(0:nk+1,1,nseal) , &
1388  wadats(imod)%DCXDX(1,nseal) , &
1389  wadats(imod)%DCYDX(1,nseal) , &
1390  wadats(imod)%DCXDY(1,nseal) , &
1391  wadats(imod)%DCYDY(1,nseal) , &
1392  stat=istat )
1393  ENDIF
1394  check_alloc_status( istat )
1395  wadats(imod)%DDDX = 0.
1396  wadats(imod)%DDDY = 0.
1397  wadats(imod)%DCDX = 0.
1398  wadats(imod)%DCDY = 0.
1399  wadats(imod)%DCXDX = 0.
1400  wadats(imod)%DCYDX = 0.
1401  wadats(imod)%DCXDY = 0.
1402  wadats(imod)%DCYDY = 0.
1403  !
1404 #ifdef W3_SMC
1405  ALLOCATE ( wadats(imod)%DHDX(nsea) , &
1406  wadats(imod)%DHDY(nsea) , &
1407  wadats(imod)%DHLMT(nth,nsea) , stat=istat )
1408  check_alloc_status( istat )
1409 #endif
1410  !
1411  ALLOCATE ( wadats(imod)%ALPHA(nk,nseal) , stat=istat )
1412  check_alloc_status( istat )
1413  !
1414 #ifdef W3_PR1
1415  ALLOCATE ( wadats(imod)%IS0(nspec) , &
1416  wadats(imod)%IS2(nspec) , &
1417  wadats(imod)%FACVX(ny*nx) , &
1418  wadats(imod)%FACVY(ny*nx) , stat=istat )
1419  check_alloc_status( istat )
1420 #endif
1421  !
1422 #ifdef W3_PR2
1423  ALLOCATE ( wadats(imod)%MAPX2(ny*nx) , &
1424  wadats(imod)%MAPY2(ny*nx) , &
1425  wadats(imod)%MAPAXY(ny*nx) , &
1426  wadats(imod)%MAPXY(nsea) , &
1427  wadats(imod)%MAPTH2((nk+2)*nth) , &
1428  wadats(imod)%MAPWN2(nspec+nth) , stat=istat )
1429  check_alloc_status( istat )
1430  wadats(imod)%MAPTH2 = 0
1431 #endif
1432  !
1433  IF (gtype .EQ. ungtype) THEN
1434  ALLOCATE( wadats(imod)%ITER(nk,nth), stat=istat )
1435  check_alloc_status( istat )
1436  END IF
1437  !
1438 #ifdef W3_PR3
1439  ALLOCATE ( wadats(imod)%MAPX2(ny*nx) , &
1440  wadats(imod)%MAPY2(ny*nx) , &
1441  wadats(imod)%MAPAXY(ny*nx) , &
1442  wadats(imod)%MAPCXY(nsea) , &
1443  wadats(imod)%MAPTH2((nk+2)*nth) , &
1444  wadats(imod)%MAPWN2(nspec+nth) , &
1445  wadats(imod)%MAPTRN(ny*nx) , stat=istat )
1446  check_alloc_status( istat )
1447  wadats(imod)%MAPTH2 = 0
1448 #endif
1449  !
1450  ALLOCATE ( wadats(imod)%IAPPRO(nspec) , &
1451  wadats(imod)%SPPNT(nth,nk,4), stat=istat )
1452  check_alloc_status( istat )
1453  !
1454  END IF
1455  !
1456  wadats(imod)%AINIT = .true.
1457 
1458  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 11')
1459  !
1460 #ifdef W3_T
1461  WRITE (ndst,9001)
1462 #endif
1463  !
1464  ! -------------------------------------------------------------------- /
1465  ! 3. Point to allocated arrays
1466  !
1467  CALL w3seta ( imod, ndse, ndst )
1468 
1469  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA 12')
1470  !
1471 #ifdef W3_T
1472  WRITE (ndst,9002)
1473 #endif
1474  !
1475  ! -------------------------------------------------------------------- /
1476  ! 4. Update counters in grid
1477  !
1478 #ifdef W3_T
1479  WRITE (ndst,9003)
1480 #endif
1481  !
1482  ! -------------------------------------------------------------------- /
1483  ! 5. Restore previous grid setting if necessary
1484  !
1485  IF ( jgrid .NE. imod ) CALL w3setg ( jgrid, ndse, ndst )
1486 
1487  call print_memcheck(memunit, 'memcheck_____:'//' W3DIMA END')
1488  !
1489  RETURN
1490  !
1491  ! Formats
1492  !
1493 1001 FORMAT (/' *** ERROR W3DIMA : GRIDS NOT INITIALIZED *** '/ &
1494  ' RUN W3NMOD FIRST '/)
1495 1002 FORMAT (/' *** ERROR W3DIMA : ILLEGAL MODEL NUMBER *** '/ &
1496  ' IMOD = ',i10/ &
1497  ' NADATA = ',i10/)
1498 1003 FORMAT (/' *** ERROR W3DIMA : ARRAY(S) ALREADY ALLOCATED *** ')
1499  !
1500 #ifdef W3_T
1501 9000 FORMAT (' TEST W3DIMA : MODEL ',i4)
1502 9001 FORMAT (' TEST W3DIMA : ARRAYS ALLOCATED')
1503 9002 FORMAT (' TEST W3DIMA : POINTERS RESET')
1504 9003 FORMAT (' TEST W3DIMA : DIMENSIONS STORED')
1505 #endif
1506  !/
1507  !/ End of W3DIMA ----------------------------------------------------- /
1508  !/
1509  END SUBROUTINE w3dima
1510  !/ ------------------------------------------------------------------- /
1522  SUBROUTINE w3xdma ( IMOD, NDSE, NDST, OUTFLAGS )
1523  !/
1524  !/ +-----------------------------------+
1525  !/ | WAVEWATCH III NOAA/NCEP |
1526  !/ | H. L. Tolman |
1527  !/ | FORTRAN 90 |
1528  !/ | Last update : 22-Mar-2021 !
1529  !/ +-----------------------------------+
1530  !/
1531  !/ 26-Dec-2012 : Origination. ( version 3.06 )
1532  !/ 10-Dec-2014 : Add checks for allocate status ( version 5.04 )
1533  !/ 22-Mar-2021 : Adds WNMEAN, TAUOC parameters ( version 7.13 )
1534  !/
1535  ! 1. Purpose :
1536  !
1537  ! Version of W3DIMX for extended ouput arrays only.
1538  !
1539  !
1540  ! 10. Source code :
1541  !
1542  !/ ------------------------------------------------------------------- /
1543  USE w3gdatmd, ONLY: ngrids, igrid, w3setg, nk, nx, ny, nsea, &
1544  nseal, nspec, nth, e3df, p2msf, us3df, &
1545  usspf, gtype, ungtype
1546  USE w3odatmd, ONLY: iaproc, naproc, ntproc, napfld, &
1547  noswll, noextr, undef, flogrd, flogr2, &
1548  nogrp, ngrpp
1549  USE w3servmd, ONLY: extcde
1550 #ifdef W3_S
1551  USE w3servmd, ONLY: strace
1552 #endif
1553  !
1554  !/
1555  !/ ------------------------------------------------------------------- /
1556  !/ Parameter list
1557  !/
1558  INTEGER, INTENT(IN) :: IMOD, NDSE, NDST
1559  LOGICAL, INTENT(IN) :: OUTFLAGS(NOGRP,NGRPP)
1560  !/
1561  !/ ------------------------------------------------------------------- /
1562  !/ Local parameters
1563  !/
1564  INTEGER :: JGRID, NXXX, I
1565 #ifdef W3_S
1566  INTEGER, SAVE :: IENT = 0
1567  CALL strace (ient, 'W3XDMA')
1568 #endif
1569  integer :: memunit
1570  !
1571  ! -------------------------------------------------------------------- /
1572  ! 1. Test input and module status
1573  !
1574  memunit = 30000+iaproc
1575  IF ( ngrids .EQ. -1 ) THEN
1576  WRITE (ndse,1001)
1577  CALL extcde (1)
1578  END IF
1579  !
1580  IF ( imod.LT.1 .OR. imod.GT.nadata ) THEN
1581  WRITE (ndse,1002) imod, nadata
1582  CALL extcde (2)
1583  END IF
1584  !
1585  IF ( wadats(imod)%AINIT2 ) THEN
1586  WRITE (ndse,1003)
1587  CALL extcde (3)
1588  END IF
1589  !
1590 #ifdef W3_T
1591  WRITE (ndst,9000) imod
1592 #endif
1593  !
1594  jgrid = igrid
1595  IF ( jgrid .NE. imod ) CALL w3setg ( imod, ndse, ndst )
1596  !
1597  ! -------------------------------------------------------------------- /
1598  ! 2. Allocate arrays
1599  !
1600  nxxx = nsealm * naproc
1601  !
1602  IF ( outflags( 2, 1) ) THEN
1603  ALLOCATE ( wadats(imod)%XHS(nxxx), stat=istat )
1604  check_alloc_status( istat )
1605  ELSE
1606  ALLOCATE ( wadats(imod)%XHS(1), stat=istat )
1607  check_alloc_status( istat )
1608  END IF
1609  !
1610  IF ( outflags( 2, 2) ) THEN
1611  ALLOCATE ( wadats(imod)%XWLM(nxxx), stat=istat )
1612  check_alloc_status( istat )
1613  ELSE
1614  ALLOCATE ( wadats(imod)%XWLM(1), stat=istat )
1615  check_alloc_status( istat )
1616  END IF
1617  !
1618  IF ( outflags( 2, 3) ) THEN
1619  ALLOCATE ( wadats(imod)%XT02(nxxx), stat=istat )
1620  check_alloc_status( istat )
1621  ELSE
1622  ALLOCATE ( wadats(imod)%XT02(1), stat=istat )
1623  check_alloc_status( istat )
1624  END IF
1625  !
1626  IF ( outflags( 2, 4) ) THEN
1627  ALLOCATE ( wadats(imod)%XT0M1(nxxx), stat=istat )
1628  check_alloc_status( istat )
1629  ELSE
1630  ALLOCATE ( wadats(imod)%XT0M1(1), stat=istat )
1631  check_alloc_status( istat )
1632  END IF
1633  !
1634  IF ( outflags( 2, 5) ) THEN
1635  ALLOCATE ( wadats(imod)%XT01 (nxxx), stat=istat )
1636  check_alloc_status( istat )
1637  ELSE
1638  ALLOCATE ( wadats(imod)%XT01 (1), stat=istat )
1639  check_alloc_status( istat )
1640  END IF
1641  !
1642  IF ( outflags( 2, 6) .OR. outflags( 2,18) ) THEN
1643  ! TP output shares FP0 internal field with FP
1644  ALLOCATE ( wadats(imod)%XFP0(nxxx), stat=istat )
1645  check_alloc_status( istat )
1646  ELSE
1647  ALLOCATE ( wadats(imod)%XFP0(1), stat=istat )
1648  check_alloc_status( istat )
1649  END IF
1650  !
1651  IF ( outflags( 2, 7) ) THEN
1652  ALLOCATE ( wadats(imod)%XTHM(nxxx), stat=istat )
1653  check_alloc_status( istat )
1654  ELSE
1655  ALLOCATE ( wadats(imod)%XTHM(1), stat=istat )
1656  check_alloc_status( istat )
1657  END IF
1658  !
1659  IF ( outflags( 2, 8) ) THEN
1660  ALLOCATE ( wadats(imod)%XTHS(nxxx), stat=istat )
1661  check_alloc_status( istat )
1662  ELSE
1663  ALLOCATE ( wadats(imod)%XTHS(1), stat=istat )
1664  check_alloc_status( istat )
1665  END IF
1666  !
1667  IF ( outflags( 2, 9) ) THEN
1668  ALLOCATE ( wadats(imod)%XTHP0(nxxx), stat=istat )
1669  check_alloc_status( istat )
1670  ELSE
1671  ALLOCATE ( wadats(imod)%XTHP0(1), stat=istat )
1672  check_alloc_status( istat )
1673  END IF
1674  !
1675  IF ( outflags( 2, 10) ) THEN
1676  ALLOCATE ( wadats(imod)%XHSIG(nxxx), stat=istat )
1677  check_alloc_status( istat )
1678  ELSE
1679  ALLOCATE ( wadats(imod)%XHSIG(1), stat=istat )
1680  check_alloc_status( istat )
1681  END IF
1682  !
1683  IF ( outflags( 2, 11) ) THEN
1684  ALLOCATE ( wadats(imod)%XSTMAXE(nxxx) )
1685  ELSE
1686  ALLOCATE ( wadats(imod)%XSTMAXE(1) )
1687  END IF
1688  !
1689  IF ( outflags( 2, 12) ) THEN
1690  ALLOCATE ( wadats(imod)%XSTMAXD(nxxx) )
1691  ELSE
1692  ALLOCATE ( wadats(imod)%XSTMAXD(1) )
1693  END IF
1694  !
1695  IF ( outflags( 2, 13) ) THEN
1696  ALLOCATE ( wadats(imod)%XHMAXE(nxxx) )
1697  ELSE
1698  ALLOCATE ( wadats(imod)%XHMAXE(1) )
1699  END IF
1700  !
1701  IF ( outflags( 2, 14) ) THEN
1702  ALLOCATE ( wadats(imod)%XHCMAXE(nxxx) )
1703  ELSE
1704  ALLOCATE ( wadats(imod)%XHCMAXE(1) )
1705  END IF
1706  !
1707  !
1708  IF ( outflags( 2, 15) ) THEN
1709  ALLOCATE ( wadats(imod)%XHMAXD(nxxx) )
1710  ELSE
1711  ALLOCATE ( wadats(imod)%XHMAXD(1) )
1712  END IF
1713  !
1714  IF ( outflags( 2, 16) ) THEN
1715  ALLOCATE ( wadats(imod)%XHCMAXD(nxxx) )
1716  ELSE
1717  ALLOCATE ( wadats(imod)%XHCMAXD(1) )
1718  END IF
1719  !
1720  IF ( outflags( 2, 17) ) THEN
1721  ALLOCATE ( wadats(imod)%XWBT (nxxx), stat=istat )
1722  check_alloc_status( istat )
1723  ELSE
1724  ALLOCATE ( wadats(imod)%XWBT (1), stat=istat )
1725  check_alloc_status( istat )
1726  END IF
1727  !
1728  IF ( outflags( 2, 19) ) THEN
1729  ALLOCATE ( wadats(imod)%XWNMEAN(nxxx), stat=istat )
1730  check_alloc_status( istat )
1731  ELSE
1732  ALLOCATE ( wadats(imod)%XWNMEAN(1), stat=istat )
1733  check_alloc_status( istat )
1734  END IF
1735  !
1736  wadats(imod)%XHS = undef
1737  wadats(imod)%XWLM = undef
1738  wadats(imod)%XT02 = undef
1739  wadats(imod)%XT0M1 = undef
1740  wadats(imod)%XT01 = undef
1741  wadats(imod)%XFP0 = undef
1742  wadats(imod)%XTHM = undef
1743  wadats(imod)%XTHS = undef
1744  wadats(imod)%XTHP0 = undef
1745  wadats(imod)%XHSIG = undef
1746  wadats(imod)%XSTMAXE= undef
1747  wadats(imod)%XSTMAXD= undef
1748  wadats(imod)%XHMAXE = undef
1749  wadats(imod)%XHMAXD = undef
1750  wadats(imod)%XHCMAXE= undef
1751  wadats(imod)%XHCMAXD= undef
1752  wadats(imod)%XWBT = undef
1753  wadats(imod)%XWNMEAN= undef
1754  !
1755  IF ( outflags( 3, 1) ) THEN
1756  ALLOCATE ( wadats(imod)%XEF(nxxx,e3df(2,1):e3df(3,1)), stat=istat )
1757  check_alloc_status( istat )
1758  ELSE
1759  ALLOCATE ( wadats(imod)%XEF(1,1), stat=istat )
1760  check_alloc_status( istat )
1761  END IF
1762 
1763  IF ( outflags( 3, 2) ) THEN
1764  ALLOCATE ( wadats(imod)%XTH1M(nxxx,e3df(2,2):e3df(3,2)), stat=istat )
1765  check_alloc_status( istat )
1766  ELSE
1767  ALLOCATE ( wadats(imod)%XTH1M(1,1), stat=istat )
1768  check_alloc_status( istat )
1769  END IF
1770 
1771  IF ( outflags( 3, 3) ) THEN
1772  ALLOCATE ( wadats(imod)%XSTH1M(nxxx,e3df(2,3):e3df(3,3)), stat=istat )
1773  check_alloc_status( istat )
1774  ELSE
1775  ALLOCATE ( wadats(imod)%XSTH1M(1,1), stat=istat )
1776  check_alloc_status( istat )
1777  END IF
1778 
1779  IF ( outflags( 3, 4) ) THEN
1780  ALLOCATE ( wadats(imod)%XTH2M(nxxx,e3df(2,4):e3df(3,4)), stat=istat )
1781  check_alloc_status( istat )
1782  ELSE
1783  ALLOCATE ( wadats(imod)%XTH2M(1,1), stat=istat )
1784  check_alloc_status( istat )
1785  END IF
1786  !
1787  IF ( outflags( 3, 5) ) THEN
1788  ALLOCATE ( wadats(imod)%XSTH2M(nxxx,e3df(2,5):e3df(3,5)), stat=istat )
1789  check_alloc_status( istat )
1790  ELSE
1791  ALLOCATE ( wadats(imod)%XSTH2M(1,1), stat=istat )
1792  check_alloc_status( istat )
1793  END IF
1794  !
1795  wadats(imod)%XEF = undef
1796  wadats(imod)%XTH1M = undef
1797  wadats(imod)%XSTH1M = undef
1798  wadats(imod)%XTH2M = undef
1799  wadats(imod)%XSTH2M = undef
1800  !
1801  IF ( outflags( 4, 1) ) THEN
1802  ALLOCATE ( wadats(imod)%XPHS(nxxx,0:noswll), stat=istat )
1803  check_alloc_status( istat )
1804  ELSE
1805  ALLOCATE ( wadats(imod)%XPHS(1,1), stat=istat )
1806  check_alloc_status( istat )
1807  END IF
1808  !
1809  IF ( outflags( 4, 2) ) THEN
1810  ALLOCATE ( wadats(imod)%XPTP(nxxx,0:noswll), stat=istat )
1811  check_alloc_status( istat )
1812  ELSE
1813  ALLOCATE ( wadats(imod)%XPTP(1,1), stat=istat )
1814  check_alloc_status( istat )
1815  END IF
1816  !
1817  IF ( outflags( 4, 3) ) THEN
1818  ALLOCATE ( wadats(imod)%XPLP(nxxx,0:noswll), stat=istat )
1819  check_alloc_status( istat )
1820  ELSE
1821  ALLOCATE ( wadats(imod)%XPLP(1,1), stat=istat )
1822  check_alloc_status( istat )
1823  END IF
1824  !
1825  IF ( outflags( 4, 4) ) THEN
1826  ALLOCATE ( wadats(imod)%XPDIR(nxxx,0:noswll), stat=istat )
1827  check_alloc_status( istat )
1828  ELSE
1829  ALLOCATE ( wadats(imod)%XPDIR(1,1), stat=istat )
1830  check_alloc_status( istat )
1831  END IF
1832  !
1833  IF ( outflags( 4, 5) ) THEN
1834  ALLOCATE ( wadats(imod)%XPSI(nxxx,0:noswll), stat=istat )
1835  check_alloc_status( istat )
1836  ELSE
1837  ALLOCATE ( wadats(imod)%XPSI(1,1), stat=istat )
1838  check_alloc_status( istat )
1839  END IF
1840  !
1841  IF ( outflags( 4, 6) ) THEN
1842  ALLOCATE ( wadats(imod)%XPWS(nxxx,0:noswll), stat=istat )
1843  check_alloc_status( istat )
1844  ELSE
1845  ALLOCATE ( wadats(imod)%XPWS(1,1), stat=istat )
1846  check_alloc_status( istat )
1847  END IF
1848  !
1849  IF ( outflags( 4, 7) ) THEN
1850  ALLOCATE ( wadats(imod)%XPTHP0(nxxx,0:noswll), stat=istat )
1851  check_alloc_status( istat )
1852  ELSE
1853  ALLOCATE ( wadats(imod)%XPTHP0(1,1), stat=istat )
1854  check_alloc_status( istat )
1855  END IF
1856  !
1857  IF ( outflags( 4, 8) ) THEN
1858  ALLOCATE ( wadats(imod)%XPQP(nxxx,0:noswll), stat=istat )
1859  check_alloc_status( istat )
1860  ELSE
1861  ALLOCATE ( wadats(imod)%XPQP(1,1), stat=istat )
1862  check_alloc_status( istat )
1863  END IF
1864  !
1865  IF ( outflags( 4, 9) ) THEN
1866  ALLOCATE ( wadats(imod)%XPPE(nxxx,0:noswll), stat=istat )
1867  check_alloc_status( istat )
1868  ELSE
1869  ALLOCATE ( wadats(imod)%XPPE(1,1), stat=istat )
1870  check_alloc_status( istat )
1871  END IF
1872  !
1873  IF ( outflags( 4,10) ) THEN
1874  ALLOCATE ( wadats(imod)%XPGW(nxxx,0:noswll), stat=istat )
1875  check_alloc_status( istat )
1876  ELSE
1877  ALLOCATE ( wadats(imod)%XPGW(1,1), stat=istat )
1878  check_alloc_status( istat )
1879  END IF
1880  !
1881  IF ( outflags( 4,11) ) THEN
1882  ALLOCATE ( wadats(imod)%XPSW(nxxx,0:noswll), stat=istat )
1883  check_alloc_status( istat )
1884  ELSE
1885  ALLOCATE ( wadats(imod)%XPSW(1,1), stat=istat )
1886  check_alloc_status( istat )
1887  END IF
1888  !
1889  IF ( outflags( 4,12) ) THEN
1890  ALLOCATE ( wadats(imod)%XPTM1(nxxx,0:noswll), stat=istat )
1891  check_alloc_status( istat )
1892  ELSE
1893  ALLOCATE ( wadats(imod)%XPTM1(1,1), stat=istat )
1894  check_alloc_status( istat )
1895  END IF
1896  !
1897  IF ( outflags( 4,13) ) THEN
1898  ALLOCATE ( wadats(imod)%XPT1(nxxx,0:noswll), stat=istat )
1899  check_alloc_status( istat )
1900  ELSE
1901  ALLOCATE ( wadats(imod)%XPT1(1,1), stat=istat )
1902  check_alloc_status( istat )
1903  END IF
1904  !
1905  IF ( outflags( 4,14) ) THEN
1906  ALLOCATE ( wadats(imod)%XPT2(nxxx,0:noswll), stat=istat )
1907  check_alloc_status( istat )
1908  ELSE
1909  ALLOCATE ( wadats(imod)%XPT2(1,1), stat=istat )
1910  check_alloc_status( istat )
1911  END IF
1912  !
1913  IF ( outflags( 4,15) ) THEN
1914  ALLOCATE ( wadats(imod)%XPEP(nxxx,0:noswll), stat=istat )
1915  check_alloc_status( istat )
1916  ELSE
1917  ALLOCATE ( wadats(imod)%XPEP(1,1), stat=istat )
1918  check_alloc_status( istat )
1919  END IF
1920  !
1921  IF ( outflags( 4,16) ) THEN
1922  ALLOCATE ( wadats(imod)%XPWST(nxxx), stat=istat )
1923  check_alloc_status( istat )
1924  ELSE
1925  ALLOCATE ( wadats(imod)%XPWST(1), stat=istat )
1926  check_alloc_status( istat )
1927  END IF
1928  !
1929  IF ( outflags( 4,17) ) THEN
1930  ALLOCATE ( wadats(imod)%XPNR(nxxx), stat=istat )
1931  check_alloc_status( istat )
1932  ELSE
1933  ALLOCATE ( wadats(imod)%XPNR(1), stat=istat )
1934  check_alloc_status( istat )
1935  END IF
1936  !
1937  wadats(imod)%XPHS = undef
1938  wadats(imod)%XPTP = undef
1939  wadats(imod)%XPLP = undef
1940  wadats(imod)%XPDIR = undef
1941  wadats(imod)%XPSI = undef
1942  wadats(imod)%XPWS = undef
1943  wadats(imod)%XPWST = undef
1944  wadats(imod)%XPNR = undef
1945  wadats(imod)%XPTHP0 = undef
1946  wadats(imod)%XPQP = undef
1947  wadats(imod)%XPPE = undef
1948  wadats(imod)%XPGW = undef
1949  wadats(imod)%XPSW = undef
1950  wadats(imod)%XPTM1 = undef
1951  wadats(imod)%XPT1 = undef
1952  wadats(imod)%XPT2 = undef
1953  wadats(imod)%XPEP = undef
1954  !
1955  IF ( outflags( 5, 2) ) THEN
1956  ALLOCATE ( wadats(imod)%XCHARN(nxxx), stat=istat )
1957  check_alloc_status( istat )
1958  ELSE
1959  ALLOCATE ( wadats(imod)%XCHARN(1), stat=istat )
1960  check_alloc_status( istat )
1961  END IF
1962  !
1963  IF ( outflags( 5, 3) ) THEN
1964  ALLOCATE ( wadats(imod)%XCGE(nxxx), stat=istat )
1965  check_alloc_status( istat )
1966  ELSE
1967  ALLOCATE ( wadats(imod)%XCGE(1), stat=istat )
1968  check_alloc_status( istat )
1969  END IF
1970  !
1971  IF ( outflags( 5, 4) ) THEN
1972  ALLOCATE ( wadats(imod)%XPHIAW(nxxx), stat=istat )
1973  check_alloc_status( istat )
1974  ELSE
1975  ALLOCATE ( wadats(imod)%XPHIAW(1), stat=istat )
1976  check_alloc_status( istat )
1977  END IF
1978  !
1979  IF ( outflags( 5, 5) ) THEN
1980  ALLOCATE ( wadats(imod)%XTAUWIX(nxxx), stat=istat )
1981  check_alloc_status( istat )
1982  ALLOCATE ( wadats(imod)%XTAUWIY(nxxx), stat=istat )
1983  check_alloc_status( istat )
1984  ELSE
1985  ALLOCATE ( wadats(imod)%XTAUWIX(1), stat=istat )
1986  check_alloc_status( istat )
1987  ALLOCATE ( wadats(imod)%XTAUWIY(1) )
1988  check_alloc_status( istat )
1989  END IF
1990  !
1991  IF ( outflags( 5, 6) ) THEN
1992  ALLOCATE ( wadats(imod)%XTAUWNX(nxxx), stat=istat )
1993  check_alloc_status( istat )
1994  ALLOCATE ( wadats(imod)%XTAUWNY(nxxx), stat=istat )
1995  check_alloc_status( istat )
1996  ELSE
1997  ALLOCATE ( wadats(imod)%XTAUWNX(1), stat=istat )
1998  check_alloc_status( istat )
1999  ALLOCATE ( wadats(imod)%XTAUWNY(1), stat=istat )
2000  check_alloc_status( istat )
2001  END IF
2002  !
2003  IF ( outflags( 5, 7) .OR. outflags( 5, 8) .OR. &
2004  outflags( 5, 9) .OR. outflags( 5,10)) THEN
2005  ALLOCATE ( wadats(imod)%XWHITECAP(nxxx,4), stat=istat )
2006  check_alloc_status( istat )
2007  ELSE
2008  ALLOCATE ( wadats(imod)%XWHITECAP(1,4), stat=istat )
2009  check_alloc_status( istat )
2010  END IF
2011  !
2012  IF ( outflags( 5, 11) ) THEN
2013  ALLOCATE ( wadats(imod)%XTWS(nxxx), stat=istat )
2014  check_alloc_status( istat )
2015  ELSE
2016  ALLOCATE ( wadats(imod)%XTWS(1), stat=istat )
2017  check_alloc_status( istat )
2018  END IF
2019  !
2020  wadats(imod)%XCHARN = undef
2021  wadats(imod)%XTWS = undef
2022  wadats(imod)%XCGE = undef
2023  wadats(imod)%XPHIAW = undef
2024  wadats(imod)%XTAUWIX = undef
2025  wadats(imod)%XTAUWIY = undef
2026  wadats(imod)%XTAUWNX = undef
2027  wadats(imod)%XTAUWNY = undef
2028  wadats(imod)%XWHITECAP = undef
2029  !
2030  IF ( outflags( 6, 1) ) THEN
2031  ALLOCATE ( wadats(imod)%XSXX(nxxx), stat=istat )
2032  check_alloc_status( istat )
2033  ALLOCATE ( wadats(imod)%XSYY(nxxx), stat=istat )
2034  check_alloc_status( istat )
2035  ALLOCATE ( wadats(imod)%XSXY(nxxx), stat=istat )
2036  check_alloc_status( istat )
2037  ELSE
2038  ALLOCATE ( wadats(imod)%XSXX(1), stat=istat )
2039  check_alloc_status( istat )
2040  ALLOCATE ( wadats(imod)%XSYY(1), stat=istat )
2041  check_alloc_status( istat )
2042  ALLOCATE ( wadats(imod)%XSXY(1), stat=istat )
2043  check_alloc_status( istat )
2044  END IF
2045  !
2046  IF ( outflags( 6, 2) ) THEN
2047  ALLOCATE ( wadats(imod)%XTAUOX(nxxx), stat=istat )
2048  check_alloc_status( istat )
2049  ALLOCATE ( wadats(imod)%XTAUOY(nxxx), stat=istat )
2050  check_alloc_status( istat )
2051  ELSE
2052  ALLOCATE ( wadats(imod)%XTAUOX(1), stat=istat )
2053  check_alloc_status( istat )
2054  ALLOCATE ( wadats(imod)%XTAUOY(1), stat=istat )
2055  check_alloc_status( istat )
2056  END IF
2057  !
2058  IF ( outflags( 6, 3) ) THEN
2059  ALLOCATE ( wadats(imod)%XBHD(nxxx), stat=istat )
2060  check_alloc_status( istat )
2061  ELSE
2062  ALLOCATE ( wadats(imod)%XBHD(1), stat=istat )
2063  check_alloc_status( istat )
2064  END IF
2065  !
2066  IF ( outflags( 6, 4) ) THEN
2067  ALLOCATE ( wadats(imod)%XPHIOC(nxxx), stat=istat )
2068  check_alloc_status( istat )
2069  ELSE
2070  ALLOCATE ( wadats(imod)%XPHIOC(1), stat=istat )
2071  check_alloc_status( istat )
2072  END IF
2073  !
2074  IF ( outflags( 6, 5) ) THEN
2075  ALLOCATE ( wadats(imod)%XTUSX(nxxx), stat=istat )
2076  check_alloc_status( istat )
2077  ALLOCATE ( wadats(imod)%XTUSY(nxxx), stat=istat )
2078  check_alloc_status( istat )
2079  ELSE
2080  ALLOCATE ( wadats(imod)%XTUSX(1), stat=istat )
2081  check_alloc_status( istat )
2082  ALLOCATE ( wadats(imod)%XTUSY(1), stat=istat )
2083  check_alloc_status( istat )
2084  END IF
2085  !
2086  IF ( outflags( 6, 6) ) THEN
2087  ALLOCATE ( wadats(imod)%XUSSX(nxxx), stat=istat )
2088  check_alloc_status( istat )
2089  ALLOCATE ( wadats(imod)%XUSSY(nxxx), stat=istat )
2090  check_alloc_status( istat )
2091  ELSE
2092  ALLOCATE ( wadats(imod)%XUSSX(1), stat=istat )
2093  check_alloc_status( istat )
2094  ALLOCATE ( wadats(imod)%XUSSY(1), stat=istat )
2095  check_alloc_status( istat )
2096  END IF
2097  !
2098  IF ( outflags( 6, 7) ) THEN
2099  ALLOCATE ( wadats(imod)%XPRMS(nxxx), stat=istat )
2100  check_alloc_status( istat )
2101  ALLOCATE ( wadats(imod)%XTPMS(nxxx), stat=istat )
2102  check_alloc_status( istat )
2103  ELSE
2104  ALLOCATE ( wadats(imod)%XPRMS(1), stat=istat )
2105  check_alloc_status( istat )
2106  ALLOCATE ( wadats(imod)%XTPMS(1), stat=istat )
2107  check_alloc_status( istat )
2108  END IF
2109  !
2110  IF ( outflags( 6, 8) ) THEN
2111  ALLOCATE ( wadats(imod)%XUS3D(nxxx,2*nk), stat=istat )
2112  check_alloc_status( istat )
2113  ELSE
2114  ALLOCATE ( wadats(imod)%XUS3D(1,1), stat=istat )
2115  check_alloc_status( istat )
2116  END IF
2117  !
2118  IF ( outflags( 6, 9) ) THEN
2119  ALLOCATE ( wadats(imod)%XP2SMS(nxxx,p2msf(2):p2msf(3)), stat=istat )
2120  check_alloc_status( istat )
2121  ELSE
2122  ALLOCATE ( wadats(imod)%XP2SMS(1,1), stat=istat )
2123  check_alloc_status( istat )
2124  END IF
2125  !
2126  IF ( outflags( 6,10) ) THEN
2127  ALLOCATE ( wadats(imod)%XTAUICE(nxxx,2), stat=istat )
2128  check_alloc_status( istat )
2129  ELSE
2130  ALLOCATE ( wadats(imod)%XTAUICE(1,2), stat=istat )
2131  check_alloc_status( istat )
2132  END IF
2133  !
2134  IF ( outflags( 6,11) ) THEN
2135  ALLOCATE ( wadats(imod)%XPHICE(nxxx), stat=istat )
2136  check_alloc_status( istat )
2137  ELSE
2138  ALLOCATE ( wadats(imod)%XPHICE(1), stat=istat )
2139  check_alloc_status( istat )
2140  END IF
2141  !
2142  IF ( outflags( 6, 12) ) THEN
2143  ALLOCATE ( wadats(imod)%XUSSP(nxxx,2*nk), stat=istat )
2144  check_alloc_status( istat )
2145  ELSE
2146  ALLOCATE ( wadats(imod)%XUSSP(1,1), stat=istat )
2147  check_alloc_status( istat )
2148  END IF
2149  !
2150  IF ( outflags( 6, 13) ) THEN
2151  ALLOCATE ( wadats(imod)%XTAUOCX(nxxx), stat=istat )
2152  check_alloc_status( istat )
2153  ALLOCATE ( wadats(imod)%XTAUOCY(nxxx), stat=istat )
2154  check_alloc_status( istat )
2155  ELSE
2156  ALLOCATE ( wadats(imod)%XTAUOCX(1), stat=istat )
2157  check_alloc_status( istat )
2158  ALLOCATE ( wadats(imod)%XTAUOCY(1), stat=istat )
2159  check_alloc_status( istat )
2160  END IF
2161  !
2162  wadats(imod)%XSXX = undef
2163  wadats(imod)%XSYY = undef
2164  wadats(imod)%XSXY = undef
2165  wadats(imod)%XTAUOX = undef
2166  wadats(imod)%XTAUOY = undef
2167  wadats(imod)%XBHD = undef
2168  wadats(imod)%XPHIOC = undef
2169  wadats(imod)%XTUSX = undef
2170  wadats(imod)%XTUSY = undef
2171  wadats(imod)%XUSSX = undef
2172  wadats(imod)%XUSSY = undef
2173  wadats(imod)%XPRMS = undef
2174  wadats(imod)%XTPMS = undef
2175  wadats(imod)%XUS3D = undef
2176  wadats(imod)%XP2SMS = undef
2177  wadats(imod)%XPHICE = undef
2178  wadats(imod)%XTAUICE = undef
2179  wadats(imod)%XUSSP = undef
2180  wadats(imod)%XTAUOCX = undef
2181  wadats(imod)%XTAUOCY = undef
2182  !
2183  IF ( outflags( 7, 1) ) THEN
2184  ALLOCATE ( wadats(imod)%XABA(nxxx), stat=istat )
2185  check_alloc_status( istat )
2186  ALLOCATE ( wadats(imod)%XABD(nxxx), stat=istat )
2187  check_alloc_status( istat )
2188  ELSE
2189  ALLOCATE ( wadats(imod)%XABA(1), stat=istat )
2190  check_alloc_status( istat )
2191  ALLOCATE ( wadats(imod)%XABD(1), stat=istat )
2192  check_alloc_status( istat )
2193  END IF
2194  !
2195  IF ( outflags( 7, 2) ) THEN
2196  ALLOCATE ( wadats(imod)%XUBA(nxxx), stat=istat )
2197  check_alloc_status( istat )
2198  ALLOCATE ( wadats(imod)%XUBD(nxxx), stat=istat )
2199  check_alloc_status( istat )
2200  ELSE
2201  ALLOCATE ( wadats(imod)%XUBA(1), stat=istat )
2202  check_alloc_status( istat )
2203  ALLOCATE ( wadats(imod)%XUBD(1), stat=istat )
2204  check_alloc_status( istat )
2205  END IF
2206  !
2207  IF ( outflags( 7, 3) ) THEN
2208  ALLOCATE ( wadats(imod)%XBEDFORMS(nxxx,3), stat=istat )
2209  check_alloc_status( istat )
2210  ELSE
2211  ALLOCATE ( wadats(imod)%XBEDFORMS(1,3), stat=istat )
2212  check_alloc_status( istat )
2213  END IF
2214  !
2215  IF ( outflags( 7, 4) ) THEN
2216  ALLOCATE ( wadats(imod)%XPHIBBL(nxxx), stat=istat )
2217  check_alloc_status( istat )
2218  ELSE
2219  ALLOCATE ( wadats(imod)%XPHIBBL(1), stat=istat )
2220  check_alloc_status( istat )
2221  END IF
2222  !
2223  IF ( outflags( 7, 5) ) THEN
2224  ALLOCATE ( wadats(imod)%XTAUBBL(nxxx,2), stat=istat )
2225  check_alloc_status( istat )
2226  ELSE
2227  ALLOCATE ( wadats(imod)%XTAUBBL(1,2), stat=istat )
2228  check_alloc_status( istat )
2229  END IF
2230  !
2231  wadats(imod)%XABA = undef
2232  wadats(imod)%XABD = undef
2233  wadats(imod)%XUBA = undef
2234  wadats(imod)%XUBD = undef
2235  wadats(imod)%XBEDFORMS = undef
2236  wadats(imod)%XPHIBBL = undef
2237  wadats(imod)%XTAUBBL = undef
2238  !
2239  IF ( outflags( 8, 1) ) THEN
2240  ALLOCATE ( wadats(imod)%XMSSX(nxxx), stat=istat )
2241  check_alloc_status( istat )
2242  ALLOCATE ( wadats(imod)%XMSSY(nxxx), stat=istat )
2243  check_alloc_status( istat )
2244  ELSE
2245  ALLOCATE ( wadats(imod)%XMSSX(1), stat=istat )
2246  check_alloc_status( istat )
2247  ALLOCATE ( wadats(imod)%XMSSY(1), stat=istat )
2248  check_alloc_status( istat )
2249  END IF
2250  !
2251  IF ( outflags( 8, 2) ) THEN
2252  ALLOCATE ( wadats(imod)%XMSCX(nxxx), stat=istat )
2253  check_alloc_status( istat )
2254  ALLOCATE ( wadats(imod)%XMSCY(nxxx), stat=istat )
2255  check_alloc_status( istat )
2256  ELSE
2257  ALLOCATE ( wadats(imod)%XMSCX(1), stat=istat )
2258  check_alloc_status( istat )
2259  ALLOCATE ( wadats(imod)%XMSCY(1), stat=istat )
2260  check_alloc_status( istat )
2261  END IF
2262  !
2263  IF ( outflags( 8, 3) ) THEN
2264  ALLOCATE ( wadats(imod)%XMSSD(nxxx), stat=istat )
2265  check_alloc_status( istat )
2266  ELSE
2267  ALLOCATE ( wadats(imod)%XMSSD(1), stat=istat )
2268  check_alloc_status( istat )
2269  END IF
2270  !
2271  IF ( outflags( 8, 4) ) THEN
2272  ALLOCATE ( wadats(imod)%XMSCD(nxxx), stat=istat )
2273  check_alloc_status( istat )
2274  ELSE
2275  ALLOCATE ( wadats(imod)%XMSCD(1), stat=istat )
2276  check_alloc_status( istat )
2277  END IF
2278  !
2279  IF ( outflags( 8, 5) ) THEN
2280  ALLOCATE ( wadats(imod)%XQP(nxxx) )
2281  ELSE
2282  ALLOCATE ( wadats(imod)%XQP(1) )
2283  END IF
2284  !
2285  IF ( outflags( 8, 6) ) THEN
2286  ALLOCATE ( wadats(imod)%XQKK(nxxx) )
2287  ELSE
2288  ALLOCATE ( wadats(imod)%XQKK(1) )
2289  END IF
2290  !
2291  IF ( outflags( 8, 7) ) THEN
2292  ALLOCATE ( wadats(imod)%XSKEW(nxxx) )
2293  ELSE
2294  ALLOCATE ( wadats(imod)%XSKEW(1) )
2295  END IF
2296  !
2297  IF ( outflags( 8, 8) ) THEN
2298  ALLOCATE ( wadats(imod)%XEMBIA1(nxxx) )
2299  ELSE
2300  ALLOCATE ( wadats(imod)%XEMBIA1(1) )
2301  END IF
2302  !
2303  IF ( outflags( 8, 9) ) THEN
2304  ALLOCATE ( wadats(imod)%XEMBIA2(nxxx) )
2305  ELSE
2306  ALLOCATE ( wadats(imod)%XEMBIA2(1) )
2307  END IF
2308  !
2309  wadats(imod)%XMSSX = undef
2310  wadats(imod)%XMSSY = undef
2311  wadats(imod)%XMSSD = undef
2312  wadats(imod)%XMSCX = undef
2313  wadats(imod)%XMSCY = undef
2314  wadats(imod)%XMSCD = undef
2315  wadats(imod)%XQP(1) = undef
2316  wadats(imod)%XQKK = undef
2317  wadats(imod)%XSKEW = undef
2318  wadats(imod)%XEMBIA1 = undef
2319  wadats(imod)%XEMBIA2 = undef
2320  !
2321  IF ( outflags( 9, 1) ) THEN
2322  ALLOCATE ( wadats(imod)%XDTDYN(nxxx), stat=istat )
2323  check_alloc_status( istat )
2324  ELSE
2325  ALLOCATE ( wadats(imod)%XDTDYN(1), stat=istat )
2326  check_alloc_status( istat )
2327  END IF
2328  !
2329  IF ( outflags( 9, 2) ) THEN
2330  ALLOCATE ( wadats(imod)%XFCUT(nxxx), stat=istat )
2331  check_alloc_status( istat )
2332  ELSE
2333  ALLOCATE ( wadats(imod)%XFCUT(1), stat=istat )
2334  check_alloc_status( istat )
2335  END IF
2336  !
2337  IF ( outflags( 9, 3) ) THEN
2338  ALLOCATE ( wadats(imod)%XCFLXYMAX(nxxx), stat=istat )
2339  check_alloc_status( istat )
2340  ELSE
2341  ALLOCATE ( wadats(imod)%XCFLXYMAX(1), stat=istat )
2342  check_alloc_status( istat )
2343  END IF
2344  !
2345  IF ( outflags( 9, 4) ) THEN
2346  ALLOCATE ( wadats(imod)%XCFLTHMAX(nxxx), stat=istat )
2347  check_alloc_status( istat )
2348  ELSE
2349  ALLOCATE ( wadats(imod)%XCFLTHMAX(1), stat=istat )
2350  check_alloc_status( istat )
2351  END IF
2352  !
2353  IF ( outflags( 9, 5) ) THEN
2354  ALLOCATE ( wadats(imod)%XCFLKMAX(nxxx), stat=istat )
2355  check_alloc_status( istat )
2356  ELSE
2357  ALLOCATE ( wadats(imod)%XCFLKMAX(1), stat=istat )
2358  check_alloc_status( istat )
2359  END IF
2360  !
2361  wadats(imod)%XDTDYN = undef
2362  wadats(imod)%XFCUT = undef
2363  wadats(imod)%XCFLXYMAX = undef
2364  wadats(imod)%XCFLTHMAX = undef
2365  wadats(imod)%XCFLKMAX = undef
2366  !
2367  DO i=1, noextr
2368  IF ( outflags(10, i) ) THEN
2369  ALLOCATE ( wadats(imod)%XUSERO(nxxx,i), stat=istat )
2370  check_alloc_status( istat )
2371  ELSE
2372  ALLOCATE ( wadats(imod)%XUSERO(1,i), stat=istat )
2373  check_alloc_status( istat )
2374  END IF
2375  END DO
2376  !
2377  wadats(imod)%XUSERO = undef
2378  !
2379  wadats(imod)%AINIT2 = .true.
2380  !
2381 #ifdef W3_T
2382  WRITE (ndst,9001)
2383  WRITE (ndst,9001)
2384 #endif
2385  !
2386  ! -------------------------------------------------------------------- /
2387  ! 5. Restore previous grid setting if necessary
2388  !
2389  IF ( jgrid .NE. imod ) CALL w3setg ( jgrid, ndse, ndst )
2390 
2391  call print_memcheck(memunit, 'memcheck_____:'//' W3XDMA')
2392  !
2393  RETURN
2394  !
2395  ! Formats
2396  !
2397 1001 FORMAT (/' *** ERROR W3XDMA : GRIDS NOT INITIALIZED *** '/ &
2398  ' RUN W3NMOD FIRST '/)
2399 1002 FORMAT (/' *** ERROR W3XDMA : ILLEGAL MODEL NUMBER *** '/ &
2400  ' IMOD = ',i10/ &
2401  ' NADATA = ',i10/)
2402 1003 FORMAT (/' *** ERROR W3XDMA : ARRAY(S) ALREADY ALLOCATED *** ')
2403  !
2404 #ifdef W3_T
2405 9000 FORMAT (' TEST W3XDMA : MODEL ',i4)
2406 9001 FORMAT (' TEST W3XDMA : ARRAYS ALLOCATED')
2407 9002 FORMAT (' TEST W3XDMA : POINTERS RESET')
2408 9003 FORMAT (' TEST W3XDMA : DIMENSIONS STORED')
2409 #endif
2410  !/
2411  !/ End of W3XDMA ----------------------------------------------------- /
2412  !/
2413  END SUBROUTINE w3xdma
2414  !/ ------------------------------------------------------------------- /
2430  SUBROUTINE w3dmnl ( IMOD, NDSE, NDST, NSP, NSPX )
2431  !/
2432  !/ +-----------------------------------+
2433  !/ | WAVEWATCH III NOAA/NCEP |
2434  !/ | H. L. Tolman |
2435  !/ | FORTRAN 90 |
2436  !/ | Last update : 10-Dec-2014 |
2437  !/ +-----------------------------------+
2438  !/
2439  !/ 24-Dec-2004 : Origination. ( version 3.06 )
2440  !/ 04-Oct-2006 : Add filter to array pointers. ( version 3.10 )
2441  !/ 10-Dec-2014 : Add checks for allocate status ( version 5.04 )
2442  !/
2443  ! 1. Purpose :
2444  !
2445  ! Initialize an individual data grid at the proper dimensions (DIA).
2446  !
2447  ! 2. Method :
2448  !
2449  ! Allocate directly into the structure array. Note that
2450  ! this cannot be done through the pointer alias!
2451  !
2452  ! 3. Parameters :
2453  !
2454  ! Parameter list
2455  ! ----------------------------------------------------------------
2456  ! IMOD Int. I Model number to point to.
2457  ! NDSE Int. I Error output unit number.
2458  ! NDST Int. I Test output unit number.
2459  ! NSP(X) Int. I Array dimensions.
2460  ! ----------------------------------------------------------------
2461  !
2462  ! 4. Subroutines used :
2463  !
2464  ! See module documentation.
2465  !
2466  ! 5. Called by :
2467  !
2468  ! Name Type Module Description
2469  ! ----------------------------------------------------------------
2470  ! INSNL1 Subr. W3SNL1MD Traditional DIA approach to Snl.
2471  ! ----------------------------------------------------------------
2472  !
2473  ! 6. Error messages :
2474  !
2475  ! - Check on input parameters.
2476  ! - Check on previous allocation.
2477  !
2478  ! 7. Remarks :
2479  !
2480  ! - W3SETA needs to be called after allocation to point to
2481  ! proper allocated arrays.
2482  !
2483  ! 8. Structure :
2484  !
2485  ! See source code.
2486  !
2487  ! 9. Switches :
2488  !
2489  ! !/S Enable subroutine tracing.
2490  ! !/T Enable test output
2491  !
2492  ! 10. Source code :
2493  !
2494  !/ ------------------------------------------------------------------- /
2495  USE w3gdatmd, ONLY: ngrids, igrid, nk, nx, ny, nsea, nseal, &
2496  nspec, nth, gtype, ungtype
2497  USE w3odatmd, ONLY: naproc
2498  USE w3servmd, ONLY: extcde
2499 #ifdef W3_S
2500  USE w3servmd, ONLY: strace
2501 #endif
2502  !
2503  !/
2504  !/ ------------------------------------------------------------------- /
2505  !/ Parameter list
2506  !/
2507  INTEGER, INTENT(IN) :: IMOD, NDSE, NDST, NSP, NSPX
2508  !/
2509  !/ ------------------------------------------------------------------- /
2510  !/ Local parameters
2511  !/
2512 #ifdef W3_S
2513  INTEGER, SAVE :: IENT = 0
2514  CALL strace (ient, 'W3DMNL')
2515 #endif
2516  !
2517  ! -------------------------------------------------------------------- /
2518  ! 1. Test input and module status
2519  !
2520  IF ( ngrids .EQ. -1 ) THEN
2521  WRITE (ndse,1001)
2522  CALL extcde (1)
2523  END IF
2524  !
2525  IF ( imod.LT.1 .OR. imod.GT.nadata ) THEN
2526  WRITE (ndse,1002) imod, nadata
2527  CALL extcde (2)
2528  END IF
2529  !
2530 #ifdef W3_NL1
2531  IF ( wadats(imod)%NLINIT ) THEN
2532  WRITE (ndse,1003)
2533  CALL extcde (3)
2534  END IF
2535 #endif
2536  !
2537 #ifdef W3_T
2538  WRITE (ndst,9000) imod
2539 #endif
2540  !
2541  ! -------------------------------------------------------------------- /
2542  ! 2. Allocate arrays
2543  !
2544 #ifdef W3_NL1
2545  ALLOCATE ( wadats(imod)%IP11(nspx), &
2546  wadats(imod)%IP12(nspx), &
2547  wadats(imod)%IP13(nspx), &
2548  wadats(imod)%IP14(nspx), &
2549  wadats(imod)%IM11(nspx), &
2550  wadats(imod)%IM12(nspx), &
2551  wadats(imod)%IM13(nspx), &
2552  wadats(imod)%IM14(nspx), &
2553  wadats(imod)%IP21(nspx), &
2554  wadats(imod)%IP22(nspx), &
2555  wadats(imod)%IP23(nspx), &
2556  wadats(imod)%IP24(nspx), &
2557  wadats(imod)%IM21(nspx), &
2558  wadats(imod)%IM22(nspx), &
2559  wadats(imod)%IM23(nspx), &
2560  wadats(imod)%IM24(nspx), &
2561  wadats(imod)%IC11(nsp) , &
2562  wadats(imod)%IC12(nsp) , &
2563  wadats(imod)%IC21(nsp) , &
2564  wadats(imod)%IC22(nsp) , &
2565  wadats(imod)%IC31(nsp) , &
2566  wadats(imod)%IC32(nsp) , &
2567  wadats(imod)%IC41(nsp) , &
2568  wadats(imod)%IC42(nsp) , &
2569  wadats(imod)%IC51(nsp) , &
2570  wadats(imod)%IC52(nsp) , &
2571  wadats(imod)%IC61(nsp) , &
2572  wadats(imod)%IC62(nsp) , &
2573  wadats(imod)%IC71(nsp) , &
2574  wadats(imod)%IC72(nsp) , &
2575  wadats(imod)%IC81(nsp) , &
2576  wadats(imod)%IC82(nsp) , &
2577  wadats(imod)%AF11(nspx), &
2578  stat=istat )
2579  check_alloc_status( istat )
2580  wadats(imod)%NLINIT = .true.
2581 #endif
2582  !
2583 #ifdef W3_T
2584  WRITE (ndst,9001)
2585 #endif
2586  !
2587  ! -------------------------------------------------------------------- /
2588  ! 3. Point to allocated arrays
2589  !
2590  CALL w3seta ( imod, ndse, ndst )
2591  !
2592 #ifdef W3_T
2593  WRITE (ndst,9002)
2594 #endif
2595  !
2596  ! -------------------------------------------------------------------- /
2597  ! 4. Update counters in grid
2598  !
2599 #ifdef W3_NL1
2600  nspecx = nspx
2601 #endif
2602  !
2603 #ifdef W3_T
2604  WRITE (ndst,9003)
2605 #endif
2606  !
2607  RETURN
2608  !
2609  ! Formats
2610  !
2611 1001 FORMAT (/' *** ERROR W3DMNL : GRIDS NOT INITIALIZED *** '/ &
2612  ' RUN W3NMOD FIRST '/)
2613 1002 FORMAT (/' *** ERROR W3DMNL : ILLEGAL MODEL NUMBER *** '/ &
2614  ' IMOD = ',i10/ &
2615  ' NADATA = ',i10/)
2616 #ifdef W3_NL1
2617 1003 FORMAT (/' *** ERROR W3DMNL : ARRAY(S) ALREADY ALLOCATED *** ')
2618 #endif
2619  !
2620 #ifdef W3_T
2621 9000 FORMAT (' TEST W3DMNL : MODEL ',i4)
2622 9001 FORMAT (' TEST W3DMNL : ARRAYS ALLOCATED')
2623 9002 FORMAT (' TEST W3DMNL : POINTERS RESET')
2624 9003 FORMAT (' TEST W3DMNL : DIMENSIONS STORED')
2625 #endif
2626  !/
2627  !/ End of W3DMNL ----------------------------------------------------- /
2628  !/
2629  END SUBROUTINE w3dmnl
2630  !/ ------------------------------------------------------------------- /
2644  SUBROUTINE w3seta ( IMOD, NDSE, NDST )
2645  !/
2646  !/ +-----------------------------------+
2647  !/ | WAVEWATCH III NOAA/NCEP |
2648  !/ | H. L. Tolman |
2649  !/ | FORTRAN 90 |
2650  !/ | Last update : 22-Mar-2021 |
2651  !/ +-----------------------------------+
2652  !/
2653  !/ 28-Dec-2004 : Origination. ( version 3.06 )
2654  !/ 04-May-2005 : Adding MPI_COMM_WAVE. ( version 3.07 )
2655  !/ 20-Jul-2005 : Adding output fields. ( version 3.07 )
2656  !/ 09-Nov-2005 : Removing soft boundary option. ( version 3.08 )
2657  !/ 13-Jun-2006 : Splitting STORE in G/SSTORE. ( version 3.09 )
2658  !/ 04-Oct-2006 : Add filter to array pointers. ( version 3.10 )
2659  !/ 28_Mar-2007 : Add partitioned data arrays. ( version 3.11 )
2660  !/ Add aditional undefined arrays.
2661  !/ 22-Mar-2021 : Adds TAUA, WNMEAN, TAUOC parameters ( version 7.13 )
2662  !/
2663  ! 1. Purpose :
2664  !
2665  ! Select one of the WAVEWATCH III grids / models.
2666  !
2667  ! 2. Method :
2668  !
2669  ! Point pointers to the proper variables in the proper element of
2670  ! the GRIDS array.
2671  !
2672  ! 3. Parameters :
2673  !
2674  ! Parameter list
2675  ! ----------------------------------------------------------------
2676  ! IMOD Int. I Model number to point to.
2677  ! NDSE Int. I Error output unit number.
2678  ! NDST Int. I Test output unit number.
2679  ! ----------------------------------------------------------------
2680  !
2681  ! 4. Subroutines used :
2682  !
2683  ! See module documentation below.
2684  !
2685  ! 5. Called by :
2686  !
2687  ! Many subroutines in the WAVEWATCH system.
2688  !
2689  ! 6. Error messages :
2690  !
2691  ! Checks on parameter list IMOD.
2692  !
2693  ! 7. Remarks :
2694  !
2695  ! 8. Structure :
2696  !
2697  ! 9. Switches :
2698  !
2699  ! !/MPI Paralllel model environment.
2700  !
2701  ! !/PRn Propagation scheme selection.
2702  !
2703  ! !/S Enable subroutine tracing.
2704  ! !/T Enable test output
2705  !
2706  ! 10. Source code :
2707  !
2708  !/ ------------------------------------------------------------------- /
2709  !
2710  USE w3idatmd, ONLY: inputs
2711  USE w3gdatmd, ONLY: e3df, p2msf, us3df, usspf, gtype, ungtype
2712  !
2713  USE w3servmd, ONLY: extcde
2714 #ifdef W3_S
2715  USE w3servmd, ONLY: strace
2716 #endif
2717  !
2718  !/
2719  !/ ------------------------------------------------------------------- /
2720  !/ Parameter list
2721  !/
2722  INTEGER, INTENT(IN) :: IMOD, NDSE, NDST
2723  !/
2724  !/ ------------------------------------------------------------------- /
2725  !/ Local parameters
2726  !/
2727 #ifdef W3_S
2728  INTEGER, SAVE :: IENT = 0
2729  CALL strace (ient, 'W3SETA')
2730 #endif
2731  !
2732  ! -------------------------------------------------------------------- /
2733  ! 1. Test input and module status
2734  !
2735  IF ( nadata .EQ. -1 ) THEN
2736  WRITE (ndse,1001)
2737  CALL extcde (1)
2738  END IF
2739  !
2740  IF ( imod.LT.1 .OR. imod.GT.nadata ) THEN
2741  WRITE (ndse,1002) imod, nadata
2742  CALL extcde (2)
2743  END IF
2744  !
2745 #ifdef W3_T
2746  WRITE (ndst,9000) imod
2747 #endif
2748  !
2749  ! -------------------------------------------------------------------- /
2750  ! 2. Set model numbers
2751  !
2752  iadata = imod
2753  !
2754  ! -------------------------------------------------------------------- /
2755  ! 3. Set pointers
2756  !
2757  itime => wadats(imod)%ITIME
2758  ipass => wadats(imod)%IPASS
2759  idlast => wadats(imod)%IDLAST
2760  nsealm => wadats(imod)%NSEALM
2761  flcold => wadats(imod)%FLCOLD
2762  fliwnd => wadats(imod)%FLIWND
2763  ainit => wadats(imod)%AINIT
2764  ainit2 => wadats(imod)%AINIT2
2765  fl_all => wadats(imod)%FL_ALL
2766  !
2767 #ifdef W3_PR2
2768  nmx0 => wadats(imod)%NMX0
2769  nmx1 => wadats(imod)%NMX1
2770  nmx2 => wadats(imod)%NMX2
2771  nmy0 => wadats(imod)%NMY0
2772  nmy1 => wadats(imod)%NMY1
2773  nmy2 => wadats(imod)%NMY2
2774  nact => wadats(imod)%NACT
2775  nmxy => wadats(imod)%NMXY
2776 #endif
2777  !
2778 #ifdef W3_PR3
2779  nmx0 => wadats(imod)%NMX0
2780  nmx1 => wadats(imod)%NMX1
2781  nmx2 => wadats(imod)%NMX2
2782  nmy0 => wadats(imod)%NMY0
2783  nmy1 => wadats(imod)%NMY1
2784  nmy2 => wadats(imod)%NMY2
2785  nact => wadats(imod)%NACT
2786  ncent => wadats(imod)%NCENT
2787 #endif
2788  !
2789 #ifdef W3_NL1
2790  nfr => wadats(imod)%NFR
2791  nfrhgh => wadats(imod)%NFRHGH
2792  nfrchg => wadats(imod)%NFRCHG
2793  nspecx => wadats(imod)%NSPECX
2794  nspecy => wadats(imod)%NSPECY
2795  dal1 => wadats(imod)%DAL1
2796  dal2 => wadats(imod)%DAL2
2797  dal3 => wadats(imod)%DAL3
2798  awg1 => wadats(imod)%AWG1
2799  awg2 => wadats(imod)%AWG2
2800  awg3 => wadats(imod)%AWG3
2801  awg4 => wadats(imod)%AWG4
2802  awg5 => wadats(imod)%AWG5
2803  awg6 => wadats(imod)%AWG6
2804  awg7 => wadats(imod)%AWG7
2805  awg8 => wadats(imod)%AWG8
2806  swg1 => wadats(imod)%SWG1
2807  swg2 => wadats(imod)%SWG2
2808  swg3 => wadats(imod)%SWG3
2809  swg4 => wadats(imod)%SWG4
2810  swg5 => wadats(imod)%SWG5
2811  swg6 => wadats(imod)%SWG6
2812  swg7 => wadats(imod)%SWG7
2813  swg8 => wadats(imod)%SWG8
2814  nlinit => wadats(imod)%NLINIT
2815 #endif
2816  !
2817 #ifdef W3_MPI
2818  mpi_comm_wave => wadats(imod)%MPI_COMM_WAVE
2819  mpi_comm_wcmp => wadats(imod)%MPI_COMM_WCMP
2820  ww3_field_vec => wadats(imod)%WW3_FIELD_VEC
2821  ww3_spec_vec => wadats(imod)%WW3_SPEC_VEC
2822  nrqsg1 => wadats(imod)%NRQSG1
2823  nrqsg2 => wadats(imod)%NRQSG2
2824  ibfloc => wadats(imod)%IBFLOC
2825  isploc => wadats(imod)%ISPLOC
2826  nsploc => wadats(imod)%NSPLOC
2827  bstat => wadats(imod)%BSTAT
2828  bispl => wadats(imod)%BISPL
2829 #endif
2830  !
2831  IF ( ainit ) THEN
2832  !
2833  dw => wadats(imod)%DW
2834  ua => wadats(imod)%UA
2835  ud => wadats(imod)%UD
2836  u10 => wadats(imod)%U10
2837  u10d => wadats(imod)%U10D
2838  as => wadats(imod)%AS
2839  cx => wadats(imod)%CX
2840  cy => wadats(imod)%CY
2841  taua => wadats(imod)%TAUA
2842  tauadir=> wadats(imod)%TAUADIR
2843  !
2844  hs => wadats(imod)%HS
2845  wlm => wadats(imod)%WLM
2846  t02 => wadats(imod)%T02
2847  t0m1 => wadats(imod)%T0M1
2848  t01 => wadats(imod)%T01
2849  fp0 => wadats(imod)%FP0
2850  thm => wadats(imod)%THM
2851  ths => wadats(imod)%THS
2852  thp0 => wadats(imod)%THP0
2853  hsig => wadats(imod)%HSIG
2854  stmaxe => wadats(imod)%STMAXE
2855  stmaxd => wadats(imod)%STMAXD
2856  hmaxe => wadats(imod)%HMAXE
2857  hmaxd => wadats(imod)%HMAXD
2858  hcmaxe => wadats(imod)%HCMAXE
2859  hcmaxd => wadats(imod)%HCMAXD
2860  qp => wadats(imod)%QP
2861  wbt => wadats(imod)%WBT
2862  wnmean => wadats(imod)%WNMEAN
2863  !
2864  ef => wadats(imod)%EF
2865  th1m => wadats(imod)%TH1M
2866  sth1m => wadats(imod)%STH1M
2867  th2m => wadats(imod)%TH2M
2868  sth2m => wadats(imod)%STH2M
2869  !
2870  phs => wadats(imod)%PHS
2871  ptp => wadats(imod)%PTP
2872  plp => wadats(imod)%PLP
2873  pdir => wadats(imod)%PDIR
2874  psi => wadats(imod)%PSI
2875  pws => wadats(imod)%PWS
2876  pwst => wadats(imod)%PWST
2877  pnr => wadats(imod)%PNR
2878  pthp0 => wadats(imod)%PTHP0
2879  pqp => wadats(imod)%PQP
2880  ppe => wadats(imod)%PPE
2881  pgw => wadats(imod)%PGW
2882  psw => wadats(imod)%PSW
2883  ptm1 => wadats(imod)%PTM1
2884  pt1 => wadats(imod)%PT1
2885  pt2 => wadats(imod)%PT2
2886  pep => wadats(imod)%PEP
2887  !
2888  charn => wadats(imod)%CHARN
2889  tws => wadats(imod)%TWS
2890  cge => wadats(imod)%CGE
2891  phiaw => wadats(imod)%PHIAW
2892  tauwix => wadats(imod)%TAUWIX
2893  tauwiy => wadats(imod)%TAUWIY
2894  tauwnx => wadats(imod)%TAUWNX
2895  tauwny => wadats(imod)%TAUWNY
2896  whitecap => wadats(imod)%WHITECAP
2897  !
2898  sxx => wadats(imod)%SXX
2899  syy => wadats(imod)%SYY
2900  sxy => wadats(imod)%SXY
2901  tauox => wadats(imod)%TAUOX
2902  tauoy => wadats(imod)%TAUOY
2903  bhd => wadats(imod)%BHD
2904  phioc => wadats(imod)%PHIOC
2905  tusx => wadats(imod)%TUSX
2906  tusy => wadats(imod)%TUSY
2907  ussx => wadats(imod)%USSX
2908  ussy => wadats(imod)%USSY
2909  prms => wadats(imod)%PRMS
2910  tpms => wadats(imod)%TPMS
2911  p2sms => wadats(imod)%P2SMS
2912  us3d => wadats(imod)%US3D
2913  phice => wadats(imod)%PHICE
2914  tauice => wadats(imod)%TAUICE
2915  ussp => wadats(imod)%USSP
2916  tauocx => wadats(imod)%TAUOCX
2917  tauocy => wadats(imod)%TAUOCY
2918  !
2919  aba => wadats(imod)%ABA
2920  abd => wadats(imod)%ABD
2921  uba => wadats(imod)%UBA
2922  ubd => wadats(imod)%UBD
2923  bedforms=> wadats(imod)%BEDFORMS
2924  phibbl => wadats(imod)%PHIBBL
2925  taubbl => wadats(imod)%TAUBBL
2926  !
2927  mssx => wadats(imod)%MSSX
2928  mssy => wadats(imod)%MSSY
2929  mssd => wadats(imod)%MSSD
2930  mscx => wadats(imod)%MSCX
2931  mscy => wadats(imod)%MSCY
2932  mscd => wadats(imod)%MSCD
2933  qkk => wadats(imod)%QKK
2934  skew => wadats(imod)%SKEW
2935  embia1 => wadats(imod)%EMBIA1
2936  embia2 => wadats(imod)%EMBIA2
2937  !
2938  dtdyn => wadats(imod)%DTDYN
2939  fcut => wadats(imod)%FCUT
2940  cflxymax => wadats(imod)%CFLXYMAX
2941  cflthmax => wadats(imod)%CFLTHMAX
2942  cflkmax => wadats(imod)%CFLKMAX
2943  !
2944  usero => wadats(imod)%USERO
2945  !
2946  wn => wadats(imod)%WN
2947 #ifdef W3_IC3
2948  ic3wn_r=> wadats(imod)%IC3WN_R
2949  ic3wn_i=> wadats(imod)%IC3WN_I
2950 #endif
2951  !
2952  IF ( fl_all ) THEN
2953  !
2954  cg => wadats(imod)%CG
2955 #ifdef W3_IC3
2956  ic3cg => wadats(imod)%IC3CG
2957 #endif
2958  !
2959  atrnx => wadats(imod)%ATRNX
2960  atrny => wadats(imod)%ATRNY
2961  !
2962  dddx => wadats(imod)%DDDX
2963  dddy => wadats(imod)%DDDY
2964  dcdx => wadats(imod)%DCDX
2965  dcdy => wadats(imod)%DCDY
2966  dcxdx => wadats(imod)%DCXDX
2967  dcydx => wadats(imod)%DCYDX
2968  dcxdy => wadats(imod)%DCXDY
2969  dcydy => wadats(imod)%DCYDY
2970  !
2971 #ifdef W3_SMC
2972  dhdx => wadats(imod)%DHDX
2973  dhdy => wadats(imod)%DHDY
2974  dhlmt => wadats(imod)%DHLMT
2975 #endif
2976  !
2977  alpha => wadats(imod)%ALPHA
2978  !
2979  IF ( inputs(imod)%INFLAGS1(2) ) THEN
2980  ca0 => wadats(imod)%CA0
2981  cai => wadats(imod)%CAI
2982  cd0 => wadats(imod)%CD0
2983  cdi => wadats(imod)%CDI
2984  END IF
2985  !
2986  IF ( inputs(imod)%INFLAGS1(3) ) THEN
2987  ua0 => wadats(imod)%UA0
2988  uai => wadats(imod)%UAI
2989  ud0 => wadats(imod)%UD0
2990  udi => wadats(imod)%UDI
2991  as0 => wadats(imod)%AS0
2992  asi => wadats(imod)%ASI
2993  END IF
2994  !
2995  IF ( inputs(imod)%INFLAGS1(5) ) THEN
2996  ma0 => wadats(imod)%MA0
2997  mai => wadats(imod)%MAI
2998  md0 => wadats(imod)%MD0
2999  mdi => wadats(imod)%MDI
3000  END IF
3001  !
3002  IF ( inputs(imod)%INFLAGS1(6) ) THEN
3003  ra0 => wadats(imod)%RA0
3004  rai => wadats(imod)%RAI
3005  END IF
3006  !
3007 #ifdef W3_PR1
3008  is0 => wadats(imod)%IS0
3009  is2 => wadats(imod)%IS2
3010  facvx => wadats(imod)%FACVX
3011  facvy => wadats(imod)%FACVY
3012 #endif
3013  !
3014 #ifdef W3_PR2
3015  mapx2 => wadats(imod)%MAPX2
3016  mapy2 => wadats(imod)%MAPY2
3017  mapaxy => wadats(imod)%MAPAXY
3018  mapxy => wadats(imod)%MAPXY
3019  mapth2 => wadats(imod)%MAPTH2
3020  mapwn2 => wadats(imod)%MAPWN2
3021 #endif
3022  !
3023 #ifdef W3_PR3
3024  mapx2 => wadats(imod)%MAPX2
3025  mapy2 => wadats(imod)%MAPY2
3026  mapaxy => wadats(imod)%MAPAXY
3027  mapcxy => wadats(imod)%MAPCXY
3028  mapth2 => wadats(imod)%MAPTH2
3029  mapwn2 => wadats(imod)%MAPWN2
3030  maptrn => wadats(imod)%MAPTRN
3031 #endif
3032  !
3033  IF (gtype .EQ. ungtype) iter => wadats(imod)%ITER
3034  !
3035  iappro => wadats(imod)%IAPPRO
3036  sppnt => wadats(imod)%SPPNT
3037  !
3038  END IF
3039  !
3040  END IF
3041  !
3042 #ifdef W3_NL1
3043  IF ( nlinit ) THEN
3044  ip11 => wadats(imod)%IP11
3045  ip12 => wadats(imod)%IP12
3046  ip13 => wadats(imod)%IP13
3047  ip14 => wadats(imod)%IP14
3048  im11 => wadats(imod)%IM11
3049  im12 => wadats(imod)%IM12
3050  im13 => wadats(imod)%IM13
3051  im14 => wadats(imod)%IM14
3052  ip21 => wadats(imod)%IP21
3053  ip22 => wadats(imod)%IP22
3054  ip23 => wadats(imod)%IP23
3055  ip24 => wadats(imod)%IP24
3056  im21 => wadats(imod)%IM21
3057  im22 => wadats(imod)%IM22
3058  im23 => wadats(imod)%IM23
3059  im24 => wadats(imod)%IM24
3060  ic11 => wadats(imod)%IC11
3061  ic12 => wadats(imod)%IC12
3062  ic21 => wadats(imod)%IC21
3063  ic22 => wadats(imod)%IC22
3064  ic31 => wadats(imod)%IC31
3065  ic32 => wadats(imod)%IC32
3066  ic41 => wadats(imod)%IC41
3067  ic42 => wadats(imod)%IC42
3068  ic51 => wadats(imod)%IC51
3069  ic52 => wadats(imod)%IC52
3070  ic61 => wadats(imod)%IC61
3071  ic62 => wadats(imod)%IC62
3072  ic71 => wadats(imod)%IC71
3073  ic72 => wadats(imod)%IC72
3074  ic81 => wadats(imod)%IC81
3075  ic82 => wadats(imod)%IC82
3076  af11 => wadats(imod)%AF11
3077  END IF
3078 #endif
3079 
3080 #ifdef W3_MPI
3081  IF ( nrqsg1 .NE. 0 ) THEN
3082  irqsg1 => wadats(imod)%IRQSG1
3083  irqsg2 => wadats(imod)%IRQSG2
3084  END IF
3085  gstore => wadats(imod)%GSTORE
3086  sstore => wadats(imod)%SSTORE
3087 #endif
3088  !
3089  RETURN
3090  !
3091  ! Formats
3092  !
3093 1001 FORMAT (/' *** ERROR W3SETA : GRIDS NOT INITIALIZED *** '/ &
3094  ' RUN W3NMOD FIRST '/)
3095 1002 FORMAT (/' *** ERROR W3SETA : ILLEGAL MODEL NUMBER *** '/ &
3096  ' IMOD = ',i10/ &
3097  ' NADATA = ',i10/)
3098  !
3099 #ifdef W3_T
3100 9000 FORMAT (' TEST W3SETA : MODEL ',i4,' SELECTED')
3101 #endif
3102  !/
3103  !/ End of W3SETA ----------------------------------------------------- /
3104  !/
3105  END SUBROUTINE w3seta
3106  !/ ------------------------------------------------------------------- /
3117  SUBROUTINE w3xeta ( IMOD, NDSE, NDST )
3118  !/
3119  !/ +-----------------------------------+
3120  !/ | WAVEWATCH III NOAA/NCEP |
3121  !/ | H. L. Tolman |
3122  !/ | FORTRAN 90 |
3123  !/ | Last update : 22-Mar-2021 |
3124  !/ +-----------------------------------+
3125  !/
3126  !/ 25-Dec-2012 : Origination. ( version 4.11 )
3127  !/ 30-Apr-2014 : Add s/th1-2m ( version 5.01 )
3128  !/ 22-Mar-2021 : Adds WNMEAN, TAUOC parameters ( version 7.13 )
3129  !/
3130  ! 1. Purpose :
3131  !
3132  ! Reduced version of W3SETA to point t expended output arrays.
3133  !
3134  ! 10. Source code :
3135  !
3136  !/ ------------------------------------------------------------------- /
3137  !
3138  USE w3idatmd, ONLY: inputs
3139  USE w3gdatmd, ONLY: e3df, p2msf, us3df, usspf, gtype, ungtype
3140  !
3141  USE w3servmd, ONLY: extcde
3142 #ifdef W3_S
3143  USE w3servmd, ONLY: strace
3144 #endif
3145  !
3146  !/
3147  !/ ------------------------------------------------------------------- /
3148  !/ Parameter list
3149  !/
3150  INTEGER, INTENT(IN) :: IMOD, NDSE, NDST
3151  !/
3152  !/ ------------------------------------------------------------------- /
3153  !/ Local parameters
3154  !/
3155 #ifdef W3_S
3156  INTEGER, SAVE :: IENT = 0
3157  CALL strace (ient, 'W3XETA')
3158 #endif
3159  !
3160  ! -------------------------------------------------------------------- /
3161  ! 1. Test input and module status
3162  !
3163  IF ( nadata .EQ. -1 ) THEN
3164  WRITE (ndse,1001)
3165  CALL extcde (1)
3166  END IF
3167  !
3168  IF ( imod.LT.1 .OR. imod.GT.nadata ) THEN
3169  WRITE (ndse,1002) imod, nadata
3170  CALL extcde (2)
3171  END IF
3172  !
3173 #ifdef W3_T
3174  WRITE (ndst,9000) imod
3175 #endif
3176  !
3177  ! -------------------------------------------------------------------- /
3178  ! 2. Set model numbers
3179  !
3180  iadata = imod
3181  !
3182  ! -------------------------------------------------------------------- /
3183  ! 3. Set pointers
3184  !
3185  IF ( ainit2 ) THEN
3186  !
3187  hs => wadats(imod)%XHS
3188  wlm => wadats(imod)%XWLM
3189  t02 => wadats(imod)%XT02
3190  t0m1 => wadats(imod)%XT0M1
3191  t01 => wadats(imod)%XT01
3192  fp0 => wadats(imod)%XFP0
3193  thm => wadats(imod)%XTHM
3194  ths => wadats(imod)%XTHS
3195  thp0 => wadats(imod)%XTHP0
3196  hsig => wadats(imod)%XHSIG
3197  stmaxe => wadats(imod)%XSTMAXE
3198  stmaxd => wadats(imod)%XSTMAXD
3199  hmaxe => wadats(imod)%XHMAXE
3200  hmaxd => wadats(imod)%XHMAXD
3201  hcmaxe => wadats(imod)%XHCMAXE
3202  hcmaxd => wadats(imod)%XHCMAXD
3203  qp => wadats(imod)%XQP
3204  wbt => wadats(imod)%XWBT
3205  wnmean => wadats(imod)%XWNMEAN
3206  !
3207  ef => wadats(imod)%XEF
3208  th1m => wadats(imod)%XTH1M
3209  sth1m => wadats(imod)%XSTH1M
3210  th2m => wadats(imod)%XTH2M
3211  sth2m => wadats(imod)%XSTH2M
3212  !
3213  phs => wadats(imod)%XPHS
3214  ptp => wadats(imod)%XPTP
3215  plp => wadats(imod)%XPLP
3216  pdir => wadats(imod)%XPDIR
3217  psi => wadats(imod)%XPSI
3218  pws => wadats(imod)%XPWS
3219  pwst => wadats(imod)%XPWST
3220  pnr => wadats(imod)%XPNR
3221  pthp0 => wadats(imod)%XPTHP0
3222  pqp => wadats(imod)%XPQP
3223  ppe => wadats(imod)%XPPE
3224  pgw => wadats(imod)%XPGW
3225  psw => wadats(imod)%XPSW
3226  ptm1 => wadats(imod)%XPTM1
3227  pt1 => wadats(imod)%XPT1
3228  pt2 => wadats(imod)%XPT2
3229  pep => wadats(imod)%XPEP
3230  !
3231  charn => wadats(imod)%XCHARN
3232  tws => wadats(imod)%XTWS
3233  cge => wadats(imod)%XCGE
3234  phiaw => wadats(imod)%XPHIAW
3235  tauwix => wadats(imod)%XTAUWIX
3236  tauwiy => wadats(imod)%XTAUWIY
3237  tauwnx => wadats(imod)%XTAUWNX
3238  tauwny => wadats(imod)%XTAUWNY
3239  whitecap => wadats(imod)%XWHITECAP
3240  !
3241  sxx => wadats(imod)%XSXX
3242  syy => wadats(imod)%XSYY
3243  sxy => wadats(imod)%XSXY
3244  tauox => wadats(imod)%XTAUOX
3245  tauoy => wadats(imod)%XTAUOY
3246  bhd => wadats(imod)%XBHD
3247  phioc => wadats(imod)%XPHIOC
3248  tusx => wadats(imod)%XTUSX
3249  tusy => wadats(imod)%XTUSY
3250  ussx => wadats(imod)%XUSSX
3251  ussy => wadats(imod)%XUSSY
3252  prms => wadats(imod)%XPRMS
3253  tpms => wadats(imod)%XTPMS
3254  p2sms => wadats(imod)%XP2SMS
3255  us3d => wadats(imod)%XUS3D
3256  phice => wadats(imod)%XPHICE
3257  tauice => wadats(imod)%XTAUICE
3258  ussp => wadats(imod)%XUSSP
3259  tauocx => wadats(imod)%XTAUOCX
3260  tauocy => wadats(imod)%XTAUOCY
3261  aba => wadats(imod)%XABA
3262  abd => wadats(imod)%XABD
3263  uba => wadats(imod)%XUBA
3264  ubd => wadats(imod)%XUBD
3265  bedforms=> wadats(imod)%XBEDFORMS
3266  phibbl => wadats(imod)%XPHIBBL
3267  taubbl => wadats(imod)%XTAUBBL
3268  !
3269  mssx => wadats(imod)%XMSSX
3270  mssy => wadats(imod)%XMSSY
3271  mssd => wadats(imod)%XMSSD
3272  mscx => wadats(imod)%XMSCX
3273  mscy => wadats(imod)%XMSCY
3274  mscd => wadats(imod)%XMSCD
3275  qkk => wadats(imod)%XQKK
3276  skew => wadats(imod)%XSKEW
3277  embia1 => wadats(imod)%XEMBIA1
3278  embia2 => wadats(imod)%XEMBIA2
3279  !
3280  dtdyn => wadats(imod)%XDTDYN
3281  fcut => wadats(imod)%XFCUT
3282  cflxymax => wadats(imod)%XCFLXYMAX
3283  cflthmax => wadats(imod)%XCFLTHMAX
3284  cflkmax => wadats(imod)%XCFLKMAX
3285  !
3286  usero => wadats(imod)%XUSERO
3287  !
3288  END IF
3289  !
3290  RETURN
3291  !
3292  ! Formats
3293  !
3294 1001 FORMAT (/' *** ERROR W3XETA : GRIDS NOT INITIALIZED *** '/ &
3295  ' RUN W3NMOD FIRST '/)
3296 1002 FORMAT (/' *** ERROR W3XETA : ILLEGAL MODEL NUMBER *** '/ &
3297  ' IMOD = ',i10/ &
3298  ' NADATA = ',i10/)
3299  !
3300 #ifdef W3_T
3301 9000 FORMAT (' TEST W3XETA : MODEL ',i4,' SELECTED')
3302 #endif
3303  !/
3304  !/ End of W3XETA ----------------------------------------------------- /
3305  !/
3306  END SUBROUTINE w3xeta
3307  !/
3308  !/ End of module W3ADATMD -------------------------------------------- /
3309  !/
3310 END MODULE w3adatmd
w3adatmd::is0
integer, dimension(:), pointer is0
Definition: w3adatmd.F90:635
w3adatmd::pt2
real, dimension(:,:), pointer pt2
Definition: w3adatmd.F90:597
w3adatmd::psw
real, dimension(:,:), pointer psw
Definition: w3adatmd.F90:597
w3gdatmd::nk
integer, pointer nk
Definition: w3gdatmd.F90:1230
w3adatmd::hcmaxe
real, dimension(:), pointer hcmaxe
Definition: w3adatmd.F90:587
w3adatmd::awg2
real, pointer awg2
Definition: w3adatmd.F90:666
w3gdatmd::nseal
integer, pointer nseal
Definition: w3gdatmd.F90:1097
w3adatmd::ic3cg
real, dimension(:,:), pointer ic3cg
Definition: w3adatmd.F90:576
w3adatmd::ic3wn_i
real, dimension(:,:), pointer ic3wn_i
Definition: w3adatmd.F90:576
w3adatmd::phice
real, dimension(:), pointer phice
Definition: w3adatmd.F90:607
w3adatmd::th2m
real, dimension(:,:), pointer th2m
Definition: w3adatmd.F90:594
w3adatmd::swg6
real, pointer swg6
Definition: w3adatmd.F90:666
w3adatmd::charn
real, dimension(:), pointer charn
Definition: w3adatmd.F90:603
w3adatmd::dtdyn
real, dimension(:), pointer dtdyn
Definition: w3adatmd.F90:620
w3adatmd::nsealm
integer, pointer nsealm
Definition: w3adatmd.F90:686
w3adatmd::as
real, dimension(:), pointer as
Definition: w3adatmd.F90:584
w3adatmd::ic52
integer, dimension(:), pointer ic52
Definition: w3adatmd.F90:658
w3adatmd::awg8
real, pointer awg8
Definition: w3adatmd.F90:666
w3adatmd
Define data structures to set up wave model auxiliary data for several models simultaneously.
Definition: w3adatmd.F90:26
w3adatmd::hcmaxd
real, dimension(:), pointer hcmaxd
Definition: w3adatmd.F90:587
w3gdatmd::nspec
integer, pointer nspec
Definition: w3gdatmd.F90:1230
w3adatmd::ip14
integer, dimension(:), pointer ip14
Definition: w3adatmd.F90:658
w3adatmd::sth1m
real, dimension(:,:), pointer sth1m
Definition: w3adatmd.F90:594
w3adatmd::ussy
real, dimension(:), pointer ussy
Definition: w3adatmd.F90:607
w3adatmd::mapy2
integer, dimension(:), pointer mapy2
Definition: w3adatmd.F90:642
w3adatmd::im13
integer, dimension(:), pointer im13
Definition: w3adatmd.F90:658
w3adatmd::im21
integer, dimension(:), pointer im21
Definition: w3adatmd.F90:658
w3adatmd::dcdx
real, dimension(:,:,:), pointer dcdx
Definition: w3adatmd.F90:629
w3adatmd::nrqsg2
integer, pointer nrqsg2
Definition: w3adatmd.F90:676
w3adatmd::pep
real, dimension(:,:), pointer pep
Definition: w3adatmd.F90:597
w3adatmd::mscd
real, dimension(:), pointer mscd
Definition: w3adatmd.F90:617
w3adatmd::ic11
integer, dimension(:), pointer ic11
Definition: w3adatmd.F90:658
w3adatmd::abd
real, dimension(:), pointer abd
Definition: w3adatmd.F90:614
w3gdatmd::ungtype
integer, parameter ungtype
Definition: w3gdatmd.F90:626
w3adatmd::ip13
integer, dimension(:), pointer ip13
Definition: w3adatmd.F90:658
w3adatmd::ic31
integer, dimension(:), pointer ic31
Definition: w3adatmd.F90:658
w3gdatmd::p2msf
integer, dimension(:), pointer p2msf
Definition: w3gdatmd.F90:1098
w3adatmd::atrnx
real, dimension(:,:), pointer atrnx
Definition: w3adatmd.F90:578
w3adatmd::dhdy
real, dimension(:), pointer dhdy
Definition: w3adatmd.F90:631
w3adatmd::stmaxe
real, dimension(:), pointer stmaxe
Definition: w3adatmd.F90:587
w3adatmd::dhlmt
real, dimension(:,:), pointer dhlmt
Definition: w3adatmd.F90:631
w3adatmd::tauice
real, dimension(:,:), pointer tauice
Definition: w3adatmd.F90:607
w3adatmd::t02
real, dimension(:), pointer t02
Definition: w3adatmd.F90:587
w3adatmd::us3d
real, dimension(:,:), pointer us3d
Definition: w3adatmd.F90:612
w3adatmd::wadats
type(wadat), dimension(:), allocatable, target wadats
Definition: w3adatmd.F90:571
w3adatmd::cflxymax
real, dimension(:), pointer cflxymax
Definition: w3adatmd.F90:620
w3adatmd::cg
real, dimension(:,:), pointer cg
Definition: w3adatmd.F90:575
w3adatmd::tws
real, dimension(:), pointer tws
Definition: w3adatmd.F90:603
w3adatmd::tusx
real, dimension(:), pointer tusx
Definition: w3adatmd.F90:607
w3adatmd::mapxy
integer, dimension(:), pointer mapxy
Definition: w3adatmd.F90:642
w3adatmd::fcut
real, dimension(:), pointer fcut
Definition: w3adatmd.F90:620
w3adatmd::nsploc
integer, pointer nsploc
Definition: w3adatmd.F90:676
w3odatmd::flogr2
logical, dimension(:,:), pointer flogr2
Definition: w3odatmd.F90:478
w3adatmd::bispl
integer, dimension(:), pointer bispl
Definition: w3adatmd.F90:680
w3idatmd::inputs
type(input), dimension(:), allocatable, target inputs
Definition: w3idatmd.F90:232
w3adatmd::nmx0
integer, pointer nmx0
Definition: w3adatmd.F90:640
w3adatmd::ic21
integer, dimension(:), pointer ic21
Definition: w3adatmd.F90:658
w3adatmd::im11
integer, dimension(:), pointer im11
Definition: w3adatmd.F90:658
w3adatmd::tusy
real, dimension(:), pointer tusy
Definition: w3adatmd.F90:607
w3adatmd::ptp
real, dimension(:,:), pointer ptp
Definition: w3adatmd.F90:597
w3adatmd::dw
real, dimension(:), pointer dw
Definition: w3adatmd.F90:584
w3adatmd::u10d
real, dimension(:), pointer u10d
Definition: w3adatmd.F90:584
w3adatmd::atrny
real, dimension(:,:), pointer atrny
Definition: w3adatmd.F90:578
w3adatmd::swg5
real, pointer swg5
Definition: w3adatmd.F90:666
w3odatmd::ntproc
integer, pointer ntproc
Definition: w3odatmd.F90:457
w3adatmd::awg6
real, pointer awg6
Definition: w3adatmd.F90:666
w3adatmd::af11
real, dimension(:), pointer af11
Definition: w3adatmd.F90:670
w3adatmd::th1m
real, dimension(:,:), pointer th1m
Definition: w3adatmd.F90:594
w3adatmd::t01
real, dimension(:), pointer t01
Definition: w3adatmd.F90:587
w3adatmd::cge
real, dimension(:), pointer cge
Definition: w3adatmd.F90:603
w3adatmd::ic32
integer, dimension(:), pointer ic32
Definition: w3adatmd.F90:658
w3adatmd::mapwn2
integer, dimension(:), pointer mapwn2
Definition: w3adatmd.F90:642
w3adatmd::ic82
integer, dimension(:), pointer ic82
Definition: w3adatmd.F90:658
w3adatmd::ipass
integer, pointer ipass
Definition: w3adatmd.F90:686
w3odatmd::iaproc
integer, pointer iaproc
Definition: w3odatmd.F90:457
w3odatmd::ngrpp
integer, parameter ngrpp
Definition: w3odatmd.F90:324
w3adatmd::iadata
integer iadata
Definition: w3adatmd.F90:374
w3adatmd::pdir
real, dimension(:,:), pointer pdir
Definition: w3adatmd.F90:597
w3adatmd::thp0
real, dimension(:), pointer thp0
Definition: w3adatmd.F90:587
w3adatmd::tauocy
real, dimension(:), pointer tauocy
Definition: w3adatmd.F90:607
w3adatmd::ic3wn_r
real, dimension(:,:), pointer ic3wn_r
Definition: w3adatmd.F90:576
w3adatmd::phs
real, dimension(:,:), pointer phs
Definition: w3adatmd.F90:597
w3adatmd::hs
real, dimension(:), pointer hs
Definition: w3adatmd.F90:587
w3adatmd::im22
integer, dimension(:), pointer im22
Definition: w3adatmd.F90:658
w3adatmd::fl_all
logical, pointer fl_all
Definition: w3adatmd.F90:688
w3adatmd::ainit2
logical, pointer ainit2
Definition: w3adatmd.F90:688
w3gdatmd::ny
integer, pointer ny
Definition: w3gdatmd.F90:1097
w3adatmd::nlinit
logical, pointer nlinit
Definition: w3adatmd.F90:671
w3adatmd::w3dmnl
subroutine w3dmnl(IMOD, NDSE, NDST, NSP, NSPX)
Initialize an individual data grid at the proper dimensions (DIA).
Definition: w3adatmd.F90:2431
w3adatmd::swg4
real, pointer swg4
Definition: w3adatmd.F90:666
w3adatmd::asi
real, dimension(:), pointer asi
Definition: w3adatmd.F90:578
w3adatmd::uba
real, dimension(:), pointer uba
Definition: w3adatmd.F90:614
w3adatmd::iappro
integer, dimension(:), pointer iappro
Definition: w3adatmd.F90:674
w3adatmd::pqp
real, dimension(:,:), pointer pqp
Definition: w3adatmd.F90:597
w3adatmd::cdi
real, dimension(:), pointer cdi
Definition: w3adatmd.F90:578
w3adatmd::ip21
integer, dimension(:), pointer ip21
Definition: w3adatmd.F90:658
w3adatmd::facvy
real, dimension(:), pointer facvy
Definition: w3adatmd.F90:636
w3idatmd::flcur
logical, pointer flcur
Definition: w3idatmd.F90:261
w3adatmd::im12
integer, dimension(:), pointer im12
Definition: w3adatmd.F90:658
w3adatmd::tauwix
real, dimension(:), pointer tauwix
Definition: w3adatmd.F90:603
w3adatmd::w3dima
subroutine w3dima(IMOD, NDSE, NDST, D_ONLY)
Initialize an individual data grid at the proper dimensions.
Definition: w3adatmd.F90:846
w3adatmd::pthp0
real, dimension(:,:), pointer pthp0
Definition: w3adatmd.F90:597
w3adatmd::dcydx
real, dimension(:,:), pointer dcydx
Definition: w3adatmd.F90:627
w3gdatmd::w3setg
subroutine w3setg(IMOD, NDSE, NDST)
Definition: w3gdatmd.F90:2152
w3adatmd::swg7
real, pointer swg7
Definition: w3adatmd.F90:666
w3adatmd::w3xeta
subroutine w3xeta(IMOD, NDSE, NDST)
Reduced version of W3SETA to point to expended output arrays.
Definition: w3adatmd.F90:3118
w3adatmd::tauwiy
real, dimension(:), pointer tauwiy
Definition: w3adatmd.F90:603
w3adatmd::ca0
real, dimension(:), pointer ca0
Definition: w3adatmd.F90:578
w3adatmd::taubbl
real, dimension(:,:), pointer taubbl
Definition: w3adatmd.F90:614
w3adatmd::nfrhgh
integer, pointer nfrhgh
Definition: w3adatmd.F90:657
w3adatmd::ainit
logical, pointer ainit
Definition: w3adatmd.F90:688
w3adatmd::ef
real, dimension(:,:), pointer ef
Definition: w3adatmd.F90:594
w3adatmd::w3seta
subroutine w3seta(IMOD, NDSE, NDST)
Select one of the WAVEWATCH III grids / models.
Definition: w3adatmd.F90:2645
w3adatmd::tauwny
real, dimension(:), pointer tauwny
Definition: w3adatmd.F90:603
w3adatmd::ra0
real, dimension(:), pointer ra0
Definition: w3adatmd.F90:578
w3adatmd::dcdy
real, dimension(:,:,:), pointer dcdy
Definition: w3adatmd.F90:629
w3adatmd::phiaw
real, dimension(:), pointer phiaw
Definition: w3adatmd.F90:603
w3adatmd::pws
real, dimension(:,:), pointer pws
Definition: w3adatmd.F90:597
w3adatmd::plp
real, dimension(:,:), pointer plp
Definition: w3adatmd.F90:597
w3adatmd::phibbl
real, dimension(:), pointer phibbl
Definition: w3adatmd.F90:614
w3adatmd::mpibuf
integer, parameter mpibuf
Definition: w3adatmd.F90:376
w3adatmd::awg7
real, pointer awg7
Definition: w3adatmd.F90:666
constants::lpdlib
logical lpdlib
LPDLIB Logical for using the PDLIB library.
Definition: constants.F90:101
w3adatmd::im23
integer, dimension(:), pointer im23
Definition: w3adatmd.F90:658
w3adatmd::swg2
real, pointer swg2
Definition: w3adatmd.F90:666
w3adatmd::cflthmax
real, dimension(:), pointer cflthmax
Definition: w3adatmd.F90:620
w3gdatmd::nsea
integer, pointer nsea
Definition: w3gdatmd.F90:1097
w3adatmd::ic72
integer, dimension(:), pointer ic72
Definition: w3adatmd.F90:658
w3adatmd::skew
real, dimension(:), pointer skew
Definition: w3adatmd.F90:617
w3adatmd::md0
real, dimension(:), pointer md0
Definition: w3adatmd.F90:578
w3adatmd::psi
real, dimension(:,:), pointer psi
Definition: w3adatmd.F90:597
w3adatmd::sth2m
real, dimension(:,:), pointer sth2m
Definition: w3adatmd.F90:594
w3adatmd::tpms
real, dimension(:), pointer tpms
Definition: w3adatmd.F90:607
w3servmd
Definition: w3servmd.F90:3
w3adatmd::dcxdx
real, dimension(:,:), pointer dcxdx
Definition: w3adatmd.F90:627
w3adatmd::embia1
real, dimension(:), pointer embia1
Definition: w3adatmd.F90:617
w3odatmd::flogrd
logical, dimension(:,:), pointer flogrd
Definition: w3odatmd.F90:478
w3adatmd::ths
real, dimension(:), pointer ths
Definition: w3adatmd.F90:587
w3adatmd::bedforms
real, dimension(:,:), pointer bedforms
Definition: w3adatmd.F90:614
w3adatmd::ud
real, dimension(:), pointer ud
Definition: w3adatmd.F90:584
w3adatmd::hmaxd
real, dimension(:), pointer hmaxd
Definition: w3adatmd.F90:587
w3adatmd::pwst
real, dimension(:), pointer pwst
Definition: w3adatmd.F90:597
w3adatmd::qkk
real, dimension(:), pointer qkk
Definition: w3adatmd.F90:617
w3adatmd::irqsg1
integer, dimension(:,:), pointer irqsg1
Definition: w3adatmd.F90:681
w3adatmd::swg1
real, pointer swg1
Definition: w3adatmd.F90:666
w3adatmd::nspecx
integer, pointer nspecx
Definition: w3adatmd.F90:657
w3adatmd::cai
real, dimension(:), pointer cai
Definition: w3adatmd.F90:578
w3adatmd::sxy
real, dimension(:), pointer sxy
Definition: w3adatmd.F90:607
w3adatmd::ud0
real, dimension(:), pointer ud0
Definition: w3adatmd.F90:578
w3servmd::print_memcheck
subroutine print_memcheck(iun, msg)
Write memory statistics if requested.
Definition: w3servmd.F90:2033
w3adatmd::gstore
real, dimension(:,:), pointer gstore
Definition: w3adatmd.F90:682
w3adatmd::tauwnx
real, dimension(:), pointer tauwnx
Definition: w3adatmd.F90:603
w3adatmd::alpha
real, dimension(:,:), pointer alpha
Definition: w3adatmd.F90:687
w3adatmd::mapcxy
integer, dimension(:), pointer mapcxy
Definition: w3adatmd.F90:649
w3gdatmd::nth
integer, pointer nth
Definition: w3gdatmd.F90:1230
w3odatmd
Definition: w3odatmd.F90:3
w3adatmd::dddy
real, dimension(:,:), pointer dddy
Definition: w3adatmd.F90:627
w3adatmd::mdi
real, dimension(:), pointer mdi
Definition: w3adatmd.F90:578
w3adatmd::wbt
real, dimension(:), pointer wbt
Definition: w3adatmd.F90:587
w3adatmd::bhd
real, dimension(:), pointer bhd
Definition: w3adatmd.F90:607
w3adatmd::w3naux
subroutine w3naux(NDSE, NDST)
Set up the number of grids to be used.
Definition: w3adatmd.F90:704
w3adatmd::cy
real, dimension(:), pointer cy
Definition: w3adatmd.F90:584
w3adatmd::pnr
real, dimension(:), pointer pnr
Definition: w3adatmd.F90:597
w3adatmd::ip11
integer, dimension(:), pointer ip11
Definition: w3adatmd.F90:658
w3adatmd::idlast
integer, pointer idlast
Definition: w3adatmd.F90:686
w3adatmd::taua
real, dimension(:), pointer taua
Definition: w3adatmd.F90:584
w3adatmd::hmaxe
real, dimension(:), pointer hmaxe
Definition: w3adatmd.F90:587
w3adatmd::pt1
real, dimension(:,:), pointer pt1
Definition: w3adatmd.F90:597
w3adatmd::wlm
real, dimension(:), pointer wlm
Definition: w3adatmd.F90:587
w3odatmd::naproc
integer, pointer naproc
Definition: w3odatmd.F90:457
w3adatmd::ma0
real, dimension(:), pointer ma0
Definition: w3adatmd.F90:578
w3adatmd::ic22
integer, dimension(:), pointer ic22
Definition: w3adatmd.F90:658
w3adatmd::ic12
integer, dimension(:), pointer ic12
Definition: w3adatmd.F90:658
w3gdatmd::us3df
integer, dimension(:), pointer us3df
Definition: w3gdatmd.F90:1098
w3gdatmd::igrid
integer igrid
Definition: w3gdatmd.F90:618
w3adatmd::mapth2
integer, dimension(:), pointer mapth2
Definition: w3adatmd.F90:642
w3adatmd::uai
real, dimension(:), pointer uai
Definition: w3adatmd.F90:578
w3adatmd::wnmean
real, dimension(:), pointer wnmean
Definition: w3adatmd.F90:587
w3adatmd::cflkmax
real, dimension(:), pointer cflkmax
Definition: w3adatmd.F90:620
w3adatmd::sstore
real, dimension(:,:), pointer sstore
Definition: w3adatmd.F90:682
w3adatmd::ic41
integer, dimension(:), pointer ic41
Definition: w3adatmd.F90:658
w3adatmd::ua0
real, dimension(:), pointer ua0
Definition: w3adatmd.F90:578
w3adatmd::phioc
real, dimension(:), pointer phioc
Definition: w3adatmd.F90:607
w3adatmd::ic62
integer, dimension(:), pointer ic62
Definition: w3adatmd.F90:658
w3adatmd::wn
real, dimension(:,:), pointer wn
Definition: w3adatmd.F90:575
w3adatmd::mai
real, dimension(:), pointer mai
Definition: w3adatmd.F90:578
w3adatmd::u10
real, dimension(:), pointer u10
Definition: w3adatmd.F90:584
w3adatmd::nmx2
integer, pointer nmx2
Definition: w3adatmd.F90:640
w3adatmd::fliwnd
logical, pointer fliwnd
Definition: w3adatmd.F90:688
w3adatmd::dcydy
real, dimension(:,:), pointer dcydy
Definition: w3adatmd.F90:627
w3adatmd::ww3_spec_vec
integer, pointer ww3_spec_vec
Definition: w3adatmd.F90:676
w3adatmd::qp
real, dimension(:), pointer qp
Definition: w3adatmd.F90:587
w3adatmd::awg5
real, pointer awg5
Definition: w3adatmd.F90:666
w3idatmd::flwind
logical, pointer flwind
Definition: w3idatmd.F90:261
w3adatmd::swg8
real, pointer swg8
Definition: w3adatmd.F90:666
w3adatmd::stmaxd
real, dimension(:), pointer stmaxd
Definition: w3adatmd.F90:587
w3servmd::strace
subroutine strace(IENT, SNAME)
Definition: w3servmd.F90:148
w3adatmd::ic71
integer, dimension(:), pointer ic71
Definition: w3adatmd.F90:658
w3gdatmd::gtype
integer, pointer gtype
Definition: w3gdatmd.F90:1094
w3idatmd
Define data structures to set up wave model input data for several models simultaneously.
Definition: w3idatmd.F90:16
w3adatmd::nfr
integer, pointer nfr
Definition: w3adatmd.F90:657
w3adatmd::ip22
integer, dimension(:), pointer ip22
Definition: w3adatmd.F90:658
w3adatmd::as0
real, dimension(:), pointer as0
Definition: w3adatmd.F90:578
w3adatmd::whitecap
real, dimension(:,:), pointer whitecap
Definition: w3adatmd.F90:603
w3adatmd::awg3
real, pointer awg3
Definition: w3adatmd.F90:666
w3adatmd::sppnt
real, dimension(:,:,:), pointer sppnt
Definition: w3adatmd.F90:684
w3adatmd::ic51
integer, dimension(:), pointer ic51
Definition: w3adatmd.F90:658
w3adatmd::facvx
real, dimension(:), pointer facvx
Definition: w3adatmd.F90:636
w3adatmd::ip12
integer, dimension(:), pointer ip12
Definition: w3adatmd.F90:658
w3adatmd::tauox
real, dimension(:), pointer tauox
Definition: w3adatmd.F90:607
w3adatmd::mapx2
integer, dimension(:), pointer mapx2
Definition: w3adatmd.F90:642
w3odatmd::noextr
integer, parameter noextr
Definition: w3odatmd.F90:328
w3adatmd::p2sms
real, dimension(:,:), pointer p2sms
Definition: w3adatmd.F90:612
w3adatmd::prms
real, dimension(:), pointer prms
Definition: w3adatmd.F90:607
w3odatmd::napfld
integer, pointer napfld
Definition: w3odatmd.F90:457
w3adatmd::usero
real, dimension(:,:), pointer usero
Definition: w3adatmd.F90:623
w3adatmd::nspecy
integer, pointer nspecy
Definition: w3adatmd.F90:657
w3adatmd::mssy
real, dimension(:), pointer mssy
Definition: w3adatmd.F90:617
w3idatmd::fltaua
logical, pointer fltaua
Definition: w3idatmd.F90:261
w3adatmd::im14
integer, dimension(:), pointer im14
Definition: w3adatmd.F90:658
w3adatmd::mpi_comm_wave
integer, pointer mpi_comm_wave
Definition: w3adatmd.F90:676
w3adatmd::udi
real, dimension(:), pointer udi
Definition: w3adatmd.F90:578
w3adatmd::maptrn
logical, dimension(:), pointer maptrn
Definition: w3adatmd.F90:651
w3adatmd::ic42
integer, dimension(:), pointer ic42
Definition: w3adatmd.F90:658
constants
Define some much-used constants for global use (all defined as PARAMETER).
Definition: constants.F90:20
w3adatmd::ua
real, dimension(:), pointer ua
Definition: w3adatmd.F90:584
w3adatmd::tauoy
real, dimension(:), pointer tauoy
Definition: w3adatmd.F90:607
w3odatmd::nogrp
integer, parameter nogrp
Definition: w3odatmd.F90:323
w3adatmd::fp0
real, dimension(:), pointer fp0
Definition: w3adatmd.F90:587
w3adatmd::iter
integer, dimension(:,:), pointer iter
Definition: w3adatmd.F90:654
w3gdatmd
Definition: w3gdatmd.F90:16
w3adatmd::nact
integer, pointer nact
Definition: w3adatmd.F90:640
w3adatmd::dal2
real, pointer dal2
Definition: w3adatmd.F90:666
w3adatmd::flcold
logical, pointer flcold
Definition: w3adatmd.F90:688
w3adatmd::hsig
real, dimension(:), pointer hsig
Definition: w3adatmd.F90:587
w3idatmd::flrhoa
logical, pointer flrhoa
Definition: w3idatmd.F90:261
w3adatmd::dal1
real, pointer dal1
Definition: w3adatmd.F90:666
w3adatmd::dcxdy
real, dimension(:,:), pointer dcxdy
Definition: w3adatmd.F90:627
w3adatmd::im24
integer, dimension(:), pointer im24
Definition: w3adatmd.F90:658
w3adatmd::ussx
real, dimension(:), pointer ussx
Definition: w3adatmd.F90:607
w3adatmd::embia2
real, dimension(:), pointer embia2
Definition: w3adatmd.F90:617
w3adatmd::mssx
real, dimension(:), pointer mssx
Definition: w3adatmd.F90:617
w3adatmd::ip24
integer, dimension(:), pointer ip24
Definition: w3adatmd.F90:658
w3adatmd::dal3
real, pointer dal3
Definition: w3adatmd.F90:666
w3adatmd::tauadir
real, dimension(:), pointer tauadir
Definition: w3adatmd.F90:584
w3servmd::extcde
subroutine extcde(IEXIT, UNIT, MSG, FILE, LINE, COMM)
Definition: w3servmd.F90:736
w3adatmd::mscy
real, dimension(:), pointer mscy
Definition: w3adatmd.F90:617
w3adatmd::w3xdma
subroutine w3xdma(IMOD, NDSE, NDST, OUTFLAGS)
Version of W3DIMX for extended ouput arrays only.
Definition: w3adatmd.F90:1523
w3adatmd::ibfloc
integer, pointer ibfloc
Definition: w3adatmd.F90:676
w3adatmd::bstat
integer, dimension(:), pointer bstat
Definition: w3adatmd.F90:680
w3adatmd::nrqsg1
integer, pointer nrqsg1
Definition: w3adatmd.F90:676
w3adatmd::is2
integer, dimension(:), pointer is2
Definition: w3adatmd.F90:635
w3adatmd::mapaxy
integer, dimension(:), pointer mapaxy
Definition: w3adatmd.F90:642
w3adatmd::nadata
integer nadata
Definition: w3adatmd.F90:374
w3adatmd::ww3_field_vec
integer, pointer ww3_field_vec
Definition: w3adatmd.F90:676
w3adatmd::itime
integer, pointer itime
Definition: w3adatmd.F90:686
w3adatmd::ptm1
real, dimension(:,:), pointer ptm1
Definition: w3adatmd.F90:597
w3adatmd::dddx
real, dimension(:,:), pointer dddx
Definition: w3adatmd.F90:627
w3adatmd::wadat
Definition: w3adatmd.F90:381
w3odatmd::noswll
integer, pointer noswll
Definition: w3odatmd.F90:460
w3adatmd::nfrchg
integer, pointer nfrchg
Definition: w3adatmd.F90:657
w3adatmd::irqsg2
integer, dimension(:,:), pointer irqsg2
Definition: w3adatmd.F90:681
w3adatmd::nmy1
integer, pointer nmy1
Definition: w3adatmd.F90:640
w3adatmd::thm
real, dimension(:), pointer thm
Definition: w3adatmd.F90:587
w3adatmd::rai
real, dimension(:), pointer rai
Definition: w3adatmd.F90:578
w3adatmd::mscx
real, dimension(:), pointer mscx
Definition: w3adatmd.F90:617
w3adatmd::ppe
real, dimension(:,:), pointer ppe
Definition: w3adatmd.F90:597
w3adatmd::isploc
integer, pointer isploc
Definition: w3adatmd.F90:676
w3adatmd::nmy2
integer, pointer nmy2
Definition: w3adatmd.F90:640
w3adatmd::cx
real, dimension(:), pointer cx
Definition: w3adatmd.F90:584
w3adatmd::cd0
real, dimension(:), pointer cd0
Definition: w3adatmd.F90:578
w3gdatmd::nx
integer, pointer nx
Definition: w3gdatmd.F90:1097
w3adatmd::ic81
integer, dimension(:), pointer ic81
Definition: w3adatmd.F90:658
w3gdatmd::ngrids
integer ngrids
Definition: w3gdatmd.F90:618
w3adatmd::pgw
real, dimension(:,:), pointer pgw
Definition: w3adatmd.F90:597
w3adatmd::awg1
real, pointer awg1
Definition: w3adatmd.F90:666
w3adatmd::ussp
real, dimension(:,:), pointer ussp
Definition: w3adatmd.F90:612
w3adatmd::nmx1
integer, pointer nmx1
Definition: w3adatmd.F90:640
w3adatmd::nmy0
integer, pointer nmy0
Definition: w3adatmd.F90:640
w3adatmd::ncent
integer, pointer ncent
Definition: w3adatmd.F90:647
w3gdatmd::usspf
integer, dimension(:), pointer usspf
Definition: w3gdatmd.F90:1098
w3adatmd::swg3
real, pointer swg3
Definition: w3adatmd.F90:666
w3adatmd::mssd
real, dimension(:), pointer mssd
Definition: w3adatmd.F90:617
w3gdatmd::e3df
integer, dimension(:,:), pointer e3df
Definition: w3gdatmd.F90:1098
w3adatmd::t0m1
real, dimension(:), pointer t0m1
Definition: w3adatmd.F90:587
w3adatmd::nmxy
integer, pointer nmxy
Definition: w3adatmd.F90:640
w3adatmd::tauocx
real, dimension(:), pointer tauocx
Definition: w3adatmd.F90:607
w3adatmd::aba
real, dimension(:), pointer aba
Definition: w3adatmd.F90:614
w3adatmd::syy
real, dimension(:), pointer syy
Definition: w3adatmd.F90:607
w3adatmd::ubd
real, dimension(:), pointer ubd
Definition: w3adatmd.F90:614
w3adatmd::sxx
real, dimension(:), pointer sxx
Definition: w3adatmd.F90:607
w3adatmd::dhdx
real, dimension(:), pointer dhdx
Definition: w3adatmd.F90:631
w3adatmd::awg4
real, pointer awg4
Definition: w3adatmd.F90:666
w3adatmd::ip23
integer, dimension(:), pointer ip23
Definition: w3adatmd.F90:658
w3adatmd::ic61
integer, dimension(:), pointer ic61
Definition: w3adatmd.F90:658
w3adatmd::mpi_comm_wcmp
integer, pointer mpi_comm_wcmp
Definition: w3adatmd.F90:676