New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
trcsub.F90 in branches/2011/dev_NEMO_MERGE_2011/NEMOGCM/NEMO/TOP_SRC – NEMO

source: branches/2011/dev_NEMO_MERGE_2011/NEMOGCM/NEMO/TOP_SRC/trcsub.F90 @ 3153

Last change on this file since 3153 was 3153, checked in by cetlod, 13 years ago

Add timing on TOP routines

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