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

source: trunk/NEMOGCM/NEMO/TOP_SRC/trcsub.F90 @ 5836

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

merge the simplification branch onto the trunk, see ticket #1612

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