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

source: branches/UKMO/dev_r5518_GO6_under_ice_relax_dr_hook/NEMOGCM/NEMO/TOP_SRC/trcsub.F90

Last change on this file was 11738, checked in by marc, 5 years ago

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