source: branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90 @ 3488

Last change on this file since 3488 was 3488, checked in by acc, 9 years ago

Branch: dev_r3385_NOCS04_HAMF; #665. Stage 3 of 2012 development: Rationalisation of code. Added LIM3 changes, corrected coupled changes and highlighted areas of concern in CICE interface

  • Property svn:keywords set to Id
File size: 13.4 KB
Line 
1MODULE trasbc
2   !!==============================================================================
3   !!                       ***  MODULE  trasbc  ***
4   !! Ocean active tracers:  surface boundary condition
5   !!==============================================================================
6   !! History :  OPA  !  1998-10  (G. Madec, G. Roullet, M. Imbard)  Original code
7   !!            8.2  !  2001-02  (D. Ludicone)  sea ice and free surface
8   !!  NEMO      1.0  !  2002-06  (G. Madec)  F90: Free form and module
9   !!            3.3  !  2010-04  (M. Leclair, G. Madec)  Forcing averaged over 2 time steps
10   !!             -   !  2010-09  (C. Ethe, G. Madec) Merge TRA-TRC
11   !!----------------------------------------------------------------------
12
13   !!----------------------------------------------------------------------
14   !!   tra_sbc      : update the tracer trend at ocean surface
15   !!----------------------------------------------------------------------
16   USE oce             ! ocean dynamics and active tracers
17   USE sbc_oce         ! surface boundary condition: ocean
18   USE dom_oce         ! ocean space domain variables
19   USE phycst          ! physical constant
20   USE traqsr          ! solar radiation penetration
21   USE trdmod_oce      ! ocean trends
22   USE trdtra          ! ocean trends
23   USE in_out_manager  ! I/O manager
24   USE prtctl          ! Print control
25   USE restart         ! ocean restart
26   USE sbcrnf          ! River runoff 
27   USE sbcmod          ! ln_rnf 
28   USE iom
29   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
30   USE wrk_nemo        ! Memory Allocation
31   USE timing          ! Timing
32
33   IMPLICIT NONE
34   PRIVATE
35
36   PUBLIC   tra_sbc    ! routine called by step.F90
37
38   !! * Substitutions
39#  include "domzgr_substitute.h90"
40#  include "vectopt_loop_substitute.h90"
41   !!----------------------------------------------------------------------
42   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
43   !! $Id$
44   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
45   !!----------------------------------------------------------------------
46CONTAINS
47
48   SUBROUTINE tra_sbc ( kt )
49      !!----------------------------------------------------------------------
50      !!                  ***  ROUTINE tra_sbc  ***
51      !!                   
52      !! ** Purpose :   Compute the tracer surface boundary condition trend of
53      !!      (flux through the interface, concentration/dilution effect)
54      !!      and add it to the general trend of tracer equations.
55      !!
56      !! ** Method :
57      !!      Following Roullet and Madec (2000), the air-sea flux can be divided
58      !!      into three effects: (1) Fext, external forcing;
59      !!      (2) Fwi, concentration/dilution effect due to water exchanged
60      !!         at the surface by evaporation, precipitations and runoff (E-P-R);
61      !!      (3) Fwe, tracer carried with the water that is exchanged.
62      !!            - salinity    : salt flux only due to freezing/melting
63      !!            sa = sa +  sfx / rau0 / e3t  for k=1
64      !!
65      !!      Fext, flux through the air-sea interface for temperature and salt:
66      !!            - temperature : heat flux q (w/m2). If penetrative solar
67      !!         radiation q is only the non solar part of the heat flux, the
68      !!         solar part is added in traqsr.F routine.
69      !!            ta = ta + q /(rau0 rcp e3t)  for k=1
70      !!            - salinity    : no salt flux
71      !!
72      !!      The formulation for Fwb and Fwi vary according to the free
73      !!      surface formulation (linear or variable volume).
74      !!      * Linear free surface
75      !!            The surface freshwater flux modifies the ocean volume
76      !!         and thus the concentration of a tracer and the temperature.
77      !!         First order of the effect of surface freshwater exchange
78      !!         for salinity, it can be neglected on temperature (especially
79      !!         as the temperature of precipitations and runoffs is usually
80      !!         unknown).
81      !!            - temperature : we assume that the temperature of both
82      !!         precipitations and runoffs is equal to the SST, thus there
83      !!         is no additional flux since in this case, the concentration
84      !!         dilution effect is balanced by the net heat flux associated
85      !!         to the freshwater exchange (Fwe+Fwi=0):
86      !!            (Tp P - Te E) + SST (P-E) = 0 when Tp=Te=SST
87      !!            - salinity    : evaporation, precipitation and runoff
88      !!         water has a zero salinity  but there is a salt flux due to
89      !!         freezing/melting, thus:
90      !!            sa = sa + emp * sn / rau0 / e3t   for k=1
91      !!                    + sfx    / rau0 / e3t
92      !!         where emp, the surface freshwater budget (evaporation minus
93      !!         precipitation minus runoff) given in kg/m2/s is divided
94      !!         by rau0 = 1020 kg/m3 (density of sea water) to obtain m/s.   
95      !!         Note: even though Fwe does not appear explicitly for
96      !!         temperature in this routine, the heat carried by the water
97      !!         exchanged through the surface is part of the total heat flux
98      !!         forcing and must be taken into account in the global heat
99      !!         balance).
100      !!      * nonlinear free surface (variable volume, lk_vvl)
101      !!         contrary to the linear free surface case, Fwi is properly
102      !!         taken into account by using the true layer thicknesses to       
103      !!         calculate tracer content and advection. There is no need to
104      !!         deal with it in this routine.
105      !!           - temperature: Fwe=SST (P-E+R) is added to Fext.
106      !!           - salinity:  Fwe = 0, there is no surface flux of salt.
107      !!
108      !! ** Action  : - Update the 1st level of (ta,sa) with the trend associated
109      !!                with the tracer surface boundary condition
110      !!              - save the trend it in ttrd ('key_trdtra')
111      !!----------------------------------------------------------------------
112      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
113      !!
114      INTEGER  ::   ji, jj, jk, jn           ! dummy loop indices 
115      REAL(wp) ::   zfact, z1_e3t, zdep
116      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrdt, ztrds
117      !!----------------------------------------------------------------------
118      !
119      IF( nn_timing == 1 )  CALL timing_start('tra_sbc')
120      !
121      IF( kt == nit000 ) THEN
122         IF(lwp) WRITE(numout,*)
123         IF(lwp) WRITE(numout,*) 'tra_sbc : TRAcer Surface Boundary Condition'
124         IF(lwp) WRITE(numout,*) '~~~~~~~ '
125      ENDIF
126
127      IF( l_trdtra )   THEN                    !* Save ta and sa trends
128         CALL wrk_alloc( jpi, jpj, jpk, ztrdt, ztrds ) 
129         ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
130         ztrds(:,:,:) = tsa(:,:,:,jp_sal)
131      ENDIF
132
133!!gm      IF( .NOT.ln_traqsr )   qsr(:,:) = 0.e0   ! no solar radiation penetration
134      IF( .NOT.ln_traqsr ) THEN     ! no solar radiation penetration
135         qns(:,:) = qns(:,:) + qsr(:,:)      ! total heat flux in qns
136         qsr(:,:) = 0.e0                     ! qsr set to zero
137      ENDIF
138
139      !----------------------------------------
140      !        EMP, EMPS and QNS effects
141      !----------------------------------------
142      !                                          Set before sbc tracer content fields
143      !                                          ************************************
144      IF( kt == nit000 ) THEN                      ! Set the forcing field at nit000 - 1
145         !                                         ! -----------------------------------
146         IF( ln_rstart .AND.    &                     ! Restart: read in restart file
147              & iom_varid( numror, 'sbc_hc_b', ldstop = .FALSE. ) > 0 ) THEN
148            IF(lwp) WRITE(numout,*) '          nit000-1 surface tracer content forcing fields red in the restart file'
149            zfact = 0.5e0
150            CALL iom_get( numror, jpdom_autoglo, 'sbc_hc_b', sbc_tsc_b(:,:,jp_tem) )   ! before heat content sbc trend
151            CALL iom_get( numror, jpdom_autoglo, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal) )   ! before salt content sbc trend
152         ELSE                                         ! No restart or restart not found: Euler forward time stepping
153            zfact = 1.e0
154            sbc_tsc_b(:,:,:) = 0.e0
155         ENDIF
156      ELSE                                         ! Swap of forcing fields
157         !                                         ! ----------------------
158         zfact = 0.5e0
159         sbc_tsc_b(:,:,:) = sbc_tsc(:,:,:)
160      ENDIF
161      !                                          Compute now sbc tracer content fields
162      !                                          *************************************
163
164                                                   ! Concentration dilution effect on (t,s) due to 
165                                                   ! evaporation, precipitation and qns, but not river runoff
166                                               
167      IF( lk_vvl ) THEN                            ! Variable Volume case  ==>> heat content of mass flux is in qns
168         DO jj = 1, jpj
169            DO ji = 1, jpi 
170               sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj)                              ! non solar heat flux
171               sbc_tsc(ji,jj,jp_sal) = r1_rau0     * sfx(ji,jj)                              ! salt flux due to freezing/melting
172            END DO
173         END DO
174      ELSE                                         ! Constant Volume case ==>> Concentration dilution effect
175         DO jj = 2, jpj
176            DO ji = fs_2, fs_jpim1   ! vector opt.
177               ! temperature : heat flux
178               sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj)                          &   ! non solar heat flux
179                  &                  + r1_rau0     * emp(ji,jj)  * tsn(ji,jj,1,jp_tem)       ! concent./dilut. effect
180               ! salinity    : salt flux + concent./dilut. effect (both in sfx)
181               sbc_tsc(ji,jj,jp_sal) = r1_rau0  * (  sfx(ji,jj)                          &   ! salt flux (freezing/melting)
182                  &                                + emp(ji,jj) * tsn(ji,jj,1,jp_sal) )      ! concent./dilut. effect
183            END DO
184         END DO
185         CALL iom_put( "emp_x_sst", emp (:,:) * tsn(:,:,1,jp_tem) )                          ! c/d term on sst
186         CALL iom_put( "emp_x_sss", emp (:,:) * tsn(:,:,1,jp_sal) )                          ! c/d term on sss
187      ENDIF
188      ! Concentration dilution effect on (t,s) due to evapouration, precipitation and qns, but not river runoff 
189      DO jn = 1, jpts
190         DO jj = 2, jpj
191            DO ji = fs_2, fs_jpim1   ! vector opt.
192               z1_e3t = zfact / fse3t(ji,jj,1)
193               tsa(ji,jj,1,jn) = tsa(ji,jj,1,jn) + ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) ) * z1_e3t
194            END DO
195         END DO
196      END DO
197      !                                          Write in the ocean restart file
198      !                                          *******************************
199      IF( lrst_oce ) THEN
200         IF(lwp) WRITE(numout,*)
201         IF(lwp) WRITE(numout,*) 'sbc : ocean surface tracer content forcing fields written in ocean restart file ',   &
202            &                    'at it= ', kt,' date= ', ndastp
203         IF(lwp) WRITE(numout,*) '~~~~'
204         CALL iom_rstput( kt, nitrst, numrow, 'sbc_hc_b', sbc_tsc(:,:,jp_tem) )
205         CALL iom_rstput( kt, nitrst, numrow, 'sbc_sc_b', sbc_tsc(:,:,jp_sal) )
206      ENDIF
207      !
208      !----------------------------------------
209      !        River Runoff effects
210      !----------------------------------------
211      !
212      zfact = 0.5e0
213
214      ! Effect on (t,s) due to river runoff (dilution effect automatically applied via vertical tracer advection)
215      IF( ln_rnf ) THEN 
216         DO jj = 2, jpj 
217            DO ji = fs_2, fs_jpim1
218               zdep = 1. / h_rnf(ji,jj)
219               zdep = zfact * zdep 
220               IF ( rnf(ji,jj) /= 0._wp ) THEN
221                  DO jk = 1, nk_rnf(ji,jj)
222                                        tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)   &
223                                          &               +  ( rnf_tsc_b(ji,jj,jp_tem) + rnf_tsc(ji,jj,jp_tem) ) * zdep
224                     IF( ln_rnf_sal )   tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal)   &
225                                          &               +  ( rnf_tsc_b(ji,jj,jp_sal) + rnf_tsc(ji,jj,jp_sal) ) * zdep 
226                  END DO
227               ENDIF
228            END DO 
229         END DO 
230      ENDIF 
231!!gm  It should be useless
232      CALL lbc_lnk( tsa(:,:,:,jp_tem), 'T', 1. )    ;    CALL lbc_lnk( tsa(:,:,:,jp_sal), 'T', 1. )
233
234      IF( l_trdtra )   THEN                      ! save the horizontal diffusive trends for further diagnostics
235         ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
236         ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
237         CALL trd_tra( kt, 'TRA', jp_tem, jptra_trd_nsr, ztrdt )
238         CALL trd_tra( kt, 'TRA', jp_sal, jptra_trd_nsr, ztrds )
239         CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds ) 
240      ENDIF
241      !
242      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' sbc  - Ta: ', mask1=tmask,   &
243         &                       tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
244      !
245      IF( nn_timing == 1 )  CALL timing_stop('tra_sbc')
246      !
247   END SUBROUTINE tra_sbc
248
249   !!======================================================================
250END MODULE trasbc
Note: See TracBrowser for help on using the repository browser.