source: NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcdmp.F90 @ 10946

Last change on this file since 10946 was 10946, checked in by acc, 19 months ago

2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps : Convert STO, TRD and USR modules and all knock on effects of these conversions. Note change to USR module may have implications for the TEST CASES (not tested yet). Standard SETTE tested only

  • Property svn:keywords set to Id
File size: 18.6 KB
Line 
1MODULE trcdmp
2   !!======================================================================
3   !!                       ***  MODULE  trcdmp  ***
4   !! Ocean physics: internal restoring trend on passive tracers
5   !!======================================================================
6   !! History :  OPA  !  1991-03  (O. Marti, G. Madec)  Original code
7   !!                 !  1996-01  (G. Madec) statement function for e3
8   !!                 !  1997-05  (H. Loukos)  adapted for passive tracers
9   !!    NEMO    9.0  !  2004-03  (C. Ethe)    free form + modules
10   !!            3.2  !  2007-02  (C. Deltel)  Diagnose ML trends for passive tracers
11   !!            3.3  !  2010-06  (C. Ethe, G. Madec) merge TRA-TRC
12   !!----------------------------------------------------------------------
13#if  defined key_top 
14   !!----------------------------------------------------------------------
15   !!   trc_dmp      : update the tracer trend with the internal damping
16   !!   trc_dmp_init : initialization, namlist read, parameters control
17   !!----------------------------------------------------------------------
18   USE oce_trc         ! ocean dynamics and tracers variables
19   USE trc             ! ocean passive tracers variables
20   USE trcdta
21   USE tradmp
22   USE trdtra
23   USE trd_oce
24   !
25   USE iom
26   USE prtctl_trc      ! Print control for debbuging
27
28   IMPLICIT NONE
29   PRIVATE
30
31   PUBLIC trc_dmp     
32   PUBLIC trc_dmp_clo   
33   PUBLIC trc_dmp_alloc 
34   PUBLIC trc_dmp_ini   
35
36   INTEGER            , PUBLIC ::   nn_zdmp_tr    !: = 0/1/2 flag for damping in the mixed layer
37   CHARACTER(LEN=200) , PUBLIC ::   cn_resto_tr   !: File containing restoration coefficient
38
39   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   restotr   ! restoring coeff. on tracers (s-1)
40
41   INTEGER, PARAMETER         ::   npncts = 8       ! number of closed sea
42   INTEGER, DIMENSION(npncts) ::   nctsi1, nctsj1   ! south-west closed sea limits (i,j)
43   INTEGER, DIMENSION(npncts) ::   nctsi2, nctsj2   ! north-east closed sea limits (i,j)
44
45   !! * Substitutions
46#  include "vectopt_loop_substitute.h90"
47   !!----------------------------------------------------------------------
48   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
49   !! $Id$
50   !! Software governed by the CeCILL license (see ./LICENSE)
51   !!----------------------------------------------------------------------
52CONTAINS
53
54   INTEGER FUNCTION trc_dmp_alloc()
55      !!----------------------------------------------------------------------
56      !!                   ***  ROUTINE trc_dmp_alloc  ***
57      !!----------------------------------------------------------------------
58      ALLOCATE( restotr(jpi,jpj,jpk) , STAT=trc_dmp_alloc )
59      !
60      IF( trc_dmp_alloc /= 0 )   CALL ctl_warn('trc_dmp_alloc: failed to allocate array')
61      !
62   END FUNCTION trc_dmp_alloc
63
64
65   SUBROUTINE trc_dmp( kt, Kmm, Krhs )
66      !!----------------------------------------------------------------------
67      !!                   ***  ROUTINE trc_dmp  ***
68      !!                 
69      !! ** Purpose :   Compute the passive tracer trend due to a newtonian damping
70      !!      of the tracer field towards given data field and add it to the
71      !!      general tracer trends.
72      !!
73      !! ** Method  :   Newtonian damping towards trdta computed
74      !!      and add to the general tracer trends:
75      !!                     trn = tra + restotr * (trdta - trb)
76      !!         The trend is computed either throughout the water column
77      !!      (nlmdmptr=0) or in area of weak vertical mixing (nlmdmptr=1) or
78      !!      below the well mixed layer (nlmdmptr=2)
79      !!
80      !! ** Action  : - update the tracer trends tra with the newtonian
81      !!                damping trends.
82      !!              - save the trends ('key_trdmxl_trc')
83      !!----------------------------------------------------------------------
84      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
85      INTEGER, INTENT(in) ::   Kmm, Krhs  ! time level indices
86      !
87      INTEGER ::   ji, jj, jk, jn, jl   ! dummy loop indices
88      CHARACTER (len=22) ::   charout
89      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrtrd
90      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   ztrcdta   ! 3D  workspace
91      !!----------------------------------------------------------------------
92      !
93      IF( ln_timing )   CALL timing_start('trc_dmp')
94      !
95      IF( l_trdtrc )   ALLOCATE( ztrtrd(jpi,jpj,jpk) )   ! temporary save of trends
96      !
97      IF( nb_trcdta > 0 ) THEN  ! Initialisation of tracer from a file that may also be used for damping
98         !
99         ALLOCATE( ztrcdta(jpi,jpj,jpk) )    ! Memory allocation
100         !                                                          ! ===========
101         DO jn = 1, jptra                                           ! tracer loop
102            !                                                       ! ===========
103            IF( l_trdtrc ) ztrtrd(:,:,:) = tra(:,:,:,jn)    ! save trends
104            !
105            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file
106               !
107               jl = n_trc_index(jn) 
108               CALL trc_dta( kt, sf_trcdta(jl), rf_trfac(jl), ztrcdta )   ! read tracer data at nit000
109               !
110               SELECT CASE ( nn_zdmp_tr )
111               !
112               CASE( 0 )                !==  newtonian damping throughout the water column  ==!
113                  DO jk = 1, jpkm1
114                     DO jj = 2, jpjm1
115                        DO ji = fs_2, fs_jpim1   ! vector opt.
116                           tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) )
117                        END DO
118                     END DO
119                  END DO
120                  !
121               CASE ( 1 )                !==  no damping in the turbocline (avt > 5 cm2/s)  ==!
122                  DO jk = 1, jpkm1
123                     DO jj = 2, jpjm1
124                        DO ji = fs_2, fs_jpim1   ! vector opt.
125                           IF( avt(ji,jj,jk) <= avt_c )  THEN
126                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) )
127                           ENDIF
128                        END DO
129                     END DO
130                  END DO
131                  !
132               CASE ( 2 )               !==  no damping in the mixed layer   ==!
133                  DO jk = 1, jpkm1
134                     DO jj = 2, jpjm1
135                        DO ji = fs_2, fs_jpim1   ! vector opt.
136                           IF( gdept_n(ji,jj,jk) >= hmlp (ji,jj) ) THEN
137                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) )
138                           END IF
139                        END DO
140                     END DO
141                  END DO
142                 
143               END SELECT
144               !
145            ENDIF
146            !
147            IF( l_trdtrc ) THEN
148               ztrtrd(:,:,:) = tra(:,:,:,jn) -  ztrtrd(:,:,:)
149               CALL trd_tra( kt, Kmm, Krhs, 'TRC', jn, jptra_dmp, ztrtrd )
150            END IF
151            !                                                       ! ===========
152         END DO                                                     ! tracer loop
153         !                                                          ! ===========
154         DEALLOCATE( ztrcdta )
155      ENDIF
156      !
157      IF( l_trdtrc )  DEALLOCATE( ztrtrd )
158      !                                          ! print mean trends (used for debugging)
159      IF( ln_ctl ) THEN
160         WRITE(charout, FMT="('dmp ')")
161         CALL prt_ctl_trc_info(charout)
162         CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
163      ENDIF
164      !
165      IF( ln_timing )   CALL timing_stop('trc_dmp')
166      !
167   END SUBROUTINE trc_dmp
168
169
170   SUBROUTINE trc_dmp_ini
171      !!----------------------------------------------------------------------
172      !!                  ***  ROUTINE trc_dmp_ini  ***
173      !!
174      !! ** Purpose :   Initialization for the newtonian damping
175      !!
176      !! ** Method  :   read the nammbf namelist and check the parameters
177      !!              called by trc_dmp at the first timestep (nittrc000)
178      !!----------------------------------------------------------------------
179      INTEGER ::   ios, imask  ! local integers
180      !!
181      NAMELIST/namtrc_dmp/ nn_zdmp_tr , cn_resto_tr
182      !!----------------------------------------------------------------------
183      !
184      REWIND( numnat_ref )              ! Namelist namtrc_dmp in reference namelist : Passive tracers newtonian damping
185      READ  ( numnat_ref, namtrc_dmp, IOSTAT = ios, ERR = 909)
186909   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtrc_dmp in reference namelist', lwp )
187      REWIND( numnat_cfg )              ! Namelist namtrc_dmp in configuration namelist : Passive tracers newtonian damping
188      READ  ( numnat_cfg, namtrc_dmp, IOSTAT = ios, ERR = 910)
189910   IF( ios >  0 )   CALL ctl_nam ( ios , 'namtrc_dmp in configuration namelist', lwp )
190      IF(lwm) WRITE ( numont, namtrc_dmp )
191
192      IF(lwp) THEN                       ! Namelist print
193         WRITE(numout,*)
194         WRITE(numout,*) 'trc_dmp : Passive tracers newtonian damping'
195         WRITE(numout,*) '~~~~~~~'
196         WRITE(numout,*) '   Namelist namtrc_dmp : set damping parameter'
197         WRITE(numout,*) '      mixed layer damping option     nn_zdmp_tr  = ', nn_zdmp_tr, '(zoom: forced to 0)'
198         WRITE(numout,*) '      Restoration coeff file         cn_resto_tr = ', cn_resto_tr
199      ENDIF
200      !                          ! Allocate arrays
201      IF( trc_dmp_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'trc_dmp_ini: unable to allocate arrays' )
202      !
203      SELECT CASE ( nn_zdmp_tr )
204      CASE ( 0 )   ;   IF(lwp) WRITE(numout,*) '      ===>>   tracer damping throughout the water column'
205      CASE ( 1 )   ;   IF(lwp) WRITE(numout,*) '      ===>>   no tracer damping in the turbocline (avt > 5 cm2/s)'
206      CASE ( 2 )   ;   IF(lwp) WRITE(numout,*) '      ===>>   no tracer damping in the mixed layer'
207      CASE DEFAULT
208         WRITE(ctmp1,*) 'bad flag value for nn_zdmp_tr = ', nn_zdmp_tr
209         CALL ctl_stop(ctmp1)
210      END SELECT
211
212      IF( .NOT.lk_c1d ) THEN
213         IF( .NOT.ln_tradmp )   &
214            &   CALL ctl_stop( 'passive tracer damping need ln_tradmp to compute damping coef.' )
215         !
216         !                          ! Read damping coefficients from file
217         !Read in mask from file
218         CALL iom_open ( cn_resto_tr, imask)
219         CALL iom_get  ( imask, jpdom_autoglo, 'resto', restotr)
220         CALL iom_close( imask )
221         !
222      ENDIF
223      !
224   END SUBROUTINE trc_dmp_ini
225
226
227   SUBROUTINE trc_dmp_clo( kt )
228      !!---------------------------------------------------------------------
229      !!                  ***  ROUTINE trc_dmp_clo  ***
230      !!
231      !! ** Purpose :   Closed sea domain initialization
232      !!
233      !! ** Method  :   if a closed sea is located only in a model grid point
234      !!                we restore to initial data
235      !!
236      !! ** Action  :   nctsi1(), nctsj1() : south-west closed sea limits (i,j)
237      !!                nctsi2(), nctsj2() : north-east Closed sea limits (i,j)
238      !!----------------------------------------------------------------------
239      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
240      !
241      INTEGER :: ji , jj, jk, jn, jl, jc                    ! dummy loop indicesa
242      INTEGER :: isrow                                      ! local index
243      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrcdta       ! 3D  workspace
244      !!----------------------------------------------------------------------
245
246      IF( kt == nit000 ) THEN
247         ! initial values
248         nctsi1(:) = 1  ;  nctsi2(:) = 1
249         nctsj1(:) = 1  ;  nctsj2(:) = 1
250
251         ! set the closed seas (in data domain indices)
252         ! -------------------
253
254         IF( cn_cfg == "orca" .OR. cn_cfg == "ORCA") THEN
255            !
256            SELECT CASE ( nn_cfg )
257            !                                           ! =======================
258            CASE ( 1 )                                  ! eORCA_R1 configuration
259            !                                           ! =======================
260            isrow = 332 - jpjglo
261            !
262            nctsi1(1)   = 333  ; nctsj1(1)   = 243 - isrow   ! Caspian Sea
263            nctsi2(1)   = 342  ; nctsj2(1)   = 274 - isrow
264            !                                       
265            nctsi1(2)   = 198  ; nctsj1(2)   = 258 - isrow   ! Lake Superior
266            nctsi2(2)   = 204  ; nctsj2(2)   = 262 - isrow
267            !                                         
268            nctsi1(3)   = 201  ; nctsj1(3)   = 250 - isrow   ! Lake Michigan
269            nctsi2(3)   = 203  ; nctsj2(3)   = 256 - isrow
270            !                                       
271            nctsi1(4)   = 204  ; nctsj1(4)   = 252 - isrow   ! Lake Huron
272            nctsi2(4)   = 209  ; nctsj2(4)   = 256 - isrow
273            !                                       
274            nctsi1(5)   = 206  ; nctsj1(5)   = 249 - isrow   ! Lake Erie
275            nctsi2(5)   = 209  ; nctsj2(5)   = 251 - isrow
276            !                                       
277            nctsi1(6)   = 210  ; nctsj1(6)   = 252 - isrow   ! Lake Ontario
278            nctsi2(6)   = 212  ; nctsj2(6)   = 252 - isrow
279            !                                       
280            nctsi1(7)   = 321  ; nctsj1(7)   = 180 - isrow   ! Victoria Lake
281            nctsi2(7)   = 322  ; nctsj2(7)   = 189 - isrow
282            !                                       
283            nctsi1(8)   = 297  ; nctsj1(8)   = 270 - isrow   ! Baltic Sea
284            nctsi2(8)   = 308  ; nctsj2(8)   = 293 - isrow
285            !                                       
286            !                                           ! =======================
287            CASE ( 2 )                                  !  ORCA_R2 configuration
288               !                                        ! =======================
289               !                                     
290               nctsi1(1)   =  11  ;  nctsj1(1)   = 103       ! Caspian Sea
291               nctsi2(1)   =  17  ;  nctsj2(1)   = 112
292               !                                     
293               nctsi1(2)   =  97  ;  nctsj1(2)   = 107       ! Great North American Lakes
294               nctsi2(2)   = 103  ;  nctsj2(2)   = 111
295               !                                     
296               nctsi1(3)   = 174  ;  nctsj1(3)   = 107       ! Black Sea 1 : west part of the Black Sea
297               nctsi2(3)   = 181  ;  nctsj2(3)   = 112
298              !                                     
299               nctsi1(4)   =   2  ;  nctsj1(4)   = 107      ! Black Sea 2 : est part of the Black Sea
300               nctsi2(4)   =   6  ;  nctsj2(4)   = 112
301               !                                     
302               nctsi1(5)   =  145 ;  nctsj1(5)   = 116       ! Baltic Sea
303               nctsi2(5)   =  150 ;  nctsj2(5)   = 126
304               !                                        ! =======================
305            CASE ( 4 )                                  !  ORCA_R4 configuration
306               !                                        ! =======================
307               !                                   
308               nctsi1(1)   =  4  ;  nctsj1(1)   = 53         ! Caspian Sea
309               nctsi2(1)   =  4  ;  nctsj2(1)   = 56
310               !                                   
311               nctsi1(2)   = 49  ;  nctsj1(2)   = 55         ! Great North American Lakes
312               nctsi2(2)   = 51  ;  nctsj2(2)   = 56
313               !                                   
314               nctsi1(3)   = 88  ;  nctsj1(3)   = 55         ! Black Sea
315               nctsi2(3)   = 91  ;  nctsj2(3)   = 56
316               !                                   
317               nctsi1(4)   = 75  ;  nctsj1(4)   = 59         ! Baltic Sea
318               nctsi2(4)   = 76  ;  nctsj2(4)   = 61
319               !                                        ! =======================
320            CASE ( 025 )                                ! ORCA_R025 configuration
321               !                                        ! =======================
322               !                                   
323               nctsi1(1)   = 1330 ; nctsj1(1)   = 645        ! Caspian + Aral sea
324               nctsi2(1)   = 1400 ; nctsj2(1)   = 795
325               !                                   
326               nctsi1(2)   = 1284 ; nctsj1(2)   = 722        ! Azov Sea
327               nctsi2(2)   = 1304 ; nctsj2(2)   = 747
328               !
329            END SELECT
330            !
331         ENDIF
332         !
333         ! convert the position in local domain indices
334         ! --------------------------------------------
335         DO jc = 1, npncts
336            nctsi1(jc)   = mi0( nctsi1(jc) )
337            nctsj1(jc)   = mj0( nctsj1(jc) )
338            !
339            nctsi2(jc)   = mi1( nctsi2(jc) )
340            nctsj2(jc)   = mj1( nctsj2(jc) )
341         END DO
342         !
343      ENDIF
344
345      ! Restore close seas values to initial data
346      IF( ln_trcdta .AND. nb_trcdta > 0 )  THEN   ! Initialisation of tracer from a file that may also be used for damping
347         !
348         IF(lwp)  WRITE(numout,*)
349         IF(lwp)  WRITE(numout,*) ' trc_dmp_clo : Restoring of nutrients on close seas at time-step kt = ', kt
350         IF(lwp)  WRITE(numout,*)
351         !
352         ALLOCATE( ztrcdta(jpi,jpj,jpk) )   ! Memory allocation
353         !
354         DO jn = 1, jptra
355            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file
356                jl = n_trc_index(jn)
357                CALL trc_dta( kt, sf_trcdta(jl), rf_trfac(jl), ztrcdta )   ! read tracer data at nit000
358                DO jc = 1, npncts
359                   DO jk = 1, jpkm1
360                      DO jj = nctsj1(jc), nctsj2(jc)
361                         DO ji = nctsi1(jc), nctsi2(jc)
362                            trn(ji,jj,jk,jn) = ztrcdta(ji,jj,jk)
363                            trb(ji,jj,jk,jn) = trn(ji,jj,jk,jn)
364                         END DO
365                      END DO
366                   END DO
367                END DO
368             ENDIF
369          END DO
370          DEALLOCATE( ztrcdta )
371      ENDIF
372      !
373   END SUBROUTINE trc_dmp_clo
374 
375#else
376   !!----------------------------------------------------------------------
377   !!  Dummy module :                                     No passive tracer
378   !!----------------------------------------------------------------------
379CONTAINS
380   SUBROUTINE trc_dmp( kt )        ! Empty routine
381      INTEGER, INTENT(in) :: kt
382      WRITE(*,*) 'trc_dmp: You should not have seen this print! error?', kt
383   END SUBROUTINE trc_dmp
384#endif
385
386   !!======================================================================
387END MODULE trcdmp
Note: See TracBrowser for help on using the repository browser.