ROMS
Loading...
Searching...
No Matches
checkdefs.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 SUBROUTINE checkdefs
3!
4!git $Id$
5!================================================== Hernan G. Arango ===
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md !
9!=======================================================================
10! !
11! This subroutine checks activated C-preprocessing options for !
12! consistency. !
13! !
14!=======================================================================
15!
16 USE mod_param
17 USE mod_parallel
18 USE mod_iounits
19 USE mod_scalars
20 USE mod_strings
21!
22 USE strings_mod, ONLY : uppercase
23!
24 implicit none
25!
26! Local variable declarations.
27!
28 integer :: iatms, ibbl, ibiology, idriver
29 integer :: ivelHadv, ivelVadv, ivmix, nearshore
30#ifdef SEDIMENT
31 integer :: sbed_load, sbed_type
32#endif
33 integer :: is, lstr, ng
34!
35!-----------------------------------------------------------------------
36! Initialze counters
37!-----------------------------------------------------------------------
38!
39 iatms = 0
40 ibbl = 0
41 ibiology = 0
42 idriver = 0
43 ivelhadv = 0
44 ivelvadv = 0
45 ivmix = 0
46 nearshore = 0
47#ifdef SEDIMENT
48 sbed_load = 0
49 sbed_type = 0
50#endif
51!
52!-----------------------------------------------------------------------
53! Report activated C-preprocessing options.
54!-----------------------------------------------------------------------
55!
56 coptions=' '
57 IF (master) WRITE (stdout,10)
58 10 FORMAT (/,' Activated C-preprocessing Options:',/)
59 20 FORMAT (1x,a,t27,a)
60!
61 IF (master) THEN
62 WRITE (stdout,20) trim(adjustl(myappcpp)), trim(adjustl(title))
63 END IF
64 is=len_trim(coptions)+1
65 lstr=len_trim(myappcpp)
66 coptions(is:is+lstr)=trim(adjustl(myappcpp))
67 is=len_trim(coptions)+1
68 coptions(is:is)=','
69
70#if defined AD_AVERAGES && defined ADJOINT
71!
72 IF (master) WRITE (stdout,20) 'AD_AVERAGES', &
73 & 'Writing out time-averaged adjoint model fields'
74 is=len_trim(coptions)+1
75 coptions(is:is+13)=' AD_AVERAGES,'
76#endif
77#if defined AD_OUTPUT_STATE && defined ADJOINT
78!
79 IF (master) WRITE (stdout,20) 'AD_OUTPUT_STATE', &
80 & 'Writing full adjoint state solution by adding time levels'
81 is=len_trim(coptions)+1
82 coptions(is:is+17)=' AD_OUTPUT_STATE,'
83#endif
84#if defined ADD_FSOBC && defined SSH_TIDES
85!
86 IF (master) WRITE (stdout,20) 'ADD_FSOBC', &
87 & 'Adding tidal elevation to processed OBC data'
88 is=len_trim(coptions)+1
89 coptions(is:is+11)=' ADD_FSOBC,'
90#endif
91#if defined ADD_M2OBC && defined UV_TIDES
92!
93 IF (master) WRITE (stdout,20) 'ADD_M2OBC', &
94 & 'Adding tidal currents to processed OBC data'
95 is=len_trim(coptions)+1
96 coptions(is:is+11)=' ADD_M2OBC,'
97#endif
98#ifdef ADJOINT
99!
100 IF (master) WRITE (stdout,20) 'ADJOINT', &
101 & 'Adjoint Model'
102 is=len_trim(coptions)+1
103 coptions(is:is+9)=' ADJOINT,'
104#endif
105#if defined ADJUST_BOUNDARY && \
106 (defined correlation || \
107 defined i4dvar || \
108 defined i4dvar_ana_sensitivity || \
109 defined rbl4dvar || \
110 defined r4dvar || \
111 defined sensitivity_4dvar || \
112 defined sp4dvar)
113!
114 IF (master) WRITE (stdout,20) 'ADJUST_BOUNDARY', &
115 & 'Including boundary conditions in 4D-Var state estimation'
116 is=len_trim(coptions)+1
117 coptions(is:is+17)=' ADJUST_BOUNDARY,'
118#endif
119#if defined ADJUST_STFLUX && \
120 (defined correlation || \
121 defined i4dvar || \
122 defined i4dvar_ana_sensitivity || \
123 defined rbl4dvar || \
124 defined r4dvar || \
125 defined sensitivity_4dvar || \
126 defined sp4dvar)
127!
128 IF (master) WRITE (stdout,20) 'ADJUST_STFLUX', &
129 & 'Including surface tracer flux in 4D-Var state estimation'
130 is=len_trim(coptions)+1
131 coptions(is:is+15)=' ADJUST_STFLUX,'
132#endif
133#if defined ADJUST_WSTRESS && \
134 (defined correlation || \
135 defined i4dvar || \
136 defined i4dvar_ana_sensitivity || \
137 defined rbl4dvar || \
138 defined r4dvar || \
139 defined sensitivity_4dvar || \
140 defined sp4dvar)
141!
142 IF (master) WRITE (stdout,20) 'ADJUST_WSTRESS', &
143 & 'Including surface wind stress in 4D-Var state estimation'
144 is=len_trim(coptions)+1
145 coptions(is:is+16)=' ADJUST_WSTRESS,'
146#endif
147#if defined AD_IMPULSE
148!
149 IF (master) WRITE (stdout,20) 'AD_IMPULSE', &
150 & 'Force adjoint model with intermittent impulses'
151 is=len_trim(coptions)+1
152 coptions(is:is+12)=' AD_IMPULSE,'
153#endif
154#ifdef ADM_DRIVER
155!
156 IF (master) WRITE (stdout,20) 'ADM_DRIVER', &
157 & 'Generic adjoint model driver'
158 is=len_trim(coptions)+1
159 coptions(is:is+12)=' ADM_DRIVER,'
160 idriver=idriver+1
161#endif
162#ifdef AD_SENSITIVITY
163!
164 IF (master) WRITE (stdout,20) 'AD_SENSITIVITY', &
165 & 'Adjoint Sensitivity Analysis'
166 is=len_trim(coptions)+1
167 coptions(is:is+16)=' AD_SENSITIVITY,'
168 idriver=idriver+1
169#endif
170#ifdef AFT_EIGENMODES
171!
172 IF (master) WRITE (stdout,20) 'AFT_EIGENMODES', &
173 & 'Adjoint Finite Time Eigenvalues'
174 is=len_trim(coptions)+1
175 coptions(is:is+16)=' AFT_EIGENMODES,'
176 idriver=idriver+1
177#endif
178#if defined AGE_MEAN && defined T_PASSIVE && defined SOLVE3D
179!
180 IF (master) WRITE (stdout,20) 'AGE_MEAN', &
181 & 'Computing Mean Age of inert passive tracer'
182 is=len_trim(coptions)+1
183 coptions(is:is+10)=' AGE_MEAN,'
184 idriver=idriver+1
185#endif
186#if defined ALBEDO && defined ANA_SRFLUX
187!
188 IF (master) WRITE (stdout,20) 'ALBEDO', &
189 & 'Shortwave radiation from albedo equation'
190 is=len_trim(coptions)+1
191 coptions(is:is+8)=' ALBEDO,'
192#endif
193#if defined ALLOW_BOTTOM_OBS && \
194 (defined four_dvar || \
195 defined verification) && defined observations
196!
197 IF (master) WRITE (stdout,20) 'ALLOW_BOTTOM_OBS', &
198 & 'Allow observations that lie in the lower bottom grid cell'
199 is=len_trim(coptions)+1
200 coptions(is:is+18)=' ALLOW_BOTTOM_OBS,'
201#endif
202#if defined BIOLOGY && defined ANA_BIOLOGY
203!
204 IF (master) WRITE (stdout,20) 'ANA_BIOLOGY', &
205 & 'Analytical biology initial conditions'
206 is=len_trim(coptions)+1
207 coptions(is:is+13)=' ANA_BIOLOGY,'
208#endif
209#if defined BIOLOGY || defined SEDIMENT || defined T_PASSIVE
210# ifdef ANA_BPFLUX
211!
212 IF (master) WRITE (stdout,20) 'ANA_BPFLUX', &
213 & 'Analytical bottom passive tracers fluxes'
214 is=len_trim(coptions)+1
215 coptions(is:is+12)=' ANA_BPFLUX,'
216# endif
217#endif
218#ifdef ANA_BSFLUX
219!
220 IF (master) WRITE (stdout,20) 'ANA_BSFLUX', &
221 & 'Analytical kinematic bottom salinity flux'
222 is=len_trim(coptions)+1
223 coptions(is:is+12)=' ANA_BSFLUX,'
224#endif
225#ifdef ANA_BTFLUX
226!
227 IF (master) WRITE (stdout,20) 'ANA_BTFLUX', &
228 & 'Analytical kinematic bottom temperature flux'
229 is=len_trim(coptions)+1
230 coptions(is:is+12)=' ANA_BTFLUX,'
231#endif
232#if defined ANA_CLOUD && defined BULK_FLUXES
233!
234 IF (master) WRITE (stdout,20) 'ANA_CLOUD', &
235 & 'Analytical cloud fraction'
236 is=len_trim(coptions)+1
237 coptions(is:is+11)=' ANA_CLOUD,'
238#endif
239#if defined ANA_DRAG && defined UV_DRAG_GRID
240!
241 IF (master) WRITE (stdout,20) 'ANA_DRAG_GRID', &
242# if defined UV_LOGDRAG || defined BBL_MODEL
243 & 'Analytical spatially varying bottom roughness length'
244# elif defined UV_LDRAG
245 & 'Analytical spatially varying linear drag coefficient'
246# elif defined UV_QDRAG
247 & 'Analytical spatially varying quadratic drag coefficient'
248# endif
249 is=len_trim(coptions)+1
250 coptions(is:is+10)=' ANA_DRAG,'
251#endif
252#ifdef ANA_DIAG
253!
254 IF (master) WRITE (stdout,20) 'ANA_DIAG', &
255 & 'Customized diagnostics'
256 is=len_trim(coptions)+1
257 coptions(is:is+10)=' ANA_DIAG,'
258#endif
259#if defined ANA_DQDSST && defined QCORRECTION
260!
261 IF (master) WRITE (stdout,20) 'ANA_DQDSST', &
262 & 'Analytical surface heat flux sensitivity to SST'
263 is=len_trim(coptions)+1
264 coptions(is:is+12)=' ANA_DQDSST,'
265#endif
266#ifdef ANA_FSOBC
267!
268 IF (master) WRITE (stdout,20) 'ANA_FSOBC', &
269 & 'Analytical free-surface boundary conditions'
270 is=len_trim(coptions)+1
271 coptions(is:is+11)=' ANA_FSOBC,'
272#endif
273#ifdef ANA_GRID
274!
275 IF (master) WRITE (stdout,20) 'ANA_GRID', &
276 & 'Analytical grid set-up'
277 is=len_trim(coptions)+1
278 coptions(is:is+10)=' ANA_GRID,'
279#endif
280#if defined ANA_HUMIDITY && defined BULK_FLUXES
281!
282 IF (master) WRITE (stdout,20) 'ANA_HUMIDITY', &
283 & 'Analytical surface air humidity'
284 is=len_trim(coptions)+1
285 coptions(is:is+14)=' ANA_HUMIDITY,'
286#endif
287#ifdef ANA_INITIAL
288!
289 IF (master) WRITE (stdout,20) 'ANA_INITIAL', &
290 & 'Analytical initial conditions'
291 is=len_trim(coptions)+1
292 coptions(is:is+13)=' ANA_INITIAL,'
293#endif
294#ifdef ANA_M2CLIMA
295!
296 IF (master) WRITE (stdout,20) 'ANA_M2CLIMA', &
297 & 'Analytical 2D momentum climatology'
298 is=len_trim(coptions)+1
299 coptions(is:is+13)=' ANA_M2CLIMA,'
300#endif
301#ifdef ANA_M2OBC
302!
303 IF (master) WRITE (stdout,20) 'ANA_M2OBC', &
304 & 'Analytical 2D momentum boundary conditions'
305 is=len_trim(coptions)+1
306 coptions(is:is+11)=' ANA_M2OBC,'
307#endif
308#if defined ANA_M3CLIMA && defined SOLVE3D
309!
310 IF (master) WRITE (stdout,20) 'ANA_M3CLIMA', &
311 & 'Analytical 3D momentum climatology'
312 is=len_trim(coptions)+1
313 coptions(is:is+13)=' ANA_M3CLIMA,'
314#endif
315#if defined ANA_M3OBC && defined SOLVE3D
316!
317 IF (master) WRITE (stdout,20) 'ANA_M3OBC', &
318 & 'Analytical 3D momentum boundary conditions'
319 is=len_trim(coptions)+1
320 coptions(is:is+11)=' ANA_M3OBC,'
321#endif
322#ifdef ANA_MASK
323!
324 IF (master) WRITE (stdout,20) 'ANA_MASK', &
325 & 'Analytical Land/Sea Masking'
326 is=len_trim(coptions)+1
327 coptions(is:is+10)=' ANA_MASK,'
328#endif
329#ifdef ANA_NUDGCOEF
330!
331 IF (any(lnudging(:))) THEN
332 IF (master) WRITE (stdout,20) 'ANA_NUDGCOEF', &
333 & 'Analytical spatially varying nudging time-scales'
334 is=len_trim(coptions)+1
335 coptions(is:is+14)=' ANA_NUDGCOEF,'
336 END IF
337#endif
338#if defined ANA_PAIR && defined BULK_FLUXES
339!
340 IF (master) WRITE (stdout,20) 'ANA_PAIR', &
341 & 'Analytical surface air pressure'
342 is=len_trim(coptions)+1
343 coptions(is:is+10)=' ANA_PAIR,'
344#endif
345#if defined ANA_PASSIVE && defined T_PASSIVE && defined SOLVE3D
346!
347 IF (master) WRITE (stdout,20) 'ANA_PASSIVE', &
348 & 'Analytical initial conditions for inert tracers'
349 is=len_trim(coptions)+1
350 coptions(is:is+13)=' ANA_PASSIVE,'
351#endif
352#if defined ANA_PERTURB
353!
354 IF (master) WRITE (stdout,20) 'ANA_PERTURB', &
355 & 'Perturb initial conditions'
356 is=len_trim(coptions)+1
357 coptions(is:is+13)=' ANA_PERTURB,'
358#endif
359#ifdef ANA_PSOURCE
360!
361 IF (master) WRITE (stdout,20) 'ANA_PSOURCE', &
362 & 'Analytical point sources and sinks'
363 is=len_trim(coptions)+1
364 coptions(is:is+13)=' ANA_PSOURCE,'
365#endif
366#if defined ANA_RAIN && defined BULK_FLUXES
367!
368 IF (master) WRITE (stdout,20) 'ANA_RAIN', &
369 & 'Analytical rain fall rate'
370 is=len_trim(coptions)+1
371 coptions(is:is+10)=' ANA_RAIN,'
372#endif
373#if defined ANA_RESPIRATION && defined HYPOXIA_SRM
374!
375 IF (master) WRITE (stdout,20) 'ANA_RESPIRATION', &
376 & 'Analytical respiration rate'
377 is=len_trim(coptions)+1
378 coptions(is:is+17)=' ANA_RESPIRATION,'
379#endif
380#if defined ANA_SCOPE && \
381 (defined ad_sensitivity || \
382 defined i4dvar_ana_sensitivity || \
383 defined opt_observations || \
384 defined sensitivity_4dvar || \
385 defined so_semi)
386!
387 IF (master) WRITE (stdout,20) 'ANA_SCOPE', &
388 & 'Analytical spatial scope masking arrays'
389 is=len_trim(coptions)+1
390 coptions(is:is+12)=' ANA_SCOPE,'
391#endif
392#if defined ANA_SEDIMENT && (defined SEDIMENT || defined BBL_MODEL)
393!
394 IF (master) WRITE (stdout,20) 'ANA_SEDIMENT', &
395 & 'Analytical sediment initial conditions'
396 is=len_trim(coptions)+1
397 coptions(is:is+14)=' ANA_SEDIMENT,'
398#endif
399#ifdef ANA_SMFLUX
400!
401 IF (master) WRITE (stdout,20) 'ANA_SMFLUX', &
402 & 'Analytical kinematic surface momentum flux'
403 is=len_trim(coptions)+1
404 coptions(is:is+12)=' ANA_SMFLUX,'
405#endif
406#if defined BIOLOGY || defined SEDIMENT || defined T_PASSIVE
407# ifdef ANA_SPFLUX
408!
409 IF (master) WRITE (stdout,20) 'ANA_SPFLUX', &
410 & 'Analytical surface passive tracer fluxes'
411 is=len_trim(coptions)+1
412 coptions(is:is+12)=' ANA_SPFLUX,'
413# endif
414#endif
415#ifdef ANA_SPINNING
416!
417 IF (master) WRITE (stdout,20) 'ANA_SPINNING', &
418 & 'Analytical time-varying rotation force'
419 is=len_trim(coptions)+1
420 coptions(is:is+14)=' ANA_SPINNING,'
421#endif
422#ifdef ANA_SPONGE
423!
424 IF (any(lsponge)) THEN
425 IF (master) WRITE (stdout,20) 'ANA_SPONGE', &
426 & 'Analytical enhanced viscosity/diffusion sponge'
427 is=len_trim(coptions)+1
428 coptions(is:is+12)=' ANA_SPONGE,'
429 END IF
430#endif
431#ifdef ANA_SRFLUX
432!
433 IF (master) WRITE (stdout,20) 'ANA_SRFLUX', &
434 & 'Analytical kinematic shortwave radiation flux'
435 is=len_trim(coptions)+1
436 coptions(is:is+12)=' ANA_SRFLUX,'
437#endif
438#ifdef ANA_SSH
439!
440 IF (master) WRITE (stdout,20) 'ANA_SSH', &
441 & 'Analytical sea surface height climatology'
442 is=len_trim(coptions)+1
443 coptions(is:is+9)=' ANA_SSH,'
444#endif
445#ifdef ANA_SSS
446!
447 IF (master) WRITE (stdout,20) 'ANA_SSS', &
448 & 'Analytical sea surface salinity'
449 is=len_trim(coptions)+1
450 coptions(is:is+9)=' ANA_SSS,'
451#endif
452#if defined ANA_SST && defined QCORRECTION
453!
454 IF (master) WRITE (stdout,20) 'ANA_SST', &
455 & 'Analytical sea surface temperature, SST'
456 is=len_trim(coptions)+1
457 coptions(is:is+9)=' ANA_SST,'
458#endif
459#ifdef ANA_SSFLUX
460!
461 IF (master) WRITE (stdout,20) 'ANA_SSFLUX', &
462 & 'Analytical kinematic surface salinity flux'
463 is=len_trim(coptions)+1
464 coptions(is:is+12)=' ANA_SSFLUX,'
465#endif
466#ifdef ANA_STFLUX
467!
468 IF (master) WRITE (stdout,20) 'ANA_STFLUX', &
469 & 'Analytical kinematic surface temperature flux'
470 is=len_trim(coptions)+1
471 coptions(is:is+12)=' ANA_STFLUX,'
472#endif
473#if defined ANA_TCLIMA && defined SOLVE3D
474!
475 IF (master) WRITE (stdout,20) 'ANA_TCLIMA', &
476 & 'Analytical tracer climatology'
477 is=len_trim(coptions)+1
478 coptions(is:is+12)=' ANA_TCLIMA,'
479#endif
480#if defined ANA_TOBC && defined SOLVE3D
481!
482 IF (master) WRITE (stdout,20) 'ANA_TOBC', &
483 & 'Analytical tracers boundary conditions'
484 is=len_trim(coptions)+1
485 coptions(is:is+10)=' ANA_TOBC,'
486#endif
487#if defined ANA_VMIX && defined SOLVE3D
488!
489 IF (master) WRITE (stdout,20) 'ANA_VMIX', &
490 & 'Analytical vertical mixing coefficients'
491 is=len_trim(coptions)+1
492 coptions(is:is+10)=' ANA_VMIX,'
493 ivmix=ivmix+1
494#endif
495#if defined ANA_WINDS && defined BULK_FLUXES
496!
497 IF (master) WRITE (stdout,20) 'ANA_WINDS', &
498 & 'Analytical surface wind components'
499 is=len_trim(coptions)+1
500 coptions(is:is+11)=' ANA_WINDS,'
501#endif
502#if defined ANA_WTYPE && defined WTYPE_GRID && \
503 (defined lmd_skpp || defined solar_source)
504!
505 IF (master) WRITE (stdout,20) 'ANA_WTYPE', &
506 & 'Analytical spatially varying Jerlov water type index'
507 is=len_trim(coptions)+1
508 coptions(is:is+11)=' ANA_WTYPE,'
509#endif
510#ifdef ANA_WWAVE
511!
512 IF (master) WRITE (stdout,20) 'ANA_WWAVE', &
513 & 'Analytical wind induced waves'
514 is=len_trim(coptions)+1
515 coptions(is:is+11)=' ANA_WWAVE,'
516#endif
517#ifdef ARRAY_MODES
518!
519 IF (master) WRITE (stdout,20) 'ARRAY_MODES', &
520 & 'Stabilized Representer Matrix Array Modes'
521 is=len_trim(coptions)+1
522 coptions(is:is+13)=' ARRAY_MODES,'
523 idriver=idriver+1
524# ifdef ARRAY_MODES_SPLIT
525!
526 IF (master) WRITE (stdout,20) 'ARRAY_MODES_SPLIT', &
527 & 'Separate analysis due to IC, surface forcing, and OBC'
528 is=len_trim(coptions)+1
529 coptions(is:is+19)=' ARRAY_MODES_SPLIT,'
530# endif
531#endif
532#if defined DISTRIBUTE && defined NESTING
533# if defined ASSEMBLE_ALLGATHER
534!
535 IF (master) WRITE (stdout,20) 'ASSEMBLE_ALLGATHER', &
536 & 'Using mpi_allgather in mp_assemble routine'
537 is=len_trim(coptions)+1
538 coptions(is:is+20)=' ASSEMBLE_ALLGATHER,'
539# elif defined ASSEMBLE_ALLREDUCE
540!
541 IF (master) WRITE (stdout,20) 'ASSEMBLE_ALLREDUCE', &
542 & 'Using mpi_allreduce in mp_assemble routine'
543 is=len_trim(coptions)+1
544 coptions(is:is+20)=' ASSEMBLE_ALLGATHER,'
545# elif defined ASSEMBLE_SENDRECV
546!
547 IF (master) WRITE (stdout,20) 'ASSEMBLE_SENDRECV', &
548 & 'Using mpi_isend/mpi_recv in mp_assemble routine'
549 is=len_trim(coptions)+1
550 coptions(is:is+19)=' ASSEMBLE_SENDRECV,'
551# endif
552#endif
553#ifdef ASSUMED_SHAPE
554!
555 IF (master) WRITE (stdout,20) 'ASSUMED_SHAPE', &
556 & 'Using assumed-shape arrays'
557 is=len_trim(coptions)+1
558 coptions(is:is+15)=' ASSUMED_SHAPE,'
559#endif
560#if defined ASYNCHRONOUS_PIO && defined PIO_LIB && defined DISTRIBUTE
561!
562 IF (master) WRITE (stdout,20) 'ASYNCHRONOUS_PIO', &
563 & 'PIO library with disjointed dedicated I/O processes'
564 is=len_trim(coptions)+1
565 coptions(is:is+18)=' ASYNCHRONOUS_PIO,'
566#endif
567#if defined ASYNCHRONOUS_SCORPIO && defined PIO_LIB && defined DISTRIBUTE
568!
569 IF (master) WRITE (stdout,20) 'ASYNCHRONOUS_SCORPIO', &
570 & 'SCORPIO library with disjointed dedicated I/O processes'
571 is=len_trim(coptions)+1
572 coptions(is:is+22)=' ASYNCHRONOUS_SCORPIO,'
573#endif
574#ifdef ATM_PRESS
575!
576 IF (master) WRITE (stdout,20) 'ATM_PRESS', &
577 & 'Impose inverse barometer effect in pressure gradient term'
578 is=len_trim(coptions)+1
579 coptions(is:is+11)=' ATM_PRESS,'
580#endif
581#ifdef AVERAGES
582!
583 IF (master) WRITE (stdout,20) 'AVERAGES', &
584 & 'Writing out time-averaged nonlinear model fields'
585 is=len_trim(coptions)+1
586 coptions(is:is+10)=' AVERAGES,'
587# if defined AVERAGES_DETIDE && (defined SSH_TIDES || defined UV_TIDES)
588!
589 IF (master) WRITE (stdout,20) 'AVERAGES_DETIDE', &
590 & 'Writing out time-averaged nonlinear model detided fields'
591 is=len_trim(coptions)+1
592 coptions(is:is+17)=' AVERAGES_DETIDE,'
593# endif
594# ifdef AVERAGES_UV_KIRBY
595 IF (master) WRITE (stdout,20) 'AVERAGES_UV_KIRBY', &
596 & 'Writing out time-averaged uwave and vwave Kirby velocities'
597 is=len_trim(coptions)+1
598 coptions(is:is+17)=' AVERAGES_UV_KIRBY,'
599# endif
600# ifdef AVERAGES_WEC
601 IF (master) WRITE (stdout,20) 'AVERAGES_WEC', &
602 & 'Writing out time-averaged nearshore radiation fields'
603 is=len_trim(coptions)+1
604 coptions(is:is+20)=' AVERAGES_WEC,'
605# endif
606#endif
607#if defined BACKGROUND && defined FOUR_DVAR
608!
609 IF (master) WRITE (stdout,20) 'BACKGROUND', &
610 & 'Include background cost function'
611 is=len_trim(coptions)+1
612 coptions(is:is+12)=' BACKGROUND,'
613#endif
614#ifdef BALANCE_OPERATOR
615!
616 IF (master) WRITE (stdout,20) 'BALANCE_OPERATOR', &
617 & 'Error Covariance Multivariate Balance Operator'
618 is=len_trim(coptions)+1
619 coptions(is:is+18)=' BALANCE_OPERATOR,'
620#endif
621#if defined SEDIMENT && defined BEDLOAD_MPM
622!
623 IF (master) WRITE (stdout,20) 'BEDLOAD_MPM', &
624 & 'Activate bed load sediment transport Meyer-Peter-Mueller'
625 is=len_trim(coptions)+1
626 coptions(is:is+13)=' BEDLOAD_MPM,'
627 sbed_load=sbed_load+1
628#endif
629#if defined SEDIMENT && defined BEDLOAD_SOULSBY
630!
631 IF (master) WRITE (stdout,20) 'BEDLOAD_SOULSBY', &
632 & 'Activate bed load sediment transport Soulsby formula'
633 is=len_trim(coptions)+1
634 coptions(is:is+17)=' BEDLOAD_SOULSBY,'
635#endif
636#if defined SEDIMENT && defined BEDLOAD_VANDERA
637!
638 IF (master) WRITE (stdout,20) 'BEDLOAD_VANDERA', &
639 & 'Activate bed load sediment transport vanDerA formula'
640 is=len_trim(coptions)+1
641 coptions(is:is+17)=' BEDLOAD_VANDERA,'
642#endif
643#if defined SEDIMENT && defined BEDLOAD_VANDERA_ASYM_LIMITS
644!
645 IF (master) WRITE (stdout,20) 'BEDLOAD_VANDERA_ASYM_LIMITS', &
646 & 'Activate bed load vanDerA asymmetry limits'
647 is=len_trim(coptions)+1
648 coptions(is:is+29)=' BEDLOAD_VANDERA_ASYM_LIMITS,'
649#endif
650#if defined SEDIMENT && defined BEDLOAD_VANDERA_SURFACE_WAVE
651!
652 IF (master) WRITE (stdout,20) 'BEDLOAD_VANDERA_SURFACE_WAVE', &
653 & 'Activate bed load vanDerA surface wave'
654 is=len_trim(coptions)+1
655 coptions(is:is+30)=' BEDLOAD_VANDERA_SURFACE_WAVE,'
656#endif
657#if defined SEDIMENT && defined BEDLOAD_VANDERA_WAVE_AVGD_STRESS
658!
659 IF (master) WRITE (stdout,20) 'BEDLOAD_VANDERA_WAVE_AVGD_STRESS', &
660 & 'Activate bed load vanDerA wave avg stress'
661 is=len_trim(coptions)+1
662 coptions(is:is+34)=' BEDLOAD_VANDERA_WAVE_AVGD_STRESS,'
663#endif
664#if defined SEDIMENT && defined BEDLOAD_VANDERA_MADSEN_UDELTA
665!
666 IF (master) WRITE (stdout,20) 'BEDLOAD_VANDERA_MADSEN_UDELTA', &
667 & 'Activate bed load madsen current calculation'
668 is=len_trim(coptions)+1
669 coptions(is:is+31)=' BEDLOAD_VANDERA_MADSEN_UDELTA,'
670#endif
671#if defined SEDIMENT && defined BEDLOAD_VANDERA_DIRECT_UDELTA
672!
673 IF (master) WRITE (stdout,20) 'BEDLOAD_VANDERA_DIRECT_UDELTA', &
674 & 'Activate bed load direct elevation current calculation'
675 is=len_trim(coptions)+1
676 coptions(is:is+30)=' BEDLOAD_VANDERA_DIRECT_UDELTA,'
677#endif
678#if defined BGQC && defined FOUR_DVAR
679!
680 IF (master) WRITE (stdout,20) 'BGQC', &
681 & 'Background quality control of observations'
682 is=len_trim(coptions)+1
683 coptions(is:is+6)=' BGQC,'
684#endif
685#if defined BEOFS_ONLY && defined WEAK_CONSTRAINT
686!
687 IF (master) WRITE (stdout,20) 'BEOFS_ONLY', &
688 & 'Compute EOFs of the background error covariance matrix'
689 is=len_trim(coptions)+1
690 coptions(is:is+12)=' BEOFS_ONLY,'
691#endif
692#ifdef BIO_FENNEL
693!
694 IF (master) WRITE (stdout,20) 'BIO_FENNEL', &
695 & 'Fennel et al. (2006) nitrogen-based model'
696 is=len_trim(coptions)+1
697 coptions(is:is+12)=' BIO_FENNEL,'
698 ibiology=ibiology+1
699#endif
700#if defined BNORM && defined HESSIAN_SV
701!
702 IF (master) WRITE (stdout,20) 'BNORM', &
703 & 'Background normalization of Hessian singular vectors'
704 is=len_trim(coptions)+1
705 coptions(is:is+7)=' BNORM,'
706#endif
707#if defined BIO_OPTIC && defined ECOSIM
708!
709 IF (master) WRITE (stdout,20) 'BIO_OPTIC', &
710 & 'Underwater light field calculations'
711 is=len_trim(coptions)+1
712 coptions(is:is+11)=' BIO_OPTIC,'
713#endif
714#if defined BIO_SEDIMENT && \
715 (defined bio_fennel || defined nemuro || defined ecosim)
716!
717 IF (master) WRITE (stdout,20) 'BIO_SEDIMENT', &
718 & 'Restore fallen particulate material to the nutrient pool'
719 is=len_trim(coptions)+1
720 coptions(is:is+14)=' BIO_SEDIMENT,'
721#endif
722#ifdef BODYFORCE
723!
724 IF (master) WRITE (stdout,20) 'BODYFORCE', &
725 & 'Momentum stresses as body-forces'
726 is=len_trim(coptions)+1
727 coptions(is:is+11)=' BODYFORCE,'
728#endif
729#ifdef DISTRIBUTE
730# if defined BOUNDARY_ALLGATHER
731!
732 IF (master) WRITE (stdout,20) 'BOUNDARY_ALLGATHER', &
733 & 'Using mpi_allgather in mp_boundary routine'
734 is=len_trim(coptions)+1
735 coptions(is:is+20)=' BOUNDARY_ALLGATHER,'
736# elif defined BOUNDARY_ALLREDUCE
737 IF (master) WRITE (stdout,20) 'BOUNDARY_ALLREDUCE', &
738 & 'Using mpi_allreduce in mp_boundary routine'
739 is=len_trim(coptions)+1
740 coptions(is:is+20)=' BOUNDARY_ALLREDUCE,'
741# endif
742#endif
743#if defined BOTTOM_STREAMING && defined WEC_VF
744 IF (master) WRITE (stdout,20) 'BOTTOM_STREAMING', &
745 & 'Wave current streaming term'
746 is=len_trim(coptions)+1
747 coptions(is:is+18)=' BOTTOM_STREAMING,'
748#endif
749#ifdef BULK_FLUXES
750!
751 IF (master) WRITE (stdout,20) 'BULK_FLUXES', &
752 & 'Surface bulk fluxes parameterization'
753 is=len_trim(coptions)+1
754 coptions(is:is+13)=' BULK_FLUXES,'
755#endif
756#ifdef BVF_MIXING
757!
758 IF (master) WRITE (stdout,20) 'BVF_MIXING', &
759 & 'Brunt-Vaisala frequency based vertical mixing'
760 is=len_trim(coptions)+1
761 coptions(is:is+12)=' BVF_MIXING,'
762 ivmix=ivmix+1
763#endif
764#if defined CANUTO_A && defined GLS_MIXING
765!
766 IF (master) WRITE (stdout,20) 'CANUTO_A', &
767 & 'Canuto A-stability function formulation'
768 is=len_trim(coptions)+1
769 coptions(is:is+8)=' CANUTO_A,'
770#endif
771#if defined CANUTO_B && defined GLS_MIXING
772!
773 IF (master) WRITE (stdout,20) 'CANUTO_B', &
774 & 'Canuto B-stability function formulation'
775 is=len_trim(coptions)+1
776 coptions(is:is+8)=' CANUTO_B,'
777#endif
778#if defined CARBON && defined BIO_FENNEL
779!
780 IF (master) WRITE (stdout,20) 'CARBON', &
781 & 'Add Carbon constituents'
782 is=len_trim(coptions)+1
783 coptions(is:is+8)=' CARBON,'
784#endif
785#if defined CELERITY_READ && defined FORWARD_READ
786!
787 IF (master) WRITE (stdout,20) 'CELERITY_READ', &
788 & 'Read in radiation celerity from forward file'
789 is=len_trim(coptions)+1
790 coptions(is:is+15)=' CELERETY_READ,'
791#endif
792#if defined CELERITY_WRITE && defined FORWARD_WRITE
793!
794 IF (master) WRITE (stdout,20) 'CELERITY_WRITE', &
795 & 'Write out radiation celerity in forward file'
796 is=len_trim(coptions)+1
797 coptions(is:is+16)=' CELERITY_WRITE,'
798#endif
799#ifdef CHECKSUM
800!
801 IF (master) WRITE (stdout,20) 'CHECKSUM', &
802 & 'Report order-invariant checksum (hash) when processing I/O'
803 is=len_trim(coptions)+1
804 coptions(is:is+10)=' CHECKSUM,'
805#endif
806#ifdef CHECK_OPEN_FILES
807!
808 IF (master) WRITE (stdout,20) 'CHECK_OPEN_FILES', &
809 & 'Monitor the number opened, closed, and created NetCDF files'
810 is=len_trim(coptions)+1
811 coptions(is:is+18)=' CHECK_OPEN_FILES,'
812#endif
813#ifdef CLIPPING
814!
815 IF (master) WRITE (stdout,20) 'CLIPPING', &
816 & 'Stabilized Representer Matrix clipping spectrum analysis'
817 is=len_trim(coptions)+1
818 coptions(is:is+10)=' CLIPPING,'
819 idriver=idriver+1
820# ifdef CLIPPING_SPLIT
821!
822 IF (master) WRITE (stdout,20) 'CLIPPING_SPLIT', &
823 & 'Separate analysis due to IC, surface forcing, and OBC'
824 is=len_trim(coptions)+1
825 coptions(is:is+16)=' CLIPPING_SPLIT,'
826# endif
827#endif
828#if defined CMEPS && defined ESMF_LIB
829!
830 IF (master) WRITE (stdout,20) 'CMEPS', &
831 & 'Using Community Mediator for Earth Prediction Systems'
832 is=len_trim(coptions)+1
833 coptions(is:is+7)=' CMEPS,'
834#endif
835#ifdef COAMPS_COUPLING
836!
837 IF (master) WRITE (stdout,20) 'COAMPS_COUPLING', &
838 & 'Atmosphere coupling with COAMPS'
839 is=len_trim(coptions)+1
840 coptions(is:is+17)=' COAMPS_COUPLING,'
841 iatms=iatms+1
842#endif
843#if defined DISTRIBUTE
844# if defined COLLECT_ALLGATHER
845!
846 IF (master) WRITE (stdout,20) 'COLLECT_ALLGATHER', &
847 & 'Using mpi_allgather in mp_collect routine'
848 is=len_trim(coptions)+1
849 coptions(is:is+19)=' COLLECT_ALLGATHER,'
850# elif defined COLLECT_ALLREDUCE
851!
852 IF (master) WRITE (stdout,20) 'COLLECT_ALLREDUCE', &
853 & 'Using mpi_allreduce in mp_collect routine'
854 is=len_trim(coptions)+1
855 coptions(is:is+19)=' COLLECT_ALLGATHER,'
856# elif defined COLLECT_SENDRECV
857!
858 IF (master) WRITE (stdout,20) 'COLLECT_SENDRECV', &
859 & 'Using mpi_isend/mpi_recv in mp_collect routine'
860 is=len_trim(coptions)+1
861 coptions(is:is+18)=' COLLECT_SENDRECV,'
862# endif
863#endif
864#if defined COMPUTE_MLD && defined STD_MODEL && defined FOUR_DVAR
865!
866 IF (master) WRITE (stdout,20) 'COMPUTE_MLD', &
867 & 'Compute Mixed Layer Depth for background error modeling'
868 is=len_trim(coptions)+1
869 coptions(is:is+13)=' COMPUTE_MLD,'
870#endif
871#if defined CONCURRENT_KERNEL && defined DISJOINTED && \
872 defined distribute && defined sp4dvar
873!
874 IF (master) WRITE (stdout,20) 'CONCURRENT_KERNEL', &
875 & 'Concurrent TLM and ADM in saddle-point 4D-Var'
876 is=len_trim(coptions)+1
877 coptions(is:is+19)=' CONCURRENT_KERNEL,'
878#endif
879#ifdef CORRELATION
880!
881 IF (master) WRITE (stdout,20) 'CORRELATION', &
882 & 'Background-error Correlation Model'
883 is=len_trim(coptions)+1
884 coptions(is:is+13)=' CORRELATION,'
885 idriver=idriver+1
886#endif
887#ifdef DATA_COUPLING
888!
889 IF (master) WRITE (stdout,20) 'DATA_COUPLING', &
890 & 'DATA model included in coupling'
891 is=len_trim(coptions)+1
892 coptions(is:is+15)=' DATA_COUPLING,'
893#endif
894#if defined DATALESS_LOOPS && defined R4DVAR
895!
896 IF (master) WRITE (stdout,20) 'DATALESS_LOOPS', &
897 & 'Testing convergence of RPM Picard iterations'
898 is=len_trim(coptions)+1
899 coptions(is:is+16)=' DATALESS_LOOPS,'
900#endif
901#if defined GLS_MIXING && defined CHARNOK
902!
903 IF (master) WRITE (stdout,20) 'CHARNOK', &
904 & 'Charnock surface roughness from wind stress'
905 is=len_trim(coptions)+1
906 coptions(is:is+9)=' CHARNOK,'
907#endif
908#if defined PROPAGATOR && defined CHECKPOINTING
909!
910 IF (master) WRITE (stdout,20) 'CHECKPOINTING', &
911 & 'Processing checkpointing NetCDF for GST analysis'
912 is=len_trim(coptions)+1
913 coptions(is:is+15)=' CHECKPOINTING,'
914#endif
915#if defined GLS_MIXING && defined CRAIG_BANNER
916!
917 IF (master) WRITE (stdout,20) 'CRAIG_BANNER', &
918 & 'Craig and Banner wave breaking surface flux'
919 is=len_trim(coptions)+1
920 coptions(is:is+14)=' CRAIG_BANNER,'
921#endif
922#if defined COARE_OOST
923!
924 IF (master) WRITE (stdout,20) 'COARE_OOST', &
925 & 'Oost et al (2002) relation for ZoW in bulk fluxes'
926 is=len_trim(coptions)+1
927 coptions(is:is+12)=' COARE_OOST,'
928#endif
929#if defined COARE_TAYLOR_YELLAND
930!
931 IF (master) WRITE (stdout,20) 'COARE_TAYLOR_YELLAND', &
932 & 'Taylor and Yelland (2001) relation for ZoW in bulk fluxes'
933 is=len_trim(coptions)+1
934 coptions(is:is+22)=' COARE_TAYLOR_YELLAND,'
935#endif
936#if defined COOL_SKIN && defined BULK_FLUXES
937!
938 IF (master) WRITE (stdout,20) 'COOL_SKIN', &
939 & 'Surface cool skin correction'
940 is=len_trim(coptions)+1
941 coptions(is:is+11)=' COOL_SKIN,'
942#endif
943#if defined COSINE2 && defined SOLVE3D
944!
945 IF (master) WRITE (stdout,20) 'COSINE2', &
946 & 'Cosine-squared shape time-averaging barotropic filter'
947 is=len_trim(coptions)+1
948 coptions(is:is+9)=' COSINE2,'
949#endif
950#ifdef CURVGRID
951!
952 IF (master) WRITE (stdout,20) 'CURVGRID', &
953 & 'Orthogonal curvilinear grid'
954 is=len_trim(coptions)+1
955 coptions(is:is+10)=' CURVGRID,'
956#endif
957#if defined DAILY_SHORTWAVE && defined RED_TIDE
958!
959 IF (master) WRITE (stdout,20) 'DAILY_SHORTWAVE', &
960 & 'Daily-averaged surface solar shortwave radiation flux'
961 is=len_trim(coptions)+1
962 coptions(is:is+17)=' DAILY_SHORTWAVE,'
963#endif
964#ifdef DEBUGGING
965!
966 IF (master) WRITE (stdout,20) 'DEBUGGING', &
967 & 'Internal debugging switch activated'
968 is=len_trim(coptions)+1
969 coptions(is:is+11)=' DEBUGGING,'
970#endif
971#ifdef DEEPWATER_WAVES
972!
973 IF (master) WRITE (stdout,20) 'DEEPWATER_WAVES', &
974 & 'Deep water waves approx in bulk fluxes'
975 is=len_trim(coptions)+1
976 coptions(is:is+17)=' DEEPWATER_WAVES,'
977#endif
978#if defined DEFLATE && defined HDF5
979!
980 IF (master) WRITE (stdout,20) 'DEFLATE', &
981 & 'Setting compression in output NetCDF-4/HDF5 files'
982 is=len_trim(coptions)+1
983 coptions(is:is+9)=' DEFLATE,'
984#endif
985#if defined DENITRIFICATION && defined BIO_FENNEL
986!
987 IF (master) WRITE (stdout,20) 'DENITRIFICATION', &
988 & 'Add denitrification processes'
989 is=len_trim(coptions)+1
990 coptions(is:is+17)=' DENITRIFICATION,'
991#endif
992#if defined DIAGNOSTICS_BIO && defined BIO_FENNEL
993!
994 IF (master) WRITE (stdout,20) 'DIAGNOSTICS_BIO', &
995 & 'Computing and writing biological diagnostic terms'
996 is=len_trim(coptions)+1
997 coptions(is:is+17)=' DIAGNOSTICS_BIO,'
998#endif
999#ifdef DIAGNOSTICS_TS
1000!
1001 IF (master) WRITE (stdout,20) 'DIAGNOSTICS_TS', &
1002 & 'Computing and writing tracer diagnostic terms'
1003 is=len_trim(coptions)+1
1004 coptions(is:is+16)=' DIAGNOSTICS_TS,'
1005#endif
1006#ifdef DIAGNOSTICS_UV
1007!
1008 IF (master) WRITE (stdout,20) 'DIAGNOSTICS_UV', &
1009 & 'Computing and writing momentum diagnostic terms'
1010 is=len_trim(coptions)+1
1011 coptions(is:is+16)=' DIAGNOSTICS_UV,'
1012#endif
1013#if defined DIFF_3DCOEF
1014!
1015 IF (master) WRITE (stdout,20) 'DIFF_3DCOEF', &
1016 & 'Horizontal, time-dependent 3D diffusion coefficient'
1017 is=len_trim(coptions)+1
1018 coptions(is:is+13)=' DIFF_3DCOEF,'
1019#endif
1020#if defined TS_DIF2 || defined TS_DIF4
1021# ifdef DIFF_GRID
1022!
1023 IF (master) WRITE (stdout,20) 'DIFF_GRID', &
1024 & 'Horizontal diffusion coefficient scaled by grid size'
1025 is=len_trim(coptions)+1
1026 coptions(is:is+11)=' DIFF_GRID,'
1027# endif
1028#endif
1029#if defined DISJOINTED && defined DISTRIBUTE
1030!
1031 IF (master) WRITE (stdout,20) 'DISJOINTED', &
1032 & 'Split mpi-communicator into disjointed subgroups'
1033 is=len_trim(coptions)+1
1034 coptions(is:is+12)=' DISJOINTED,'
1035#endif
1036#ifdef DIURNAL_SRFLUX
1037!
1038 IF (master) WRITE (stdout,20) 'DIURNAL_SRFLUX', &
1039 & 'Modulate shortwave radiation by the local diurnal cycle'
1040 is=len_trim(coptions)+1
1041 coptions(is:is+16)=' DIURNAL_SRFLUX,'
1042#endif
1043#ifdef DJ_GRADPS
1044!
1045 IF (master) WRITE (stdout,20) 'DJ_GRADPS', &
1046 & 'Parabolic Splines density Jacobian (Shchepetkin, 2002)'
1047 is=len_trim(coptions)+1
1048 coptions(is:is+11)=' DJ_GRADPS,'
1049#endif
1050#ifdef DOUBLE_PRECISION
1051!
1052 IF (master) WRITE (stdout,20) 'DOUBLE_PRECISION', &
1053 & 'Double precision arithmetic numerical kernel'
1054 is=len_trim(coptions)+1
1055 coptions(is:is+18)=' DOUBLE_PRECISION,'
1056#endif
1057#if defined DRENNAN
1058 IF (master) WRITE (stdout,20) 'DRENNAN', &
1059 & 'Drennan et al (2003) relation for ZoW in bulk fluxes'
1060 is=len_trim(coptions)+1
1061 coptions(is:is+9)=' DRENNAN,'
1062#endif
1063#ifdef ECOSIM
1064!
1065 IF (master) WRITE (stdout,20) 'ECOSIM', &
1066 & 'Bio-Optical EcoSim model'
1067 is=len_trim(coptions)+1
1068 coptions(is:is+8)=' ECOSIM,'
1069 ibiology=ibiology+1
1070#endif
1071#if defined EMINUSP && defined BULK_FLUXES
1072!
1073 IF (master) WRITE (stdout,20) 'EMINUSP', &
1074 & 'Compute Salt Flux using E-P'
1075 is=len_trim(coptions)+1
1076 coptions(is:is+9)=' EMINUSP,'
1077#endif
1078#if defined ESMF_LIB && defined MODEL_COUPLING
1079!
1080 IF (master) WRITE (stdout,20) 'ESMF_LIB', &
1081 & 'Using Earth System Modeling Framework library'
1082 is=len_trim(coptions)+1
1083 coptions(is:is+10)=' ESMF_LIB,'
1084#endif
1085#ifdef ENKF_RESTART
1086!
1087 IF (master) WRITE (stdout,20) 'ENKF_RESTART', &
1088 & 'Writting fields for offline Ensemble Kalman Filter'
1089 is=len_trim(coptions)+1
1090 coptions(is:is+14)=' ENKF_RESTART,'
1091#endif
1092#ifdef ENSEMBLE
1093!
1094 IF (master) WRITE (stdout,20) 'ENSEMBLE', &
1095 & 'Ensemble Forecasting Propagator'
1096 is=len_trim(coptions)+1
1097 coptions(is:is+10)=' ENSEMBLE,'
1098#endif
1099#if defined ESM_SETRUNCLOCK && defined MODEL_COUPLING && \
1100 defined esmf_lib
1101!
1102 IF (master) WRITE (stdout,20) 'ESM_SETRUNCLOCK', &
1103 & 'Sets coupled component specialized SetRunClock method'
1104 is=len_trim(coptions)+1
1105 coptions(is:is+17)=' ESM_SETRUNCLOCK,'
1106#endif
1107#if defined EVOLVED_LCZ && defined FOUR_DVAR
1108!
1109 IF (master) WRITE (stdout,20) 'EVOLVED_LCZ', &
1110 & 'Compute 4D-Var evolved Hessian singular vectors'
1111 is=len_trim(coptions)+1
1112 coptions(is:is+13)=' EVOLVED_LCZ,'
1113#endif
1114# if defined EXCLUDE_SPONGE && defined MODEL_COUPLING && \
1115 (defined data_coupling && !defined ANA_SPONGE)
1116!
1117 IF (master) WRITE (stdout,20) 'EXCLUDE_SPONGE', &
1118 & 'Exclude sponge points in export fields during coupling'
1119 is=len_trim(coptions)+1
1120 coptions(is:is+16)=' EXCLUDE_SPONGE,'
1121#endif
1122#ifdef FLOATS
1123!
1124 IF (master) WRITE (stdout,20) 'FLOATS', &
1125 & 'Simulated Lagrangian drifters'
1126 is=len_trim(coptions)+1
1127 coptions(is:is+8)=' FLOATS,'
1128#endif
1129#if defined FLOAT_BIOLOGY && defined FLOATS
1130!
1131 IF (master) WRITE (stdout,20) 'FLOAT_BIOLOGY', &
1132 & 'Biological behavior on Lagrangian drifters'
1133 is=len_trim(coptions)+1
1134 coptions(is:is+15)=' FLOAT_BIOLOGY,'
1135#endif
1136#if defined FLOAT_OYSTER && defined FLOATS
1137!
1138 IF (master) WRITE (stdout,20) 'FLOAT_OYSTER', &
1139 & 'Oyster model behavior on Lagrangian drifters'
1140 is=len_trim(coptions)+1
1141 coptions(is:is+14)=' FLOAT_OYSTER,'
1142#endif
1143#if defined FLOAT_STICKY && defined FLOATS && defined SOLVE3D
1144!
1145 IF (master) WRITE (stdout,20) 'FLOAT_STICKY', &
1146 & 'Reflect (stick) floats that hit the surface (bottom)'
1147 is=len_trim(coptions)+1
1148 coptions(is:is+14)=' FLOAT_STICKY,'
1149#endif
1150#if defined FLOATS && defined FLOAT_VWALK && defined SOLVE3D
1151!
1152 IF (master) WRITE (stdout,20) 'FLOAT_VWALK', &
1153 & 'Lagrangian drifters with vertical diffusion, random walk'
1154 is=len_trim(coptions)+1
1155 coptions(is:is+13)=' FLOAT_VWALK,'
1156#endif
1157#ifdef FORCING_SV
1158!
1159 IF (master) WRITE (stdout,20) 'FORCING_SV', &
1160 & 'Forcing Singular Vectors Propagator'
1161 is=len_trim(coptions)+1
1162 coptions(is:is+12)=' FORCING_SV,'
1163 idriver=idriver+1
1164#endif
1165#if defined FORWARD_FLUXES
1166!
1167 IF (master) WRITE (stdout,20) 'FORWARD_FLUXES', &
1168 & 'Using nonlinear model Forward fluxes for Tangent/Adjoint'
1169 is=len_trim(coptions)+1
1170 coptions(is:is+16)=' FORWARD_FLUXES,'
1171#endif
1172#if defined FORWARD_MIXING && defined SOLVE3D
1173!
1174 IF (master) WRITE (stdout,20) 'FORWARD_MIXING', &
1175 & 'Read in Forward vertical mixing for Tangent/Adjoint'
1176 is=len_trim(coptions)+1
1177 coptions(is:is+16)=' FORWARD_MIXING,'
1178#endif
1179#ifdef FORWARD_READ
1180!
1181 IF (master) WRITE (stdout,20) 'FORWARD_READ', &
1182 & 'Read in Forward solution for Tangent/Adjoint'
1183 is=len_trim(coptions)+1
1184 coptions(is:is+14)=' FORWARD_READ,'
1185#endif
1186#ifdef FORWARD_RHS
1187!
1188 IF (master) WRITE (stdout,20) 'FORWARD_RHS', &
1189 & 'Process Forward RHS terms for Tangent/Adjoint'
1190 is=len_trim(coptions)+1
1191 coptions(is:is+13)=' FORWARD_RHS,'
1192#endif
1193#ifdef FORWARD_WRITE
1194!
1195 IF (master) WRITE (stdout,20) 'FORWARD_WRITE', &
1196 & 'Write out Forward solution for Tangent/Adjoint'
1197 is=len_trim(coptions)+1
1198 coptions(is:is+15)=' FORWARD_WRITE,'
1199#endif
1200#if defined FRC_COUPLING && defined MODEL_COUPLING
1201!
1202 IF (master) WRITE (stdout,20) 'FRC_COUPLING', &
1203 & 'Atmospheric forcing fields are from ESM coupling'
1204 is=len_trim(coptions)+1
1205 coptions(is:is+14)=' FRC_COUPLING,'
1206#endif
1207#ifdef FSOBC_REDUCED
1208!
1209 IF (master) WRITE (stdout,20) 'FSOBC_REDUCED', &
1210 & 'Using free-surface data in reduced physics conditions'
1211 is=len_trim(coptions)+1
1212 coptions(is:is+15)=' FSOBC_REDUCED,'
1213#endif
1214#ifdef FT_EIGENMODES
1215!
1216 IF (master) WRITE (stdout,20) 'FT_EIGENMODES', &
1217 & 'Finite Time Eigenmodes: Normal Modes'
1218 is=len_trim(coptions)+1
1219 coptions(is:is+15)=' FT_EIGENMODES,'
1220 idriver=idriver+1
1221#endif
1222#if defined FOUR_DVAR || defined PROPAGATOR || defined VERIFICATION
1223# ifdef FULL_GRID
1224!
1225 IF (master) WRITE (stdout,20) 'FULL_GRID', &
1226# if (defined FOUR_DVAR || defined VERIFICATION)
1227 & 'Considering observations at interior and boundary points'
1228# elif defined PROPAGATOR
1229 & 'State vector includes interior and boundary points'
1230# endif
1231 is=len_trim(coptions)+1
1232 coptions(is:is+11)=' FULL_GRID,'
1233# else
1234!
1235 IF (master) WRITE (stdout,20) '!FULL_GRID', &
1236# if (defined FOUR_DVAR || defined VERIFICATION)
1237 & 'Considering observations at interior points only'
1238# elif defined PROPAGATOR
1239 & 'State vector includes interior points only'
1240# endif
1241 is=len_trim(coptions)+1
1242 coptions(is:is+12)=' !FULL_GRID,'
1243# endif
1244#endif
1245#ifdef DISTRIBUTE
1246# ifdef GATHER_SENDRECV
1247 IF (master) WRITE (stdout,20) 'GATHER_SENDRECV', &
1248 & 'Using mpi_isend/mpi_irecv in mp_gather2d/mp_gather3d routines'
1249 is=len_trim(coptions)+1
1250 coptions(is:is+17)=' GATHER_SENDRECV,'
1251# else
1252 IF (master) WRITE (stdout,20) '!GATHER_SENDRECV', &
1253 & 'Using mpi_gatherv in mp_gather2d/mp_gather3d routines'
1254 is=len_trim(coptions)+1
1255 coptions(is:is+18)=' !GATHER_SENDRECV,'
1256# endif
1257#endif
1258#ifdef GENERIC_DSTART
1259!
1260 IF (master) WRITE (stdout,20) 'GENERIC_DSTART', &
1261 & 'Avoid using ''dstart'' as initial time in ADM, TLM, and RPM'
1262 is=len_trim(coptions)+1
1263 coptions(is:is+16)=' GENERIC_DSTART,'
1264#endif
1265#if defined GEOPOTENTIAL_HCONV && defined FOUR_DVAR && defined SOLVE3D
1266!
1267 IF (master) WRITE (stdout,20) 'GEOPOTENTIAL_HCONV', &
1268 & 'Horizontal convolutions along geopetential surfaces'
1269 is=len_trim(coptions)+1
1270 coptions(is:is+20)=' GEOPOTENTIAL_HCONV,'
1271#endif
1272#ifdef GLS_MIXING
1273!
1274 IF (master) WRITE (stdout,20) 'GLS_MIXING', &
1275 & 'Generic Length-Scale turbulence closure'
1276 is=len_trim(coptions)+1
1277 coptions(is:is+12)=' GLS_MIXING,'
1278 ivmix=ivmix+1
1279#endif
1280#ifdef GRID_EXTRACT
1281!
1282 IF (master) WRITE (stdout,20) 'GRID_EXTRACT', &
1283 & 'Writing output extraction history file'
1284 is=len_trim(coptions)+1
1285 coptions(is:is+14)=' GRID_EXTRACT,'
1286#endif
1287#ifdef HDF5
1288!
1289 IF (master) WRITE (stdout,20) 'HDF5', &
1290 & 'Creating NetCDF-4/HDF5 format files'
1291 is=len_trim(coptions)+1
1292 coptions(is:is+6)=' HDF5,'
1293#endif
1294#ifdef HESSIAN_FSV
1295!
1296 IF (master) WRITE (stdout,20) 'HESSIAN_FSV', &
1297 & 'Hessian (4D-Var) Forcing Singular Vectors Propagator'
1298 is=len_trim(coptions)+1
1299 coptions(is:is+13)=' HESSIAN_FSV,'
1300#endif
1301#ifdef HESSIAN_SO
1302!
1303 IF (master) WRITE (stdout,20) 'HESSIAN_SO', &
1304 & 'Hessian (4D-Var) Stochastic Optimals Propagator'
1305 is=len_trim(coptions)+1
1306 coptions(is:is+12)=' HESSIAN_SO,'
1307#endif
1308#ifdef HESSIAN_SV
1309!
1310 IF (master) WRITE (stdout,20) 'HESSIAN_SV', &
1311 & 'Hessian (4D-Var) Singular Vectors Propagator'
1312 is=len_trim(coptions)+1
1313 coptions(is:is+12)=' HESSIAN_SV,'
1314#endif
1315#if defined HOLLING_GRAZING && defined NEMURO
1316!
1317 IF (master) WRITE (stdout,20) 'HOLLING_GRAZING', &
1318 & 'Holling-type s-shaped grazing algorithm, Implicit'
1319 is=len_trim(coptions)+1
1320 coptions(is:is+17)=' HOLLING_GRAZING,'
1321#endif
1322#ifdef HYPOXIA_SRM
1323!
1324 IF (master) WRITE (stdout,20) 'HYPOXIA_SRM', &
1325 & 'Hypoxia Simple Respiration Model'
1326 is=len_trim(coptions)+1
1327 coptions(is:is+13)=' HYPOXIA_SRM,'
1328 ibiology=ibiology+1
1329#endif
1330#if defined ICE_MODEL && defined SOLVE3D
1331!
1332 IF (master) WRITE (stdout,20) 'ICE_MODEL', &
1333 & 'Native Sea-Ice Model (Budgell; Durski and Kurapov)'
1334 is=len_trim(coptions)+1
1335 coptions(is:is+11)=' ICE_MODEL,'
1336# ifdef ICE_ADVECT
1337!
1338 IF (master) WRITE (stdout,20) 'ICE_ADVECT', &
1339 & 'Including advection of ice tracers'
1340 is=len_trim(coptions)+1
1341 coptions(is:is+12)=' ICE_ADVECT,'
1342# endif
1343# ifdef ICE_ALB_EC92
1344!
1345 IF (master) WRITE (stdout,20) 'ICE_ALB_EC92', &
1346 & 'Albedo computation from Ebert and Curry'
1347 is=len_trim(coptions)+1
1348 coptions(is:is+14)=' ICE_ALB_EC92,'
1349# endif
1350# ifdef ICE_BULK_FLUXES
1351!
1352 IF (master) WRITE (stdout,20) 'ICE_BULK_FLUXES', &
1353 & 'Including ice terms in the bulk flux computation'
1354 is=len_trim(coptions)+1
1355 coptions(is:is+17)=' ICE_BULK_FLUXES,'
1356# endif
1357# ifdef ICE_CONVSNOW
1358!
1359 IF (master) WRITE (stdout,20) 'ICE_CONVSNOW', &
1360 & 'Conversion of flooded snow to ice'
1361 is=len_trim(coptions)+1
1362 coptions(is:is+14)=' ICE_CONVSNOW,'
1363# endif
1364# if defined ICE_EVP && defined ICE_MOMENTUM
1365!
1366 IF (master) WRITE (stdout,20) 'ICE_EVP', &
1367 & 'Elastic-viscous-plastic rheology'
1368 is=len_trim(coptions)+1
1369 coptions(is:is+9)=' ICE_EVP,'
1370# endif
1371# ifdef ICE_MK
1372!
1373 IF (master) WRITE (stdout,20) 'ICE_MK', &
1374 & 'Mellor-Kantha thermodynamics'
1375 is=len_trim(coptions)+1
1376 coptions(is:is+8)=' ICE_MK,'
1377# endif
1378# ifdef ICE_MOMENTUM
1379!
1380 IF (master) WRITE (stdout,20) 'ICE_MOMENTUM', &
1381 & 'Including momentum term'
1382 is=len_trim(coptions)+1
1383 coptions(is:is+14)=' ICE_MOMENTUM,'
1384# endif
1385# ifdef ICE_MOM_BULK
1386!
1387 IF (master) WRITE (stdout,20) 'ICE_MOM_BULK', &
1388 & 'Alternated ice-water stress computation'
1389 is=len_trim(coptions)+1
1390 coptions(is:is+14)=' ICE_MOM_BULK,'
1391# endif
1392# ifdef ICE_STRENGTH_QUAD
1393!
1394 IF (master) WRITE (stdout,20) 'ICE_STRENGTH_QUAD', &
1395 & 'Quadratic ice strength, a function of thickness'
1396 is=len_trim(coptions)+1
1397 coptions(is:is+19)=' ICE_STRENGTH_QUAD,'
1398# endif
1399# if defined ICE_SMOLAR && defined ICE_ADVECT
1400!
1401 IF (master) WRITE (stdout,20) 'ICE_SMOLAR', &
1402 & 'MPDATA advection scheme'
1403 is=len_trim(coptions)+1
1404 coptions(is:is+12)=' ICE_SMOLAR,'
1405# endif
1406# if defined ICE_UPWIND && defined ICE_ADVECT
1407!
1408 IF (master) WRITE (stdout,20) 'ICE_UPWIND', &
1409 & 'Upwind advection scheme'
1410 is=len_trim(coptions)+1
1411 coptions(is:is+12)=' ICE_UPWIND,'
1412# endif
1413#endif
1414#ifdef ICESHELF
1415!
1416 IF (master) WRITE (stdout,20) 'ICESHELF', &
1417 & 'Include Ice Shelf Cavities'
1418 is=len_trim(coptions)+1
1419 coptions(is:is+10)=' ICESHELF,'
1420#endif
1421#if defined IMPACT_INNER && \
1422 defined obs_impact && \
1423 (defined i4dvar_ana_sensitivity || \
1424 defined rbl4dvar_ana_sensitivity || \
1425 defined rbl4dvar_fct_sensitivity || \
1426 defined r4dvar_ana_sensitivity)
1427!
1428 IF (master) WRITE (stdout,20) 'IMPACT_INNER', &
1429 & 'Writing observations impacts for each inner loop'
1430 is=len_trim(coptions)+1
1431 coptions(is:is+14)=' IMPACT_INNER,'
1432#endif
1433#ifdef IMPLICIT_NUDGING
1434!
1435 IF (master) WRITE (stdout,20) 'IMPLICIT_NUDGING', &
1436 & 'Implicit nudging term in momentum radiation conditions'
1437 is=len_trim(coptions)+1
1438 coptions(is:is+18)=' IMPLICIT_NUDGING,'
1439#endif
1440#if defined IMPLICIT_VCONV && defined FOUR_DVAR && defined VCONVOLUTION
1441!
1442 IF (master) WRITE (stdout,20) 'IMPLICIT_VCONV', &
1443 & 'Implicit Vertical Convolution Algorithm'
1444 is=len_trim(coptions)+1
1445 coptions(is:is+16)=' IMPLICIT_VCONV,'
1446#endif
1447#ifdef INITIALIZE_AUTOMATIC
1448!
1449 IF (master) WRITE (stdout,20) 'INITIALIZE_AUTOMATIC', &
1450 & 'Initialize automatic arrays in step2d for debugging'
1451 is=len_trim(coptions)+1
1452 coptions(is:is+22)=' INITIALIZE_AUTOMATIC,'
1453#endif
1454#ifdef JEDI
1455 IF (master) WRITE (stdout,20) 'JEDI', &
1456 & 'Joint Effort for Data Assimilation Integration'
1457 is=len_trim(coptions)+1
1458 coptions(is:is+6)=' JEDI,'
1459#endif
1460#if defined LIMIT_BSTRESS && defined SOLVE3D && !defined BBL_MODEL
1461!
1462 IF (master) WRITE (stdout,20) 'LIMIT_BSTRESS', &
1463 & 'Limit bottom stress to maintain bottom velocity direction'
1464 is=len_trim(coptions)+1
1465 coptions(is:is+15)=' LIMIT_BSTRESS,'
1466#endif
1467#if defined LIMIT_STFLX_COOLING && defined SOLVE3D
1468!
1469 IF (master) WRITE (stdout,20) 'LIMIT_STFLX_COOLING', &
1470 & 'Suppress further cooling if SST is at freezing point'
1471 is=len_trim(coptions)+1
1472 coptions(is:is+21)=' LIMIT_STFLX_COOLING,'
1473#endif
1474#if defined NPZD_IRON && defined IRON_LIMIT
1475!
1476 IF (master) WRITE (stdout,20) 'IRON_LIMIT', &
1477 & 'Include dissolved iron limitation on phytoplankton growth'
1478 is=len_trim(coptions)+1
1479 coptions(is:is+12)=' IRON_LIMIT,'
1480#endif
1481#if defined NPZD_IRON && defined IRON_LIMIT && defined IRON_RELAX
1482!
1483 IF (master) WRITE (stdout,20) 'IRON_RELAX', &
1484 & 'Nudge dissolved iron over the shelf, h <= FeHmin'
1485 is=len_trim(coptions)+1
1486 coptions(is:is+12)=' IRON_RELAX,'
1487#endif
1488#ifdef IMPULSE
1489!
1490 IF (master) WRITE (stdout,20) 'IMPULSE', &
1491 & 'Processing Adjoint Impulse forcing'
1492 is=len_trim(coptions)+1
1493 coptions(is:is+9)=' IMPULSE,'
1494#endif
1495#ifdef INNER_PRODUCT
1496!
1497 IF (master) WRITE (stdout,20) 'INNER_PRODUCT', &
1498 & 'Tangent/Adjoint State Matrices Inner Product Test'
1499 is=len_trim(coptions)+1
1500 coptions(is:is+15)=' INNER_PRODUCT,'
1501 idriver=idriver+1
1502#endif
1503#if defined I4DVAR && !defined SPLIT_4DVAR
1504!
1505 IF (master) WRITE (stdout,20) 'I4DVAR', &
1506 & 'Incremental, strong constraint 4D-Var data assimilation'
1507 is=len_trim(coptions)+1
1508 coptions(is:is+8)=' I4DVAR,'
1509 idriver=idriver+1
1510#endif
1511#ifdef I4DVAR_ANA_SENSITIVITY
1512!
1513 IF (master) WRITE (stdout,20) 'I4DVAR_ANA_SENSITIVITY', &
1514 & 'I4D-Var Observations Sensitivity Analysis'
1515 is=len_trim(coptions)+1
1516 coptions(is:is+24)=' I4DVAR_ANA_SENSITIVITY,'
1517 idriver=idriver+1
1518#endif
1519#if defined INLINE_2DIO && defined DISTRIBUTE
1520!
1521 IF (master) WRITE (stdout,20) 'INLINE_2DIO', &
1522 & 'Processing 3D IO level by level to reduce memory needs'
1523 is=len_trim(coptions)+1
1524 coptions(is:is+13)=' INLINE_2DIO,'
1525#endif
1526#if defined IVLEV_EXPLICIT && defined NEMURO
1527!
1528 IF (master) WRITE (stdout,20) 'IVLEV_EXPLICIT', &
1529 & 'Ivlev Grazing algorithm, Explicit'
1530 is=len_trim(coptions)+1
1531 coptions(is:is+16)=' IVLEV_EXPLICIT,'
1532#endif
1533#if defined KANTHA_CLAYSON && (defined GLS_MIXING || defined MY25_MIXING)
1534!
1535 IF (master) WRITE (stdout,20) 'KANTHA_CLAYSON', &
1536 & 'Kantha and Clayson stability function formulation'
1537 is=len_trim(coptions)+1
1538 coptions(is:is+16)=' KANTHA_CLAYSON,'
1539#endif
1540#if defined LCZ_FINAL && defined FOUR_DVAR
1541!
1542 IF (master) WRITE (stdout,20) 'LCZ_FINAL', &
1543 & 'Compute 4D-Var Hessian singular vectors'
1544 is=len_trim(coptions)+1
1545 coptions(is:is+11)=' LCZ_FINAL,'
1546#endif
1547#if defined LIMIT_VDIFF && \
1548 (defined gls_mixing || defined lmd_mixing || defined my25_mixing)
1549!
1550 IF (master) WRITE (stdout,20) 'LIMIT_VDIFF', &
1551 & 'Impose an upper limit on vertical diffusion coefficient'
1552 is=len_trim(coptions)+1
1553 coptions(is:is+13)=' LIMIT_VDIFF,'
1554#endif
1555#if defined LIMIT_VVISC && \
1556 (defined gls_mixing || defined lmd_mixing || defined my25_mixing)
1557!
1558 IF (master) WRITE (stdout,20) 'LIMIT_VVISC', &
1559 & 'Impose an upper limit on vertical viscosity coefficient'
1560 is=len_trim(coptions)+1
1561 coptions(is:is+13)=' LIMIT_VVISC,'
1562#endif
1563#ifdef LMD_BKPP
1564!
1565 IF (master) WRITE (stdout,20) 'LMD_BKPP', &
1566 & 'KPP bottom boundary layer mixing'
1567 is=len_trim(coptions)+1
1568 coptions(is:is+10)=' LMD_BKPP,'
1569#endif
1570#ifdef LMD_CONVEC
1571!
1572 IF (master) WRITE (stdout,20) 'LMD_CONVEC', &
1573 & 'LMD convective mixing due to shear instability'
1574 is=len_trim(coptions)+1
1575 coptions(is:is+12)=' LMD_CONVEC,'
1576#endif
1577#ifdef LMD_DDMIX
1578!
1579 IF (master) WRITE (stdout,20) 'LMD_DDMIX', &
1580 & 'LMD double-diffusive mixing'
1581 is=len_trim(coptions)+1
1582 coptions(is:is+11)=' LMD_DDMIX,'
1583#endif
1584#ifdef LMD_MIXING
1585!
1586 IF (master) WRITE (stdout,20) 'LMD_MIXING', &
1587 & 'Large/McWilliams/Doney interior mixing'
1588 is=len_trim(coptions)+1
1589 coptions(is:is+12)=' LMD_MIXING,'
1590 ivmix=ivmix+1
1591#endif
1592#ifdef LMD_NONLOCAL
1593!
1594 IF (master) WRITE (stdout,20) 'LMD_NONLOCAL', &
1595 & 'LMD convective nonlocal transport'
1596 is=len_trim(coptions)+1
1597 coptions(is:is+14)=' LMD_NONLOCAL,'
1598#endif
1599#ifdef LMD_RIMIX
1600!
1601 IF (master) WRITE (stdout,20) 'LMD_RIMIX', &
1602 & 'LMD diffusivity due to shear instability'
1603 is=len_trim(coptions)+1
1604 coptions(is:is+11)=' LMD_RIMIX,'
1605#endif
1606#ifdef LMD_SHAPIRO
1607!
1608 IF (master) WRITE (stdout,20) 'LMD_SHAPIRO', &
1609 & 'Shapiro filtering boundary layer depth'
1610 is=len_trim(coptions)+1
1611 coptions(is:is+13)=' LMD_SHAPIRO,'
1612#endif
1613#ifdef LMD_SKPP
1614!
1615 IF (master) WRITE (stdout,20) 'LMD_SKPP', &
1616 & 'KPP surface boundary layer mixing'
1617 is=len_trim(coptions)+1
1618 coptions(is:is+10)=' LMD_SKPP,'
1619#endif
1620#if defined LONGWAVE && defined BULK_FLUXES
1621!
1622 IF (master) WRITE (stdout,20) 'LONGWAVE', &
1623 & 'Compute net longwave radiation internally'
1624 is=len_trim(coptions)+1
1625 coptions(is:is+10)=' LONGWAVE,'
1626#endif
1627#if defined LONGWAVE_OUT && defined BULK_FLUXES
1628!
1629 IF (master) WRITE (stdout,20) 'LONGWAVE_OUT', &
1630 & 'Compute outgoing longwave radiation internally'
1631 is=len_trim(coptions)+1
1632 coptions(is:is+14)=' LONGWAVE_OUT,'
1633#endif
1634#ifdef MASKING
1635!
1636 IF (master) WRITE (stdout,20) 'MASKING', &
1637 & 'Land/Sea masking'
1638 is=len_trim(coptions)+1
1639 coptions(is:is+9)=' MASKING,'
1640#endif
1641#ifdef MB_BBL
1642!
1643 IF (master) WRITE (stdout,20) 'MB_BBL', &
1644 & 'Blaas Bottom Boundary Layer'
1645 is=len_trim(coptions)+1
1646 coptions(is:is+8)=' MB_BBL,'
1647 ibbl=ibbl+1
1648# ifdef MB_CALC_UB
1649!
1650 IF (master) WRITE (stdout,20) 'MB_CALC_UB', &
1651 & 'Internal computation of bottom orbital velocity'
1652 is=len_trim(coptions)+1
1653 coptions(is:is+12)=' MB_CALC_UB,'
1654# endif
1655# ifdef MB_CALC_ZNOT
1656!
1657 IF (master) WRITE (stdout,20) 'MB_CALC_ZNOT', &
1658 & 'Internal computation of bottom roughness'
1659 is=len_trim(coptions)+1
1660 coptions(is:is+14)=' MB_CALC_ZNOT,'
1661# endif
1662# ifdef MB_Z0BIO
1663!
1664 IF (master) WRITE (stdout,20) 'MB_Z0BIO', &
1665 & 'Biogenic bedform roughness ripple height and length'
1666 is=len_trim(coptions)+1
1667 coptions(is:is+10)=' MB_Z0BIO,'
1668# endif
1669# ifdef MB_Z0BL
1670!
1671 IF (master) WRITE (stdout,20) 'MB_Z0BL', &
1672 & 'Bedload roughness for ripples'
1673 is=len_trim(coptions)+1
1674 coptions(is:is+9)=' MB_Z0BL,'
1675# endif
1676# ifdef MB_Z0RIP
1677!
1678 IF (master) WRITE (stdout,20) 'MB_Z0RIP', &
1679 & 'Bedform roughness ripple height and length'
1680 is=len_trim(coptions)+1
1681 coptions(is:is+10)=' MB_Z0RIP,'
1682# endif
1683#endif
1684#if defined MCT_LIB && defined MODEL_COUPLING
1685!
1686 IF (master) WRITE (stdout,20) 'MCT_LIB', &
1687 & 'Using Model Coupling Toolkit library'
1688 is=len_trim(coptions)+1
1689 coptions(is:is+9)=' MCT_LIB,'
1690#endif
1691#ifdef METADATA_REPORT
1692!
1693 IF (master) WRITE (stdout,20) 'METADATA_REPORT', &
1694 & 'Reports/dumps YAML metadata dictionary'
1695 is=len_trim(coptions)+1
1696 coptions(is:is+10)=' METADATA_REPORT,'
1697#endif
1698#if defined MINRES && defined FOUR_DVAR
1699!
1700 IF (master) WRITE (stdout,20) 'MINRES', &
1701 & 'Minimal Residual Method for 4D-Var minimization'
1702 is=len_trim(coptions)+1
1703 coptions(is:is+8)=' MINRES,'
1704#endif
1705#if (defined TS_DIF2 || defined TS_DIF4) && defined SOLVE3D
1706# ifdef MIX_GEO_TS
1707!
1708 IF (master) WRITE (stdout,20) 'MIX_GEO_TS', &
1709 & 'Mixing of tracers along geopotential surfaces'
1710 is=len_trim(coptions)+1
1711 coptions(is:is+12)=' MIX_GEO_TS,'
1712# endif
1713# ifdef MIX_ISO_TS
1714!
1715 IF (master) WRITE (stdout,20) 'MIX_ISO_TS', &
1716 & 'Mixing of tracers along isopycnal surfaces'
1717 is=len_trim(coptions)+1
1718 coptions(is:is+12)=' MIX_ISO_TS,'
1719# endif
1720# ifdef MIX_S_TS
1721!
1722 IF (master) WRITE (stdout,20) 'MIX_S_TS', &
1723 & 'Mixing of tracers along constant S-surfaces'
1724 is=len_trim(coptions)+1
1725 coptions(is:is+10)=' MIX_S_TS,'
1726# endif
1727#endif
1728#if (defined UV_VIS2 || defined UV_VIS4) && defined SOLVE3D
1729# ifdef MIX_GEO_UV
1730!
1731 IF (master) WRITE (stdout,20) 'MIX_GEO_UV', &
1732 & 'Mixing of momentum along geopotential surfaces'
1733 is=len_trim(coptions)+1
1734 coptions(is:is+11)=' MIX_GEO_UV,'
1735# endif
1736# ifdef MIX_S_UV
1737!
1738 IF (master) WRITE (stdout,20) 'MIX_S_UV', &
1739 & 'Mixing of momentum along constant S-surfaces'
1740 is=len_trim(coptions)+1
1741 coptions(is:is+10)=' MIX_S_UV,'
1742# endif
1743#endif
1744#ifdef MPI
1745!
1746 IF (master) WRITE (stdout,20) 'MPI', &
1747 & 'MPI distributed-memory configuration'
1748 is=len_trim(coptions)+1
1749 coptions(is:is+4)=' MPI,'
1750#endif
1751#if defined DISTRIBUTE
1752# if defined MULTI_THREAD
1753!
1754 IF (master) WRITE (stdout,20) 'MULTI_THREAD', &
1755 & 'Using mpi_init_thread with multiple thread level processing'
1756 is=len_trim(coptions)+1
1757 coptions(is:is+17)=' MULTIPLE_THREAD,'
1758# else
1759!
1760 IF (master) WRITE (stdout,20) '!MULTI_THREAD', &
1761 & 'Using mpi_init with single thread level processing'
1762 is=len_trim(coptions)+1
1763 coptions(is:is+18)=' !MULTIPLE_THREAD,'
1764# endif
1765#endif
1766#if defined MULTIPLE_TLM && defined TANGENT
1767!
1768 IF (master) WRITE (stdout,20) 'MULTIPLE_TLM', &
1769 & 'Create multiple TLM history files'
1770 is=len_trim(coptions)+1
1771 coptions(is:is+14)=' MULTIPLE_TLM,'
1772#endif
1773#ifdef MY25_MIXING
1774!
1775 IF (master) WRITE (stdout,20) 'MY25_MIXING', &
1776 & 'Mellor/Yamada Level-2.5 mixing closure'
1777 is=len_trim(coptions)+1
1778 coptions(is:is+13)=' MY25_MIXING,'
1779 ivmix=ivmix+1
1780#endif
1781#ifdef NEMURO
1782!
1783 IF (master) WRITE (stdout,20) 'NEMURO', &
1784 & 'Nemuro Lower Trophic Level Ecosystem Model'
1785 is=len_trim(coptions)+1
1786 coptions(is:is+8)=' NEMURO,'
1787 ibiology=ibiology+1
1788#endif
1789#ifdef NESTING
1790!
1791 IF (master) WRITE (stdout,20) 'NESTING', &
1792 & 'Nesting grids: Composite and Refinement'
1793 is=len_trim(coptions)+1
1794 coptions(is:is+9)=' NESTING,'
1795# ifdef NESTING_DEBUG
1796!
1797 IF (master) WRITE (stdout,20) 'NESTING_DEBUG', &
1798 & 'Checking mass fluxes conservation in refinement'
1799 is=len_trim(coptions)+1
1800 coptions(is:is+15)=' NESTING_DEBUG,'
1801# endif
1802#endif
1803#if defined NLM_OUTER && defined WEAK_CONSTRAINT
1804!
1805 IF (master) WRITE (stdout,20) 'NLM_OUTER', &
1806 & 'Using the Nonlinear model as Basic State in Outer Loop'
1807 is=len_trim(coptions)+1
1808 coptions(is:is+11)=' NLM_OUTER,'
1809#endif
1810#ifdef NONLINEAR
1811!
1812 IF (master) WRITE (stdout,20) 'NONLINEAR', &
1813 & 'Nonlinear Model'
1814 is=len_trim(coptions)+1
1815 coptions(is:is+10)=' NONLINEAR,'
1816#endif
1817#if defined NO_SCORRECTION_ICE && defined ICE_MODEL && defined SOLVE3D
1818!
1819 IF (master) WRITE (stdout,20) 'NO_SCORRECTION_ICE', &
1820 & 'No salinity correction under the ice'
1821 is=len_trim(coptions)+1
1822 coptions(is:is+20)=' NO_SCORRECTION_ICE,'
1823#endif
1824#ifdef SOLVE3D
1825# if defined NONLIN_EOS
1826!
1827 IF (master) WRITE (stdout,20) 'NONLIN_EOS', &
1828 & 'Nonlinear Equation of State for seawater'
1829 is=len_trim(coptions)+1
1830 coptions(is:is+12)=' NONLIN_EOS,'
1831# else
1832!
1833 IF (master) WRITE (stdout,20) '!NONLIN_EOS', &
1834 & 'Linear Equation of State for seawater'
1835 is=len_trim(coptions)+1
1836 coptions(is:is+13)=' !NONLIN_EOS,'
1837# endif
1838#endif
1839#if defined NO_CORRECT_TRACER && defined NESTING && defined SOLVE3D
1840!
1841 IF (master) WRITE (stdout,20) 'NO_CORRECT_TRACER', &
1842 & 'Avoid two-way correction of coaser grid boundary tracer'
1843 is=len_trim(coptions)+1
1844 coptions(is:is+19)=' NO_CORRECT_TRACER,'
1845#endif
1846#ifdef NO_LBC_ATT
1847!
1848 IF (master) WRITE (stdout,20) 'NO_LBC_ATT', &
1849 & 'Not checking NetCDF global attribute NLM_LBC during restart'
1850 is=len_trim(coptions)+1
1851 coptions(is:is+12)=' NO_LBC_ATT,'
1852#endif
1853#if defined NO_READ_GHOST && defined DISTRIBUTE
1854!
1855 IF (master) WRITE (stdout,20) 'NO_READ_GHOST', &
1856 & 'Not including ghost points during scattering read data'
1857 is=len_trim(coptions)+1
1858 coptions(is:is+15)=' NO_READ_GHOST,'
1859#endif
1860#ifdef NO_WRITE_GRID
1861!
1862 IF (master) WRITE (stdout,20) 'NO_WRITE_GRID', &
1863 & 'Not Writing grid arrays into NetCDF ouput files'
1864 is=len_trim(coptions)+1
1865 coptions(is:is+15)=' NO_WRITE_GRID,'
1866#endif
1867#ifdef NPZD_FRANKS
1868!
1869 IF (master) WRITE (stdout,20) 'NPZD_FRANKS', &
1870 & 'NPZD Biological Model, Franks et al. fomulation'
1871 is=len_trim(coptions)+1
1872 coptions(is:is+13)=' NPZD_FRANKS,'
1873 ibiology=ibiology+1
1874#endif
1875#ifdef NPZD_IRON
1876!
1877 IF (master) WRITE (stdout,20) 'NPZD_IRON', &
1878 & 'NPZD Biological Model (Powell et al.) with iron limitation'
1879 is=len_trim(coptions)+1
1880 coptions(is:is+11)=' NPZD_IRON,'
1881 ibiology=ibiology+1
1882#endif
1883#ifdef NPZD_POWELL
1884!
1885 IF (master) WRITE (stdout,20) 'NPZD_POWELL', &
1886 & 'NPZD Biological Model, Powell et al. fomulation'
1887 is=len_trim(coptions)+1
1888 coptions(is:is+13)=' NPZD_POWELL,'
1889 ibiology=ibiology+1
1890#endif
1891#if defined N2S2_HORAVG && (defined GLS_MIXING || defined MY25_MIXING)
1892!
1893 IF (master) WRITE (stdout,20) 'N2S2_HORAVG', &
1894 & 'Horizontal smoothing of buoyancy and shear'
1895 is=len_trim(coptions)+1
1896 coptions(is:is+13)=' N2S2_HORAVG,'
1897#endif
1898#ifdef OBSERVATIONS
1899!
1900 IF (master) WRITE (stdout,20) 'OBSERVATIONS', &
1901 & 'Processing 4D-Var observations'
1902 is=len_trim(coptions)+1
1903 coptions(is:is+13)=' OBSERVATIONS,'
1904#endif
1905#ifdef OBS_IMPACT
1906!
1907 IF (master) WRITE (stdout,20) 'OBS_IMPACT', &
1908 & 'Compute observation impact to the 4D-Var system'
1909 is=len_trim(coptions)+1
1910 coptions(is:is+12)=' OBS_IMPACT,'
1911#endif
1912#if defined OBS_IMPACT && defined OBS_IMPACT_SPLIT
1913!
1914 IF (master) WRITE (stdout,20) 'OBS_IMPACT_SPLIT', &
1915 & 'Separate impact due to IC, surface forcing, and OBC'
1916 is=len_trim(coptions)+1
1917 coptions(is:is+18)=' OBS_IMPACT_SPLIT,'
1918#endif
1919#if defined OBS_SPACE && defined RBL4DVAR_FCT_SENSITIVITY
1920!
1921 IF (master) WRITE (stdout,20) 'OBS_SPACE', &
1922 & 'Forecast metric in observation-space instead of state-space'
1923 is=len_trim(coptions)+1
1924 coptions(is:is+11)=' OBS_SPACE,'
1925#endif
1926#if defined OMEGA_IMPLICIT && defined SOLVE3D
1927!
1928 IF (master) WRITE (stdout,20) 'OMEGA_IMPLICIT', &
1929 & 'Adaptive, Courant-number based implicit vertical advection'
1930 is=len_trim(coptions)+1
1931 coptions(is:is+16)=' OMEGA_IMPLICIT,'
1932#endif
1933#ifdef NESTING
1934# ifdef ONE_WAY
1935!
1936 IF (master) WRITE (stdout,20) 'ONE_WAY', &
1937 & 'One-way nesting in refinement grids'
1938 is=len_trim(coptions)+1
1939 coptions(is:is+9)=' ONE_WAY,'
1940# else
1941!
1942 IF (master) WRITE (stdout,20) '!ONE_WAY', &
1943 & 'Two-way nesting in refinement grids'
1944 is=len_trim(coptions)+1
1945 coptions(is:is+10)=' !ONE_WAY,'
1946# endif
1947#endif
1948#ifdef OUT_DOUBLE
1949!
1950 IF (master) WRITE (stdout,20) 'OUT_DOUBLE', &
1951 & 'Double precision output fields in NetCDF files'
1952 is=len_trim(coptions)+1
1953 coptions(is:is+12)=' OUT_DOUBLE,'
1954#endif
1955#ifdef OUTPUT_STATS
1956!
1957 IF (master) WRITE (stdout,20) 'OUTPUT_STATS', &
1958 & 'Reporting NetCDF output fields statistics'
1959 is=len_trim(coptions)+1
1960 coptions(is:is+14)=' OUTPUT_STATS,'
1961#endif
1962#ifdef _OPENMP
1963!
1964 IF (master) WRITE (stdout,20) '_OPENMP', &
1965 & 'OpenMP parallel shared-memory directives'
1966 is=len_trim(coptions)+1
1967 coptions(is:is+9)=' _OPENMP,'
1968#endif
1969#ifdef OPT_OBSERVATIONS
1970!
1971 IF (master) WRITE (stdout,20) 'OPT_OBSERVATIONS', &
1972 & 'Optimal Observations Driver'
1973 is=len_trim(coptions)+1
1974 coptions(is:is+18)=' OPT_OBSERVATIONS,'
1975 idriver=idriver+1
1976#endif
1977#ifdef OPT_PERTURBATION
1978!
1979 IF (master) WRITE (stdout,20) 'OPT_PERTURBATION', &
1980 & 'Optimal Perturbation Propagator'
1981 is=len_trim(coptions)+1
1982 coptions(is:is+18)=' OPT_PERTURBATION,'
1983 idriver=idriver+1
1984#endif
1985#if defined OXYGEN && defined BIO_FENNEL
1986# ifdef OCMIP_OXYGEN_SC
1987!
1988 IF (master) WRITE (stdout,20) 'OCMIP_OXYGEN_SC', &
1989 & 'Schmidt number from Keeling et al. (1998)'
1990 is=len_trim(coptions)+1
1991 coptions(is:is+17)=' OCMIP_OXYGEN_SC,'
1992# else
1993!
1994 IF (master) WRITE (stdout,20) '!OCMIP_OXYGEN_SC', &
1995 & 'Schmidt number from Wanninkhof (1992)'
1996 is=len_trim(coptions)+1
1997 coptions(is:is+18)=' !OCMIP_OXYGEN_SC,'
1998# endif
1999!
2000 IF (master) WRITE (stdout,20) 'OXYGEN', &
2001 & 'Add oxygen dynamics'
2002 is=len_trim(coptions)+1
2003 coptions(is:is+8)=' OXYGEN,'
2004#endif
2005#if defined HYPOXIA_SRM && !defined SURFACE_DO_SATURATION
2006# ifdef OCMIP_OXYGEN_SC
2007!
2008 IF (master) WRITE (stdout,20) 'OCMIP_OXYGEN_SC', &
2009 & 'Schmidt number from Keeling et al. (1998)'
2010 is=len_trim(coptions)+1
2011 coptions(is:is+17)=' OCMIP_OXYGEN_SC,'
2012# else
2013!
2014 IF (master) WRITE (stdout,20) '!OCMIP_OXYGEN_SC', &
2015 & 'Schmidt number from Wanninkhof (1992)'
2016 is=len_trim(coptions)+1
2017 coptions(is:is+18)=' !OCMIP_OXYGEN_SC,'
2018# endif
2019#endif
2020#if defined PARALLEL_IO && defined DISTRIBUTE
2021!
2022 IF (master) WRITE (stdout,20) 'PARALLEL_IO', &
2023 & 'Parallel I/O processing'
2024 is=len_trim(coptions)+1
2025 coptions(is:is+13)=' PARALLEL_IO,'
2026#endif
2027#if defined CARBON && defined BIO_FENNEL
2028# if defined PCO2AIR_DATA
2029!
2030 IF (master) WRITE (stdout,20) 'PCO2AIR_DATA', &
2031 & 'CO2 gas exchange from Laurent et al. (2017)'
2032 is=len_trim(coptions)+1
2033 coptions(is:is+14)=' PCO2AIR_DATA,'
2034# elif defined PCO2AIR_SECULAR
2035!
2036 IF (master) WRITE (stdout,20) 'PCO2AIR_SECULAR', &
2037 & 'CO2 gas exchange time-dependent evolution'
2038 is=len_trim(coptions)+1
2039 coptions(is:is+17)=' PCO2AIR_SECULAR,'
2040# endif
2041#endif
2042#ifdef PERFECT_RESTART
2043!
2044 IF (master) WRITE (stdout,20) 'PERFECT_RESTART', &
2045 & 'Processing perfect restart variables'
2046 is=len_trim(coptions)+1
2047 coptions(is:is+17)=' PERFECT_RESTART,'
2048#endif
2049#ifdef PICARD_TEST
2050!
2051 IF (master) WRITE (stdout,20) 'PICARD_TEST', &
2052 & 'Representers tangent linear model Picard iterations test'
2053 is=len_trim(coptions)+1
2054 coptions(is:is+13)=' PICARD TEST,'
2055 idriver=idriver+1
2056#endif
2057#if defined PIO_LIB && defined DISTRIBUTE
2058!
2059 IF (master) WRITE (stdout,20) 'PIO_LIB', &
2060 & 'Using Parallel-IO from the PIO library'
2061 is=len_trim(coptions)+1
2062 coptions(is:is+9)=' PIO_LIB,'
2063#endif
2064#ifdef PJ_GRADP
2065!
2066 IF (master) WRITE (stdout,20) 'PJ_GRADP', &
2067 & 'Finite volume Pressure Jacobian (Lin, 1997)'
2068 is=len_trim(coptions)+1
2069 coptions(is:is+10)=' PJ_GRADP,'
2070#endif
2071#ifdef PJ_GRADPQ2
2072!
2073 IF (master) WRITE (stdout,20) 'PJ_GRADPQ2', &
2074 & 'Quartic 2 polynomial Pressure Jacobian (Shchepetkin, 2002)'
2075 is=len_trim(coptions)+1
2076 coptions(is:is+12)=' PJ_GRADPQ2,'
2077#endif
2078#ifdef PJ_GRADPQ4
2079!
2080 IF (master) WRITE (stdout,20) 'PJ_GRADPQ4', &
2081 & 'Quartic 4 polynomial Pressure Jacobian (Shchepetkin, 2002)'
2082 is=len_trim(coptions)+1
2083 coptions(is:is+12)=' PJ_GRADPQ4,'
2084#endif
2085#if defined PRESS_COMPENSATE && defined ATM_PRESS
2086!
2087 IF (master) WRITE (stdout,20) 'PRESS_COMPENSATE', &
2088 & 'Compensate for boundary data without ATM pressure effect'
2089 is=len_trim(coptions)+1
2090 coptions(is:is+18)=' PRESS_COMPENSATE,'
2091#endif
2092#if defined PO4 && defined BIO_FENNEL
2093!
2094 IF (master) WRITE (stdout,20) 'PO4', &
2095 & 'Phytoplanckton growth limited by Phosphorus'
2096 is=len_trim(coptions)+1
2097 coptions(is:is+5)=' P04,'
2098#endif
2099#if defined PRIOR_BULK_FLUXES && defined BULK_FLUXES && \
2100 defined four_dvar
2101!
2102 IF (master) WRITE (stdout,20) 'PRIOR_BULK_FLUXES', &
2103 & 'Impose prior surface fluxes in other outer loops and analysis'
2104 is=len_trim(coptions)+1
2105 coptions(is:is+19)=' PRIOR_BULK_FLUXES,'
2106#endif
2107#ifdef POSITIVE_ZERO
2108!
2109 IF (master) WRITE (stdout,20) 'POSITIVE_ZERO', &
2110 & 'Impose positive zero in output data'
2111 is=len_trim(coptions)+1
2112 coptions(is:is+15)=' POSITIVE_ZERO,'
2113#endif
2114#if defined POSTERIOR_EOFS && defined WEAK_CONSTRAINT
2115!
2116 IF (master) WRITE (stdout,20) 'POSTERIOR_EOFS', &
2117 & 'Estimating posterior analysis error covariace matrix EOFs'
2118 is=len_trim(coptions)+1
2119 coptions(is:is+16)=' POSTERIOR_EOFS,'
2120#endif
2121#if defined POSTERIOR_ERROR_F && defined WEAK_CONSTRAINT
2122!
2123 IF (master) WRITE (stdout,20) 'POSTERIOR_ERROR_F', &
2124 & 'Estimating final posterior analysis error covariace matrix'
2125 is=len_trim(coptions)+1
2126 coptions(is:is+16)=' POSTERIOR_ERROR_F,'
2127#endif
2128#if defined POSTERIOR_ERROR_I && defined WEAK_CONSTRAINT
2129!
2130 IF (master) WRITE (stdout,20) 'POSTERIOR_ERROR_I', &
2131 & 'Estimating initial posterior analysis error covariace matrix'
2132 is=len_trim(coptions)+1
2133 coptions(is:is+16)=' POSTERIOR_ERROR_F,'
2134#endif
2135#if defined POWER_LAW && defined SOLVE3D
2136!
2137 IF (master) WRITE (stdout,20) 'POWER_LAW', &
2138 & 'Power-law shape time-averaging barotropic filter'
2139 is=len_trim(coptions)+1
2140 coptions(is:is+11)=' POWER_LAW,'
2141#endif
2142#if !(defined PJ_GRADPQ4 || defined PJ_GRADPQ2 || defined PJ_GRADP || \
2143 defined dj_gradps) && defined solve3d
2144!
2145 IF (master) WRITE (stdout,20) 'PRSGRD31', &
2146 & 'Standard density Jacobian formulation (Song, 1998)'
2147 is=len_trim(coptions)+1
2148 coptions(is:is+10)=' PRSGRD31,'
2149#endif
2150#ifdef PROFILE
2151!
2152 IF (master) WRITE (stdout,20) 'PROFILE', &
2153 & 'Time profiling activated'
2154 is=len_trim(coptions)+1
2155 coptions(is:is+9)=' PROFILE,'
2156#endif
2157#ifdef PSEUDOSPECTRA
2158!
2159 IF (master) WRITE (stdout,20) 'PSEUDOSPECTRA', &
2160 & 'Pseudospectra Propagator'
2161 is=len_trim(coptions)+1
2162 coptions(is:is+15)=' PSEUDOSPECTRA,'
2163#endif
2164#ifdef QCORRECTION
2165!
2166 IF (master) WRITE (stdout,20) 'QCORRECTION', &
2167 & 'Surface net heat flux correction'
2168 is=len_trim(coptions)+1
2169 coptions(is:is+13)=' QCORRECTION,'
2170#endif
2171#ifdef OBSERVATIONS
2172# ifdef RADIAL_ANGLE_CCW_EAST
2173!
2174 IF (master) WRITE (stdout,20) 'RADIAL_ANGLE_CCW_EAST', &
2175 & 'Radial observations as an azimuth counterclockwise from EAST'
2176 is=len_trim(coptions)+1
2177 coptions(is:is+23)=' RADIAL_ANGLE_CCW_EAST,'
2178# else
2179 IF (master) WRITE (stdout,20) '!RADIAL_ANGLE_CCW_EAST', &
2180 & 'Radial observations as an azimuth clockwise from NORTH'
2181 is=len_trim(coptions)+1
2182 coptions(is:is+24)=' !RADIAL_ANGLE_CCW_EAST,'
2183# endif
2184#endif
2185#if defined NESTING && defined REFINE_BOUNDARY
2186!
2187 IF (master) WRITE (stdout,20) 'REFINE_BOUNDARY', &
2188 & 'Fine-to-Coarse momentum averaging at coarse grid boundary'
2189 is=len_trim(coptions)+1
2190 coptions(is:is+17)=' REFINE_BOUNDARY,'
2191#endif
2192#if defined REMOVE_LAPACK_GOTOS && defined FOUR_DVAR
2193!
2194 IF (master) WRITE (stdout,20) 'REMOVE_LAPACK_GOTOS', &
2195 & 'Modernizing adapted LAPACK routines by replacing GOTOs'
2196 is=len_trim(coptions)+1
2197 coptions(is:is+21)=' REMOVE_LAPACK_GOTOS,'
2198#endif
2199#if defined OXYGEN && defined BIO_FENNEL
2200# if defined OCMIP_OXYGEN_SC
2201!
2202 IF (master) WRITE (stdout,20) 'OCMIP_OXYGEN_SC', &
2203 & 'Schmidt number from Keeling et al. (1998)'
2204 is=len_trim(coptions)+1
2205 coptions(is:is+17)=' OCMIP_OXYGEN_SC,'
2206# elif defined RW14_OXYGEN_SC
2207!
2208 IF (master) WRITE (stdout,20) 'RW14_OXYGEN_SC', &
2209 & 'O2 Schmidt number from Wanninkhof (2014)'
2210 is=len_trim(coptions)+1
2211 coptions(is:is+16)=' RW14_OXYGEN_SC,'
2212# else
2213!
2214 IF (master) WRITE (stdout,20) '!RW14_OXYGEN_SC', &
2215 & 'Schmidt number from Wanninkhof (1992)'
2216 is=len_trim(coptions)+1
2217 coptions(is:is+17)=' !RW14_OXYGEN_SC,'
2218# endif
2219#endif
2220#if defined GLS_MIXING || defined MY25_MIXING
2221# if defined K_C2ADVECTION
2222!
2223 IF (master) WRITE (stdout,20) 'K_C2ADVECTION', &
2224 & 'Second-order centered differences advection of TKE fields'
2225 is=len_trim(coptions)+1
2226 coptions(is:is+15)=' K_C2ADVECTION,'
2227# elif defined K_C4ADVECTION
2228!
2229 IF (master) WRITE (stdout,20) 'K_C4ADVECTION', &
2230 & 'Fourth-order centered differences advection of TKE fields'
2231 is=len_trim(coptions)+1
2232 coptions(is:is+15)=' K_C4ADVECTION,'
2233# else
2234!
2235 IF (master) WRITE (stdout,20) 'K_GSCHEME', &
2236 & 'Third-order upstream advection of TKE fields'
2237 is=len_trim(coptions)+1
2238 coptions(is:is+11)=' K_GSCHEME,'
2239# endif
2240# ifdef TKE_DIF2
2241!
2242 IF (master) WRITE (stdout,20) 'TKE_DIF2', &
2243 & 'Harmonic mixing of TKE fields'
2244 is=len_trim(coptions)+1
2245 coptions(is:is+10)=' TKE_DIF2,'
2246# endif
2247# ifdef TKE_DIF4
2248!
2249 IF (master) WRITE (stdout,20) 'TKE_DIF4', &
2250 & 'Biharmonic mixing of TKE fields'
2251 is=len_trim(coptions)+1
2252 coptions(is:is+10)=' TKE_DIF4,'
2253# endif
2254#endif
2255#ifdef RADIATION_2D
2256!
2257 IF (master) WRITE (stdout,20) 'RADIATION_2D', &
2258 & 'Use tangential phase speed in radiation conditions'
2259 is=len_trim(coptions)+1
2260 coptions(is:is+14)=' RADIATION_2D,'
2261#endif
2262#if defined RAMP_TIDES && (defined SSH_TIDES || defined UV_TIDES)
2263!
2264 IF (master) WRITE (stdout,20) 'RAMP_TIDES', &
2265 & 'Ramping tidal forcing for one day'
2266 is=len_trim(coptions)+1
2267 coptions(is:is+14)=' RAMP_TIDES,'
2268#endif
2269#if defined RBL4DVAR && !defined SPLIT_RBL4DVAR
2270!
2271 IF (master) WRITE (stdout,20) 'RBL4DVAR', &
2272 & 'Strong/Weak constraint RBL4D-VAR data assimilation'
2273 is=len_trim(coptions)+1
2274 coptions(is:is+10)=' RBL4DVAR,'
2275 idriver=idriver+1
2276#endif
2277#ifdef RBL4DVAR_ANA_SENSITIVITY
2278!
2279 IF (master) WRITE (stdout,20) 'RBL4DVAR_ANA_SENSITIVITY', &
2280 & 'RBL4D-VAR analysis observation sensitivity'
2281 is=len_trim(coptions)+1
2282 coptions(is:is+26)=' RBL4DVAR_ANA_SENSITIVITY,'
2283 idriver=idriver+1
2284#endif
2285#ifdef RBL4DVAR_FCT_SENSITIVITY
2286!
2287 IF (master) WRITE (stdout,20) 'RBL4DVAR_FCT_SENSITIVITY', &
2288 & 'RBL4D-VAR forecast observation sensitivity'
2289 is=len_trim(coptions)+1
2290 coptions(is:is+26)=' RBL4DVAR_FCT_SENSITIVITY,'
2291 idriver=idriver+1
2292#endif
2293#if defined READ_WATER && defined MASKING
2294!
2295 IF (master) WRITE (stdout,20) 'READ_WATER', &
2296 & 'Reading data at water points only'
2297 is=len_trim(coptions)+1
2298 coptions(is:is+12)=' READ_WATER,'
2299#endif
2300#if defined RECOMPUTE_4DVAR && defined FOUR_DVAR
2301!
2302 IF (master) WRITE (stdout,20) 'RECOMPUTE_4DVAR', &
2303 & 'Recomputing 4D-Var during data assimilation analysis'
2304 is=len_trim(coptions)+1
2305 coptions(is:is+17)=' RECOMPUTE_4DVAR,'
2306#endif
2307#if defined DISTRIBUTE
2308# if defined REDUCE_ALLGATHER
2309!
2310 IF (master) WRITE (stdout,20) 'REDUCE_ALLGATHER', &
2311 & 'Using mpi_allgather in mp_reduce routine'
2312 is=len_trim(coptions)+1
2313 coptions(is:is+18)=' REDUCE_ALLGATHER,'
2314# elif defined REDUCE_ALLREDUCE
2315!
2316 IF (master) WRITE (stdout,20) 'REDUCE_ALLREDUCE', &
2317 & 'Using mpi_allreduce in mp_reduce routine'
2318 is=len_trim(coptions)+1
2319 coptions(is:is+18)=' REDUCE_ALLREDUCE,'
2320# elif defined REDUCE_SENDRECV
2321!
2322 IF (master) WRITE (stdout,20) 'REDUCE_SENDRECV', &
2323 & 'Using mpi_isend/mpi_recv in mp_reduce routine'
2324 is=len_trim(coptions)+1
2325 coptions(is:is+17)=' REDUCE_SENDRECV,'
2326# endif
2327#endif
2328#ifdef REGCM_COUPLING
2329!
2330 IF (master) WRITE (stdout,20) 'COAMPS_COUPLING', &
2331 & 'Atmosphere coupling with RegCM'
2332 is=len_trim(coptions)+1
2333 coptions(is:is+16)=' REGCM_COUPLING,'
2334 iatms=iatms+1
2335#endif
2336#ifdef RED_TIDE
2337!
2338 IF (master) WRITE (stdout,20) 'RED_TIDE', &
2339 & 'Red tide biological model'
2340 is=len_trim(coptions)+1
2341 coptions(is:is+10)=' RED_TIDE,'
2342 ibiology=ibiology+1
2343#endif
2344#if defined REGRESS_STARTCLOCK && defined MODEL_COUPLING && \
2345 defined esmf_lib
2346!
2347 IF (master) WRITE (stdout,20) 'REGRESS_STARTCLOCK', &
2348 & 'Regressing start clock by one coupling interval'
2349 is=len_trim(coptions)+1
2350 coptions(is:is+20)=' REGRESS_STARTCLOCK,'
2351#endif
2352#ifdef REGRID_SHAPIRO
2353!
2354 IF (master) WRITE (stdout,20) 'REGRID_SHAPIRO', &
2355 & 'Applying Shapiro filter to regridded data'
2356 is=len_trim(coptions)+1
2357 coptions(is:is+16)=' REGRID_SHAPIRO,'
2358#endif
2359#if defined RIVER_DON && defined BIO_FENNEL
2360!
2361 IF (master) WRITE (stdout,20) 'RIVER_DON', &
2362 & 'DON non-skinking source from rivers'
2363 is=len_trim(coptions)+1
2364 coptions(is:is+5)=' RIVER_DON,'
2365#endif
2366#ifdef ROMS_STDOUT
2367!
2368 IF (master) WRITE (stdout,20) 'ROMS_STDOUT', &
2369 & "Writes ROMS standard output into 'log.roms' file"
2370 is=len_trim(coptions)+1
2371 coptions(is:is+13)=' ROMS_STDOUT,'
2372#endif
2373#ifdef R_SYMMETRY
2374!
2375 IF (master) WRITE (stdout,20) 'REP_MATRIX', &
2376 & 'Representers Matrix Symmetry Test'
2377 is=len_trim(coptions)+1
2378 coptions(is:is+12)=' REP_MATRIX,'
2379 idriver=idriver+1
2380#endif
2381#if defined LMD_MIXING
2382# ifdef RI_HORAVG
2383!
2384 IF (master) WRITE (stdout,20) 'RI_HORAVG', &
2385 & 'Smooth Richardson number horizontally'
2386 is=len_trim(coptions)+1
2387 coptions(is:is+11)=' RI_HORAVG,'
2388# endif
2389# ifdef RI_VERAVG
2390!
2391 IF (master) WRITE (stdout,20) 'RI_VERAVG', &
2392 & 'Smooth Richardson number vertically'
2393 is=len_trim(coptions)+1
2394 coptions(is:is+11)=' RI_VERAVG,'
2395# endif
2396#endif
2397#if defined SOLVE3D && \
2398 (defined gls_mixing || defined lmd_mixing || defined my25_mixing)
2399# ifdef RI_SPLINES
2400!
2401 IF (master) WRITE (stdout,20) 'RI_SPLINES', &
2402 & 'Parabolic Spline Reconstruction for Richardson Number'
2403 is=len_trim(coptions)+1
2404 coptions(is:is+12)=' RI_SPLINES,'
2405# endif
2406#endif
2407#if !(defined PJ_GRADPQ4 || defined PJ_GRADPQ2 || defined PJ_GRADP || \
2408 defined dj_gradps) && defined rho_surf && defined solve3d
2409!
2410 IF (master) WRITE (stdout,20) 'RHO_SURF', &
2411 & 'Include difference between rho0 and surface density'
2412 is=len_trim(coptions)+1
2413 coptions(is:is+10)=' RHO_SURF,'
2414#endif
2415#if defined ROLLER_MONO && defined WEC_ROLLER && defined WEC_VF
2416!
2417 IF (master) WRITE (stdout,20) 'ROLLER_MONO', &
2418 & 'Wave energy roller from monochromatic waves'
2419 is=len_trim(coptions)+1
2420 coptions(is:is+13)=' ROLLER_MONO,'
2421#endif
2422#ifdef ROLLER_RENIERS && defined WEC_ROLLER && defined WEC_VF
2423!
2424 IF (master) WRITE (stdout,20) 'ROLLER_RENIERS', &
2425 & 'Wave energy roller based on Reniers 2004'
2426 is=len_trim(coptions)+1
2427 coptions(is:is+16)=' ROLLER_RENIERS,'
2428#endif
2429#ifdef ROLLER_SVENDSEN && defined WEC_ROLLER && defined WEC_VF
2430!
2431 IF (master) WRITE (stdout,20) 'ROLLER_SVENDSEN', &
2432 & 'Wave energy roller based on Svendsen 1984'
2433 is=len_trim(coptions)+1
2434 coptions(is:is+17)=' ROLLER_SVENDSEN,'
2435#endif
2436#if defined RPCG && defined WEAK_CONSTRAINT
2437!
2438 IF (master) WRITE (stdout,20) 'RPCG', &
2439 & 'Restricted B-preconditioned Lanczos minimization solver'
2440 is=len_trim(coptions)+1
2441 coptions(is:is+6)=' RPCG,'
2442#endif
2443#if defined RP_AVERAGES && defined TL_IOMS
2444!
2445 IF (master) WRITE (stdout,20) 'RP_AVERAGES', &
2446 & 'Writing out time-averaged representer model fields'
2447 is=len_trim(coptions)+1
2448 coptions(is:is+13)=' RP_AVERAGES,'
2449#endif
2450#ifdef RPM_DRIVER
2451!
2452 IF (master) WRITE (stdout,20) 'RPM_DRIVER', &
2453 & 'Generic representers tangent linear model driver'
2454 is=len_trim(coptions)+1
2455 coptions(is:is+12)=' RPM_DRIVER,'
2456 idriver=idriver+1
2457#endif
2458#ifdef RPM_RELAXATION
2459!
2460 IF (master) WRITE (stdout,20) 'RPM_RELAXATION', &
2461 & 'Diffusive Relaxation of RPM in Picard Iterations'
2462 is=len_trim(coptions)+1
2463 coptions(is:is+16)=' RPM_RELAXATION,'
2464#endif
2465#ifdef RST_SINGLE
2466!
2467 IF (master) WRITE (stdout,20) 'RST_SINGLE', &
2468 & 'Single precision fields in restart NetCDF file'
2469 is=len_trim(coptions)+1
2470 coptions(is:is+12)=' RST_SINGLE,'
2471#else
2472!
2473 IF (master) WRITE (stdout,20) '!RST_SINGLE', &
2474 & 'Double precision fields in restart NetCDF file'
2475 is=len_trim(coptions)+1
2476 coptions(is:is+13)=' !RST_SINGLE,'
2477#endif
2478#if defined R4DVAR && !defined SPLIT_R4DVAR
2479!
2480 IF (master) WRITE (stdout,20) 'R4DVAR', &
2481 & 'Strong/Weak constraint, R4D-Var data assimilation'
2482 is=len_trim(coptions)+1
2483 coptions(is:is+8)=' R4DVAR,'
2484 idriver=idriver+1
2485#endif
2486#ifdef R4DVAR_ANA_SENSITIVITY
2487!
2488 IF (master) WRITE (stdout,20) 'R4DVAR_ANA_SENSITIVITY', &
2489 & 'R4D-Var analysis observation sensitivity'
2490 is=len_trim(coptions)+1
2491 coptions(is:is+24)=' R4DVAR_ANA_SENSITIVITY,'
2492 idriver=idriver+1
2493#endif
2494#ifdef SALINITY
2495!
2496 IF (master) WRITE (stdout,20) 'SALINITY', &
2497 & 'Using salinity'
2498 is=len_trim(coptions)+1
2499 coptions(is:is+10)=' SALINITY,'
2500#endif
2501#ifdef SANITY_CHECK
2502!
2503 IF (master) WRITE (stdout,20) 'SANITY_CHECK', &
2504 & 'Tangent Linear and Adjoint Models Correctness Check'
2505 is=len_trim(coptions)+1
2506 coptions(is:is+14)=' SANITY_CHECK,'
2507 idriver=idriver+1
2508#endif
2509#ifdef DISTRIBUTE
2510# ifdef SCATTER_BCAST
2511 IF (master) WRITE (stdout,20) 'SCATTER_BCAST', &
2512 & 'Using mpi_bcast in mp_scatter2d/mp_scatter3d routines'
2513 is=len_trim(coptions)+1
2514 coptions(is:is+15)=' SCATTER_BCAST,'
2515# else
2516 IF (master) WRITE (stdout,20) '!SCATTER_BCAST', &
2517 & 'Using mpi_scatterv in mp_scatter2d/mp_scatter3d routines'
2518 is=len_trim(coptions)+1
2519 coptions(is:is+16)=' !SCATTER_BCAST,'
2520# endif
2521#endif
2522#ifdef SCORRECTION
2523!
2524 IF (master) WRITE (stdout,20) 'SCORRECTION', &
2525 & 'Surface salinity flux correction'
2526 is=len_trim(coptions)+1
2527 coptions(is:is+13)=' SCORRECTION,'
2528#endif
2529#ifdef SEDIMENT
2530!
2531 IF (master) WRITE (stdout,20) 'SEDIMENT', &
2532 & 'Cohesive and noncohesive sediments'
2533 is=len_trim(coptions)+1
2534 coptions(is:is+10)=' SEDIMENT,'
2535# ifdef SED_DENS
2536!
2537 IF (master) WRITE (stdout,20) 'SED_DENS', &
2538 & 'Include density of suspended sediment in equation of state'
2539 is=len_trim(coptions)+1
2540 coptions(is:is+10)=' SED_DENS,'
2541# endif
2542# ifdef SED_MORPH
2543!
2544 IF (master) WRITE (stdout,20) 'SED_MORPH', &
2545 & 'Allow bottom model elevation to evolve'
2546 is=len_trim(coptions)+1
2547 coptions(is:is+11)=' SED_MORPH,'
2548# endif
2549# ifdef SUSPLOAD
2550!
2551 IF (master) WRITE (stdout,20) 'SUSPLOAD', &
2552 & 'Activate suspended sediment transport'
2553 is=len_trim(coptions)+1
2554 coptions(is:is+10)=' SUSPLOAD,'
2555# endif
2556#endif
2557#if defined STD_MODEL && defined FOUR_DVAR
2558!
2559 IF (master) WRITE (stdout,20) 'STD_MODEL', &
2560 & 'Background error covariance Standard Deviation model'
2561 is=len_trim(coptions)+1
2562 coptions(is:is+11)=' STD_MODEL,'
2563#endif
2564#if defined STEP2D_FB_AB3_AM4
2565!
2566 IF (master) WRITE (stdout,20) 'STEP2D_FB_AB3_AM4', &
2567 & 'Generalized Forward-Backward AB3-AM4 stepping scheme'
2568 is=len_trim(coptions)+1
2569 coptions(is:is+19)=' STEP2D_FB_AB3_AM4,'
2570#elif defined STEP2D_FB_LF_AM3
2571!
2572 IF (master) WRITE (stdout,20) 'STEP2D_FB_LF_AM3', &
2573 & 'Predictor/Corrector Forward-Backward LF-AM3 stepping scheme'
2574 is=len_trim(coptions)+1
2575 coptions(is:is+18)=' STEP2D_FB_LF_AM3,'
2576#else
2577!
2578 IF (master) WRITE (stdout,20) 'STEP2D_LF_AM3', &
2579 & 'Predictor/Corrector LF-AM3 stepping scheme'
2580 is=len_trim(coptions)+1
2581 coptions(is:is+15)=' STEP2D_LF_AM3,'
2582#endif
2583#ifdef SG_BBL
2584!
2585 IF (master) WRITE (stdout,20) 'SG_BBL', &
2586 & 'Styles and Glenn Bottom Boundary Layer'
2587 is=len_trim(coptions)+1
2588 coptions(is:is+8)=' SG_BBL,'
2589 ibbl=ibbl+1
2590# ifdef SG_CALC_UB
2591!
2592 IF (master) WRITE (stdout,20) 'SG_CALC_UB', &
2593 & 'Internal computation of bottom orbital velocity'
2594 is=len_trim(coptions)+1
2595 coptions(is:is+12)=' SG_CALC_UB,'
2596# endif
2597# ifdef SG_CALC_ZNOT
2598!
2599 IF (master) WRITE (stdout,20) 'SG_CALC_ZNOT', &
2600 & 'Internal computation of bottom roughness'
2601 is=len_trim(coptions)+1
2602 coptions(is:is+14)=' SG_CALC_ZNOT,'
2603# endif
2604# ifdef SG_LOGINT
2605!
2606 IF (master) WRITE (stdout,20) 'SG_LOGINT', &
2607 & 'Bottom currents logarithmic interpolation'
2608 is=len_trim(coptions)+1
2609 coptions(is:is+11)=' SG_LOGINT,'
2610# endif
2611#endif
2612#ifdef SINGLE_PRECISION
2613!
2614 IF (master) WRITE (stdout,20) 'SINGLE_PRECISION', &
2615 & 'Single precision arithmetic numerical kernel'
2616 is=len_trim(coptions)+1
2617 coptions(is:is+18)=' SINGLE_PRECISION,'
2618#endif
2619#if defined SKIP_NLM && defined SENSITIVITY_4DVAR
2620!
2621 IF (master) WRITE (stdout,20) 'SKIP_NLM', &
2622 & 'Skipping running NLM, reading NLM state trajectory'
2623 is=len_trim(coptions)+1
2624 coptions(is:is+10)=' SKIP_NLM,'
2625#endif
2626#ifdef SPECTRUM_STOKES
2627!
2628 IF (master) WRITE (stdout,20) 'SPECTRUM_STOKES', &
2629 & 'Stokes drift estimated from wave spectrum'
2630 is=len_trim(coptions)+1
2631 coptions(is:is+17)=' SPECTRUM_STOKES,'
2632#endif
2633#ifdef SRELAXATION
2634!
2635 IF (master) WRITE (stdout,20) 'SRELAXATION', &
2636 & 'Surface salinity relaxation as surface flux'
2637 is=len_trim(coptions)+1
2638 coptions(is:is+13)=' SRELAXATION,'
2639#endif
2640#ifdef SOLAR_SOURCE
2641!
2642 IF (master) WRITE (stdout,20) 'SOLAR_SOURCE', &
2643 & 'Solar Radiation Source Term'
2644 is=len_trim(coptions)+1
2645 coptions(is:is+14)=' SOLAR_SOURCE,'
2646#endif
2647#ifdef SOLVE3D
2648!
2649 IF (master) WRITE (stdout,20) 'SOLVE3D', &
2650 & 'Solving 3D Primitive Equations'
2651 is=len_trim(coptions)+1
2652 coptions(is:is+9)=' SOLVE3D,'
2653#endif
2654#ifdef SO_SEMI
2655!
2656 IF (master) WRITE (stdout,20) 'SO_SEMI', &
2657 & 'Stochastic Optimals Propagator, Semi-norm Estimation'
2658 is=len_trim(coptions)+1
2659 coptions(is:is+9)=' SO_SEMI,'
2660 idriver=idriver+1
2661#endif
2662#ifdef SO_TRACE
2663!
2664 IF (master) WRITE (stdout,20) 'SO_TRACE', &
2665 & 'Stochastic Optimals Propagator, Randomized Trace Estimation'
2666 is=len_trim(coptions)+1
2667 coptions(is:is+10)=' SO_TRACE,'
2668 idriver=idriver+1
2669#endif
2670#ifdef SO_SEMI
2671# ifdef SO_SEMI_WHITE
2672!
2673 IF (master) WRITE (stdout,20) 'SO_SEMI_WHITE', &
2674 & 'Stochastic Optimals, white noise processes'
2675 is=len_trim(coptions)+1
2676 coptions(is:is+15)=' SO_SEMI_WHITE,'
2677# else
2678!
2679 IF (master) WRITE (stdout,20) '!SO_WHITE', &
2680 & 'Stochastic Optimals, red noise processes'
2681 is=len_trim(coptions)+1
2682 coptions(is:is+16)=' !SO_SEMI_WHITE,'
2683# endif
2684#endif
2685#ifdef SPLIT_I4DVAR
2686!
2687 IF (master) WRITE (stdout,20) 'SPLIT_I4DVAR', &
2688 & 'Incremental, strong constraint Split 4D-Var data assimilation'
2689 is=len_trim(coptions)+1
2690 coptions(is:is+14)=' SPLIT_I4DVAR,'
2691 idriver=idriver+1
2692#endif
2693#ifdef SPLIT_RBL4DVAR
2694!
2695 IF (master) WRITE (stdout,20) 'SPLIT_RBL4DVAR', &
2696 & 'Strong/Weak constraint Split RBL4D-VAR data assimilation'
2697 is=len_trim(coptions)+1
2698 coptions(is:is+16)=' SPLIT_RBL4DVAR,'
2699 idriver=idriver+1
2700#endif
2701#ifdef SPLIT_R4DVAR
2702!
2703 IF (master) WRITE (stdout,20) 'SPLIT_R4DVAR', &
2704 & 'Strong/Weak constraint, Split R4D-Var data assimilation'
2705 is=len_trim(coptions)+1
2706 coptions(is:is+14)=' SPLIT_R4DVAR,'
2707 idriver=idriver+1
2708#endif
2709#ifdef SPLIT_SP4DVAR
2710!
2711 IF (master) WRITE (stdout,20) 'SPLIT_SP4DVAR', &
2712 & 'Weak Constraint, Saddle-Point Split 4D-Var Data Assimilation'
2713 is=len_trim(coptions)+1
2714 coptions(is:is+15)=' SPLIT_SP4DVAR,'
2715 idriver=idriver+1
2716#endif
2717#if defined SP4DVAR && !defined SPLIT_4DVAR
2718!
2719 IF (master) WRITE (stdout,20) 'SP4DVAR', &
2720 & 'Weak Constraint, Saddle-Point 4D-Var Data Assimilation'
2721 is=len_trim(coptions)+1
2722 coptions(is:is+9)=' SP4DVAR,'
2723 idriver=idriver+1
2724#endif
2725#ifdef SPHERICAL
2726!
2727 IF (master) WRITE (stdout,20) 'SPHERICAL', &
2728 & 'Spherical grid configuration'
2729 is=len_trim(coptions)+1
2730 coptions(is:is+11)=' SPHERICAL,'
2731#endif
2732#if defined SPLINES_VCONV && defined FOUR_DVAR && defined VCONVOLUTION
2733!
2734 IF (master) WRITE (stdout,20) 'SPLINES_VCONV', &
2735 & 'Implicit, Parabolic Splines Vertical Convolution Algorithm'
2736 is=len_trim(coptions)+1
2737 coptions(is:is+16)=' SPLINES_VCONV,'
2738#endif
2739#if defined SPLINES_VDIFF && defined SOLVE3D
2740!
2741 IF (master) WRITE (stdout,20) 'SPLINES_VDIFF', &
2742 & 'Parabolic Spline Reconstruction for Vertical Diffusion'
2743 is=len_trim(coptions)+1
2744 coptions(is:is+15)=' SPLINES_VDIFF,'
2745#endif
2746#if defined SPLINES_VVISC && defined SOLVE3D
2747!
2748 IF (master) WRITE (stdout,20) 'SPLINES_VVISC', &
2749 & 'Parabolic Spline Reconstruction for Vertical Viscosity'
2750 is=len_trim(coptions)+1
2751 coptions(is:is+15)=' SPLINES_VVISC,'
2752#endif
2753#if defined SSH_TIDES
2754!
2755 IF (master) WRITE (stdout,20) 'SSH_TIDES', &
2756 & 'Adding tidal elevation to SSH climatology'
2757 is=len_trim(coptions)+1
2758 coptions(is:is+11)=' SSH_TIDES,'
2759#endif
2760#if defined STEP2D_CORILIS && defined SOLVE3D && \
2761 (defined step2d_fb_ab3_am4 || defined step2d_fb_lf_am3)
2762!
2763 IF (master) WRITE (stdout,20) 'STEP2D_CORILIS', &
2764 & 'Adding Coriolis force to barotropic kernel in 3D set-up'
2765 is=len_trim(coptions)+1
2766 coptions(is:is+16)=' STEP2D_CORILIS,'
2767#endif
2768#ifdef SSW_BBL
2769!
2770 IF (master) WRITE (stdout,20) 'SSW_BBL', &
2771 & 'Styles and Glenn Bottom Boundary Layer - modified'
2772 is=len_trim(coptions)+1
2773 coptions(is:is+8)=' SSW_BBL,'
2774 ibbl=ibbl+1
2775# ifdef SSW_CALC_UB
2776!
2777 IF (master) WRITE (stdout,20) 'SSW_CALC_UB', &
2778 & 'Internal computation of bottom orbital velocity'
2779 is=len_trim(coptions)+1
2780 coptions(is:is+14)=' SSW_CALC_UB,'
2781# endif
2782# ifdef SSW_CALC_ZNOT
2783!
2784 IF (master) WRITE (stdout,20) 'SSW_CALC_ZNOT', &
2785 & 'Internal computation of bottom roughness'
2786 is=len_trim(coptions)+1
2787 coptions(is:is+14)=' SSW_CALC_ZNOT,'
2788# endif
2789# ifdef SSW_FORM_DRAG_COR
2790!
2791 IF (master) WRITE (stdout,20) 'SSW_FORM_DRAG_COR', &
2792 & 'Form drag coefficient'
2793 is=len_trim(coptions)+1
2794 coptions(is:is+19)=' SSW_FORM_DRAG_COR,'
2795# endif
2796# ifdef SSW_ZOBIO
2797!
2798 IF (master) WRITE (stdout,20) 'SSW_Z0BIO', &
2799 & 'Biogenic bedfrom roughness from ripples'
2800 is=len_trim(coptions)+1
2801 coptions(is:is+11)=' SSW_Z0BIO,'
2802# endif
2803# ifdef SSW_ZOBL
2804!
2805 IF (master) WRITE (stdout,20) 'SSW_Z0BL', &
2806 & 'Bedload roughness from ripples'
2807 is=len_trim(coptions)+1
2808 coptions(is:is+10)=' SSW_Z0BL,'
2809# endif
2810# ifdef SSW_ZORIP
2811!
2812 IF (master) WRITE (stdout,20) 'SSW_Z0RIP', &
2813 & 'Bedfrom roughness from ripples'
2814 is=len_trim(coptions)+1
2815 coptions(is:is+11)=' SSW_Z0RIP,'
2816# endif
2817# ifdef SSW_LOGINT
2818!
2819 IF (master) WRITE (stdout,20) 'SSW_LOGINT', &
2820 & 'Bottom currents logarithmic interpolation'
2821 is=len_trim(coptions)+1
2822 coptions(is:is+12)=' SSW_LOGINT,'
2823# endif
2824# ifdef SSW_LOGINT_WBL
2825!
2826 IF (master) WRITE (stdout,20) 'SSW_LOGINT_WBL', &
2827 & 'Bottom currents logarithmic interpolation at wbl'
2828 is=len_trim(coptions)+1
2829 coptions(is:is+16)=' SSW_LOGINT_WBL,'
2830# endif
2831# ifdef SSW_LOGINT_DIRECT
2832!
2833 IF (master) WRITE (stdout,20) 'SSW_LOGINT_DIRECT', &
2834 & 'Bottom currents logarithmic interpolation at user elev'
2835 is=len_trim(coptions)+1
2836 coptions(is:is+19)=' SSW_LOGINT_DIRECT,'
2837# endif
2838#endif
2839#ifdef STATIONS
2840!
2841 IF (master) WRITE (stdout,20) 'STATIONS', &
2842 & 'Writing out station data'
2843 is=len_trim(coptions)+1
2844 coptions(is:is+10)=' STATIONS,'
2845#endif
2846#if defined STATIONS_CGRID && defined STATIONS
2847!
2848 IF (master) WRITE (stdout,20) 'STATIONS_CGRID', &
2849 & 'Extracting station data at native C-grid locations'
2850 is=len_trim(coptions)+1
2851 coptions(is:is+16)=' STATIONS_CGRID,'
2852#endif
2853#ifdef STOCHASTIC_OPT
2854!
2855 IF (master) WRITE (stdout,20) 'STOCHASTIC_OPT', &
2856 & 'Stochastic Optimals Propagator'
2857 is=len_trim(coptions)+1
2858 coptions(is:is+16)=' STOCHASTIC_OPT,'
2859 idriver=idriver+1
2860# ifdef STOCH_OPT_WHITE
2861!
2862 IF (master) WRITE (stdout,20) 'STOCH_OPT_WHITE', &
2863 & 'Stochastic Optimals, white noise processes'
2864 is=len_trim(coptions)+1
2865 coptions(is:is+17)=' STOCH_OPT_WHITE,'
2866# else
2867!
2868 IF (master) WRITE (stdout,20) '!STOCH_OPT_WHITE', &
2869 & 'Stochastic Optimals, red noise processes'
2870 is=len_trim(coptions)+1
2871 coptions(is:is+18)=' !STOCH_OPT_WHITE,'
2872# endif
2873#endif
2874#if defined SURFACE_DO_SATURATION && defined HYPOXIA_SRM
2875!
2876 IF (master) WRITE (stdout,20) 'SURFACE_DO_SATURATION', &
2877 & 'Use surface dissolved oxygen saturation'
2878 is=len_trim(coptions)+1
2879 coptions(is:is+23)=' SURFACE_DO_SATURATION,'
2880#endif
2881#if defined SUBOBJECT_DEALLOCATION
2882!
2883 IF (master) WRITE (stdout,20) 'SUBOBJECT_DEALLOCATION', &
2884 & 'Deallocate each variable in structure objects separately'
2885 is=len_trim(coptions)+1
2886 coptions(is:is+24)=' SUBOBJECT_DEALLOCATION,'
2887#endif
2888#if defined SUPPRESS_REPORT && defined SPLIT_4DVAR
2889!
2890 IF (master) WRITE (stdout,20) 'SUPPRESS_REPORT', &
2891 & 'Suppress report of repetitive configuration information'
2892 is=len_trim(coptions)+1
2893 coptions(is:is+17)=' SUPPRESS_REPORT,'
2894#endif
2895#if defined SURFACE_STREAMING && defined WEC_VF
2896!
2897 IF (master) WRITE (stdout,20) 'SURFACE_STREAMING', &
2898 & 'Wave current streaming term'
2899 is=len_trim(coptions)+1
2900 coptions(is:is+18)=' SURFACE_STREAMING,'
2901#endif
2902#ifdef SWAN_COUPLING
2903!
2904 IF (master) WRITE (stdout,20) 'SWAN_COUPLING', &
2905 & 'Two-way SWAN/ROMS coupling'
2906 is=len_trim(coptions)+1
2907 coptions(is:is+15)=' SWAN_COUPLING,'
2908#endif
2909#if defined CARBON && defined BIO_FENNEL
2910# ifdef TALK_NONCONSERV
2911!
2912 IF (master) WRITE (stdout,20) 'TALK_NONCONSERV', &
2913 & 'Alkalinity is affected by changes in nitrate or ammonium'
2914 is=len_trim(coptions)+1
2915 coptions(is:is+17)=' TALK_NONCONSERV,'
2916# else
2917!
2918 IF (master) WRITE (stdout,20) '!TALK_NONCONSERV', &
2919 & 'Alkalinity is passive and unaffected by nitrate or ammonium'
2920 is=len_trim(coptions)+1
2921 coptions(is:is+18)=' !TALK_NONCONSERV,'
2922# endif
2923#endif
2924#ifdef TANGENT
2925!
2926 IF (master) WRITE (stdout,20) 'TANGENT', &
2927 & 'Tangent Linear Model'
2928 is=len_trim(coptions)+1
2929 coptions(is:is+9)=' TANGENT,'
2930#endif
2931#ifdef TIDE_GENERATING_FORCES
2932!
2933 IF (master) WRITE (stdout,20) 'TIDE_GENERATING_FORCES', &
2934 & 'Adding astronomical TGF term to pressure gradient'
2935 is=len_trim(coptions)+1
2936 coptions(is:is+24)=' TIDE_GENERATING_FORCES,'
2937#endif
2938#if defined TIME_CONV && defined WEAK_CONSTRAINT
2939!
2940 IF (master) WRITE (stdout,20) 'TIME_CONV', &
2941 & 'Weak-constraint 4D-Var time convolution'
2942 is=len_trim(coptions)+1
2943 coptions(is:is+11)=' TIME_CONV,'
2944#endif
2945#if defined TIME_INTERP && defined MODEL_COUPLING
2946!
2947 IF (master) WRITE (stdout,20) 'TIME_INTERP', &
2948 & 'Time interpolation from imported coupling snapshots data'
2949 is=len_trim(coptions)+1
2950 coptions(is:is+13)=' TIME_INTERP,'
2951#endif
2952#if defined TIME_INTERP_FLUX && defined NESTING && defined SOLVE3D
2953!
2954 IF (master) WRITE (stdout,20) 'TIME_INTERP_FLUX', &
2955 & 'Time interpolate coaser grid mass flux in refinement'
2956 is=len_trim(coptions)+1
2957 coptions(is:is+18)=' TIME_INTERP_FLUX,'
2958#endif
2959#if defined TL_AVERAGES && defined TANGENT
2960!
2961 IF (master) WRITE (stdout,20) 'TL_AVERAGES', &
2962 & 'Writing out time-averaged tangent linear model fields'
2963 is=len_trim(coptions)+1
2964 coptions(is:is+13)=' TL_AVERAGES,'
2965#endif
2966#ifdef TLM_DRIVER
2967!
2968 IF (master) WRITE (stdout,20) 'TLM_DRIVER', &
2969 & 'Generic tangent linear model driver'
2970 is=len_trim(coptions)+1
2971 coptions(is:is+12)=' TLM_DRIVER,'
2972 idriver=idriver+1
2973#endif
2974#if defined GLS_MIXING && defined TKE_WAVEDISS
2975!
2976 IF (master) WRITE (stdout,20) 'TKE_WAVEDISS', &
2977 & 'Wave breaking surface tke flux based on wave amplitude'
2978 is=len_trim(coptions)+1
2979 coptions(is:is+14)=' TKE_WAVEDISS,'
2980#endif
2981#ifdef TL_IOMS
2982!
2983 IF (master) WRITE (stdout,20) 'TL_IOMS', &
2984 & 'Representers Tangent Linear Model'
2985 is=len_trim(coptions)+1
2986 coptions(is:is+9)=' TL_IOMS,'
2987#endif
2988#ifdef TLM_CHECK
2989!
2990 IF (master) WRITE (stdout,20) 'TLM_CHECK', &
2991 & 'Tangent Linear Model Linearization Check'
2992 is=len_trim(coptions)+1
2993 coptions(is:is+11)=' TLM_CHECK,'
2994#endif
2995#if defined TRACING && defined PROFILE
2996!
2997 IF (master) WRITE (stdout,20) 'TRACING', &
2998 & 'Tracing sequence of routine calls for debugging'
2999 is=len_trim(coptions)+1
3000 coptions(is:is+9)=' TRACING,'
3001#endif
3002#ifdef SOLVE3D
3003# if defined T_PASSIVE
3004!
3005 IF (master) WRITE (stdout,20) 'T_PASSIVE', &
3006 & 'Advecting and diffusing inert passive tracer'
3007 is=len_trim(coptions)+1
3008 coptions(is:is+11)=' T_PASSIVE,'
3009# endif
3010# if defined TS_MIX_CLIMA && \
3011 (defined ts_dif2 || defined ts_dif4)
3012!
3013 IF (master) WRITE (stdout,20) 'TS_MIX_CLIMA', &
3014 & 'Horizontal diffusion of tracer perturbation (T-Tclm)'
3015 is=len_trim(coptions)+1
3016 coptions(is:is+14)=' TS_MIX_CLIMA,'
3017# endif
3018# if defined TS_MIX_MAX_SLOPE && defined MIX_ISO_TS && \
3019 (defined ts_dif2 || defined ts_dif4)
3020!
3021 IF (master) WRITE (stdout,20) 'TS_MIX_MAX_SLOPE', &
3022 & 'Maximum density slope in isopycnal diffusion'
3023 is=len_trim(coptions)+1
3024 coptions(is:is+18)=' TS_MIX_MAX_SLOPE,'
3025# endif
3026# if defined TS_MIX_MAX_SLOPE && defined MIX_ISO_TS && \
3027 (defined ts_dif2 || defined ts_dif4)
3028!
3029 IF (master) WRITE (stdout,20) 'TS_MIX_MAX_STRAT', &
3030 & 'Minimum density stratification in isopycnal diffusion'
3031 is=len_trim(coptions)+1
3032 coptions(is:is+18)=' TS_MIX_MAX_STRAT,'
3033# endif
3034# if defined TS_MIX_STABILITY && \
3035 (defined ts_dif2 || defined ts_dif4)
3036!
3037 IF (master) WRITE (stdout,20) 'TS_MIX_STABILITY', &
3038 & 'Weighting diffusion between two time levels'
3039 is=len_trim(coptions)+1
3040 coptions(is:is+18)=' TS_MIX_STABILITY,'
3041# endif
3042# ifdef TS_MPDATA_LIMIT
3043!
3044 IF (master) WRITE (stdout,20) 'TS_MPDATA_LIMIT', &
3045 & 'Further limiter in upwind corrector fluxes for stability'
3046 is=len_trim(coptions)+1
3047 coptions(is:is+17)=' TS_MPDATA_LIMIT,'
3048# endif
3049#endif
3050#if defined TS_DIF2 && defined SOLVE3D
3051!
3052 IF (master) WRITE (stdout,20) 'TS_DIF2', &
3053 & 'Harmonic mixing of tracers'
3054 is=len_trim(coptions)+1
3055 coptions(is:is+9)=' TS_DIF2,'
3056#endif
3057#if defined TS_DIF4 && defined SOLVE3D
3058!
3059 IF (master) WRITE (stdout,20) 'TS_DIF4', &
3060 & 'Biharmonic mixing of tracers'
3061 is=len_trim(coptions)+1
3062 coptions(is:is+9)=' TS_DIF4,'
3063#endif
3064#ifdef TS_FIXED
3065!
3066 IF (master) WRITE (stdout,20) 'TS_FIXED', &
3067 & 'Diagnostic configuration, no evolution of tracer'
3068 is=len_trim(coptions)+1
3069 coptions(is:is+10)=' TS_FIXED,'
3070#endif
3071#if defined TS_SMAGORINSKY && (defined TS_DIF2 || defined TS_DIF4)
3072!
3073 IF (master) WRITE (stdout,20) 'TS_SMAGORINSKY', &
3074 & 'Smagorinksy-like time-dependent diffusion coefficients'
3075 is=len_trim(coptions)+1
3076 coptions(is:is+16)=' TS_SMAGORINSKY,'
3077#endif
3078#ifdef UV_ADV
3079!
3080 IF (master) WRITE (stdout,20) 'UV_ADV', &
3081 & 'Advection of momentum'
3082 is=len_trim(coptions)+1
3083 coptions(is:is+8)=' UV_ADV,'
3084#endif
3085#ifdef UV_COR
3086!
3087 IF (master) WRITE (stdout,20) 'UV_COR', &
3088 & 'Coriolis term'
3089 is=len_trim(coptions)+1
3090 coptions(is:is+8)=' UV_COR,'
3091#endif
3092#ifdef UV_DESTAGGERED
3093!
3094 IF (master) WRITE (stdout,20) 'UV_DESTAGGERED', &
3095 & 'Currents destaggering to cell-center for data assimilation'
3096 is=len_trim(coptions)+1
3097 coptions(is:is+16)=' UV_DESTAGGERED,'
3098#endif
3099#ifdef UV_ADV
3100# ifdef UV_C2ADVECTION
3101!
3102 IF (master) WRITE (stdout,20) 'UV_C2ADVECTION', &
3103 & 'Second-order centered differences advection of momentum'
3104 is=len_trim(coptions)+1
3105 coptions(is:is+16)=' UV_C2ADVECTION,'
3106 ivelhadv=ivelhadv+1
3107 ivelvadv=ivelvadv+1
3108# endif
3109# if defined UV_C4ADVECTION || defined UV_U3ADV_SPLIT
3110!
3111 IF (master) WRITE (stdout,20) 'UV_C4ADVECTION', &
3112 & 'Fourth-order centered differences advection of momentum'
3113 is=len_trim(coptions)+1
3114 coptions(is:is+16)=' UV_C4ADVECTION,'
3115 ivelhadv=ivelhadv+1
3116 ivelvadv=ivelvadv+1
3117# endif
3118# ifdef SOLVE3D
3119# if !(defined UV_C2ADVECTION || defined UV_C4ADVECTION)
3120!
3121 IF (master) WRITE (stdout,20) 'UV_U3HADVECTION', &
3122 & 'Third-order upstream horizontal advection of 3D momentum'
3123 is=len_trim(coptions)+1
3124 coptions(is:is+17)=' UV_U3HADVECTION,'
3125 ivelhadv=ivelhadv+1
3126# endif
3127# if !(defined UV_C2ADVECTION || defined UV_C4ADVECTION || \
3128 defined uv_sadvection)
3129!
3130 IF (master) WRITE (stdout,20) 'UV_C4VADVECTION', &
3131 & 'Fourth-order centered vertical advection of momentum'
3132 is=len_trim(coptions)+1
3133 coptions(is:is+16)=' UV_C4VADVECTION,'
3134 ivelvadv=ivelvadv+1
3135# endif
3136# ifdef UV_KIRBY
3137 IF (master) WRITE (stdout,20) 'UV_KIRBY', &
3138 & 'Compute uwave and vwave Kirby avg velocities'
3139 is=len_trim(coptions)+1
3140 coptions(is:is+15)=' UV_KIRBY,'
3141# endif
3142# ifdef UV_SADVECTION
3143!
3144 IF (master) WRITE (stdout,20) 'UV_SADVECTION', &
3145 & 'Parabolic splines vertical advection of momentum'
3146 is=len_trim(coptions)+1
3147 coptions(is:is+15)=' UV_SADVECTION,'
3148 ivelvadv=ivelvadv+1
3149# endif
3150# ifdef UV_U3ADV_SPLIT
3151!
3152 IF (master) WRITE (stdout,20) 'UV_U3ADV_SPLIT', &
3153 & 'Split third-order upstream advection of momentum'
3154 is=len_trim(coptions)+1
3155 coptions(is:is+16)=' UV_U3ADV_SPLIT,'
3156# endif
3157# endif
3158#endif
3159#ifdef UV_DRAG_GRID
3160!
3161 IF (master) WRITE (stdout,20) 'UV_DRAG_GRID', &
3162# if defined UV_LOGDRAG || defined BBL_MODEL
3163 & 'Spatially varying bottom roughness length'
3164# elif defined UV_LDRAG
3165 & 'Spatially varying linear drag coefficient'
3166# elif defined UV_QDRAG
3167 & 'Spatially varying quadratic drag coefficient'
3168# endif
3169 is=len_trim(coptions)+1
3170 coptions(is:is+14)=' UV_DRAG_GRID,'
3171#endif
3172#ifdef UV_LDRAG
3173!
3174 IF (master) WRITE (stdout,20) 'UV_LDRAG', &
3175 & 'Linear bottom stress'
3176 is=len_trim(coptions)+1
3177 coptions(is:is+10)=' UV_LDRAG,'
3178 ibbl=ibbl+1
3179#endif
3180#ifdef UV_LOGDRAG
3181!
3182 IF (master) WRITE (stdout,20) 'UV_LOGDRAG', &
3183 & 'Logarithmic bottom stress'
3184 is=len_trim(coptions)+1
3185 coptions(is:is+12)=' UV_LOGDRAG,'
3186 ibbl=ibbl+1
3187#endif
3188#ifdef UV_QDRAG
3189!
3190 IF (master) WRITE (stdout,20) 'UV_QDRAG', &
3191 & 'Quadratic bottom stress'
3192 is=len_trim(coptions)+1
3193 coptions(is:is+10)=' UV_QDRAG,'
3194 ibbl=ibbl+1
3195#endif
3196#if defined UV_TIDES
3197!
3198 IF (master) WRITE (stdout,20) 'UV_TIDES', &
3199 & 'Add tidal currents to 2D momentum climatologies'
3200 is=len_trim(coptions)+1
3201 coptions(is:is+10)=' UV_TIDES,'
3202#endif
3203#ifdef UV_VIS2
3204!
3205 IF (master) WRITE (stdout,20) 'UV_VIS2', &
3206 & 'Harmonic mixing of momentum'
3207 is=len_trim(coptions)+1
3208 coptions(is:is+9)=' UV_VIS2,'
3209#endif
3210#ifdef UV_VIS4
3211!
3212 IF (master) WRITE (stdout,20) 'UV_VIS4', &
3213 & 'Biharmonic mixing of momentum'
3214 is=len_trim(coptions)+1
3215 coptions(is:is+9)=' UV_VIS4,'
3216#endif
3217#if defined UV_SMAGORINSKY && (defined UV_VIS2 || defined UV_VIS4)
3218!
3219 IF (master) WRITE (stdout,20) 'UV_SMAGORINSKY', &
3220 & 'Smagorinksy-like time-dependent viscosity coefficients'
3221 is=len_trim(coptions)+1
3222 coptions(is:is+16)=' UV_SMAGORINSKY,'
3223#endif
3224#if defined VAR_RHO_2D && defined SOLVE3D
3225!
3226 IF (master) WRITE (stdout,20) 'VAR_RHO_2D', &
3227 & 'Variable density barotropic mode'
3228 is=len_trim(coptions)+1
3229 coptions(is:is+12)=' VAR_RHO_2D,'
3230#endif
3231#if defined UV_VIS2 || defined UV_VIS4
3232# ifdef VISC_GRID
3233!
3234 IF (master) WRITE (stdout,20) 'VISC_GRID', &
3235 & 'Horizontal viscosity coefficient scaled by grid size'
3236 is=len_trim(coptions)+1
3237 coptions(is:is+11)=' VISC_GRID,'
3238# endif
3239#endif
3240#if defined VCONVOLUTION && defined FOUR_DVAR && defined SOLVE3D
3241!
3242 IF (master) WRITE (stdout,20) 'VCONVOLUTION', &
3243 & 'Include vertical correlations in convolutions'
3244 is=len_trim(coptions)+1
3245 coptions(is:is+14)=' VCONVOLUTION,'
3246#endif
3247#ifdef VERIFICATION
3248!
3249 IF (master) WRITE (stdout,20) 'VERIFICATION', &
3250 & 'Proccess model solution at observation locations'
3251 is=len_trim(coptions)+1
3252 coptions(is:is+14)=' VERIFICATION,'
3253#endif
3254#if defined FLOATS && defined FLOAT_VWALK && defined SOLVE3D
3255# ifdef VWALK_FORWARD
3256!
3257 IF (master) WRITE (stdout,20) 'VWALK_FORWARD', &
3258 & 'Lagrangian drifters, forward stepping random walk'
3259 is=len_trim(coptions)+1
3260 coptions(is:is+15)=' VWALK_FORWARD,'
3261# else
3262!
3263 IF (master) WRITE (stdout,20) '!VWALK_FORWARD', &
3264 & 'Lagrangian drifters, predictor/corrector random walk'
3265 is=len_trim(coptions)+1
3266 coptions(is:is+16)=' !VWALK_FORWARD,'
3267# endif
3268#endif
3269#if defined VISC_3DCOEF
3270!
3271 IF (master) WRITE (stdout,20) 'VISC_3DCOEF', &
3272 & 'Horizontal, time-dependent 3D viscosity coefficient'
3273 is=len_trim(coptions)+1
3274 coptions(is:is+13)=' VISC_3DCOEF,'
3275#endif
3276#ifdef WAVE_MIXING
3277!
3278 IF (master) WRITE (stdout,20) 'WAVE_MIXING', &
3279 & 'Enhanced veritcal mixing from waves term'
3280 is=len_trim(coptions)+1
3281 coptions(is:is+13)=' WAVE_MIXING,'
3282#endif
3283#if defined WDISS_THORGUZA && defined WEC_VF
3284!
3285 IF (master) WRITE (stdout,20) 'WDISS_THORGUZA', &
3286 & 'Wave energy dissipation computed from Thorton/Guza 1986'
3287 is=len_trim(coptions)+1
3288 coptions(is:is+16)=' WDISS_THORGUZA,'
3289#endif
3290#if defined WDISS_CHURTHOR && defined WEC_VF
3291!
3292 IF (master) WRITE (stdout,20) 'WDISS_CHURTHOR', &
3293 & 'Wave energy dissipation computed from Church/Thorton 1986'
3294 is=len_trim(coptions)+1
3295 coptions(is:is+16)=' WDISS_CHURTHOR,'
3296#endif
3297#ifdef WDISS_WAVEMOD
3298!
3299 IF (master) WRITE (stdout,20) 'WDISS_WAVEMOD', &
3300 & 'Wave energy dissipation acquired from coupled wave model'
3301 is=len_trim(coptions)+1
3302 coptions(is:is+15)=' WDISS_WAVEMOD,'
3303#endif
3304#ifdef WAV_COUPLING
3305!
3306 IF (master) WRITE (stdout,20) 'WAV_COUPLING', &
3307 & 'Two-way wave-ocean models coupling'
3308 is=len_trim(coptions)+1
3309 coptions(is:is+13)=' WAV_COUPLING,'
3310#endif
3311#ifdef WEAK_CONSTRAINT
3312!
3313 IF (master) WRITE (stdout,20) 'WEAK_CONSTRAINT', &
3314 & 'Activated weak constraint assimilation set-up'
3315 is=len_trim(coptions)+1
3316 coptions(is:is+17)=' WEAK_CONSTRAINT,'
3317#endif
3318#if defined WEAK_NOINTERP && defined WEAK_CONSTRAINT
3319!
3320 IF (master) WRITE (stdout,20) 'WEAK_NOINTERP', &
3321 & 'Avoid time interpolation in weak constraint forcing'
3322 is=len_trim(coptions)+1
3323 coptions(is:is+15)=' WEAK_NOINTERP,'
3324#endif
3325#ifdef WEC_VF
3326!
3327 IF (master) WRITE (stdout,20) 'WEC_VF', &
3328 & 'Vortex Force wave current interaction'
3329 is=len_trim(coptions)+1
3330 coptions(is:is+8)=' WEC_VF,'
3331 nearshore=nearshore+1
3332#endif
3333#if defined WRITE_WATER && defined MASKING
3334!
3335 IF (master) WRITE (stdout,20) 'WRITE_WATER', &
3336 & 'Writing data at water points only'
3337 is=len_trim(coptions)+1
3338 coptions(is:is+13)=' WRITE_WATER,'
3339#endif
3340#ifdef WET_DRY
3341!
3342 IF (master) WRITE (stdout,20) 'WET_DRY', &
3343 & 'Wetting and drying activated'
3344 is=len_trim(coptions)+1
3345 coptions(is:is+9)=' WET_DRY,'
3346#endif
3347#ifdef WIND_MINUS_CURRENT && defined BULK_FLUXES && defined SOLVE3D
3348!
3349 IF (master) WRITE (stdout,20) 'WIND_MINUS_CURRENT', &
3350 & 'Compute effective wind by removing surface ocean current'
3351 is=len_trim(coptions)+1
3352 coptions(is:is+20)=' WIND_MINUS_CURRENT,'
3353#endif
3354#ifdef WIND_WAVES
3355!
3356 IF (master) WRITE (stdout,20) 'WIND_WAVES', &
3357 & 'Wind-induced wave data available'
3358 is=len_trim(coptions)+1
3359 coptions(is:is+11)=' WIND_WAVES,'
3360#endif
3361#if !(defined PJ_GRADPQ4 || defined PJ_GRADPQ2 || defined PJ_GRADP || \
3362 defined dj_gradps) && defined wj_gradp
3363!
3364 IF (master) WRITE (stdout,20) 'WJ_GRADP', &
3365 & 'Weighted Jacobians pressure gradient adjustment'
3366 is=len_trim(coptions)+1
3367 coptions(is:is+10)=' WJ_GRADP,'
3368#endif
3369#ifdef WRF_COUPLING
3370!
3371 IF (master) WRITE (stdout,20) 'WRF_COUPLING', &
3372 & 'Atmosphere coupling with WRF'
3373 is=len_trim(coptions)+1
3374 coptions(is:is+14)=' WRF_COUPLING,'
3375 iatms=iatms+1
3376# ifdef WRF_TIMEAVG
3377!
3378 IF (master) WRITE (stdout,20) 'WRF_TIMEAVG', &
3379 & 'WRF Exports time-averaged fields over the coupling interval'
3380 is=len_trim(coptions)+1
3381 coptions(is:is+13)=' WRF_TIMEAVG,'
3382# endif
3383#endif
3384#if defined WTYPE_GRID && \
3385 (defined lmd_skpp || defined solar_source) && \
3386 !defined ANA_WTYPE
3387!
3388 IF (master) WRITE (stdout,20) 'WTYPE_GRID', &
3389 & 'Spatially varying Jerlov water type index'
3390 is=len_trim(coptions)+1
3391 coptions(is:is+12)=' WTYPE_GRID,'
3392#endif
3393#if defined ZETA_ELLIPTIC && defined BALANCE_OPERATOR && \
3394 defined four_dvar
3395!
3396 IF (master) WRITE (stdout,20) 'ZETA_ELLIPTIC', &
3397 & 'Solving SSH elliptic equation in balance operator'
3398 is=len_trim(coptions)+1
3399 coptions(is:is+15)=' ZETA_ELLIPTIC,'
3400#endif
3401#ifdef ZOS_HSIG
3402!
3403 IF (master) WRITE (stdout,20) 'ZOS_HSIG', &
3404 & 'Wave amplitude used for Zos calculation'
3405 is=len_trim(coptions)+1
3406 coptions(is:is+10)=' ZOS_HSIG,'
3407#endif
3408!
3409!-----------------------------------------------------------------------
3410! Stop if unsupported C-preprocessing options or report issues with
3411! particular options.
3412!-----------------------------------------------------------------------
3413!
3414 CALL checkadj
3415!
3416!-----------------------------------------------------------------------
3417! Check C-preprocessing options.
3418!-----------------------------------------------------------------------
3419!
3420! Stop if more than one vertical closure scheme is selected.
3421!
3422 IF (master.and.(ivmix.gt.1)) THEN
3423 WRITE (stdout,30)
3424 30 FORMAT (/,' CHECKDEFS - only one vertical closure scheme', &
3425 & ' is allowed.')
3426 exit_flag=5
3427 END IF
3428!
3429! Stop if more that one bottom stress formulation is selected.
3430!
3431 IF (master.and.(ibbl.gt.1)) THEN
3432 WRITE (stdout,40)
3433 40 FORMAT (/,' CHECKDEFS - only one bottom stress formulation is', &
3434 & ' allowed.')
3435 exit_flag=5
3436 END IF
3437!
3438! Stop if no bottom stress formulation is selected.
3439!
3440 IF (master.and.(ibbl.eq.0)) THEN
3441 WRITE (stdout,50)
3442 50 FORMAT (/,' CHECKDEFS - no bottom stress formulation is', &
3443 & ' selected.')
3444 exit_flag=5
3445 END IF
3446!
3447! Stop if more than one biological module is selected.
3448!
3449 IF (master.and.(ibiology.gt.1)) THEN
3450 WRITE (stdout,60)
3451 60 FORMAT (/,' CHECKDEFS - only one biology MODULE is allowed.')
3452 exit_flag=5
3453 END IF
3454!
3455! Stop if more that one model driver is selected.
3456!
3457 IF (master.and.(idriver.gt.1)) THEN
3458 WRITE (stdout,70)
3459 70 FORMAT (/,' CHECKDEFS - only one model example is allowed.')
3460 exit_flag=5
3461 END IF
3462
3463#ifndef SOLVE3D
3464!
3465! Stop it explicit time-step splitting on shallow water set-up.
3466!
3467 DO ng=1,ngrids
3468 IF (master.and.(ndtfast(ng).gt.1)) THEN
3469 WRITE (stdout,80)
3470 80 FORMAT (/,' CHECKDEFS - explicit time-step splitting is ', &
3471 & ' inconsistent.', &
3472 & /,13x,'Change parameter NDTFAST to unity.')
3473 exit_flag=5
3474 END IF
3475 END DO
3476#endif
3477#if defined ADJOINT && defined _OPENMP
3478 IF (master) THEN
3479 WRITE (stdout,90)
3480 90 FORMAT (/,' CHECKDEFS - cannot run the adjoint model in', &
3481 & ' shared-memory mode.')
3482 exit_flag=5
3483 END IF
3484#endif
3485#if defined NEMURO && defined HOLLING_GRAZING && defined IVLEV_EXPLICIT
3486!
3487! Stop if using more than one grazing algorithm.
3488!
3489 IF (master) THEN
3490 WRITE (stdout,100)
3491 100 FORMAT (/,' CHECKDEFS - More than one grazing algorithm', &
3492 & ' selected for Nemuro model.')
3493 exit_flag=5
3494 END IF
3495#endif
3496#if defined FOUR_DVAR && defined VCONVOLUTION
3497# if defined IMPLICIT_VCONV && defined SPLINES_VCONV
3498!
3499! Stop if using more than vertical convolution algorithm.
3500!
3501 IF (master) THEN
3502 WRITE (stdout,110)
3503 110 FORMAT (/,' CHECKDEFS - More than one vertical convolution', &
3504 & ' algorithm selected.')
3505 exit_flag=5
3506 END IF
3507# endif
3508#endif
3509#if defined MODEL_COUPLING && (defined ESMF_LIB && defined MCT_LIB)
3510 IF (master) THEN
3511 WRITE (stdout,120)
3512 120 FORMAT (/,' CHECKDEFS - More than one coupling library', &
3513 & ' selected.')
3514 exit_flag=5
3515 END IF
3516#endif
3517#if defined POSTERIOR_ERROR_F && defined POSTERIOR_ERROR_I && \
3518 defined weak_constraint
3519 IF (master) THEN
3520 WRITE (stdout,130)
3521 130 FORMAT (/,' CHECKDEFS - Both initial and final posterior', &
3522 & ' error covariance analysis selected, choose one.')
3523 exit_flag=5
3524 END IF
3525#endif
3526#ifdef UV_ADV
3527!
3528! Stop if more than one advection scheme has been activated.
3529!
3530 IF (master.and.(ivelhadv.gt.1)) THEN
3531 WRITE (stdout,140) 'horizontal','momentum','ivelHadv =',ivelhadv
3532 exit_flag=5
3533 END IF
3534# ifdef SOLVE3D
3535 IF (master.and.(ivelvadv.gt.1)) THEN
3536 WRITE (stdout,140) 'vertical','momentum','ivelVadv =',ivelvadv
3537 exit_flag=5
3538 END IF
3539# endif
3540#endif
3541 140 FORMAT (/,' CHECKDEFS - only one ',a,' advection scheme', &
3542 & /,13x,'is allowed for ',a,', ',a,1x,i1)
3543#if defined SOLVE3D && defined NEMURO
3544# if !(defined IVLEV_EXPLICIT || defined HOLLING_GRAZING)
3545 IF (master) THEN
3546 WRITE (stdout,150) uppercase('ivlev_explicit'), &
3547 & uppercase('holling_grazing')
3548 150 FORMAT (/,' CHECKDEFS - Need to choose a zooplankton grazing', &
3549 & ' option;',/,13x,'use either ',a,' or implicit ',a, &
3550 & '.',/,13x,'The default implicit IVLEV algorithm', &
3551 & ' does not work well yet.')
3552 exit_flag=5
3553 END IF
3554# endif
3555#endif
3556#ifdef WEC
3557!
3558! Stop if more that one radiation stress formulation is activated.
3559!
3560 IF (master.and.(nearshore.gt.1)) THEN
3561 WRITE (stdout,160)
3562 160 FORMAT (/,' CHECKDEFS - only one wave current formulation' &
3563 & ' is allowed.')
3564 exit_flag=5
3565 END IF
3566#endif
3567#ifdef SEDIMENT
3568 IF (master.and.(sbed_type.gt.1)) THEN
3569 WRITE (stdout,170)
3570 170 FORMAT (/,' CHECKDEFS - More than one sediment bed type ', &
3571 & ' selected.')
3572 exit_flag=5
3573 END IF
3574 IF (master.and.(sbed_load.gt.1)) THEN
3575 WRITE (stdout,180)
3576 180 FORMAT (/,' CHECKDEFS - More than one sediment bed load ', &
3577 & ' selected.')
3578 exit_flag=5
3579 END IF
3580#endif
3581#ifdef PARALLEL_IO
3582# if !(defined HDF5 || defined PNETCDF)
3583 IF (master) THEN
3584 WRITE (stdout,190) uppercase('hdf5'), &
3585 & uppercase('pnetcdf')
3586 190 FORMAT (/,' CHECKDEFS - Need to activate either ',a,' or ',a, &
3587 & ' options',/,13x,'for parallel I/O processing.')
3588 exit_flag=5
3589 END IF
3590# endif
3591# if defined HDF5 && defined PNETCDF
3592 IF (master) THEN
3593 WRITE (stdout,200) uppercase('hdf5'), &
3594 & uppercase('pnetcdf')
3595 200 FORMAT (/,' CHECKDEFS - Both ',a,' and ',a,' options', &
3596 & ' selected',/,13x,'for parallel I/O processing.', &
3597 & /,13x,'Choose only one.')
3598 exit_flag=5
3599 END IF
3600# endif
3601#endif
3602#if defined WRITE_WATER && defined WET_DRY
3603!
3604! Stop if writting only water points when wetting and drying is
3605! activated. This is not possible since the Land/Sea masking is
3606! time dependent.
3607!
3608 IF (master) THEN
3609 WRITE (stdout,210) uppercase('write_water'), &
3610 & uppercase('wet_dry')
3611 210 FORMAT (/,' CHECKDEFS - cannot activate ',a,' and ',a, &
3612 & ' together',/,13x,'because of time dependent ', &
3613 & ' Land/Sea masking.')
3614 exit_flag=5
3615 END IF
3616#endif
3617#if ((defined AD_AVERAGES && defined ADJOINT) && defined AVERAGES) || \
3618 ((defined rp_averages && defined tl_ioms) && defined averages) || \
3619 ((defined tl_averages && defined tangent) && defined averages) || \
3620 (defined ad_averages && defined tl_averages) || \
3621 (defined ad_averages && defined rp_averages) || \
3622 (defined rp_averages && defined tl_averages)
3623!
3624! Stop if activating more that one time-averaged option. This is not
3625! possible in the same executable because the same internal arrays
3626! are used to store the cummulative sums.
3627!
3628 IF (master) THEN
3629 WRITE (stdout,220) uppercase('averages'), &
3630 & uppercase('ad_averages'), &
3631 & uppercase('rp_averages'), &
3632 & uppercase('tl_averages')
3633 220 FORMAT (/,' CHECKDEFS - cannot activate ',a,' ',a,' ',a, &
3634 & ' or',a,/,13x,'at the same time. Choose one!')
3635 exit_flag=5
3636 END IF
3637#endif
3638#if !defined DISTRIBUTE && defined TS_MPDATA
3639!
3640! Stop if activating MPDATA in serial with partitions or shared-
3641! memory.
3642!
3643 IF (master) THEN
3644 DO ng=1,ngrids
3645 IF (ntilex(ng)*ntilee(ng).gt.1) THEN
3646 WRITE (stdout,230) uppercase('ts_mpdata')
3647 exit_flag=5
3648 EXIT
3649 END IF
3650 END DO
3651 230 FORMAT (/,' CHECKDEFS - cannot activate option: ',a, &
3652 & /,13x,'in serial with partitions or shared-memory...', &
3653 & /,13x,'Use distributed-memory (MPI) in parallel runs.')
3654 END IF
3655#endif
3656#ifdef Q_PSOURCE
3657!
3658! Stop if activating obsolete mass point Sources/Sinks option. This
3659! capability is now activated with standard input switch "LwSrc".
3660!
3661 IF (master) THEN
3662 WRITE (stdout,240) uppercase('q_psource'), 'LwSrc'
3663 240 FORMAT (/,' CHECKDEFS - cannot use obsolete option: ',a, &
3664 & /,13x,'Use instead standard input switch: ',a, &
3665 & /,13x,'Edit header file or build script!')
3666 exit_flag=5
3667 END IF
3668#endif
3669#if defined TS_PSOURCE && defined SOLVE3D
3670!
3671! Stop if activating obsolete tracers point Sources/Sinks option. This
3672! capability is now activated with standard input switch "LtracerSrc".
3673!
3674 IF (master) THEN
3675 WRITE (stdout,250) uppercase('ts_psource'), 'LtracerSrc'
3676 250 FORMAT (/,' CHECKDEFS - cannot use obsolete option: ',a, &
3677 & /,13x,'Use instead standard input switch: ',a, &
3678 & /,13x,'Edit header file or build script!')
3679 exit_flag=5
3680 END IF
3681#endif
3682#ifdef UV_PSOURCE
3683!
3684! Stop if activating obsolete mass point Sources/Sinks option. This
3685! capability is now activated with standard input switch "LwSrc".
3686!
3687 IF (master) THEN
3688 WRITE (stdout,260) uppercase('uv_psource'), 'LuvSrc'
3689 260 FORMAT (/,' CHECKDEFS - cannot use obsolete option: ',a, &
3690 & /,13x,'Use instead standard input switch: ',a, &
3691 & /,13x,'Edit header file or build script!')
3692 exit_flag=5
3693 END IF
3694#endif
3695#if defined SPLINES && defined SOLVE3D
3696!
3697! Stop if activating parabolic splines reconstruction with deprecated
3698! option. See https://www.myroms.org/projects/src/ticket/681 for more
3699! details.
3700!
3701 IF (master) THEN
3702 WRITE (stdout,270) uppercase('splines'), &
3703 & 'https://www.myroms.org/projects/src/ticket/681'
3704 270 FORMAT (/,' CHECKDEFS - cannot use obsolete option: ',a, &
3705 & /,13x,'See following trac ticket for details: ', &
3706 & /,13x,a,/,13x,'please change header file accordingly.')
3707 exit_flag=5
3708 END IF
3709#endif
3710#if !defined DISTRIBUTE && defined NESTING
3711!
3712! Stop if activating nesting in serial with partitions or shared-
3713! memory.
3714!
3715 IF (master) THEN
3716 DO ng=1,ngrids
3717 IF (ntilex(ng)*ntilee(ng).gt.1) THEN
3718 WRITE (stdout,280) uppercase('nesting')
3719 exit_flag=5
3720 EXIT
3721 END IF
3722 END DO
3723 280 FORMAT (/,' CHECKDEFS - cannot activate option: ',a, &
3724 & /,13x,'in serial with partitions or shared-memory...', &
3725 & /,13x,'We have bugs in shared-memory that need fixing.',&
3726 & /,13x,'Use distributed-memory (MPI) in parallel runs.')
3727 END IF
3728#endif
3729#ifdef ATM_COUPLING
3730!
3731! Stop if more than one atmosphere model is selected for coupling.
3732!
3733 IF (master.and.(iatms.gt.1)) THEN
3734 WRITE (stdout,290)
3735 290 FORMAT (/,' CHECKDEFS - only one atmosphere model coupling', &
3736 & ' is allowed.')
3737 exit_flag=5
3738 END IF
3739#endif
3740#if defined TIDE_GENERATING_FORCES && \
3741 (defined ana_grid && !defined SPHERICAL)
3742!
3743! Stop if activating tide generating forces on analytical Cartesian
3744! grids. This forcing is only for spherical grids since the (lon,lat)
3745! is required to compute the equilibrium tide harmonics that depend
3746! on the geographical location with respect to the Moon and Sun.
3747!
3748 IF (master) THEN
3749 WRITE (stdout,300) uppercase('tide_generating_forces')
3750 exit_flag=5
3751 EXIT
3752 END IF
3753 END DO
3754 300 FORMAT (/,' CHECKDEFS - cannot activate option: ',a, &
3755 & /,13x,'in analytical Cartesian grids applications...', &
3756 & /,13x,'The grid (lon,lat) is needed to compute the', &
3757 & /,13x,'equilibrium tide harmonics.')
3758 END IF
3759#endif
3760!
3761 RETURN
3762 END SUBROUTINE checkdefs
subroutine ana_grid(ng, tile, model)
Definition ana_grid.h:3
subroutine checkadj
Definition checkadj.F:3
subroutine checkdefs
Definition checkdefs.F:3
subroutine lmd_skpp(ng, tile)
Definition lmd_skpp.F:45
character(len=256) myappcpp
integer stdout
logical master
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
integer ngrids
Definition mod_param.F:113
integer, dimension(:), allocatable ntilee
Definition mod_param.F:686
logical, dimension(:), allocatable lsponge
logical, dimension(:), allocatable lnudging
integer, dimension(:), allocatable ndtfast
integer exit_flag
character(len=80) title
character(len=2048) coptions
character(len(sinp)) function, public uppercase(sinp)
Definition strings.F:582