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

source: branches/2016/dev_r7012_ROBUST5_CMCC/NEMOGCM/NEMO/TOP_SRC/trcsub.F90 @ 7079

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

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