2
3
4
5
6
7
8
9
10
11
12
13
14
20
22
23 implicit none
24
25
26
27 logical, intent(in) :: Lwrite
28 integer, intent(in) :: model, inp, out
29
30
31
32 integer :: Npts, Nval
33 integer :: iTrcStr, iTrcEnd
34 integer :: i, ifield, igrid, itracer, itrc, ng, nline, status
35
36 logical, dimension(NBT,Ngrids) :: Ltrc
37
38 real(r8), dimension(NBT,Ngrids) :: Rbio
39
40 real(dp), dimension(nRval) :: Rval
41
42 character (len=40 ) :: KeyWord
43 character (len=256) :: line
44 character (len=256), dimension(nCval) :: Cval
45
46
47
48
49
50 igrid=1
51 itracer=0
52 itrcstr=1
54 nline=0
55
56
57
58
59
60#ifdef ANA_BIOLOGY
62#endif
63 DO WHILE (.true.)
64 READ (inp,'(a)',err=10,END=20) line
66 IF (status.gt.0) THEN
67 SELECT CASE (trim(keyword))
68 CASE ('Lbiology')
70 CASE ('BioIter')
72#ifdef ANA_BIOLOGY
73 CASE ('BioIni(iNO3_)')
75 CASE ('BioIni(iPhyt)')
77 CASE ('BioIni(iZoop)')
79 CASE ('BioIni(iSDet)')
81#endif
82 CASE ('K_ext')
84 CASE ('K_NO3')
86 CASE ('K_Phy')
88 CASE ('Vm_NO3')
90 CASE ('PhyMR')
92 CASE ('ZooGR')
94 CASE ('ZooMR')
96 CASE ('ZooMD')
98 CASE ('ZooGA')
100 CASE ('ZooEC')
102 CASE ('DetRR')
104 CASE ('wDet')
106 CASE ('TNU2')
112 END DO
113 END DO
114 CASE ('TNU4')
120 END DO
121 END DO
122 CASE ('ad_TNU2')
129 END DO
130 END DO
131 CASE ('ad_TNU4')
138 END DO
139 END DO
140 CASE ('LtracerSponge')
146 END DO
147 END DO
148 CASE ('AKT_BAK')
154 END DO
155 END DO
156 CASE ('ad_AKT_fac')
163 END DO
164 END DO
165 CASE ('TNUDG')
170 tnudg(i,ng)=rbio(itrc,ng)
171 END DO
172 END DO
173 CASE ('Hadvection')
174 IF (itracer.lt.
nbt)
THEN
175 itracer=itracer+1
176 ELSE
177 itracer=1
178 END IF
180 npts=
load_tadv(nval, cval, line, nline, itrc, igrid, &
184 CASE ('Vadvection')
185 IF (itracer.lt.
nbt)
THEN
186 itracer=itracer+1
187 ELSE
188 itracer=1
189 END IF
191 npts=
load_tadv(nval, cval, line, nline, itrc, igrid, &
195#if defined ADJOINT || defined TANGENT || defined TL_IOMS
196 CASE ('ad_Hadvection')
197 IF (itracer.lt.
nbt)
THEN
198 itracer=itracer+1
199 ELSE
200 itracer=1
201 END IF
203 npts=
load_tadv(nval, cval, line, nline, itrc, igrid, &
207 CASE ('Vadvection')
208 IF (itracer.lt.(
nbt)
THEN
209 itracer=itracer+1
210 ELSE
211 itracer=1
212 END IF
214 npts=
load_tadv(nval, cval, line, nline, itrc, igrid, &
218#endif
219 CASE ('LBC(isTvar)')
220 IF (itracer.lt.
nbt)
THEN
221 itracer=itracer+1
222 ELSE
223 itracer=1
224 END IF
226 npts=
load_lbc(nval, cval, line, nline, ifield, igrid, &
229#if defined ADJOINT || defined TANGENT || defined TL_IOMS
230 CASE ('ad_LBC(isTvar)')
231 IF (itracer.lt.
nbt)
THEN
232 itracer=itracer+1
233 ELSE
234 itracer=1
235 END IF
237 npts=
load_lbc(nval, cval, line, nline, ifield, igrid, &
240#endif
241 CASE ('LtracerSrc')
247 END DO
248 END DO
249 CASE ('LtracerCLM')
255 END DO
256 END DO
257 CASE ('LnudgeTCLM')
263 END DO
264 END DO
265 CASE ('Hout(idTvar)')
270 IF (i.eq.0) THEN
271 IF (
master)
WRITE (out,30) &
272 & 'idTvar(idbio(', itrc, '))'
274 RETURN
275 END IF
276 hout(i,ng)=ltrc(itrc,ng)
277 END DO
278 END DO
279 CASE ('Hout(idTsur)')
284 IF (i.eq.0) THEN
285 IF (
master)
WRITE (out,30) &
286 & 'idTsur(idbio(', itrc, '))'
288 RETURN
289 END IF
290 hout(i,ng)=ltrc(itrc,ng)
291 END DO
292 END DO
293 CASE ('Qout(idTvar)')
298 qout(i,ng)=ltrc(itrc,ng)
299 END DO
300 END DO
301 CASE ('Qout(idsurT)')
306 IF (i.eq.0) THEN
307 IF (
master)
WRITE (out,30) &
308 & 'idsurT(idbio(', itrc, '))'
310 RETURN
311 END IF
312 qout(i,ng)=ltrc(itrc,ng)
313 END DO
314 END DO
315 CASE ('Qout(idTsur)')
320 qout(i,ng)=ltrc(itrc,ng)
321 END DO
322 END DO
323#if defined AVERAGES || \
324 (defined ad_averages && defined adjoint) || \
325 (defined rp_averages && defined tl_ioms) || \
326 (defined tl_averages && defined tangent)
327 CASE ('Aout(idTvar)')
332 aout(i,ng)=ltrc(itrc,ng)
333 END DO
334 END DO
335 CASE ('Aout(idTTav)')
340 aout(i,ng)=ltrc(itrc,ng)
341 END DO
342 END DO
343 CASE ('Aout(idUTav)')
348 aout(i,ng)=ltrc(itrc,ng)
349 END DO
350 END DO
351 CASE ('Aout(idVTav)')
356 aout(i,ng)=ltrc(itrc,ng)
357 END DO
358 END DO
359 CASE ('Aout(iHUTav)')
364 aout(i,ng)=ltrc(itrc,ng)
365 END DO
366 END DO
367 CASE ('Aout(iHVTav)')
372 aout(i,ng)=ltrc(itrc,ng)
373 END DO
374 END DO
375#endif
376#ifdef DIAGNOSTICS_TS
377 CASE ('Dout(iTrate)')
383 END DO
384 END DO
385 CASE ('Dout(iThadv)')
391 END DO
392 END DO
393 CASE ('Dout(iTxadv)')
399 END DO
400 END DO
401 CASE ('Dout(iTyadv)')
407 END DO
408 END DO
409 CASE ('Dout(iTvadv)')
415 END DO
416 END DO
417# if defined TS_DIF2 || defined TS_DIF4
418 CASE ('Dout(iThdif)')
424 END DO
425 END DO
426 CASE ('Dout(iTxdif)')
432 END DO
433 END DO
434 CASE ('Dout(iTydif)')
440 END DO
441 END DO
442# if defined MIX_GEO_TS || defined MIX_ISO_TS
443 CASE ('Dout(iTsdif)')
449 END DO
450 END DO
451# endif
452# endif
453 CASE ('Dout(iTvdif)')
459 END DO
460 END DO
461#endif
462 END SELECT
463 END IF
464 END DO
465 10
IF (
master)
WRITE (out,40) line
467 RETURN
468 20 CONTINUE
469
470
471
472
473
474 IF (
master.and.lwrite)
THEN
477 WRITE (out,50) ng
478 WRITE (out,60)
bioiter(ng),
'BioIter', &
479 & 'Number of iterations for nonlinear convergence.'
480#ifdef ANA_BIOLOGY
481 WRITE (out,70)
bioini(
ino3_,ng),
'BioIni(iNO3_)', &
482 & 'Nitrate initial concentration (mmol/m3).'
483 WRITE (out,70)
bioini(
iphyt,ng),
'BioIni(iPhyt)', &
484 & 'Phytoplankton initial concentration (mmol/m3).'
485 WRITE (out,70)
bioini(
izoop,ng),
'BioIni(iZoop)', &
486 & 'Zooplankton initial concentration (mmol/m3).'
487 WRITE (out,70)
bioini(
isdet,ng),
'BioIni(iSDet)', &
488 & 'Small detritus initial concentration (mmol/m3).'
489#endif
490 WRITE (out,70)
k_ext(ng),
'K_ext', &
491 & 'Light extinction coefficient (m-1).'
492 WRITE (out,80)
k_no3(ng),
'K_NO3', &
493 & 'Inverse half-saturation for phytoplankton NO3', &
494 & 'uptake (1/(mmol m-3)).'
495 WRITE (out,80)
k_phy(ng),
'K_Phy', &
496 & 'Phytoplankton saturation coefficient', &
497 & '(mmol/m3)^2.'
498 WRITE (out,70)
vm_no3(ng),
'Vm_NO3', &
499 & 'Nitrate upatake rate (day-1).'
500 WRITE (out,70)
phymr(ng),
'PhyMR', &
501 & 'Phytoplankton senescence/mortality rate (day-1)'
502 WRITE (out,70)
zoogr(ng),
'ZooGR', &
503 & 'Zooplankton maximum growth rate (day-1).'
504 WRITE (out,70)
zoomr(ng),
'ZooMR', &
505 & 'Zooplankton mortality rate (day-1).'
506 WRITE (out,70)
zoomd(ng),
'ZooMD', &
507 & 'Zooplankton death bits rate (day-1).'
508 WRITE (out,70)
zooga(ng),
'ZooGA', &
509 & 'Zooplankton grazing inefficiency (nondimensional).'
510 WRITE (out,70)
zooec(ng),
'ZooEC', &
511 & 'Zooplankton excreted fraction (nondimensional).'
512 WRITE (out,70)
detrr(ng),
'DetRR', &
513 & 'Detritus remineralization rate (day-1).'
514 WRITE (out,70)
wdet(ng),
'wDet', &
515 & 'Detrital sinking rate (m/day).'
516#ifdef TS_DIF2
519 WRITE (out,90)
nl_tnu2(i,ng),
'nl_tnu2', i, &
520 & 'NLM Horizontal, harmonic mixing coefficient', &
521 &
'(m2/s) for tracer ', i, trim(
vname(1,
idtvar(i)))
522# ifdef ADJOINT
523 WRITE (out,90)
ad_tnu2(i,ng),
'ad_tnu2', i, &
524 & 'ADM Horizontal, harmonic mixing coefficient', &
525 &
'(m2/s) for tracer ', i, trim(
vname(1,
idtvar(i)))
526# endif
527# if defined TANGENT || defined TL_IOMS
528 WRITE (out,90)
tl_tnu2(i,ng),
'tl_tnu2', i, &
529 & 'TLM Horizontal, harmonic mixing coefficient', &
530 &
'(m2/s) for tracer ', i, trim(
vname(1,
idtvar(i)))
531# endif
532 END DO
533#endif
534#ifdef TS_DIF4
537 WRITE (out,90)
nl_tnu4(i,ng),
'nl_tnu4', i, &
538 & 'NLM Horizontal, biharmonic mixing coefficient', &
539 &
'(m4/s) for tracer ', i, trim(
vname(1,
idtvar(i)))
540# ifdef ADJOINT
541 WRITE (out,90)
ad_tnu4(i,ng),
'ad_tnu4', i, &
542 & 'ADM Horizontal, biharmonic mixing coefficient', &
543 &
'(m4/s) for tracer ', i, trim(
vname(1,
idtvar(i)))
544# endif
545# if defined TANGENT || defined TL_IOMS
546 WRITE (out,90)
tl_tnu4(i,ng),
'tl_tnu4', i, &
547 & 'TLM Horizontal, biharmonic mixing coefficient', &
548 &
'(m4/s) for tracer ', i, trim(
vname(1,
idtvar(i)))
549# endif
550 END DO
551#endif
556 & i, 'Turning ON sponge on tracer ', i, &
558 ELSE
560 & i, 'Turning OFF sponge on tracer ', i, &
562 END IF
563 END DO
566 WRITE(out,90)
akt_bak(i,ng),
'Akt_bak', i, &
567 & 'Background vertical mixing coefficient (m2/s)', &
569 END DO
570#ifdef FORWARD_MIXING
573# ifdef ADJOINT
574 WRITE (out,90)
ad_akt_fac(i,ng),
'ad_Akt_fac', i, &
575 & 'ADM basic state vertical mixing scale factor', &
577# endif
578# if defined TANGENT || defined TL_IOMS
579 WRITE (out,90)
tl_akt_fac(i,ng),
'tl_Akt_fac', i, &
580 & 'TLM basic state vertical mixing scale factor', &
582# endif
583 END DO
584#endif
587 WRITE (out,90)
tnudg(i,ng),
'Tnudg', i, &
588 & 'Nudging/relaxation time scale (days)', &
590 END DO
594 WRITE (out,100)
ltracersrc(i,ng),
'LtracerSrc', &
595 & i, 'Turning ON point sources/Sink on tracer ', i, &
597 ELSE
598 WRITE (out,100)
ltracersrc(i,ng),
'LtracerSrc', &
599 & i, 'Turning OFF point sources/Sink on tracer ', i, &
601 END IF
602 END DO
606 WRITE (out,100)
ltracerclm(i,ng),
'LtracerCLM', i, &
607 & 'Turning ON processing of climatology tracer ', i, &
609 ELSE
610 WRITE (out,100)
ltracerclm(i,ng),
'LtracerCLM', i, &
611 & 'Turning OFF processing of climatology tracer ', i, &
613 END IF
614 END DO
618 WRITE (out,100)
lnudgetclm(i,ng),
'LnudgeTCLM', i, &
619 & 'Turning ON nudging of climatology tracer ', i, &
621 ELSE
622 WRITE (out,100)
lnudgetclm(i,ng),
'LnudgeTCLM', i, &
623 & 'Turning OFF nudging of climatology tracer ', i, &
625 END IF
626 END DO
627 IF ((
nhis(ng).gt.0).and.any(
hout(:,ng)))
THEN
628 WRITE (out,'(1x)')
634 END DO
639 & 'Write out tracer flux ', i, &
641 END DO
642 END IF
643 IF ((
nqck(ng).gt.0).and.any(
qout(:,ng)))
THEN
644 WRITE (out,'(1x)')
650 END DO
655 & 'Write out surface tracer ', i, &
657 END DO
662 & 'Write out tracer flux ', i, &
664 END DO
665 END IF
666#if defined AVERAGES || \
667 (defined ad_averages && defined adjoint) || \
668 (defined rp_averages && defined tl_ioms) || \
669 (defined tl_averages && defined tangent)
670 IF ((
navg(ng).gt.0).and.any(
aout(:,ng)))
THEN
671 WRITE (out,'(1x)')
676 & 'Write out averaged tracer ', i, &
678 END DO
683 & 'Write out averaged <t*t> for tracer ', i, &
685 END DO
690 & 'Write out averaged <u*t> for tracer ', i, &
692 END DO
697 & 'Write out averaged <v*t> for tracer ', i, &
699 END DO
704 & 'Write out averaged <Huon*t> for tracer ', i, &
706 END DO
711 & 'Write out averaged <Hvom*t> for tracer ', i, &
713 END DO
714 END IF
715#endif
716#ifdef DIAGNOSTICS_TS
717 IF ((
ndia(ng).gt.0).and.any(
dout(:,ng)))
THEN
718 WRITE (out,'(1x)')
722 & WRITE (out,110) .true., 'Dout(iTrate)', &
723 & 'Write out rate of change of tracer ', itrc, &
725 END DO
729 & WRITE (out,110) .true., 'Dout(iThadv)', &
730 & 'Write out horizontal advection, tracer ', itrc, &
732 END DO
736 & WRITE (out,110) .true., 'Dout(iTxadv)', &
737 & 'Write out horizontal X-advection, tracer ', itrc, &
739 END DO
743 & WRITE (out,110) .true., 'Dout(iTyadv)', &
744 & 'Write out horizontal Y-advection, tracer ', itrc, &
746 END DO
750 & WRITE (out,110) .true., 'Dout(iTvadv)', &
751 & 'Write out vertical advection, tracer ', itrc, &
753 END DO
754# if defined TS_DIF2 || defined TS_DIF4
758 & WRITE (out,110) .true., 'Dout(iThdif)', &
759 & 'Write out horizontal diffusion, tracer ', itrc, &
761 END DO
765 & WRITE (out,110) .true., 'Dout(iTxdif)', &
766 & 'Write out horizontal X-diffusion, tracer ', itrc, &
768 END DO
772 & WRITE (out,110) .true., 'Dout(iTydif)', &
773 & 'Write out horizontal Y-diffusion, tracer ', itrc, &
775 END DO
776# if defined MIX_GEO_TS || defined MIX_ISO_TS
780 & WRITE (out,110) .true., 'Dout(iTsdif)', &
781 & 'Write out horizontal S-diffusion, tracer ', itrc, &
783 END DO
784# endif
785# endif
789 & WRITE (out,110) .true., 'Dout(iTvdif)', &
790 & 'Write out vertical diffusion, tracer ', itrc, &
792 END DO
793 END IF
794#endif
795 END IF
796 END DO
797 END IF
798
799
800
801
802
803
804
805
810#ifdef ADJOINT
812#endif
813#if defined TANGENT || defined TL_IOMS
815#endif
816
817
818
819 IF (
tnudg(i,ng).gt.0.0_r8)
THEN
821 ELSE
823 END IF
824 END DO
825 END DO
826
827 30 FORMAT (/,' read_BioPar - variable info not yet loaded, ', &
828 & a,i2.2,a)
829 40 FORMAT (/,' read_BioPar - Error while processing line: ',/,a)
830 50 FORMAT (/,/,' NPZD Model Parameters, Grid: ',i2.2, &
831 & /, ' ===============================',/)
832 60 FORMAT (1x,i10,2x,a,t32,a)
833 70 FORMAT (1p,e11.4,2x,a,t32,a)
834 80 FORMAT (1p,e11.4,2x,a,t32,a,/,t34,a)
835 90 FORMAT (1p,e11.4,2x,a,'(',i2.2,')',t32,a,/,t34,a,i2.2,':',1x,a)
836 100 FORMAT (10x,l1,2x,a,'(',i2.2,')',t32,a,i2.2,':',1x,a)
837 110 FORMAT (10x,l1,2x,a,t32,a,i2.2,':',1x,a)
838
839 RETURN
integer function decode_line(line_text, keyword, nval, cval, rval)
integer function load_lbc(ninp, vinp, line, nline, ifield, igrid, itrcstr, itrcend, svname, s)
integer function load_tadv(ninp, vinp, line, nline, itrc, igrid, itracer, itrcstr, itrcend, svname, s)
real(r8), dimension(:), allocatable zooga
real(r8), dimension(:), allocatable zoomr
real(r8), dimension(:), allocatable k_phy
real(r8), dimension(:), allocatable detrr
real(r8), dimension(:), allocatable wdet
real(r8), dimension(:), allocatable zoomd
integer, dimension(:), allocatable bioiter
real(r8), dimension(:,:), allocatable bioini
real(r8), dimension(:), allocatable zoogr
real(r8), dimension(:), allocatable k_no3
real(r8), dimension(:), allocatable phymr
real(r8), dimension(:), allocatable zooec
integer, dimension(:), allocatable idbio
real(r8), dimension(:), allocatable vm_no3
real(r8), dimension(:), allocatable k_ext
integer, dimension(:), allocatable idttav
logical, dimension(:,:), allocatable hout
integer, dimension(:), allocatable ihvtav
integer, dimension(:), allocatable ihutav
integer, dimension(:), allocatable idutav
integer, dimension(:), allocatable idsurt
integer, dimension(:), allocatable idtsur
integer, dimension(:), allocatable idtvar
integer, dimension(:), allocatable istvar
logical, dimension(:,:), allocatable qout
character(len=maxlen), dimension(6, 0:nv) vname
integer, dimension(:), allocatable idvtav
logical, dimension(:,:), allocatable dout
logical, dimension(:,:), allocatable aout
integer, dimension(:,:), allocatable iddtrc
type(t_adv), dimension(:,:), allocatable hadvection
type(t_lbc), dimension(:,:,:), allocatable ad_lbc
type(t_adv), dimension(:,:), allocatable ad_hadvection
type(t_lbc), dimension(:,:,:), allocatable lbc
type(t_adv), dimension(:,:), allocatable vadvection
type(t_adv), dimension(:,:), allocatable ad_vadvection
real(r8), dimension(:,:), allocatable tl_tnu2
logical, dimension(:,:), allocatable ltracersrc
real(r8), dimension(:,:), allocatable nl_tnu2
integer, dimension(:), allocatable nqck
real(dp), dimension(:,:), allocatable tnudg
logical, dimension(:,:), allocatable ltracersponge
real(r8), dimension(:,:), allocatable tl_tnu4
integer, dimension(:), allocatable navg
real(r8), dimension(:,:), allocatable ad_akt_fac
real(r8), dimension(:,:), allocatable ad_tnu2
real(r8), dimension(:,:), allocatable akt_bak
integer, dimension(:), allocatable nhis
real(r8), dimension(:,:), allocatable ad_tnu4
real(r8), dimension(:,:), allocatable nl_tnu4
logical, dimension(:), allocatable lbiology
logical, dimension(:,:), allocatable ltracerclm
integer, dimension(:), allocatable ndia
logical, dimension(:,:), allocatable lnudgetclm
real(r8), dimension(:,:), allocatable tl_akt_fac