source: branches/2012/dev_r3385_NOCS04_HAMF/NEMOGCM/NEMO/TOP_SRC/trcsub.F90 @ 3403

Last change on this file since 3403 was 3403, checked in by acc, 8 years ago

Branch: dev_r3385_NOCS04_HAMF; #665. Stage 2a of 2012 development: further suppression of emps array in TOP routines (to be checked) and inclusion of new namelist parameter in reference configuration namelists

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