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 @ 5385

Last change on this file since 5385 was 5215, checked in by nicolasmartin, 9 years ago

SVN keywords tag inconsistencies resolution

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