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/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC – NEMO

source: branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/TOP_SRC/trcsub.F90 @ 5105

Last change on this file since 5105 was 5105, checked in by cbricaud, 9 years ago

bug correction
File size: 26.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, ONLY : jpnf90
16   USE in_out_manager, ONLY : jprstlib
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, ONLY : div_cur           ! hor. divergence and curl      (div & cur routines)
28   USE sbcrnf, ONLY: h_rnf, nk_rnf   ! River runoff
29   USE bdy_oce
30#if defined key_agrif
31   USE agrif_opa_update
32   USE agrif_opa_interp
33#endif
34
35   IMPLICIT NONE
36
37   PUBLIC   trc_sub_stp      ! called by trc_stp
38   PUBLIC   trc_sub_ini      ! called by trc_ini to initialize substepping arrays.
39   PUBLIC   trc_sub_reset    ! called by trc_stp to reset physics variables
40   PUBLIC   trc_sub_ssh      ! called by trc_stp to reset physics variables
41
42   !!* Module variables
43   REAL(wp)  :: r1_ndttrc     !    1 /  nn_dttrc
44   REAL(wp)  :: r1_ndttrcp1   !    1 / (nn_dttrc+1)
45
46   !!* Substitution
47#  include "top_substitute.h90"
48   !!----------------------------------------------------------------------
49   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
50   !! $Id: trcstp.F90 2528 2010-12-27 17:33:53Z rblod $
51   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
52   !!----------------------------------------------------------------------
53CONTAINS
54
55   SUBROUTINE trc_sub_stp( kt )
56      !!-------------------------------------------------------------------
57      !!                     ***  ROUTINE trc_stp  ***
58      !!                     
59      !! ** Purpose : Average variables needed for sub-stepping passive tracers
60      !!
61      !! ** Method  : Called every timestep to increment _tm (time mean) variables
62      !!              on TOP steps, calculate averages.
63      !!-------------------------------------------------------------------
64      INTEGER, INTENT( in ) ::  kt        ! ocean time-step index
65      INTEGER               ::  ji,jj,jk  ! dummy loop indices
66      REAL(wp)              ::  z1_ne3t, z1_ne3u, z1_ne3v, z1_ne3w
67      !!-------------------------------------------------------------------
68      !
69      IF( nn_timing == 1 )  CALL timing_start('trc_sub_stp')
70      !
71      IF( kt == nit000 ) THEN
72           IF(lwp) WRITE(numout,*)
73           IF(lwp) WRITE(numout,*) 'trc_sub_stp : substepping of the passive tracers'
74           IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
75           !
76           sshb_hold  (:,:) = sshn  (:,:)
77           emp_b_hold (:,:) = emp_b (:,:)
78           !
79           r1_ndttrc        = 1._wp / REAL( nn_dttrc    , wp ) 
80           r1_ndttrcp1      = 1._wp / REAL( nn_dttrc + 1, wp )
81           !
82      ENDIF 
83
84       IF( MOD( kt , nn_dttrc ) /= 0 ) THEN
85          !
86          un_tm   (:,:,:)        = un_tm   (:,:,:)        + un   (:,:,:)        * fse3u(:,:,:) 
87          vn_tm   (:,:,:)        = vn_tm   (:,:,:)        + vn   (:,:,:)        * fse3v(:,:,:) 
88          tsn_tm  (:,:,:,jp_tem) = tsn_tm  (:,:,:,jp_tem) + tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
89          tsn_tm  (:,:,:,jp_sal) = tsn_tm  (:,:,:,jp_sal) + tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
90          rhop_tm (:,:,:)        = rhop_tm (:,:,:)        + rhop (:,:,:)        * 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(:,:,:)
97          wslpj_tm(:,:,:)        = wslpj_tm(:,:,:)        + wslpj(:,:,:)
98          uslp_tm (:,:,:)        = uslp_tm (:,:,:)        + uslp (:,:,:)
99          vslp_tm (:,:,:)        = vslp_tm (:,:,:)        + vslp (:,:,:)
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          rnf_tm   (:,:)         = rnf_tm   (:,:)         + rnf   (:,:) 
114          h_rnf_tm (:,:)         = h_rnf_tm (:,:)         + h_rnf (:,:) 
115          hmld_tm  (:,:)         = hmld_tm  (:,:)         + hmld  (:,:)
116          fr_i_tm  (:,:)         = fr_i_tm  (:,:)         + fr_i  (:,:)
117          emp_tm   (:,:)         = emp_tm   (:,:)         + emp   (:,:) 
118          fmmflx_tm(:,:)         = fmmflx_tm(:,:)         + fmmflx(:,:)
119          qsr_tm   (:,:)         = qsr_tm   (:,:)         + qsr   (:,:)
120          wndm_tm  (:,:)         = wndm_tm  (:,:)         + wndm  (:,:)
121
122      ELSE                           !  It is time to substep
123         !   1. set temporary arrays to hold physics variables
124         un_temp    (:,:,:)      = un    (:,:,:)
125         vn_temp    (:,:,:)      = vn    (:,:,:)
126         wn_temp    (:,:,:)      = wn    (:,:,:)
127         tsn_temp   (:,:,:,:)    = tsn   (:,:,:,:)
128         rhop_temp  (:,:,:)      = rhop  (:,:,:)   
129         avt_temp   (:,:,:)      = avt   (:,:,:)
130# if defined key_zdfddm
131         avs_temp   (:,:,:)      = avs   (:,:,:)
132# endif
133#if defined key_ldfslp
134         wslpi_temp (:,:,:)      = wslpi (:,:,:)
135         wslpj_temp (:,:,:)      = wslpj (:,:,:)
136         uslp_temp  (:,:,:)      = uslp  (:,:,:)
137         vslp_temp  (:,:,:)      = vslp  (:,:,:)
138#endif
139# if defined key_trabbl
140          IF( nn_bbl_ldf == 1 ) THEN
141             ahu_bbl_temp(:,:)   = ahu_bbl(:,:) 
142             ahv_bbl_temp(:,:)   = ahv_bbl(:,:) 
143          ENDIF
144          IF( nn_bbl_adv == 1 ) THEN
145             utr_bbl_temp(:,:)   = utr_bbl(:,:) 
146             vtr_bbl_temp(:,:)   = vtr_bbl(:,:) 
147          ENDIF
148# endif
149         sshn_temp  (:,:)        = sshn  (:,:)
150         sshb_temp  (:,:)        = sshb  (:,:)
151         ssha_temp  (:,:)        = ssha  (:,:)
152         rnf_temp   (:,:)        = rnf   (:,:)
153         h_rnf_temp (:,:)        = h_rnf (:,:)
154         hmld_temp  (:,:)        = hmld  (:,:)
155         fr_i_temp  (:,:)        = fr_i  (:,:)
156         emp_temp   (:,:)        = emp   (:,:)
157         emp_b_temp (:,:)        = emp_b (:,:)
158         fmmflx_temp(:,:)        = fmmflx(:,:)
159         qsr_temp   (:,:)        = qsr   (:,:)
160         wndm_temp  (:,:)        = wndm  (:,:)
161         !                                    !  Variables reset in trc_sub_ssh
162#if ! defined key_crs
163         rotn_temp  (:,:,:)      = rotn  (:,:,:)
164# endif
165         hdivn_temp (:,:,:)      = hdivn (:,:,:)
166#if ! defined key_crs
167         rotb_temp  (:,:,:)      = rotb  (:,:,:)
168# endif
169         hdivb_temp (:,:,:)      = hdivb (:,:,:)
170         !
171         ! 2. Create averages and reassign variables
172         un_tm    (:,:,:)        = un_tm   (:,:,:)        + un   (:,:,:)        * fse3u(:,:,:) 
173         vn_tm    (:,:,:)        = vn_tm   (:,:,:)        + vn   (:,:,:)        * fse3v(:,:,:) 
174         tsn_tm   (:,:,:,jp_tem) = tsn_tm  (:,:,:,jp_tem) + tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
175         tsn_tm   (:,:,:,jp_sal) = tsn_tm  (:,:,:,jp_sal) + tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
176         rhop_tm (:,:,:)         = rhop_tm (:,:,:)        + rhop (:,:,:)        * fse3t(:,:,:) 
177         avt_tm   (:,:,:)        = avt_tm  (:,:,:)        + avt  (:,:,:)        * fse3w(:,:,:) 
178# if defined key_zdfddm
179         avs_tm   (:,:,:)        = avs_tm  (:,:,:)        + avs  (:,:,:)        * fse3w(:,:,:) 
180# endif
181#if defined key_ldfslp
182         wslpi_tm (:,:,:)        = wslpi_tm(:,:,:)        + wslpi(:,:,:) 
183         wslpj_tm (:,:,:)        = wslpj_tm(:,:,:)        + wslpj(:,:,:) 
184         uslp_tm  (:,:,:)        = uslp_tm (:,:,:)        + uslp (:,:,:) 
185         vslp_tm  (:,:,:)        = vslp_tm (:,:,:)        + vslp (:,:,:)
186#endif
187# if defined key_trabbl
188          IF( nn_bbl_ldf == 1 ) THEN
189             ahu_bbl_tm(:,:)     = ahu_bbl_tm(:,:)        + ahu_bbl(:,:) 
190             ahv_bbl_tm(:,:)     = ahv_bbl_tm(:,:)        + ahv_bbl(:,:) 
191          ENDIF
192          IF( nn_bbl_adv == 1 ) THEN
193             utr_bbl_tm(:,:)     = utr_bbl_tm(:,:)        + utr_bbl(:,:) 
194             vtr_bbl_tm(:,:)     = vtr_bbl_tm(:,:)        + vtr_bbl(:,:) 
195          ENDIF
196# endif
197         sshn_tm  (:,:)          = sshn_tm    (:,:)       + sshn  (:,:) 
198         rnf_tm   (:,:)          = rnf_tm     (:,:)       + rnf   (:,:) 
199         h_rnf_tm (:,:)          = h_rnf_tm   (:,:)       + h_rnf (:,:) 
200         hmld_tm  (:,:)          = hmld_tm    (:,:)       + hmld  (:,:)
201         fr_i_tm  (:,:)          = fr_i_tm    (:,:)       + fr_i  (:,:)
202         emp_tm   (:,:)          = emp_tm     (:,:)       + emp   (:,:) 
203         fmmflx_tm(:,:)          = fmmflx_tm  (:,:)       + fmmflx(:,:)
204         qsr_tm   (:,:)          = qsr_tm     (:,:)       + qsr   (:,:)
205         wndm_tm  (:,:)          = wndm_tm    (:,:)       + wndm  (:,:)
206         !
207         sshn     (:,:)          = sshn_tm    (:,:) * r1_ndttrcp1 
208         sshb     (:,:)          = sshb_hold  (:,:)
209         rnf      (:,:)          = rnf_tm     (:,:) * r1_ndttrcp1 
210         h_rnf    (:,:)          = h_rnf_tm   (:,:) * r1_ndttrcp1 
211         hmld     (:,:)          = hmld_tm    (:,:) * r1_ndttrcp1 
212         !  variables that are initialized after averages
213         emp_b    (:,:) = emp_b_hold (:,:)
214         IF( kt == nittrc000 ) THEN
215            wndm  (:,:)          = wndm_tm    (:,:) * r1_ndttrc 
216            qsr   (:,:)          = qsr_tm     (:,:) * r1_ndttrc 
217            emp   (:,:)          = emp_tm     (:,:) * r1_ndttrc 
218            fmmflx(:,:)          = fmmflx_tm  (:,:) * r1_ndttrc 
219            fr_i  (:,:)          = fr_i_tm    (:,:) * r1_ndttrc
220# if defined key_trabbl
221            IF( nn_bbl_ldf == 1 ) THEN
222               ahu_bbl(:,:)      = ahu_bbl_tm (:,:) * r1_ndttrc 
223               ahv_bbl(:,:)      = ahv_bbl_tm (:,:) * r1_ndttrc 
224            ENDIF
225            IF( nn_bbl_adv == 1 ) THEN
226               utr_bbl(:,:)      = utr_bbl_tm (:,:) * r1_ndttrc 
227               vtr_bbl(:,:)      = vtr_bbl_tm (:,:) * r1_ndttrc 
228            ENDIF
229# endif
230         ELSE
231            wndm  (:,:)          = wndm_tm    (:,:) * r1_ndttrcp1 
232            qsr   (:,:)          = qsr_tm     (:,:) * r1_ndttrcp1 
233            emp   (:,:)          = emp_tm     (:,:) * r1_ndttrcp1 
234            fmmflx(:,:)          = fmmflx_tm  (:,:) * r1_ndttrcp1 
235            fr_i  (:,:)          = fr_i_tm    (:,:) * r1_ndttrcp1 
236# if defined key_trabbl
237            IF( nn_bbl_ldf == 1 ) THEN
238               ahu_bbl(:,:)      = ahu_bbl_tm (:,:) * r1_ndttrcp1 
239               ahv_bbl(:,:)      = ahv_bbl_tm (:,:) * r1_ndttrcp1 
240            ENDIF
241            IF( nn_bbl_adv == 1 ) THEN
242               utr_bbl(:,:)      = utr_bbl_tm (:,:) * r1_ndttrcp1 
243               vtr_bbl(:,:)      = vtr_bbl_tm (:,:) * r1_ndttrcp1 
244            ENDIF
245# endif
246         ENDIF
247         !
248         DO jk = 1, jpk
249            DO jj = 1, jpj
250               DO ji = 1, jpi
251                  z1_ne3t = r1_ndttrcp1  / fse3t(ji,jj,jk)
252                  z1_ne3u = r1_ndttrcp1  / fse3u(ji,jj,jk)
253                  z1_ne3v = r1_ndttrcp1  / fse3v(ji,jj,jk)
254                  z1_ne3w = r1_ndttrcp1  / fse3w(ji,jj,jk)
255                  !
256                  un   (ji,jj,jk)        = un_tm   (ji,jj,jk)        * z1_ne3u
257                  vn   (ji,jj,jk)        = vn_tm   (ji,jj,jk)        * z1_ne3v
258                  tsn  (ji,jj,jk,jp_tem) = tsn_tm  (ji,jj,jk,jp_tem) * z1_ne3t
259                  tsn  (ji,jj,jk,jp_sal) = tsn_tm  (ji,jj,jk,jp_sal) * z1_ne3t
260                  rhop (ji,jj,jk)        = rhop_tm (ji,jj,jk)        * z1_ne3t
261                  avt  (ji,jj,jk)        = avt_tm  (ji,jj,jk)        * z1_ne3w
262# if defined key_zdfddm
263                  avs  (ji,jj,jk)        = avs_tm  (ji,jj,jk)        * z1_ne3w
264# endif
265#if defined key_ldfslp
266                  wslpi(ji,jj,jk)        = wslpi_tm(ji,jj,jk) 
267                  wslpj(ji,jj,jk)        = wslpj_tm(ji,jj,jk)
268                  uslp (ji,jj,jk)        = uslp_tm (ji,jj,jk)
269                  vslp (ji,jj,jk)        = vslp_tm (ji,jj,jk)
270#endif
271               ENDDO
272            ENDDO
273         ENDDO
274         !
275         CALL trc_sub_ssh( kt )         ! after ssh & vertical velocity
276         !
277      ENDIF
278      !
279      IF( nn_timing == 1 )  CALL timing_start('trc_sub_stp')
280      !
281   END SUBROUTINE trc_sub_stp
282
283   SUBROUTINE trc_sub_ini
284      !!-------------------------------------------------------------------
285      !!                     ***  ROUTINE trc_sub_ini  ***
286      !!                     
287      !! ** Purpose : Initialize variables needed for sub-stepping passive tracers
288      !!
289      !! ** Method  :
290      !!              Compute the averages for sub-stepping
291      !!-------------------------------------------------------------------
292      INTEGER ::   ierr
293      !!-------------------------------------------------------------------
294      !
295      IF( nn_timing == 1 )  CALL timing_start('trc_sub_ini')
296      !
297      IF(lwp) WRITE(numout,*)
298      IF(lwp) WRITE(numout,*) 'trc_sub_ini : initial set up of the passive tracers substepping'
299      IF(lwp) WRITE(numout,*) '~~~~~~~'
300
301      ierr =  trc_sub_alloc    ()
302      IF( lk_mpp    )   CALL mpp_sum( ierr )
303      IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'top_sub_alloc : unable to allocate standard ocean arrays' )
304
305      un_tm   (:,:,:)        = un   (:,:,:)        * fse3u(:,:,:) 
306      vn_tm   (:,:,:)        = vn   (:,:,:)        * fse3v(:,:,:) 
307      tsn_tm  (:,:,:,jp_tem) = tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
308      tsn_tm  (:,:,:,jp_sal) = tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
309      rhop_tm (:,:,:)        = rhop (:,:,:)        * fse3t(:,:,:) 
310      avt_tm  (:,:,:)        = avt  (:,:,:)        * fse3w(:,:,:) 
311# if defined key_zdfddm
312      avs_tm  (:,:,:)        = avs  (:,:,:)        * fse3w(:,:,:) 
313# endif
314#if defined key_ldfslp
315      wslpi_tm(:,:,:)        = wslpi(:,:,:)
316      wslpj_tm(:,:,:)        = wslpj(:,:,:)
317      uslp_tm (:,:,:)        = uslp (:,:,:)
318      vslp_tm (:,:,:)        = vslp (:,:,:)
319#endif
320      sshn_tm  (:,:) = sshn  (:,:) 
321      rnf_tm   (:,:) = rnf   (:,:) 
322      h_rnf_tm (:,:) = h_rnf (:,:) 
323      hmld_tm  (:,:) = hmld  (:,:)
324
325      ! Physics variables that are set after initialization:
326      fr_i_tm(:,:) = 0._wp
327      emp_tm (:,:) = 0._wp
328      fmmflx_tm(:,:)  = 0._wp
329      qsr_tm (:,:) = 0._wp
330      wndm_tm(:,:) = 0._wp
331# if defined key_trabbl
332      IF( nn_bbl_ldf == 1 ) THEN
333         ahu_bbl_tm(:,:) = 0._wp
334         ahv_bbl_tm(:,:) = 0._wp
335      ENDIF
336      IF( nn_bbl_adv == 1 ) THEN
337         utr_bbl_tm(:,:) = 0._wp
338         vtr_bbl_tm(:,:) = 0._wp
339      ENDIF
340# endif
341      !
342      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_ini')
343      !
344   END SUBROUTINE trc_sub_ini
345
346   SUBROUTINE trc_sub_reset( kt )
347      !!-------------------------------------------------------------------
348      !!                     ***  ROUTINE trc_sub_reset  ***
349      !!                     
350      !! ** Purpose : Reset physics variables averaged for substepping
351      !!
352      !! ** Method  :
353      !!              Compute the averages for sub-stepping
354      !!-------------------------------------------------------------------
355      INTEGER, INTENT( in ) ::  kt  ! ocean time-step index
356      INTEGER :: jk                 ! dummy loop indices
357      !!-------------------------------------------------------------------
358      !
359      IF( nn_timing == 1 )  CALL timing_start('trc_sub_reset')
360      !
361      !   restore physics variables
362      un    (:,:,:)   =  un_temp    (:,:,:)
363      vn    (:,:,:)   =  vn_temp    (:,:,:)
364      wn    (:,:,:)   =  wn_temp    (:,:,:)
365      tsn   (:,:,:,:) =  tsn_temp   (:,:,:,:)
366      rhop  (:,:,:)   =  rhop_temp  (:,:,:)
367      avt   (:,:,:)   =  avt_temp   (:,:,:)
368# if defined key_zdfddm
369      avs   (:,:,:)   =  avs_temp   (:,:,:)
370# endif
371#if defined key_ldfslp
372      wslpi (:,:,:)   =  wslpi_temp (:,:,:)
373      wslpj (:,:,:)   =  wslpj_temp (:,:,:)
374      uslp  (:,:,:)   =  uslp_temp  (:,:,:)
375      vslp  (:,:,:)   =  vslp_temp  (:,:,:)
376#endif
377      sshn  (:,:)     =  sshn_temp  (:,:)
378      sshb  (:,:)     =  sshb_temp  (:,:)
379      ssha  (:,:)     =  ssha_temp  (:,:)
380      rnf   (:,:)     =  rnf_temp   (:,:)
381      h_rnf (:,:)     =  h_rnf_temp (:,:)
382      !
383      hmld  (:,:)     =  hmld_temp  (:,:)
384      fr_i  (:,:)     =  fr_i_temp  (:,:)
385      emp   (:,:)     =  emp_temp   (:,:)
386      fmmflx(:,:)     =  fmmflx_temp(:,:)
387      emp_b (:,:)     =  emp_b_temp (:,:)
388      qsr   (:,:)     =  qsr_temp   (:,:)
389      wndm  (:,:)     =  wndm_temp  (:,:)
390# if defined key_trabbl
391      IF( nn_bbl_ldf == 1 ) THEN
392         ahu_bbl(:,:) = ahu_bbl_temp(:,:) 
393         ahv_bbl(:,:) = ahv_bbl_temp(:,:) 
394      ENDIF
395      IF( nn_bbl_adv == 1 ) THEN
396         utr_bbl(:,:) = utr_bbl_temp(:,:) 
397         vtr_bbl(:,:) = vtr_bbl_temp(:,:) 
398      ENDIF
399# endif
400      !
401      hdivn (:,:,:)   =  hdivn_temp (:,:,:)
402      hdivb (:,:,:)   =  hdivb_temp (:,:,:)
403#if ! defined key_crs
404      rotn  (:,:,:)   =  rotn_temp  (:,:,:)
405      rotb  (:,:,:)   =  rotb_temp  (:,:,:)
406#endif
407      !                                     
408
409      ! Start new averages
410         un_tm   (:,:,:)        = un   (:,:,:)        * fse3u(:,:,:) 
411         vn_tm   (:,:,:)        = vn   (:,:,:)        * fse3v(:,:,:) 
412         tsn_tm  (:,:,:,jp_tem) = tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
413         tsn_tm  (:,:,:,jp_sal) = tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
414         rhop_tm (:,:,:)        = rhop (:,:,:)        * fse3t(:,:,:) 
415         avt_tm  (:,:,:)        = avt  (:,:,:)        * fse3w(:,:,:) 
416# if defined key_zdfddm
417         avs_tm  (:,:,:)        = avs  (:,:,:)        * fse3w(:,:,:) 
418# endif
419#if defined key_ldfslp
420         wslpi_tm(:,:,:)        = wslpi(:,:,:) 
421         wslpj_tm(:,:,:)        = wslpj(:,:,:)
422         uslp_tm (:,:,:)        = uslp (:,:,:)
423         vslp_tm (:,:,:)        = vslp (:,:,:)
424#endif
425      !
426      sshb_hold  (:,:) = sshn  (:,:)
427      emp_b_hold (:,:) = emp   (:,:)
428      sshn_tm    (:,:) = sshn  (:,:) 
429      rnf_tm     (:,:) = rnf   (:,:) 
430      h_rnf_tm   (:,:) = h_rnf (:,:) 
431      hmld_tm    (:,:) = hmld  (:,:)
432      fr_i_tm    (:,:) = fr_i  (:,:)
433      emp_tm     (:,:) = emp   (:,:)
434      fmmflx_tm  (:,:) = fmmflx(:,:)
435      qsr_tm     (:,:) = qsr   (:,:)
436      wndm_tm    (:,:) = wndm  (:,:)
437# if defined key_trabbl
438      IF( nn_bbl_ldf == 1 ) THEN
439         ahu_bbl_tm(:,:) = ahu_bbl(:,:) 
440         ahv_bbl_tm(:,:) = ahv_bbl(:,:) 
441      ENDIF
442      IF( nn_bbl_adv == 1 ) THEN
443         utr_bbl_tm(:,:) = utr_bbl(:,:) 
444         vtr_bbl_tm(:,:) = vtr_bbl(:,:) 
445      ENDIF
446# endif
447      !
448      !
449      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_reset')
450      !
451   END SUBROUTINE trc_sub_reset
452
453
454   SUBROUTINE trc_sub_ssh( kt ) 
455      !!----------------------------------------------------------------------
456      !!                ***  ROUTINE trc_sub_ssh  ***
457      !!                   
458      !! ** Purpose :   compute the after ssh (ssha), the now vertical velocity
459      !!              and update the now vertical coordinate (lk_vvl=T).
460      !!
461      !! ** Method  : - Using the incompressibility hypothesis, the vertical
462      !!      velocity is computed by integrating the horizontal divergence 
463      !!      from the bottom to the surface minus the scale factor evolution.
464      !!        The boundary conditions are w=0 at the bottom (no flux) and.
465      !!
466      !! ** action  :   ssha    : after sea surface height
467      !!                wn      : now vertical velocity
468      !!                sshu_a, sshv_a, sshf_a  : after sea surface height (lk_vvl=T)
469      !!
470      !! Reference  : Leclair, M., and G. Madec, 2009, Ocean Modelling.
471      !!----------------------------------------------------------------------
472      !
473      INTEGER, INTENT(in) ::   kt   ! time step
474      !
475      INTEGER  ::   ji, jj, jk   ! dummy loop indices
476      REAL(wp) ::   zcoefu, zcoefv, zcoeff, z2dt, z1_2dt, z1_rau0   ! local scalars
477      REAL(wp), POINTER, DIMENSION(:,:) :: zhdiv
478      !!---------------------------------------------------------------------
479      !
480      IF( nn_timing == 1 )  CALL timing_start('trc_sub_ssh')
481      !
482      ! Allocate temporary workspace
483      CALL wrk_alloc( jpi, jpj, zhdiv )
484
485      IF( kt == nittrc000 ) THEN
486         !
487         IF(lwp) WRITE(numout,*)
488         IF(lwp) WRITE(numout,*) 'trc_sub_ssh : after sea surface height and now vertical velocity '
489         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
490         !
491         wn(:,:,jpk) = 0._wp                  ! bottom boundary condition: w=0 (set once for all)
492         !
493      ENDIF
494      !
495      CALL div_cur( kt )                              ! Horizontal divergence & Relative vorticity
496      !
497      z2dt = 2._wp * rdt                              ! set time step size (Euler/Leapfrog)
498      IF( neuler == 0 .AND. kt == nittrc000 )   z2dt = rdt
499
500      !                                           !------------------------------!
501      !                                           !   After Sea Surface Height   !
502      !                                           !------------------------------!
503      zhdiv(:,:) = 0._wp
504      DO jk = 1, jpkm1                                 ! Horizontal divergence of barotropic transports
505        zhdiv(:,:) = zhdiv(:,:) + fse3t(:,:,jk) * hdivn(:,:,jk)
506      END DO
507      !                                                ! Sea surface elevation time stepping
508      ! In forward Euler time stepping case, the same formulation as in the leap-frog case can be used
509      ! because emp_b field is initialized with the vlaues of emp field. Hence, 0.5 * ( emp + emp_b ) = emp
510      z1_rau0 = 0.5 / rau0
511      ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * tmask(:,:,1)
512#if ! defined key_dynspg_ts
513      ! These lines are not necessary with time splitting since
514      ! boundary condition on sea level is set during ts loop
515#if defined key_agrif
516      CALL agrif_ssh( kt )
517#endif
518#if defined key_bdy
519      ssha(:,:) = ssha(:,:) * bdytmask(:,:)
520      CALL lbc_lnk( ssha, 'T', 1. ) 
521#endif
522#endif
523
524
525      !                                           !------------------------------!
526      !                                           !     Now Vertical Velocity    !
527      !                                           !------------------------------!
528      z1_2dt = 1.e0 / z2dt
529      DO jk = jpkm1, 1, -1                             ! integrate from the bottom the hor. divergence
530         ! - ML - need 3 lines here because replacement of fse3t by its expression yields too long lines otherwise
531         wn(:,:,jk) = wn(:,:,jk+1) -   fse3t_n(:,:,jk) * hdivn(:,:,jk)        &
532            &                      - ( fse3t_a(:,:,jk) - fse3t_b(:,:,jk) )    &
533            &                         * tmask(:,:,jk) * z1_2dt
534#if defined key_bdy
535         wn(:,:,jk) = wn(:,:,jk) * bdytmask(:,:)
536#endif
537      END DO
538
539      !
540      CALL wrk_dealloc( jpi, jpj, zhdiv )
541      !
542      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_ssh')
543      !
544   END SUBROUTINE trc_sub_ssh
545
546   INTEGER FUNCTION trc_sub_alloc()
547      !!-------------------------------------------------------------------
548      !!                    *** ROUTINE trc_sub_alloc ***
549      !!-------------------------------------------------------------------
550      USE lib_mpp, ONLY: ctl_warn
551      INTEGER ::  ierr
552      !!-------------------------------------------------------------------
553      !
554      ALLOCATE( un_temp(jpi,jpj,jpk)        ,  vn_temp(jpi,jpj,jpk)  ,   &
555         &      wn_temp(jpi,jpj,jpk)        ,  avt_temp(jpi,jpj,jpk) ,   &
556         &      rhop_temp(jpi,jpj,jpk)      ,  rhop_tm(jpi,jpj,jpk) ,   &
557         &      sshn_temp(jpi,jpj)          ,  sshb_temp(jpi,jpj) ,      &
558         &      ssha_temp(jpi,jpj)          ,                           &
559#if defined key_ldfslp
560         &      wslpi_temp(jpi,jpj,jpk)     ,  wslpj_temp(jpi,jpj,jpk),  &
561         &      uslp_temp(jpi,jpj,jpk)      ,  vslp_temp(jpi,jpj,jpk),   &
562#endif
563#if defined key_trabbl
564         &      ahu_bbl_temp(jpi,jpj)       ,  ahv_bbl_temp(jpi,jpj),    &
565         &      utr_bbl_temp(jpi,jpj)       ,  vtr_bbl_temp(jpi,jpj),    &
566#endif
567         &      rnf_temp(jpi,jpj)           ,  h_rnf_temp(jpi,jpj) ,     &
568         &      tsn_temp(jpi,jpj,jpk,2)     ,  emp_b_temp(jpi,jpj),      &
569         &      emp_temp(jpi,jpj)           ,  fmmflx_temp(jpi,jpj),     &
570         &      hmld_temp(jpi,jpj)          ,  qsr_temp(jpi,jpj) ,       &
571         &      fr_i_temp(jpi,jpj)          ,  fr_i_tm(jpi,jpj) ,        &
572         &      wndm_temp(jpi,jpj)          ,  wndm_tm(jpi,jpj) ,        &
573# if defined key_zdfddm
574         &      avs_tm(jpi,jpj,jpk)         ,  avs_temp(jpi,jpj,jpk) ,   &
575# endif
576         &      hdivn_temp(jpi,jpj,jpk)     ,  hdivb_temp(jpi,jpj,jpk),  &
577         &      rotn_temp(jpi,jpj,jpk)      ,  rotb_temp(jpi,jpj,jpk),   &
578         &      un_tm(jpi,jpj,jpk)          ,  vn_tm(jpi,jpj,jpk)  ,     &
579         &      avt_tm(jpi,jpj,jpk)                                ,     &
580         &      sshn_tm(jpi,jpj)            ,  sshb_hold(jpi,jpj) ,      &
581         &      tsn_tm(jpi,jpj,jpk,2)       ,                            &
582         &      emp_tm(jpi,jpj)             ,  fmmflx_tm(jpi,jpj)  ,     &
583         &      emp_b_hold(jpi,jpj)         ,                            &
584         &      hmld_tm(jpi,jpj)            ,  qsr_tm(jpi,jpj) ,         &
585#if defined key_ldfslp
586         &      wslpi_tm(jpi,jpj,jpk)       ,  wslpj_tm(jpi,jpj,jpk),    &
587         &      uslp_tm(jpi,jpj,jpk)        ,  vslp_tm(jpi,jpj,jpk),     &
588#endif
589#if defined key_trabbl
590         &      ahu_bbl_tm(jpi,jpj)         ,  ahv_bbl_tm(jpi,jpj),      &
591         &      utr_bbl_tm(jpi,jpj)         ,  vtr_bbl_tm(jpi,jpj),      &
592#endif
593         &      rnf_tm(jpi,jpj)             ,  h_rnf_tm(jpi,jpj) ,       &
594         &                                    STAT=trc_sub_alloc ) 
595      IF( trc_sub_alloc /= 0 )   CALL ctl_warn('trc_sub_alloc: failed to allocate arrays')
596
597      !
598   END FUNCTION trc_sub_alloc
599
600#else
601   !!----------------------------------------------------------------------
602   !!   Default key                                     NO passive tracers
603   !!----------------------------------------------------------------------
604CONTAINS
605   SUBROUTINE trc_sub_stp( kt )        ! Empty routine
606      WRITE(*,*) 'trc_sub_stp: You should not have seen this print! error?', kt
607   END SUBROUTINE trc_sub_stp
608   SUBROUTINE trc_sub_ini        ! Empty routine
609      WRITE(*,*) 'trc_sub_ini: You should not have seen this print! error?', kt
610   END SUBROUTINE trc_sub_ini
611
612#endif
613
614   !!======================================================================
615END MODULE trcsub
Note: See TracBrowser for help on using the repository browser.