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/UKMO/ROMS_WAD_7832/NEMOGCM/NEMO/OPA_SRC/TRA – NEMO

source: branches/UKMO/ROMS_WAD_7832/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90 @ 8865

Last change on this file since 8865 was 8865, checked in by deazer, 6 years ago

Moving Changes from CS15mini config into NEMO main src
Updating TEST configs to run wit this version of the code
all sette tests pass at this revision other than AGRIF
Includes updates to dynnxt and tranxt required for 3D rives in WAD case to be conservative.

Next commit will update naming conventions and tidy the code.
  • Property svn:keywords set to Id
File size: 12.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   !!            3.6  !  2014-11  (P. Mathiot) isf melting forcing
12   !!----------------------------------------------------------------------
13
14   !!----------------------------------------------------------------------
15   !!   tra_sbc       : update the tracer trend at ocean surface
16   !!----------------------------------------------------------------------
17   USE oce            ! ocean dynamics and active tracers
18   USE sbc_oce        ! surface boundary condition: ocean
19   USE dom_oce        ! ocean space domain variables
20   USE phycst         ! physical constant
21   USE eosbn2         ! Equation Of State
22   USE sbcmod         ! ln_rnf 
23   USE sbcrnf         ! River runoff 
24   USE sbcisf         ! Ice shelf   
25   USE iscplini       ! Ice sheet coupling
26   USE traqsr         ! solar radiation penetration
27   USE trd_oce        ! trends: ocean variables
28   USE trdtra         ! trends manager: tracers
29   !
30   USE in_out_manager ! I/O manager
31   USE prtctl         ! Print control
32   USE iom            ! xIOS server
33   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
34   USE wrk_nemo       ! Memory Allocation
35   USE timing         ! Timing
36   USE wet_dry
37
38   IMPLICIT NONE
39   PRIVATE
40
41   PUBLIC   tra_sbc   ! routine called by step.F90
42
43   !! * Substitutions
44#  include "vectopt_loop_substitute.h90"
45   !!----------------------------------------------------------------------
46   !! NEMO/OPA 3.7 , NEMO Consortium (2014)
47   !! $Id$
48   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
49   !!----------------------------------------------------------------------
50CONTAINS
51
52   SUBROUTINE tra_sbc ( kt )
53      !!----------------------------------------------------------------------
54      !!                  ***  ROUTINE tra_sbc  ***
55      !!                   
56      !! ** Purpose :   Compute the tracer surface boundary condition trend of
57      !!      (flux through the interface, concentration/dilution effect)
58      !!      and add it to the general trend of tracer equations.
59      !!
60      !! ** Method :   The (air+ice)-sea flux has two components:
61      !!      (1) Fext, external forcing (i.e. flux through the (air+ice)-sea interface);
62      !!      (2) Fwe , tracer carried with the water that is exchanged with air+ice.
63      !!               The input forcing fields (emp, rnf, sfx, isf) contain Fext+Fwe,
64      !!             they are simply added to the tracer trend (tsa).
65      !!               In linear free surface case (ln_linssh=T), the volume of the
66      !!             ocean does not change with the water exchanges at the (air+ice)-sea
67      !!             interface. Therefore another term has to be added, to mimic the
68      !!             concentration/dilution effect associated with water exchanges.
69      !!
70      !! ** Action  : - Update tsa with the surface boundary condition trend
71      !!              - send trends to trdtra module for further diagnostics(l_trdtra=T)
72      !!----------------------------------------------------------------------
73      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
74      !
75      INTEGER  ::   ji, jj, jk, jn        ! dummy loop indices 
76      INTEGER  ::   ikt, ikb              ! local integers
77      REAL(wp) ::   zfact, z1_e3t, zdep   ! local scalar
78      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrdt, ztrds
79      !!----------------------------------------------------------------------
80      !
81      IF( nn_timing == 1 )  CALL timing_start('tra_sbc')
82      !
83      IF( kt == nit000 ) THEN
84         IF(lwp) WRITE(numout,*)
85         IF(lwp) WRITE(numout,*) 'tra_sbc : TRAcer Surface Boundary Condition'
86         IF(lwp) WRITE(numout,*) '~~~~~~~ '
87      ENDIF
88      !
89      IF( l_trdtra ) THEN                    !* Save ta and sa trends
90         CALL wrk_alloc( jpi, jpj, jpk, ztrdt, ztrds ) 
91         ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
92         ztrds(:,:,:) = tsa(:,:,:,jp_sal)
93      ENDIF
94      !
95!!gm  This should be moved into sbcmod.F90 module ? (especially now that ln_traqsr is read in namsbc namelist)
96      IF( .NOT.ln_traqsr ) THEN     ! no solar radiation penetration
97         qns(:,:) = qns(:,:) + qsr(:,:)      ! total heat flux in qns
98         qsr(:,:) = 0._wp                     ! qsr set to zero
99      ENDIF
100
101      !----------------------------------------
102      !        EMP, SFX and QNS effects
103      !----------------------------------------
104      !                             !==  Set before sbc tracer content fields  ==!
105      IF( kt == nit000 ) THEN             !* 1st time-step
106         IF( ln_rstart .AND.    &               ! Restart: read in restart file
107              & iom_varid( numror, 'sbc_hc_b', ldstop = .FALSE. ) > 0 ) THEN
108            IF(lwp) WRITE(numout,*) '          nit000-1 sbc tracer content field read in the restart file'
109            zfact = 0.5_wp
110            sbc_tsc(:,:,:) = 0._wp
111            CALL iom_get( numror, jpdom_autoglo, 'sbc_hc_b', sbc_tsc_b(:,:,jp_tem) )   ! before heat content sbc trend
112            CALL iom_get( numror, jpdom_autoglo, 'sbc_sc_b', sbc_tsc_b(:,:,jp_sal) )   ! before salt content sbc trend
113         ELSE                                   ! No restart or restart not found: Euler forward time stepping
114            zfact = 1._wp
115            sbc_tsc(:,:,:) = 0._wp
116            sbc_tsc_b(:,:,:) = 0._wp
117         ENDIF
118      ELSE                                !* other time-steps: swap of forcing fields
119         zfact = 0.5_wp
120         sbc_tsc_b(:,:,:) = sbc_tsc(:,:,:)
121      ENDIF
122      !                             !==  Now sbc tracer content fields  ==!
123      DO jj = 2, jpj
124         DO ji = fs_2, fs_jpim1   ! vector opt.
125            IF ( ln_rwd_rmp ) THEN    ! If near WAD point limite the flux for now
126               IF ( sshn(ji,jj) + ht_0(ji,jj) >  2._wp * rn_wdmin1 ) THEN
127                  sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj)   ! non solar heat flux
128               ELSE IF ( sshn(ji,jj) + ht_0(ji,jj) >  rn_wdmin1 ) THEN
129                  sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj) * (tanh(5._wp*( ( sshn(ji,jj) + ht_0(ji,jj) -  rn_wdmin1 )/rn_wdmin1)) )
130               ELSE
131                  sbc_tsc(ji,jj,jp_tem) = 0._wp
132               ENDIF
133            ELSE
134               sbc_tsc(ji,jj,jp_tem) = r1_rau0_rcp * qns(ji,jj)   ! non solar heat flux
135            ENDIF
136
137            sbc_tsc(ji,jj,jp_sal) = r1_rau0     * sfx(ji,jj)   ! salt flux due to freezing/melting
138         END DO
139      END DO
140      IF( ln_linssh ) THEN                !* linear free surface 
141         DO jj = 2, jpj                         !==>> add concentration/dilution effect due to constant volume cell
142            DO ji = fs_2, fs_jpim1   ! vector opt.
143               sbc_tsc(ji,jj,jp_tem) = sbc_tsc(ji,jj,jp_tem) + r1_rau0 * emp(ji,jj) * tsn(ji,jj,1,jp_tem)
144               sbc_tsc(ji,jj,jp_sal) = sbc_tsc(ji,jj,jp_sal) + r1_rau0 * emp(ji,jj) * tsn(ji,jj,1,jp_sal)
145            END DO
146         END DO                                 !==>> output c./d. term
147         IF( iom_use('emp_x_sst') )   CALL iom_put( "emp_x_sst", emp (:,:) * tsn(:,:,1,jp_tem) )
148         IF( iom_use('emp_x_sss') )   CALL iom_put( "emp_x_sss", emp (:,:) * tsn(:,:,1,jp_sal) )
149      ENDIF
150      !
151      DO jn = 1, jpts               !==  update tracer trend  ==!
152         DO jj = 2, jpj
153            DO ji = fs_2, fs_jpim1   ! vector opt. 
154               tsa(ji,jj,1,jn) = tsa(ji,jj,1,jn) + zfact * ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) ) / e3t_n(ji,jj,1)
155            END DO
156         END DO
157      END DO
158      !                 
159      IF( lrst_oce ) THEN           !==  write sbc_tsc in the ocean restart file  ==!
160         CALL iom_rstput( kt, nitrst, numrow, 'sbc_hc_b', sbc_tsc(:,:,jp_tem) )
161         CALL iom_rstput( kt, nitrst, numrow, 'sbc_sc_b', sbc_tsc(:,:,jp_sal) )
162      ENDIF
163      !
164      !----------------------------------------
165      !       Ice Shelf effects (ISF)
166      !     tbl treated as in Losh (2008) JGR
167      !----------------------------------------
168      !
169!!gm BUG ?   Why no differences between non-linear and linear free surface ?
170!!gm         probably taken into account in r1_hisf_tbl : to be verified
171      IF( ln_isf ) THEN
172         zfact = 0.5_wp
173         DO jj = 2, jpj
174            DO ji = fs_2, fs_jpim1
175               !
176               ikt = misfkt(ji,jj)
177               ikb = misfkb(ji,jj)
178               !
179               ! level fully include in the ice shelf boundary layer
180               ! sign - because fwf sign of evapo (rnf sign of precip)
181               DO jk = ikt, ikb - 1
182               ! compute trend
183                  tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)                                                &
184                     &           + zfact * ( risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem) )             &
185                     &           * r1_hisf_tbl(ji,jj)
186               END DO
187   
188               ! level partially include in ice shelf boundary layer
189               ! compute trend
190               tsa(ji,jj,ikb,jp_tem) = tsa(ji,jj,ikb,jp_tem)                                                 &
191                  &              + zfact * ( risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem) )             &
192                  &              * r1_hisf_tbl(ji,jj) * ralpha(ji,jj)
193
194            END DO
195         END DO
196      END IF
197      !
198      !----------------------------------------
199      !        River Runoff effects
200      !----------------------------------------
201      !
202      IF( ln_rnf ) THEN         ! input of heat and salt due to river runoff
203         zfact = 0.5_wp
204         DO jj = 2, jpj 
205            DO ji = fs_2, fs_jpim1
206               IF( rnf(ji,jj) /= 0._wp ) THEN
207                  zdep = zfact / h_rnf(ji,jj)
208                  DO jk = 1, nk_rnf(ji,jj)
209                                        tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)                                 &
210                                           &                 +  ( rnf_tsc_b(ji,jj,jp_tem) + rnf_tsc(ji,jj,jp_tem) ) * zdep
211                     IF( ln_rnf_sal )   tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal)                                 &
212                                           &                 +  ( rnf_tsc_b(ji,jj,jp_sal) + rnf_tsc(ji,jj,jp_sal) ) * zdep 
213                  END DO
214               ENDIF
215            END DO 
216         END DO 
217      ENDIF
218
219      IF( iom_use('rnf_x_sst') )   CALL iom_put( "rnf_x_sst", rnf*tsn(:,:,1,jp_tem) )   ! runoff term on sst
220      IF( iom_use('rnf_x_sss') )   CALL iom_put( "rnf_x_sss", rnf*tsn(:,:,1,jp_sal) )   ! runoff term on sss
221
222      !
223      !----------------------------------------
224      !        Ice Sheet coupling imbalance correction to have conservation
225      !----------------------------------------
226      !
227      IF( ln_iscpl .AND. ln_hsb) THEN         ! input of heat and salt due to river runoff
228         DO jk = 1,jpk
229            DO jj = 2, jpj 
230               DO ji = fs_2, fs_jpim1
231                  zdep = 1._wp / e3t_n(ji,jj,jk) 
232                  tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) - htsc_iscpl(ji,jj,jk,jp_tem)                       &
233                      &                                         * zdep
234                  tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) - htsc_iscpl(ji,jj,jk,jp_sal)                       &
235                      &                                         * zdep 
236               END DO 
237            END DO 
238         END DO
239      ENDIF
240
241      IF( l_trdtra )   THEN                      ! save the horizontal diffusive trends for further diagnostics
242         ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
243         ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
244         CALL trd_tra( kt, 'TRA', jp_tem, jptra_nsr, ztrdt )
245         CALL trd_tra( kt, 'TRA', jp_sal, jptra_nsr, ztrds )
246         CALL wrk_dealloc( jpi, jpj, jpk, ztrdt, ztrds ) 
247      ENDIF
248      !
249      IF(ln_ctl)   CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' sbc  - Ta: ', mask1=tmask,   &
250         &                       tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
251      !
252      IF( nn_timing == 1 )  CALL timing_stop('tra_sbc')
253      !
254   END SUBROUTINE tra_sbc
255
256   !!======================================================================
257END MODULE trasbc
Note: See TracBrowser for help on using the repository browser.