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 branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/TRA – NEMO

source: branches/DEV_r2106_LOCEAN2010/NEMO/OPA_SRC/TRA/trasbc.F90 @ 2224

Last change on this file since 2224 was 2148, checked in by cetlod, 14 years ago

merge LOCEAN 2010 developments branches
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 15.6 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
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC   tra_sbc    ! routine called by step.F90
34
35   !! * Substitutions
36#  include "domzgr_substitute.h90"
37#  include "vectopt_loop_substitute.h90"
38   !!----------------------------------------------------------------------
39   !! NEMO/OPA 3.3 , LOCEAN-IPSL (2010)
40   !! $Id$
41   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
42   !!----------------------------------------------------------------------
43
44CONTAINS
45
46   SUBROUTINE tra_sbc ( kt )
47      !!----------------------------------------------------------------------
48      !!                  ***  ROUTINE tra_sbc  ***
49      !!                   
50      !! ** Purpose :   Compute the tracer surface boundary condition trend of
51      !!      (flux through the interface, concentration/dilution effect)
52      !!      and add it to the general trend of tracer equations.
53      !!
54      !! ** Method :
55      !!      Following Roullet and Madec (2000), the air-sea flux can be divided
56      !!      into three effects: (1) Fext, external forcing;
57      !!      (2) Fwi, concentration/dilution effect due to water exchanged
58      !!         at the surface by evaporation, precipitations and runoff (E-P-R);
59      !!      (3) Fwe, tracer carried with the water that is exchanged.
60      !!
61      !!      Fext, flux through the air-sea interface for temperature and salt:
62      !!            - temperature : heat flux q (w/m2). If penetrative solar
63      !!         radiation q is only the non solar part of the heat flux, the
64      !!         solar part is added in traqsr.F routine.
65      !!            ta = ta + q /(rau0 rcp e3t)  for k=1
66      !!            - salinity    : no salt flux
67      !!
68      !!      The formulation for Fwb and Fwi vary according to the free
69      !!      surface formulation (linear or variable volume).
70      !!      * Linear free surface
71      !!            The surface freshwater flux modifies the ocean volume
72      !!         and thus the concentration of a tracer and the temperature.
73      !!         First order of the effect of surface freshwater exchange
74      !!         for salinity, it can be neglected on temperature (especially
75      !!         as the temperature of precipitations and runoffs is usually
76      !!         unknown).
77      !!            - temperature : we assume that the temperature of both
78      !!         precipitations and runoffs is equal to the SST, thus there
79      !!         is no additional flux since in this case, the concentration
80      !!         dilution effect is balanced by the net heat flux associated
81      !!         to the freshwater exchange (Fwe+Fwi=0):
82      !!            (Tp P - Te E) + SST (P-E) = 0 when Tp=Te=SST
83      !!            - salinity    : evaporation, precipitation and runoff
84      !!         water has a zero salinity (Fwe=0), thus only Fwi remains:
85      !!            sa = sa + emp * sn / e3t   for k=1
86      !!         where emp, the surface freshwater budget (evaporation minus
87      !!         precipitation minus runoff) given in kg/m2/s is divided
88      !!         by 1035 kg/m3 (density of ocena water) to obtain m/s.   
89      !!         Note: even though Fwe does not appear explicitly for
90      !!         temperature in this routine, the heat carried by the water
91      !!         exchanged through the surface is part of the total heat flux
92      !!         forcing and must be taken into account in the global heat
93      !!         balance).
94      !!      * nonlinear free surface (variable volume, lk_vvl)
95      !!         contrary to the linear free surface case, Fwi is properly
96      !!         taken into account by using the true layer thicknesses to       
97      !!         calculate tracer content and advection. There is no need to
98      !!         deal with it in this routine.
99      !!           - temperature: Fwe=SST (P-E+R) is added to Fext.
100      !!           - salinity:  Fwe = 0, there is no surface flux of salt.
101      !!
102      !! ** Action  : - Update the 1st level of (ta,sa) with the trend associated
103      !!                with the tracer surface boundary condition
104      !!              - save the trend it in ttrd ('key_trdtra')
105      !!----------------------------------------------------------------------
106      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
107      !!
108      INTEGER  ::   ji, jj, jk           ! dummy loop indices 
109      REAL(wp) ::   zta, zsa, zrnf       ! local scalars, adjustment to temperature and salinity 
110      REAL(wp) ::   zsrau, zse3t, zdep   ! local scalars, 1/density, 1/height of box, 1/height of effected water column 
111      REAL(wp) ::   zdheat, zdsalt       ! total change of temperature and salinity 
112      REAL(wp) ::   zfact, z1_e3t        !
113      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  ztrdt, ztrds
114      !!----------------------------------------------------------------------
115
116      IF( kt == nit000 ) THEN
117         IF(lwp) WRITE(numout,*)
118         IF(lwp) WRITE(numout,*) 'tra_sbc : TRAcer Surface Boundary Condition'
119         IF(lwp) WRITE(numout,*) '~~~~~~~ '
120      ENDIF
121
122      zsrau = 1. / rau0             ! initialization
123
124      IF( l_trdtra )   THEN                    !* Save ta and sa trends
125         ALLOCATE( ztrdt(jpi,jpj,jpk) )   ;    ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
126         ALLOCATE( ztrds(jpi,jpj,jpk) )   ;    ztrds(:,:,:) = tsa(:,:,:,jp_sal)
127      ENDIF
128
129!!gm      IF( .NOT.ln_traqsr )   qsr(:,:) = 0.e0   ! no solar radiation penetration
130      IF( .NOT.ln_traqsr ) THEN     ! no solar radiation penetration
131         qns(:,:) = qns(:,:) + qsr(:,:)      ! total heat flux in qns
132         qsr(:,:) = 0.e0                     ! qsr set to zero
133      ENDIF
134
135      !                                          Set before sbc tracer content fields
136      !                                          ************************************
137      IF( kt == nit000 ) THEN                      ! Set the forcing field at nit000 - 1
138         !                                         ! -----------------------------------
139         IF( ln_rstart .AND.    &                     ! Restart: read in restart file
140              & iom_varid( numror, 'sbc_hc_b', ldstop = .FALSE. ) > 0 ) THEN
141            IF(lwp) WRITE(numout,*) '          nit000-1 surface tracer content forcing fields red in the restart file'
142            zfact = 0.5e0
143            CALL iom_get( numror, jpdom_autoglo, 'sbc_hc_b', sbc_hc_b )   ! before heat content sbc trend
144            CALL iom_get( numror, jpdom_autoglo, 'sbc_sc_b', sbc_sc_b )   ! before salt content sbc trend
145         ELSE                                         ! No restart or restart not found: Euler forward time stepping
146            zfact = 1.e0
147            sbc_hc_b(:,:) = 0.e0
148            sbc_sc_b(:,:) = 0.e0
149         ENDIF
150      ELSE                                         ! Swap of forcing fields
151         !                                         ! ----------------------
152         zfact = 0.5e0
153         sbc_hc_b(:,:) = sbc_hc_n(:,:)
154         sbc_sc_b(:,:) = sbc_sc_n(:,:)
155      ENDIF
156      !                                          Compute now sbc tracer content fields
157      !                                          *************************************
158
159                                                   ! Concentration dilution effect on (t,s) due to 
160                                                   ! evaporation, precipitation and qns, but not river runoff
161                                               
162      IF( lk_vvl ) THEN                            ! Variable Volume case
163         DO jj = 2, jpj
164            DO ji = fs_2, fs_jpim1   ! vector opt.
165               ! temperature : heat flux + cooling/heating effet of EMP flux
166               sbc_hc_n(ji,jj) = ro0cpr * qns(ji,jj) - zsrau * emp(ji,jj) * tsn(ji,jj,1,jp_tem)
167               ! concent./dilut. effect due to sea-ice melt/formation and (possibly) SSS restoration
168               sbc_sc_n(ji,jj) = ( emps(ji,jj) - emp(ji,jj) ) * zsrau * tsn(ji,jj,1,jp_sal)
169            END DO
170         END DO
171      ELSE                                         ! Constant Volume case
172         DO jj = 2, jpj
173            DO ji = fs_2, fs_jpim1   ! vector opt.
174               ! temperature : heat flux
175               sbc_hc_n(ji,jj) = ro0cpr * qns(ji,jj)
176               ! salinity    : salt flux + concent./dilut. effect (both in emps)
177               sbc_sc_n(ji,jj) = zsrau * emps(ji,jj) * tsn(ji,jj,1,jp_sal)
178            END DO
179         END DO
180      ENDIF
181                                                   ! Concentration dilution effect on (t,s) due to 
182                                                   ! river runoff without T, S and depth attributes
183      IF( ln_rnf ) THEN 
184         !
185         IF( lk_vvl ) THEN                            ! Variable Volume case
186            !                                         ! cooling/heating effect of runoff & No salinity concent./dilut. effect
187            DO jj = 2, jpj
188               DO ji = fs_2, fs_jpim1   ! vector opt.
189                  sbc_hc_n(ji,jj) = sbc_hc_n(ji,jj) + zsrau * rnf(ji,jj) * tsn(ji,jj,1,jp_tem)
190               END DO
191            END DO
192         ELSE                                         ! Constant Volume case
193            !                                         ! concent./dilut. effect only
194            DO jj = 2, jpj 
195               DO ji = fs_2, fs_jpim1   ! vector opt.
196                  sbc_sc_n(ji,jj) = sbc_sc_n(ji,jj) - zsrau * rnf(ji,jj) * tsn(ji,jj,1,jp_sal)
197               END DO
198            END DO
199         ENDIF
200         ! 
201      ENDIF
202      !                                          Add to the general trend
203      !                                          ************************
204      DO jj = 2, jpj
205         DO ji = fs_2, fs_jpim1   ! vector opt.
206            z1_e3t = zfact / fse3t(ji,jj,1)
207            tsa(ji,jj,1,jp_tem) = tsa(ji,jj,1,jp_tem) + ( sbc_hc_b(ji,jj) + sbc_hc_n(ji,jj) ) * z1_e3t
208            tsa(ji,jj,1,jp_sal) = tsa(ji,jj,1,jp_sal) + ( sbc_sc_b(ji,jj) + sbc_sc_n(ji,jj) ) * z1_e3t
209         END DO
210      END DO
211      !                                          Write in the ocean restart file
212      !                                          *******************************
213      IF( lrst_oce ) THEN
214         IF(lwp) WRITE(numout,*)
215         IF(lwp) WRITE(numout,*) 'sbc : ocean surface tracer content forcing fields written in ocean restart file ',   &
216            &                    'at it= ', kt,' date= ', ndastp
217         IF(lwp) WRITE(numout,*) '~~~~'
218         CALL iom_rstput( kt, nitrst, numrow, 'sbc_hc_b', sbc_hc_n )
219         CALL iom_rstput( kt, nitrst, numrow, 'sbc_sc_b', sbc_sc_n )
220      ENDIF
221
222      IF( ln_rnf .AND. ln_rnf_att ) THEN        ! Concentration / dilution effect on (t,s) due to river runoff 
223        !
224        DO jj = 1, jpj 
225           DO ji = 1, jpi 
226              zdep  = 1. / rnf_dep(ji,jj) 
227              zse3t = 1. / fse3t(ji,jj,1) 
228              rnf_dep(ji,jj) = 0.e0
229              DO jk = 1, rnf_mod_dep(ji,jj)                          ! recalculates rnf_dep to be the depth 
230                rnf_dep(ji,jj) = rnf_dep(ji,jj) + fse3t(ji,jj,jk)    ! in metres to the bottom of the relevant grid box 
231              END DO
232              IF( rnf_tmp(ji,jj) == -999 )   rnf_tmp(ji,jj) = tsn(ji,jj,1,jp_tem)   ! if not specified set runoff temp to be sst 
233 
234              IF( rnf(ji,jj) > 0.e0 ) THEN 
235                !
236                zrnf = rnf(ji,jj) * zsrau * zdep
237                IF( lk_vvl ) THEN 
238                  ! indirect flux, concentration or dilution effect : force a dilution effect in all levels
239                  zdheat = 0.e0 
240                  zdsalt = 0.e0 
241                  DO jk = 1, rnf_mod_dep(ji,jj) 
242                    zta = -tsn(ji,jj,jk,jp_tem) * zrnf 
243                    zsa = -tsn(ji,jj,jk,jp_sal) * zrnf 
244                    tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + zta        ! add the trend to the general tracer trend
245                    tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) + zsa
246                    zdheat = zdheat + zta * fse3t(ji,jj,jk) 
247                    zdsalt = zdsalt + zsa * fse3t(ji,jj,jk) 
248                  END DO 
249                  ! negate this total change in heat and salt content from top level    !!gm I don't understand this
250                  tsa(ji,jj,1,jp_tem) = tsa(ji,jj,1,jp_tem) - zdheat * zse3t
251                  tsa(ji,jj,1,jp_sal) = tsa(ji,jj,1,jp_sal) - zdsalt * zse3t
252   
253                  DO jk = 1, rnf_mod_dep(ji,jj) 
254                    tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + rnf_tmp(ji,jj) * zrnf
255                    tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) + rnf_sal(ji,jj) * zrnf
256                  END DO 
257                ELSE 
258                  DO jk = 1, rnf_mod_dep(ji,jj) 
259                    tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) + ( rnf_tmp(ji,jj) - tsn(ji,jj,jk,jp_tem) ) * zrnf
260                    tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) + ( rnf_sal(ji,jj) - tsn(ji,jj,jk,jp_sal) ) * zrnf
261                  END DO 
262                ENDIF 
263 
264              ELSEIF( rnf(ji,jj) < 0.e0 ) THEN   ! for use in baltic when flow is out of domain, want no change in temp and sal 
265 
266                IF( lk_vvl ) THEN 
267                  ! calculate automatic adjustment to sal and temp due to dilution/concentraion effect   
268                  zrnf = rnf(ji,jj) * zsrau * zse3t
269                  tsa(ji,jj,1,jp_tem) = tsa(ji,jj,1,jp_tem) + tsn(ji,jj,1,jp_tem) * zrnf 
270                  tsa(ji,jj,1,jp_sal) = tsa(ji,jj,1,jp_sal) + tsn(ji,jj,1,jp_sal) * zrnf
271                ENDIF 
272 
273              ENDIF 
274 
275           END DO 
276        END DO 
277        !
278      ENDIF 
279
280      IF( l_trdtra )   THEN                      ! save the horizontal diffusive trends for further diagnostics
281         ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
282         ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
283         CALL trd_tra( kt, 'TRA', jp_tem, jptra_trd_nsr, ztrdt )
284         CALL trd_tra( kt, 'TRA', jp_sal, jptra_trd_nsr, ztrds )
285         DEALLOCATE( ztrdt )      ;     DEALLOCATE( ztrds )
286      ENDIF
287      !
288      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' sbc  - Ta: ', mask1=tmask,   &
289         &                       tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
290      !
291   END SUBROUTINE tra_sbc
292
293   !!======================================================================
294END MODULE trasbc
Note: See TracBrowser for help on using the repository browser.