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 @ 1533

Last change on this file since 1533 was 1257, checked in by cetlod, 16 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
RevLine 
[941]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
[1175]16   USE trp_trc             ! ocean passive tracers variables
17   USE trdmld_trc
18   USE trdmld_trc_oce
[941]19   USE lib_mpp
20   USE prtctl_trc          ! Print control for debbuging
21
22   IMPLICIT NONE
23   PRIVATE
24
[1257]25   PUBLIC trc_rad         ! routine called by trcstp.F90
[941]26
27   !! * Substitutions
28#  include "top_substitute.h90"
29   !!----------------------------------------------------------------------
30   !! NEMO/TOP 1.0 , LOCEAN-IPSL (2005)
[1146]31   !! $Id$
[941]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
[1257]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
[941]67
[1003]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
[1175]78   SUBROUTINE trc_rad_sms( kt, ptrb, ptrn, jp_sms0, jp_sms1, cpreserv )
[1003]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
[1175]95      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
[1003]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   !    "         "
[1175]110      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrtrdb  ! workspace arrays
111      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrtrdn  ! workspace arrays
[1257]112      REAL(wp) :: zs2rdt
[1175]113      LOGICAL ::   lldebug = .FALSE.
[1003]114
115      !!----------------------------------------------------------------------
116
[1175]117      IF( l_trdtrc ) THEN
118        !
119        ALLOCATE( ztrtrdb(jpi,jpj,jpk) )
120        ALLOCATE( ztrtrdn(jpi,jpj,jpk) )
121        !
122      ENDIF
[1003]123     
124      IF( PRESENT( cpreserv )  ) THEN   !  total tracer concentration is preserved
125     
126         DO jn = jp_sms0, jp_sms1
[1175]127         !                                                        ! ===========
[1193]128            ztrcorb = 0.e0   ;   ztrmasb = 0.e0
129            ztrcorn = 0.e0   ;   ztrmasn = 0.e0
[1003]130
[1193]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
[1175]135
136
[941]137            DO jk = 1, jpkm1
138               DO jj = 1, jpj
139                  DO ji = 1, jpi
[1257]140                     zvolk  = cvol(ji,jj,jk)
[941]141# if defined key_off_degrad
[1257]142                     zvolk  = zvolk * facvol(ji,jj,jk)
[941]143# endif
[1257]144                     ztrcorb = ztrcorb + MIN( 0., ptrb(ji,jj,jk,jn) ) * zvolk
145                     ztrcorn = ztrcorn + MIN( 0., ptrn(ji,jj,jk,jn) ) * zvolk
[941]146
[1257]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) )
[941]149
[1003]150                     ztrmasb = ztrmasb + ptrb(ji,jj,jk,jn) * zvolk
151                     ztrmasn = ztrmasn + ptrn(ji,jj,jk,jn) * zvolk
[941]152                  END DO
153               END DO
154            END DO
[1003]155
[941]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
[1003]166                  ptrb(:,:,jk,jn) = ptrb(:,:,jk,jn) * zcoef * tmask(:,:,jk)
[941]167               END DO
168            ENDIF
169
170            IF( ztrcorn /= 0 ) THEN
171               zcoef = 1. + ztrcorn / ztrmasn
172               DO jk = 1, jpkm1
[1003]173                  ptrn(:,:,jk,jn) = ptrn(:,:,jk,jn) * zcoef * tmask(:,:,jk)
[941]174               END DO
175            ENDIF
176            !
[1175]177            IF( l_trdtrc ) THEN
178               !
[1257]179               zs2rdt = 1. / ( 2. * rdt )
180               ztrtrdb(:,:,:) = ( ptrb(:,:,:,jn) - ztrtrdb(:,:,:) ) * zs2rdt
181               ztrtrdn(:,:,:) = ( ptrn(:,:,:,jn) - ztrtrdn(:,:,:) ) * zs2rdt 
[1175]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
[941]187         END DO
188         !
[1003]189         !
190      ELSE  ! total CFC content is not strictly preserved
191
192         DO jn = jp_sms0, jp_sms1 
[1257]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
[1003]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
[1257]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
[1003]219
[941]220      ENDIF
221
[1175]222      IF( l_trdtrc )   DEALLOCATE( ztrtrdb, ztrtrdn )
223
[1003]224   END SUBROUTINE trc_rad_sms
[941]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.