Opened 11 months ago

Closed 11 months ago

Last modified 11 months ago

#839 closed upgrade (Done)

VERY IMPORTANT: Tracer Advection Revisited

Reported by: arango Owned by:
Priority: major Milestone: Release ROMS/TOMS 3.9
Component: Nonlinear Version: 3.9
Keywords: Cc:

Description (last modified by arango)

All the tracer advection CPP options were removed and replaced with logical switches to facilitate applying a desired horizontal and vertical scheme to each tracer. One may have a particular scheme for temperature and a different for salinity, inert, biological, or sediment tracer. It is usually desirable to have a positive-definite or monotonic algorithm for salinity and passive tracers.

  • WARNING: ROMS standard input files for physics, biology, and sediment were modified to include the tracer advection switches. In roms.in we have:
    ! Set horizontal and vertical advection schemes for active and inert
    ! tracers. A different advection scheme is allowed for each tracer.
    ! For example, a positive-definite (monotonic) algorithm can be activated
    ! for salinity and inert tracers, while a different one is set for
    ! temperature. [1:NAT+NPT,Ngrids] values are expected.
    !
    !   Keyword    Advection Algorithm
    !
    !   A4         4th-order Akima (horizontal/vertical)
    !   C2         2nd-order centered differences (horizontal/vertical)
    !   C4         4th-order centered differences (horizontal/vertical)
    !   HSIMT      3th-order HSIMT-TVD (horizontal/vertical)
    !   MPDATA     recursive flux corrected MPDATA (horizontal/vertical)
    !   SPLINES    parabolic splines (only vertical)
    !   SU3        split third-order upstream (horizontal/vertical)
    !   U3         3rd-order upstream-biased (only horizontal)
    !
    ! The user has the option of specifying the full Keyword or the first
    ! two letters, regardless if using uppercase or lowercase. If nested
    ! grids, specify values for each grid (see glossary below).
    
       Hadvection == U3       \                     ! temperature
                     MPDATA                         ! salinity
    
       Vadvection == C4       \                     ! temperature
                     MPDATA                         ! salinity
    
    ! Adjoint-based algorithms can have different horizontal and schemes
    ! for active and inert tracers.
    
    ad_Hadvection == U3       \                     ! temperature
                     U3                             ! salinity
    
    ad_Vadvection == C4       \                     ! temperature
                     C4                             ! salinity
    

Notice that similar switches are defined for the adjoint-based algorithms.

For the Fennel biological model, we specify the tracer advection switches in bio_fennel.in :

! Set horizontal and vertical advection schemes for biological tracers.
! A different advection scheme is allowed for each tracer. For example,
! a positive-definite (monotonic) algorithm can be activated for
! salinity and biological tracers, while a different one is set for
! temperature. [1:NAT+NPT,Ngrids] values are expected.
!
!   Keyword    Advection Algorithm
!
!   A4         4th-order Akima (horizontal/vertical)
!   C2         2nd-order centered differences (horizontal/vertical)
!   C4         4th-order centered differences (horizontal/vertical)
!   HSIMT      3th-order HSIMT-TVD (horizontal/vertical)
!   MPDATA     recursive flux corrected MPDATA (horizontal/vertical)
!   SPLINES    parabolic splines (only vertical)
!   SU3        split third-order upstream (horizontal/vertical)
!   U3         3rd-order upstream-biased (only horizontal)
!
! The user has the option of specifying the full Keyword or the first
! two letters, regardless if using uppercase or lowercase. If nested
! grids, specify values for each grid.

   Hadvection == HSIMT    \                     ! idbio( 1), NO3
                 HSIMT    \                     ! idbio( 2), NH4
                 HSIMT    \                     ! idbio( 3), chlorophyll
                 HSIMT    \                     ! idbio( 4), phytoplankton
                 HSIMT    \                     ! idbio( 5), zooplankton
                 HSIMT    \                     ! idbio( 6), LdetritusN
                 HSIMT    \                     ! idbio( 7), SdetritusN
                 HSIMT    \                     ! idbio( 8), LdetritusC
                 HSIMT    \                     ! idbio( 9), SdetritusC
                 HSIMT    \                     ! idbio(10), TIC
                 HSIMT    \                     ! idbio(11), alkalinity
                 HSIMT                          ! idbio(12), oxygen

   Vadvection == HSIMT    \                     ! idbio( 1), NO3
                 HSIMT    \                     ! idbio( 2), NH4
                 HSIMT    \                     ! idbio( 3), chlorophyll
                 HSIMT    \                     ! idbio( 4), phytoplankton
                 HSIMT    \                     ! idbio( 5), zooplankton
                 HSIMT    \                     ! idbio( 6), LdetritusN
                 HSIMT    \                     ! idbio( 7), SdetritusN
                 HSIMT    \                     ! idbio( 8), LdetritusC
                 HSIMT    \                     ! idbio( 9), SdetritusC
                 HSIMT    \                     ! idbio(10), TIC
                 HSIMT    \                     ! idbio(11), alkalinity
                 HSIMT                          ! idbio(12), oxygen

Alternatively, it advantageous to specify the tracer advection switches in a compact form. For example, in ecosim.in we have:

   Hadvection == HSIMT                          ! idbio(:), compact

   Vadvection == HSIMT                          ! idbio(:), compact

! Adjoint-based algorithms can have different horizontal and schemes
! for active and inert tracers.

ad_Hadvection == U3                             ! idbio(:), compact

ad_Vadvection == C4

That is all EcoSim tracers have the same tracer advection algorithm. A similar strategy can be used for the sediment tracers in sediment.in.

In a nested application, the syntax is as follows:

   Hadvection == A4       \                     ! temperature, Grid 1
                 A4       \                     ! temperature, Grid 2
                 A4       \                     ! temperature, Grid 3
                 HSIMT    \                     ! salinity,    Grid 1
                 HSIMT    \                     ! salinity,    Grid 2
                 HSIMT                          ! salinity,    Grid 3
 
   Vadvection == A4       \                     ! temperature, Grid 1
                 A4       \                     ! temperature, Grid 2
                 A4       \                     ! temperature, Grid 3
                 HSIMT    \                     ! salinity,    Grid 1
                 HSIMT    \                     ! salinity,    Grid 2
                 HSIMT                          ! salinity,    Grid 3

  • A new tracer scheme HSIMT (Wu and Zhu, 2010) that uses a Total Variation Diminishing (TVD) limiter is included. It is adapted from Hui Wu code and implemented in COAWST by Tarandeep Kalra and John Warner. The HSIMT algorithm can be used for salinity and passive tracer to achieve monotonicity. For more information, check the following publication:
        Hui Wu and Jianrong Zhu, 2010: Advection scheme with 3rd          
          high-order spatial interpolation at the middle temporal
          level and its application to saltwater intrusion in 
          Changjiang Estuary, Ocean Modelling 33, 33-51,                  
          doi:10.1016/j.ocemod.2009.12.001
    
    The HSIMT algorithm is inlined into ROMS code in step3d_t.F instead of calling a subroutine inside of DO-loops to avoid computational inefficiency and facilitate the tangent linear and adjoint transformation.

It is recommended to use HSIMT or MPDATA for passive tracers since the have a positive-definite dynamical range. The HSIMT is more modern and efficient than the recursive MPDATA algorithm. When using HSIMT or MPDATA, the user needs to specify the same scheme for both horizontal and vertical tracer advection.

  • All the tracer advection algorithms are available at run time. For example, in the step3d_t.F the horizontal advective fluxes are computed as:
             HADV_FLUX : IF (Hadvection(itrc,ng)%CENTERED2) THEN
    
               FX(i,j) = ...
               FE(i,j) = ...
    
             ELSE IF (Hadvection(itrc,ng)%MPDATA) THEN
    
               FX(i,j) = ...
               FE(i,j) = ...
    
             ELSE IF (Hadvection(itrc,ng)%HSIMT) THEN
    
               FX(i,j) = ...
               FE(i,j) = ...
    
             ELSE IF ((Hadvection(itrc,ng)%AKIMA4).or.                     &
        &             (Hadvection(itrc,ng)%CENTERED4).or.                  &
        &             (Hadvection(itrc,ng)%SPLIT_U3).or.                   &
        &             (Hadvection(itrc,ng)%UPSTREAM3)) THEN
    
               FX(i,j) = ...
               FE(i,j) = ...
    
             END IF HADV_FLUX
    
    and for vertical advective flux
              VADV_FLUX : IF (Vadvection(itrc,ng)%SPLINES) THEN
    
                FC(i,k) = ..
    
              ELSE IF (Vadvection(itrc,ng)%AKIMA4) THEN
    
                FC(i,k) = ..
    
              ELSE IF (Vadvection(itrc,ng)%CENTERED2) THEN
    
                FC(i,k) = ..
    
              ELSE IF (Vadvection(itrc,ng)%MPDATA) THEN
    
                FC(i,k) = ..
    
              ELSE IF (Vadvection(itrc,ng)%HSIMT) THEN
    
                FC(i,k) = ..
    
              ELSE IF ((Vadvection(itrc,ng)%CENTERED4).or.                  &
         &             (Vadvection(itrc,ng)%SPLIT_U3)) THEN
    
                FC(i,k) = ..
    
              END IF VADV_FLUX
    

A similar structure is found in pre_step3d.F.

  • The new tracer advection strategy was tested to get identical results regardless of the parallel partitions in distributed-memory (MPI), shared-memory (OpenMP), and serial with tiled partitions. It is only possible in double precision computations. In serial precision with double precision IO, the solution are not identical because of round-off. The user needs to be aware of the computational sensitivity in parallel single-precision applications. The solution differs as a function of the tile partition.
  • The HSIMT and MPDATA algorithms are not implemented in the tangent linear and adjoint codes. MPDATA is very tricky because of the recursive algorithm in the anti-diffusivity correction.
  • The tracer advection switches are reported in the output NetCDF files global attribute NLM_TADV:
                   :NLM_TADV = "\n",
                           "ADVECTION:      HORIZONTAL   VERTICAL     \n",
                           "temp:           Upstream3    Centered4    \n",
                           "salt:           HSIMT        HSIMT        \n",
                           "NO3:            HSIMT        HSIMT        \n",
                           "NH4:            HSIMT        HSIMT        \n",
                           "chlorophyll:    HSIMT        HSIMT        \n",
                           "phytoplankton:  HSIMT        HSIMT        \n",
                           "zooplankton:    HSIMT        HSIMT        \n",
                           "LdetritusN:     HSIMT        HSIMT        \n",
                           "SdetritusN:     HSIMT        HSIMT        \n",
                           "LdetritusC:     HSIMT        HSIMT        \n",
                           "SdetritusC:     HSIMT        HSIMT        \n",
                           "TIC:            HSIMT        HSIMT        \n",
                           "Talk:           HSIMT        HSIMT        \n",
                           "oxyg:           HSIMT        HSIMT" ;
    
  • The old tracer advection CPP options will be ignored in ROMS. Therefore, Users may want to still keep then in the application header file to run with older versions of ROMS.
  • We strongly discourage you from using the parabolic splines vertical advection (SPLINES) scheme in realistic applications! It was intended for use in idealized toy problems. It is still available for historical reasons and backward compatibility.

Change History (3)

comment:1 Changed 11 months ago by arango

  • Resolution set to Done
  • Status changed from new to closed

comment:2 Changed 11 months ago by arango

  • Description modified (diff)

comment:3 Changed 11 months ago by arango

  • Description modified (diff)
Note: See TracTickets for help on using tickets.