ROMS
Loading...
Searching...
No Matches
tadv_mod Module Reference

Data Types

interface  tadv_putatt
 

Functions/Subroutines

subroutine tadv_putatt_nf90 (ng, ncid, ncname, aname, hadv, vadv, status)
 
subroutine tadv_putatt_pio (ng, piofile, ncname, aname, hadv, vadv, status)
 
subroutine tadv_report (iunit, model, hadv, vadv, lwrite)
 

Function/Subroutine Documentation

◆ tadv_putatt_nf90()

subroutine tadv_mod::tadv_putatt_nf90 ( integer, intent(in) ng,
integer, intent(in) ncid,
character (*), intent(in) ncname,
character (*), intent(in) aname,
type(t_adv), dimension(maxval(nt),ngrids), intent(in) hadv,
type(t_adv), dimension(maxval(nt),ngrids), intent(in) vadv,
integer, intent(out) status )

Definition at line 35 of file tadv.F.

37!***********************************************************************
38! !
39! This routine writes tracer advection scheme keywords strings into !
40! specified NetCDF file global attribute when using the standard !
41! NetCDF-3 or NetCDF-4 library. !
42! !
43! On Input: !
44! !
45! ng Nested grid number (integer) !
46! model Calling model identifier (integer) !
47! ncid NetCDF file ID (integer) !
48! ncname NetCDF filename (character) !
49! aname NetCDF global attribute name (character) !
50! Hadv Horizontal advection type structure, TYPE(T_ADV) !
51! Vadv Vertical advection type structure, TYPE(T_ADV) !
52! !
53! On Output: !
54! !
55! exit_flag Error flag (integer) stored in MOD_SCALARS !
56! ioerror NetCDF return code (integer) stored in MOD_IOUNITS !
57! status NetCDF return code (integer) !
58! !
59!***********************************************************************
60!
61 USE mod_param
62 USE mod_parallel
63 USE mod_iounits
64 USE mod_ncparam
65 USE mod_netcdf
66 USE mod_scalars
67!
68 USE strings_mod, ONLY : founderror
69!
70 implicit none
71!
72! Imported variable declarations.
73!
74 integer, intent(in) :: ng, ncid
75 integer, intent(out) :: status
76!
77 character (*), intent(in) :: ncname
78 character (*), intent(in) :: aname
79!
80 TYPE(T_ADV), intent(in) :: Hadv(MAXVAL(NT),Ngrids)
81 TYPE(T_ADV), intent(in) :: Vadv(MAXVAL(NT),Ngrids)
82!
83! Local variable declarations
84!
85 integer :: i, ie, is, itrc, lvar, lstr, nTvar
86!
87 character (len= 1) :: newline
88 character (len= 13) :: Hstring, Vstring
89 character (len= 17) :: frmt
90 character (len= 70) :: line
91 character (len=2816) :: tadv_att
92
93 character (len=*), parameter :: MyFile = &
94 & __FILE__//", tadv_putatt_nf90"
95!
96!-----------------------------------------------------------------------
97! Write lateral boundary conditions global attribute.
98!-----------------------------------------------------------------------
99!
100! Determine maximum length of state variable length.
101!
102 ntvar=maxval(nt)
103 lvar=0
104 DO itrc=1,ntvar
105 lvar=max(lvar, len_trim(vname(1,idtvar(itrc))))
106 END DO
107 WRITE (frmt,10) max(10,lvar)+4
108 10 FORMAT ("(a,':',t",i2.2,',a,a,a)')
109!
110! Initialize attribute.
111!
112 newline=char(10) ! Line Feed (LF) character for
113 lstr=len_trim(newline) ! attribute clarity with "ncdump"
114 DO i=1,len(tadv_att)
115 tadv_att(i:i)=' '
116 END DO
117 tadv_att(1:lstr)=newline(1:lstr)
118 is=lstr+1
119 WRITE (line,frmt) 'ADVECTION', &
120 & 'HORIZONTAL ', &
121 & 'VERTICAL ', &
122 & newline(1:lstr)
123 lstr=len_trim(line)
124 ie=is+lstr
125 tadv_att(is:ie)=line(1:lstr)
126 is=ie
127!
128! Check if the local string "tadv_att" is big enough to store the
129! tracer advection scheme global attribute.
130!
131 lstr=(ntvar+1)*(26+lvar+4)+1
132 IF (len(tadv_att).lt.lstr) THEN
133 IF (master) THEN
134 WRITE (stdout,20) len(tadv_att), lstr
135 20 FORMAT (/,' TADV_PUTATT_NF90 - Length of local string ', &
136 & ' tadv_att too small',/,20x,'Current = ',i5, &
137 & ' Needed = ',i5)
138 END IF
139 exit_flag=5
140 RETURN
141 END IF
142!
143! Build attribute string.
144!
145 DO itrc=1,ntvar
146 IF (hadv(itrc,ng)%AKIMA4) THEN
147 hstring='Akima4'
148 ELSE IF (hadv(itrc,ng)%CENTERED2) THEN
149 hstring='Centered2'
150 ELSE IF (hadv(itrc,ng)%CENTERED4) THEN
151 hstring='Centered4'
152 ELSE IF (hadv(itrc,ng)%HSIMT) THEN
153 hstring='HSIMT'
154 ELSE IF (hadv(itrc,ng)%MPDATA) THEN
155 hstring='MPDATA'
156 ELSE IF (hadv(itrc,ng)%SPLINES) THEN
157 hstring='Splines'
158 ELSE IF (hadv(itrc,ng)%SPLIT_U3) THEN
159 hstring='Split_U3'
160 ELSE IF (hadv(itrc,ng)%UPSTREAM3) THEN
161 hstring='Upstream3'
162 END IF
163!
164 IF (vadv(itrc,ng)%AKIMA4) THEN
165 vstring='Akima4'
166 ELSE IF (vadv(itrc,ng)%CENTERED2) THEN
167 vstring='Centered2'
168 ELSE IF (vadv(itrc,ng)%CENTERED4) THEN
169 vstring='Centered4'
170 ELSE IF (vadv(itrc,ng)%HSIMT) THEN
171 vstring='HSIMT'
172 ELSE IF (vadv(itrc,ng)%MPDATA) THEN
173 vstring='MPDATA'
174 ELSE IF (vadv(itrc,ng)%SPLINES) THEN
175 vstring='Splines'
176 ELSE IF (vadv(itrc,ng)%SPLIT_U3) THEN
177 vstring='Split_U3'
178 ELSE IF (vadv(itrc,ng)%UPSTREAM3) THEN
179 vstring='Upstream3'
180 END IF
181 IF (itrc.eq.ntvar) newline=' '
182 WRITE (line,frmt) trim(vname(1,idtvar(itrc))), &
183 & hstring, vstring, &
184 & newline
185 lstr=len_trim(line)
186 ie=is+lstr
187 tadv_att(is:ie)=line(1:lstr)
188 is=ie
189 END DO
190!
191! Write attribute to NetCDF file.
192!
193 status=nf90_put_att(ncid, nf90_global, trim(aname), &
194 & trim(tadv_att))
195 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
196 exit_flag=3
197 ioerror=status
198 RETURN
199 END IF
200!
201 RETURN
integer ioerror
integer stdout
integer, dimension(:), allocatable idtvar
character(len=maxlen), dimension(6, 0:nv) vname
logical master
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer exit_flag
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52

◆ tadv_putatt_pio()

subroutine tadv_mod::tadv_putatt_pio ( integer, intent(in) ng,
type (file_desc_t), intent(in) piofile,
character (*), intent(in) ncname,
character (*), intent(in) aname,
type(t_adv), dimension(maxval(nt),ngrids), intent(in) hadv,
type(t_adv), dimension(maxval(nt),ngrids), intent(in) vadv,
integer, intent(out) status )

Definition at line 207 of file tadv.F.

209!***********************************************************************
210! !
211! This routine writes tracer advection scheme keywords strings into !
212! specified NetCDF file global attribute when using the NCAR !
213! Parallel-IO library. !
214! !
215! On Input: !
216! !
217! ng Nested grid number (integer) !
218! model Calling model identifier (integer) !
219! pioFile PIO file descriptor structure, TYPE(File_desc_t) !
220! pioFile%fh file handler !
221! pioFile%iosystem IO system descriptor (struct) !
222! ncname PIO filename (character) !
223! aname PIO global attribute name (character) !
224! Hadv Horizontal advection type structure, TYPE(T_ADV) !
225! Vadv Vertical advection type structure, TYPE(T_ADV) !
226! !
227! On Output: !
228! !
229! exit_flag Error flag (integer) stored in MOD_SCALARS !
230! ioerror NetCDF return code (integer) stored in MOD_IOUNITS !
231! status NetCDF return code (integer) !
232! !
233!***********************************************************************
234!
235 USE mod_param
236 USE mod_parallel
237 USE mod_iounits
238 USE mod_ncparam
239 USE mod_scalars
240 USE pio
241!
242 USE strings_mod, ONLY : founderror
243!
244 implicit none
245!
246! Imported variable declarations.
247!
248 integer, intent(in) :: ng
249 integer, intent(out) :: status
250!
251 character (*), intent(in) :: ncname
252 character (*), intent(in) :: aname
253!
254 TYPE (File_desc_t), intent(in) :: pioFile
255 TYPE(T_ADV), intent(in) :: Hadv(MAXVAL(NT),Ngrids)
256 TYPE(T_ADV), intent(in) :: Vadv(MAXVAL(NT),Ngrids)
257!
258! Local variable declarations
259!
260 integer :: i, ie, is, itrc, lvar, lstr, nTvar
261!
262 character (len= 1) :: newline
263 character (len= 13) :: Hstring, Vstring
264 character (len= 17) :: frmt
265 character (len= 70) :: line
266 character (len=2816) :: tadv_att
267
268 character (len=*), parameter :: MyFile = &
269 & __FILE__//", tadv_putatt_pio"
270!
271!-----------------------------------------------------------------------
272! Write lateral boundary conditions global attribute.
273!-----------------------------------------------------------------------
274!
275! Determine maximum length of state variable length.
276!
277 ntvar=maxval(nt)
278 lvar=0
279 DO itrc=1,ntvar
280 lvar=max(lvar, len_trim(vname(1,idtvar(itrc))))
281 END DO
282 WRITE (frmt,10) max(10,lvar)+4
283 10 FORMAT ("(a,':',t",i2.2,',a,a,a)')
284!
285! Initialize attribute.
286!
287 newline=char(10) ! Line Feed (LF) character for
288 lstr=len_trim(newline) ! attribute clarity with "ncdump"
289 DO i=1,len(tadv_att)
290 tadv_att(i:i)=' '
291 END DO
292 tadv_att(1:lstr)=newline(1:lstr)
293 is=lstr+1
294 WRITE (line,frmt) 'ADVECTION', &
295 & 'HORIZONTAL ', &
296 & 'VERTICAL ', &
297 & newline(1:lstr)
298 lstr=len_trim(line)
299 ie=is+lstr
300 tadv_att(is:ie)=line(1:lstr)
301 is=ie
302!
303! Check if the local string "tadv_att" is big enough to store the
304! tracer advection scheme global attribute.
305!
306 lstr=(ntvar+1)*(26+lvar+4)+1
307 IF (len(tadv_att).lt.lstr) THEN
308 IF (master) THEN
309 WRITE (stdout,20) len(tadv_att), lstr
310 20 FORMAT (/,' TADV_PUTATT_PIO - Length of local string ', &
311 & 'tadv_att too small',/,19x,'Current = ',i5, &
312 & ' Needed = ',i5)
313 END IF
314 exit_flag=5
315 RETURN
316 END IF
317!
318! Build attribute string.
319!
320 DO itrc=1,ntvar
321 IF (hadv(itrc,ng)%AKIMA4) THEN
322 hstring='Akima4'
323 ELSE IF (hadv(itrc,ng)%CENTERED2) THEN
324 hstring='Centered2'
325 ELSE IF (hadv(itrc,ng)%CENTERED4) THEN
326 hstring='Centered4'
327 ELSE IF (hadv(itrc,ng)%HSIMT) THEN
328 hstring='HSIMT'
329 ELSE IF (hadv(itrc,ng)%MPDATA) THEN
330 hstring='MPDATA'
331 ELSE IF (hadv(itrc,ng)%SPLINES) THEN
332 hstring='Splines'
333 ELSE IF (hadv(itrc,ng)%SPLIT_U3) THEN
334 hstring='Split_U3'
335 ELSE IF (hadv(itrc,ng)%UPSTREAM3) THEN
336 hstring='Upstream3'
337 END IF
338!
339 IF (vadv(itrc,ng)%AKIMA4) THEN
340 vstring='Akima4'
341 ELSE IF (vadv(itrc,ng)%CENTERED2) THEN
342 vstring='Centered2'
343 ELSE IF (vadv(itrc,ng)%CENTERED4) THEN
344 vstring='Centered4'
345 ELSE IF (vadv(itrc,ng)%HSIMT) THEN
346 vstring='HSIMT'
347 ELSE IF (vadv(itrc,ng)%MPDATA) THEN
348 vstring='MPDATA'
349 ELSE IF (vadv(itrc,ng)%SPLINES) THEN
350 vstring='Splines'
351 ELSE IF (vadv(itrc,ng)%SPLIT_U3) THEN
352 vstring='Split_U3'
353 ELSE IF (vadv(itrc,ng)%UPSTREAM3) THEN
354 vstring='Upstream3'
355 END IF
356 IF (itrc.eq.ntvar) newline=' '
357 WRITE (line,frmt) trim(vname(1,idtvar(itrc))), &
358 & hstring, vstring, &
359 & newline
360 lstr=len_trim(line)
361 ie=is+lstr
362 tadv_att(is:ie)=line(1:lstr)
363 is=ie
364 END DO
365!
366! Write attribute to NetCDF file.
367!
368 status=pio_put_att(piofile, pio_global, trim(aname), &
369 & trim(tadv_att))
370 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
371 exit_flag=3
372 ioerror=status
373 RETURN
374 END IF
375!
376 RETURN

◆ tadv_report()

subroutine tadv_mod::tadv_report ( integer, intent(in) iunit,
integer, intent(in) model,
type(t_adv), dimension(maxval(nt),ngrids), intent(in) hadv,
type(t_adv), dimension(maxval(nt),ngrids), intent(in) vadv,
logical, intent(in) lwrite )

Definition at line 382 of file tadv.F.

383!***********************************************************************
384! !
385! This routine reports horizontal and vertical advection scheme for !
386! each tracer variable. !
387! !
388! On Input: !
389! !
390! iunit Output logical unit (integer) !
391! model Calling model identifier (integer) !
392! Hadv Horizontal advection type structure, TYPE(T_ADV) !
393! Vadv Vertical advection type structure, TYPE(T_ADV) !
394! Lwrite Switch to report information to standard output !
395! unit or file (logical) !
396! !
397!***********************************************************************
398!
399 USE mod_param
400 USE mod_parallel
401 USE mod_ncparam
402 USE mod_scalars
403!
404 implicit none
405!
406! Imported variable declarations.
407!
408 logical, intent(in) :: Lwrite
409!
410 integer, intent(in) :: iunit, model
411!
412 TYPE(T_ADV), intent(in) :: Hadv(MAXVAL(NT),Ngrids)
413 TYPE(T_ADV), intent(in) :: Vadv(MAXVAL(NT),Ngrids)
414!
415! Local variable declarations
416!
417 integer :: i, itrc, ng
418!
419 character (len=11) :: Hstring(MAXVAL(NT),Ngrids)
420 character (len=11) :: Vstring(MAXVAL(NT),Ngrids)
421!
422!-----------------------------------------------------------------------
423! Tracer horizontal and vertical switches.
424!-----------------------------------------------------------------------
425!
426 DO ng=1,ngrids
427 DO itrc=1,nt(ng)
428 IF (hadv(itrc,ng)%AKIMA4) THEN
429 hstring(itrc,ng)='Akima4'
430 ELSE IF (hadv(itrc,ng)%CENTERED2) THEN
431 hstring(itrc,ng)='Centered2'
432 ELSE IF (hadv(itrc,ng)%CENTERED4) THEN
433 hstring(itrc,ng)='Centered4'
434 ELSE IF (hadv(itrc,ng)%HSIMT) THEN
435 hstring(itrc,ng)='HSIMT'
436 ELSE IF (hadv(itrc,ng)%MPDATA) THEN
437 hstring(itrc,ng)='MPDATA'
438 ELSE IF (hadv(itrc,ng)%SPLINES) THEN
439 hstring(itrc,ng)='Splines'
440 ELSE IF (hadv(itrc,ng)%SPLIT_U3) THEN
441 hstring(itrc,ng)='Split_U3'
442 ELSE IF (hadv(itrc,ng)%UPSTREAM3) THEN
443 hstring(itrc,ng)='Upstream3'
444 END IF
445!
446 IF (vadv(itrc,ng)%AKIMA4) THEN
447 vstring(itrc,ng)='Akima4'
448 ELSE IF (vadv(itrc,ng)%CENTERED2) THEN
449 vstring(itrc,ng)='Centered2'
450 ELSE IF (vadv(itrc,ng)%CENTERED4) THEN
451 vstring(itrc,ng)='Centered4'
452 ELSE IF (vadv(itrc,ng)%HSIMT) THEN
453 vstring(itrc,ng)='HSIMT'
454 ELSE IF (vadv(itrc,ng)%MPDATA) THEN
455 vstring(itrc,ng)='MPDATA'
456 ELSE IF (vadv(itrc,ng)%SPLINES) THEN
457 vstring(itrc,ng)='Splines'
458 ELSE IF (vadv(itrc,ng)%SPLIT_U3) THEN
459 vstring(itrc,ng)='Split_U3'
460 ELSE IF (vadv(itrc,ng)%UPSTREAM3) THEN
461 vstring(itrc,ng)='Upstream3'
462 END IF
463 END DO
464 END DO
465!
466! Report.
467!
468 IF (master.and.lwrite) THEN
469 DO itrc=1,maxval(nt)
470 DO ng=1,ngrids
471 IF (ng.eq.1) THEN
472 WRITE (iunit,10) trim(vname(1,idtvar(itrc))), ng, &
473 & trim(hstring(itrc,ng)), &
474 & trim(vstring(itrc,ng))
475 ELSE
476 WRITE (iunit,20) ng, &
477 & trim(hstring(itrc,ng)), &
478 & trim(vstring(itrc,ng))
479 END IF
480 END DO
481 END DO
482!
483 IF (model.eq.inlm) THEN
484 WRITE (iunit,'(1x)')
485 WRITE (iunit,30) 'Akima4', &
486 & 'Fourth-order Akima advection'
487 WRITE (iunit,30) 'Centered2', &
488 & 'Second-order centered differences advection'
489 WRITE (iunit,30) 'Centered4', &
490 & 'Fourth-order centered differences advection'
491 WRITE (iunit,30) 'HSIMT', &
492 & 'Third High-order Spatial Inteporlation at Middle '// &
493 & 'Time Advection with TVD limiter'
494 WRITE (iunit,30) 'MPDATA', &
495 & 'Multidimensional Positive Definite Advection '// &
496 & 'Algorithm, recursive method'
497 WRITE (iunit,30) 'Splines', &
498 & 'Conservative Parabolic Splines Reconstruction '// &
499 & 'Advection (only vertical; not recommended)'
500 WRITE (iunit,30) 'Split_U3', &
501 & 'Split third-order Upstream Advection'
502 WRITE (iunit,30) 'Upstream3', &
503 & 'Third-order Upstream-biased Advection '// &
504 & '(only horizontal)'
505 WRITE (iunit,'(1x)')
506 END IF
507 END IF
508!
509! Check switches for consistency.
510!
511 IF (model.eq.inlm) THEN
512 DO ng=1,ngrids
513 DO i=1,nt(ng)
514 IF (.not.vadv(i,ng)%MPDATA.and.hadv(i,ng)%MPDATA) THEN
515 IF (master) THEN
516 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
517 & trim(hstring(i,ng)), &
518 & 'must be specified for both advective terms'
519 END IF
520 exit_flag=5
521 RETURN
522 ELSE IF (.not.hadv(i,ng)%MPDATA.and.vadv(i,ng)%MPDATA) THEN
523 IF (master) THEN
524 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
525 & trim(vstring(i,ng)), &
526 & 'must be specified for both advective terms'
527 END IF
528 exit_flag=5
529 RETURN
530 ELSE IF (.not.vadv(i,ng)%HSIMT.and.hadv(i,ng)%HSIMT) THEN
531 IF (master) THEN
532 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
533 & trim(hstring(i,ng)), &
534 & 'must be specified for both advective terms'
535 END IF
536 exit_flag=5
537 RETURN
538 ELSE IF (.not.hadv(i,ng)%HSIMT.and.vadv(i,ng)%HSIMT) THEN
539 IF (master) THEN
540 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
541 & trim(vstring(i,ng)), &
542 & 'must be specified for both advective terms'
543 END IF
544 exit_flag=5
545 RETURN
546 ELSE IF (hadv(i,ng)%SPLINES) THEN
547 IF (master) THEN
548 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
549 & trim(hstring(i,ng)), &
550 & 'is only available for the vertical term'
551
552 END IF
553 exit_flag=5
554 RETURN
555 ELSE IF (vadv(i,ng)%UPSTREAM3) THEN
556 IF (master) THEN
557 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
558 & trim(vstring(i,ng)), &
559 & 'is only available for the horizontal term'
560
561 END IF
562 exit_flag=5
563 RETURN
564 END IF
565 END DO
566 END DO
567 ELSE IF (model.eq.iadm) THEN
568 DO ng=1,ngrids
569 DO i=1,nt(ng)
570 IF (hadv(i,ng)%MPDATA.or.hadv(i,ng)%HSIMT) THEN
571 IF (master) THEN
572 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
573 & trim(hstring(i,ng)), &
574 & 'is not supported in adjoint-based algorithms'
575 END IF
576 exit_flag=5
577 RETURN
578 ELSE IF (vadv(i,ng)%MPDATA.or.vadv(i,ng)%HSIMT) THEN
579 IF (master) THEN
580 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
581 & trim(vstring(i,ng)), &
582 & 'is not supported in adjoint-based algorithms'
583 END IF
584 exit_flag=5
585 RETURN
586 ELSE IF (hadv(i,ng)%SPLINES) THEN
587 IF (master) THEN
588 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
589 & trim(hstring(i,ng)), &
590 & 'is only available for the vertical term'
591
592 END IF
593 exit_flag=5
594 RETURN
595 ELSE IF (vadv(i,ng)%UPSTREAM3) THEN
596 IF (master) THEN
597 WRITE (iunit,40) trim(vname(1,idtvar(i))), ng, &
598 & trim(vstring(i,ng)), &
599 & 'is only available for the horizontal term'
600
601 END IF
602 exit_flag=5
603 RETURN
604 END IF
605 END DO
606 END DO
607 END IF
608!
609 10 FORMAT (/,1x,a,t26,i2,t31,a,t50,a)
610 20 FORMAT (t26,i2,t31,a,t50,a)
611 30 FORMAT (1x,a,t13,a)
612 40 FORMAT (/,'TADV_REPORT - Illegal tracer advection scheme for ''', &
613 & a,''' in grid: ',i0,/,14x,'''',a,'''',1x,a,'.',/)
614!
615 RETURN
integer, parameter inlm
Definition mod_param.F:662
integer, parameter iadm
Definition mod_param.F:665
integer ngrids
Definition mod_param.F:113

References mod_scalars::exit_flag, mod_param::iadm, mod_ncparam::idtvar, mod_param::inlm, mod_parallel::master, and mod_ncparam::vname.

Referenced by inp_par_mod::inp_par().

Here is the caller graph for this function: