source: branches/2012/dev_MERCATOR_2012_rev3555/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn3d.F90 @ 3583

Last change on this file since 3583 was 3583, checked in by cbricaud, 8 years ago

add modification from dev_r3327_MERCATOR1_BDY branch in dev_MERCATOR_2012_rev3555 branch

File size: 11.8 KB
Line 
1MODULE bdydyn3d
2   !!======================================================================
3   !!                       ***  MODULE  bdydyn3d  ***
4   !! Unstructured Open Boundary Cond. :   Flow relaxation scheme on baroclinic velocities
5   !!======================================================================
6   !! History :  3.4  !  2011     (D. Storkey) new module as part of BDY rewrite
7   !!----------------------------------------------------------------------
8#if defined key_bdy 
9   !!----------------------------------------------------------------------
10   !!   'key_bdy' :                    Unstructured Open Boundary Condition
11   !!----------------------------------------------------------------------
12   !!   bdy_dyn3d        : apply open boundary conditions to baroclinic velocities
13   !!   bdy_dyn3d_frs    : apply Flow Relaxation Scheme
14   !!----------------------------------------------------------------------
15   USE timing          ! Timing
16   USE wrk_nemo        ! Memory Allocation
17   USE oce             ! ocean dynamics and tracers
18   USE dom_oce         ! ocean space and time domain
19   USE bdy_oce         ! ocean open boundary conditions
20   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
21   USE in_out_manager  !
22   Use phycst
23
24   IMPLICIT NONE
25   PRIVATE
26
27   PUBLIC   bdy_dyn3d     ! routine called by bdy_dyn
28   PUBLIC   bdy_dyn3d_dmp ! routine called by step
29
30   !! * Substitutions
31#  include "domzgr_substitute.h90"
32   !!----------------------------------------------------------------------
33   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
34   !! $Id: bdydyn.F90 2528 2010-12-27 17:33:53Z rblod $
35   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
36   !!----------------------------------------------------------------------
37CONTAINS
38
39   SUBROUTINE bdy_dyn3d( kt )
40      !!----------------------------------------------------------------------
41      !!                  ***  SUBROUTINE bdy_dyn3d  ***
42      !!
43      !! ** Purpose : - Apply open boundary conditions for baroclinic velocities
44      !!
45      !!----------------------------------------------------------------------
46      INTEGER, INTENT( in ) :: kt     ! Main time step counter
47      !!
48      INTEGER               :: ib_bdy ! loop index
49      !!
50
51      DO ib_bdy=1, nb_bdy
52
53!!$         IF ( using Orlanski radiation conditions ) THEN
54!!$            CALL bdy_rad( kt,  bdyidx(ib_bdy) )
55!!$         ENDIF
56
57         SELECT CASE( nn_dyn3d(ib_bdy) )
58         CASE(jp_none)
59            CYCLE
60         CASE(jp_frs)
61            CALL bdy_dyn3d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt )
62         CASE(2)
63            CALL bdy_dyn3d_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt )
64         CASE(3)
65            CALL bdy_dyn3d_zro( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt )
66         CASE DEFAULT
67            CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' )
68         END SELECT
69      ENDDO
70
71   END SUBROUTINE bdy_dyn3d
72
73   SUBROUTINE bdy_dyn3d_spe( idx, dta, kt )
74      !!----------------------------------------------------------------------
75      !!                  ***  SUBROUTINE bdy_dyn3d_spe  ***
76      !!
77      !! ** Purpose : - Apply a specified value for baroclinic velocities
78      !!                at open boundaries.
79      !!
80      !!----------------------------------------------------------------------
81      INTEGER                     ::   kt
82      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
83      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
84      !!
85      INTEGER  ::   jb, jk         ! dummy loop indices
86      INTEGER  ::   ii, ij, igrd   ! local integers
87      REAL(wp) ::   zwgt           ! boundary weight
88      !!----------------------------------------------------------------------
89      !
90      IF( nn_timing == 1 ) CALL timing_start('bdy_dyn3d_spe')
91      !
92      igrd = 2                      ! Relaxation of zonal 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            ua(ii,ij,jk) = dta%u3d(jb,jk) * umask(ii,ij,jk)
98         END DO
99      END DO
100      !
101      igrd = 3                      ! Relaxation of meridional velocity
102      DO jb = 1, idx%nblenrim(igrd)
103         DO jk = 1, jpkm1
104            ii   = idx%nbi(jb,igrd)
105            ij   = idx%nbj(jb,igrd)
106            va(ii,ij,jk) = dta%v3d(jb,jk) * vmask(ii,ij,jk)
107         END DO
108      END DO
109      CALL lbc_lnk( ua, 'U', -1. )   ;   CALL lbc_lnk( va, 'V', -1. )   ! Boundary points should be updated
110      !
111      IF( kt .eq. nit000 ) CLOSE( unit = 102 )
112
113      IF( nn_timing == 1 ) CALL timing_stop('bdy_dyn3d_spe')
114
115   END SUBROUTINE bdy_dyn3d_spe
116
117   SUBROUTINE bdy_dyn3d_zro( idx, dta, kt )
118      !!----------------------------------------------------------------------
119      !!                  ***  SUBROUTINE bdy_dyn3d_zro  ***
120      !!
121      !! ** Purpose : - baroclinic velocities = 0. at open boundaries.
122      !!
123      !!----------------------------------------------------------------------
124      INTEGER                     ::   kt
125      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
126      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
127      !!
128      INTEGER  ::   ib, ik         ! dummy loop indices
129      INTEGER  ::   ii, ij, igrd, zcoef   ! local integers
130      REAL(wp) ::   zwgt           ! boundary weight
131      !!----------------------------------------------------------------------
132      !
133      IF( nn_timing == 1 ) CALL timing_start('bdy_dyn3d_zro')
134      !
135      igrd = 2                       ! Everything is at T-points here
136      DO ib = 1, idx%nblenrim(igrd)
137         ii = idx%nbi(ib,igrd)
138         ij = idx%nbj(ib,igrd)
139         DO ik = 1, jpkm1
140            ua(ii,ij,ik) = 0._wp
141         END DO
142      END DO
143
144      igrd = 3                       ! Everything is at T-points here
145      DO ib = 1, idx%nblenrim(igrd)
146         ii = idx%nbi(ib,igrd)
147         ij = idx%nbj(ib,igrd)
148         DO ik = 1, jpkm1
149            va(ii,ij,ik) = 0._wp
150         END DO
151      END DO
152      !
153      CALL lbc_lnk( ua, 'U', -1. )   ;   CALL lbc_lnk( va, 'V', -1. )   ! Boundary points should be updated
154      !
155      IF( kt .eq. nit000 ) CLOSE( unit = 102 )
156
157      IF( nn_timing == 1 ) CALL timing_stop('bdy_dyn3d_zro')
158
159   END SUBROUTINE bdy_dyn3d_zro
160
161   SUBROUTINE bdy_dyn3d_frs( idx, dta, kt )
162      !!----------------------------------------------------------------------
163      !!                  ***  SUBROUTINE bdy_dyn3d_frs  ***
164      !!
165      !! ** Purpose : - Apply the Flow Relaxation Scheme for baroclinic velocities
166      !!                at open boundaries.
167      !!
168      !! References :- Engedahl H., 1995: Use of the flow relaxation scheme in
169      !!               a three-dimensional baroclinic ocean model with realistic
170      !!               topography. Tellus, 365-382.
171      !!----------------------------------------------------------------------
172      INTEGER                     ::   kt
173      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
174      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
175      !!
176      INTEGER  ::   jb, jk         ! dummy loop indices
177      INTEGER  ::   ii, ij, igrd   ! local integers
178      REAL(wp) ::   zwgt           ! boundary weight
179      !!----------------------------------------------------------------------
180      !
181      IF( nn_timing == 1 ) CALL timing_start('bdy_dyn3d_frs')
182      !
183      igrd = 2                      ! Relaxation of zonal velocity
184      DO jb = 1, idx%nblen(igrd)
185         DO jk = 1, jpkm1
186            ii   = idx%nbi(jb,igrd)
187            ij   = idx%nbj(jb,igrd)
188            zwgt = idx%nbw(jb,igrd)
189            ua(ii,ij,jk) = ( ua(ii,ij,jk) + zwgt * ( dta%u3d(jb,jk) - ua(ii,ij,jk) ) ) * umask(ii,ij,jk)
190         END DO
191      END DO
192      !
193      igrd = 3                      ! Relaxation of meridional velocity
194      DO jb = 1, idx%nblen(igrd)
195         DO jk = 1, jpkm1
196            ii   = idx%nbi(jb,igrd)
197            ij   = idx%nbj(jb,igrd)
198            zwgt = idx%nbw(jb,igrd)
199            va(ii,ij,jk) = ( va(ii,ij,jk) + zwgt * ( dta%v3d(jb,jk) - va(ii,ij,jk) ) ) * vmask(ii,ij,jk)
200         END DO
201      END DO
202      CALL lbc_lnk( ua, 'U', -1. )   ;   CALL lbc_lnk( va, 'V', -1. )   ! Boundary points should be updated
203      !
204      IF( kt .eq. nit000 ) CLOSE( unit = 102 )
205
206      IF( nn_timing == 1 ) CALL timing_stop('bdy_dyn3d_frs')
207
208   END SUBROUTINE bdy_dyn3d_frs
209
210   SUBROUTINE bdy_dyn3d_dmp( kt )
211      !!----------------------------------------------------------------------
212      !!                  ***  SUBROUTINE bdy_dyn3d_dmp  ***
213      !!
214      !! ** Purpose : Apply damping for baroclinic velocities at open boundaries.
215      !!
216      !!----------------------------------------------------------------------
217      INTEGER                     ::   kt
218      !!
219      INTEGER  ::   jb, jk         ! dummy loop indices
220      INTEGER  ::   ii, ij, igrd   ! local integers
221      REAL(wp) ::   zwgt           ! boundary weight
222      INTEGER  ::  ib_bdy          ! loop index
223      !!----------------------------------------------------------------------
224      !
225      IF( nn_timing == 1 ) CALL timing_start('bdy_dyn3d_dmp')
226      !
227      !-------------------------------------------------------
228      ! Remove barotropic part from before velocity
229      !-------------------------------------------------------
230      CALL wrk_alloc(jpi,jpj,pu2d,pv2d) 
231
232      pu2d(:,:) = 0.e0
233      pv2d(:,:) = 0.e0
234
235      DO jk = 1, jpkm1
236#if defined key_vvl
237         pu2d(:,:) = pu2d(:,:) + fse3u_b(:,:,jk)* ub(:,:,jk)   *umask(:,:,jk) 
238         pv2d(:,:) = pv2d(:,:) + fse3v_b(:,:,jk)* vb(:,:,jk)   *vmask(:,:,jk)
239#else
240         pu2d(:,:) = pu2d(:,:) + fse3u_0(:,:,jk) * ub(:,:,jk)  * umask(:,:,jk)
241         pv2d(:,:) = pv2d(:,:) + fse3v_0(:,:,jk) * vb(:,:,jk)  * vmask(:,:,jk)
242#endif
243      END DO
244
245      IF( lk_vvl ) THEN
246         pu2d(:,:) = pu2d(:,:) * umask(:,:,1) / ( hu_0(:,:) + sshu_b(:,:) + 1._wp - umask(:,:,1) )
247         pv2d(:,:) = pv2d(:,:) * vmask(:,:,1) / ( hv_0(:,:) + sshv_b(:,:) + 1._wp - vmask(:,:,1) )
248      ELSE
249         pu2d(:,:) = pv2d(:,:) * hur(:,:)
250         pv2d(:,:) = pu2d(:,:) * hvr(:,:)
251      ENDIF
252
253      DO ib_bdy=1, nb_bdy
254         IF ( ln_dyn3d_dmp(ib_bdy).and.nn_dyn3d(ib_bdy).gt.0 ) THEN
255            igrd = 2                      ! Relaxation of zonal velocity
256            DO jb = 1, idx_bdy(ib_bdy)%nblen(igrd)
257               ii   = idx_bdy(ib_bdy)%nbi(jb,igrd)
258               ij   = idx_bdy(ib_bdy)%nbj(jb,igrd)
259               zwgt = idx_bdy(ib_bdy)%nbd(jb,igrd)
260               DO jk = 1, jpkm1
261                  ua(ii,ij,jk) = ( ua(ii,ij,jk) + zwgt * ( dta_bdy(ib_bdy)%u3d(jb,jk) - ub(ii,ij,jk) + pu2d(ii,ij)) ) * umask(ii,ij,jk)
262               END DO
263            END DO
264            !
265            igrd = 3                      ! Relaxation of meridional velocity
266            DO jb = 1, idx_bdy(ib_bdy)%nblen(igrd)
267               ii   = idx_bdy(ib_bdy)%nbi(jb,igrd)
268               ij   = idx_bdy(ib_bdy)%nbj(jb,igrd)
269               zwgt = idx_bdy(ib_bdy)%nbd(jb,igrd)
270               DO jk = 1, jpkm1
271                  va(ii,ij,jk) = ( va(ii,ij,jk) + zwgt * ( dta_bdy(ib_bdy)%v3d(jb,jk) - vb(ii,ij,jk) + pv2d(ii,ij)) ) * vmask(ii,ij,jk)
272               END DO
273            END DO
274         ENDIF
275      ENDDO
276      !
277      CALL wrk_dealloc(jpi,jpj,pu2d,pv2d) 
278      !
279      CALL lbc_lnk( ua, 'U', -1. )   ;   CALL lbc_lnk( va, 'V', -1. )   ! Boundary points should be updated
280      !
281      IF( nn_timing == 1 ) CALL timing_stop('bdy_dyn3d_dmp')
282
283   END SUBROUTINE bdy_dyn3d_dmp
284
285#else
286   !!----------------------------------------------------------------------
287   !!   Dummy module                   NO Unstruct Open Boundary Conditions
288   !!----------------------------------------------------------------------
289CONTAINS
290   SUBROUTINE bdy_dyn3d( kt )      ! Empty routine
291      WRITE(*,*) 'bdy_dyn_frs: You should not have seen this print! error?', kt
292   END SUBROUTINE bdy_dyn3d
293#endif
294
295   !!======================================================================
296END MODULE bdydyn3d
Note: See TracBrowser for help on using the repository browser.