source: NEMO/branches/2020/dev_r12563_ASINTER-06_ABL_improvement/src/TOP/C14/trcsms_c14.F90 @ 12808

Last change on this file since 12808 was 12489, checked in by davestorkey, 11 months ago

Preparation for new timestepping scheme #2390.
Main changes:

  1. Initial euler timestep now handled in stp and not in TRA/DYN routines.
  2. Renaming of all timestep parameters. In summary, the namelist parameter is now rn_Dt and the current timestep is rDt (and rDt_ice, rDt_trc etc).
  3. Renaming of a few miscellaneous parameters, eg. atfp → rn_atfp (namelist parameter used everywhere) and rau0 → rho0.

This version gives bit-comparable results to the previous version of the trunk.

  • Property svn:keywords set to Id
File size: 7.1 KB
Line 
1MODULE trcsms_c14
2   !!======================================================================
3   !!                      ***  MODULE trcsms_c14  ***
4   !! TOP : Bomb C14 main module
5   !!======================================================================
6   !! History     -   ! 1994-05 ( J. Orr ) original code
7   !!            1.0  ! 2006-02 ( J.M. Molines )  Free form + modularity
8   !!            2.0  ! 2008-12 ( C. Ethe ) reorganisation
9   !!            4.0  ! 2011-02 ( A.R. Porter, STFC Daresbury ) Dynamic memory
10   !!                 ! 2015    (A. Mouchet) general C14 + update formulas
11   !!----------------------------------------------------------------------
12   !!   trc_sms_c14 :  compute and add C14 suface forcing to C14 trends
13   !!----------------------------------------------------------------------
14   USE oce_trc       ! Ocean variables
15   USE par_trc       ! TOP parameters
16   USE trc           ! TOP variables
17   USE trd_oce    ! trends
18   USE trdtrc    ! trends
19   USE sms_c14    ! atmospheric forcing
20   USE trcatm_c14    ! atmospheric forcing
21   USE iom
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC   trc_sms_c14       ! called in trcsms.F90
27
28   !! * Substitutions
29#  include "do_loop_substitute.h90"
30   !!----------------------------------------------------------------------
31   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
32   !! $Id$
33   !! Software governed by the CeCILL license (see ./LICENSE)
34   !!----------------------------------------------------------------------
35CONTAINS
36
37   SUBROUTINE trc_sms_c14( kt, Kbb, Kmm, Krhs )
38      !!----------------------------------------------------------------------
39      !!                  ***  ROUTINE trc_sms_c14  ***
40      !!
41      !! ** Purpose :   Compute the surface boundary contition on C14
42      !!      passive tracer associated with air-sea fluxes and add it to
43      !!      the general trend of tracers equations.
44      !
45      !    Method:
46      !          - transport the ratio C14/C as in Toggweiler et al. (JGR,1989)
47      !          - if on-line a passive tracer (jpcref; NO sms) allows compensating for
48      !            freshwater fluxes which should not impact the C14/C ratio
49      !
50      !        =>   Delta-C14= ( tr(...jp_c14,Kmm) -1)*1000.
51      !!
52      !!----------------------------------------------------------------------
53      !
54      INTEGER, INTENT(in) ::   kt               ! ocean time-step index
55      INTEGER, INTENT(in) ::   Kbb, Kmm, Krhs   ! ocean time level
56      !
57      INTEGER  :: ji, jj, jk        ! dummy loop indices
58      REAL(wp) :: zt, ztp, zsk      ! dummy variables
59      REAL(wp) :: zsol              ! solubility
60      REAL(wp) :: zsch              ! schmidt number
61      REAL(wp) :: zv2               ! wind speed ( square)
62      REAL(wp) :: zpv               ! piston velocity
63      !!----------------------------------------------------------------------
64      !
65      IF( ln_timing )   CALL timing_start('trc_sms_c14')
66      !
67      IF( kt ==  nittrc000 ) THEN
68         IF(lwp) WRITE(numout,*)
69         IF(lwp) WRITE(numout,*) ' trc_sms_c14:  C14 model'
70         IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~'
71      ENDIF
72      !
73      ! Get co2sbc & c14sbc(ji,jj): at 1st iter for all, at each time step for transient
74      IF( kc14typ >= 1 .OR.  kt ==  nittrc000 )   CALL trc_atm_c14( kt, co2sbc, c14sbc ) 
75
76      ! -------------------------------------------------------------------
77      !  Gas exchange coefficient (Wanninkhof, 1992, JGR, 97,7373-7382)
78      !  Schmidt number of CO2 in seawater (Wanninkhof, 1992 & 2014)
79      !  CO2 solubility (Weiss, 1974; Wanninkhof, 2014)
80      ! -------------------------------------------------------------------
81
82      DO_2D_11_11
83         IF( tmask(ji,jj,1) >  0. ) THEN
84            !
85            zt   = MIN( 40. , ts(ji,jj,1,jp_tem,Kmm) )
86            !
87            !  Computation of solubility zsol in [mol/(L * atm)]
88            !   after Wanninkhof (2014) referencing Weiss (1974)
89            ztp  = ( zt + 273.16 ) * 0.01
90            zsk  = 0.027766 + ztp * ( -0.025888 + 0.0050578 * ztp )   ! [mol/(L * atm)]
91            zsol = EXP( -58.0931 + 90.5069 / ztp  + 22.2940 * LOG( ztp ) + zsk * ts(ji,jj,1,jp_sal,Kmm) )
92            ! convert solubilities [mol/(L * atm)] -> [mol/(m^3 * ppm)]
93            zsol = zsol * 1.e-03
94
95            ! Computes the Schmidt number of CO2 in seawater
96            !               Wanninkhof-2014
97            zsch = 2116.8 + zt * ( -136.25 + zt * (4.7353 + zt * (-0.092307 + 0.0007555 * zt ) ) )
98
99            ! Wanninkhof Piston velocity: zpv in units [m/s]
100            zv2 = xkwind * (wndm(ji,jj) * wndm(ji,jj))              ! wind speed module at T points
101            ! chemical enhancement (Wanninkhof & Knox, 1996)
102            IF( ln_chemh ) zv2 = zv2 + 2.5 * ( 0.5246 + zt * (0.016256 + 0.00049946  * zt ) )
103            zv2 = zv2/360000._wp                                    ! conversion cm/h -> m/s
104            !
105            zpv  = ( zv2 * SQRT( 660./ zsch ) ) * ( 1. - fr_i(ji,jj) ) * tmask(ji,jj,1)
106
107            ! CO2 piston velocity (m/s)
108            exch_co2(ji,jj)= zpv
109            ! CO2 invasion rate (mol/ppm/m2/s) = 1st part of 14C/C exchange velocity
110            exch_c14(ji,jj)= zpv * zsol
111         ELSE
112            exch_co2(ji,jj) = 0._wp
113            exch_c14(ji,jj) = 0._wp
114         ENDIF
115      END_2D
116
117      ! Exchange velocity for 14C/C ratio (m/s)
118      zt = co2sbc / xdicsur
119      exch_c14(:,:) = zt * exch_c14(:,:)
120      !
121      ! Flux of C-14 from air-to-sea; units: (C14/C ratio) x m/s
122      !                               already masked
123      qtr_c14(:,:) = exch_c14(:,:) * ( c14sbc(:,:) - tr(:,:,1,jp_c14,Kbb) )
124           
125      ! cumulation of air-to-sea flux at each time step
126      qint_c14(:,:) = qint_c14(:,:) + qtr_c14(:,:) * rn_Dt
127      !
128      ! Add the surface flux to the trend of jp_c14
129      DO_2D_11_11
130         tr(ji,jj,1,jp_c14,Krhs) = tr(ji,jj,1,jp_c14,Krhs) + qtr_c14(ji,jj) / e3t(ji,jj,1,Kmm) 
131      END_2D
132      !
133      ! Computation of decay effects on jp_c14
134      DO_3D_11_11( 1, jpk )
135         !
136         tr(ji,jj,jk,jp_c14,Krhs) = tr(ji,jj,jk,jp_c14,Krhs) - rlam14 * tr(ji,jj,jk,jp_c14,Kbb) * tmask(ji,jj,jk) 
137         !
138      END_3D
139      !
140      IF( lrst_trc ) THEN
141         IF(lwp) WRITE(numout,*)
142         IF(lwp) WRITE(numout,*) ' trc_rst_wri_c14 : Write specific variables from c14 model '
143         IF(lwp) WRITE(numout,*) ' ~~~~~~~~~~~~~~'
144         !
145         CALL iom_rstput( kt, nitrst, numrtw, 'co2sbc', co2sbc )       ! These five need      &
146         CALL iom_rstput( kt, nitrst, numrtw, 'c14sbc', c14sbc )     ! &    to be written   &
147         CALL iom_rstput( kt, nitrst, numrtw, 'exch_co2', exch_co2 ) ! &    for temporal    &
148         CALL iom_rstput( kt, nitrst, numrtw, 'exch_c14', exch_c14 ) ! &    averages        &
149         CALL iom_rstput( kt, nitrst, numrtw, 'qtr_c14', qtr_c14 )   ! &    to be coherent.
150         CALL iom_rstput( kt, nitrst, numrtw, 'qint_c14', qint_c14 ) ! Cumulative
151         !
152      ENDIF
153
154      IF( l_trdtrc )  CALL trd_trc( tr(:,:,:,jp_c14,Krhs), 1, jptra_sms, kt, Kmm )   ! save trends
155      !
156      IF( ln_timing )   CALL timing_stop('trc_sms_c14')
157      !
158   END SUBROUTINE trc_sms_c14
159
160  !!======================================================================
161END MODULE trcsms_c14
Note: See TracBrowser for help on using the repository browser.