source: branches/NERC/dev_r5518_NOC_MEDUSA_Stable/NEMOGCM/NEMO/TOP_SRC/CFC/trcsms_cfc.F90 @ 6201

Last change on this file since 6201 was 6201, checked in by jpalmier, 6 years ago

JPALM -- 04-01-2016 -- add debugg prints under debugg_key

File size: 13.6 KB
Line 
1MODULE trcsms_cfc
2   !!======================================================================
3   !!                      ***  MODULE trcsms_cfc  ***
4   !! TOP : CFC main model
5   !!======================================================================
6   !! History :  OPA  !  1999-10  (JC. Dutay)  original code
7   !!  NEMO      1.0  !  2004-03  (C. Ethe) free form + modularity
8   !!            2.0  !  2007-12  (C. Ethe, G. Madec)  reorganisation
9   !!----------------------------------------------------------------------
10#if defined key_cfc
11   !!----------------------------------------------------------------------
12   !!   'key_cfc'                                               CFC tracers
13   !!----------------------------------------------------------------------
14   !!   trc_sms_cfc  :  compute and add CFC suface forcing to CFC trends
15   !!   cfc_init     :  sets constants for CFC surface forcing computation
16   !!----------------------------------------------------------------------
17   USE dom_oce       ! ocean space and time domain
18   USE oce_trc       ! Ocean variables
19   USE par_trc       ! TOP parameters
20   USE trc           ! TOP variables
21   USE trd_oce
22   USE trdtrc
23   USE iom           ! I/O library
24
25   IMPLICIT NONE
26   PRIVATE
27
28   PUBLIC   trc_sms_cfc         ! called in ???   
29   PUBLIC   trc_sms_cfc_alloc   ! called in trcini_cfc.F90
30
31   INTEGER , PUBLIC, PARAMETER ::   jphem  =   2   ! parameter for the 2 hemispheres
32   INTEGER , PUBLIC            ::   jpyear         ! Number of years read in CFC1112 file
33   INTEGER , PUBLIC            ::   ndate_beg      ! initial calendar date (aammjj) for CFC
34   INTEGER , PUBLIC            ::   nyear_res      ! restoring time constant (year)
35   INTEGER , PUBLIC            ::   nyear_beg      ! initial year (aa)
36   
37   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   p_cfc    ! partial hemispheric pressure for CFC
38   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   xphem    ! spatial interpolation factor for patm
39   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qtr_cfc  ! flux at surface
40   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qint_cfc ! cumulative flux
41   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   patm     ! atmospheric function
42
43   REAL(wp), DIMENSION(4,2) ::   soa   ! coefficient for solubility of CFC [mol/l/atm]
44   REAL(wp), DIMENSION(3,2) ::   sob   !    "               "
45   REAL(wp), DIMENSION(4,2) ::   sca   ! coefficients for schmidt number in degre Celcius
46     
47   !                          ! coefficients for conversion
48   REAL(wp) ::   xconv1 = 1.0          ! conversion from to
49   REAL(wp) ::   xconv2 = 0.01/3600.   ! conversion from cm/h to m/s:
50   REAL(wp) ::   xconv3 = 1.0e+3       ! conversion from mol/l/atm to mol/m3/atm
51   REAL(wp) ::   xconv4 = 1.0e-12      ! conversion from mol/m3/atm to mol/m3/pptv
52
53   !! * Substitutions
54#  include "top_substitute.h90"
55   !!----------------------------------------------------------------------
56   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
57   !! $Id$
58   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
59   !!----------------------------------------------------------------------
60CONTAINS
61
62   SUBROUTINE trc_sms_cfc( kt )
63      !!----------------------------------------------------------------------
64      !!                     ***  ROUTINE trc_sms_cfc  ***
65      !!
66      !! ** Purpose :   Compute the surface boundary contition on CFC 11
67      !!             passive tracer associated with air-mer fluxes and add it
68      !!             to the general trend of tracers equations.
69      !!
70      !! ** Method  : - get the atmospheric partial pressure - given in pico -
71      !!              - computation of solubility ( in 1.e-12 mol/l then in 1.e-9 mol/m3)
72      !!              - computation of transfert speed ( given in cm/hour ----> cm/s )
73      !!              - the input function is given by :
74      !!                speed * ( concentration at equilibrium - concentration at surface )
75      !!              - the input function is in pico-mol/m3/s and the
76      !!                CFC concentration in pico-mol/m3
77      !!----------------------------------------------------------------------
78      !
79      INTEGER, INTENT(in) ::   kt    ! ocean time-step index
80      !
81      INTEGER  ::   ji, jj, jn, jl, jm, js
82      INTEGER  ::   iyear_beg, iyear_end
83      INTEGER  ::   im1, im2, ierr
84      REAL(wp) ::   ztap, zdtap       
85      REAL(wp) ::   zt1, zt2, zt3, zv2
86      REAL(wp) ::   zsol      ! solubility
87      REAL(wp) ::   zsch      ! schmidt number
88      REAL(wp) ::   zpp_cfc   ! atmospheric partial pressure of CFC
89      REAL(wp) ::   zca_cfc   ! concentration at equilibrium
90      REAL(wp) ::   zak_cfc   ! transfert coefficients
91      REAL(wp), ALLOCATABLE, DIMENSION(:,:)  ::   zpatm     ! atmospheric function
92      !!----------------------------------------------------------------------
93      !
94      !
95      IF( nn_timing == 1 )  CALL timing_start('trc_sms_cfc')
96      !
97      ALLOCATE( zpatm(jphem,jp_cfc), STAT=ierr )
98      IF( ierr > 0 ) THEN
99         CALL ctl_stop( 'trc_sms_cfc: unable to allocate zpatm array' )   ;   RETURN
100      ENDIF
101
102      IF( kt == nittrc000 )   CALL cfc_init
103
104      ! Temporal interpolation
105      ! ----------------------
106      iyear_beg = nyear - 1900
107      IF ( nmonth <= 6 ) THEN
108         iyear_beg = iyear_beg - 1
109         im1       =  6 - nmonth + 1
110         im2       =  6 + nmonth - 1
111      ELSE
112         im1       = 12 - nmonth + 7
113         im2       =      nmonth - 7
114      ENDIF
115      iyear_end = iyear_beg + 1
116
117      !                                                  !------------!
118      DO jl = 1, jp_cfc                                  !  CFC loop  !
119         !                                               !------------!
120         jn = jp_cfc0 + jl - 1
121         ! time interpolation at time kt
122         DO jm = 1, jphem
123            zpatm(jm,jl) = (  p_cfc(iyear_beg, jm, jl) * FLOAT (im1)  &
124               &           +  p_cfc(iyear_end, jm, jl) * FLOAT (im2) ) / 12.
125         END DO
126         
127         !                                                         !------------!
128         DO jj = 1, jpj                                            !  i-j loop  !
129            DO ji = 1, jpi                                         !------------!
130 
131               ! space interpolation
132               zpp_cfc  =       xphem(ji,jj)   * zpatm(1,jl)   &
133                  &     + ( 1.- xphem(ji,jj) ) * zpatm(2,jl)
134
135               ! Computation of concentration at equilibrium : in picomol/l
136               ! coefficient for solubility for CFC-11/12 in  mol/l/atm
137               IF( tmask(ji,jj,1) .GE. 0.5 ) THEN
138                  ztap  = ( tsn(ji,jj,1,jp_tem) + 273.16 ) * 0.01
139                  zdtap = sob(1,jl) + ztap * ( sob(2,jl) + ztap * sob(3,jl) ) 
140                  zsol  =  EXP( soa(1,jl) + soa(2,jl) / ztap + soa(3,jl) * LOG( ztap )   &
141                     &                    + soa(4,jl) * ztap * ztap + tsn(ji,jj,1,jp_sal) * zdtap ) 
142               ELSE
143                  zsol  = 0.e0
144               ENDIF
145               ! conversion from mol/l/atm to mol/m3/atm and from mol/m3/atm to mol/m3/pptv   
146               zsol = xconv4 * xconv3 * zsol * tmask(ji,jj,1) 
147               ! concentration at equilibrium
148               zca_cfc = xconv1 * zpp_cfc * zsol * tmask(ji,jj,1)             
149 
150               ! Computation of speed transfert
151               !    Schmidt number
152               zt1  = tsn(ji,jj,1,jp_tem)
153               zt2  = zt1 * zt1 
154               zt3  = zt1 * zt2
155               zsch = sca(1,jl) + sca(2,jl) * zt1 + sca(3,jl) * zt2 + sca(4,jl) * zt3
156
157               !    speed transfert : formulae of wanninkhof 1992
158               zv2     = wndm(ji,jj) * wndm(ji,jj)
159               zsch    = zsch / 660.
160               zak_cfc = ( 0.39 * xconv2 * zv2 / SQRT(zsch) ) * tmask(ji,jj,1)
161
162               ! Input function  : speed *( conc. at equil - concen at surface )
163               ! trn in pico-mol/l idem qtr; ak in en m/a
164               qtr_cfc(ji,jj,jl) = -zak_cfc * ( trb(ji,jj,1,jn) - zca_cfc )   &
165#if defined key_degrad
166                  &                         * facvol(ji,jj,1)                           &
167#endif
168                  &                         * tmask(ji,jj,1) * ( 1. - fr_i(ji,jj) )
169               ! Add the surface flux to the trend
170               tra(ji,jj,1,jn) = tra(ji,jj,1,jn) + qtr_cfc(ji,jj,jl) / fse3t(ji,jj,1) 
171
172               ! cumulation of surface flux at each time step
173               qint_cfc(ji,jj,jl) = qint_cfc(ji,jj,jl) + qtr_cfc(ji,jj,jl) * rdt
174               !                                               !----------------!
175            END DO                                             !  end i-j loop  !
176         END DO                                                !----------------!
177         !                                                  !----------------!
178      END DO                                                !  end CFC loop  !
179         !
180      IF( kt == nittrc000 ) THEN
181         DO jl = 1, jp_cfc   
182             WRITE(NUMOUT,*) ' '
183             WRITE(NUMOUT,*) 'CFC interpolation verification '  !! Jpalm 
184             WRITE(NUMOUT,*) '################################## '
185             WRITE(NUMOUT,*) ' '
186               if (jl.EQ.1) then
187                   WRITE(NUMOUT,*) 'Traceur = CFC11: '
188               elseif (jl.EQ.2) then
189                   WRITE(NUMOUT,*) 'Traceur = CFC12: '
190               endif
191             WRITE(NUMOUT,*) 'nyear    = ', nyear
192             WRITE(NUMOUT,*) 'nmonth   = ', nmonth
193             WRITE(NUMOUT,*) 'iyear_beg= ', iyear_beg
194             WRITE(NUMOUT,*) 'iyear_end= ', iyear_end
195             WRITE(NUMOUT,*) 'p_cfc(iyear_beg)= ',p_cfc(iyear_beg, 1, jl)
196             WRITE(NUMOUT,*) 'p_cfc(iyear_end)= ',p_cfc(iyear_end, 1, jl)
197             WRITE(NUMOUT,*) 'Im1= ',im1
198             WRITE(NUMOUT,*) 'Im2= ',im2
199             WRITE(NUMOUT,*) 'zpp_cfc = ',zpp_cfc
200             WRITE(NUMOUT,*) ' '
201         END DO 
202# if defined key_debug_medusa
203         CALL flush(numout)
204# endif
205      ENDIF
206        !
207      IF( lrst_trc ) THEN
208         IF(lwp) WRITE(numout,*)
209         IF(lwp) WRITE(numout,*) 'trc_sms_cfc : cumulated input function fields written in ocean restart file ',   &
210            &                    'at it= ', kt,' date= ', ndastp
211         IF(lwp) WRITE(numout,*) '~~~~'
212         DO jn = jp_cfc0, jp_cfc1
213            CALL iom_rstput( kt, nitrst, numrtw, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jn) )
214         END DO
215      ENDIF                                           
216      !
217      IF( lk_iomput ) THEN
218         CALL iom_put( "qtrCFC11"  , qtr_cfc (:,:,1) )
219         CALL iom_put( "qintCFC11" , qint_cfc(:,:,1) )
220      ELSE
221         IF( ln_diatrc ) THEN
222            trc2d(:,:,jp_cfc0_2d    ) = qtr_cfc (:,:,1)
223            trc2d(:,:,jp_cfc0_2d + 1) = qint_cfc(:,:,1)
224         END IF
225      END IF
226      !
227      IF( l_trdtrc ) THEN
228          DO jn = jp_cfc0, jp_cfc1
229            CALL trd_trc( tra(:,:,:,jn), jn, jptra_sms, kt )   ! save trends
230          END DO
231      END IF
232      !
233      IF( nn_timing == 1 )  CALL timing_stop('trc_sms_cfc')
234      !
235   END SUBROUTINE trc_sms_cfc
236
237
238   SUBROUTINE cfc_init
239      !!---------------------------------------------------------------------
240      !!                     ***  cfc_init  *** 
241      !!
242      !! ** Purpose : sets constants for CFC model
243      !!---------------------------------------------------------------------
244      INTEGER :: jn
245
246      ! coefficient for CFC11
247      !----------------------
248
249      ! Solubility
250      soa(1,1) = -229.9261 
251      soa(2,1) =  319.6552
252      soa(3,1) =  119.4471
253      soa(4,1) =  -1.39165
254
255      sob(1,1) =  -0.142382
256      sob(2,1) =   0.091459
257      sob(3,1) =  -0.0157274
258
259      ! Schmidt number
260      sca(1,1) = 3501.8
261      sca(2,1) = -210.31
262      sca(3,1) =  6.1851
263      sca(4,1) = -0.07513
264
265      ! coefficient for CFC12
266      !----------------------
267
268      ! Solubility
269      soa(1,2) = -218.0971
270      soa(2,2) =  298.9702
271      soa(3,2) =  113.8049
272      soa(4,2) =  -1.39165
273
274      sob(1,2) =  -0.143566
275      sob(2,2) =   0.091015
276      sob(3,2) =  -0.0153924
277
278      ! schmidt number
279      sca(1,2) =  3845.4 
280      sca(2,2) =  -228.95
281      sca(3,2) =  6.1908 
282      sca(4,2) =  -0.067430
283
284      IF( ln_rsttr ) THEN
285         IF(lwp) WRITE(numout,*)
286         IF(lwp) WRITE(numout,*) ' Read specific variables from CFC model '
287         IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~'
288         !
289         DO jn = jp_cfc0, jp_cfc1
290            CALL iom_get( numrtr, jpdom_autoglo, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jn) ) 
291         END DO
292      ENDIF
293      IF(lwp) WRITE(numout,*)
294      !
295   END SUBROUTINE cfc_init
296
297
298   INTEGER FUNCTION trc_sms_cfc_alloc()
299      !!----------------------------------------------------------------------
300      !!                     ***  ROUTINE trc_sms_cfc_alloc  ***
301      !!----------------------------------------------------------------------
302      ALLOCATE( xphem   (jpi,jpj)        ,     &
303         &      qtr_cfc (jpi,jpj,jp_cfc) ,     &
304         &      qint_cfc(jpi,jpj,jp_cfc) , STAT=trc_sms_cfc_alloc )
305         !
306      IF( trc_sms_cfc_alloc /= 0 ) CALL ctl_warn('trc_sms_cfc_alloc : failed to allocate arrays.')
307      !
308   END FUNCTION trc_sms_cfc_alloc
309
310#else
311   !!----------------------------------------------------------------------
312   !!   Dummy module                                         No CFC tracers
313   !!----------------------------------------------------------------------
314CONTAINS
315   SUBROUTINE trc_sms_cfc( kt )       ! Empty routine
316      WRITE(*,*) 'trc_sms_cfc: You should not have seen this print! error?', kt
317   END SUBROUTINE trc_sms_cfc
318#endif
319
320   !!======================================================================
321END MODULE trcsms_cfc
Note: See TracBrowser for help on using the repository browser.