ROMS
Loading...
Searching...
No Matches
def_gst.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined PROPAGATOR && defined CHECKPOINTING
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This module creates the checkpointing restart file for the GST !
13! (Generalized Stability Theory) analysis drivers using either the !
14! standard NetCDF library or the Parallel-IO (PIO) library. !
15! !
16! The NetCDF file contains all the necessary eigenproblem fields !
17! for restating ARPACK. !
18! !
19!=======================================================================
20!
21 USE mod_param
22 USE mod_parallel
23 USE mod_iounits
24 USE mod_ncparam
25 USE mod_scalars
26 USE mod_storage
27!
28 USE def_dim_mod, ONLY : def_dim
29 USE def_info_mod, ONLY : def_info
30 USE def_var_mod, ONLY : def_var
31 USE strings_mod, ONLY : founderror
32!
33 implicit none
34!
35 PUBLIC :: def_gst
36 PRIVATE :: def_gst_nf90
37# if defined PIO_LIB && defined DISTRIBUTE
38 PRIVATE :: def_gst_pio
39# endif
40!
41 CONTAINS
42!
43!***********************************************************************
44 SUBROUTINE def_gst (ng, model)
45!***********************************************************************
46!
47! Imported variable declarations.
48!
49 integer, intent(in) :: ng, model
50!
51! Local variable declarations.
52!
53 character (len=*), parameter :: myfile = &
54 & __FILE__
55!
56!-----------------------------------------------------------------------
57! Create a new history file according to IO type.
58!-----------------------------------------------------------------------
59!
60 SELECT CASE (gst(ng)%IOtype)
61 CASE (io_nf90)
62 CALL def_gst_nf90 (ng, model)
63
64# if defined PIO_LIB && defined DISTRIBUTE
65 CASE (io_pio)
66 CALL def_gst_pio (ng, model)
67# endif
68 CASE DEFAULT
69 IF (master) THEN
70 WRITE (stdout,10) gst(ng)%IOtype
71 END IF
72 exit_flag=3
73 END SELECT
74 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
75!
76 10 FORMAT (' DEF_GST - Illegal output file type, io_type = ',i0, &
77 & /,15x,'Check KeyWord ''OUT_LIB'' in ''roms.in''.')
78!
79 RETURN
80 END SUBROUTINE def_gst
81!
82!***********************************************************************
83 SUBROUTINE def_gst_nf90 (ng, model)
84!***********************************************************************
85!
86 USE mod_netcdf
87!
88! Imported variable declarations.
89!
90 integer, intent(in) :: ng, model
91!
92! Local variable declarations.
93!
94 integer, parameter :: natt = 25
95!
96 integer :: mstatedim, ncvdim, sworkddim, sworkldim
97 integer :: char1dim, char2dim, iaitrdim, iaupddim, iaup2dim
98 integer :: iparamdim, ipntrdim, laitrdim, laup2dim, raitrdim
99 integer :: raup2dim
100 integer :: i, j, status, varid
101 integer :: dimids(ndimid), vardim(2)
102!
103 real(r8) :: aval(6)
104!
105 character (len=256) :: ncname
106 character (len=MaxLen) :: vinfo(natt)
107
108 character (len=*), parameter :: myfile = &
109 & __FILE__//", def_gst_nf90"
110!
111 sourcefile=myfile
112!
113!-----------------------------------------------------------------------
114! Set and report file name.
115!-----------------------------------------------------------------------
116!
117 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
118 ncname=gst(ng)%name
119!
120 IF (master) THEN
121 IF (.not.lrstgst) THEN
122 WRITE (stdout,10) ng, trim(ncname)
123 ELSE
124 WRITE (stdout,20) ng, trim(ncname)
125 END IF
126 END IF
127!
128!=======================================================================
129! Create a new averages NetCDF file.
130!=======================================================================
131!
132 define : IF (.not.lrstgst) THEN
133 CALL netcdf_create (ng, model, trim(ncname), gst(ng)%ncid)
134 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
135 IF (master) WRITE (stdout,30) trim(ncname)
136 RETURN
137 END IF
138!
139!-----------------------------------------------------------------------
140! Define the dimensions of staggered fields.
141!-----------------------------------------------------------------------
142!
143 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'Mstate', &
144 & mstate(ng), mstatedim)
145 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
146
147 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'NCV', &
148 & ncv, ncvdim)
149 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
150
151 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'LworkD', &
152 & 3*mstate(ng), sworkddim)
153 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
154
155 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'LworkL', &
156 & lworkl, sworkldim)
157 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
158
159 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'iparam', &
160 & SIZE(iparam), iparamdim)
161 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
162
163 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'ipntr', &
164 & SIZE(ipntr), ipntrdim)
165 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
166
167 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'iaupd', &
168 & SIZE(iaupd), iaupddim)
169 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
170
171 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'laitr', &
172 & SIZE(laitr), laitrdim)
173 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
174
175 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'iaitr', &
176 & SIZE(iaitr), iaitrdim)
177 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
178
179 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'raitr', &
180 & SIZE(raitr), raitrdim)
181 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
182
183 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'laup2', &
184 & SIZE(laup2), laup2dim)
185 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
186
187 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'iaup2', &
188 & SIZE(iaup2), iaup2dim)
189 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
190
191 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'raup2', &
192 & SIZE(raup2), raup2dim)
193 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
194
195 status=def_dim(ng, model, gst(ng)%ncid, ncname, 'char2', &
196 & 2, char2dim)
197 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
198!
199!-----------------------------------------------------------------------
200! Define global attributes.
201!-----------------------------------------------------------------------
202!
203 CALL def_info (ng, model, gst(ng)%ncid, ncname, dimids)
204 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
205!
206!-----------------------------------------------------------------------
207! Define variables and their attributes.
208!-----------------------------------------------------------------------
209!
210! Initialize local information variable arrays.
211!
212 DO i=1,natt
213 DO j=1,len(vinfo(1))
214 vinfo(i)(j:j)=' '
215 END DO
216 END DO
217 DO i=1,6
218 aval(i)=0.0_r8
219 END DO
220!
221! Define number of eigenvalues to compute.
222!
223 vinfo( 1)='NEV'
224 vinfo( 2)='number of eigenvalues to compute'
225 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
226 & 1, (/0/), aval, vinfo, ncname, &
227 & setparaccess = .false.)
228 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
229!
230! Define number of Lanczos vectors to compute.
231!
232 vinfo( 1)='NCV'
233 vinfo( 2)='number of Lanczos vectors to compute'
234 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
235 & 1, (/0/), aval, vinfo, ncname, &
236 & setparaccess = .false.)
237 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
238!
239! Define size of eigenvalue problem.
240!
241 vinfo( 1)='Mstate'
242 vinfo( 2)='total size of eigenvalue problem'
243 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
244 & 1, (/0/), aval, vinfo, ncname, &
245 & setparaccess = .false.)
246 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
247!
248! Define iteration number.
249!
250 vinfo( 1)='iter'
251 vinfo( 2)='iteration number'
252 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
253 & 1, (/0/), aval, vinfo, ncname, &
254 & setparaccess = .false.)
255 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
256!
257! Define reverse communications flag.
258!
259 vinfo( 1)='ido'
260 vinfo( 2)='reverse communications flag'
261 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
262 & 1, (/0/), aval, vinfo, ncname, &
263 & setparaccess = .false.)
264 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
265!
266! Define information and error flag.
267!
268 vinfo( 1)='info'
269 vinfo( 2)='information and error flag'
270 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
271 & 1, (/0/), aval, vinfo, ncname, &
272 & setparaccess = .false.)
273 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
274!
275! Define eigenvalue problem type.
276!
277 vinfo( 1)='bmat'
278 vinfo( 2)='eigenvalue problem type'
279 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_char, &
280 & 1, (/0/), aval, vinfo, ncname, &
281 & setparaccess = .false.)
282 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
283!
284! Define Ritz eigenvalues to compute.
285!
286 vinfo( 1)='which'
287 vinfo( 2)='Ritz eigenvalues to compute'
288 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_char, &
289 & 1, (/char2dim/), aval, vinfo, ncname, &
290 & setparaccess = .false.)
291 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
292!
293! Define form of basis function.
294!
295 vinfo( 1)='howmany'
296 vinfo( 2)='form of basis function'
297 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_char, &
298 & 1, (/0/), aval, vinfo, ncname, &
299 & setparaccess = .false.)
300 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
301!
302! Define relative accuracy of computed Ritz values.
303!
304 vinfo( 1)='Ritz_tol'
305 vinfo( 2)='relative accuracy of computed Ritz values'
306 status=def_var(ng, model, gst(ng)%ncid, varid, nf_frst, &
307 & 1, (/0/), aval, vinfo, ncname, &
308 & setparaccess = .false.)
309 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
310!
311! Define eigenproblem parameters.
312!
313 vinfo( 1)='iparam'
314 vinfo( 2)='eigenproblem parameters'
315 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
316 & 1, (/iparamdim/), aval, vinfo, ncname, &
317 & setparaccess = .false.)
318 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
319!
320! Define pointers to mark starting location in work arrays.
321!
322 vinfo( 1)='ipntr'
323 vinfo( 2)='pointers to mark starting location in work arrays'
324 status=def_var(ng, model, gst(ng)%ncid,varid, nf90_int, &
325 & 1, (/ipntrdim/), aval, vinfo, ncname, &
326 & setparaccess = .false.)
327 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
328!
329! Define ARPACK internal integer parameters to _aupd routines.
330!
331 vinfo( 1)='iaupd'
332 vinfo( 2)='ARPACK internal integer parameters to _aupd routines'
333 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
334 & 1, (/iaupddim/), aval, vinfo, ncname, &
335 & setparaccess = .false.)
336 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
337!
338! Define ARPACK internal integer parameters to _aitr routines.
339!
340 vinfo( 1)='iaitr'
341 vinfo( 2)='ARPACK internal integer parameters to _aitr routines'
342 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
343 & 1, (/iaitrdim/), aval, vinfo, ncname, &
344 & setparaccess = .false.)
345 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
346!
347! Define ARPACK internal integer parameters to _aup2 routines.
348!
349 vinfo( 1)='iaup2'
350 vinfo( 2)='ARPACK internal integer parameters to _aup2 routines'
351 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
352 & 1, (/iaup2dim/), aval, vinfo, ncname, &
353 & setparaccess = .false.)
354 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
355!
356! Define ARPACK internal logical parameters to _aitr routines.
357!
358 vinfo( 1)='laitr'
359 vinfo( 2)='ARPACK internal logical parameters to _aitr routines'
360 vinfo( 9)='.FALSE.'
361 vinfo(10)='.TRUE.'
362 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
363 & 1, (/laitrdim/), aval, vinfo, ncname, &
364 & setparaccess = .false.)
365 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
366!
367! Define ARPACK internal logical parameters to _aup2 routines.
368!
369 vinfo( 1)='laup2'
370 vinfo( 2)='ARPACK internal logical parameters to _aup2 routines'
371 vinfo( 9)='.FALSE.'
372 vinfo(10)='.TRUE.'
373 status=def_var(ng, model, gst(ng)%ncid, varid, nf90_int, &
374 & 1, (/laup2dim/), aval, vinfo, ncname, &
375 & setparaccess = .false.)
376 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
377!
378! Define ARPACK internal real parameters to _aitr routines.
379!
380 vinfo( 1)='raitr'
381 vinfo( 2)='ARPACK internal real parameters to _aitr routines'
382 status=def_var(ng, model, gst(ng)%ncid, varid, nf_frst, &
383 & 1, (/raitrdim/), aval, vinfo, ncname, &
384 & setparaccess = .false.)
385 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
386!
387! Define ARPACK internal real parameters to _aup2 routines.
388!
389 vinfo( 1)='raup2'
390 vinfo( 2)='ARPACK internal real parameters to _aup2 routines'
391 status=def_var(ng, model, gst(ng)%ncid, varid, nf_frst, &
392 & 1, (/raup2dim/), aval, vinfo, ncname, &
393 & setparaccess = .false.)
394 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
395!
396! Define Lanczos/Arnoldi basis vectors.
397!
398 vinfo( 1)='Bvec'
399 vinfo( 2)='Lanczos/Arnoldi basis vectors'
400 vardim(1)=mstatedim
401 vardim(2)=ncvdim
402 status=def_var(ng, model, gst(ng)%ncid, varid, nf_frst, &
403 & 2, vardim, aval, vinfo, ncname)
404 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
405!
406! Define eigenproblem residual vector.
407!
408 vinfo( 1)='resid'
409 vinfo( 2)='eigenproblem residual vector'
410 status=def_var(ng, model, gst(ng)%ncid, varid, nf_frst, &
411 & 1, (/mstatedim/), aval, vinfo, ncname)
412 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
413!
414! Define state reverse communications work array.
415!
416 vinfo( 1)='SworkD'
417 vinfo( 2)='reverse communications state array'
418 status=def_var(ng, model, gst(ng)%ncid, varid, nf_frst, &
419 & 1, (/sworkddim/), aval, vinfo, ncname)
420 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
421!
422! Define eigenproblem work array.
423!
424 vinfo( 1)='SworkL'
425 vinfo( 2)='eigenproblem work array'
426 status=def_var(ng, model, gst(ng)%ncid, varid, nf_frst, &
427 & 1, (/sworkldim/), aval, vinfo, ncname)
428 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
429!
430!-----------------------------------------------------------------------
431! Leave definition mode.
432!-----------------------------------------------------------------------
433!
434 CALL netcdf_enddef (ng, model, ncname, gst(ng)%ncid)
435 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
436
437 END IF define
438!
439 10 FORMAT (/,2x,'DEF_GST_NF90 - creating checkpointing file,', &
440 & t56,'Grid ',i2.2,': ',a)
441 20 FORMAT (/,2x,'DEF_GST_NF90 - inquiring checkpointing file,', &
442 & t56,'Grid ',i2.2,': ',a)
443 30 FORMAT (/,' DEF_GST_NF90 - unable to create checkpointing', &
444 & ' NetCDF file: ',a)
445!
446 RETURN
447 END SUBROUTINE def_gst_nf90
448
449# if defined PIO_LIB && defined DISTRIBUTE
450!
451!***********************************************************************
452 SUBROUTINE def_gst_pio (ng, model)
453!***********************************************************************
454!
456!
457! Imported variable declarations.
458!
459 integer, intent(in) :: ng, model
460!
461! Local variable declarations.
462!
463 integer, parameter :: natt = 25
464!
465 integer :: mstatedim, ncvdim, sworkddim, sworkldim
466 integer :: char1dim, char2dim, iaitrdim, iaupddim, iaup2dim
467 integer :: iparamdim, ipntrdim, laitrdim, laup2dim, raitrdim
468 integer :: raup2dim
469 integer :: i, j, status
470 integer :: dimids(ndimid), vardim(2)
471!
472 real(r8) :: aval(6)
473!
474 character (len=256) :: ncname
475 character (len=MaxLen) :: vinfo(natt)
476
477 character (len=*), parameter :: myfile = &
478 & __FILE__//", def_gst_pio"
479!
480 TYPE (var_desc_t) :: vardesc
481!
482 sourcefile=myfile
483!
484!-----------------------------------------------------------------------
485! Set and report file name.
486!-----------------------------------------------------------------------
487!
488 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
489 ncname=gst(ng)%name
490!
491 IF (master) THEN
492 IF (.not.lrstgst) THEN
493 WRITE (stdout,10) ng, trim(ncname)
494 ELSE
495 WRITE (stdout,20) ng, trim(ncname)
496 END IF
497 END IF
498!
499!=======================================================================
500! Create a new averages NetCDF file.
501!=======================================================================
502!
503 define : IF (.not.lrstgst) THEN
504 CALL pio_netcdf_create (ng, model, trim(ncname), &
505 & gst(ng)%pioFile)
506 IF (founderror(exit_flag, noerror, __line__, myfile)) THEN
507 IF (master) WRITE (stdout,30) trim(ncname)
508 RETURN
509 END IF
510!
511!-----------------------------------------------------------------------
512! Define the dimensions of staggered fields.
513!-----------------------------------------------------------------------
514!
515 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'Mstate', &
516 & mstate(ng), mstatedim)
517 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
518
519 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'NCV', &
520 & ncv, ncvdim)
521 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
522
523 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'LworkD', &
524 & 3*mstate(ng), sworkddim)
525 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
526
527 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'LworkL', &
528 & lworkl, sworkldim)
529 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
530
531 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'iparam', &
532 & SIZE(iparam), iparamdim)
533 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
534
535 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'ipntr', &
536 & SIZE(ipntr), ipntrdim)
537 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
538
539 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'iaupd', &
540 & SIZE(iaupd), iaupddim)
541 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
542
543 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'laitr', &
544 & SIZE(laitr), laitrdim)
545 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
546
547 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'iaitr', &
548 & SIZE(iaitr), iaitrdim)
549 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
550
551 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'raitr', &
552 & SIZE(raitr), raitrdim)
553 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
554
555 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'laup2', &
556 & SIZE(laup2), laup2dim)
557 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
558
559 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'iaup2', &
560 & SIZE(iaup2), iaup2dim)
561 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
562
563 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'raup2', &
564 & SIZE(raup2), raup2dim)
565 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
566
567 status=def_dim(ng, model, gst(ng)%pioFile, ncname, 'char2', &
568 & 2, char2dim)
569 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
570!
571!-----------------------------------------------------------------------
572! Define global attributes.
573!-----------------------------------------------------------------------
574!
575 CALL def_info (ng, model, gst(ng)%pioFile, ncname, dimids)
576 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
577!
578!-----------------------------------------------------------------------
579! Define variables and their attributes.
580!-----------------------------------------------------------------------
581!
582! Initialize local information variable arrays.
583!
584 DO i=1,natt
585 DO j=1,len(vinfo(1))
586 vinfo(i)(j:j)=' '
587 END DO
588 END DO
589 DO i=1,6
590 aval(i)=0.0_r8
591 END DO
592!
593! Define number of eigenvalues to compute.
594!
595 vinfo( 1)='NEV'
596 vinfo( 2)='number of eigenvalues to compute'
597 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
598 & 1, (/0/), aval, vinfo, ncname, &
599 & setparaccess = .false.)
600 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
601!
602! Define number of Lanczos vectors to compute.
603!
604 vinfo( 1)='NCV'
605 vinfo( 2)='number of Lanczos vectors to compute'
606 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
607 & 1, (/0/), aval, vinfo, ncname, &
608 & setparaccess = .false.)
609 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
610!
611! Define size of eigenvalue problem.
612!
613 vinfo( 1)='Mstate'
614 vinfo( 2)='total size of eigenvalue problem'
615 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
616 & 1, (/0/), aval, vinfo, ncname, &
617 & setparaccess = .false.)
618 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
619!
620! Define iteration number.
621!
622 vinfo( 1)='iter'
623 vinfo( 2)='iteration number'
624 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
625 & 1, (/0/), aval, vinfo, ncname, &
626 & setparaccess = .false.)
627 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
628!
629! Define reverse communications flag.
630!
631 vinfo( 1)='ido'
632 vinfo( 2)='reverse communications flag'
633 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
634 & 1, (/0/), aval, vinfo, ncname, &
635 & setparaccess = .false.)
636 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
637!
638! Define information and error flag.
639!
640 vinfo( 1)='info'
641 vinfo( 2)='information and error flag'
642 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
643 & 1, (/0/), aval, vinfo, ncname, &
644 & setparaccess = .false.)
645 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
646!
647! Define eigenvalue problem type.
648!
649 vinfo( 1)='bmat'
650 vinfo( 2)='eigenvalue problem type'
651 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_char, &
652 & 1, (/0/), aval, vinfo, ncname, &
653 & setparaccess = .false.)
654 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
655!
656! Define Ritz eigenvalues to compute.
657!
658 vinfo( 1)='which'
659 vinfo( 2)='Ritz eigenvalues to compute'
660 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_char, &
661 & 1, (/char2dim/), aval, vinfo, ncname, &
662 & setparaccess = .false.)
663 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
664!
665! Define form of basis function.
666!
667 vinfo( 1)='howmany'
668 vinfo( 2)='form of basis function'
669 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_char, &
670 & 1, (/0/), aval, vinfo, ncname, &
671 & setparaccess = .false.)
672 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
673!
674! Define relative accuracy of computed Ritz values.
675!
676 vinfo( 1)='Ritz_tol'
677 vinfo( 2)='relative accuracy of computed Ritz values'
678 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_frst, &
679 & 1, (/0/), aval, vinfo, ncname, &
680 & setparaccess = .false.)
681 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
682!
683! Define eigenproblem parameters.
684!
685 vinfo( 1)='iparam'
686 vinfo( 2)='eigenproblem parameters'
687 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
688 & 1, (/iparamdim/), aval, vinfo, ncname, &
689 & setparaccess = .false.)
690 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
691!
692! Define pointers to mark starting location in work arrays.
693!
694 vinfo( 1)='ipntr'
695 vinfo( 2)='pointers to mark starting location in work arrays'
696 status=def_var(ng, model, gst(ng)%pioFile,vardesc, pio_int, &
697 & 1, (/ipntrdim/), aval, vinfo, ncname, &
698 & setparaccess = .false.)
699 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
700!
701! Define ARPACK internal integer parameters to _aupd routines.
702!
703 vinfo( 1)='iaupd'
704 vinfo( 2)='ARPACK internal integer parameters to _aupd routines'
705 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
706 & 1, (/iaupddim/), aval, vinfo, ncname, &
707 & setparaccess = .false.)
708 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
709!
710! Define ARPACK internal integer parameters to _aitr routines.
711!
712 vinfo( 1)='iaitr'
713 vinfo( 2)='ARPACK internal integer parameters to _aitr routines'
714 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
715 & 1, (/iaitrdim/), aval, vinfo, ncname, &
716 & setparaccess = .false.)
717 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
718!
719! Define ARPACK internal integer parameters to _aup2 routines.
720!
721 vinfo( 1)='iaup2'
722 vinfo( 2)='ARPACK internal integer parameters to _aup2 routines'
723 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
724 & 1, (/iaup2dim/), aval, vinfo, ncname, &
725 & setparaccess = .false.)
726 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
727!
728! Define ARPACK internal logical parameters to _aitr routines.
729!
730 vinfo( 1)='laitr'
731 vinfo( 2)='ARPACK internal logical parameters to _aitr routines'
732 vinfo( 9)='.FALSE.'
733 vinfo(10)='.TRUE.'
734 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
735 & 1, (/laitrdim/), aval, vinfo, ncname, &
736 & setparaccess = .false.)
737 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
738!
739! Define ARPACK internal logical parameters to _aup2 routines.
740!
741 vinfo( 1)='laup2'
742 vinfo( 2)='ARPACK internal logical parameters to _aup2 routines'
743 vinfo( 9)='.FALSE.'
744 vinfo(10)='.TRUE.'
745 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_int, &
746 & 1, (/laup2dim/), aval, vinfo, ncname, &
747 & setparaccess = .false.)
748 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
749!
750! Define ARPACK internal real parameters to _aitr routines.
751!
752 vinfo( 1)='raitr'
753 vinfo( 2)='ARPACK internal real parameters to _aitr routines'
754 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_frst, &
755 & 1, (/raitrdim/), aval, vinfo, ncname, &
756 & setparaccess = .false.)
757 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
758!
759! Define ARPACK internal real parameters to _aup2 routines.
760!
761 vinfo( 1)='raup2'
762 vinfo( 2)='ARPACK internal real parameters to _aup2 routines'
763 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_frst, &
764 & 1, (/raup2dim/), aval, vinfo, ncname, &
765 & setparaccess = .false.)
766 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
767!
768! Define Lanczos/Arnoldi basis vectors.
769!
770 vinfo( 1)='Bvec'
771 vinfo( 2)='Lanczos/Arnoldi basis vectors'
772 vardim(1)=mstatedim
773 vardim(2)=ncvdim
774 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_frst, &
775 & 2, vardim, aval, vinfo, ncname)
776 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
777!
778! Define eigenproblem residual vector.
779!
780 vinfo( 1)='resid'
781 vinfo( 2)='eigenproblem residual vector'
782 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_frst, &
783 & 1, (/mstatedim/), aval, vinfo, ncname)
784 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
785!
786! Define state reverse communications work array.
787!
788 vinfo( 1)='SworkD'
789 vinfo( 2)='reverse communications state array'
790 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_frst, &
791 & 1, (/sworkddim/), aval, vinfo, ncname)
792 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
793!
794! Define eigenproblem work array.
795!
796 vinfo( 1)='SworkL'
797 vinfo( 2)='eigenproblem work array'
798 status=def_var(ng, model, gst(ng)%pioFile, vardesc, pio_frst, &
799 & 1, (/sworkldim/), aval, vinfo, ncname)
800 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
801!
802!-----------------------------------------------------------------------
803! Leave definition mode.
804!-----------------------------------------------------------------------
805!
806 CALL pio_netcdf_enddef (ng, model, ncname, gst(ng)%pioFile)
807 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
808
809 END IF define
810!
811 10 FORMAT (/,2x,'DEF_GST_PIO - creating checkpointing file,', &
812 & t56,'Grid ',i2.2,': ',a)
813 20 FORMAT (/,2x,'DEF_GST_PIO - inquiring checkpointing file,', &
814 & t56,'Grid ',i2.2,': ',a)
815 30 FORMAT (/,' DEF_GST_PIO - unable to create checkpointing', &
816 & ' NetCDF file: ',a)
817!
818 RETURN
819 END SUBROUTINE def_gst_pio
820# endif
821#endif
822 END MODULE def_gst_mod
subroutine, private def_gst_nf90(ng, model)
Definition def_gst.F:84
subroutine, private def_gst_pio(ng, model)
Definition def_gst.F:453
subroutine, public def_gst(ng, model)
Definition def_gst.F:45
type(t_io), dimension(:), allocatable gst
integer stdout
character(len=256) sourcefile
integer, parameter io_nf90
Definition mod_ncparam.F:95
integer, parameter io_pio
Definition mod_ncparam.F:96
integer, parameter ndimid
subroutine, public netcdf_enddef(ng, model, ncname, ncid)
subroutine, public netcdf_create(ng, model, ncname, ncid)
integer, parameter nf_frst
Definition mod_netcdf.F:193
logical master
integer, dimension(:), allocatable mstate
Definition mod_param.F:644
integer, parameter pio_frst
subroutine, public pio_netcdf_create(ng, model, ncname, piofile)
subroutine, public pio_netcdf_enddef(ng, model, ncname, piofile)
logical lrstgst
integer exit_flag
integer noerror
real(r8), dimension(2) raup2
integer, dimension(8) iaup2
integer, dimension(20) iaupd
integer, dimension(:,:), allocatable ipntr
integer ncv
real(r8), dimension(8) raitr
logical, dimension(5) laup2
integer lworkl
integer, dimension(:,:), allocatable iparam
integer, dimension(8) iaitr
logical, dimension(5) laitr
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52