New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 7058 for branches/2016/dev_r6522_SIMPLIF_3/NEMOGCM/NEMO/OPA_SRC/BDY/bdytra.F90 – NEMO

Ignore:
Timestamp:
2016-10-20T15:19:01+02:00 (8 years ago)
Author:
lovato
Message:

#1783 - trunk: Generalize the open boundary schemes and revise TRA and TRC BDY wrappers

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2016/dev_r6522_SIMPLIF_3/NEMOGCM/NEMO/OPA_SRC/BDY/bdytra.F90

    r6862 r7058  
    88   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge 
    99   !!            3.5  !  2012     (S. Mocavero, I. Epicoco) Optimization of BDY communications 
     10   !!            4.0  !  2016     (T. Lovato) Generalize OBC structure 
    1011   !!---------------------------------------------------------------------- 
    11    !!   bdy_tra            : Apply open boundary conditions to T and S 
    12    !!   bdy_tra_frs        : Apply Flow Relaxation Scheme 
     12   !!   bdy_tra       : Apply open boundary conditions & damping to T and S 
    1313   !!---------------------------------------------------------------------- 
    1414   USE oce            ! ocean dynamics and tracers variables 
     
    1616   USE bdy_oce        ! ocean open boundary conditions 
    1717   USE bdylib         ! for orlanski library routines 
    18    USE bdydta   , ONLY:   bf   !  
    1918   ! 
    2019   USE in_out_manager ! I/O manager 
     
    2524   PRIVATE 
    2625 
     26   ! Local structure to rearrange tracers data 
     27   TYPE, PUBLIC ::   ztrabdy 
     28      REAL(wp), POINTER, DIMENSION(:,:) ::  tra 
     29   END TYPE 
     30 
    2731   PUBLIC   bdy_tra      ! called in tranxt.F90  
    2832   PUBLIC   bdy_tra_dmp  ! called in step.F90  
    2933 
    3034   !!---------------------------------------------------------------------- 
    31    !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
     35   !! NEMO/OPA 4.0, NEMO Consortium (2016) 
    3236   !! $Id$  
    3337   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     
    4448      INTEGER, INTENT(in) ::   kt   ! Main time step counter 
    4549      ! 
    46       INTEGER ::   ib_bdy   ! Loop index 
     50      INTEGER                        :: ib_bdy, jn   ! Loop indeces 
     51      TYPE(ztrabdy), DIMENSION(jpts) :: zdta         ! Temporary data structure 
    4752      !!---------------------------------------------------------------------- 
    4853 
    4954      DO ib_bdy=1, nb_bdy 
    5055         ! 
    51          SELECT CASE( cn_tra(ib_bdy) ) 
    52          CASE('none'        )   ;   CYCLE 
    53          CASE('frs'         )   ;   CALL bdy_tra_frs     ( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt ) 
    54          CASE('specified'   )   ;   CALL bdy_tra_spe     ( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt ) 
    55          CASE('neumann'     )   ;   CALL bdy_tra_nmn     ( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt ) 
    56          CASE('orlanski'    )   ;   CALL bdy_tra_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ll_npo=.false. ) 
    57          CASE('orlanski_npo')   ;   CALL bdy_tra_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ll_npo=.true. ) 
    58          CASE('runoff'      )   ;   CALL bdy_tra_rnf     ( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt ) 
    59          CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' ) 
    60          END SELECT 
    61          ! Boundary points should be updated 
    62          CALL lbc_bdy_lnk( tsa(:,:,:,jp_tem), 'T', 1., ib_bdy ) 
    63          CALL lbc_bdy_lnk( tsa(:,:,:,jp_sal), 'T', 1., ib_bdy ) 
     56         zdta(1)%tra => dta_bdy(ib_bdy)%tem 
     57         zdta(2)%tra => dta_bdy(ib_bdy)%sal 
     58         ! 
     59         DO jn = 1, jpts 
     60            ! 
     61            SELECT CASE( cn_tra(ib_bdy) ) 
     62            CASE('none'        )   ;   CYCLE 
     63            CASE('frs'         )   ;   CALL bdy_frs ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra ) 
     64            CASE('specified'   )   ;   CALL bdy_spe ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra ) 
     65            CASE('neumann'     )   ;   CALL bdy_nmn ( idx_bdy(ib_bdy),                tsa(:,:,:,jn)               ) 
     66            CASE('orlanski'    )   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.false. ) 
     67            CASE('orlanski_npo')   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.true. ) 
     68            CASE('runoff'      )   ;   CALL bdy_rnf ( idx_bdy(ib_bdy),                tsa(:,:,:,jn),               jn ) 
     69            CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' ) 
     70            END SELECT 
     71            ! Boundary points should be updated 
     72            CALL lbc_bdy_lnk( tsa(:,:,:,jn), 'T', 1., ib_bdy ) 
     73            !  
     74         END DO 
    6475      END DO 
    6576      ! 
    6677   END SUBROUTINE bdy_tra 
    6778 
    68  
    69    SUBROUTINE bdy_tra_frs( idx, dta, kt ) 
     79   SUBROUTINE bdy_rnf( idx, pta, jpa ) 
    7080      !!---------------------------------------------------------------------- 
    71       !!                 ***  SUBROUTINE bdy_tra_frs  *** 
     81      !!                 ***  SUBROUTINE bdy_rnf  *** 
    7282      !!                     
    73       !! ** Purpose : Apply the Flow Relaxation Scheme for tracers at open boundaries. 
    74       !!  
    75       !! Reference : Engedahl H., 1995, Tellus, 365-382. 
    76       !!---------------------------------------------------------------------- 
    77       INTEGER,         INTENT(in) ::   kt    ! 
    78       TYPE(OBC_INDEX), INTENT(in) ::   idx   ! OBC indices 
    79       TYPE(OBC_DATA),  INTENT(in) ::   dta   ! OBC external data 
    80       ! 
    81       REAL(wp) ::   zwgt           ! boundary weight 
    82       INTEGER  ::   ib, ik, igrd   ! dummy loop indices 
    83       INTEGER  ::   ii, ij         ! 2D addresses 
    84       !!---------------------------------------------------------------------- 
    85       ! 
    86       IF( nn_timing == 1 )   CALL timing_start('bdy_tra_frs') 
    87       ! 
    88       igrd = 1                       ! Everything is at T-points here 
    89       DO ib = 1, idx%nblen(igrd) 
    90          DO ik = 1, jpkm1 
    91             ii = idx%nbi(ib,igrd) 
    92             ij = idx%nbj(ib,igrd) 
    93             zwgt = idx%nbw(ib,igrd) 
    94             tsa(ii,ij,ik,jp_tem) = ( tsa(ii,ij,ik,jp_tem) + zwgt * ( dta%tem(ib,ik) - tsa(ii,ij,ik,jp_tem) ) ) * tmask(ii,ij,ik)          
    95             tsa(ii,ij,ik,jp_sal) = ( tsa(ii,ij,ik,jp_sal) + zwgt * ( dta%sal(ib,ik) - tsa(ii,ij,ik,jp_sal) ) ) * tmask(ii,ij,ik) 
    96          END DO 
    97       END DO  
    98       ! 
    99       IF( kt .eq. nit000 )   CLOSE( unit = 102 ) 
    100       ! 
    101       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_frs') 
    102       ! 
    103    END SUBROUTINE bdy_tra_frs 
    104  
    105  
    106    SUBROUTINE bdy_tra_spe( idx, dta, kt ) 
    107       !!---------------------------------------------------------------------- 
    108       !!                 ***  SUBROUTINE bdy_tra_frs  *** 
    109       !!                     
    110       !! ** Purpose : Apply a specified value for tracers at open boundaries. 
     83      !! ** Purpose : Specialized routine to apply TRA runoff values at OBs: 
     84      !!                  - duplicate the neighbour value for the temperature 
     85      !!                  - specified to 0.1 PSU for the salinity 
    11186      !!  
    11287      !!---------------------------------------------------------------------- 
    113       INTEGER,         INTENT(in) ::   kt    ! 
    114       TYPE(OBC_INDEX), INTENT(in) ::   idx   ! OBC indices 
    115       TYPE(OBC_DATA),  INTENT(in) ::   dta   ! OBC external data 
    116       ! 
    117       REAL(wp) ::   zwgt           ! boundary weight 
    118       INTEGER  ::   ib, ik, igrd   ! dummy loop indices 
    119       INTEGER  ::   ii, ij         ! 2D addresses 
    120       !!---------------------------------------------------------------------- 
    121       ! 
    122       IF( nn_timing == 1 ) CALL timing_start('bdy_tra_spe') 
    123       ! 
    124       igrd = 1                       ! Everything is at T-points here 
    125       DO ib = 1, idx%nblenrim(igrd) 
    126          ii = idx%nbi(ib,igrd) 
    127          ij = idx%nbj(ib,igrd) 
    128          DO ik = 1, jpkm1 
    129             tsa(ii,ij,ik,jp_tem) = dta%tem(ib,ik) * tmask(ii,ij,ik) 
    130             tsa(ii,ij,ik,jp_sal) = dta%sal(ib,ik) * tmask(ii,ij,ik) 
    131          END DO 
    132       END DO 
    133       ! 
    134       IF( kt == nit000 )   CLOSE( unit = 102 ) 
    135       ! 
    136       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_spe') 
    137       ! 
    138    END SUBROUTINE bdy_tra_spe 
    139  
    140  
    141    SUBROUTINE bdy_tra_nmn( idx, dta, kt ) 
    142       !!---------------------------------------------------------------------- 
    143       !!                 ***  SUBROUTINE bdy_tra_nmn  *** 
    144       !!                     
    145       !! ** Purpose : Duplicate the value for tracers at open boundaries. 
    146       !!  
    147       !!---------------------------------------------------------------------- 
    148       INTEGER,         INTENT(in) ::   kt    !  
    149       TYPE(OBC_INDEX), INTENT(in) ::   idx   ! OBC indices 
    150       TYPE(OBC_DATA) , INTENT(in) ::   dta   ! OBC external data 
    151       ! 
    152       REAL(wp) ::   zwgt           ! boundary weight 
    153       INTEGER  ::   ib, ik, igrd   ! dummy loop indices 
    154       INTEGER  ::   ii, ij,zcoef, zcoef1,zcoef2, ip, jp   ! 2D addresses 
    155       !!---------------------------------------------------------------------- 
    156       ! 
    157       IF( nn_timing == 1 )   CALL timing_start('bdy_tra_nmn') 
    158       ! 
    159       igrd = 1                       ! Everything is at T-points here 
    160       DO ib = 1, idx%nblenrim(igrd) 
    161          ii = idx%nbi(ib,igrd) 
    162          ij = idx%nbj(ib,igrd) 
    163          DO ik = 1, jpkm1 
    164             ! search the sense of the gradient 
    165             zcoef1 = bdytmask(ii-1,ij  ) +  bdytmask(ii+1,ij  ) 
    166             zcoef2 = bdytmask(ii  ,ij-1) +  bdytmask(ii  ,ij+1) 
    167             IF ( zcoef1+zcoef2 == 0) THEN 
    168                ! corner 
    169                zcoef = tmask(ii-1,ij,ik) + tmask(ii+1,ij,ik) +  tmask(ii,ij-1,ik) +  tmask(ii,ij+1,ik) 
    170                tsa(ii,ij,ik,jp_tem) = tsa(ii-1,ij  ,ik,jp_tem) * tmask(ii-1,ij  ,ik) + & 
    171                  &                    tsa(ii+1,ij  ,ik,jp_tem) * tmask(ii+1,ij  ,ik) + & 
    172                  &                    tsa(ii  ,ij-1,ik,jp_tem) * tmask(ii  ,ij-1,ik) + & 
    173                  &                    tsa(ii  ,ij+1,ik,jp_tem) * tmask(ii  ,ij+1,ik) 
    174                tsa(ii,ij,ik,jp_tem) = ( tsa(ii,ij,ik,jp_tem) / MAX( 1, zcoef) ) * tmask(ii,ij,ik) 
    175                tsa(ii,ij,ik,jp_sal) = tsa(ii-1,ij  ,ik,jp_sal) * tmask(ii-1,ij  ,ik) + & 
    176                  &                    tsa(ii+1,ij  ,ik,jp_sal) * tmask(ii+1,ij  ,ik) + & 
    177                  &                    tsa(ii  ,ij-1,ik,jp_sal) * tmask(ii  ,ij-1,ik) + & 
    178                  &                    tsa(ii  ,ij+1,ik,jp_sal) * tmask(ii  ,ij+1,ik) 
    179                tsa(ii,ij,ik,jp_sal) = ( tsa(ii,ij,ik,jp_sal) / MAX( 1, zcoef) ) * tmask(ii,ij,ik) 
    180             ELSE 
    181                ip = bdytmask(ii+1,ij  ) - bdytmask(ii-1,ij  ) 
    182                jp = bdytmask(ii  ,ij+1) - bdytmask(ii  ,ij-1) 
    183                tsa(ii,ij,ik,jp_tem) = tsa(ii+ip,ij+jp,ik,jp_tem) * tmask(ii+ip,ij+jp,ik) 
    184                tsa(ii,ij,ik,jp_sal) = tsa(ii+ip,ij+jp,ik,jp_sal) * tmask(ii+ip,ij+jp,ik) 
    185             ENDIF 
    186          END DO 
    187       END DO 
    188       ! 
    189       IF( kt == nit000 )   CLOSE( unit = 102 ) 
    190       ! 
    191       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_nmn') 
    192       ! 
    193    END SUBROUTINE bdy_tra_nmn 
    194   
    195  
    196    SUBROUTINE bdy_tra_orlanski( idx, dta, ll_npo ) 
    197       !!---------------------------------------------------------------------- 
    198       !!                 ***  SUBROUTINE bdy_tra_orlanski  *** 
    199       !!              
    200       !!              - Apply Orlanski radiation to temperature and salinity.  
    201       !!              - Wrapper routine for bdy_orlanski_3d 
    202       !!  
    203       !! 
    204       !! References:  Marchesiello, McWilliams and Shchepetkin, Ocean Modelling vol. 3 (2001)     
    205       !!---------------------------------------------------------------------- 
    206       TYPE(OBC_INDEX), INTENT(in) ::   idx     ! OBC indices 
    207       TYPE(OBC_DATA) , INTENT(in) ::   dta     ! OBC external data 
    208       LOGICAL        , INTENT(in) ::   ll_npo  ! switch for NPO version 
    209       ! 
    210       INTEGER  ::   igrd                                    ! grid index 
    211       !!---------------------------------------------------------------------- 
    212       ! 
    213       IF( nn_timing == 1 ) CALL timing_start('bdy_tra_orlanski') 
    214       ! 
    215       igrd = 1      ! Orlanski bc on temperature;  
    216       !             
    217       CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_tem), tsa(:,:,:,jp_tem), dta%tem, ll_npo ) 
    218  
    219       igrd = 1      ! Orlanski bc on salinity; 
    220       !   
    221       CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_sal), tsa(:,:,:,jp_sal), dta%sal, ll_npo ) 
    222       ! 
    223       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_orlanski') 
    224       ! 
    225    END SUBROUTINE bdy_tra_orlanski 
    226  
    227  
    228    SUBROUTINE bdy_tra_rnf( idx, dta, kt ) 
    229       !!---------------------------------------------------------------------- 
    230       !!                 ***  SUBROUTINE bdy_tra_rnf  *** 
    231       !!                     
    232       !! ** Purpose : Apply the runoff values for tracers at open boundaries: 
    233       !!                  - specified to 0.1 PSU for the salinity 
    234       !!                  - duplicate the value for the temperature 
    235       !!  
    236       !!---------------------------------------------------------------------- 
    237       INTEGER        , INTENT(in) ::   kt    !  
    238       TYPE(OBC_INDEX), INTENT(in) ::   idx   ! OBC indices 
    239       TYPE(OBC_DATA) , INTENT(in) ::   dta   ! OBC external data 
     88      TYPE(OBC_INDEX),                     INTENT(in) ::   idx  ! OBC indices 
     89      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend 
     90      INTEGER,                             INTENT(in) ::   jpa  ! TRA index 
    24091      ! 
    24192      REAL(wp) ::   zwgt           ! boundary weight 
     
    24495      !!---------------------------------------------------------------------- 
    24596      ! 
    246       IF( nn_timing == 1 )   CALL timing_start('bdy_tra_rnf') 
     97      IF( nn_timing == 1 )   CALL timing_start('bdy_rnf') 
    24798      ! 
    24899      igrd = 1                       ! Everything is at T-points here 
     
    253104            ip = bdytmask(ii+1,ij  ) - bdytmask(ii-1,ij  ) 
    254105            jp = bdytmask(ii  ,ij+1) - bdytmask(ii  ,ij-1) 
    255             tsa(ii,ij,ik,jp_tem) = tsa(ii+ip,ij+jp,ik,jp_tem) * tmask(ii,ij,ik) 
    256             tsa(ii,ij,ik,jp_sal) =                        0.1 * tmask(ii,ij,ik) 
     106            if (jpa == jp_tem) pta(ii,ij,ik) = pta(ii+ip,ij+jp,ik) * tmask(ii,ij,ik) 
     107            if (jpa == jp_sal) pta(ii,ij,ik) =                 0.1 * tmask(ii,ij,ik) 
    257108         END DO 
    258109      END DO 
    259110      ! 
    260       IF( kt == nit000 )   CLOSE( unit = 102 ) 
     111      IF( nn_timing == 1 )   CALL timing_stop('bdy_rnf') 
    261112      ! 
    262       IF( nn_timing == 1 )   CALL timing_stop('bdy_tra_rnf') 
    263       ! 
    264    END SUBROUTINE bdy_tra_rnf 
    265  
     113   END SUBROUTINE bdy_rnf 
    266114 
    267115   SUBROUTINE bdy_tra_dmp( kt ) 
Note: See TracChangeset for help on using the changeset viewer.