source: trunk/NEMOGCM/NEMO/OPA_SRC/USR/usrdef_closea.F90 @ 7646

Last change on this file since 7646 was 7427, checked in by timgraham, 5 years ago

#1811 Another lk_mpp_rep

  • Property svn:keywords set to Id
File size: 18.4 KB
Line 
1MODULE usrdef_closea
2   !!======================================================================
3   !!                   ***  MODULE  usrdef_closea  ***
4   !!
5   !!                      ===  ORCA configuration  ===
6   !!                         (2, 1 and 1/4 degrees)
7   !!
8   !! User define : specific treatments associated with closed seas
9   !!======================================================================
10   !! History :   8.2  !  2000-05  (O. Marti)  Original code
11   !!   NEMO      1.0  !  2002-06  (E. Durand, G. Madec)  F90
12   !!             3.0  !  2006-07  (G. Madec)  add clo_rnf, clo_ups, clo_bat
13   !!             3.4  !  2014-12  (P.G. Fogli) sbc_clo bug fix & mpp reproducibility
14   !!             4.0  !  2016-06  (G. Madec)  move to usrdef_closea, remove clo_ups
15   !!----------------------------------------------------------------------
16
17   !!----------------------------------------------------------------------
18   !!   dom_clo    : modification of the ocean domain for closed seas cases
19   !!   sbc_clo    : Special handling of closed seas
20   !!   clo_rnf    : set close sea outflows as river mouths (see sbcrnf)
21   !!   clo_bat    : set to zero a field over closed sea (see domzrg)
22   !!----------------------------------------------------------------------
23   USE oce             ! dynamics and tracers
24   USE dom_oce         ! ocean space and time domain
25   USE phycst          ! physical constants
26   USE sbc_oce         ! ocean surface boundary conditions
27   !
28   USE in_out_manager  ! I/O manager
29   USE lib_fortran,    ONLY: glob_sum, DDPDD
30   USE lbclnk          ! lateral boundary condition - MPP exchanges
31   USE lib_mpp         ! MPP library
32   USE timing
33
34   IMPLICIT NONE
35   PRIVATE
36
37   PUBLIC dom_clo      ! called by domain module
38   PUBLIC sbc_clo      ! called by step module
39   PUBLIC clo_rnf      ! called by sbcrnf module
40   PUBLIC clo_bat      ! called in domzgr module
41
42   INTEGER, PUBLIC, PARAMETER          ::   jpncs   = 4      !: number of closed sea
43   INTEGER, PUBLIC, DIMENSION(jpncs)   ::   ncstt            !: Type of closed sea
44   INTEGER, PUBLIC, DIMENSION(jpncs)   ::   ncsi1, ncsj1     !: south-west closed sea limits (i,j)
45   INTEGER, PUBLIC, DIMENSION(jpncs)   ::   ncsi2, ncsj2     !: north-east closed sea limits (i,j)
46   INTEGER, PUBLIC, DIMENSION(jpncs)   ::   ncsnr            !: number of point where run-off pours
47   INTEGER, PUBLIC, DIMENSION(jpncs,4) ::   ncsir, ncsjr     !: Location of runoff
48
49   REAL(wp), DIMENSION (jpncs+1)       ::   surf             ! closed sea surface
50
51   !! * Substitutions
52#  include "vectopt_loop_substitute.h90"
53   !!----------------------------------------------------------------------
54   !! NEMO/OPA 4.0 , NEMO Consortium (2016)
55   !! $Id$
56   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
57   !!----------------------------------------------------------------------
58CONTAINS
59
60   SUBROUTINE dom_clo( cd_cfg, kcfg )
61      !!---------------------------------------------------------------------
62      !!                  ***  ROUTINE dom_clo  ***
63      !!       
64      !! ** Purpose :   Closed sea domain initialization
65      !!
66      !! ** Method  :   if a closed sea is located only in a model grid point
67      !!                just the thermodynamic processes are applied.
68      !!
69      !! ** Action  :   ncsi1(), ncsj1() : south-west closed sea limits (i,j)
70      !!                ncsi2(), ncsj2() : north-east Closed sea limits (i,j)
71      !!                ncsir(), ncsjr() : Location of runoff
72      !!                ncsnr            : number of point where run-off pours
73      !!                ncstt            : Type of closed sea
74      !!                                   =0 spread over the world ocean
75      !!                                   =2 put at location runoff
76      !!----------------------------------------------------------------------
77      CHARACTER(len=*), INTENT(in   ) ::   cd_cfg   ! configuration name
78      INTEGER         , INTENT(in   ) ::   kcfg     ! configuration identifier
79      !
80      INTEGER ::   jc      ! dummy loop indices
81      INTEGER ::   isrow   ! local index
82      !!----------------------------------------------------------------------
83      !
84      IF(lwp) WRITE(numout,*)
85      IF(lwp) WRITE(numout,*)'dom_clo : closed seas '
86      IF(lwp) WRITE(numout,*)'~~~~~~~'
87      !
88      ! initial values
89      ncsnr(:) = 1  ;  ncsi1(:) = 1  ;  ncsi2(:) = 1  ;  ncsir(:,:) = 1
90      ncstt(:) = 0  ;  ncsj1(:) = 1  ;  ncsj2(:) = 1  ;  ncsjr(:,:) = 1
91      !
92      ! set the closed seas (in data domain indices)
93      ! -------------------
94      !
95      IF( cd_cfg == "orca" ) THEN      !==  ORCA configuration  ==!
96         !
97         SELECT CASE ( kcfg )
98         !                                           ! =======================
99         CASE ( 1 )                                  !  ORCA_R1 configuration
100            !                                        ! =======================
101            IF(lwp) WRITE(numout,*)'   ORCA_R1 closed seas :  only the Caspian Sea'
102            ! This dirty section will be suppressed by simplification process:
103            ! all this will come back in input files
104            ! Currently these hard-wired indices relate to configuration with
105            ! extend grid (jpjglo=332)
106            isrow = 332 - jpjglo
107            !
108            ncsnr(1)   = 1    ; ncstt(1)   = 0           ! Caspian Sea  (spread over the globe)
109            ncsi1(1)   = 332  ; ncsj1(1)   = 243 - isrow
110            ncsi2(1)   = 344  ; ncsj2(1)   = 275 - isrow
111            ncsir(1,1) = 1    ; ncsjr(1,1) = 1
112            !                                       
113            !                                        ! =======================
114         CASE ( 2 )                                  !  ORCA_R2 configuration
115            !                                        ! =======================
116            IF(lwp) WRITE(numout,*)'   ORCA_R2 closed seas and lakes : '
117            !                                            ! Caspian Sea
118            IF(lwp) WRITE(numout,*)'      Caspian Sea  '
119            ncsnr(1)   =   1  ;  ncstt(1)   =   0           ! spread over the globe
120            ncsi1(1)   =  11  ;  ncsj1(1)   = 103
121            ncsi2(1)   =  17  ;  ncsj2(1)   = 112
122            ncsir(1,1) =   1  ;  ncsjr(1,1) =   1 
123            !                                            ! Great North American Lakes
124            IF(lwp) WRITE(numout,*)'      Great North American Lakes  '
125            ncsnr(2)   =   1  ;  ncstt(2)   =   2           ! put at St Laurent mouth
126            ncsi1(2)   =  97  ;  ncsj1(2)   = 107
127            ncsi2(2)   = 103  ;  ncsj2(2)   = 111
128            ncsir(2,1) = 110  ;  ncsjr(2,1) = 111           
129            !                                            ! Black Sea (crossed by the cyclic boundary condition)
130            IF(lwp) WRITE(numout,*)'      Black Sea  '
131            ncsnr(3:4) =   4  ;  ncstt(3:4) =   2           ! put in Med Sea (north of Aegean Sea)
132            ncsir(3:4,1) = 171;  ncsjr(3:4,1) = 106         !
133            ncsir(3:4,2) = 170;  ncsjr(3:4,2) = 106 
134            ncsir(3:4,3) = 171;  ncsjr(3:4,3) = 105 
135            ncsir(3:4,4) = 170;  ncsjr(3:4,4) = 105 
136            ncsi1(3)   = 174  ;  ncsj1(3)   = 107           ! 1 : west part of the Black Sea     
137            ncsi2(3)   = 181  ;  ncsj2(3)   = 112           !            (ie west of the cyclic b.c.)
138            ncsi1(4)   =   2  ;  ncsj1(4)   = 107           ! 2 : east part of the Black Sea
139            ncsi2(4)   =   6  ;  ncsj2(4)   = 112           !           (ie east of the cyclic b.c.)
140            !
141            !                                        ! =========================
142         CASE ( 025 )                                !  ORCA_R025 configuration
143            !                                        ! =========================
144            IF(lwp) WRITE(numout,*)'   ORCA_R025 closed seas : '
145            !                                            ! Caspian Sea
146            IF(lwp) WRITE(numout,*)'      Caspian Sea  '
147            ncsnr(1)   = 1    ; ncstt(1)   = 0               ! Caspian + Aral sea
148            ncsi1(1)   = 1330 ; ncsj1(1)   = 645
149            ncsi2(1)   = 1400 ; ncsj2(1)   = 795
150            ncsir(1,1) = 1    ; ncsjr(1,1) = 1
151            !                                       
152            IF(lwp) WRITE(numout,*)'      Azov Sea  '
153            ncsnr(2)   = 1    ; ncstt(2)   = 0               ! Azov Sea
154            ncsi1(2)   = 1284 ; ncsj1(2)   = 722
155            ncsi2(2)   = 1304 ; ncsj2(2)   = 747
156            ncsir(2,1) = 1    ; ncsjr(2,1) = 1
157            !
158         END SELECT
159         !
160      ELSE                             !==  No closed sea in the configuration  ==!
161         !
162         IF(lwp) WRITE(numout,*)'   No closed seas or lakes in the configuration '
163         !
164      ENDIF
165
166      ! convert the position in local domain indices
167      ! --------------------------------------------
168      DO jc = 1, jpncs
169         ncsi1(jc)   = mi0( ncsi1(jc) )
170         ncsj1(jc)   = mj0( ncsj1(jc) )
171         !
172         ncsi2(jc)   = mi1( ncsi2(jc) )   
173         ncsj2(jc)   = mj1( ncsj2(jc) ) 
174      END DO
175      !
176   END SUBROUTINE dom_clo
177
178
179   SUBROUTINE sbc_clo( kt, cd_cfg, kcfg )
180      !!---------------------------------------------------------------------
181      !!                  ***  ROUTINE sbc_clo  ***
182      !!                   
183      !! ** Purpose :   Special handling of closed seas
184      !!
185      !! ** Method  :   Water flux is forced to zero over closed sea
186      !!      Excess is shared between remaining ocean, or
187      !!      put as run-off in open ocean.
188      !!
189      !! ** Action  :   emp updated surface freshwater fluxes and associated heat content at kt
190      !!----------------------------------------------------------------------
191      INTEGER         , INTENT(in   ) ::   kt       ! ocean model time step
192      CHARACTER(len=*), INTENT(in   ) ::   cd_cfg   ! configuration name
193      INTEGER         , INTENT(in   ) ::   kcfg     ! configuration identifier
194      !
195      INTEGER             ::   ji, jj, jc, jn   ! dummy loop indices
196      REAL(wp), PARAMETER ::   rsmall = 1.e-20_wp    ! Closed sea correction epsilon
197      REAL(wp)            ::   zze2, ztmp, zcorr     !
198      REAL(wp)            ::   zcoef, zcoef1         !
199      COMPLEX(wp)         ::   ctmp 
200      REAL(wp), DIMENSION(jpncs) ::   zfwf   ! 1D workspace
201      !!----------------------------------------------------------------------
202      !
203      IF( nn_timing == 1 )  CALL timing_start('sbc_clo')
204      !                                                   !------------------!
205      IF( kt == nit000 ) THEN                             !  Initialisation  !
206         !                                                !------------------!
207         IF(lwp) WRITE(numout,*)
208         IF(lwp) WRITE(numout,*)'sbc_clo : closed seas '
209         IF(lwp) WRITE(numout,*)'~~~~~~~'
210
211         surf(:) = 0._wp
212         !
213         surf(jpncs+1) = glob_sum( e1e2t(:,:) )   ! surface of the global ocean
214         !
215         !                                        ! surface of closed seas
216         DO jc = 1, jpncs
217            ctmp = CMPLX( 0.e0, 0.e0, wp )
218            DO jj = ncsj1(jc), ncsj2(jc)
219               DO ji = ncsi1(jc), ncsi2(jc)
220                  ztmp = e1e2t(ji,jj) * tmask_i(ji,jj)
221                  CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp )
222               END DO
223            END DO
224            IF( lk_mpp )   CALL mpp_sum( ctmp )
225            surf(jc) = REAL(ctmp,wp)
226         END DO
227
228         IF(lwp) WRITE(numout,*)'     Closed sea surfaces'
229         DO jc = 1, jpncs
230            IF(lwp)WRITE(numout,FMT='(1I3,4I4,5X,F16.2)') jc, ncsi1(jc), ncsi2(jc), ncsj1(jc), ncsj2(jc), surf(jc)
231         END DO
232
233         ! jpncs+1 : surface of sea, closed seas excluded
234         DO jc = 1, jpncs
235            surf(jpncs+1) = surf(jpncs+1) - surf(jc)
236         END DO           
237         !
238      ENDIF
239      !                                                   !--------------------!
240      !                                                   !  update emp        !
241      zfwf = 0.e0_wp                                      !--------------------!
242      DO jc = 1, jpncs
243         ctmp = CMPLX( 0.e0, 0.e0, wp )
244         DO jj = ncsj1(jc), ncsj2(jc)
245            DO ji = ncsi1(jc), ncsi2(jc)
246               ztmp = e1e2t(ji,jj) * ( emp(ji,jj)-rnf(ji,jj) ) * tmask_i(ji,jj)
247               CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp )
248            END DO 
249         END DO
250         IF( lk_mpp )   CALL mpp_sum( ctmp )
251         zfwf(jc) = REAL(ctmp,wp)
252      END DO
253
254      IF( cd_cfg == "orca" .AND. kcfg == 2 ) THEN      ! Black Sea case for ORCA_R2 configuration
255         zze2    = ( zfwf(3) + zfwf(4) ) * 0.5_wp
256         zfwf(3) = zze2
257         zfwf(4) = zze2
258      ENDIF
259
260      zcorr = 0._wp
261
262      DO jc = 1, jpncs
263         !
264         ! The following if avoids the redistribution of the round off
265         IF ( ABS(zfwf(jc) / surf(jpncs+1) ) > rsmall) THEN
266            !
267            IF( ncstt(jc) == 0 ) THEN           ! water/evap excess is shared by all open ocean
268               zcoef    = zfwf(jc) / surf(jpncs+1)
269               zcoef1   = rcp * zcoef
270               emp(:,:) = emp(:,:) + zcoef
271               qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:)
272               ! accumulate closed seas correction
273               zcorr    = zcorr    + zcoef
274               !
275            ELSEIF( ncstt(jc) == 1 ) THEN       ! Excess water in open sea, at outflow location, excess evap shared
276               IF ( zfwf(jc) <= 0.e0_wp ) THEN
277                   DO jn = 1, ncsnr(jc)
278                     ji = mi0(ncsir(jc,jn))
279                     jj = mj0(ncsjr(jc,jn)) ! Location of outflow in open ocean
280                     IF (      ji > 1 .AND. ji < jpi   &
281                         .AND. jj > 1 .AND. jj < jpj ) THEN
282                         zcoef      = zfwf(jc) / ( REAL(ncsnr(jc)) * e1e2t(ji,jj) )
283                         zcoef1     = rcp * zcoef
284                         emp(ji,jj) = emp(ji,jj) + zcoef
285                         qns(ji,jj) = qns(ji,jj) - zcoef1 * sst_m(ji,jj)
286                     ENDIF
287                   END DO
288               ELSE
289                   zcoef    = zfwf(jc) / surf(jpncs+1)
290                   zcoef1   = rcp * zcoef
291                   emp(:,:) = emp(:,:) + zcoef
292                   qns(:,:) = qns(:,:) - zcoef1 * sst_m(:,:)
293                   ! accumulate closed seas correction
294                   zcorr    = zcorr    + zcoef
295               ENDIF
296            ELSEIF( ncstt(jc) == 2 ) THEN       ! Excess e-p-r (either sign) goes to open ocean, at outflow location
297               DO jn = 1, ncsnr(jc)
298                  ji = mi0(ncsir(jc,jn))
299                  jj = mj0(ncsjr(jc,jn)) ! Location of outflow in open ocean
300                  IF(      ji > 1 .AND. ji < jpi    &
301                     .AND. jj > 1 .AND. jj < jpj ) THEN
302                     zcoef      = zfwf(jc) / ( REAL(ncsnr(jc)) *  e1e2t(ji,jj) )
303                     zcoef1     = rcp * zcoef
304                     emp(ji,jj) = emp(ji,jj) + zcoef
305                     qns(ji,jj) = qns(ji,jj) - zcoef1 * sst_m(ji,jj)
306                  ENDIF
307               END DO
308            ENDIF 
309            !
310            DO jj = ncsj1(jc), ncsj2(jc)
311               DO ji = ncsi1(jc), ncsi2(jc)
312                  zcoef      = zfwf(jc) / surf(jc)
313                  zcoef1     = rcp * zcoef
314                  emp(ji,jj) = emp(ji,jj) - zcoef
315                  qns(ji,jj) = qns(ji,jj) + zcoef1 * sst_m(ji,jj)
316               END DO 
317            END DO 
318            !
319         END IF
320      END DO
321
322      IF ( ABS(zcorr) > rsmall ) THEN      ! remove the global correction from the closed seas
323         DO jc = 1, jpncs                  ! only if it is large enough
324            DO jj = ncsj1(jc), ncsj2(jc)
325               DO ji = ncsi1(jc), ncsi2(jc)
326                  emp(ji,jj) = emp(ji,jj) - zcorr
327                  qns(ji,jj) = qns(ji,jj) + rcp * zcorr * sst_m(ji,jj)
328               END DO 
329             END DO
330          END DO
331      ENDIF
332      !
333      emp (:,:) = emp (:,:) * tmask(:,:,1)
334      !
335      CALL lbc_lnk( emp , 'T', 1._wp )
336      !
337      IF( nn_timing == 1 )  CALL timing_stop('sbc_clo')
338      !
339   END SUBROUTINE sbc_clo
340
341
342   SUBROUTINE clo_rnf( p_rnfmsk )
343      !!---------------------------------------------------------------------
344      !!                  ***  ROUTINE sbc_rnf  ***
345      !!                   
346      !! ** Purpose :   allow the treatment of closed sea outflow grid-points
347      !!                to be the same as river mouth grid-points
348      !!
349      !! ** Method  :   set to 1 the runoff mask (mskrnf, see sbcrnf module)
350      !!                at the closed sea outflow grid-point.
351      !!
352      !! ** Action  :   update (p_)mskrnf (set 1 at closed sea outflow)
353      !!----------------------------------------------------------------------
354      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   p_rnfmsk   ! river runoff mask (rnfmsk array)
355      !
356      INTEGER  ::   jc, jn, ji, jj      ! dummy loop indices
357      !!----------------------------------------------------------------------
358      !
359      DO jc = 1, jpncs
360         IF( ncstt(jc) >= 1 ) THEN            ! runoff mask set to 1 at closed sea outflows
361             DO jn = 1, 4
362                DO jj =    mj0( ncsjr(jc,jn) ), mj1( ncsjr(jc,jn) )
363                   DO ji = mi0( ncsir(jc,jn) ), mi1( ncsir(jc,jn) )
364                      p_rnfmsk(ji,jj) = MAX( p_rnfmsk(ji,jj), 1.0_wp )
365                   END DO
366                END DO
367            END DO
368         ENDIF
369      END DO 
370      !
371   END SUBROUTINE clo_rnf
372   
373     
374   SUBROUTINE clo_bat( k_top, k_bot )
375      !!---------------------------------------------------------------------
376      !!                  ***  ROUTINE clo_bat  ***
377      !!                   
378      !! ** Purpose :   suppress closed sea from the domain
379      !!
380      !! ** Method  :   set first and last ocean level to 0 over the closed seas.
381      !!
382      !! ** Action  :   set pbat=0 and kbat=0 over closed seas
383      !!----------------------------------------------------------------------
384      INTEGER, DIMENSION(:,:), INTENT(inout) ::   k_top, k_bot   ! ocean first and last level indices
385      !
386      INTEGER  ::   jc, ji, jj      ! dummy loop indices
387      !!----------------------------------------------------------------------
388      !
389      DO jc = 1, jpncs
390         DO jj = ncsj1(jc), ncsj2(jc)
391            DO ji = ncsi1(jc), ncsi2(jc)
392               k_top(ji,jj) = 0   
393               k_bot(ji,jj) = 0   
394            END DO
395         END DO
396       END DO 
397       !
398   END SUBROUTINE clo_bat
399
400   !!======================================================================
401END MODULE usrdef_closea
402
Note: See TracBrowser for help on using the repository browser.