New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
trcrad.F90 in branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/TOP_SRC/TRP – NEMO

source: branches/2015/dev_r5803_NOC_WAD/NEMOGCM/NEMO/TOP_SRC/TRP/trcrad.F90 @ 5870

Last change on this file since 5870 was 5870, checked in by acc, 8 years ago

Branch 2015/dev_r5803_NOC_WAD. Merge in trunk changes from 5803 to 5869 in preparation for merge. Also tidied and reorganised some wetting and drying code. Renamed wadlmt.F90 to wetdry.F90. Wetting drying code changes restricted to domzgr.F90, domvvl.F90 nemogcm.F90 sshwzv.F90, dynspg_ts.F90, wetdry.F90 and dynhpg.F90. Code passes full SETTE tests with ln_wd=.false.. Still awaiting test case for checking with ln_wd=.false.
  • Property svn:keywords set to Id
File size: 10.8 KB
Line 
1MODULE trcrad
2   !!======================================================================
3   !!                       ***  MODULE  trcrad  ***
4   !! Ocean passive tracers:  correction of negative concentrations
5   !!======================================================================
6   !! History :   -   !  01-01  (O. Aumont & E. Kestenare)  Original code
7   !!            1.0  !  04-03  (C. Ethe)  free form F90
8   !!----------------------------------------------------------------------
9#if defined key_top
10   !!----------------------------------------------------------------------
11   !!   'key_top'                                                TOP models
12   !!----------------------------------------------------------------------
13   !!   trc_rad    : correction of negative concentrations
14   !!----------------------------------------------------------------------
15   USE oce_trc             ! ocean dynamics and tracers variables
16   USE trc                 ! ocean passive tracers variables
17   USE trd_oce
18   USE trdtra
19   USE prtctl_trc          ! Print control for debbuging
20
21   IMPLICIT NONE
22   PRIVATE
23
24   PUBLIC trc_rad     
25   PUBLIC trc_rad_ini 
26
27   LOGICAL , PUBLIC ::   ln_trcrad           !: flag to artificially correct negative concentrations
28
29   !!----------------------------------------------------------------------
30   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
31   !! $Id$
32   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
33   !!----------------------------------------------------------------------
34   
35CONTAINS
36
37   SUBROUTINE trc_rad( kt )
38      !!----------------------------------------------------------------------
39      !!                  ***  ROUTINE trc_rad  ***
40      !!
41      !! ** Purpose :   "crappy" routine to correct artificial negative
42      !!              concentrations due to isopycnal scheme
43      !!
44      !! ** Method  : - PISCES or LOBSTER: Set negative concentrations to zero
45      !!                while computing the corresponding tracer content that
46      !!                is added to the tracers. Then, adjust the tracer
47      !!                concentration using a multiplicative factor so that
48      !!                the total tracer concentration is preserved.
49      !!              - CFC: simply set to zero the negative CFC concentration
50      !!                (the total CFC content is not strictly preserved)
51      !!----------------------------------------------------------------------
52      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index     
53      CHARACTER (len=22) :: charout
54      !!----------------------------------------------------------------------
55      !
56      IF( nn_timing == 1 )  CALL timing_start('trc_rad')
57      !
58      IF( kt == nittrc000 ) THEN
59         IF(lwp) WRITE(numout,*)
60         IF(lwp) WRITE(numout,*) 'trc_rad : Correct artificial negative concentrations '
61         IF(lwp) WRITE(numout,*) '~~~~~~~ '
62      ENDIF
63
64      IF( lk_cfc     )   CALL trc_rad_sms( kt, trb, trn, jp_cfc0 , jp_cfc1               )  ! CFC model
65      IF( lk_c14b    )   CALL trc_rad_sms( kt, trb, trn, jp_c14b0, jp_c14b1              )  ! bomb C14
66      IF( lk_pisces  )   CALL trc_rad_sms( kt, trb, trn, jp_pcs0 , jp_pcs1, cpreserv='Y' )  ! PISCES model
67      IF( lk_my_trc  )   CALL trc_rad_sms( kt, trb, trn, jp_myt0 , jp_myt1               )  ! MY_TRC model
68
69      !
70      IF(ln_ctl) THEN      ! print mean trends (used for debugging)
71         WRITE(charout, FMT="('rad')")
72         CALL prt_ctl_trc_info( charout )
73         CALL prt_ctl_trc( tab4d=trn, mask=tmask, clinfo=ctrcnm )
74      ENDIF
75      !
76      IF( nn_timing == 1 )  CALL timing_stop('trc_rad')
77      !
78   END SUBROUTINE trc_rad
79
80   SUBROUTINE trc_rad_ini
81      !!---------------------------------------------------------------------
82      !!                  ***  ROUTINE trc _rad_ini ***
83      !!
84      !! ** Purpose : read  namelist options
85      !!----------------------------------------------------------------------
86      INTEGER ::  ios                 ! Local integer output status for namelist read
87      NAMELIST/namtrc_rad/ ln_trcrad
88      !!----------------------------------------------------------------------
89
90      !
91      REWIND( numnat_ref )              ! namtrc_rad in reference namelist
92      READ  ( numnat_ref, namtrc_rad, IOSTAT = ios, ERR = 907)
93907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_rad in reference namelist', lwp )
94
95      REWIND( numnat_cfg )              ! namtrc_rad in configuration namelist
96      READ  ( numnat_cfg, namtrc_rad, IOSTAT = ios, ERR = 908 )
97908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namtrc_rad in configuration namelist', lwp )
98      IF(lwm) WRITE ( numont, namtrc_rad )
99
100      IF(lwp) THEN                     !   ! Control print
101         WRITE(numout,*)
102         WRITE(numout,*) '   Namelist namtrc_rad : treatment of negative concentrations'
103         WRITE(numout,*) '      correct artificially negative concen. or not ln_trcrad = ', ln_trcrad
104      ENDIF
105      !
106   END SUBROUTINE trc_rad_ini
107
108   SUBROUTINE trc_rad_sms( kt, ptrb, ptrn, jp_sms0, jp_sms1, cpreserv )
109      !!-----------------------------------------------------------------------------
110      !!                  ***  ROUTINE trc_rad_sms  ***
111      !!
112      !! ** Purpose :   "crappy" routine to correct artificial negative
113      !!              concentrations due to isopycnal scheme
114      !!
115      !! ** Method  : 2 cases :
116      !!                - Set negative concentrations to zero while computing
117      !!                  the corresponding tracer content that is added to the
118      !!                  tracers. Then, adjust the tracer concentration using
119      !!                  a multiplicative factor so that the total tracer
120      !!                  concentration is preserved.
121      !!                - simply set to zero the negative CFC concentration
122      !!                  (the total content of concentration is not strictly preserved)
123      !!--------------------------------------------------------------------------------
124      !! Arguments
125      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
126      INTEGER  , INTENT( in ) ::  &
127         jp_sms0, &       !: First index of the passive tracer model
128         jp_sms1          !: Last  index of  the passive tracer model
129
130      REAL(wp), DIMENSION (jpi,jpj,jpk,jptra), INTENT( inout )  :: &
131         ptrb, ptrn       !: before and now traceur concentration
132
133      CHARACTER( len = 1) , INTENT(in), OPTIONAL  :: &
134         cpreserv          !: flag to preserve content or not
135     
136      ! Local declarations
137      INTEGER  :: ji, jj, jk, jn     ! dummy loop indices
138      REAL(wp) :: ztrcorb, ztrmasb   ! temporary scalars
139      REAL(wp) :: zcoef, ztrcorn, ztrmasn   !    "         "
140      REAL(wp), POINTER, DIMENSION(:,:,:) ::   ztrtrdb, ztrtrdn   ! workspace arrays
141      REAL(wp) :: zs2rdt
142      LOGICAL ::   lldebug = .FALSE.
143      !!----------------------------------------------------------------------
144
145 
146      IF( l_trdtrc )  CALL wrk_alloc( jpi, jpj, jpk, ztrtrdb, ztrtrdn )
147     
148      IF( PRESENT( cpreserv )  ) THEN   !  total tracer concentration is preserved
149     
150         DO jn = jp_sms0, jp_sms1
151            !                                                        ! ===========
152            ztrcorb = 0.e0   ;   ztrmasb = 0.e0
153            ztrcorn = 0.e0   ;   ztrmasn = 0.e0
154
155            IF( l_trdtrc ) THEN
156               ztrtrdb(:,:,:) = ptrb(:,:,:,jn)                        ! save input trb for trend computation
157               ztrtrdn(:,:,:) = ptrn(:,:,:,jn)                        ! save input trn for trend computation
158            ENDIF
159            !                                                         ! sum over the global domain
160            ztrcorb = glob_sum( MIN( 0., ptrb(:,:,:,jn) ) * cvol(:,:,:) )
161            ztrcorn = glob_sum( MIN( 0., ptrn(:,:,:,jn) ) * cvol(:,:,:) )
162
163            ztrmasb = glob_sum( MAX( 0., ptrb(:,:,:,jn) ) * cvol(:,:,:) )
164            ztrmasn = glob_sum( MAX( 0., ptrn(:,:,:,jn) ) * cvol(:,:,:) )
165
166            IF( ztrcorb /= 0 ) THEN
167               zcoef = 1. + ztrcorb / ztrmasb
168               DO jk = 1, jpkm1
169                  ptrb(:,:,jk,jn) = MAX( 0., ptrb(:,:,jk,jn) )
170                  ptrb(:,:,jk,jn) = ptrb(:,:,jk,jn) * zcoef * tmask(:,:,jk)
171               END DO
172            ENDIF
173
174            IF( ztrcorn /= 0 ) THEN
175               zcoef = 1. + ztrcorn / ztrmasn
176               DO jk = 1, jpkm1
177                  ptrn(:,:,jk,jn) = MAX( 0., ptrn(:,:,jk,jn) )
178                  ptrn(:,:,jk,jn) = ptrn(:,:,jk,jn) * zcoef * tmask(:,:,jk)
179               END DO
180            ENDIF
181            !
182            IF( l_trdtrc ) THEN
183               !
184               zs2rdt = 1. / ( 2. * rdt )
185               ztrtrdb(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrdb(:,:,:) ) * zs2rdt
186               ztrtrdn(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrdn(:,:,:) ) * zs2rdt 
187               CALL trd_tra( kt, 'TRC', jn, jptra_radb, ztrtrdb )       ! Asselin-like trend handling
188               CALL trd_tra( kt, 'TRC', jn, jptra_radn, ztrtrdn )       ! standard     trend handling
189              !
190            ENDIF
191
192         END DO
193         !
194         !
195      ELSE  ! total CFC content is not strictly preserved
196
197         DO jn = jp_sms0, jp_sms1 
198
199           IF( l_trdtrc ) THEN
200              ztrtrdb(:,:,:) = ptrb(:,:,:,jn)                        ! save input trb for trend computation
201              ztrtrdn(:,:,:) = ptrn(:,:,:,jn)                        ! save input trn for trend computation
202           ENDIF
203
204            DO jk = 1, jpkm1
205               DO jj = 1, jpj
206                  DO ji = 1, jpi
207                     ptrn(ji,jj,jk,jn) = MAX( 0. , ptrn(ji,jj,jk,jn) )
208                     ptrb(ji,jj,jk,jn) = MAX( 0. , ptrb(ji,jj,jk,jn) )
209                  END DO
210               END DO
211            END DO
212         
213            IF( l_trdtrc ) THEN
214               !
215               zs2rdt = 1. / ( 2. * rdt * FLOAT( nn_dttrc ) )
216               ztrtrdb(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrdb(:,:,:) ) * zs2rdt
217               ztrtrdn(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrdn(:,:,:) ) * zs2rdt 
218               CALL trd_tra( kt, 'TRC', jn, jptra_radb, ztrtrdb )       ! Asselin-like trend handling
219               CALL trd_tra( kt, 'TRC', jn, jptra_radn, ztrtrdn )       ! standard     trend handling
220              !
221            ENDIF
222            !
223         ENDDO
224
225      ENDIF
226
227      IF( l_trdtrc )  CALL wrk_dealloc( jpi, jpj, jpk, ztrtrdb, ztrtrdn )
228
229   END SUBROUTINE trc_rad_sms
230#else
231   !!----------------------------------------------------------------------
232   !!   Dummy module :                                         NO TOP model
233   !!----------------------------------------------------------------------
234CONTAINS
235   SUBROUTINE trc_rad( kt )              ! Empty routine
236      INTEGER, INTENT(in) ::   kt
237      WRITE(*,*) 'trc_rad: You should not have seen this print! error?', kt
238   END SUBROUTINE trc_rad
239#endif
240   
241   !!======================================================================
242END MODULE trcrad
Note: See TracBrowser for help on using the repository browser.