ana_biology.h

Go to the documentation of this file.

!!
      SUBROUTINE ana_biology (ng, tile, model)
!
!! git $Id$
!!======================================================================
!! Copyright (c) 2002-2025 The ROMS Group                              !
!!   Licensed under a MIT/X style license                              !
!!   See License_ROMS.md                                               !
!=======================================================================
!                                                                      !
!  This routine sets initial conditions for biological tracer fields   !
!  using analytical expressions.                                       !
!                                                                      !
!=======================================================================
!
      USE mod_param
      USE mod_ncparam
      USE mod_ocean
!
! Imported variable declarations.
!
      integer, intent(in) :: ng, tile, model
!
! Local variable declarations.
!
      character (len=*), parameter :: MyFile =                          &
     &  __FILE__
!
#include "tile.h"
!
      CALL ana_biology_tile (ng, tile, model,                           &
     &                       lbi, ubi, lbj, ubj,                        &
     &                       imins, imaxs, jmins, jmaxs,                &
     &                       ocean(ng) % t)
!
! Set analytical header file name used.
!
#ifdef DISTRIBUTE
      IF (lanafile) THEN
#else
      IF (lanafile.and.(tile.eq.0)) THEN
#endif
        ananame( 1)=myfile
      END IF
!
      RETURN
      END SUBROUTINE ana_biology
!
!***********************************************************************
      SUBROUTINE ana_biology_tile (ng, tile, model,                     &
     &                             LBi, UBi, LBj, UBj,                  &
     &                             IminS, ImaxS, JminS, JmaxS,          &
     &                             t)
!***********************************************************************
!
      USE mod_param
      USE mod_parallel
      USE mod_biology
      USE mod_ncparam
      USE mod_iounits
      USE mod_scalars
!
      USE stats_mod, ONLY : stats_3dfld
!
!  Imported variable declarations.
!
      integer, intent(in) :: ng, tile, model
      integer, intent(in) :: LBi, UBi, LBj, UBj
      integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
!
#ifdef ASSUMED_SHAPE
      real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
#else
      real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
#endif
!
!  Local variable declarations.
!
      logical, save :: first = .true.
 
      integer :: i, is, itrc, j, k
 
#if defined BIO_FENNEL || defined NEMURO
      real(r8) :: SiO4, cff1, cff2, temp
#elif defined ECOSIM
      real(r8) :: cff1, cff2, cff3, cff4, cff5, cff6, cff7, cff8, cff9
      real(r8) :: cff10, cff11, cff12, cff13, cff14, cff15
      real(r8) :: salt, sftm, temp
#endif
!
!   Maximum 80 biological tracers consider for field statistics.
!
      TYPE (T_STATS), save :: Stats(80)
 
#include "set_bounds.h"
!
!-----------------------------------------------------------------------
!  Initialize field statistics structure.
!-----------------------------------------------------------------------
!
      IF (first) THEN
        first=.false.
        DO i=1,SIZE(stats,1)
          stats(i) % checksum=0_i8b
          stats(i) % count=0.0_r8
          stats(i) % min=large
          stats(i) % max=-large
          stats(i) % avg=0.0_r8
          stats(i) % rms=0.0_r8
        END DO
      END IF
 
#if defined BIO_FENNEL
!
!-----------------------------------------------------------------------
!  Fennel et al. (2006), nitrogen-based biology model.
!-----------------------------------------------------------------------
!
      cff1=20.0_r8/3.0_r8
      cff2= 2.0_r8/3.0_r8
      DO k=1,n(ng)
        DO j=jstrt,jendt
          DO i=istrt,iendt
            temp=t(i,j,k,1,itemp)
            IF (temp.lt.8.0_r8) THEN
              sio4=30.0_r8
            ELSE IF ((temp.ge.8.0_r8).and.(temp.le.11.0_r8)) THEN
              sio4=30.0_r8-((temp-8.0_r8)*cff1)
            ELSE IF ((temp.gt.11.0_r8).and.(temp.le.13.0_r8)) THEN
              sio4=10.0_r8-((temp-11.0_r8)*4.0_r8)
            ELSE IF ((temp.gt.13.0_r8).and.(temp.le.16.0_r8)) THEN
              sio4=2.0_r8-((temp-13.0_r8)*cff2)
            ELSE IF (temp.gt.16.0_r8) THEN
              sio4=0.0_r8
            END IF
            t(i,j,k,1,ino3_)=1.67_r8+0.5873_r8*sio4+                    &
     &                               0.0144_r8*sio4**2+                 &
     &                               0.0003099_r8*sio4**3
            t(i,j,k,1,iphyt)=0.08_r8
            t(i,j,k,1,izoop)=0.06_r8
            t(i,j,k,1,inh4_)=0.1_r8
            t(i,j,k,1,ilden)=0.02_r8
            t(i,j,k,1,isden)=0.04_r8
            t(i,j,k,1,ichlo)=0.02_r8
#ifdef CARBON
            t(i,j,k,1,itic_)=2100.0_r8
            t(i,j,k,1,italk)=2350.0_r8
            t(i,j,k,1,ildec)=0.002_r8
            t(i,j,k,1,isdec)=0.06_r8
#endif
#ifdef OXYGEN
            t(i,j,k,1,ioxyg)=10.0_r8/0.02241_r8
#endif
          END DO
        END DO
      END DO
 
#elif defined NEMURO
!
!-----------------------------------------------------------------------
!  Nemuro lower trophic level ecosystem model.
!-----------------------------------------------------------------------
!
      cff1=20.0_r8/3.0_r8
      cff2= 2.0_r8/3.0_r8
      DO k=1,n(ng)
        DO j=jstrt,jendt
          DO i=istrt,iendt
            temp=t(i,j,k,1,itemp)
            IF (temp.lt.8.0_r8) THEN
              sio4=30.0_r8
            ELSE IF ((temp.ge.8.0_r8).and.(temp.le.11.0_r8)) THEN
              sio4=30.0_r8-((temp-8.0_r8)*cff1)
            ELSE IF ((temp.gt.11.0_r8).and.(temp.le.13.0_r8)) THEN
              sio4=10.0_r8-((temp-11.0_r8)*4.0_r8)
            ELSE IF ((temp.gt.13.0_r8).and.(temp.le.16.0_r8)) THEN
              sio4=2.0_r8-((temp-13.0_r8)*cff2)
            ELSE IF (temp.gt.16.0_r8) THEN
              sio4=0.0_r8
            END IF
            t(i,j,k,1,ino3_)=1.67_r8+0.5873_r8*sio4+                    &
     &                               0.0144_r8*sio4**2+                 &
     &                               0.0003099_r8*sio4**3
            t(i,j,k,1,isphy)=0.06_r8
            t(i,j,k,1,ilphy)=0.06_r8
            t(i,j,k,1,iszoo)=0.05_r8
            t(i,j,k,1,ilzoo)=0.05_r8
            t(i,j,k,1,ipzoo)=0.05_r8
            t(i,j,k,1,inh4_)=0.1_r8
            t(i,j,k,1,ipon_)=0.001_r8
            t(i,j,k,1,idon_)=0.001_r8
            t(i,j,k,1,isioh)=sio4
            t(i,j,k,1,iopal)=0.001_r8
          END DO
        END DO
      END DO
 
#elif defined NPZD_FRANKS || defined NPZD_POWELL
!
!-----------------------------------------------------------------------
!  NPZD biology model.
!-----------------------------------------------------------------------
!
      DO k=1,n(ng)
        DO j=jstrt,jendt
          DO i=istrt,iendt
            t(i,j,k,1,ino3_)=bioini(ino3_,ng)
            t(i,j,k,1,iphyt)=bioini(iphyt,ng)
            t(i,j,k,1,izoop)=bioini(izoop,ng)
            t(i,j,k,1,isdet)=bioini(isdet,ng)
          END DO
        END DO
      END DO
 
#elif defined NPZD_IRON
!
!-----------------------------------------------------------------------
!  NPZD biology model with or without iron limitation on phytoplankton
!  growth.
!-----------------------------------------------------------------------
!
      DO k=1,n(ng)
        DO j=jstrt,jendt
          DO i=istrt,iendt
            t(i,j,k,1,ino3_)=bioini(ino3_,ng)
            t(i,j,k,1,iphyt)=bioini(iphyt,ng)
            t(i,j,k,1,izoop)=bioini(izoop,ng)
            t(i,j,k,1,isdet)=bioini(isdet,ng)
# ifdef IRON_LIMIT
            t(i,j,k,1,ifphy)=bioini(ifphy,ng)
            t(i,j,k,1,ifdis)=bioini(ifdis,ng)
# endif
          END DO
        END DO
      END DO
 
#elif defined ECOSIM
!
!---------------------------------------------------------------------
!  EcoSim initial fields.
!---------------------------------------------------------------------
!
! Assumed maximum temperature gradient.
!
      cff3=1.0_r8/14.0_r8
      cff4=1.0_r8/16.0_r8
      cff5=32.0_r8
      cff7=1.0_r8/0.0157_r8
      cff8=1.0_r8/6.625_r8
      cff9=1.0_r8/16.0_r8
      cff10=1.0_r8/15.0_r8
      cff11=1.0_r8/8.0_r8
      cff12=1.0_r8/128.0_r8
      cff13=1.0_r8/1000.0_r8
      cff14=1.0_r8/12.0_r8
      cff15=cff5*cff8*cff14                  ! mole N : gram Chl
 
      DO k=n(ng),1,-1
        DO j=jstrt,jendt
          DO i=istrt,iendt
!
! Initialization of surface chlorophyll.
!
            sftm=t(i,j,n(ng),1,itemp)
            temp=t(i,j,k,1,itemp)
# ifdef SALINITY
            salt=t(i,j,k,1,isalt)
# endif
            cff1=-0.0827_r8*sftm+2.6386_r8
            cff2=max(0.00001_r8,cff1*(1.0_r8-(sftm-temp)*cff3))
!
! Initialization of nutrients.
!
            t(i,j,k,1,inh4_)=0.053_r8*temp+0.7990_r8
            t(i,j,k,1,ino3_)=8.5_r8-cff2*cff15-t(i,j,k,1,inh4_)
            t(i,j,k,1,ipo4_)=(t(i,j,k,1,inh4_)+t(i,j,k,1,ino3_))*cff4
            t(i,j,k,1,ifeo_)=1.0_r8
!
! Assuming diatoms are 75% of initialized chlorophyll.
!
            t(i,j,k,1,isio_)=5.5_r8-(cff2*0.75_r8)*cff15*1.20_r8
            t(i,j,k,1,idic_)=2000.0_r8
!
! Bacteria Initialization.
!
            DO is=1,nbac
              t(i,j,k,1,ibacc(is))=0.85_r8
              t(i,j,k,1,ibacn(is))=t(i,j,k,1,ibacc(is))*n2cbac(ng)
              t(i,j,k,1,ibacp(is))=t(i,j,k,1,ibacc(is))*p2cbac(ng)
              t(i,j,k,1,ibacf(is))=t(i,j,k,1,ibacc(is))*fe2cbac(ng)
            END DO
!
! Initialize phytoplankton populations.
!
            t(i,j,k,1,iphyc(1))=max(0.02_r8,                            &
     &                              0.75_r8*0.75_r8*cff5*cff2*cff14)
            t(i,j,k,1,iphyc(2))=max(0.02_r8,                            &
     &                              0.75_r8*0.25_r8*cff5*cff2*cff14)
            t(i,j,k,1,iphyc(3))=max(0.02_r8,                            &
     &                              0.125_r8*cff5*cff2*cff14)
            t(i,j,k,1,iphyc(4))=t(i,j,k,1,iphyc(3))
            DO is=1,nphy
              t(i,j,k,1,iphyn(is))=t(i,j,k,1,iphyc(is))*cff8
              t(i,j,k,1,iphyp(is))=t(i,j,k,1,iphyn(is))*cff4
              t(i,j,k,1,iphyf(is))=t(i,j,k,1,iphyc(is))*cff13
              IF (iphys(is).gt.0) THEN
                t(i,j,k,1,iphys(is))=t(i,j,k,1,iphyn(is))*1.20_r8
              END IF
!
!  Initialize Pigments in ugrams/liter (not umole/liter).
!  Chlorophyll-a
!
              cff6=12.0_r8/cff5
              t(i,j,k,1,ipigs(is,1))=cff6*t(i,j,k,1,iphyc(is))
!
!  Chlorophyll-b.
!
              cff6=cff5-b_c2cl(is,ng)
              IF (ipigs(is,2).gt.0) THEN
                 t(i,j,k,1,ipigs(is,2))=t(i,j,k,1,ipigs(is,1))*         &
     &                                  (b_chlb(is,ng)+                 &
     &                                   mxchlb(is,ng)*cff6)
              END IF
!
!  Chlorophyll-c.
!
              IF (ipigs(is,3).gt.0) THEN
                 t(i,j,k,1,ipigs(is,3))=t(i,j,k,1,ipigs(is,1))*         &
     &                                  (b_chlc(is,ng)+                 &
     &                                   mxchlc(is,ng)*cff6)
              END IF
!
!  Photosynthetic Carotenoids.
!
              IF (ipigs(is,4).gt.0) THEN
                 t(i,j,k,1,ipigs(is,4))=t(i,j,k,1,ipigs(is,1))*         &
     &                                  (b_psc(is,ng)+                  &
     &                                   mxpsc(is,ng)*cff6)
              END IF
!
!  Photoprotective Carotenoids.
!
              IF (ipigs(is,5).gt.0) THEN
                 t(i,j,k,1,ipigs(is,5))=t(i,j,k,1,ipigs(is,1))*         &
     &                                  (b_ppc(is,ng)+                  &
     &                                   mxppc(is,ng)*cff6)
              END IF
!
!  Low Urobilin Phycoeurythin Carotenoids.
!
              IF (ipigs(is,6).gt.0) THEN
                 t(i,j,k,1,ipigs(is,6))=t(i,j,k,1,ipigs(is,1))*         &
     &                                  (b_lpub(is,ng)+                 &
     &                                   mxlpub(is,ng)*cff6)
              END IF
!
!  High Urobilin Phycoeurythin Carotenoids.
!
              IF (ipigs(is,7).gt.0) THEN
                 t(i,j,k,1,ipigs(is,7))=t(i,j,k,1,ipigs(is,1))*         &
     &                                  (b_hpub(is,ng)+                 &
     &                                   mxhpub(is,ng)*cff6)
              END IF
            END DO
!
! DOC initialization.
!
# ifdef SALINITY
            cff6=max(0.001_r8,-0.9833_r8*salt+33.411_r8)
# else
            cff6=0.0_r8
# endif
            t(i,j,k,1,idomc(1))=0.1_r8
            t(i,j,k,1,idomn(1))=t(i,j,k,1,idomc(1))*cff8
            t(i,j,k,1,idomp(1))=t(i,j,k,1,idomn(1))*cff9
            t(i,j,k,1,icdmc(1))=t(i,j,k,1,idomc(1))*cdocfrac_c(1,ng)
            t(i,j,k,1,idomc(2))=15.254_r8*cff6+70.0_r8
            t(i,j,k,1,idomn(2))=t(i,j,k,1,idomc(2))*cff10
            t(i,j,k,1,idomp(2))=0.0_r8
            t(i,j,k,1,icdmc(2))=(0.243_r8*cff6+0.055_r8)*cff7
!
! Fecal Initialization.
!
            DO is=1,nfec
              t(i,j,k,1,ifecc(is))=0.002_r8
              t(i,j,k,1,ifecn(is))=t(i,j,k,1,ifecc(is))*cff11
              t(i,j,k,1,ifecp(is))=t(i,j,k,1,ifecc(is))*cff12
              t(i,j,k,1,ifecf(is))=t(i,j,k,1,ifecc(is))*cff13
              t(i,j,k,1,ifecs(is))=t(i,j,k,1,ifecc(is))*cff11
            END DO
          END DO
        END DO
      END DO
#endif
!
!  Report statistics.
!
      DO itrc=1,nbt
        i=idbio(itrc)
        CALL stats_3dfld (ng, tile, inlm, r3dvar, stats(itrc), 0,       &
     &                    lbi, ubi, lbj, ubj, 1, n(ng), t(:,:,:,1,i))
        IF (domain(ng)%NorthEast_Corner(tile)) THEN
          WRITE (stdout,10) trim(vname(2,idtvar(i)))//': '//            &
     &                      trim(vname(1,idtvar(i))),                   &
     &                      ng, stats(itrc)%min, stats(itrc)%max
        END IF
      END DO
!
  10  FORMAT (3x,' ANA_BIOLOGY - ',a,/,19x,                             &
     &        '(Grid = ',i2.2,', Min = ',1p,e15.8,0p,                   &
     &                         ' Max = ',1p,e15.8,0p,')')
!
      RETURN
      END SUBROUTINE ana_biology_tile