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.
bdydyn3d.F90 in NEMO/trunk/src/OCE/BDY – NEMO

source: NEMO/trunk/src/OCE/BDY/bdydyn3d.F90 @ 10529

Last change on this file since 10529 was 10529, checked in by smasson, 5 years ago

trunk: bugfix in BDY routines for non-MPP compilation, see #2211

  • Property svn:keywords set to Id
File size: 15.6 KB
RevLine 
[3117]1MODULE bdydyn3d
2   !!======================================================================
[3182]3   !!                       ***  MODULE  bdydyn3d  ***
[3191]4   !! Unstructured Open Boundary Cond. :   Flow relaxation scheme on baroclinic velocities
[3117]5   !!======================================================================
[3191]6   !! History :  3.4  !  2011     (D. Storkey) new module as part of BDY rewrite
[3680]7   !!            3.5  !  2012     (S. Mocavero, I. Epicoco) Optimization of BDY communications
[3117]8   !!----------------------------------------------------------------------
9   !!   bdy_dyn3d        : apply open boundary conditions to baroclinic velocities
10   !!   bdy_dyn3d_frs    : apply Flow Relaxation Scheme
11   !!----------------------------------------------------------------------
[3182]12   USE timing          ! Timing
[3117]13   USE oce             ! ocean dynamics and tracers
14   USE dom_oce         ! ocean space and time domain
15   USE bdy_oce         ! ocean open boundary conditions
[4292]16   USE bdylib          ! for orlanski library routines
[3117]17   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
18   USE in_out_manager  !
[10529]19   USE lib_mpp, ONLY: ctl_stop
[3651]20   Use phycst
[3117]21
22   IMPLICIT NONE
23   PRIVATE
24
[3191]25   PUBLIC   bdy_dyn3d     ! routine called by bdy_dyn
[3651]26   PUBLIC   bdy_dyn3d_dmp ! routine called by step
[3117]27
28   !!----------------------------------------------------------------------
[9598]29   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
[5215]30   !! $Id$
[10068]31   !! Software governed by the CeCILL license (see ./LICENSE)
[3117]32   !!----------------------------------------------------------------------
33CONTAINS
34
35   SUBROUTINE bdy_dyn3d( kt )
36      !!----------------------------------------------------------------------
37      !!                  ***  SUBROUTINE bdy_dyn3d  ***
38      !!
39      !! ** Purpose : - Apply open boundary conditions for baroclinic velocities
40      !!
41      !!----------------------------------------------------------------------
[6140]42      INTEGER, INTENT(in) ::   kt   ! Main time step counter
43      !
44      INTEGER ::   ib_bdy   ! loop index
45      !!----------------------------------------------------------------------
46      !
[3117]47      DO ib_bdy=1, nb_bdy
[6140]48         !
[4292]49         SELECT CASE( cn_dyn3d(ib_bdy) )
[6140]50         CASE('none')        ;   CYCLE
51         CASE('frs' )        ;   CALL bdy_dyn3d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
52         CASE('specified')   ;   CALL bdy_dyn3d_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
53         CASE('zero')        ;   CALL bdy_dyn3d_zro( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
54         CASE('orlanski' )   ;   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.false. )
55         CASE('orlanski_npo');   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.true. )
[7646]56         CASE('zerograd')    ;   CALL bdy_dyn3d_zgrad( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
57         CASE('neumann')     ;   CALL bdy_dyn3d_nmn( idx_bdy(ib_bdy), ib_bdy )
[6140]58         CASE DEFAULT        ;   CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' )
[3117]59         END SELECT
[6140]60      END DO
61      !
[3117]62   END SUBROUTINE bdy_dyn3d
63
[6140]64
[3680]65   SUBROUTINE bdy_dyn3d_spe( idx, dta, kt , ib_bdy )
[3651]66      !!----------------------------------------------------------------------
67      !!                  ***  SUBROUTINE bdy_dyn3d_spe  ***
68      !!
69      !! ** Purpose : - Apply a specified value for baroclinic velocities
70      !!                at open boundaries.
71      !!
72      !!----------------------------------------------------------------------
[6140]73      INTEGER        , INTENT(in) ::   kt      ! time step index
74      TYPE(OBC_INDEX), INTENT(in) ::   idx     ! OBC indices
75      TYPE(OBC_DATA) , INTENT(in) ::   dta     ! OBC external data
76      INTEGER        , INTENT(in) ::   ib_bdy  ! BDY set index
77      !
[3651]78      INTEGER  ::   jb, jk         ! dummy loop indices
79      INTEGER  ::   ii, ij, igrd   ! local integers
80      REAL(wp) ::   zwgt           ! boundary weight
81      !!----------------------------------------------------------------------
82      !
83      igrd = 2                      ! Relaxation of zonal velocity
84      DO jb = 1, idx%nblenrim(igrd)
85         DO jk = 1, jpkm1
86            ii   = idx%nbi(jb,igrd)
87            ij   = idx%nbj(jb,igrd)
88            ua(ii,ij,jk) = dta%u3d(jb,jk) * umask(ii,ij,jk)
89         END DO
90      END DO
91      !
92      igrd = 3                      ! Relaxation of meridional velocity
93      DO jb = 1, idx%nblenrim(igrd)
94         DO jk = 1, jpkm1
95            ii   = idx%nbi(jb,igrd)
96            ij   = idx%nbj(jb,igrd)
97            va(ii,ij,jk) = dta%v3d(jb,jk) * vmask(ii,ij,jk)
98         END DO
99      END DO
[10425]100      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )   ! Boundary points should be updated 
101      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )   
[3651]102      !
[6140]103      IF( kt == nit000 )   CLOSE( unit = 102 )
104      !
[3651]105   END SUBROUTINE bdy_dyn3d_spe
106
[9124]107
[7646]108   SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt , ib_bdy )
109      !!----------------------------------------------------------------------
110      !!                  ***  SUBROUTINE bdy_dyn3d_zgrad  ***
111      !!
112      !! ** Purpose : - Enforce a zero gradient of normal velocity
113      !!
114      !!----------------------------------------------------------------------
115      INTEGER                     ::   kt
116      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
117      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
118      INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index
119      !!
120      INTEGER  ::   jb, jk         ! dummy loop indices
121      INTEGER  ::   ii, ij, igrd   ! local integers
122      REAL(wp) ::   zwgt           ! boundary weight
123      INTEGER  ::   fu, fv
124      !!----------------------------------------------------------------------
125      !
126      igrd = 2                      ! Copying tangential velocity into bdy points
127      DO jb = 1, idx%nblenrim(igrd)
128         DO jk = 1, jpkm1
129            ii   = idx%nbi(jb,igrd)
130            ij   = idx%nbj(jb,igrd)
131            fu   = ABS( ABS (NINT( idx%flagu(jb,igrd) ) ) - 1 )
132            ua(ii,ij,jk) = ua(ii,ij,jk) * REAL( 1 - fu) + ( ua(ii,ij+fu,jk) * umask(ii,ij+fu,jk) &
133                        &+ ua(ii,ij-fu,jk) * umask(ii,ij-fu,jk) ) * umask(ii,ij,jk) * REAL( fu )
134         END DO
135      END DO
136      !
137      igrd = 3                      ! Copying tangential velocity into bdy points
138      DO jb = 1, idx%nblenrim(igrd)
139         DO jk = 1, jpkm1
140            ii   = idx%nbi(jb,igrd)
141            ij   = idx%nbj(jb,igrd)
142            fv   = ABS( ABS (NINT( idx%flagv(jb,igrd) ) ) - 1 )
143            va(ii,ij,jk) = va(ii,ij,jk) * REAL( 1 - fv ) + ( va(ii+fv,ij,jk) * vmask(ii+fv,ij,jk) &
144                        &+ va(ii-fv,ij,jk) * vmask(ii-fv,ij,jk) ) * vmask(ii,ij,jk) * REAL( fv )
145         END DO
146      END DO
[10425]147      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )   ! Boundary points should be updated 
148      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )   
[7646]149      !
[9124]150      IF( kt == nit000 )   CLOSE( unit = 102 )
151      !
152   END SUBROUTINE bdy_dyn3d_zgrad
[6140]153
[7646]154
[3680]155   SUBROUTINE bdy_dyn3d_zro( idx, dta, kt, ib_bdy )
[3651]156      !!----------------------------------------------------------------------
157      !!                  ***  SUBROUTINE bdy_dyn3d_zro  ***
158      !!
159      !! ** Purpose : - baroclinic velocities = 0. at open boundaries.
160      !!
161      !!----------------------------------------------------------------------
[6140]162      INTEGER        , INTENT(in) ::   kt      ! time step index
163      TYPE(OBC_INDEX), INTENT(in) ::   idx     ! OBC indices
164      TYPE(OBC_DATA) , INTENT(in) ::   dta     ! OBC external data
[3680]165      INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index
[6140]166      !
[3651]167      INTEGER  ::   ib, ik         ! dummy loop indices
[6140]168      INTEGER  ::   ii, ij, igrd   ! local integers
[3651]169      REAL(wp) ::   zwgt           ! boundary weight
170      !!----------------------------------------------------------------------
171      !
172      igrd = 2                       ! Everything is at T-points here
173      DO ib = 1, idx%nblenrim(igrd)
174         ii = idx%nbi(ib,igrd)
175         ij = idx%nbj(ib,igrd)
176         DO ik = 1, jpkm1
177            ua(ii,ij,ik) = 0._wp
178         END DO
179      END DO
180
181      igrd = 3                       ! Everything is at T-points here
182      DO ib = 1, idx%nblenrim(igrd)
183         ii = idx%nbi(ib,igrd)
184         ij = idx%nbj(ib,igrd)
185         DO ik = 1, jpkm1
186            va(ii,ij,ik) = 0._wp
187         END DO
188      END DO
189      !
[10425]190      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )   ;   CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1.,ib_bdy )   ! Boundary points should be updated
[3651]191      !
[6140]192      IF( kt == nit000 )   CLOSE( unit = 102 )
193      !
194   END SUBROUTINE bdy_dyn3d_zro
[3651]195
196
[3680]197   SUBROUTINE bdy_dyn3d_frs( idx, dta, kt, ib_bdy )
[3117]198      !!----------------------------------------------------------------------
199      !!                  ***  SUBROUTINE bdy_dyn3d_frs  ***
200      !!
201      !! ** Purpose : - Apply the Flow Relaxation Scheme for baroclinic velocities
202      !!                at open boundaries.
203      !!
204      !! References :- Engedahl H., 1995: Use of the flow relaxation scheme in
205      !!               a three-dimensional baroclinic ocean model with realistic
206      !!               topography. Tellus, 365-382.
207      !!----------------------------------------------------------------------
[6140]208      INTEGER        , INTENT(in) ::   kt      ! time step index
209      TYPE(OBC_INDEX), INTENT(in) ::   idx     ! OBC indices
210      TYPE(OBC_DATA) , INTENT(in) ::   dta     ! OBC external data
[3680]211      INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index
[6140]212      !
[3117]213      INTEGER  ::   jb, jk         ! dummy loop indices
214      INTEGER  ::   ii, ij, igrd   ! local integers
215      REAL(wp) ::   zwgt           ! boundary weight
216      !!----------------------------------------------------------------------
217      !
218      igrd = 2                      ! Relaxation of zonal velocity
219      DO jb = 1, idx%nblen(igrd)
220         DO jk = 1, jpkm1
221            ii   = idx%nbi(jb,igrd)
222            ij   = idx%nbj(jb,igrd)
223            zwgt = idx%nbw(jb,igrd)
224            ua(ii,ij,jk) = ( ua(ii,ij,jk) + zwgt * ( dta%u3d(jb,jk) - ua(ii,ij,jk) ) ) * umask(ii,ij,jk)
225         END DO
226      END DO
227      !
228      igrd = 3                      ! Relaxation of meridional velocity
229      DO jb = 1, idx%nblen(igrd)
230         DO jk = 1, jpkm1
231            ii   = idx%nbi(jb,igrd)
232            ij   = idx%nbj(jb,igrd)
233            zwgt = idx%nbw(jb,igrd)
234            va(ii,ij,jk) = ( va(ii,ij,jk) + zwgt * ( dta%v3d(jb,jk) - va(ii,ij,jk) ) ) * vmask(ii,ij,jk)
235         END DO
236      END DO
[10425]237      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )    ! Boundary points should be updated
238      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )   
[3117]239      !
[6140]240      IF( kt == nit000 )   CLOSE( unit = 102 )
241      !
242   END SUBROUTINE bdy_dyn3d_frs
[3117]243
[3182]244
[4292]245   SUBROUTINE bdy_dyn3d_orlanski( idx, dta, ib_bdy, ll_npo )
246      !!----------------------------------------------------------------------
247      !!                 ***  SUBROUTINE bdy_dyn3d_orlanski  ***
248      !!             
249      !!              - Apply Orlanski radiation to baroclinic velocities.
250      !!              - Wrapper routine for bdy_orlanski_3d
251      !!
252      !!
253      !! References:  Marchesiello, McWilliams and Shchepetkin, Ocean Modelling vol. 3 (2001)   
254      !!----------------------------------------------------------------------
255      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices
256      TYPE(OBC_DATA),               INTENT(in) ::   dta  ! OBC external data
257      INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index
258      LOGICAL,                      INTENT(in) ::   ll_npo  ! switch for NPO version
259
260      INTEGER  ::   jb, igrd                               ! dummy loop indices
261      !!----------------------------------------------------------------------
262      !
263      !! Note that at this stage the ub and ua arrays contain the baroclinic velocities.
264      !
265      igrd = 2      ! Orlanski bc on u-velocity;
266      !           
267      CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, ll_npo )
268
269      igrd = 3      ! Orlanski bc on v-velocity
270     
271      CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo )
272      !
[10425]273      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )    ! Boundary points should be updated
274      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )   
[4292]275      !
276   END SUBROUTINE bdy_dyn3d_orlanski
277
278
[3651]279   SUBROUTINE bdy_dyn3d_dmp( kt )
280      !!----------------------------------------------------------------------
281      !!                  ***  SUBROUTINE bdy_dyn3d_dmp  ***
282      !!
283      !! ** Purpose : Apply damping for baroclinic velocities at open boundaries.
284      !!
285      !!----------------------------------------------------------------------
[6140]286      INTEGER, INTENT(in) ::   kt   ! time step index
287      !
[3651]288      INTEGER  ::   jb, jk         ! dummy loop indices
[6140]289      INTEGER  ::   ib_bdy         ! loop index
[3651]290      INTEGER  ::   ii, ij, igrd   ! local integers
291      REAL(wp) ::   zwgt           ! boundary weight
292      !!----------------------------------------------------------------------
293      !
[9124]294      IF( ln_timing )   CALL timing_start('bdy_dyn3d_dmp')
[3651]295      !
296      DO ib_bdy=1, nb_bdy
[4292]297         IF ( ln_dyn3d_dmp(ib_bdy) .and. cn_dyn3d(ib_bdy) /= 'none' ) THEN
[3651]298            igrd = 2                      ! Relaxation of zonal velocity
299            DO jb = 1, idx_bdy(ib_bdy)%nblen(igrd)
300               ii   = idx_bdy(ib_bdy)%nbi(jb,igrd)
301               ij   = idx_bdy(ib_bdy)%nbj(jb,igrd)
302               zwgt = idx_bdy(ib_bdy)%nbd(jb,igrd)
303               DO jk = 1, jpkm1
[3703]304                  ua(ii,ij,jk) = ( ua(ii,ij,jk) + zwgt * ( dta_bdy(ib_bdy)%u3d(jb,jk) - &
[4354]305                                   ub(ii,ij,jk) + ub_b(ii,ij)) ) * umask(ii,ij,jk)
[3651]306               END DO
307            END DO
308            !
309            igrd = 3                      ! Relaxation of meridional velocity
310            DO jb = 1, idx_bdy(ib_bdy)%nblen(igrd)
311               ii   = idx_bdy(ib_bdy)%nbi(jb,igrd)
312               ij   = idx_bdy(ib_bdy)%nbj(jb,igrd)
313               zwgt = idx_bdy(ib_bdy)%nbd(jb,igrd)
314               DO jk = 1, jpkm1
[3703]315                  va(ii,ij,jk) = ( va(ii,ij,jk) + zwgt * ( dta_bdy(ib_bdy)%v3d(jb,jk) -  &
[4354]316                                   vb(ii,ij,jk) + vb_b(ii,ij)) ) * vmask(ii,ij,jk)
[3651]317               END DO
318            END DO
319         ENDIF
[6140]320      END DO
[3651]321      !
[10425]322      CALL lbc_lnk_multi( 'bdydyn3d', ua, 'U', -1.,  va, 'V', -1. )   ! Boundary points should be updated
[3651]323      !
[9124]324      IF( ln_timing )   CALL timing_stop('bdy_dyn3d_dmp')
[6140]325      !
[3651]326   END SUBROUTINE bdy_dyn3d_dmp
327
[9124]328
[7646]329   SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy )
330      !!----------------------------------------------------------------------
331      !!                 ***  SUBROUTINE bdy_dyn3d_nmn  ***
332      !!             
333      !!              - Apply Neumann condition to baroclinic velocities.
334      !!              - Wrapper routine for bdy_nmn
335      !!
336      !!
337      !!----------------------------------------------------------------------
338      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices
339      INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index
[3117]340
[7646]341      INTEGER  ::   jb, igrd                               ! dummy loop indices
342      !!----------------------------------------------------------------------
343      !
344      !! Note that at this stage the ub and ua arrays contain the baroclinic velocities.
345      !
346      igrd = 2      ! Neumann bc on u-velocity;
347      !           
348      CALL bdy_nmn( idx, igrd, ua )
349
350      igrd = 3      ! Neumann bc on v-velocity
351     
352      CALL bdy_nmn( idx, igrd, va )
353      !
[10425]354      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )    ! Boundary points should be updated
355      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )
[7646]356      !
357   END SUBROUTINE bdy_dyn3d_nmn
358
[3117]359   !!======================================================================
360END MODULE bdydyn3d
Note: See TracBrowser for help on using the repository browser.