source: branches/2013/dev_LOCEAN_2013/NEMOGCM/NEMO/TOP_SRC/trcsub.F90 @ 4148

Last change on this file since 4148 was 4148, checked in by cetlod, 8 years ago

merge in trunk changes between r3853 and r3940 and commit the changes, see ticket #1169

File size: 51.2 KB
Line 
1MODULE trcsub
2   !!======================================================================
3   !!                       ***  MODULE trcsubstp  ***
4   !!TOP :   Averages physics variables for TOP substepping.
5   !!======================================================================
6   !! History :  1.0  !  2011-10  (K. Edwards)  Original
7   !!----------------------------------------------------------------------
8#if defined key_top
9   !!----------------------------------------------------------------------
10   !!   trc_sub    : passive tracer system sub-stepping
11   !!----------------------------------------------------------------------
12   USE oce_trc          ! ocean dynamics and active tracers variables
13   USE trc
14   USE prtctl_trc       ! Print control for debbuging
15   USE iom
16   USE in_out_manager
17   USE lbclnk
18#if defined key_zdftke
19   USE zdftke          ! twice TKE (en)
20#endif
21#if defined key_zdfgls
22   USE zdfgls, ONLY: en
23#endif
24   USE trabbl
25   USE zdf_oce
26   USE domvvl
27   USE divcur          ! hor. divergence and curl      (div & cur routines)
28   USE sbcrnf, ONLY: h_rnf, nk_rnf   ! River runoff
29   USE sbc_oce         ! surface boundary condition: ocean
30   USE bdy_oce
31#if defined key_obc
32   USE obc_oce, ONLY: obctmsk
33#endif
34#if defined key_agrif
35   USE agrif_opa_update
36   USE agrif_opa_interp
37#endif
38
39   IMPLICIT NONE
40
41   PUBLIC   trc_sub_stp      ! called by trc_stp
42   PUBLIC   trc_sub_ini      ! called by trc_ini to initialize substepping arrays.
43   PUBLIC   trc_sub_reset    ! called by trc_stp to reset physics variables
44   PUBLIC   trc_sub_ssh      ! called by trc_stp to reset physics variables
45
46   !!* Module variables
47   REAL(wp)  :: r1_ndttrc     !    1 /  nn_dttrc
48   REAL(wp)  :: r1_ndttrcp1   !    1 / (nn_dttrc+1)
49
50   !!* Substitution
51#  include "top_substitute.h90"
52   !!----------------------------------------------------------------------
53   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
54   !! $Id: trcstp.F90 2528 2010-12-27 17:33:53Z rblod $
55   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
56   !!----------------------------------------------------------------------
57CONTAINS
58
59   SUBROUTINE trc_sub_stp( kt )
60      !!-------------------------------------------------------------------
61      !!                     ***  ROUTINE trc_stp  ***
62      !!                     
63      !! ** Purpose : Average variables needed for sub-stepping passive tracers
64      !!
65      !! ** Method  : Called every timestep to increment _tm (time mean) variables
66      !!              on TOP steps, calculate averages.
67      !!-------------------------------------------------------------------
68      INTEGER, INTENT( in ) ::  kt        ! ocean time-step index
69      INTEGER               ::  ji,jj,jk  ! dummy loop indices
70      REAL(wp)              ::  z1_ne3t, z1_ne3u, z1_ne3v, z1_ne3w
71      !!-------------------------------------------------------------------
72      !
73      IF( nn_timing == 1 )  CALL timing_start('trc_sub_stp')
74      !
75      IF( kt == nit000 ) THEN
76           IF(lwp) WRITE(numout,*)
77           IF(lwp) WRITE(numout,*) 'trc_sub_stp : substepping of the passive tracers'
78           IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
79           !
80           sshb_hold  (:,:) = sshn  (:,:)
81           sshu_b_hold(:,:) = sshu_n(:,:)
82           sshv_b_hold(:,:) = sshv_n(:,:)
83           emp_b_hold (:,:) = emp_b (:,:)
84           !
85           r1_ndttrc        = 1._wp / REAL( nn_dttrc    , wp ) 
86           r1_ndttrcp1      = 1._wp / REAL( nn_dttrc + 1, wp )
87           !
88      ENDIF 
89
90       IF( MOD( kt , nn_dttrc ) /= 0 ) THEN
91          !
92          un_tm   (:,:,:)        = un_tm   (:,:,:)        + un   (:,:,:)        * fse3u(:,:,:) 
93          vn_tm   (:,:,:)        = vn_tm   (:,:,:)        + vn   (:,:,:)        * fse3v(:,:,:) 
94          tsn_tm  (:,:,:,jp_tem) = tsn_tm  (:,:,:,jp_tem) + tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
95          tsn_tm  (:,:,:,jp_sal) = tsn_tm  (:,:,:,jp_sal) + tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
96          rhop_tm (:,:,:)        = rhop_tm (:,:,:)        + rhop (:,:,:)        * fse3t(:,:,:) 
97          avt_tm  (:,:,:)        = avt_tm  (:,:,:)        + avt  (:,:,:)        * fse3w(:,:,:) 
98# if defined key_zdfddm
99          avs_tm  (:,:,:)        = avs_tm  (:,:,:)        + avs  (:,:,:)        * fse3w(:,:,:) 
100# endif
101#if defined key_ldfslp
102          wslpi_tm(:,:,:)        = wslpi_tm(:,:,:)        + wslpi(:,:,:)        * fse3w(:,:,:) 
103          wslpj_tm(:,:,:)        = wslpj_tm(:,:,:)        + wslpj(:,:,:)        * fse3w(:,:,:) 
104          uslp_tm (:,:,:)        = uslp_tm (:,:,:)        + uslp (:,:,:)        * fse3u(:,:,:) 
105          vslp_tm (:,:,:)        = vslp_tm (:,:,:)        + vslp (:,:,:)        * fse3v(:,:,:) 
106#endif
107# if defined key_trabbl
108          IF( nn_bbl_ldf == 1 ) THEN
109             ahu_bbl_tm(:,:)     = ahu_bbl_tm(:,:)        + ahu_bbl(:,:) 
110             ahv_bbl_tm(:,:)     = ahv_bbl_tm(:,:)        + ahv_bbl(:,:) 
111          ENDIF
112          IF( nn_bbl_adv == 1 ) THEN
113             utr_bbl_tm(:,:)     = utr_bbl_tm(:,:)        + utr_bbl(:,:) 
114             vtr_bbl_tm(:,:)     = vtr_bbl_tm(:,:)        + vtr_bbl(:,:) 
115          ENDIF
116# endif
117          !
118          sshn_tm  (:,:)         = sshn_tm  (:,:)         + sshn  (:,:) 
119          sshu_n_tm(:,:)         = sshu_n_tm(:,:)         + sshu_n(:,:) 
120          sshv_n_tm(:,:)         = sshv_n_tm(:,:)         + sshv_n(:,:) 
121          rnf_tm   (:,:)         = rnf_tm   (:,:)         + rnf   (:,:) 
122          h_rnf_tm (:,:)         = h_rnf_tm (:,:)         + h_rnf (:,:) 
123          hmld_tm  (:,:)         = hmld_tm  (:,:)         + hmld  (:,:)
124          fr_i_tm  (:,:)         = fr_i_tm  (:,:)         + fr_i  (:,:)
125          emp_tm   (:,:)         = emp_tm   (:,:)         + emp   (:,:) 
126          fmmflx_tm(:,:)         = fmmflx_tm(:,:)         + fmmflx(:,:)
127          qsr_tm   (:,:)         = qsr_tm   (:,:)         + qsr   (:,:)
128          wndm_tm  (:,:)         = wndm_tm  (:,:)         + wndm  (:,:)
129          !
130#if defined key_traldf_c3d
131          ahtt_tm  (:,:,:)       = ahtt_tm  (:,:,:)       + ahtt(:,:,:)         * fse3t(:,:,:)
132          ahtu_tm  (:,:,:)       = ahtu_tm  (:,:,:)       + ahtu(:,:,:)         * fse3u(:,:,:)
133          ahtv_tm  (:,:,:)       = ahtv_tm  (:,:,:)       + ahtv(:,:,:)         * fse3v(:,:,:)
134          ahtw_tm  (:,:,:)       = ahtw_tm  (:,:,:)       + ahtw(:,:,:)         * fse3w(:,:,:)
135#elif defined key_traldf_c2d
136          ahtt_tm  (:,:)         = ahtt_tm  (:,:)         + ahtt(:,:)
137          ahtu_tm  (:,:)         = ahtu_tm  (:,:)         + ahtu(:,:)
138          ahtv_tm  (:,:)         = ahtv_tm  (:,:)         + ahtv(:,:)
139          ahtw_tm  (:,:)         = ahtw_tm  (:,:)         + ahtw(:,:)
140#elif defined key_traldf_c1d
141          ahtt_tm  (:)           = ahtt_tm  (:,:)         + ahtt(:)
142          ahtu_tm  (:)           = ahtu_tm  (:,:)         + ahtu(:)
143          ahtv_tm  (:)           = ahtv_tm  (:,:)         + ahtv(:)
144          ahtw_tm  (:)           = ahtw_tm  (:,:)         + ahtw(:)
145#else
146          ahtt_tm                = ahtt_tm                + ahtt
147          ahtu_tm                = ahtu_tm                + ahtu
148          ahtv_tm                = ahtv_tm                + ahtv
149          ahtw_tm                = ahtw_tm                + ahtw
150#endif
151#if defined key_traldf_eiv
152#  if defined key_traldf_c3d
153          aeiu_tm  (:,:,:)       = aeiu_tm  (:,:,:)       + aeiu(:,:,:)         * fse3u(:,:,:)
154          aeiv_tm  (:,:,:)       = aeiv_tm  (:,:,:)       + aeiv(:,:,:)         * fse3v(:,:,:)
155          aeiw_tm  (:,:,:)       = aeiw_tm  (:,:,:)       + aeiw(:,:,:)         * fse3w(:,:,:)
156#  elif defined key_traldf_c2d
157          aeiu_tm  (:,:)         = aeiu_tm  (:,:)         + aeiu(:,:)
158          aeiv_tm  (:,:)         = aeiv_tm  (:,:)         + aeiv(:,:)
159          aeiw_tm  (:,:)         = aeiw_tm  (:,:)         + aeiw(:,:)
160#  elif defined key_traldf_c1d
161          aeiu_tm  (:)           = aeiu_tm  (:,:)         + aeiu(:)
162          aeiv_tm  (:)           = aeiv_tm  (:,:)         + aeiv(:)
163          aeiw_tm  (:)           = aeiw_tm  (:,:)         + aeiw(:)
164#  else
165          aeiu_tm                = aeiu_tm                + aeiu
166          aeiv_tm                = aeiv_tm                + aeiv
167          aeiw_tm                = aeiw_tm                + aeiw
168#  endif
169#endif
170
171      ELSE                           !  It is time to substep
172         !   1. set temporary arrays to hold physics variables
173         un_temp    (:,:,:)      = un    (:,:,:)
174         vn_temp    (:,:,:)      = vn    (:,:,:)
175         wn_temp    (:,:,:)      = wn    (:,:,:)
176         tsn_temp   (:,:,:,:)    = tsn   (:,:,:,:)
177         rhop_temp  (:,:,:)      = rhop  (:,:,:)   
178         avt_temp   (:,:,:)      = avt   (:,:,:)
179# if defined key_zdfddm
180         avs_temp   (:,:,:)      = avs   (:,:,:)
181# endif
182#if defined key_ldfslp
183         wslpi_temp (:,:,:)      = wslpi (:,:,:)
184         wslpj_temp (:,:,:)      = wslpj (:,:,:)
185         uslp_temp  (:,:,:)      = uslp  (:,:,:)
186         vslp_temp  (:,:,:)      = vslp  (:,:,:)
187#endif
188# if defined key_trabbl
189          IF( nn_bbl_ldf == 1 ) THEN
190             ahu_bbl_temp(:,:)   = ahu_bbl(:,:) 
191             ahv_bbl_temp(:,:)   = ahv_bbl(:,:) 
192          ENDIF
193          IF( nn_bbl_adv == 1 ) THEN
194             utr_bbl_temp(:,:)   = utr_bbl(:,:) 
195             vtr_bbl_temp(:,:)   = vtr_bbl(:,:) 
196          ENDIF
197# endif
198         sshn_temp  (:,:)        = sshn  (:,:)
199         sshu_n_temp(:,:)        = sshu_n(:,:)
200         sshv_n_temp(:,:)        = sshv_n(:,:)
201         sshf_n_temp(:,:)        = sshf_n(:,:)
202         sshb_temp  (:,:)        = sshb  (:,:)
203         sshu_b_temp(:,:)        = sshu_b(:,:)
204         sshv_b_temp(:,:)        = sshv_b(:,:)
205         ssha_temp  (:,:)        = ssha  (:,:)
206         sshu_a_temp(:,:)        = sshu_a(:,:)
207         sshv_a_temp(:,:)        = sshv_a(:,:)
208         rnf_temp   (:,:)        = rnf   (:,:)
209         h_rnf_temp (:,:)        = h_rnf (:,:)
210         hmld_temp  (:,:)        = hmld  (:,:)
211         fr_i_temp  (:,:)        = fr_i  (:,:)
212         emp_temp   (:,:)        = emp   (:,:)
213         emp_b_temp (:,:)        = emp_b (:,:)
214         fmmflx_temp(:,:)        = fmmflx(:,:)
215         qsr_temp   (:,:)        = qsr   (:,:)
216         wndm_temp  (:,:)        = wndm  (:,:)
217#if defined key_traldf_c3d
218         ahtu_temp  (:,:,:)      = ahtu  (:,:,:)
219         ahtv_temp  (:,:,:)      = ahtv  (:,:,:)
220         ahtw_temp  (:,:,:)      = ahtw  (:,:,:)
221         ahtt_temp  (:,:,:)      = ahtt  (:,:,:)
222#elif defined key_traldf_c2d
223         ahtu_temp  (:,:)        = ahtu  (:,:)
224         ahtv_temp  (:,:)        = ahtv  (:,:)
225         ahtw_temp  (:,:)        = ahtw  (:,:)
226         ahtt_temp  (:,:)        = ahtt  (:,:)
227#elif defined key_traldf_c1d
228         ahtu_temp  (:)          = ahtu  (:)
229         ahtv_temp  (:)          = ahtv  (:)
230         ahtw_temp  (:)          = ahtw  (:)
231         ahtt_temp  (:)          = ahtt  (:)
232#else
233         ahtu_temp               = ahtu
234         ahtv_temp               = ahtv
235         ahtw_temp               = ahtw
236         ahtt_temp               = ahtt
237#endif
238
239#if defined key_traldf_eiv
240# if defined key_traldf_c3d
241         aeiu_temp  (:,:,:)      = aeiu(:,:,:)
242         aeiv_temp  (:,:,:)      = aeiv(:,:,:)
243         aeiw_temp  (:,:,:)      = aeiw(:,:,:)
244# elif defined key_traldf_c2d
245         aeiu_temp  (:,:)        = aeiu(:,:)
246         aeiv_temp  (:,:)        = aeiv(:,:)
247         aeiw_temp  (:,:)        = aeiw(:,:)
248# elif defined key_traldf_c1d
249         aeiu_temp  (:)          = aeiu(:)
250         aeiv_temp  (:)          = aeiv(:)
251         aeiw_temp  (:)          = aeiw(:)
252# else
253         aeiu_temp               = aeiu
254         aeiv_temp               = aeiv
255         aeiw_temp               = aeiw
256# endif
257#endif
258         !                                    !  Variables reset in trc_sub_ssh
259         rotn_temp  (:,:,:)      = rotn  (:,:,:)
260         hdivn_temp (:,:,:)      = hdivn (:,:,:)
261         rotb_temp  (:,:,:)      = rotb  (:,:,:)
262         hdivb_temp (:,:,:)      = hdivb (:,:,:)
263         hu_temp    (:,:)        = hu    (:,:)
264         hv_temp    (:,:)        = hv    (:,:)
265         hur_temp   (:,:)        = hur   (:,:)
266         hvr_temp   (:,:)        = hvr   (:,:)
267         !
268         DO jk = 1, jpk
269            e3t_temp(:,:,jk)     = fse3t(:,:,jk)
270            e3u_temp(:,:,jk)     = fse3u(:,:,jk)
271            e3v_temp(:,:,jk)     = fse3v(:,:,jk)
272            e3w_temp(:,:,jk)     = fse3w(:,:,jk)
273         ENDDO
274         IF( lk_vvl ) THEN                      !  Update Now Vertical coord.  !   (only in vvl case)
275           !                                    !------------------------------!
276           DO jk = 1, jpk
277              fse3t (:,:,jk)     = fse3t_n (:,:,jk)   ! vertical scale factors stored in fse3. arrays
278              fse3u (:,:,jk)     = fse3u_n (:,:,jk)
279              fse3v (:,:,jk)     = fse3v_n (:,:,jk)
280              fse3w (:,:,jk)     = fse3w_n (:,:,jk)
281           END DO
282         ENDIF
283
284         ! 2. Create averages and reassign variables
285         un_tm    (:,:,:)        = un_tm   (:,:,:)        + un   (:,:,:)        * e3u_temp(:,:,:) 
286         vn_tm    (:,:,:)        = vn_tm   (:,:,:)        + vn   (:,:,:)        * e3v_temp(:,:,:) 
287         tsn_tm   (:,:,:,jp_tem) = tsn_tm  (:,:,:,jp_tem) + tsn  (:,:,:,jp_tem) * e3t_temp(:,:,:) 
288         tsn_tm   (:,:,:,jp_sal) = tsn_tm  (:,:,:,jp_sal) + tsn  (:,:,:,jp_sal) * e3t_temp(:,:,:) 
289         rhop_tm (:,:,:)         = rhop_tm (:,:,:)        + rhop (:,:,:)        * e3t_temp(:,:,:) 
290         avt_tm   (:,:,:)        = avt_tm  (:,:,:)        + avt  (:,:,:)        * e3w_temp(:,:,:) 
291# if defined key_zdfddm
292         avs_tm   (:,:,:)        = avs_tm  (:,:,:)        + avs  (:,:,:)        * e3w_temp(:,:,:) 
293# endif
294#if defined key_ldfslp
295         wslpi_tm (:,:,:)        = wslpi_tm(:,:,:)        + wslpi(:,:,:)        * e3w_temp(:,:,:) 
296         wslpj_tm (:,:,:)        = wslpj_tm(:,:,:)        + wslpj(:,:,:)        * e3w_temp(:,:,:) 
297         uslp_tm  (:,:,:)        = uslp_tm (:,:,:)        + uslp (:,:,:)        * e3u_temp(:,:,:) 
298         vslp_tm  (:,:,:)        = vslp_tm (:,:,:)        + vslp (:,:,:)        * e3v_temp(:,:,:) 
299#endif
300# if defined key_trabbl
301          IF( nn_bbl_ldf == 1 ) THEN
302             ahu_bbl_tm(:,:)     = ahu_bbl_tm(:,:)        + ahu_bbl(:,:) 
303             ahv_bbl_tm(:,:)     = ahv_bbl_tm(:,:)        + ahv_bbl(:,:) 
304          ENDIF
305          IF( nn_bbl_adv == 1 ) THEN
306             utr_bbl_tm(:,:)     = utr_bbl_tm(:,:)        + utr_bbl(:,:) 
307             vtr_bbl_tm(:,:)     = vtr_bbl_tm(:,:)        + vtr_bbl(:,:) 
308          ENDIF
309# endif
310         sshn_tm  (:,:)          = sshn_tm    (:,:)       + sshn  (:,:) 
311         sshu_n_tm(:,:)          = sshu_n_tm  (:,:)       + sshu_n(:,:) 
312         sshv_n_tm(:,:)          = sshv_n_tm  (:,:)       + sshv_n(:,:) 
313         rnf_tm   (:,:)          = rnf_tm     (:,:)       + rnf   (:,:) 
314         h_rnf_tm (:,:)          = h_rnf_tm   (:,:)       + h_rnf (:,:) 
315         hmld_tm  (:,:)          = hmld_tm    (:,:)       + hmld  (:,:)
316         fr_i_tm  (:,:)          = fr_i_tm    (:,:)       + fr_i  (:,:)
317         emp_tm   (:,:)          = emp_tm     (:,:)       + emp   (:,:) 
318         fmmflx_tm(:,:)          = fmmflx_tm  (:,:)       + fmmflx(:,:)
319         qsr_tm   (:,:)          = qsr_tm     (:,:)       + qsr   (:,:)
320         wndm_tm  (:,:)          = wndm_tm    (:,:)       + wndm  (:,:)
321         !
322         sshn     (:,:)          = sshn_tm    (:,:) * r1_ndttrcp1 
323         sshu_n   (:,:)          = sshu_n_tm  (:,:) * r1_ndttrcp1 
324         sshv_n   (:,:)          = sshv_n_tm  (:,:) * r1_ndttrcp1 
325         sshb     (:,:)          = sshb_hold  (:,:)
326         sshu_b   (:,:)          = sshu_b_hold(:,:)
327         sshv_b   (:,:)          = sshv_b_hold(:,:)
328         rnf      (:,:)          = rnf_tm     (:,:) * r1_ndttrcp1 
329         h_rnf    (:,:)          = h_rnf_tm   (:,:) * r1_ndttrcp1 
330         hmld     (:,:)          = hmld_tm    (:,:) * r1_ndttrcp1 
331         !  variables that are initialized after averages initialized
332         emp_b    (:,:) = emp_b_hold (:,:)
333         IF( kt == nittrc000 ) THEN
334            wndm  (:,:)          = wndm_tm    (:,:) * r1_ndttrc 
335            qsr   (:,:)          = qsr_tm     (:,:) * r1_ndttrc 
336            emp   (:,:)          = emp_tm     (:,:) * r1_ndttrc 
337            fmmflx(:,:)          = fmmflx_tm  (:,:) * r1_ndttrc 
338            fr_i  (:,:)          = fr_i_tm    (:,:) * r1_ndttrc
339# if defined key_trabbl
340            IF( nn_bbl_ldf == 1 ) THEN
341               ahu_bbl(:,:)      = ahu_bbl_tm (:,:) * r1_ndttrc 
342               ahv_bbl(:,:)      = ahv_bbl_tm (:,:) * r1_ndttrc 
343            ENDIF
344            IF( nn_bbl_adv == 1 ) THEN
345               utr_bbl(:,:)      = utr_bbl_tm (:,:) * r1_ndttrc 
346               vtr_bbl(:,:)      = vtr_bbl_tm (:,:) * r1_ndttrc 
347            ENDIF
348# endif
349         ELSE
350            wndm  (:,:)          = wndm_tm    (:,:) * r1_ndttrcp1 
351            qsr   (:,:)          = qsr_tm     (:,:) * r1_ndttrcp1 
352            emp   (:,:)          = emp_tm     (:,:) * r1_ndttrcp1 
353            fmmflx(:,:)          = fmmflx_tm  (:,:) * r1_ndttrcp1 
354            fr_i  (:,:)          = fr_i_tm    (:,:) * r1_ndttrcp1 
355# if defined key_trabbl
356            IF( nn_bbl_ldf == 1 ) THEN
357               ahu_bbl(:,:)      = ahu_bbl_tm (:,:) * r1_ndttrcp1 
358               ahv_bbl(:,:)      = ahv_bbl_tm (:,:) * r1_ndttrcp1 
359            ENDIF
360            IF( nn_bbl_adv == 1 ) THEN
361               utr_bbl(:,:)      = utr_bbl_tm (:,:) * r1_ndttrcp1 
362               vtr_bbl(:,:)      = vtr_bbl_tm (:,:) * r1_ndttrcp1 
363            ENDIF
364# endif
365         ENDIF
366         !
367         DO jk = 1, jpk
368            DO jj = 1, jpj
369               DO ji = 1, jpi
370                  z1_ne3t = r1_ndttrcp1  / fse3t(ji,jj,jk)
371                  z1_ne3u = r1_ndttrcp1  / fse3u(ji,jj,jk)
372                  z1_ne3v = r1_ndttrcp1  / fse3v(ji,jj,jk)
373                  z1_ne3w = r1_ndttrcp1  / fse3w(ji,jj,jk)
374                  !
375                  un   (ji,jj,jk)        = un_tm   (ji,jj,jk)        * z1_ne3u
376                  vn   (ji,jj,jk)        = vn_tm   (ji,jj,jk)        * z1_ne3v
377                  tsn  (ji,jj,jk,jp_tem) = tsn_tm  (ji,jj,jk,jp_tem) * z1_ne3t
378                  tsn  (ji,jj,jk,jp_sal) = tsn_tm  (ji,jj,jk,jp_sal) * z1_ne3t
379                  rhop (ji,jj,jk)        = rhop_tm (ji,jj,jk)        * z1_ne3t
380                  avt  (ji,jj,jk)        = avt_tm  (ji,jj,jk)        * z1_ne3w
381# if defined key_zdfddm
382                  avs  (ji,jj,jk)        = avs_tm  (ji,jj,jk)        * z1_ne3w
383# endif
384#if defined key_ldfslp
385                  wslpi(ji,jj,jk)        = wslpi_tm(ji,jj,jk)        * z1_ne3w 
386                  wslpj(ji,jj,jk)        = wslpj_tm(ji,jj,jk)        * z1_ne3w 
387                  uslp (ji,jj,jk)        = uslp_tm (ji,jj,jk)        * z1_ne3u 
388                  vslp (ji,jj,jk)        = vslp_tm (ji,jj,jk)        * z1_ne3v 
389#endif
390               ENDDO
391            ENDDO
392         ENDDO
393
394#if defined key_traldf_c3d
395          ahtt_tm  (:,:,:)       = ahtt_tm  (:,:,:)       + ahtt(:,:,:)         * e3t_temp(:,:,:)
396          ahtu_tm  (:,:,:)       = ahtu_tm  (:,:,:)       + ahtu(:,:,:)         * e3u_temp(:,:,:)
397          ahtv_tm  (:,:,:)       = ahtv_tm  (:,:,:)       + ahtv(:,:,:)         * e3v_temp(:,:,:)
398          ahtw_tm  (:,:,:)       = ahtw_tm  (:,:,:)       + ahtw(:,:,:)         * e3w_temp(:,:,:)
399          !
400          ahtt     (:,:,:)       = ahtt_tm  (:,:,:) * r1_ndttrcp1  / fse3t(:,:,:)
401          ahtu     (:,:,:)       = ahtu_tm  (:,:,:) * r1_ndttrcp1  / fse3u(:,:,:)
402          ahtv     (:,:,:)       = ahtv_tm  (:,:,:) * r1_ndttrcp1  / fse3v(:,:,:)
403          ahtw     (:,:,:)       = ahtw_tm  (:,:,:) * r1_ndttrcp1  / fse3w(:,:,:)
404#elif defined key_traldf_c2d
405          ahtt_tm  (:,:)         = ahtt_tm  (:,:)         + ahtt(:,:)
406          ahtu_tm  (:,:)         = ahtu_tm  (:,:)         + ahtu(:,:)
407          ahtv_tm  (:,:)         = ahtv_tm  (:,:)         + ahtv(:,:)
408          ahtw_tm  (:,:)         = ahtw_tm  (:,:)         + ahtw(:,:)
409          !
410          ahtt     (:,:)         = ahtt_tm  (:,:)   * r1_ndttrcp1
411          ahtu     (:,:)         = ahtu_tm  (:,:)   * r1_ndttrcp1
412          ahtv     (:,:)         = ahtv_tm  (:,:)   * r1_ndttrcp1
413          ahtw     (:,:)         = ahtw_tm  (:,:)   * r1_ndttrcp1
414#elif defined key_traldf_c1d
415          ahtt_tm  (:)           = ahtt_tm  (:,:)         + ahtt(:)
416          ahtu_tm  (:)           = ahtu_tm  (:,:)         + ahtu(:)
417          ahtv_tm  (:)           = ahtv_tm  (:,:)         + ahtv(:)
418          ahtw_tm  (:)           = ahtw_tm  (:,:)         + ahtw(:)
419          !
420          ahtt     (:)           = ahtt_tm  (:)     * r1_ndttrcp1
421          ahtu     (:)           = ahtu_tm  (:)     * r1_ndttrcp1
422          ahtv     (:)           = ahtv_tm  (:)     * r1_ndttrcp1
423          ahtw     (:)           = ahtw_tm  (:)     * r1_ndttrcp1
424#else
425          ahtt_tm                = ahtt_tm                + ahtt
426          ahtu_tm                = ahtu_tm                + ahtu
427          ahtv_tm                = ahtv_tm                + ahtv
428          ahtw_tm                = ahtw_tm                + ahtw
429          !
430          ahtt                   = ahtt_tm          * r1_ndttrcp1
431          ahtu                   = ahtu_tm          * r1_ndttrcp1
432          ahtv                   = ahtv_tm          * r1_ndttrcp1
433          ahtw                   = ahtw_tm          * r1_ndttrcp1
434#endif
435
436#if defined key_traldf_eiv
437# if defined key_traldf_c3d
438          aeiu_tm  (:,:,:)       = aeiu_tm  (:,:,:)       + aeiu(:,:,:)         * e3u_temp(:,:,:)
439          aeiv_tm  (:,:,:)       = aeiv_tm  (:,:,:)       + aeiv(:,:,:)         * e3v_temp(:,:,:)
440          aeiw_tm  (:,:,:)       = aeiw_tm  (:,:,:)       + aeiw(:,:,:)         * e3w_temp(:,:,:)
441          !
442          aeiu     (:,:,:)       = aeiu_tm  (:,:,:) * r1_ndttrcp1  / fse3u(:,:,:)
443          aeiv     (:,:,:)       = aeiv_tm  (:,:,:) * r1_ndttrcp1  / fse3v(:,:,:)
444          aeiw     (:,:,:)       = aeiw_tm  (:,:,:) * r1_ndttrcp1  / fse3w(:,:,:)
445# elif defined key_traldf_c2d
446          aeiu_tm  (:,:)         = aeiu_tm  (:,:)         + aeiu(:,:)
447          aeiv_tm  (:,:)         = aeiv_tm  (:,:)         + aeiv(:,:)
448          aeiw_tm  (:,:)         = aeiw_tm  (:,:)         + aeiw(:,:)
449          !
450          aeiu     (:,:)         = aeiu_tm  (:,:)   * r1_ndttrcp1
451          aeiv     (:,:)         = aeiv_tm  (:,:)   * r1_ndttrcp1
452          aeiw     (:,:)         = aeiw_tm  (:,:)   * r1_ndttrcp1
453# elif defined key_traldf_c1d
454          aeiu_tm  (:)           = aeiu_tm  (:,:)         + aeiu(:)
455          aeiv_tm  (:)           = aeiv_tm  (:,:)         + aeiv(:)
456          aeiw_tm  (:)           = aeiw_tm  (:,:)         + aeiw(:)
457          !
458          aeiu     (:)           = aeiu_tm  (:)     * r1_ndttrcp1
459          aeiv     (:)           = aeiv_tm  (:)     * r1_ndttrcp1
460          aeiw     (:)           = aeiw_tm  (:)     * r1_ndttrcp1
461# else
462          aeiu_tm                = aeiu_tm                + aeiu
463          aeiv_tm                = aeiv_tm                + aeiv
464          aeiw_tm                = aeiw_tm                + aeiw
465          !
466          aeiu                   = aeiu_tm          * r1_ndttrcp1
467          aeiv                   = aeiv_tm          * r1_ndttrcp1
468          aeiw                   = aeiw_tm          * r1_ndttrcp1
469# endif
470#endif
471
472         CALL lbc_lnk( un    (:,:,:)       , 'U',-1. ) 
473         CALL lbc_lnk( vn    (:,:,:)       , 'V',-1. ) 
474         CALL lbc_lnk( tsn   (:,:,:,jp_tem), 'T', 1. ) 
475         CALL lbc_lnk( tsn   (:,:,:,jp_sal), 'T', 1. ) 
476         CALL lbc_lnk( rhop  (:,:,:)       , 'T', 1. ) 
477         CALL lbc_lnk( avt   (:,:,:)       , 'W', 1. ) 
478# if defined key_zdfddm
479          CALL lbc_lnk( avs  (:,:,:)       , 'W', 1. ) 
480# endif
481#if defined key_ldfslp
482         CALL lbc_lnk( uslp  (:,:,:)       , 'U',-1. ) 
483         CALL lbc_lnk( vslp  (:,:,:)       , 'V',-1. ) 
484         CALL lbc_lnk( wslpi (:,:,:)       , 'W',-1. ) 
485         CALL lbc_lnk( wslpj (:,:,:)       , 'W',-1. ) 
486#endif
487         CALL lbc_lnk( sshn  (:,:)         , 'T', 1. ) 
488         CALL lbc_lnk( sshu_n(:,:)         , 'U', 1. ) 
489         CALL lbc_lnk( sshv_n(:,:)         , 'V', 1. ) 
490         CALL lbc_lnk( sshf_n(:,:)         , 'F', 1. ) 
491         CALL lbc_lnk( sshb  (:,:)         , 'T', 1. ) 
492         CALL lbc_lnk( sshu_b(:,:)         , 'U', 1. ) 
493         CALL lbc_lnk( sshv_b(:,:)         , 'V', 1. ) 
494         CALL lbc_lnk( ssha  (:,:)         , 'T', 1. ) 
495         CALL lbc_lnk( sshu_a(:,:)         , 'U', 1. ) 
496         CALL lbc_lnk( sshv_a(:,:)         , 'V', 1. ) 
497         CALL lbc_lnk( rnf   (:,:)         , 'T', 1. ) 
498         CALL lbc_lnk( h_rnf (:,:)         , 'T', 1. ) 
499         CALL lbc_lnk( hmld  (:,:)         , 'T', 1. ) 
500         CALL lbc_lnk( fr_i  (:,:)         , 'T', 1. ) 
501         CALL lbc_lnk( emp   (:,:)         , 'T', 1. ) 
502         CALL lbc_lnk( emp_b (:,:)         , 'T', 1. ) 
503         CALL lbc_lnk( fmmflx(:,:)         , 'T', 1. ) 
504         CALL lbc_lnk( qsr   (:,:)         , 'T', 1. ) 
505         CALL lbc_lnk( wndm  (:,:)         , 'T', 1. ) 
506# if defined key_trabbl
507         IF( nn_bbl_ldf == 1 ) THEN
508            CALL lbc_lnk( ahu_bbl(:,:)     , 'U', 1. ) 
509            CALL lbc_lnk( ahv_bbl(:,:)     , 'v', 1. ) 
510         ENDIF
511         IF( nn_bbl_adv == 1 ) THEN
512            CALL lbc_lnk( utr_bbl(:,:)     , 'U', 1. ) 
513            CALL lbc_lnk( vtr_bbl(:,:)     , 'U', 1. ) 
514         ENDIF
515# endif
516#if defined key_traldf_c3d
517         CALL lbc_lnk( ahtt  (:,:,:)       , 'T', 1. ) 
518         CALL lbc_lnk( ahtu  (:,:,:)       , 'U', 1. ) 
519         CALL lbc_lnk( ahtv  (:,:,:)       , 'V', 1. ) 
520         CALL lbc_lnk( ahtw  (:,:,:)       , 'W', 1. ) 
521#elif defined key_traldf_c2d
522         CALL lbc_lnk( ahtt  (:,:)         , 'T', 1. ) 
523         CALL lbc_lnk( ahtu  (:,:)         , 'U', 1. ) 
524         CALL lbc_lnk( ahtv  (:,:)         , 'V', 1. ) 
525         CALL lbc_lnk( ahtw  (:,:)         , 'W', 1. ) 
526#endif
527#if defined key_traldf_eiv
528#if defined key_traldf_c3d
529         CALL lbc_lnk( aeiu  (:,:,:)       , 'U', 1. ) 
530         CALL lbc_lnk( aeiv  (:,:,:)       , 'V', 1. ) 
531         CALL lbc_lnk( aeiw  (:,:,:)       , 'W', 1. ) 
532#elif defined key_traldf_c2d
533         CALL lbc_lnk( aeiu  (:,:)         , 'U', 1. ) 
534         CALL lbc_lnk( aeiv  (:,:)         , 'V', 1. ) 
535         CALL lbc_lnk( aeiw  (:,:)         , 'W', 1. ) 
536#endif
537#endif
538         !
539         CALL trc_sub_ssh( kt )         ! after ssh & vertical velocity
540         !
541         CALL lbc_lnk( wn    (:,:,:)       , 'W',-1. ) 
542         CALL lbc_lnk( rotn  (:,:,:)       , 'F', 1. ) 
543         CALL lbc_lnk( hdivn (:,:,:)       , 'T', 1. ) 
544         CALL lbc_lnk( rotb  (:,:,:)       , 'F', 1. ) 
545         CALL lbc_lnk( hdivb (:,:,:)       , 'T', 1. ) 
546         CALL lbc_lnk( hu    (:,:)         , 'U', 1. ) 
547         CALL lbc_lnk( hv    (:,:)         , 'V', 1. ) 
548         CALL lbc_lnk( hur   (:,:)         , 'U', 1. ) 
549         CALL lbc_lnk( hvr   (:,:)         , 'V', 1. ) 
550      ENDIF
551      !
552      IF( nn_timing == 1 )  CALL timing_start('trc_sub_stp')
553      !
554   END SUBROUTINE trc_sub_stp
555
556   SUBROUTINE trc_sub_ini
557      !!-------------------------------------------------------------------
558      !!                     ***  ROUTINE trc_sub_ini  ***
559      !!                     
560      !! ** Purpose : Initialize variables needed for sub-stepping passive tracers
561      !!
562      !! ** Method  :
563      !!              Compute the averages for sub-stepping
564      !!-------------------------------------------------------------------
565      INTEGER ::   ierr
566      !!-------------------------------------------------------------------
567      !
568      IF( nn_timing == 1 )  CALL timing_start('trc_sub_ini')
569      !
570      IF(lwp) WRITE(numout,*)
571      IF(lwp) WRITE(numout,*) 'trc_sub_ini : initial set up of the passive tracers substepping'
572      IF(lwp) WRITE(numout,*) '~~~~~~~'
573
574      ierr =  trc_sub_alloc    ()
575      IF( lk_mpp    )   CALL mpp_sum( ierr )
576      IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'top_sub_alloc : unable to allocate standard ocean arrays' )
577
578      un_tm   (:,:,:)        = un   (:,:,:)        * fse3u(:,:,:) 
579      vn_tm   (:,:,:)        = vn   (:,:,:)        * fse3v(:,:,:) 
580      tsn_tm  (:,:,:,jp_tem) = tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
581      tsn_tm  (:,:,:,jp_sal) = tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
582      rhop_tm (:,:,:)        = rhop (:,:,:)        * fse3t(:,:,:) 
583      avt_tm  (:,:,:)        = avt  (:,:,:)        * fse3w(:,:,:) 
584# if defined key_zdfddm
585      avs_tm  (:,:,:)        = avs  (:,:,:)        * fse3w(:,:,:) 
586# endif
587#if defined key_ldfslp
588      wslpi_tm(:,:,:)        = wslpi(:,:,:)        * fse3w(:,:,:) 
589      wslpj_tm(:,:,:)        = wslpj(:,:,:)        * fse3w(:,:,:) 
590      uslp_tm (:,:,:)        = uslp (:,:,:)        * fse3u(:,:,:) 
591      vslp_tm (:,:,:)        = vslp (:,:,:)        * fse3v(:,:,:) 
592#endif
593      sshn_tm  (:,:) = sshn  (:,:) 
594      sshu_n_tm(:,:) = sshu_n(:,:) 
595      sshv_n_tm(:,:) = sshv_n(:,:) 
596      rnf_tm   (:,:) = rnf   (:,:) 
597      h_rnf_tm (:,:) = h_rnf (:,:) 
598      hmld_tm  (:,:) = hmld  (:,:)
599
600      ! Physics variables that are set after initialization:
601      fr_i_tm(:,:) = 0._wp
602      emp_tm (:,:) = 0._wp
603      fmmflx_tm(:,:)  = 0._wp
604      qsr_tm (:,:) = 0._wp
605      wndm_tm(:,:) = 0._wp
606# if defined key_trabbl
607      IF( nn_bbl_ldf == 1 ) THEN
608         ahu_bbl_tm(:,:) = 0._wp
609         ahv_bbl_tm(:,:) = 0._wp
610      ENDIF
611      IF( nn_bbl_adv == 1 ) THEN
612         utr_bbl_tm(:,:) = 0._wp
613         vtr_bbl_tm(:,:) = 0._wp
614      ENDIF
615# endif
616      !
617#if defined key_traldf_c3d
618      ahtt_tm(:,:,:) = ahtt(:,:,:) * fse3t(:,:,:)
619      ahtu_tm(:,:,:) = ahtu(:,:,:) * fse3u(:,:,:)
620      ahtv_tm(:,:,:) = ahtv(:,:,:) * fse3v(:,:,:)
621      ahtw_tm(:,:,:) = ahtw(:,:,:) * fse3w(:,:,:)
622#elif defined key_traldf_c2d
623      ahtt_tm(:,:)   = ahtt(:,:)
624      ahtu_tm(:,:)   = ahtu(:,:)
625      ahtv_tm(:,:)   = ahtv(:,:)
626      ahtw_tm(:,:)   = ahtw(:,:)
627#elif defined key_traldf_c1d
628      ahtt_tm(:)     = ahtt(:)
629      ahtu_tm(:)     = ahtu(:)
630      ahtv_tm(:)     = ahtv(:)
631      ahtw_tm(:)     = ahtw(:)
632#else
633      ahtt_tm        = ahtt
634      ahtu_tm        = ahtu
635      ahtv_tm        = ahtv
636      ahtw_tm        = ahtw
637#endif
638      !
639#if defined key_traldf_eiv
640#  if defined key_traldf_c3d
641      aeiu_tm(:,:,:) = aeiu(:,:,:) * fse3u(:,:,:)
642      aeiv_tm(:,:,:) = aeiv(:,:,:) * fse3v(:,:,:)
643      aeiw_tm(:,:,:) = aeiw(:,:,:) * fse3w(:,:,:)
644#  elif defined key_traldf_c2d
645      aeiu_tm(:,:)   = aeiu(:,:)
646      aeiv_tm(:,:)   = aeiv(:,:)
647      aeiw_tm(:,:)   = aeiw(:,:)
648#  elif defined key_traldf_c1d
649      aeiu_tm(:)     = aeiu(:)
650      aeiv_tm(:)     = aeiv(:)
651      aeiw_tm(:)     = aeiw(:)
652#  else
653      aeiu_tm        = aeiu
654      aeiv_tm        = aeiv
655      aeiw_tm        = aeiw
656#  endif
657#endif
658      !
659      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_ini')
660      !
661   END SUBROUTINE trc_sub_ini
662
663   SUBROUTINE trc_sub_reset( kt )
664      !!-------------------------------------------------------------------
665      !!                     ***  ROUTINE trc_sub_reset  ***
666      !!                     
667      !! ** Purpose : Reset physics variables averaged for substepping
668      !!
669      !! ** Method  :
670      !!              Compute the averages for sub-stepping
671      !!-------------------------------------------------------------------
672      INTEGER, INTENT( in ) ::  kt  ! ocean time-step index
673      INTEGER :: jk                 ! dummy loop indices
674      !!-------------------------------------------------------------------
675      !
676      IF( nn_timing == 1 )  CALL timing_start('trc_sub_reset')
677      !
678      !   restore physics variables
679      un    (:,:,:)   =  un_temp    (:,:,:)
680      vn    (:,:,:)   =  vn_temp    (:,:,:)
681      wn    (:,:,:)   =  wn_temp    (:,:,:)
682      tsn   (:,:,:,:) =  tsn_temp   (:,:,:,:)
683      rhop  (:,:,:)   =  rhop_temp  (:,:,:)
684      avt   (:,:,:)   =  avt_temp   (:,:,:)
685# if defined key_zdfddm
686      avs   (:,:,:)   =  avs_temp   (:,:,:)
687# endif
688#if defined key_ldfslp
689      wslpi (:,:,:)   =  wslpi_temp (:,:,:)
690      wslpj (:,:,:)   =  wslpj_temp (:,:,:)
691      uslp  (:,:,:)   =  uslp_temp  (:,:,:)
692      vslp  (:,:,:)   =  vslp_temp  (:,:,:)
693#endif
694      sshn  (:,:)     =  sshn_temp  (:,:)
695      sshb  (:,:)     =  sshb_temp  (:,:)
696      ssha  (:,:)     =  ssha_temp  (:,:)
697      sshu_n(:,:)     =  sshu_n_temp(:,:)
698      sshu_b(:,:)     =  sshu_b_temp(:,:)
699      sshu_a(:,:)     =  sshu_a_temp(:,:)
700      sshv_n(:,:)     =  sshv_n_temp(:,:)
701      sshv_b(:,:)     =  sshv_b_temp(:,:)
702      sshv_a(:,:)     =  sshv_a_temp(:,:)
703      sshf_n(:,:)     =  sshf_n_temp(:,:)
704      rnf   (:,:)     =  rnf_temp   (:,:)
705      h_rnf (:,:)     =  h_rnf_temp (:,:)
706      !
707      hmld  (:,:)     =  hmld_temp  (:,:)
708      fr_i  (:,:)     =  fr_i_temp  (:,:)
709      emp   (:,:)     =  emp_temp   (:,:)
710      fmmflx(:,:)     =  fmmflx_temp(:,:)
711      emp_b (:,:)     =  emp_b_temp (:,:)
712      qsr   (:,:)     =  qsr_temp   (:,:)
713      wndm  (:,:)     =  wndm_temp  (:,:)
714# if defined key_trabbl
715      IF( nn_bbl_ldf == 1 ) THEN
716         ahu_bbl(:,:) = ahu_bbl_temp(:,:) 
717         ahv_bbl(:,:) = ahv_bbl_temp(:,:) 
718      ENDIF
719      IF( nn_bbl_adv == 1 ) THEN
720         utr_bbl(:,:) = utr_bbl_temp(:,:) 
721         vtr_bbl(:,:) = vtr_bbl_temp(:,:) 
722      ENDIF
723# endif
724      !
725#if defined key_traldf_c3d
726      ahtu  (:,:,:)   =  ahtu_temp  (:,:,:)
727      ahtv  (:,:,:)   =  ahtv_temp  (:,:,:)
728      ahtw  (:,:,:)   =  ahtw_temp  (:,:,:)
729      ahtt  (:,:,:)   =  ahtt_temp  (:,:,:)
730#elif defined key_traldf_c2d
731      ahtu  (:,:)     =  ahtu_temp  (:,:)
732      ahtv  (:,:)     =  ahtv_temp  (:,:)
733      ahtw  (:,:)     =  ahtw_temp  (:,:)
734      ahtt  (:,:)     =  ahtt_temp  (:,:)
735#elif defined key_traldf_c1d
736      ahtu  (:)       =  ahtu_temp  (:)
737      ahtv  (:)       =  ahtv_temp  (:)
738      ahtw  (:)       =  ahtw_temp  (:)
739      ahtt  (:)       =  ahtt_temp  (:)
740#else
741      ahtu            =  ahtu_temp
742      ahtv            =  ahtv_temp
743      ahtw            =  ahtw_temp
744      ahtt            =  ahtt_temp
745#endif
746      !
747#if defined key_traldf_eiv
748#if defined key_traldf_c3d
749      aeiu  (:,:,:)  =  aeiu_temp(:,:,:)
750      aeiv  (:,:,:)  =  aeiv_temp(:,:,:)
751      aeiw  (:,:,:)  =  aeiw_temp(:,:,:)
752#elif defined key_traldf_c2d
753      aeiu  (:,:)    =  aeiu_temp(:,:)
754      aeiv  (:,:)    =  aeiv_temp(:,:)
755      aeiw  (:,:)    =  aeiw_temp(:,:)
756#elif defined key_traldf_c1d
757      aeiu  (:)      =  aeiu_temp(:)
758      aeiv  (:)      =  aeiv_temp(:)
759      aeiw  (:)      =  aeiw_temp(:)
760#else
761      aeiu           =  aeiu_temp
762      aeiv           =  aeiv_temp
763      aeiw           =  aeiw_temp
764#endif
765#endif
766      hdivn (:,:,:)   =  hdivn_temp (:,:,:)
767      rotn  (:,:,:)   =  rotn_temp  (:,:,:)
768      hdivb (:,:,:)   =  hdivb_temp (:,:,:)
769      rotb  (:,:,:)   =  rotb_temp  (:,:,:)
770      hu    (:,:)     =  hu_temp    (:,:)
771      hv    (:,:)     =  hv_temp    (:,:)
772      hur   (:,:)     =  hur_temp   (:,:)
773      hvr   (:,:)     =  hvr_temp   (:,:)
774      !                                     
775      DO jk = 1, jpk
776         fse3t(:,:,jk)= e3t_temp(:,:,jk) 
777         fse3u(:,:,jk)= e3u_temp(:,:,jk) 
778         fse3v(:,:,jk)= e3v_temp(:,:,jk) 
779         fse3w(:,:,jk)= e3w_temp(:,:,jk) 
780      END DO
781      !                                           !------------------------------!
782      IF( lk_vvl ) THEN                           !  Update Now Vertical coord.  !   (only in vvl case)
783        !                                           !------------------------------!
784         DO jk = 1, jpkm1
785            fsdept(:,:,jk) = fsdept_n(:,:,jk)          ! now local depths stored in fsdep. arrays
786            fsdepw(:,:,jk) = fsdepw_n(:,:,jk)
787            fsde3w(:,:,jk) = fsde3w_n(:,:,jk)
788            !
789            fse3t (:,:,jk) = fse3t_n (:,:,jk)          ! vertical scale factors stored in fse3. arrays
790            fse3u (:,:,jk) = fse3u_n (:,:,jk)
791            fse3v (:,:,jk) = fse3v_n (:,:,jk)
792            fse3f (:,:,jk) = fse3f_n (:,:,jk)
793            fse3w (:,:,jk) = fse3w_n (:,:,jk)
794            fse3uw(:,:,jk) = fse3uw_n(:,:,jk)
795            fse3vw(:,:,jk) = fse3vw_n(:,:,jk)
796         END DO
797         !
798      ENDIF
799
800      ! Start new averages
801         un_tm   (:,:,:)        = un   (:,:,:)        * fse3u(:,:,:) 
802         vn_tm   (:,:,:)        = vn   (:,:,:)        * fse3v(:,:,:) 
803         tsn_tm  (:,:,:,jp_tem) = tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
804         tsn_tm  (:,:,:,jp_sal) = tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
805         rhop_tm (:,:,:)        = rhop (:,:,:)        * fse3t(:,:,:) 
806         avt_tm  (:,:,:)        = avt  (:,:,:)        * fse3w(:,:,:) 
807# if defined key_zdfddm
808         avs_tm  (:,:,:)        = avs  (:,:,:)        * fse3w(:,:,:) 
809# endif
810#if defined key_ldfslp
811         wslpi_tm(:,:,:)        = wslpi(:,:,:)        * fse3w(:,:,:) 
812         wslpj_tm(:,:,:)        = wslpj(:,:,:)        * fse3w(:,:,:) 
813         uslp_tm (:,:,:)        = uslp (:,:,:)        * fse3u(:,:,:) 
814         vslp_tm (:,:,:)        = vslp (:,:,:)        * fse3v(:,:,:) 
815#endif
816      !
817      sshb_hold  (:,:) = sshn  (:,:)
818      sshu_b_hold(:,:) = sshu_n(:,:)
819      sshv_b_hold(:,:) = sshv_n(:,:)
820      emp_b_hold (:,:) = emp   (:,:)
821      sshn_tm    (:,:) = sshn  (:,:) 
822      sshu_n_tm  (:,:) = sshu_n(:,:) 
823      sshv_n_tm  (:,:) = sshv_n(:,:) 
824      rnf_tm     (:,:) = rnf   (:,:) 
825      h_rnf_tm   (:,:) = h_rnf (:,:) 
826      hmld_tm    (:,:) = hmld  (:,:)
827      fr_i_tm    (:,:) = fr_i  (:,:)
828      emp_tm     (:,:) = emp   (:,:)
829      fmmflx_tm  (:,:) = fmmflx(:,:)
830      qsr_tm     (:,:) = qsr   (:,:)
831      wndm_tm    (:,:) = wndm  (:,:)
832# if defined key_trabbl
833      IF( nn_bbl_ldf == 1 ) THEN
834         ahu_bbl_tm(:,:) = ahu_bbl(:,:) 
835         ahv_bbl_tm(:,:) = ahv_bbl(:,:) 
836      ENDIF
837      IF( nn_bbl_adv == 1 ) THEN
838         utr_bbl_tm(:,:) = utr_bbl(:,:) 
839         vtr_bbl_tm(:,:) = vtr_bbl(:,:) 
840      ENDIF
841# endif
842      !
843#if defined key_traldf_c3d
844      DO jk = 1, jpkm1
845         ahtt_tm(:,:,jk) = ahtt(:,:,jk) * fse3t(:,:,jk)
846         ahtu_tm(:,:,jk) = ahtu(:,:,jk) * fse3u(:,:,jk)
847         ahtv_tm(:,:,jk) = ahtv(:,:,jk) * fse3v(:,:,jk)
848         ahtw_tm(:,:,jk) = ahtw(:,:,jk) * fse3w(:,:,jk)
849      END DO
850#elif defined key_traldf_c2d
851      ahtt_tm(:,:)   = ahtt(:,:)
852      ahtu_tm(:,:)   = ahtu(:,:)
853      ahtv_tm(:,:)   = ahtv(:,:)
854      ahtw_tm(:,:)   = ahtw(:,:)
855#elif defined key_traldf_c1d
856      ahtt_tm(:)     = ahtt(:)
857      ahtu_tm(:)     = ahtu(:)
858      ahtv_tm(:)     = ahtv(:)
859      ahtw_tm(:)     = ahtw(:)
860#else
861      ahtt_tm        = ahtt
862      ahtu_tm        = ahtu
863      ahtv_tm        = ahtv
864      ahtw_tm        = ahtw
865#endif
866      !
867#if defined key_traldf_eiv
868#  if defined key_traldf_c3d
869      DO jk = 1, jpk
870         aeiu_tm(:,:,jk) = aeiu(:,:,jk) * fse3u(:,:,jk)
871         aeiv_tm(:,:,jk) = aeiv(:,:,jk) * fse3v(:,:,jk)
872         aeiw_tm(:,:,jk) = aeiw(:,:,jk) * fse3w(:,:,jk)
873      END DO
874#  elif defined key_traldf_c2d
875      aeiu_tm(:,:)   = aeiu(:,:)
876      aeiv_tm(:,:)   = aeiv(:,:)
877      aeiw_tm(:,:)   = aeiw(:,:)
878#  elif defined key_traldf_c1d
879      aeiu_tm(:)     = aeiu(:)
880      aeiv_tm(:)     = aeiv(:)
881      aeiw_tm(:)     = aeiw(:)
882#  else
883      aeiu_tm        = aeiu
884      aeiv_tm        = aeiv
885      aeiw_tm        = aeiw
886#  endif
887#endif
888      !
889      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_reset')
890      !
891   END SUBROUTINE trc_sub_reset
892
893
894   SUBROUTINE trc_sub_ssh( kt ) 
895      !!----------------------------------------------------------------------
896      !!                ***  ROUTINE trc_sub_ssh  ***
897      !!                   
898      !! ** Purpose :   compute the after ssh (ssha), the now vertical velocity
899      !!              and update the now vertical coordinate (lk_vvl=T).
900      !!
901      !! ** Method  : - Using the incompressibility hypothesis, the vertical
902      !!      velocity is computed by integrating the horizontal divergence 
903      !!      from the bottom to the surface minus the scale factor evolution.
904      !!        The boundary conditions are w=0 at the bottom (no flux) and.
905      !!
906      !! ** action  :   ssha    : after sea surface height
907      !!                wn      : now vertical velocity
908      !!                sshu_a, sshv_a, sshf_a  : after sea surface height (lk_vvl=T)
909      !!                hu, hv, hur, hvr        : ocean depth and its inverse at u-,v-points
910      !!
911      !! Reference  : Leclair, M., and G. Madec, 2009, Ocean Modelling.
912      !!----------------------------------------------------------------------
913      !
914      INTEGER, INTENT(in) ::   kt   ! time step
915      !
916      INTEGER  ::   ji, jj, jk   ! dummy loop indices
917      REAL(wp) ::   zcoefu, zcoefv, zcoeff, z2dt, z1_2dt, z1_rau0   ! local scalars
918      REAL(wp), POINTER, DIMENSION(:,:) :: zhdiv
919      !!---------------------------------------------------------------------
920      !
921      IF( nn_timing == 1 )  CALL timing_start('trc_sub_ssh')
922      !
923      ! Allocate temporary workspace
924      CALL wrk_alloc( jpi, jpj, zhdiv )
925
926      IF( kt == nittrc000 ) THEN
927         !
928         IF(lwp) WRITE(numout,*)
929         IF(lwp) WRITE(numout,*) 'trc_sub_ssh : after sea surface height and now vertical velocity '
930         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
931         !
932         wn(:,:,jpk) = 0._wp                  ! bottom boundary condition: w=0 (set once for all)
933         !
934      ENDIF
935
936      !                                           !------------------------------------------!
937      IF( lk_vvl ) THEN                           !  Regridding: Update Now Vertical coord.  !   (only in vvl case)
938         !                                        !------------------------------------------!
939         DO jk = 1, jpkm1
940            fsdept(:,:,jk) = fsdept_n(:,:,jk)         ! now local depths stored in fsdep. arrays
941            fsdepw(:,:,jk) = fsdepw_n(:,:,jk)
942            fsde3w(:,:,jk) = fsde3w_n(:,:,jk)
943            !
944            fse3t (:,:,jk) = fse3t_n (:,:,jk)         ! vertical scale factors stored in fse3. arrays
945            fse3u (:,:,jk) = fse3u_n (:,:,jk)
946            fse3v (:,:,jk) = fse3v_n (:,:,jk)
947            fse3f (:,:,jk) = fse3f_n (:,:,jk)
948            fse3w (:,:,jk) = fse3w_n (:,:,jk)
949            fse3uw(:,:,jk) = fse3uw_n(:,:,jk)
950            fse3vw(:,:,jk) = fse3vw_n(:,:,jk)
951         END DO
952         !
953         hu(:,:) = hu_0(:,:) + sshu_n(:,:)            ! now ocean depth (at u- and v-points)
954         hv(:,:) = hv_0(:,:) + sshv_n(:,:)
955         !                                            ! now masked inverse of the ocean depth (at u- and v-points)
956         hur(:,:) = umask(:,:,1) / ( hu(:,:) + 1._wp - umask(:,:,1) )
957         hvr(:,:) = vmask(:,:,1) / ( hv(:,:) + 1._wp - vmask(:,:,1) )
958         !
959      ENDIF
960      !
961      CALL div_cur( kt )                              ! Horizontal divergence & Relative vorticity
962      !
963      z2dt = 2._wp * rdt                              ! set time step size (Euler/Leapfrog)
964      IF( neuler == 0 .AND. kt == nittrc000 )   z2dt = rdt
965
966      !                                           !------------------------------!
967      !                                           !   After Sea Surface Height   !
968      !                                           !------------------------------!
969      zhdiv(:,:) = 0._wp
970      DO jk = 1, jpkm1                                 ! Horizontal divergence of barotropic transports
971        zhdiv(:,:) = zhdiv(:,:) + fse3t(:,:,jk) * hdivn(:,:,jk)
972      END DO
973      !                                                ! Sea surface elevation time stepping
974      ! In forward Euler time stepping case, the same formulation as in the leap-frog case can be used
975      ! because emp_b field is initialized with the vlaues of emp field. Hence, 0.5 * ( emp + emp_b ) = emp
976      z1_rau0 = 0.5 / rau0
977      ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * tmask(:,:,1)
978
979#if defined key_agrif
980      CALL agrif_ssh( kt )
981#endif
982#if defined key_obc
983      IF( Agrif_Root() ) THEN
984         ssha(:,:) = ssha(:,:) * obctmsk(:,:)
985         CALL lbc_lnk( ssha, 'T', 1. )                 ! absolutly compulsory !! (jmm)
986      ENDIF
987#endif
988#if defined key_bdy
989      ssha(:,:) = ssha(:,:) * bdytmask(:,:)
990      CALL lbc_lnk( ssha, 'T', 1. ) 
991#endif
992
993      !                                                ! Sea Surface Height at u-,v- and f-points (vvl case only)
994      IF( lk_vvl ) THEN                                ! (required only in key_vvl case)
995         DO jj = 1, jpjm1
996            DO ji = 1, jpim1      ! NO Vector Opt.
997               sshu_a(ji,jj) = 0.5  * umask(ji,jj,1) / ( e1u(ji  ,jj) * e2u(ji  ,jj) )                   &
998                  &                                  * ( e1t(ji  ,jj) * e2t(ji  ,jj) * ssha(ji  ,jj)     &
999                  &                                    + e1t(ji+1,jj) * e2t(ji+1,jj) * ssha(ji+1,jj) )
1000               sshv_a(ji,jj) = 0.5  * vmask(ji,jj,1) / ( e1v(ji,jj  ) * e2v(ji,jj  ) )                   &
1001                  &                                  * ( e1t(ji,jj  ) * e2t(ji,jj  ) * ssha(ji,jj  )     &
1002                  &                                    + e1t(ji,jj+1) * e2t(ji,jj+1) * ssha(ji,jj+1) )
1003            END DO
1004         END DO
1005         CALL lbc_lnk( sshu_a, 'U', 1. )   ;   CALL lbc_lnk( sshv_a, 'V', 1. )      ! Boundaries conditions
1006      ENDIF
1007     
1008
1009      !                                           !------------------------------!
1010      !                                           !     Now Vertical Velocity    !
1011      !                                           !------------------------------!
1012      z1_2dt = 1.e0 / z2dt
1013      DO jk = jpkm1, 1, -1                             ! integrate from the bottom the hor. divergence
1014         ! - ML - need 3 lines here because replacement of fse3t by its expression yields too long lines otherwise
1015         wn(:,:,jk) = wn(:,:,jk+1) -   fse3t_n(:,:,jk) * hdivn(:,:,jk)        &
1016            &                      - ( fse3t_a(:,:,jk) - fse3t_b(:,:,jk) )    &
1017            &                         * tmask(:,:,jk) * z1_2dt
1018#if defined key_bdy
1019         wn(:,:,jk) = wn(:,:,jk) * bdytmask(:,:)
1020#endif
1021      END DO
1022
1023      !
1024      CALL wrk_dealloc( jpi, jpj, zhdiv )
1025      !
1026      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_ssh')
1027      !
1028   END SUBROUTINE trc_sub_ssh
1029
1030   INTEGER FUNCTION trc_sub_alloc()
1031      !!-------------------------------------------------------------------
1032      !!                    *** ROUTINE trc_sub_alloc ***
1033      !!-------------------------------------------------------------------
1034      USE lib_mpp, ONLY: ctl_warn
1035      INTEGER ::  ierr
1036      !!-------------------------------------------------------------------
1037      !
1038      ALLOCATE( un_temp(jpi,jpj,jpk)        ,  vn_temp(jpi,jpj,jpk)  ,   &
1039         &      wn_temp(jpi,jpj,jpk)        ,  avt_temp(jpi,jpj,jpk) ,   &
1040         &      rhop_temp(jpi,jpj,jpk)      ,  rhop_tm(jpi,jpj,jpk) ,   &
1041         &      sshn_temp(jpi,jpj)          ,  sshb_temp(jpi,jpj) ,      &
1042         &      ssha_temp(jpi,jpj)          ,  sshu_a_temp(jpi,jpj),     &
1043         &      sshu_n_temp(jpi,jpj)        ,  sshu_b_temp(jpi,jpj),     &
1044         &      sshv_n_temp(jpi,jpj)        ,  sshv_b_temp(jpi,jpj),     &
1045         &      sshv_a_temp(jpi,jpj)        ,  sshf_n_temp(jpi,jpj) ,   &
1046         &      e3t_temp(jpi,jpj,jpk)       ,  e3u_temp(jpi,jpj,jpk),    &
1047         &      e3v_temp(jpi,jpj,jpk)       ,  e3w_temp(jpi,jpj,jpk),    &
1048#if defined key_ldfslp
1049         &      wslpi_temp(jpi,jpj,jpk)     ,  wslpj_temp(jpi,jpj,jpk),  &
1050         &      uslp_temp(jpi,jpj,jpk)      ,  vslp_temp(jpi,jpj,jpk),   &
1051#endif
1052#if defined key_trabbl
1053         &      ahu_bbl_temp(jpi,jpj)       ,  ahv_bbl_temp(jpi,jpj),    &
1054         &      utr_bbl_temp(jpi,jpj)       ,  vtr_bbl_temp(jpi,jpj),    &
1055#endif
1056         &      rnf_temp(jpi,jpj)           ,  h_rnf_temp(jpi,jpj) ,     &
1057         &      tsn_temp(jpi,jpj,jpk,2)     ,  emp_b_temp(jpi,jpj),      &
1058         &      emp_temp(jpi,jpj)           ,  fmmflx_temp(jpi,jpj),     &
1059         &      hmld_temp(jpi,jpj)          ,  qsr_temp(jpi,jpj) ,       &
1060         &      fr_i_temp(jpi,jpj)          ,  fr_i_tm(jpi,jpj) ,        &
1061         &      wndm_temp(jpi,jpj)          ,  wndm_tm(jpi,jpj) ,        &
1062# if defined key_zdfddm
1063         &      avs_tm(jpi,jpj,jpk)         ,  avs_temp(jpi,jpj,jpk) ,   &
1064# endif
1065#if defined key_traldf_c3d
1066         &      ahtt_tm(jpi,jpj,jpk)        ,  ahtt_temp(jpi,jpj,jpk),   &
1067         &      ahtu_tm(jpi,jpj,jpk)        ,  ahtu_temp(jpi,jpj,jpk),   &
1068         &      ahtv_tm(jpi,jpj,jpk)        ,  ahtv_temp(jpi,jpj,jpk),   &
1069         &      ahtw_tm(jpi,jpj,jpk)        ,  ahtw_temp(jpi,jpj,jpk),   &
1070#elif defined key_traldf_c2d
1071         &      ahtt_tm(jpi,jpj)            ,  ahtt_temp(jpi,jpj),       &
1072         &      ahtu_tm(jpi,jpj)            ,  ahtu_temp(jpi,jpj),       &
1073         &      ahtv_tm(jpi,jpj)            ,  ahtv_temp(jpi,jpj),       &
1074         &      ahtw_tm(jpi,jpj)            ,  ahtw_temp(jpi,jpj),       &
1075#elif defined key_traldf_c1d
1076         &      ahtt_tm(jpk)                ,  ahtt_temp(jpk),           &
1077         &      ahtu_tm(jpk)                ,  ahtu_temp(jpk),           &
1078         &      ahtv_tm(jpk)                ,  ahtv_temp(jpk),           &
1079         &      ahtw_tm(jpk)                ,  ahtw_temp(jpk),           &
1080#endif
1081#if defined key_traldf_eiv
1082# if defined key_traldf_c3d
1083         &      aeiu_tm(jpi,jpj,jpk)        ,  aeiu_temp(jpi,jpj,jpk),   &
1084         &      aeiv_tm(jpi,jpj,jpk)        ,  aeiv_temp(jpi,jpj,jpk),   &
1085         &      aeiw_tm(jpi,jpj,jpk)        ,  aeiw_temp(jpi,jpj,jpk),   &
1086# elif defined key_traldf_c2d
1087         &      aeiu_tm(jpi,jpj)            ,  aeiu_temp(jpi,jpj),       &
1088         &      aeiv_tm(jpi,jpj)            ,  aeiv_temp(jpi,jpj),       &
1089         &      aeiw_tm(jpi,jpj)            ,  aeiw_temp(jpi,jpj),       &
1090# elif defined key_traldf_c1d
1091         &      aeiu_tm(jpk)                ,  aeiu_temp(jpk),           &
1092         &      aeiv_tm(jpk)                ,  aeiv_temp(jpk),           &
1093         &      aeiw_tm(jpk)                ,  aeiw_temp(jpk),           &
1094# endif
1095# endif
1096         &      hdivn_temp(jpi,jpj,jpk)     ,  hdivb_temp(jpi,jpj,jpk),  &
1097         &      rotn_temp(jpi,jpj,jpk)      ,  rotb_temp(jpi,jpj,jpk),   &
1098         &      hu_temp(jpi,jpj)            ,  hv_temp(jpi,jpj),         &
1099         &      hur_temp(jpi,jpj)           ,  hvr_temp(jpi,jpj),        &
1100         &      un_tm(jpi,jpj,jpk)          ,  vn_tm(jpi,jpj,jpk)  ,     &
1101         &      avt_tm(jpi,jpj,jpk)                                ,     &
1102         &      sshn_tm(jpi,jpj)            ,  sshb_hold(jpi,jpj) ,      &
1103         &      sshu_n_tm(jpi,jpj)          ,  sshu_b_hold(jpi,jpj),     &
1104         &      sshv_n_tm(jpi,jpj)          ,  sshv_b_hold(jpi,jpj),     &
1105         &      tsn_tm(jpi,jpj,jpk,2)       ,                            &
1106         &      emp_tm(jpi,jpj)             ,  fmmflx_tm(jpi,jpj)  ,     &
1107         &      emp_b_hold(jpi,jpj)         ,                            &
1108         &      hmld_tm(jpi,jpj)            ,  qsr_tm(jpi,jpj) ,         &
1109#if defined key_ldfslp
1110         &      wslpi_tm(jpi,jpj,jpk)       ,  wslpj_tm(jpi,jpj,jpk),    &
1111         &      uslp_tm(jpi,jpj,jpk)        ,  vslp_tm(jpi,jpj,jpk),     &
1112#endif
1113#if defined key_trabbl
1114         &      ahu_bbl_tm(jpi,jpj)         ,  ahv_bbl_tm(jpi,jpj),      &
1115         &      utr_bbl_tm(jpi,jpj)         ,  vtr_bbl_tm(jpi,jpj),      &
1116#endif
1117         &      rnf_tm(jpi,jpj)             ,  h_rnf_tm(jpi,jpj) ,       &
1118         &                                    STAT=trc_sub_alloc ) 
1119      IF( trc_sub_alloc /= 0 )   CALL ctl_warn('trc_sub_alloc: failed to allocate arrays')
1120
1121      !
1122   END FUNCTION trc_sub_alloc
1123
1124#else
1125   !!----------------------------------------------------------------------
1126   !!   Default key                                     NO passive tracers
1127   !!----------------------------------------------------------------------
1128CONTAINS
1129   SUBROUTINE trc_sub_stp( kt )        ! Empty routine
1130      WRITE(*,*) 'trc_sub_stp: You should not have seen this print! error?', kt
1131   END SUBROUTINE trc_sub_stp
1132   SUBROUTINE trc_sub_ini        ! Empty routine
1133      WRITE(*,*) 'trc_sub_ini: You should not have seen this print! error?', kt
1134   END SUBROUTINE trc_sub_ini
1135
1136#endif
1137
1138   !!======================================================================
1139END MODULE trcsub
Note: See TracBrowser for help on using the repository browser.