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

Last change on this file since 7091 was 7091, checked in by lovato, 8 years ago

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