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

source: branches/2015/dev_r5721_CNRS9_NOC3_LDF/NEMOGCM/NEMO/TOP_SRC/trcsub.F90 @ 5758

Last change on this file since 5758 was 5758, checked in by gm, 9 years ago

#1593: LDF-ADV, step II.1: phasing the improvements/simplifications of diffusive trend (see wiki)
  • Property svn:keywords set to Id
File size: 27.3 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 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   !!* Substitution
51#  include "top_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         rotn_temp  (:,:,:)      = rotn  (:,:,:)
165         hdivn_temp (:,:,:)      = hdivn (:,:,:)
166         rotb_temp  (:,:,:)      = rotb  (:,:,:)
167         hdivb_temp (:,:,:)      = hdivb (:,:,:)
168         !
169         ! 2. Create averages and reassign variables
170         un_tm    (:,:,:)        = un_tm   (:,:,:)        + un   (:,:,:)        * fse3u(:,:,:) 
171         vn_tm    (:,:,:)        = vn_tm   (:,:,:)        + vn   (:,:,:)        * fse3v(:,:,:) 
172         tsn_tm   (:,:,:,jp_tem) = tsn_tm  (:,:,:,jp_tem) + tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
173         tsn_tm   (:,:,:,jp_sal) = tsn_tm  (:,:,:,jp_sal) + tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
174         rhop_tm (:,:,:)         = rhop_tm (:,:,:)        + rhop (:,:,:)        * fse3t(:,:,:) 
175         avt_tm   (:,:,:)        = avt_tm  (:,:,:)        + avt  (:,:,:)        * fse3w(:,:,:) 
176# if defined key_zdfddm
177         avs_tm   (:,:,:)        = avs_tm  (:,:,:)        + avs  (:,:,:)        * fse3w(:,:,:) 
178# endif
179         IF( l_ldfslp ) THEN
180            uslp_tm  (:,:,:)     = uslp_tm (:,:,:)        + uslp (:,:,:) 
181            vslp_tm  (:,:,:)     = vslp_tm (:,:,:)        + vslp (:,:,:)
182            wslpi_tm (:,:,:)     = wslpi_tm(:,:,:)        + wslpi(:,:,:) 
183            wslpj_tm (:,:,:)     = wslpj_tm(:,:,:)        + wslpj(:,:,:) 
184         ENDIF
185# if defined key_trabbl
186          IF( nn_bbl_ldf == 1 ) THEN
187             ahu_bbl_tm(:,:)     = ahu_bbl_tm(:,:)        + ahu_bbl(:,:) 
188             ahv_bbl_tm(:,:)     = ahv_bbl_tm(:,:)        + ahv_bbl(:,:) 
189          ENDIF
190          IF( nn_bbl_adv == 1 ) THEN
191             utr_bbl_tm(:,:)     = utr_bbl_tm(:,:)        + utr_bbl(:,:) 
192             vtr_bbl_tm(:,:)     = vtr_bbl_tm(:,:)        + vtr_bbl(:,:) 
193          ENDIF
194# endif
195         sshn_tm  (:,:)          = sshn_tm    (:,:)       + sshn  (:,:) 
196         rnf_tm   (:,:)          = rnf_tm     (:,:)       + rnf   (:,:) 
197         h_rnf_tm (:,:)          = h_rnf_tm   (:,:)       + h_rnf (:,:) 
198         hmld_tm  (:,:)          = hmld_tm    (:,:)       + hmld  (:,:)
199         fr_i_tm  (:,:)          = fr_i_tm    (:,:)       + fr_i  (:,:)
200         emp_tm   (:,:)          = emp_tm     (:,:)       + emp   (:,:) 
201         fmmflx_tm(:,:)          = fmmflx_tm  (:,:)       + fmmflx(:,:)
202         qsr_tm   (:,:)          = qsr_tm     (:,:)       + qsr   (:,:)
203         wndm_tm  (:,:)          = wndm_tm    (:,:)       + wndm  (:,:)
204         !
205         sshn     (:,:)          = sshn_tm    (:,:) * r1_ndttrcp1 
206         sshb     (:,:)          = sshb_hold  (:,:)
207         rnf      (:,:)          = rnf_tm     (:,:) * r1_ndttrcp1 
208         h_rnf    (:,:)          = h_rnf_tm   (:,:) * r1_ndttrcp1 
209         hmld     (:,:)          = hmld_tm    (:,:) * r1_ndttrcp1 
210         !  variables that are initialized after averages
211         emp_b    (:,:) = emp_b_hold (:,:)
212         IF( kt == nittrc000 ) THEN
213            wndm  (:,:)          = wndm_tm    (:,:) * r1_ndttrc 
214            qsr   (:,:)          = qsr_tm     (:,:) * r1_ndttrc 
215            emp   (:,:)          = emp_tm     (:,:) * r1_ndttrc 
216            fmmflx(:,:)          = fmmflx_tm  (:,:) * r1_ndttrc 
217            fr_i  (:,:)          = fr_i_tm    (:,:) * r1_ndttrc
218# if defined key_trabbl
219            IF( nn_bbl_ldf == 1 ) THEN
220               ahu_bbl(:,:)      = ahu_bbl_tm (:,:) * r1_ndttrc 
221               ahv_bbl(:,:)      = ahv_bbl_tm (:,:) * r1_ndttrc 
222            ENDIF
223            IF( nn_bbl_adv == 1 ) THEN
224               utr_bbl(:,:)      = utr_bbl_tm (:,:) * r1_ndttrc 
225               vtr_bbl(:,:)      = vtr_bbl_tm (:,:) * r1_ndttrc 
226            ENDIF
227# endif
228         ELSE
229            wndm  (:,:)          = wndm_tm    (:,:) * r1_ndttrcp1 
230            qsr   (:,:)          = qsr_tm     (:,:) * r1_ndttrcp1 
231            emp   (:,:)          = emp_tm     (:,:) * r1_ndttrcp1 
232            fmmflx(:,:)          = fmmflx_tm  (:,:) * r1_ndttrcp1 
233            fr_i  (:,:)          = fr_i_tm    (:,:) * r1_ndttrcp1 
234# if defined key_trabbl
235            IF( nn_bbl_ldf == 1 ) THEN
236               ahu_bbl(:,:)      = ahu_bbl_tm (:,:) * r1_ndttrcp1 
237               ahv_bbl(:,:)      = ahv_bbl_tm (:,:) * r1_ndttrcp1 
238            ENDIF
239            IF( nn_bbl_adv == 1 ) THEN
240               utr_bbl(:,:)      = utr_bbl_tm (:,:) * r1_ndttrcp1 
241               vtr_bbl(:,:)      = vtr_bbl_tm (:,:) * r1_ndttrcp1 
242            ENDIF
243# endif
244         ENDIF
245         !
246         DO jk = 1, jpk
247            DO jj = 1, jpj
248               DO ji = 1, jpi
249                  z1_ne3t = r1_ndttrcp1  / fse3t(ji,jj,jk)
250                  z1_ne3u = r1_ndttrcp1  / fse3u(ji,jj,jk)
251                  z1_ne3v = r1_ndttrcp1  / fse3v(ji,jj,jk)
252                  z1_ne3w = r1_ndttrcp1  / fse3w(ji,jj,jk)
253                  !
254                  un   (ji,jj,jk)        = un_tm   (ji,jj,jk)        * z1_ne3u
255                  vn   (ji,jj,jk)        = vn_tm   (ji,jj,jk)        * z1_ne3v
256                  tsn  (ji,jj,jk,jp_tem) = tsn_tm  (ji,jj,jk,jp_tem) * z1_ne3t
257                  tsn  (ji,jj,jk,jp_sal) = tsn_tm  (ji,jj,jk,jp_sal) * z1_ne3t
258                  rhop (ji,jj,jk)        = rhop_tm (ji,jj,jk)        * z1_ne3t
259!!gm : BUG? ==>> for avt & avs I don't understand the division by e3w
260                  avt  (ji,jj,jk)        = avt_tm  (ji,jj,jk)        * z1_ne3w
261# if defined key_zdfddm
262                  avs  (ji,jj,jk)        = avs_tm  (ji,jj,jk)        * z1_ne3w
263# endif
264               END DO
265            END DO
266         END DO
267         IF( l_ldfslp ) THEN
268            wslpi(:,:,:)        = wslpi_tm(:,:,:) 
269            wslpj(:,:,:)        = wslpj_tm(:,:,:)
270            uslp (:,:,:)        = uslp_tm (:,:,:)
271            vslp (:,:,:)        = vslp_tm (:,:,:)
272         ENDIF
273         !
274         CALL trc_sub_ssh( kt )         ! after ssh & vertical velocity
275         !
276      ENDIF
277      !
278      IF( nn_timing == 1 )  CALL timing_start('trc_sub_stp')
279      !
280   END SUBROUTINE trc_sub_stp
281
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!!gm : BUG? ==>> for avt & avs I don't understand the division by e3w
311      avt_tm  (:,:,:)        = avt  (:,:,:)        * fse3w(:,:,:) 
312# if defined key_zdfddm
313      avs_tm  (:,:,:)        = avs  (:,:,:)        * fse3w(:,:,:) 
314# endif
315      IF( l_ldfslp ) THEN
316         wslpi_tm(:,:,:)     = wslpi(:,:,:)
317         wslpj_tm(:,:,:)     = wslpj(:,:,:)
318         uslp_tm (:,:,:)     = uslp (:,:,:)
319         vslp_tm (:,:,:)     = vslp (:,:,:)
320      ENDIF
321      sshn_tm  (:,:) = sshn  (:,:) 
322      rnf_tm   (:,:) = rnf   (:,:) 
323      h_rnf_tm (:,:) = h_rnf (:,:) 
324      hmld_tm  (:,:) = hmld  (:,:)
325
326      ! Physics variables that are set after initialization:
327      fr_i_tm(:,:) = 0._wp
328      emp_tm (:,:) = 0._wp
329      fmmflx_tm(:,:)  = 0._wp
330      qsr_tm (:,:) = 0._wp
331      wndm_tm(:,:) = 0._wp
332# if defined key_trabbl
333      IF( nn_bbl_ldf == 1 ) THEN
334         ahu_bbl_tm(:,:) = 0._wp
335         ahv_bbl_tm(:,:) = 0._wp
336      ENDIF
337      IF( nn_bbl_adv == 1 ) THEN
338         utr_bbl_tm(:,:) = 0._wp
339         vtr_bbl_tm(:,:) = 0._wp
340      ENDIF
341# endif
342      !
343      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_ini')
344      !
345   END SUBROUTINE trc_sub_ini
346
347   SUBROUTINE trc_sub_reset( kt )
348      !!-------------------------------------------------------------------
349      !!                     ***  ROUTINE trc_sub_reset  ***
350      !!                     
351      !! ** Purpose : Reset physics variables averaged for substepping
352      !!
353      !! ** Method  :
354      !!              Compute the averages for sub-stepping
355      !!-------------------------------------------------------------------
356      INTEGER, INTENT( in ) ::  kt  ! ocean time-step index
357      INTEGER :: jk                 ! dummy loop indices
358      !!-------------------------------------------------------------------
359      !
360      IF( nn_timing == 1 )  CALL timing_start('trc_sub_reset')
361      !
362      !   restore physics variables
363      un    (:,:,:)   =  un_temp    (:,:,:)
364      vn    (:,:,:)   =  vn_temp    (:,:,:)
365      wn    (:,:,:)   =  wn_temp    (:,:,:)
366      tsn   (:,:,:,:) =  tsn_temp   (:,:,:,:)
367      rhop  (:,:,:)   =  rhop_temp  (:,:,:)
368      avt   (:,:,:)   =  avt_temp   (:,:,:)
369# if defined key_zdfddm
370      avs   (:,:,:)   =  avs_temp   (:,:,:)
371# endif
372      IF( l_ldfslp ) THEN
373         wslpi (:,:,:)=  wslpi_temp (:,:,:)
374         wslpj (:,:,:)=  wslpj_temp (:,:,:)
375         uslp  (:,:,:)=  uslp_temp  (:,:,:)
376         vslp  (:,:,:)=  vslp_temp  (:,:,:)
377      ENDIF
378      sshn  (:,:)     =  sshn_temp  (:,:)
379      sshb  (:,:)     =  sshb_temp  (:,:)
380      ssha  (:,:)     =  ssha_temp  (:,:)
381      rnf   (:,:)     =  rnf_temp   (:,:)
382      h_rnf (:,:)     =  h_rnf_temp (:,:)
383      !
384      hmld  (:,:)     =  hmld_temp  (:,:)
385      fr_i  (:,:)     =  fr_i_temp  (:,:)
386      emp   (:,:)     =  emp_temp   (:,:)
387      fmmflx(:,:)     =  fmmflx_temp(:,:)
388      emp_b (:,:)     =  emp_b_temp (:,:)
389      qsr   (:,:)     =  qsr_temp   (:,:)
390      wndm  (:,:)     =  wndm_temp  (:,:)
391# if defined key_trabbl
392      IF( nn_bbl_ldf == 1 ) THEN
393         ahu_bbl(:,:) = ahu_bbl_temp(:,:) 
394         ahv_bbl(:,:) = ahv_bbl_temp(:,:) 
395      ENDIF
396      IF( nn_bbl_adv == 1 ) THEN
397         utr_bbl(:,:) = utr_bbl_temp(:,:) 
398         vtr_bbl(:,:) = vtr_bbl_temp(:,:) 
399      ENDIF
400# endif
401      !
402      hdivn (:,:,:)   =  hdivn_temp (:,:,:)
403      rotn  (:,:,:)   =  rotn_temp  (:,:,:)
404      hdivb (:,:,:)   =  hdivb_temp (:,:,:)
405      rotb  (:,:,:)   =  rotb_temp  (:,:,:)
406      !                                     
407
408      ! Start new averages
409         un_tm   (:,:,:)        = un   (:,:,:)        * fse3u(:,:,:) 
410         vn_tm   (:,:,:)        = vn   (:,:,:)        * fse3v(:,:,:) 
411         tsn_tm  (:,:,:,jp_tem) = tsn  (:,:,:,jp_tem) * fse3t(:,:,:) 
412         tsn_tm  (:,:,:,jp_sal) = tsn  (:,:,:,jp_sal) * fse3t(:,:,:) 
413         rhop_tm (:,:,:)        = rhop (:,:,:)        * fse3t(:,:,:) 
414         avt_tm  (:,:,:)        = avt  (:,:,:)        * fse3w(:,:,:) 
415# if defined key_zdfddm
416         avs_tm  (:,:,:)        = avs  (:,:,:)        * fse3w(:,:,:) 
417# endif
418      IF( l_ldfslp ) THEN
419         uslp_tm (:,:,:)        = uslp (:,:,:)
420         vslp_tm (:,:,:)        = vslp (:,:,:)
421         wslpi_tm(:,:,:)        = wslpi(:,:,:) 
422         wslpj_tm(:,:,:)        = wslpj(:,:,:)
423      ENDIF
424      !
425      sshb_hold  (:,:) = sshn  (:,:)
426      emp_b_hold (:,:) = emp   (:,:)
427      sshn_tm    (:,:) = sshn  (:,:) 
428      rnf_tm     (:,:) = rnf   (:,:) 
429      h_rnf_tm   (:,:) = h_rnf (:,:) 
430      hmld_tm    (:,:) = hmld  (:,:)
431      fr_i_tm    (:,:) = fr_i  (:,:)
432      emp_tm     (:,:) = emp   (:,:)
433      fmmflx_tm  (:,:) = fmmflx(:,:)
434      qsr_tm     (:,:) = qsr   (:,:)
435      wndm_tm    (:,:) = wndm  (:,:)
436# if defined key_trabbl
437      IF( nn_bbl_ldf == 1 ) THEN
438         ahu_bbl_tm(:,:) = ahu_bbl(:,:) 
439         ahv_bbl_tm(:,:) = ahv_bbl(:,:) 
440      ENDIF
441      IF( nn_bbl_adv == 1 ) THEN
442         utr_bbl_tm(:,:) = utr_bbl(:,:) 
443         vtr_bbl_tm(:,:) = vtr_bbl(:,:) 
444      ENDIF
445# endif
446      !
447      !
448      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_reset')
449      !
450   END SUBROUTINE trc_sub_reset
451
452
453   SUBROUTINE trc_sub_ssh( kt ) 
454      !!----------------------------------------------------------------------
455      !!                ***  ROUTINE trc_sub_ssh  ***
456      !!                   
457      !! ** Purpose :   compute the after ssh (ssha), the now vertical velocity
458      !!              and update the now vertical coordinate (lk_vvl=T).
459      !!
460      !! ** Method  : - Using the incompressibility hypothesis, the vertical
461      !!      velocity is computed by integrating the horizontal divergence 
462      !!      from the bottom to the surface minus the scale factor evolution.
463      !!        The boundary conditions are w=0 at the bottom (no flux) and.
464      !!
465      !! ** action  :   ssha    : after sea surface height
466      !!                wn      : now vertical velocity
467      !!                sshu_a, sshv_a, sshf_a  : after sea surface height (lk_vvl=T)
468      !!
469      !! Reference  : Leclair, M., and G. Madec, 2009, Ocean Modelling.
470      !!----------------------------------------------------------------------
471      !
472      INTEGER, INTENT(in) ::   kt   ! time step
473      !
474      INTEGER  ::   ji, jj, jk   ! dummy loop indices
475      REAL(wp) ::   zcoefu, zcoefv, zcoeff, z2dt, z1_2dt, z1_rau0   ! local scalars
476      REAL(wp), POINTER, DIMENSION(:,:) :: zhdiv
477      !!---------------------------------------------------------------------
478      !
479      IF( nn_timing == 1 )  CALL timing_start('trc_sub_ssh')
480      !
481      ! Allocate temporary workspace
482      CALL wrk_alloc( jpi, jpj, zhdiv )
483
484      IF( kt == nittrc000 ) THEN
485         !
486         IF(lwp) WRITE(numout,*)
487         IF(lwp) WRITE(numout,*) 'trc_sub_ssh : after sea surface height and now vertical velocity '
488         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '
489         !
490         wn(:,:,jpk) = 0._wp                  ! bottom boundary condition: w=0 (set once for all)
491         !
492      ENDIF
493      !
494      CALL div_cur( kt )                              ! Horizontal divergence & Relative vorticity
495      !
496      z2dt = 2._wp * rdt                              ! set time step size (Euler/Leapfrog)
497      IF( neuler == 0 .AND. kt == nittrc000 )   z2dt = rdt
498
499      !                                           !------------------------------!
500      !                                           !   After Sea Surface Height   !
501      !                                           !------------------------------!
502      zhdiv(:,:) = 0._wp
503      DO jk = 1, jpkm1                                 ! Horizontal divergence of barotropic transports
504        zhdiv(:,:) = zhdiv(:,:) + fse3t(:,:,jk) * hdivn(:,:,jk)
505      END DO
506      !                                                ! Sea surface elevation time stepping
507      ! In forward Euler time stepping case, the same formulation as in the leap-frog case can be used
508      ! because emp_b field is initialized with the vlaues of emp field. Hence, 0.5 * ( emp + emp_b ) = emp
509      z1_rau0 = 0.5 / rau0
510      ssha(:,:) = (  sshb(:,:) - z2dt * ( z1_rau0 * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * tmask(:,:,1)
511#if ! defined key_dynspg_ts
512      ! These lines are not necessary with time splitting since
513      ! boundary condition on sea level is set during ts loop
514#if defined key_agrif
515      CALL agrif_ssh( kt )
516#endif
517#if defined key_bdy
518      ssha(:,:) = ssha(:,:) * bdytmask(:,:)
519      CALL lbc_lnk( ssha, 'T', 1. ) 
520#endif
521#endif
522
523
524      !                                           !------------------------------!
525      !                                           !     Now Vertical Velocity    !
526      !                                           !------------------------------!
527      z1_2dt = 1.e0 / z2dt
528      DO jk = jpkm1, 1, -1                             ! integrate from the bottom the hor. divergence
529         ! - ML - need 3 lines here because replacement of fse3t by its expression yields too long lines otherwise
530         wn(:,:,jk) = wn(:,:,jk+1) -   fse3t_n(:,:,jk) * hdivn(:,:,jk)        &
531            &                      - ( fse3t_a(:,:,jk) - fse3t_b(:,:,jk) )    &
532            &                         * tmask(:,:,jk) * z1_2dt
533#if defined key_bdy
534         wn(:,:,jk) = wn(:,:,jk) * bdytmask(:,:)
535#endif
536      END DO
537
538      !
539      CALL wrk_dealloc( jpi, jpj, zhdiv )
540      !
541      IF( nn_timing == 1 )  CALL timing_stop('trc_sub_ssh')
542      !
543   END SUBROUTINE trc_sub_ssh
544
545   INTEGER FUNCTION trc_sub_alloc()
546      !!-------------------------------------------------------------------
547      !!                    *** ROUTINE trc_sub_alloc ***
548      !!-------------------------------------------------------------------
549      USE lib_mpp, ONLY: ctl_warn
550      INTEGER ::  ierr
551      !!-------------------------------------------------------------------
552      !
553      ALLOCATE( un_temp(jpi,jpj,jpk)        ,  vn_temp(jpi,jpj,jpk)  ,   &
554         &      wn_temp(jpi,jpj,jpk)        ,  avt_temp(jpi,jpj,jpk) ,   &
555         &      rhop_temp(jpi,jpj,jpk)      ,  rhop_tm(jpi,jpj,jpk) ,   &
556         &      sshn_temp(jpi,jpj)          ,  sshb_temp(jpi,jpj) ,      &
557         &      ssha_temp(jpi,jpj)          ,                           &
558#if defined key_trabbl
559         &      ahu_bbl_temp(jpi,jpj)       ,  ahv_bbl_temp(jpi,jpj),    &
560         &      utr_bbl_temp(jpi,jpj)       ,  vtr_bbl_temp(jpi,jpj),    &
561#endif
562         &      rnf_temp(jpi,jpj)           ,  h_rnf_temp(jpi,jpj) ,     &
563         &      tsn_temp(jpi,jpj,jpk,2)     ,  emp_b_temp(jpi,jpj),      &
564         &      emp_temp(jpi,jpj)           ,  fmmflx_temp(jpi,jpj),     &
565         &      hmld_temp(jpi,jpj)          ,  qsr_temp(jpi,jpj) ,       &
566         &      fr_i_temp(jpi,jpj)          ,  fr_i_tm(jpi,jpj) ,        &
567         &      wndm_temp(jpi,jpj)          ,  wndm_tm(jpi,jpj) ,        &
568# if defined key_zdfddm
569         &      avs_tm(jpi,jpj,jpk)         ,  avs_temp(jpi,jpj,jpk) ,   &
570# endif
571         &      hdivn_temp(jpi,jpj,jpk)     ,  hdivb_temp(jpi,jpj,jpk),  &
572         &      rotn_temp(jpi,jpj,jpk)      ,  rotb_temp(jpi,jpj,jpk),   &
573         &      un_tm(jpi,jpj,jpk)          ,  vn_tm(jpi,jpj,jpk)  ,     &
574         &      avt_tm(jpi,jpj,jpk)                                ,     &
575         &      sshn_tm(jpi,jpj)            ,  sshb_hold(jpi,jpj) ,      &
576         &      tsn_tm(jpi,jpj,jpk,2)       ,                            &
577         &      emp_tm(jpi,jpj)             ,  fmmflx_tm(jpi,jpj)  ,     &
578         &      emp_b_hold(jpi,jpj)         ,                            &
579         &      hmld_tm(jpi,jpj)            ,  qsr_tm(jpi,jpj) ,         &
580#if defined key_trabbl
581         &      ahu_bbl_tm(jpi,jpj)         ,  ahv_bbl_tm(jpi,jpj),      &
582         &      utr_bbl_tm(jpi,jpj)         ,  vtr_bbl_tm(jpi,jpj),      &
583#endif
584         &      rnf_tm(jpi,jpj)             ,  h_rnf_tm(jpi,jpj) , STAT=trc_sub_alloc ) 
585      !
586      IF( trc_sub_alloc /= 0 )   CALL ctl_warn('trc_sub_alloc: failed to allocate arrays')
587      !
588      IF( l_ldfslp ) THEN
589         ALLOCATE( uslp_temp(jpi,jpj,jpk)   ,  wslpi_temp(jpi,jpj,jpk),      &
590            &      vslp_temp(jpi,jpj,jpk)   ,  wslpj_temp(jpi,jpj,jpk),      &
591            &      uslp_tm  (jpi,jpj,jpk)   ,  wslpi_tm  (jpi,jpj,jpk),      &
592            &      vslp_tm  (jpi,jpj,jpk)   ,  wslpj_tm  (jpi,jpj,jpk),  STAT=trc_sub_alloc )
593      ENDIF
594      !
595      IF( trc_sub_alloc /= 0 )   CALL ctl_warn('trc_sub_alloc: failed to allocate ldf_slp arrays')
596      !
597   END FUNCTION trc_sub_alloc
598
599#else
600   !!----------------------------------------------------------------------
601   !!   Default key                                     NO passive tracers
602   !!----------------------------------------------------------------------
603CONTAINS
604   SUBROUTINE trc_sub_stp( kt )        ! Empty routine
605      WRITE(*,*) 'trc_sub_stp: You should not have seen this print! error?', kt
606   END SUBROUTINE trc_sub_stp
607   SUBROUTINE trc_sub_ini        ! Empty routine
608      WRITE(*,*) 'trc_sub_ini: You should not have seen this print! error?', kt
609   END SUBROUTINE trc_sub_ini
610
611#endif
612
613   !!======================================================================
614END MODULE trcsub
Note: See TracBrowser for help on using the repository browser.