source: branches/UKMO/r6232_collate_bgc_diagnostics/NEMOGCM/NEMO/OPA_SRC/BDY/bdytra.F90 @ 11134

Last change on this file since 11134 was 11134, checked in by jcastill, 17 months ago

Full set of changes as in the original branch

File size: 14.6 KB
Line 
1MODULE bdytra
2   !!======================================================================
3   !!                       ***  MODULE  bdytra  ***
4   !! Ocean tracers:   Apply boundary conditions for tracers
5   !!======================================================================
6   !! History :  1.0  !  2005-01  (J. Chanut, A. Sellar)  Original code
7   !!            3.0  !  2008-04  (NEMO team)  add in the reference version
8   !!            3.4  !  2011     (D. Storkey) rewrite in preparation for OBC-BDY merge
9   !!            3.5  !  2012     (S. Mocavero, I. Epicoco) Optimization of BDY communications
10   !!----------------------------------------------------------------------
11#if defined key_bdy
12   !!----------------------------------------------------------------------
13   !!   'key_bdy'                     Unstructured Open Boundary Conditions
14   !!----------------------------------------------------------------------
15   !!   bdy_tra            : Apply open boundary conditions to T and S
16   !!   bdy_tra_frs        : Apply Flow Relaxation Scheme
17   !!----------------------------------------------------------------------
18   USE timing          ! Timing
19   USE oce             ! ocean dynamics and tracers variables
20   USE dom_oce         ! ocean space and time domain variables
21   USE bdy_oce         ! ocean open boundary conditions
22   USE bdylib          ! for orlanski library routines
23   USE bdydta, ONLY:   bf
24   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
25   USE in_out_manager  ! I/O manager
26
27
28   IMPLICIT NONE
29   PRIVATE
30
31   PUBLIC bdy_tra      ! routine called in tranxt.F90
32   PUBLIC bdy_tra_dmp  ! routine called in step.F90
33
34   !!----------------------------------------------------------------------
35   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
36   !! $Id$
37   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
38   !!----------------------------------------------------------------------
39CONTAINS
40
41   SUBROUTINE bdy_tra( kt )
42      !!----------------------------------------------------------------------
43      !!                  ***  SUBROUTINE bdy_tra  ***
44      !!
45      !! ** Purpose : - Apply open boundary conditions for temperature and salinity
46      !!
47      !!----------------------------------------------------------------------
48      INTEGER, INTENT( in ) :: kt     ! Main time step counter
49      !!
50      INTEGER               :: ib_bdy ! Loop index
51
52      DO ib_bdy=1, nb_bdy
53
54         SELECT CASE( cn_tra(ib_bdy) )
55         CASE('none')
56            CYCLE
57         CASE('frs')
58            CALL bdy_tra_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt )
59         CASE('specified')
60            CALL bdy_tra_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt )
61         CASE('neumann')
62            CALL bdy_tra_nmn( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt )
63         CASE('orlanski')
64            CALL bdy_tra_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ll_npo=.false. )
65         CASE('orlanski_npo')
66            CALL bdy_tra_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ll_npo=.true. )
67         CASE('runoff')
68            CALL bdy_tra_rnf( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt )
69         CASE DEFAULT
70            CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' )
71         END SELECT
72         ! Boundary points should be updated
73         CALL lbc_bdy_lnk( tsa(:,:,:,jp_tem), 'T', 1., ib_bdy )
74         CALL lbc_bdy_lnk( tsa(:,:,:,jp_sal), 'T', 1., ib_bdy )
75      ENDDO
76      !
77
78   END SUBROUTINE bdy_tra
79
80   SUBROUTINE bdy_tra_frs( idx, dta, kt )
81      !!----------------------------------------------------------------------
82      !!                 ***  SUBROUTINE bdy_tra_frs  ***
83      !!                   
84      !! ** Purpose : Apply the Flow Relaxation Scheme for tracers at open boundaries.
85      !!
86      !! Reference : Engedahl H., 1995, Tellus, 365-382.
87      !!----------------------------------------------------------------------
88      INTEGER,         INTENT(in) ::   kt
89      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
90      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
91      !!
92      REAL(wp) ::   zwgt           ! boundary weight
93      REAL(wp) ::   zcoef, zcoef1,zcoef2 
94      INTEGER  ::   ib, ik, igrd   ! dummy loop indices
95      INTEGER  ::   ii, ij, ip, jp   ! 2D addresses
96      !!----------------------------------------------------------------------
97      !
98      IF( nn_timing == 1 ) CALL timing_start('bdy_tra_frs')
99      !
100      igrd = 1                       ! Everything is at T-points here
101      DO ib = 1, idx%nblen(igrd)
102         DO ik = 1, jpkm1
103            ii = idx%nbi(ib,igrd)
104            ij = idx%nbj(ib,igrd)
105            zwgt = idx%nbw(ib,igrd)
106            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)         
107            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)
108         END DO
109      END DO 
110      !
111      IF( kt .eq. nit000 ) CLOSE( unit = 102 )
112      !
113      IF( nn_timing == 1 ) CALL timing_stop('bdy_tra_frs')
114      !
115   END SUBROUTINE bdy_tra_frs
116 
117   SUBROUTINE bdy_tra_spe( idx, dta, kt )
118      !!----------------------------------------------------------------------
119      !!                 ***  SUBROUTINE bdy_tra_frs  ***
120      !!                   
121      !! ** Purpose : Apply a specified value for tracers at open boundaries.
122      !!
123      !!----------------------------------------------------------------------
124      INTEGER,         INTENT(in) ::   kt
125      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
126      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
127      !!
128      REAL(wp) ::   zwgt           ! boundary weight
129      INTEGER  ::   ib, ik, igrd   ! dummy loop indices
130      INTEGER  ::   ii, ij         ! 2D addresses
131      !!----------------------------------------------------------------------
132      !
133      IF( nn_timing == 1 ) CALL timing_start('bdy_tra_spe')
134      !
135      igrd = 1                       ! 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            tsa(ii,ij,ik,jp_tem) = dta%tem(ib,ik) * tmask(ii,ij,ik)
141            tsa(ii,ij,ik,jp_sal) = dta%sal(ib,ik) * tmask(ii,ij,ik)
142         END DO
143      END DO
144      !
145      IF( kt .eq. nit000 ) CLOSE( unit = 102 )
146      !
147      IF( nn_timing == 1 ) CALL timing_stop('bdy_tra_spe')
148      !
149   END SUBROUTINE bdy_tra_spe
150
151   SUBROUTINE bdy_tra_nmn( idx, dta, kt )
152      !!----------------------------------------------------------------------
153      !!                 ***  SUBROUTINE bdy_tra_nmn  ***
154      !!                   
155      !! ** Purpose : Duplicate the value for tracers at open boundaries.
156      !!
157      !!----------------------------------------------------------------------
158      INTEGER,         INTENT(in) ::   kt
159      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
160      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
161      !!
162      REAL(wp) ::   zwgt           ! boundary weight
163      REAL(wp) ::   zcoef, zcoef1,zcoef2 
164      INTEGER  ::   ib, ik, igrd   ! dummy loop indices
165      INTEGER  ::   ii, ij, ip, jp   ! 2D addresses
166      !!----------------------------------------------------------------------
167      !
168      IF( nn_timing == 1 ) CALL timing_start('bdy_tra_nmn')
169      !
170      igrd = 1                       ! Everything is at T-points here
171      DO ib = 1, idx%nblenrim(igrd)
172         ii = idx%nbi(ib,igrd)
173         ij = idx%nbj(ib,igrd)
174         DO ik = 1, jpkm1
175            ! search the sense of the gradient
176            zcoef1 = bdytmask(ii-1,ij  ) +  bdytmask(ii+1,ij  )
177            zcoef2 = bdytmask(ii  ,ij-1) +  bdytmask(ii  ,ij+1)
178            IF ( NINT(zcoef1+zcoef2) == 0) THEN 
179               ! corner
180               zcoef = tmask(ii-1,ij,ik) + tmask(ii+1,ij,ik) +  tmask(ii,ij-1,ik) +  tmask(ii,ij+1,ik)
181               tsa(ii,ij,ik,jp_tem) = tsa(ii-1,ij  ,ik,jp_tem) * tmask(ii-1,ij  ,ik) + &
182                 &                    tsa(ii+1,ij  ,ik,jp_tem) * tmask(ii+1,ij  ,ik) + &
183                 &                    tsa(ii  ,ij-1,ik,jp_tem) * tmask(ii  ,ij-1,ik) + &
184                 &                    tsa(ii  ,ij+1,ik,jp_tem) * tmask(ii  ,ij+1,ik)
185               tsa(ii,ij,ik,jp_tem) = ( tsa(ii,ij,ik,jp_tem) / MAX( 1._wp, zcoef) ) * tmask(ii,ij,ik) 
186               tsa(ii,ij,ik,jp_sal) = tsa(ii-1,ij  ,ik,jp_sal) * tmask(ii-1,ij  ,ik) + &
187                 &                    tsa(ii+1,ij  ,ik,jp_sal) * tmask(ii+1,ij  ,ik) + &
188                 &                    tsa(ii  ,ij-1,ik,jp_sal) * tmask(ii  ,ij-1,ik) + &
189                 &                    tsa(ii  ,ij+1,ik,jp_sal) * tmask(ii  ,ij+1,ik)
190               tsa(ii,ij,ik,jp_sal) = ( tsa(ii,ij,ik,jp_sal) / MAX( 1._wp, zcoef) ) * tmask(ii,ij,ik)
191            ELSE
192               ip = NINT(bdytmask(ii+1,ij  ) - bdytmask(ii-1,ij  )) 
193               jp = NINT(bdytmask(ii  ,ij+1) - bdytmask(ii  ,ij-1)) 
194               tsa(ii,ij,ik,jp_tem) = tsa(ii+ip,ij+jp,ik,jp_tem) * tmask(ii+ip,ij+jp,ik)
195               tsa(ii,ij,ik,jp_sal) = tsa(ii+ip,ij+jp,ik,jp_sal) * tmask(ii+ip,ij+jp,ik)
196            ENDIF
197         END DO
198      END DO
199      !
200      IF( kt .eq. nit000 ) CLOSE( unit = 102 )
201      !
202      IF( nn_timing == 1 ) CALL timing_stop('bdy_tra_nmn')
203      !
204   END SUBROUTINE bdy_tra_nmn
205 
206
207   SUBROUTINE bdy_tra_orlanski( idx, dta, ll_npo )
208      !!----------------------------------------------------------------------
209      !!                 ***  SUBROUTINE bdy_tra_orlanski  ***
210      !!             
211      !!              - Apply Orlanski radiation to temperature and salinity.
212      !!              - Wrapper routine for bdy_orlanski_3d
213      !!
214      !!
215      !! References:  Marchesiello, McWilliams and Shchepetkin, Ocean Modelling vol. 3 (2001)   
216      !!----------------------------------------------------------------------
217      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices
218      TYPE(OBC_DATA),               INTENT(in) ::   dta  ! OBC external data
219      LOGICAL,                      INTENT(in) ::   ll_npo  ! switch for NPO version
220
221      INTEGER  ::   igrd                                    ! grid index
222      !!----------------------------------------------------------------------
223
224      IF( nn_timing == 1 ) CALL timing_start('bdy_tra_orlanski')
225      !
226      igrd = 1      ! Orlanski bc on temperature;
227      !           
228      CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_tem), tsa(:,:,:,jp_tem), dta%tem, ll_npo )
229
230      igrd = 1      ! Orlanski bc on salinity;
231     
232      CALL bdy_orlanski_3d( idx, igrd, tsb(:,:,:,jp_sal), tsa(:,:,:,jp_sal), dta%sal, ll_npo )
233      !
234      IF( nn_timing == 1 ) CALL timing_stop('bdy_tra_orlanski')
235      !
236
237   END SUBROUTINE bdy_tra_orlanski
238
239
240   SUBROUTINE bdy_tra_rnf( idx, dta, kt )
241      !!----------------------------------------------------------------------
242      !!                 ***  SUBROUTINE bdy_tra_rnf  ***
243      !!                   
244      !! ** Purpose : Apply the runoff values for tracers at open boundaries:
245      !!                  - specified to 0.1 PSU for the salinity
246      !!                  - duplicate the value for the temperature
247      !!
248      !!----------------------------------------------------------------------
249      INTEGER,         INTENT(in) ::   kt
250      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
251      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
252      !!
253      REAL(wp) ::   zwgt           ! boundary weight
254      INTEGER  ::   ib, ik, igrd   ! dummy loop indices
255      INTEGER  ::   ii, ij, ip, jp ! 2D addresses
256      !!----------------------------------------------------------------------
257      !
258      IF( nn_timing == 1 ) CALL timing_start('bdy_tra_rnf')
259      !
260      igrd = 1                       ! Everything is at T-points here
261      DO ib = 1, idx%nblenrim(igrd)
262         ii = idx%nbi(ib,igrd)
263         ij = idx%nbj(ib,igrd)
264         DO ik = 1, jpkm1
265            ip = bdytmask(ii+1,ij  ) - bdytmask(ii-1,ij  )
266            jp = bdytmask(ii  ,ij+1) - bdytmask(ii  ,ij-1)
267            tsa(ii,ij,ik,jp_tem) = tsa(ii+ip,ij+jp,ik,jp_tem) * tmask(ii,ij,ik)
268            tsa(ii,ij,ik,jp_sal) =                        0.1 * tmask(ii,ij,ik)
269         END DO
270      END DO
271      !
272      IF( kt .eq. nit000 ) CLOSE( unit = 102 )
273      !
274      IF( nn_timing == 1 ) CALL timing_stop('bdy_tra_rnf')
275      !
276   END SUBROUTINE bdy_tra_rnf
277
278   SUBROUTINE bdy_tra_dmp( kt )
279      !!----------------------------------------------------------------------
280      !!                 ***  SUBROUTINE bdy_tra_dmp  ***
281      !!                   
282      !! ** Purpose : Apply damping for tracers at open boundaries.
283      !!
284      !!----------------------------------------------------------------------
285      INTEGER,         INTENT(in) ::   kt
286      !!
287      REAL(wp) ::   zwgt           ! boundary weight
288      REAL(wp) ::   zta, zsa, ztime
289      INTEGER  ::   ib, ik, igrd   ! dummy loop indices
290      INTEGER  ::   ii, ij         ! 2D addresses
291      INTEGER  ::   ib_bdy         ! Loop index
292      !!----------------------------------------------------------------------
293      !
294      IF( nn_timing == 1 ) CALL timing_start('bdy_tra_dmp')
295      !
296      DO ib_bdy=1, nb_bdy
297         IF ( ln_tra_dmp(ib_bdy) ) THEN
298            igrd = 1                       ! Everything is at T-points here
299            DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
300               ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
301               ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
302               zwgt = idx_bdy(ib_bdy)%nbd(ib,igrd)
303               DO ik = 1, jpkm1
304                  zta = zwgt * ( dta_bdy(ib_bdy)%tem(ib,ik) - tsb(ii,ij,ik,jp_tem) ) * tmask(ii,ij,ik)
305                  zsa = zwgt * ( dta_bdy(ib_bdy)%sal(ib,ik) - tsb(ii,ij,ik,jp_sal) ) * tmask(ii,ij,ik)
306                  tsa(ii,ij,ik,jp_tem) = tsa(ii,ij,ik,jp_tem) + zta
307                  tsa(ii,ij,ik,jp_sal) = tsa(ii,ij,ik,jp_sal) + zsa
308               END DO
309            END DO
310         ENDIF
311      ENDDO
312      !
313      IF( nn_timing == 1 ) CALL timing_stop('bdy_tra_dmp')
314      !
315   END SUBROUTINE bdy_tra_dmp
316 
317#else
318   !!----------------------------------------------------------------------
319   !!   Dummy module                   NO Unstruct Open Boundary Conditions
320   !!----------------------------------------------------------------------
321CONTAINS
322   SUBROUTINE bdy_tra(kt)      ! Empty routine
323      WRITE(*,*) 'bdy_tra: You should not have seen this print! error?', kt
324   END SUBROUTINE bdy_tra
325
326   SUBROUTINE bdy_tra_dmp(kt)      ! Empty routine
327      WRITE(*,*) 'bdy_tra_dmp: You should not have seen this print! error?', kt
328   END SUBROUTINE bdy_tra_dmp
329
330#endif
331
332   !!======================================================================
333END MODULE bdytra
Note: See TracBrowser for help on using the repository browser.