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.
trasbc.F90 in trunk/NEMOGCM/NEMO/OPA_SRC/TRA – NEMO

source: trunk/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90 @ 4832

Last change on this file since 4832 was 3764, checked in by smasson, 11 years ago

dev_MERGE_2012: report bugfixes done in the trunk from r3555 to r3763 into dev_MERGE_2012

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