ROMS
Loading...
Searching...
No Matches
lbc.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE lbc_mod
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! These routines are used to process lateral boundary condition !
12! structure: !
13! !
14! lbc_getatt Reads activated keyword strings from specified !
15! input NetCDF file and checks input script !
16! values during restart for consistency. !
17! !
18! lbc_putatt Writes activated keyword strings into specified !
19! output NetCDF file global attribute. !
20! !
21! lbc_report Reports to standard output activated keyword !
22! strings. !
23! !
24!=======================================================================
25!
26 implicit none
27!
28 INTERFACE lbc_getatt
29 MODULE PROCEDURE lbc_getatt_nf90
30#if defined PIO_LIB && defined DISTRIBUTE
31 MODULE PROCEDURE lbc_getatt_pio
32#endif
33 END INTERFACE lbc_getatt
34!
35 INTERFACE lbc_putatt
36 MODULE PROCEDURE lbc_putatt_nf90
37#if defined PIO_LIB && defined DISTRIBUTE
38 MODULE PROCEDURE lbc_putatt_pio
39#endif
40 END INTERFACE lbc_putatt
41!
42 CONTAINS
43!
44!***********************************************************************
45 SUBROUTINE lbc_getatt_nf90 (ng, model, ncid, ncname, aname, S)
46!***********************************************************************
47! !
48! If restart, this routine reads lateral boundary conditions !
49! keywords strings from restart NetCDF file global attribute. !
50! Then, it checks their values with those provided from input !
51! script for consistency. !
52! !
53! On Input: !
54! !
55! ng Nested grid number (integer). !
56! model Calling model identifier (integer). !
57! ncid NetCDF file ID (integer). !
58! ncname NetCDF file name (character). !
59! aname NetCDF global attribute name (character). !
60! S Derived type structure, TYPE(T_LBC) !
61! !
62! On Output: !
63! !
64! exit_flag Error flag (integer) stored in MOD_SCALARS !
65! ioerror NetCDF return code (integer) stored in MOD_IOUNITS !
66! !
67!***********************************************************************
68!
69 USE mod_param
70 USE mod_parallel
71 USE mod_iounits
72 USE mod_ncparam
73 USE mod_netcdf
74 USE mod_scalars
75!
76#if !defined PARALLEL_IO && defined DISTRIBUTE
78#endif
79 USE strings_mod, ONLY : founderror
80!
81 implicit none
82!
83! Imported variable declarations.
84!
85 integer, intent(in) :: ng, model, ncid
86 character (*), intent(in) :: ncname
87 character (*), intent(in) :: aname
88!
89 TYPE(t_lbc), intent(in) :: S(4,nLBCvar,Ngrids)
90!
91! Local variable declarations
92!
93 integer :: i, ibry, ie, ifield, is, ne, lstr, lvar, status
94#if !defined PARALLEL_IO && defined DISTRIBUTE
95 integer, dimension(2) :: ibuffer
96#endif
97!
98 character (len= 7) :: string(4)
99 character (len= 8) :: B(4)
100 character (len= 40) :: BryVar1, BryVar2
101 character (len= 70) :: Bstring, line
102 character (len=2816) :: lbc_att
103
104 character (len=*), parameter :: MyFile = &
105 & __FILE__//", lbc_getatt_nf90"
106!
107!-----------------------------------------------------------------------
108! If restart, read and check lateral boundary conditions global
109! attribute.
110!-----------------------------------------------------------------------
111!
112! Read in global attribute.
113!
114 IF (inpthread) THEN
115 status=nf90_get_att(ncid, nf90_global, trim(aname), lbc_att)
116 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
117 WRITE (stdout,10) trim(aname), trim(ncname), &
118 & trim(sourcefile)
119 exit_flag=3
120 ioerror=status
121 END IF
122 END IF
123
124#if !defined PARALLEL_IO && defined DISTRIBUTE
125!
126! Broadcast error flags to all processors in the group.
127!
128 ibuffer(1)=exit_flag
129 ibuffer(2)=ioerror
130 CALL mp_bcasti (ng, model, ibuffer)
131 exit_flag=ibuffer(1)
132 ioerror=ibuffer(2)
133 IF (exit_flag.eq.noerror) THEN
134 CALL mp_bcasts (ng, model, lbc_att)
135 END IF
136#endif
137!
138! Check keyword values from global attribute and compare against
139! logical switches in lateral boundary condition structure.
140!
141 b(iwest )='western '
142 b(isouth)='southern'
143 b(ieast )='eastern '
144 b(inorth)='northern'
145 DO i=1,len(bstring)
146 bstring(i:i)=' '
147 END DO
148 DO i=1,len(line)
149 line(i:i)=' '
150 END DO
151!
152 DO ifield=1,nlbcvar
153 bryvar1=trim(vname(1,idbvar(ifield))) // ':'
154 is=index(lbc_att, trim(bryvar1))
155 IF (ifield.lt.nlbcvar) THEN
156 bryvar2=trim(vname(1,idbvar(ifield+1))) // ':'
157 ie=index(lbc_att, trim(bryvar2))-1
158 ELSE
159 ie=len_trim(lbc_att)
160 END IF
161 IF ((is.gt.0).and.(ie.gt.0).and.(ie.gt.is)) THEN
162 line=lbc_att(is:ie)
163 is=index(line, ':')+1 ! find colon location
164 ie=index(line, char(10))-1 ! find Line Feed location
165 IF (ie.le.0) THEN
166 ie=len_trim(line) ! last attribute entry
167 END IF
168 bstring=trim(adjustl(line(is:ie))) ! extract boundary string
169 ne=min(len_trim(bstring), 28) ! north keyword < 7 chars
170 string(1)=bstring( 1: 7) ! west keyword
171 string(2)=bstring( 8:14) ! south keyword
172 string(3)=bstring(15:21) ! east keyword
173 string(4)=bstring(22:ne) ! north keyword
174 DO ibry=1,4
175 SELECT CASE (trim(string(ibry)))
176 CASE ('Cha')
177 IF (.not.s(ibry,ifield,ng)%Chapman_implicit) THEN
178 IF (master) THEN
179 WRITE (stdout,20) b(ibry), &
180 & trim(vname(1,idbvar(ifield))), &
181 & trim(string(ibry)), &
182 & 'S(',ibry,ifield,ng,')%Chapman_implicit', &
183 & s(ibry,ifield,ng)%Chapman_implicit, &
184 & trim(ncname)
185 END IF
186 exit_flag=5
187 END IF
188 CASE ('Che')
189 IF (.not.s(ibry,ifield,ng)%Chapman_explicit) THEN
190 IF (master) THEN
191 WRITE (stdout,20) b(ibry), &
192 & trim(vname(1,idbvar(ifield))), &
193 & trim(string(ibry)), &
194 & 'S(',ibry,ifield,ng,')%Chapman_explicit', &
195 & s(ibry,ifield,ng)%Chapman_explicit, &
196 & trim(ncname)
197 END IF
198 exit_flag=5
199 END IF
200 CASE ('Cla')
201 IF (.not.s(ibry,ifield,ng)%clamped) THEN
202 IF (master) THEN
203 WRITE (stdout,20) b(ibry), &
204 & trim(vname(1,idbvar(ifield))), &
205 & trim(string(ibry)), &
206 & 'S(',ibry,ifield,ng,')%clamped', &
207 & s(ibry,ifield,ng)%clamped, &
208 & trim(ncname)
209 END IF
210 exit_flag=5
211 END IF
212 CASE ('Clo')
213 IF (.not.s(ibry,ifield,ng)%closed) THEN
214 IF (master) THEN
215 WRITE (stdout,20) b(ibry), &
216 & trim(vname(1,idbvar(ifield))), &
217 & trim(string(ibry)), &
218 & 'S(',ibry,ifield,ng,')%closed', &
219 & s(ibry,ifield,ng)%closed, &
220 & trim(ncname)
221 END IF
222 exit_flag=5
223 END IF
224 CASE ('Fla')
225 IF (.not.s(ibry,ifield,ng)%Flather) THEN
226 IF (master) THEN
227 WRITE (stdout,20) b(ibry), &
228 & trim(vname(1,idbvar(ifield))), &
229 & trim(string(ibry)), &
230 & 'S(',ibry,ifield,ng,')%Flather', &
231 & s(ibry,ifield,ng)%Flather, &
232 & trim(ncname)
233 END IF
234 exit_flag=5
235 END IF
236 CASE ('Gra')
237 IF (.not.s(ibry,ifield,ng)%gradient) THEN
238 IF (master) THEN
239 WRITE (stdout,20) b(ibry), &
240 & trim(vname(1,idbvar(ifield))), &
241 & trim(string(ibry)), &
242 & 'S(',ibry,ifield,ng,')%gradient', &
243 & s(ibry,ifield,ng)%gradient, &
244 & trim(ncname)
245 END IF
246 exit_flag=5
247 END IF
248 CASE ('Mix')
249 IF (.not.s(ibry,ifield,ng)%mixed) THEN
250 IF (master) THEN
251 WRITE (stdout,20) b(ibry), &
252 & trim(vname(1,idbvar(ifield))), &
253 & trim(string(ibry)), &
254 & 'S(',ibry,ifield,ng,')%mixed', &
255 & s(ibry,ifield,ng)%mixed, &
256 & trim(ncname)
257 END IF
258 exit_flag=5
259 END IF
260 CASE ('Nes')
261 IF (.not.s(ibry,ifield,ng)%nested) THEN
262 IF (master) THEN
263 WRITE (stdout,20) b(ibry), &
264 & trim(vname(1,idbvar(ifield))), &
265 & trim(string(ibry)), &
266 & 'S(',ibry,ifield,ng,')%nested', &
267 & s(ibry,ifield,ng)%nested, &
268 & trim(ncname)
269 END IF
270 exit_flag=5
271 END IF
272 CASE ('Per')
273 IF (.not.s(ibry,ifield,ng)%periodic) THEN
274 IF (master) THEN
275 WRITE (stdout,20) b(ibry), &
276 & trim(vname(1,idbvar(ifield))), &
277 & trim(string(ibry)), &
278 & 'S(',ibry,ifield,ng,')%periodic', &
279 & s(ibry,ifield,ng)%periodic, &
280 & trim(ncname)
281 END IF
282 exit_flag=5
283 END IF
284 CASE ('Rad')
285 IF (.not.s(ibry,ifield,ng)%radiation) THEN
286 IF (master) THEN
287 WRITE (stdout,20) b(ibry), &
288 & trim(vname(1,idbvar(ifield))), &
289 & trim(string(ibry)), &
290 & 'S(',ibry,ifield,ng,')%radiation', &
291 & s(ibry,ifield,ng)%radiation, &
292 & trim(ncname)
293 END IF
294 exit_flag=5
295 END IF
296 CASE ('RadNud')
297 IF (.not.(s(ibry,ifield,ng)%radiation.and. &
298 & s(ibry,ifield,ng)%nudging)) THEN
299 IF (master) THEN
300 WRITE (stdout,20) b(ibry), &
301 & trim(vname(1,idbvar(ifield))), &
302 & trim(string(ibry)), &
303 & 'S(',ibry,ifield,ng,')%radiation', &
304 & s(ibry,ifield,ng)%radiation, &
305 & trim(ncname)
306 END IF
307 exit_flag=5
308 END IF
309 CASE ('Red')
310 IF (.not.s(ibry,ifield,ng)%reduced) THEN
311 IF (master) THEN
312 WRITE (stdout,20) b(ibry), &
313 & trim(vname(1,idbvar(ifield))), &
314 & trim(string(ibry)), &
315 & 'S(',ibry,ifield,ng,')%reduced', &
316 & s(ibry,ifield,ng)%reduced, &
317 & trim(ncname)
318 END IF
319 exit_flag=5
320 END IF
321 CASE ('Shc')
322 IF (.not.s(ibry,ifield,ng)%Shchepetkin) THEN
323 IF (master) THEN
324 WRITE (stdout,20) b(ibry), &
325 & trim(vname(1,idbvar(ifield))), &
326 & trim(string(ibry)), &
327 & 'S(',ibry,ifield,ng,')%Shchepetkin', &
328 & s(ibry,ifield,ng)%Shchepetkin, &
329 & trim(ncname)
330 END IF
331 exit_flag=5
332 END IF
333 CASE DEFAULT
334 IF (master) THEN
335 WRITE (stdout,30) b(ibry), &
336 & trim(vname(1,idbvar(ifield))), &
337 & trim(string(ibry)), trim(ncname)
338 END IF
339 exit_flag=5
340 END SELECT
341 END DO
342 END IF
343 END DO
344!
345 10 FORMAT (/,' LBC_GETATT_NF90 - error while reading global ', &
346 & 'attribute:',2x,a,/,19x,'in restart file:',2x,a,/, &
347 & 19x,'call from:',2x,a, &
348 & /,19x,'Probably global attribute was not found ...', &
349 & /,19x,'restart file needs to be generated by ROMS ', &
350 & 'version 3.6 or higher', &
351 & /,19x,'Alternatively, you may use NO_LBC_ATT at your ', &
352 & 'own risk!')
353 20 FORMAT (/,' LBC_GETATT_NF90 - inconsistent ',a,' lateral', &
354 & 'boundary condition for variable: ',2x,a, &
355 & /,19x,'restart file LBC keyword = ',1x,a, &
356 & /,19x,'but assigned structure switch: ', &
357 & 1x,a,i1,',',i2,',',i1,a,' = ',l1, &
358 & /,19x,'check input script LBC keyword for consitency ...',&
359 & /,19x,'restart file:',2x,a)
360 30 FORMAT (/,' LBC_GETATT_NF90 - inconsistent ',a,' boundary for ', &
361 & 'variable: ',a,2x,'Keyword = ',a,/,19x,'in input file:', &
362 & 2x,a)
363!
364 RETURN
365 END SUBROUTINE lbc_getatt_nf90
366
367#if defined PIO_LIB && defined DISTRIBUTE
368!
369!***********************************************************************
370 SUBROUTINE lbc_getatt_pio (ng, model, pioFile, ncname, aname, S)
371!***********************************************************************
372! !
373! If restart, this routine reads lateral boundary conditions !
374! keywords strings from restart NetCDF file global attribute. !
375! Then, it checks their values with those provided from input !
376! script for consistency. !
377! !
378! On Input: !
379! !
380! ng Nested grid number (integer) !
381! model Calling model identifier (integer) !
382! pioFile PIO file descriptor structure, TYPE(File_desc_t) !
383! pioFile%fh file handler !
384! pioFile%iosystem IO system descriptor (struct) !
385! ncname PIO filename (character) !
386! aname PIO global attribute name (character) !
387! S Derived type structure, TYPE(T_LBC) !
388! !
389! On Output: !
390! !
391! exit_flag Error flag (integer) stored in MOD_SCALARS !
392! ioerror NetCDF return code (integer) stored in MOD_IOUNITS !
393! !
394!***********************************************************************
395!
396 USE mod_param
397 USE mod_parallel
398 USE mod_iounits
399 USE mod_ncparam
400 USE mod_scalars
401 USE pio
402!
403 USE strings_mod, ONLY : founderror
404!
405 implicit none
406!
407! Imported variable declarations.
408!
409 integer, intent(in) :: ng, model
410 character (*), intent(in) :: ncname
411 character (*), intent(in) :: aname
412!
413 TYPE (File_desc_t), intent(in) :: pioFile
414 TYPE(t_lbc), intent(in) :: S(4,nLBCvar,Ngrids)
415!
416! Local variable declarations
417!
418 integer :: i, ibry, ie, ifield, is, ne, lstr, lvar, status
419!
420 character (len= 7) :: string(4)
421 character (len= 8) :: B(4)
422 character (len= 40) :: BryVar1, BryVar2
423 character (len= 70) :: Bstring, line
424 character (len=2816) :: lbc_att
425
426 character (len=*), parameter :: MyFile = &
427 & __FILE__//", lbc_getatt_pio"
428!
429!-----------------------------------------------------------------------
430! If restart, read and check lateral boundary conditions global
431! attribute.
432!-----------------------------------------------------------------------
433!
434! Read in global attribute.
435!
436 status=pio_get_att(piofile, pio_global, trim(aname), lbc_att)
437 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
438 WRITE (stdout,10) trim(aname), trim(ncname), &
439 & trim(sourcefile)
440 exit_flag=3
441 ioerror=status
442 END IF
443!
444! Check keyword values from global attribute and compare against
445! logical switches in lateral boundary condition structure.
446!
447 b(iwest )='western '
448 b(isouth)='southern'
449 b(ieast )='eastern '
450 b(inorth)='northern'
451 DO i=1,len(bstring)
452 bstring(i:i)=' '
453 END DO
454 DO i=1,len(line)
455 line(i:i)=' '
456 END DO
457!
458 DO ifield=1,nlbcvar
459 bryvar1=trim(vname(1,idbvar(ifield))) // ':'
460 is=index(lbc_att, trim(bryvar1))
461 IF (ifield.lt.nlbcvar) THEN
462 bryvar2=trim(vname(1,idbvar(ifield+1))) // ':'
463 ie=index(lbc_att, trim(bryvar2))-1
464 ELSE
465 ie=len_trim(lbc_att)
466 END IF
467 IF ((is.gt.0).and.(ie.gt.0).and.(ie.gt.is)) THEN
468 line=lbc_att(is:ie)
469 is=index(line, ':')+1 ! find colon location
470 ie=index(line, char(10))-1 ! find Line Feed location
471 IF (ie.le.0) THEN
472 ie=len_trim(line) ! last attribute entry
473 END IF
474 bstring=trim(adjustl(line(is:ie))) ! extract boundary string
475 ne=min(len_trim(bstring), 28) ! north keyword < 7 chars
476 string(1)=bstring( 1: 7) ! west keyword
477 string(2)=bstring( 8:14) ! south keyword
478 string(3)=bstring(15:21) ! east keyword
479 string(4)=bstring(22:ne) ! north keyword
480 DO ibry=1,4
481 SELECT CASE (trim(string(ibry)))
482 CASE ('Cha')
483 IF (.not.s(ibry,ifield,ng)%Chapman_implicit) THEN
484 IF (master) THEN
485 WRITE (stdout,20) b(ibry), &
486 & trim(vname(1,idbvar(ifield))), &
487 & trim(string(ibry)), &
488 & 'S(',ibry,ifield,ng,')%Chapman_implicit', &
489 & s(ibry,ifield,ng)%Chapman_implicit, &
490 & trim(ncname)
491 END IF
492 exit_flag=5
493 END IF
494 CASE ('Che')
495 IF (.not.s(ibry,ifield,ng)%Chapman_explicit) THEN
496 IF (master) THEN
497 WRITE (stdout,20) b(ibry), &
498 & trim(vname(1,idbvar(ifield))), &
499 & trim(string(ibry)), &
500 & 'S(',ibry,ifield,ng,')%Chapman_explicit', &
501 & s(ibry,ifield,ng)%Chapman_explicit, &
502 & trim(ncname)
503 END IF
504 exit_flag=5
505 END IF
506
507 CASE ('Cla')
508 IF (.not.s(ibry,ifield,ng)%clamped) THEN
509 IF (master) THEN
510 WRITE (stdout,20) b(ibry), &
511 & trim(vname(1,idbvar(ifield))), &
512 & trim(string(ibry)), &
513 & 'S(',ibry,ifield,ng,')%clamped', &
514 & s(ibry,ifield,ng)%clamped, &
515 & trim(ncname)
516 END IF
517 exit_flag=5
518 END IF
519 CASE ('Clo')
520 IF (.not.s(ibry,ifield,ng)%closed) THEN
521 IF (master) THEN
522 WRITE (stdout,20) b(ibry), &
523 & trim(vname(1,idbvar(ifield))), &
524 & trim(string(ibry)), &
525 & 'S(',ibry,ifield,ng,')%closed', &
526 & s(ibry,ifield,ng)%closed, &
527 & trim(ncname)
528 END IF
529 exit_flag=5
530 END IF
531 CASE ('Fla')
532 IF (.not.s(ibry,ifield,ng)%Flather) THEN
533 IF (master) THEN
534 WRITE (stdout,20) b(ibry), &
535 & trim(vname(1,idbvar(ifield))), &
536 & trim(string(ibry)), &
537 & 'S(',ibry,ifield,ng,')%Flather', &
538 & s(ibry,ifield,ng)%Flather, &
539 & trim(ncname)
540 END IF
541 exit_flag=5
542 END IF
543 CASE ('Gra')
544 IF (.not.s(ibry,ifield,ng)%gradient) THEN
545 IF (master) THEN
546 WRITE (stdout,20) b(ibry), &
547 & trim(vname(1,idbvar(ifield))), &
548 & trim(string(ibry)), &
549 & 'S(',ibry,ifield,ng,')%gradient', &
550 & s(ibry,ifield,ng)%gradient, &
551 & trim(ncname)
552 END IF
553 exit_flag=5
554 END IF
555 CASE ('Mix')
556 IF (.not.s(ibry,ifield,ng)%mixed) THEN
557 IF (master) THEN
558 WRITE (stdout,20) b(ibry), &
559 & trim(vname(1,idbvar(ifield))), &
560 & trim(string(ibry)), &
561 & 'S(',ibry,ifield,ng,')%mixed', &
562 & s(ibry,ifield,ng)%mixed, &
563 & trim(ncname)
564 END IF
565 exit_flag=5
566 END IF
567 CASE ('Nes')
568 IF (.not.s(ibry,ifield,ng)%nested) THEN
569 IF (master) THEN
570 WRITE (stdout,20) b(ibry), &
571 & trim(vname(1,idbvar(ifield))), &
572 & trim(string(ibry)), &
573 & 'S(',ibry,ifield,ng,')%nested', &
574 & s(ibry,ifield,ng)%nested, &
575 & trim(ncname)
576 END IF
577 exit_flag=5
578 END IF
579 CASE ('Per')
580 IF (.not.s(ibry,ifield,ng)%periodic) THEN
581 IF (master) THEN
582 WRITE (stdout,20) b(ibry), &
583 & trim(vname(1,idbvar(ifield))), &
584 & trim(string(ibry)), &
585 & 'S(',ibry,ifield,ng,')%periodic', &
586 & s(ibry,ifield,ng)%periodic, &
587 & trim(ncname)
588 END IF
589 exit_flag=5
590 END IF
591 CASE ('Rad')
592 IF (.not.s(ibry,ifield,ng)%radiation) THEN
593 IF (master) THEN
594 WRITE (stdout,20) b(ibry), &
595 & trim(vname(1,idbvar(ifield))), &
596 & trim(string(ibry)), &
597 & 'S(',ibry,ifield,ng,')%radiation', &
598 & s(ibry,ifield,ng)%radiation, &
599 & trim(ncname)
600 END IF
601 exit_flag=5
602 END IF
603 CASE ('RadNud')
604 IF (.not.(s(ibry,ifield,ng)%radiation.and. &
605 & s(ibry,ifield,ng)%nudging)) THEN
606 IF (master) THEN
607 WRITE (stdout,20) b(ibry), &
608 & trim(vname(1,idbvar(ifield))), &
609 & trim(string(ibry)), &
610 & 'S(',ibry,ifield,ng,')%radiation', &
611 & s(ibry,ifield,ng)%radiation, &
612 & trim(ncname)
613 END IF
614 exit_flag=5
615 END IF
616 CASE ('Red')
617 IF (.not.s(ibry,ifield,ng)%reduced) THEN
618 IF (master) THEN
619 WRITE (stdout,20) b(ibry), &
620 & trim(vname(1,idbvar(ifield))), &
621 & trim(string(ibry)), &
622 & 'S(',ibry,ifield,ng,')%reduced', &
623 & s(ibry,ifield,ng)%reduced, &
624 & trim(ncname)
625 END IF
626 exit_flag=5
627 END IF
628 CASE ('Shc')
629 IF (.not.s(ibry,ifield,ng)%Shchepetkin) THEN
630 IF (master) THEN
631 WRITE (stdout,20) b(ibry), &
632 & trim(vname(1,idbvar(ifield))), &
633 & trim(string(ibry)), &
634 & 'S(',ibry,ifield,ng,')%Shchepetkin', &
635 & s(ibry,ifield,ng)%Shchepetkin, &
636 & trim(ncname)
637 END IF
638 exit_flag=5
639 END IF
640 CASE DEFAULT
641 IF (master) THEN
642 WRITE (stdout,30) b(ibry), &
643 & trim(vname(1,idbvar(ifield))), &
644 & trim(string(ibry)), trim(ncname)
645 END IF
646 exit_flag=5
647 END SELECT
648 END DO
649 END IF
650 END DO
651!
652 10 FORMAT (/,' LBC_GETATT_PIO - error while reading global ', &
653 & 'attribute:',2x,a,/,18x,'in restart file:',2x,a,/, &
654 & 18x,'call from:',2x,a, &
655 & /,18x,'Probably global attribute was not found ...', &
656 & /,18x,'restart file needs to be generated by ROMS ', &
657 & 'version 3.6 or higher', &
658 & /,18x,'Alternatively, you may use NO_LBC_ATT at your ', &
659 & 'own risk!')
660 20 FORMAT (/,' LBC_GETATT_PIO - inconsistent ',a,' lateral', &
661 & 'boundary condition for variable: ',2x,a, &
662 & /,18x,'restart file LBC keyword = ',1x,a, &
663 & /,18x,'but assigned structure switch: ', &
664 & 1x,a,i1,',',i2,',',i1,a,' = ',l1, &
665 & /,18x,'check input script LBC keyword for consitency ...',&
666 & /,18x,'restart file:',2x,a)
667 30 FORMAT (/,' LBC_GETATT_PIO - inconsistent ',a,' boundary for ', &
668 & 'variable: ',a,2x,'Keyword = ',a,/,18x,'in input file:', &
669 & 2x,a)
670!
671 RETURN
672 END SUBROUTINE lbc_getatt_pio
673#endif
674!
675!***********************************************************************
676 SUBROUTINE lbc_putatt_nf90 (ng, ncid, ncname, aname, S, status)
677!***********************************************************************
678! !
679! This routine writes lateral boundary conditions keywords strings !
680! into specified output NetCDF file global attribute. !
681! !
682! On Input: !
683! !
684! ng Nested grid number (integer) !
685! ncid NetCDF file ID (integer) !
686! ncname NetCDF filename (character) !
687! aname NetCDF global attribute name (character) !
688! S Derived type structure, TYPE(T_LBC) !
689! !
690! On Output: !
691! !
692! exit_flag Error flag (integer) stored in MOD_SCALARS !
693! ioerror NetCDF return code (integer) stored in MOD_IOUNITS !
694! status NetCDF return code (integer) !
695! !
696!***********************************************************************
697!
698 USE mod_param
699 USE mod_parallel
700 USE mod_iounits
701 USE mod_ncparam
702 USE mod_netcdf
703 USE mod_scalars
704!
705 USE strings_mod, ONLY : founderror
706!
707 implicit none
708!
709! Imported variable declarations.
710!
711 integer, intent(in) :: ng, ncid
712 integer, intent(out) :: status
713!
714 character (*), intent(in) :: ncname
715 character (*), intent(in) :: aname
716!
717 TYPE(t_lbc), intent(in) :: S(4,nLBCvar,Ngrids)
718!
719! Local variable declarations
720!
721 integer :: i, ibry, ie, ifield, is, lstr, lvar
722!
723 character (len= 1) :: newline
724 character (len= 7) :: string(4)
725 character (len= 21) :: frmt
726 character (len= 100) :: line
727 character (len=2816) :: lbc_att
728
729 character (len=*), parameter :: MyFile = &
730 & __FILE__//", lbc_putatt_nf90"
731!
732!-----------------------------------------------------------------------
733! Write lateral boundary conditions global attribute.
734!-----------------------------------------------------------------------
735!
736! Determine maximum length of state variable length.
737!
738 lvar=0
739 DO ifield=1,nlbcvar
740 IF (idbvar(ifield).gt.0) THEN
741 lvar=max(lvar, len_trim(vname(1,idbvar(ifield))))
742 END IF
743 END DO
744 WRITE (frmt,10) "(a,':',t", lvar+4, ",5a)"
745 10 FORMAT (a,i0,a)
746!
747! Initialize attribute.
748!
749 newline=char(10) ! Line Feed (LF) character for
750 lstr=len_trim(newline) ! attribute clarity with "ncdump"
751 DO i=1,len(lbc_att)
752 lbc_att(i:i)=' '
753 END DO
754 lbc_att(1:lstr)=newline(1:lstr)
755 is=lstr+1
756 WRITE (line,frmt) 'EDGE', &
757 & 'WEST ', &
758 & 'SOUTH ', &
759 & 'EAST ', &
760 & 'NORTH ', &
761 & newline(1:lstr)
762 lstr=len_trim(line)
763 ie=is+lstr
764 lbc_att(is:ie)=line(1:lstr)
765 is=ie
766!
767! Check if the local string "lbc_att" is big enough to store the
768! lateral boundary conditions global attribute.
769!
770 lstr=(nlbcvar+1)*(29+lvar+4)+1
771 IF (len(lbc_att).lt.lstr) THEN
772 IF (master) THEN
773 WRITE (stdout,20) len(lbc_att), lstr
774 20 FORMAT (/,' LBC_PUTATT_NF90 - Length of local string lbc_att',&
775 & ' too small',/,19x,'Current = ',i5,' Needed = ',i5)
776 END IF
777 exit_flag=5
778 RETURN
779 END IF
780!
781! Build attribute string.
782!
783 DO ifield=1,nlbcvar
784 IF (idbvar(ifield).gt.0) THEN
785 DO ibry=1,4
786 IF (s(ibry,ifield,ng)%Chapman_explicit) THEN
787 string(ibry)='Che '
788 ELSE IF (s(ibry,ifield,ng)%Chapman_implicit) THEN
789 string(ibry)='Cha '
790 ELSE IF (s(ibry,ifield,ng)%clamped) THEN
791 string(ibry)='Cla '
792 ELSE IF (s(ibry,ifield,ng)%closed) THEN
793 string(ibry)='Clo '
794 ELSE IF (s(ibry,ifield,ng)%Flather) THEN
795 string(ibry)='Fla '
796 ELSE IF (s(ibry,ifield,ng)%gradient) THEN
797 string(ibry)='Gra '
798 ELSE IF (s(ibry,ifield,ng)%nested) THEN
799 string(ibry)='Nes '
800 ELSE IF (s(ibry,ifield,ng)%periodic) THEN
801 string(ibry)='Per '
802 ELSE IF (s(ibry,ifield,ng)%radiation) THEN
803 IF (s(ibry,ifield,ng)%nudging) THEN
804 string(ibry)='RadNud '
805 ELSE
806 string(ibry)='Rad '
807 END IF
808 ELSE IF (s(ibry,ifield,ng)%reduced) THEN
809 string(ibry)='Red '
810 ELSE IF (s(ibry,ifield,ng)%Shchepetkin) THEN
811 string(ibry)='Shc '
812 END IF
813 END DO
814 IF (ifield.eq.nlbcvar) newline=' '
815 WRITE (line,frmt) trim(vname(1,idbvar(ifield))), &
816 & string(iwest), &
817 & string(isouth), &
818 & string(ieast), &
819 & string(inorth), &
820 & newline
821 lstr=len_trim(line)
822 ie=is+lstr
823 lbc_att(is:ie)=line(1:lstr)
824 is=ie
825 END IF
826 END DO
827!
828! Write attribute to NetCDF file.
829!
830 status=nf90_put_att(ncid, nf90_global, trim(aname), &
831 & trim(lbc_att))
832 IF (founderror(status, nf90_noerr, __line__, myfile)) THEN
833 exit_flag=3
834 ioerror=status
835 RETURN
836 END IF
837!
838 RETURN
839 END SUBROUTINE lbc_putatt_nf90
840
841#if defined PIO_LIB && defined DISTRIBUTE
842!
843!
844!***********************************************************************
845 SUBROUTINE lbc_putatt_pio (ng, pioFile, ncname, aname, S, status)
846!***********************************************************************
847! !
848! This routine writes lateral boundary conditions keywords strings !
849! into specified output NetCDF file global attribute. !
850! !
851! On Input: !
852! !
853! ng Nested grid number (integer) !
854! pioFile PIO file descriptor structure, TYPE(File_desc_t) !
855! pioFile%fh file handler !
856! pioFile%iosystem IO system descriptor (struct) !
857! ncname PIO filename (character) !
858! aname PIO global attribute name (character) !
859! S Derived type structure, TYPE(T_LBC) !
860! !
861! On Output: !
862! !
863! exit_flag Error flag (integer) stored in MOD_SCALARS !
864! ioerror NetCDF return code (integer) stored in MOD_IOUNITS !
865! status NetCDF return code (integer) !
866! !
867!***********************************************************************
868!
869 USE mod_param
870 USE mod_parallel
871 USE mod_iounits
872 USE mod_ncparam
873 USE mod_scalars
874 USE pio
875!
876 USE strings_mod, ONLY : founderror
877!
878 implicit none
879!
880! Imported variable declarations.
881!
882 integer, intent(in) :: ng
883 integer, intent(out) :: status
884!
885 character (*), intent(in) :: ncname
886 character (*), intent(in) :: aname
887!
888 TYPE (File_desc_t), intent(in) :: pioFile
889 TYPE(t_lbc), intent(in) :: S(4,nLBCvar,Ngrids)
890!
891! Local variable declarations
892!
893 integer :: i, ibry, ie, ifield, is, lstr, lvar
894!
895 character (len= 1) :: newline
896 character (len= 7) :: string(4)
897 character (len= 21) :: frmt
898 character (len= 100) :: line
899 character (len=2816) :: lbc_att
900
901 character (len=*), parameter :: MyFile = &
902 & __FILE__//", lbc_putatt_pio"
903!
904!-----------------------------------------------------------------------
905! Write lateral boundary conditions global attribute.
906!-----------------------------------------------------------------------
907!
908! Determine maximum length of state variable length.
909!
910 lvar=0
911 DO ifield=1,nlbcvar
912 IF (idbvar(ifield).gt.0) THEN
913 lvar=max(lvar, len_trim(vname(1,idbvar(ifield))))
914 END IF
915 END DO
916 WRITE (frmt,10) "(a,':',t", lvar+4, ",5a)"
917 10 FORMAT (a,i0,a)
918!
919! Initialize attribute.
920!
921 newline=char(10) ! Line Feed (LF) character for
922 lstr=len_trim(newline) ! attribute clarity with "ncdump"
923 DO i=1,len(lbc_att)
924 lbc_att(i:i)=' '
925 END DO
926 lbc_att(1:lstr)=newline(1:lstr)
927 is=lstr+1
928 WRITE (line,frmt) 'EDGE', &
929 & 'WEST ', &
930 & 'SOUTH ', &
931 & 'EAST ', &
932 & 'NORTH ', &
933 & newline(1:lstr)
934 lstr=len_trim(line)
935 ie=is+lstr
936 lbc_att(is:ie)=line(1:lstr)
937 is=ie
938!
939! Check if the local string "lbc_att" is big enough to store the
940! lateral boundary conditions global attribute.
941!
942 lstr=(nlbcvar+1)*(29+lvar+4)+1
943 IF (len(lbc_att).lt.lstr) THEN
944 IF (master) THEN
945 WRITE (stdout,20) len(lbc_att), lstr
946 20 FORMAT (/,' LBC_PUTATT_PIO - Length of local string lbc_att', &
947 & ' too small',/,18x,'Current = ',i5,' Needed = ',i5)
948 END IF
949 exit_flag=5
950 RETURN
951 END IF
952!
953! Build attribute string.
954!
955 DO ifield=1,nlbcvar
956 IF (idbvar(ifield).gt.0) THEN
957 DO ibry=1,4
958 IF (s(ibry,ifield,ng)%Chapman_explicit) THEN
959 string(ibry)='Che '
960 ELSE IF (s(ibry,ifield,ng)%Chapman_implicit) THEN
961 string(ibry)='Cha '
962 ELSE IF (s(ibry,ifield,ng)%clamped) THEN
963 string(ibry)='Cla '
964 ELSE IF (s(ibry,ifield,ng)%closed) THEN
965 string(ibry)='Clo '
966 ELSE IF (s(ibry,ifield,ng)%Flather) THEN
967 string(ibry)='Fla '
968 ELSE IF (s(ibry,ifield,ng)%gradient) THEN
969 string(ibry)='Gra '
970 ELSE IF (s(ibry,ifield,ng)%nested) THEN
971 string(ibry)='Nes '
972 ELSE IF (s(ibry,ifield,ng)%periodic) THEN
973 string(ibry)='Per '
974 ELSE IF (s(ibry,ifield,ng)%radiation) THEN
975 IF (s(ibry,ifield,ng)%nudging) THEN
976 string(ibry)='RadNud '
977 ELSE
978 string(ibry)='Rad '
979 END IF
980 ELSE IF (s(ibry,ifield,ng)%reduced) THEN
981 string(ibry)='Red '
982 ELSE IF (s(ibry,ifield,ng)%Shchepetkin) THEN
983 string(ibry)='Shc '
984 END IF
985 END DO
986 IF (ifield.eq.nlbcvar) newline=' '
987 WRITE (line,frmt) trim(vname(1,idbvar(ifield))), &
988 & string(iwest), &
989 & string(isouth), &
990 & string(ieast), &
991 & string(inorth), &
992 & newline
993 lstr=len_trim(line)
994 ie=is+lstr
995 lbc_att(is:ie)=line(1:lstr)
996 is=ie
997 END IF
998 END DO
999!
1000! Write attribute to NetCDF file.
1001!
1002 status=pio_put_att(piofile, pio_global, trim(aname), &
1003 & trim(lbc_att))
1004 IF (founderror(status, pio_noerr, __line__, myfile)) THEN
1005 exit_flag=3
1006 ioerror=status
1007 RETURN
1008 END IF
1009!
1010 RETURN
1011 END SUBROUTINE lbc_putatt_pio
1012#endif
1013!
1014!***********************************************************************
1015 SUBROUTINE lbc_report (iunit, ifield, S)
1016!***********************************************************************
1017! !
1018! This routine reports lateral boundary conditions to requested !
1019! state variable. !
1020! !
1021! On Input: !
1022! !
1023! iunit Output logical unit (integer) !
1024! ifield Lateral boundary variable index (integer) !
1025! S Derived type structure, TYPE(T_LBC) !
1026! !
1027!***********************************************************************
1028!
1029 USE mod_param
1030 USE mod_ncparam
1031 USE mod_scalars
1032!
1033 implicit none
1034!
1035! Imported variable declarations.
1036!
1037 integer, intent(in) :: ifield, iunit
1038
1039 TYPE(t_lbc), intent(in) :: S(4,nLBCvar,Ngrids)
1040!
1041! Local variable declarations
1042!
1043 integer :: ibry, ng
1044
1045 character (len=11) :: string(4,Ngrids)
1046!
1047!-----------------------------------------------------------------------
1048! Report lateral boundary conditions.
1049!-----------------------------------------------------------------------
1050!
1051! Build information strings.
1052!
1053 DO ng=1,ngrids
1054 DO ibry=1,4
1055 IF (s(ibry,ifield,ng)%Chapman_explicit) THEN
1056 string(ibry,ng)='Chapman Exp'
1057 ELSE IF (s(ibry,ifield,ng)%Chapman_implicit) THEN
1058 string(ibry,ng)='Chapman Imp'
1059 ELSE IF (s(ibry,ifield,ng)%clamped) THEN
1060 string(ibry,ng)='Clamped '
1061 ELSE IF (s(ibry,ifield,ng)%closed) THEN
1062 string(ibry,ng)='Closed '
1063 ELSE IF (s(ibry,ifield,ng)%Flather) THEN
1064 string(ibry,ng)='Flather '
1065 ELSE IF (s(ibry,ifield,ng)%gradient) THEN
1066 string(ibry,ng)='Gradient '
1067 ELSE IF (s(ibry,ifield,ng)%nested) THEN
1068 string(ibry,ng)='Nested '
1069 ELSE IF (s(ibry,ifield,ng)%periodic) THEN
1070 string(ibry,ng)='Periodic '
1071 ELSE IF (s(ibry,ifield,ng)%radiation) THEN
1072 IF (s(ibry,ifield,ng)%nudging) THEN
1073 string(ibry,ng)='Rad + Nud '
1074 ELSE
1075 string(ibry,ng)='Radiation '
1076 END IF
1077 ELSE IF (s(ibry,ifield,ng)%reduced) THEN
1078 string(ibry,ng)='Reduced '
1079 ELSE IF (s(ibry,ifield,ng)%Shchepetkin) THEN
1080 string(ibry,ng)='Shchepetkin'
1081 END IF
1082 END DO
1083 END DO
1084!
1085! Report.
1086!
1087 DO ng=1,ngrids
1088 IF (ng.eq.1) THEN
1089 WRITE (iunit,10) trim(vname(1,idbvar(ifield))), ng, &
1090 & trim(string(1,ng)), &
1091 & trim(string(2,ng)), &
1092 & trim(string(3,ng)), &
1093 & trim(string(4,ng))
1094 ELSE
1095 WRITE (iunit,20) ng, &
1096 & trim(string(1,ng)), &
1097 & trim(string(2,ng)), &
1098 & trim(string(3,ng)), &
1099 & trim(string(4,ng))
1100 END IF
1101 END DO
1102!
1103! Check periodic bounadry conditions: both east and west or north and
1104! need to be specified.
1105!
1106 DO ng=1,ngrids
1107 IF (.not.s(iwest,ifield,ng)%periodic.and. &
1108 & s(ieast,ifield,ng)%periodic) THEN
1109 WRITE (iunit,30) 'Western Edge boundary', &
1110 & trim(vname(1,idbvar(ifield))), ng
1111 exit_flag=5
1112 ELSE IF (.not.s(ieast,ifield,ng)%periodic.and. &
1113 & s(iwest,ifield,ng)%periodic) THEN
1114 WRITE (iunit,30) 'Eastern Edge boundary', &
1115 & trim(vname(1,idbvar(ifield))), ng
1116 exit_flag=5
1117 ELSE IF (.not.s(inorth,ifield,ng)%periodic.and. &
1118 & s(isouth,ifield,ng)%periodic) THEN
1119 WRITE (iunit,30) 'Northern Edge boundary', &
1120 & trim(vname(1,idbvar(ifield))), ng
1121 exit_flag=5
1122 ELSE IF (.not.s(isouth,ifield,ng)%periodic.and. &
1123 & s(inorth,ifield,ng)%periodic) THEN
1124 WRITE (iunit,30) 'Southern Edge boundary', &
1125 & trim(vname(1,idbvar(ifield))), ng
1126 exit_flag=5
1127 END IF
1128 END DO
1129
1130 10 FORMAT (/,1x,a,t26,i2,t31,a,t44,a,t57,a,t70,a)
1131 20 FORMAT (t26,i2,t31,a,t44,a,t57,a,t70,a)
1132 30 FORMAT (/,' LBC_REPORT - illegal configuration: The ',a, &
1133 & ' needs to be periodic too!',/,14x,'Variable: ',a,3x, &
1134 & 'Grid = ',i2.2)
1135!
1136 RETURN
1137 END SUBROUTINE lbc_report
1138
1139 END MODULE lbc_mod
Definition lbc.F:2
subroutine lbc_report(iunit, ifield, s)
Definition lbc.F:1016
subroutine lbc_putatt_pio(ng, piofile, ncname, aname, s, status)
Definition lbc.F:846
subroutine lbc_putatt_nf90(ng, ncid, ncname, aname, s, status)
Definition lbc.F:677
subroutine lbc_getatt_nf90(ng, model, ncid, ncname, aname, s)
Definition lbc.F:46
subroutine lbc_getatt_pio(ng, model, piofile, ncname, aname, s)
Definition lbc.F:371
integer ioerror
integer stdout
character(len=256) sourcefile
integer, dimension(:), allocatable idbvar
character(len=maxlen), dimension(6, 0:nv) vname
logical inpthread
logical master
integer, parameter iwest
integer exit_flag
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer noerror
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52