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 trunk/NEMO/TOP_SRC/TRP – NEMO

source: trunk/NEMO/TOP_SRC/TRP/trcrad.F90 @ 1284

Last change on this file since 1284 was 1257, checked in by cetlod, 15 years ago

update trcrad to take into account the new C14 bomb tracer model, see ticket:298
  • Property svn:executable set to *
  • Property svn:keywords set to Id
File size: 10.2 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 trp_trc             ! ocean passive tracers variables
17   USE trdmld_trc
18   USE trdmld_trc_oce
19   USE lib_mpp
20   USE prtctl_trc          ! Print control for debbuging
21
22   IMPLICIT NONE
23   PRIVATE
24
25   PUBLIC trc_rad         ! routine called by trcstp.F90
26
27   !! * Substitutions
28#  include "top_substitute.h90"
29   !!----------------------------------------------------------------------
30   !! NEMO/TOP 1.0 , LOCEAN-IPSL (2005)
31   !! $Id$
32   !! Software governed by the CeCILL licence (modipsl/doc/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( kt == nittrc000 ) THEN
57         IF(lwp) WRITE(numout,*)
58         IF(lwp) WRITE(numout,*) 'trc_rad : Correct artificial negative concentrations '
59         IF(lwp) WRITE(numout,*) '~~~~~~~ '
60      ENDIF
61
62      IF( lk_cfc     )   CALL trc_rad_sms( kt, trb, trn, jp_cfc0 , jp_cfc1               )  ! CFC model
63      IF( lk_c14b    )   CALL trc_rad_sms( kt, trb, trn, jp_c14b0, jp_c14b1              )  ! bomb C14
64      IF( lk_lobster )   CALL trc_rad_sms( kt, trb, trn, jp_lob0 , jp_lob1, cpreserv='Y' )  ! LOBSTER model
65      IF( lk_pisces  )   CALL trc_rad_sms( kt, trb, trn, jp_pcs0 , jp_pcs1, cpreserv='Y' )  ! PISCES model
66      IF( lk_my_trc  )   CALL trc_rad_sms( kt, trb, trn, jp_myt0 , jp_myt1               )  ! MY_TRC model
67
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   END SUBROUTINE trc_rad
77
78   SUBROUTINE trc_rad_sms( kt, ptrb, ptrn, jp_sms0, jp_sms1, cpreserv )
79      !!-----------------------------------------------------------------------------
80      !!                  ***  ROUTINE trc_rad_sms  ***
81      !!
82      !! ** Purpose :   "crappy" routine to correct artificial negative
83      !!              concentrations due to isopycnal scheme
84      !!
85      !! ** Method  : 2 cases :
86      !!                - Set negative concentrations to zero while computing
87      !!                  the corresponding tracer content that is added to the
88      !!                  tracers. Then, adjust the tracer concentration using
89      !!                  a multiplicative factor so that the total tracer
90      !!                  concentration is preserved.
91      !!                - simply set to zero the negative CFC concentration
92      !!                  (the total content of concentration is not strictly preserved)
93      !!--------------------------------------------------------------------------------
94      !! Arguments
95      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
96      INTEGER  , INTENT( in ) ::  &
97         jp_sms0, &       !: First index of the passive tracer model
98         jp_sms1          !: Last  index of  the passive tracer model
99
100      REAL(wp), DIMENSION (jpi,jpj,jpk,jptra), INTENT( inout )  :: &
101         ptrb, ptrn       !: before and now traceur concentration
102
103      CHARACTER( len = 1) , INTENT(in), OPTIONAL  :: &
104         cpreserv          !: flag to preserve content or not
105     
106      ! Local declarations
107      INTEGER  ::  ji, jj, jk, jn     ! dummy loop indices
108      REAL(wp) :: zvolk, ztrcorb, ztrmasb   ! temporary scalars
109      REAL(wp) :: zcoef, ztrcorn, ztrmasn   !    "         "
110      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrtrdb  ! workspace arrays
111      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrtrdn  ! workspace arrays
112      REAL(wp) :: zs2rdt
113      LOGICAL ::   lldebug = .FALSE.
114
115      !!----------------------------------------------------------------------
116
117      IF( l_trdtrc ) THEN
118        !
119        ALLOCATE( ztrtrdb(jpi,jpj,jpk) )
120        ALLOCATE( ztrtrdn(jpi,jpj,jpk) )
121        !
122      ENDIF
123     
124      IF( PRESENT( cpreserv )  ) THEN   !  total tracer concentration is preserved
125     
126         DO jn = jp_sms0, jp_sms1
127         !                                                        ! ===========
128            ztrcorb = 0.e0   ;   ztrmasb = 0.e0
129            ztrcorn = 0.e0   ;   ztrmasn = 0.e0
130
131           IF( l_trdtrc ) THEN
132              ztrtrdb(:,:,:) = ptrb(:,:,:,jn)                        ! save input trb for trend computation
133              ztrtrdn(:,:,:) = ptrn(:,:,:,jn)                        ! save input trn for trend computation
134           ENDIF
135
136
137            DO jk = 1, jpkm1
138               DO jj = 1, jpj
139                  DO ji = 1, jpi
140                     zvolk  = cvol(ji,jj,jk)
141# if defined key_off_degrad
142                     zvolk  = zvolk * facvol(ji,jj,jk)
143# endif
144                     ztrcorb = ztrcorb + MIN( 0., ptrb(ji,jj,jk,jn) ) * zvolk
145                     ztrcorn = ztrcorn + MIN( 0., ptrn(ji,jj,jk,jn) ) * zvolk
146
147                     ptrb(ji,jj,jk,jn) = MAX( 0., ptrb(ji,jj,jk,jn) )
148                     ptrn(ji,jj,jk,jn) = MAX( 0., ptrn(ji,jj,jk,jn) )
149
150                     ztrmasb = ztrmasb + ptrb(ji,jj,jk,jn) * zvolk
151                     ztrmasn = ztrmasn + ptrn(ji,jj,jk,jn) * zvolk
152                  END DO
153               END DO
154            END DO
155
156            IF( lk_mpp ) THEN
157               CALL mpp_sum( ztrcorb )      ! sum over the global domain
158               CALL mpp_sum( ztrcorn )      ! sum over the global domain
159               CALL mpp_sum( ztrmasb )      ! sum over the global domain
160               CALL mpp_sum( ztrmasn )      ! sum over the global domain
161            ENDIF
162
163            IF( ztrcorb /= 0 ) THEN
164               zcoef = 1. + ztrcorb / ztrmasb
165               DO jk = 1, jpkm1
166                  ptrb(:,:,jk,jn) = ptrb(:,:,jk,jn) * zcoef * tmask(:,:,jk)
167               END DO
168            ENDIF
169
170            IF( ztrcorn /= 0 ) THEN
171               zcoef = 1. + ztrcorn / ztrmasn
172               DO jk = 1, jpkm1
173                  ptrn(:,:,jk,jn) = ptrn(:,:,jk,jn) * zcoef * tmask(:,:,jk)
174               END DO
175            ENDIF
176            !
177            IF( l_trdtrc ) THEN
178               !
179               zs2rdt = 1. / ( 2. * rdt )
180               ztrtrdb(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrdb(:,:,:) ) * zs2rdt
181               ztrtrdn(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrdn(:,:,:) ) * zs2rdt 
182               IF (luttrd(jn)) CALL trd_mod_trc( ztrtrdb, jn, jptrc_trd_radb, kt )       ! Asselin-like trend handling
183               IF (luttrd(jn)) CALL trd_mod_trc( ztrtrdn, jn, jptrc_trd_radn, kt )       ! standard     trend handling
184              !
185            ENDIF
186
187         END DO
188         !
189         !
190      ELSE  ! total CFC content is not strictly preserved
191
192         DO jn = jp_sms0, jp_sms1 
193
194           IF( l_trdtrc ) THEN
195              ztrtrdb(:,:,:) = ptrb(:,:,:,jn)                        ! save input trb for trend computation
196              ztrtrdn(:,:,:) = ptrn(:,:,:,jn)                        ! save input trn for trend computation
197           ENDIF
198
199            DO jk = 1, jpkm1
200               DO jj = 1, jpj
201                  DO ji = 1, jpi
202                     ptrn(ji,jj,jk,jn) = MAX( 0. , ptrn(ji,jj,jk,jn) )
203                     ptrb(ji,jj,jk,jn) = MAX( 0. , ptrb(ji,jj,jk,jn) )
204                  END DO
205               END DO
206            END DO
207         
208            IF( l_trdtrc ) THEN
209               !
210               zs2rdt = 1. / ( 2. * rdt * FLOAT(ndttrc) )
211               ztrtrdb(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrdb(:,:,:) ) * zs2rdt
212               ztrtrdn(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrdn(:,:,:) ) * zs2rdt 
213               IF (luttrd(jn)) CALL trd_mod_trc( ztrtrdb, jn, jptrc_trd_radb, kt )       ! Asselin-like trend handling
214               IF (luttrd(jn)) CALL trd_mod_trc( ztrtrdn, jn, jptrc_trd_radn, kt )       ! standard     trend handling
215              !
216            ENDIF
217            !
218         ENDDO
219
220      ENDIF
221
222      IF( l_trdtrc )   DEALLOCATE( ztrtrdb, ztrtrdn )
223
224   END SUBROUTINE trc_rad_sms
225#else
226   !!----------------------------------------------------------------------
227   !!   Dummy module :                                         NO TOP model
228   !!----------------------------------------------------------------------
229CONTAINS
230   SUBROUTINE trc_rad( kt )              ! Empty routine
231      INTEGER, INTENT(in) ::   kt
232      WRITE(*,*) 'trc_rad: You should not have seen this print! error?', kt
233   END SUBROUTINE trc_rad
234#endif
235   
236   !!======================================================================
237END MODULE trcrad
Note: See TracBrowser for help on using the repository browser.