source: branches/2013/dev_r3996_CMCC6_topbc/NEMOGCM/NEMO/TOP_SRC/TRP/trcdmp.F90 @ 4011

Last change on this file since 4011 was 4011, checked in by vichi, 8 years ago

Make a generic interface for trcdta when using other BGCM

This change introduces a more general trcdta structure that
is not strictly dependent on the number of tracers defined
in PISCES. The loop on the number of tracers is moved outside
trcdta and the tracer info and array is passed as an argument.
This allows to use trcdta as a library subroutine by the BFM and
other models.
NOTE: it must be tested throughly with all the PISCES configurations

This commit also updates the GYRE_BFM configuration and corrects
some minor missing cpp keys and real type definitions

  • Property svn:keywords set to Id
File size: 16.4 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 trcnam_trp      ! passive tracers transport namelist variables
21   USE trcdta
22   USE tradmp
23   USE prtctl_trc      ! Print control for debbuging
24   USE trdtra
25   USE trdmod_oce
26
27   IMPLICIT NONE
28   PRIVATE
29
30   PUBLIC trc_dmp            ! routine called by step.F90
31   PUBLIC trc_dmp_clo        ! routine called by step.F90
32   PUBLIC trc_dmp_alloc      ! routine called by nemogcm.F90
33
34   !                                !!* Namelist namtrc_dmp : passive tracer newtonian damping *
35   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   restotr   ! restoring coeff. on tracers (s-1)
36
37   INTEGER, PARAMETER           ::   npncts   = 5        ! number of closed sea
38   INTEGER, DIMENSION(npncts)   ::   nctsi1, nctsj1      ! south-west closed sea limits (i,j)
39   INTEGER, DIMENSION(npncts)   ::   nctsi2, nctsj2      ! north-east closed sea limits (i,j)
40
41   !! * Substitutions
42#  include "top_substitute.h90"
43   !!----------------------------------------------------------------------
44   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
45   !! $Header: /home/opalod/NEMOCVSROOT/NEMO/TOP_SRC/TRP/trcdmp.F90,v 1.11 2006/09/01 14:03:49 opalod Exp $
46   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
47   !!----------------------------------------------------------------------
48CONTAINS
49
50   INTEGER FUNCTION trc_dmp_alloc()
51      !!----------------------------------------------------------------------
52      !!                   ***  ROUTINE trc_dmp_alloc  ***
53      !!----------------------------------------------------------------------
54      ALLOCATE( restotr(jpi,jpj,jpk) , STAT=trc_dmp_alloc )
55      !
56      IF( trc_dmp_alloc /= 0 )   CALL ctl_warn('trc_dmp_alloc: failed to allocate array')
57      !
58   END FUNCTION trc_dmp_alloc
59
60
61   SUBROUTINE trc_dmp( kt )
62      !!----------------------------------------------------------------------
63      !!                   ***  ROUTINE trc_dmp  ***
64      !!                 
65      !! ** Purpose :   Compute the passive tracer trend due to a newtonian damping
66      !!      of the tracer field towards given data field and add it to the
67      !!      general tracer trends.
68      !!
69      !! ** Method  :   Newtonian damping towards trdta computed
70      !!      and add to the general tracer trends:
71      !!                     trn = tra + restotr * (trdta - trb)
72      !!         The trend is computed either throughout the water column
73      !!      (nlmdmptr=0) or in area of weak vertical mixing (nlmdmptr=1) or
74      !!      below the well mixed layer (nlmdmptr=2)
75      !!
76      !! ** Action  : - update the tracer trends tra with the newtonian
77      !!                damping trends.
78      !!              - save the trends ('key_trdmld_trc')
79      !!----------------------------------------------------------------------
80      !!
81      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
82      !!
83      INTEGER  ::   ji, jj, jk, jn, jl       ! dummy loop indices
84      REAL(wp) ::   ztra                 ! temporary scalars
85      CHARACTER (len=22) :: charout
86      REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztrtrd
87      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrcdta   ! 3D  workspace
88      !!----------------------------------------------------------------------
89      !
90      IF( nn_timing == 1 )  CALL timing_start('trc_dmp')
91      !
92      ! 0. Initialization (first time-step only)
93      !    --------------
94      IF( kt == nittrc000 ) CALL trc_dmp_init
95
96      IF( l_trdtrc )   CALL wrk_alloc( jpi, jpj, jpk, ztrtrd )   ! temporary save of trends
97      !
98      IF( nb_trcdta > 0 ) THEN  ! Initialisation of tracer from a file that may also be used for damping
99         !
100         CALL wrk_alloc( jpi, jpj, jpk, ztrcdta )    ! Memory allocation
101         !                                                          ! ===========
102         DO jn = 1, jptra                                           ! tracer loop
103            !                                                       ! ===========
104            IF( l_trdtrc ) ztrtrd(:,:,:) = tra(:,:,:,jn)    ! save trends
105            !
106            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file
107               
108               jl = n_trc_index(jn) 
109               CALL trc_dta( kt, sf_trcdta(jl),rf_trfac(jl) )   ! read tracer data at nit000
110               ztrcdta(:,:,:) = sf_trcdta(jl)%fnow(:,:,:)
111
112               SELECT CASE ( nn_zdmp_tr )
113               !
114               CASE( 0 )                !==  newtonian damping throughout the water column  ==!
115                  DO jk = 1, jpkm1
116                     DO jj = 2, jpjm1
117                        DO ji = fs_2, fs_jpim1   ! vector opt.
118                           ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) )
119                           tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
120                        END DO
121                     END DO
122                  END DO
123               !
124               CASE ( 1 )                !==  no damping in the turbocline (avt > 5 cm2/s)  ==!
125                  DO jk = 1, jpkm1
126                     DO jj = 2, jpjm1
127                        DO ji = fs_2, fs_jpim1   ! vector opt.
128                           IF( avt(ji,jj,jk) <= 5.e-4 )  THEN
129                              ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) )
130                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
131                           ENDIF
132                        END DO
133                     END DO
134                  END DO
135               !
136               CASE ( 2 )               !==  no damping in the mixed layer   ==!
137                  DO jk = 1, jpkm1
138                     DO jj = 2, jpjm1
139                        DO ji = fs_2, fs_jpim1   ! vector opt.
140                           IF( fsdept(ji,jj,jk) >= hmlp (ji,jj) ) THEN
141                              ztra = restotr(ji,jj,jk) * ( ztrcdta(ji,jj,jk) - trb(ji,jj,jk,jn) )
142                              tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn) + ztra
143                           END IF
144                        END DO
145                     END DO
146                  END DO
147               
148               END SELECT
149               !
150            ENDIF
151            !
152            IF( l_trdtrc ) THEN
153               ztrtrd(:,:,:) = tra(:,:,:,jn) -  ztrtrd(:,:,:)
154               CALL trd_tra( kt, 'TRC', jn, jptra_trd_dmp, ztrtrd )
155            END IF
156            !                                                       ! ===========
157         END DO                                                     ! tracer loop
158         !                                                          ! ===========
159         CALL wrk_dealloc( jpi, jpj, jpk, ztrcdta )
160      ENDIF
161      !
162      IF( l_trdtrc )  CALL wrk_dealloc( jpi, jpj, jpk, ztrtrd )
163      !                                          ! print mean trends (used for debugging)
164      IF( ln_ctl )   THEN
165         WRITE(charout, FMT="('dmp ')") ;  CALL prt_ctl_trc_info(charout)
166                                           CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=ctrcnm, clinfo2='trd' )
167      ENDIF
168      !
169      IF( nn_timing == 1 )  CALL timing_stop('trc_dmp')
170      !
171   END SUBROUTINE trc_dmp
172
173   SUBROUTINE trc_dmp_clo( kt )
174      !!---------------------------------------------------------------------
175      !!                  ***  ROUTINE trc_dmp_clo  ***
176      !!
177      !! ** Purpose :   Closed sea domain initialization
178      !!
179      !! ** Method  :   if a closed sea is located only in a model grid point
180      !!                we restore to initial data
181      !!
182      !! ** Action  :   nctsi1(), nctsj1() : south-west closed sea limits (i,j)
183      !!                nctsi2(), nctsj2() : north-east Closed sea limits (i,j)
184      !!----------------------------------------------------------------------
185      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
186      !
187      INTEGER :: ji, jj, jk, jn, jl, jc                     ! dummy loop indicesa
188      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrcdta       ! 3D  workspace
189
190      !!----------------------------------------------------------------------
191
192      IF( kt == nit000 ) THEN
193         ! initial values
194         nctsi1(:) = 1  ;  nctsi2(:) = 1
195         nctsj1(:) = 1  ;  nctsj2(:) = 1
196
197         ! set the closed seas (in data domain indices)
198         ! -------------------
199
200         IF( cp_cfg == "orca" ) THEN
201            !
202            SELECT CASE ( jp_cfg )
203            !                                           ! =======================
204            CASE ( 2 )                                  !  ORCA_R2 configuration
205               !                                        ! =======================
206               !                                            ! Caspian Sea
207               nctsi1(1)   =  11  ;  nctsj1(1)   = 103
208               nctsi2(1)   =  17  ;  nctsj2(1)   = 112
209               !                                            ! Great North American Lakes
210               nctsi1(2)   =  97  ;  nctsj1(2)   = 107
211               nctsi2(2)   = 103  ;  nctsj2(2)   = 111
212               !                                            ! Black Sea 1 : west part of the Black Sea
213               nctsi1(3)   = 174  ;  nctsj1(3)   = 107
214               nctsi2(3)   = 181  ;  nctsj2(3)   = 112
215              !                                            ! Black Sea 2 : est part of the Black Sea
216               nctsi1(4)   =   2  ;  nctsj1(4)   = 107
217               nctsi2(4)   =   6  ;  nctsj2(4)   = 112
218               !                                            ! Baltic Sea
219               nctsi1(5)   =  145 ;  nctsj1(5)   = 116
220               nctsi2(5)   =  150 ;  nctsj2(5)   = 126
221               !                                        ! =======================
222            CASE ( 4 )                                  !  ORCA_R4 configuration
223               !                                        ! =======================
224               !                                            ! Caspian Sea
225               nctsi1(1)   =  4  ;  nctsj1(1)   = 53
226               nctsi2(1)   =  4  ;  nctsj2(1)   = 56
227               !                                            ! Great North American Lakes
228               nctsi1(2)   = 49  ;  nctsj1(2)   = 55
229               nctsi2(2)   = 51  ;  nctsj2(2)   = 56
230               !                                            ! Black Sea
231               nctsi1(3)   = 88  ;  nctsj1(3)   = 55
232               nctsi2(3)   = 91  ;  nctsj2(3)   = 56
233               !                                            ! Baltic Sea
234               nctsi1(4)   = 75  ;  nctsj1(4)   = 59
235               nctsi2(4)   = 76  ;  nctsj2(4)   = 61
236               !                                        ! =======================
237            CASE ( 025 )                                ! ORCA_R025 configuration
238               !                                        ! =======================
239                                                     ! Caspian + Aral sea
240               nctsi1(1)   = 1330 ; nctsj1(1)   = 645
241               nctsi2(1)   = 1400 ; nctsj2(1)   = 795
242               !                                        ! Azov Sea
243               nctsi1(2)   = 1284 ; nctsj1(2)   = 722
244               nctsi2(2)   = 1304 ; nctsj2(2)   = 747
245               !
246            END SELECT
247            !
248         ENDIF
249         !
250
251         ! convert the position in local domain indices
252         ! --------------------------------------------
253         DO jc = 1, npncts
254            nctsi1(jc)   = mi0( nctsi1(jc) )
255            nctsj1(jc)   = mj0( nctsj1(jc) )
256
257            nctsi2(jc)   = mi1( nctsi2(jc) )
258            nctsj2(jc)   = mj1( nctsj2(jc) )
259         END DO
260         !
261      ENDIF
262
263      ! Restore close seas values to initial data
264      IF( ln_trcdta .AND. nb_trcdta > 0 )  THEN   ! Initialisation of tracer from a file that may also be used for damping
265         !
266         IF(lwp)  WRITE(numout,*)
267         IF(lwp)  WRITE(numout,*) ' trc_dmp_clo : Restoring of nutrients on close seas at time-step kt = ', kt
268         IF(lwp)  WRITE(numout,*)
269         !
270         CALL wrk_alloc( jpi, jpj, jpk, ztrcdta )   ! Memory allocation
271         !
272         DO jn = 1, jptra
273            IF( ln_trc_ini(jn) ) THEN      ! update passive tracers arrays with input data read from file
274                jl = n_trc_index(jn)
275                CALL trc_dta( kt, sf_trcdta(jl),rf_trfac(jl) )   ! read tracer data at nit000
276                ztrcdta(:,:,:) = sf_trcdta(jl)%fnow(:,:,:)
277                DO jc = 1, npncts
278                   DO jk = 1, jpkm1
279                      DO jj = nctsj1(jc), nctsj2(jc)
280                         DO ji = nctsi1(jc), nctsi2(jc)
281                            trn(ji,jj,jk,jn) = ztrcdta(ji,jj,jk) * tmask(ji,jj,jk)
282                            trb(ji,jj,jk,jn) = trn(ji,jj,jk,jn)
283                         ENDDO
284                      ENDDO
285                   ENDDO
286                ENDDO
287             ENDIF
288          ENDDO
289          CALL wrk_dealloc( jpi, jpj, jpk, ztrcdta )
290      ENDIF
291      !
292   END SUBROUTINE trc_dmp_clo
293
294
295   SUBROUTINE trc_dmp_init
296      !!----------------------------------------------------------------------
297      !!                  ***  ROUTINE trc_dmp_init  ***
298      !!
299      !! ** Purpose :   Initialization for the newtonian damping
300      !!
301      !! ** Method  :   read the nammbf namelist and check the parameters
302      !!              called by trc_dmp at the first timestep (nittrc000)
303      !!----------------------------------------------------------------------
304      !
305      IF( nn_timing == 1 )  CALL timing_start('trc_dmp_init')
306      !
307      SELECT CASE ( nn_hdmp_tr )
308      CASE (  -1  )   ;   IF(lwp) WRITE(numout,*) '   tracer damping in the Med & Red seas only'
309      CASE ( 1:90 )   ;   IF(lwp) WRITE(numout,*) '   tracer damping poleward of', nn_hdmp_tr, ' degrees'
310      CASE DEFAULT
311         WRITE(ctmp1,*) '          bad flag value for nn_hdmp_tr = ', nn_hdmp_tr
312         CALL ctl_stop(ctmp1)
313      END SELECT
314
315      SELECT CASE ( nn_zdmp_tr )
316      CASE ( 0 )   ;   IF(lwp) WRITE(numout,*) '   tracer damping throughout the water column'
317      CASE ( 1 )   ;   IF(lwp) WRITE(numout,*) '   no tracer damping in the turbocline (avt > 5 cm2/s)'
318      CASE ( 2 )   ;   IF(lwp) WRITE(numout,*) '   no tracer damping in the mixed layer'
319      CASE DEFAULT
320         WRITE(ctmp1,*) 'bad flag value for nn_zdmp_tr = ', nn_zdmp_tr
321         CALL ctl_stop(ctmp1)
322      END SELECT
323
324      IF( .NOT. ln_tradmp )   &
325         &   CALL ctl_stop( 'passive trace damping need key_tradmp to compute damping coef.' )
326      !
327      !                          ! Damping coefficients initialization
328      IF( lzoom ) THEN   ;   CALL dtacof_zoom( restotr )
329      ELSE               ;   CALL dtacof( nn_hdmp_tr, rn_surf_tr, rn_bot_tr, rn_dep_tr,  &
330                             &            nn_file_tr, 'TRC'     , restotr                )
331      ENDIF
332      !
333      IF( nn_timing == 1 )  CALL timing_stop('trc_dmp_init')
334      !
335   END SUBROUTINE trc_dmp_init
336
337#else
338   !!----------------------------------------------------------------------
339   !!  Dummy module :                                     No passive tracer
340   !!----------------------------------------------------------------------
341CONTAINS
342   SUBROUTINE trc_dmp( kt )        ! Empty routine
343      INTEGER, INTENT(in) :: kt
344      WRITE(*,*) 'trc_dmp: You should not have seen this print! error?', kt
345   END SUBROUTINE trc_dmp
346#endif
347
348
349   !!======================================================================
350END MODULE trcdmp
Note: See TracBrowser for help on using the repository browser.