source: trunk/NEMOGCM/NEMO/OPA_SRC/TRD/trdicp.F90 @ 3294

Last change on this file since 3294 was 3294, checked in by rblod, 9 years ago

Merge of 3.4beta into the trunk

  • Property svn:keywords set to Id
File size: 29.9 KB
Line 
1MODULE trdicp
2   !!======================================================================
3   !!                       ***  MODULE  trdicp  ***
4   !! Ocean diagnostics:  ocean tracers and dynamic trends
5   !!=====================================================================
6   !! History :  1.0  !  2004-08 (C. Talandier) New trends organization
7   !!----------------------------------------------------------------------
8#if  defined key_trdtra   ||   defined key_trddyn   ||   defined key_esopa
9   !!----------------------------------------------------------------------
10   !!   'key_trdtra'  or                  active tracers trends diagnostics
11   !!   'key_trddyn'                            momentum trends diagnostics
12   !!----------------------------------------------------------------------
13   !!   trd_icp          : compute the basin averaged properties for tra/dyn
14   !!   trd_dwr          : print dynmaic trends in ocean.output file
15   !!   trd_twr          : print tracers trends in ocean.output file
16   !!   trd_icp_init     : initialization step
17   !!----------------------------------------------------------------------
18   USE oce             ! ocean dynamics and tracers variables
19   USE dom_oce         ! ocean space and time domain variables
20   USE trdmod_oce      ! ocean variables trends
21   USE ldftra_oce      ! ocean active tracers: lateral physics
22   USE ldfdyn_oce      ! ocean dynamics: lateral physics
23   USE zdf_oce         ! ocean vertical physics
24   USE in_out_manager  ! I/O manager
25   USE lib_mpp         ! distibuted memory computing library
26   USE eosbn2          ! equation of state
27   USE phycst          ! physical constants
28   USE wrk_nemo        ! Memory allocation
29
30
31   IMPLICIT NONE
32   PRIVATE
33
34   INTERFACE trd_icp
35      MODULE PROCEDURE trd_2d, trd_3d
36   END INTERFACE
37
38   PUBLIC   trd_icp       ! called by trdmod.F90
39   PUBLIC   trd_dwr       ! called by step.F90
40   PUBLIC   trd_twr       ! called by step.F90
41   PUBLIC   trd_icp_init  ! called by opa.F90
42
43   !! * Substitutions
44#  include "domzgr_substitute.h90"
45#  include "vectopt_loop_substitute.h90"
46   !!----------------------------------------------------------------------
47   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
48   !! $Id$
49   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
50   !!----------------------------------------------------------------------
51CONTAINS
52
53   SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype )
54      !!---------------------------------------------------------------------
55      !!                  ***  ROUTINE trd_2d  ***
56      !!
57      !! ** Purpose : verify the basin averaged properties of tracers and/or
58      !!              momentum equations at every time step frequency nn_trd.
59      !!----------------------------------------------------------------------
60      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ptrd2dx   ! Temperature or U trend
61      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ptrd2dy   ! Salinity    or V trend
62      INTEGER                     , INTENT(in   ) ::   ktrd      ! tracer trend index
63      CHARACTER(len=3)            , INTENT(in   ) ::   ctype     ! momentum ('DYN') or tracers ('TRA') trends
64      !!
65      INTEGER ::   ji, jj   ! loop indices
66      !!----------------------------------------------------------------------
67
68      SELECT CASE( ctype )    !==  Mask trends  ==!
69      !
70      CASE( 'DYN' )                    ! Momentum
71         DO jj = 1, jpjm1
72            DO ji = 1, jpim1
73               ptrd2dx(ji,jj) = ptrd2dx(ji,jj) * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,1)
74               ptrd2dy(ji,jj) = ptrd2dy(ji,jj) * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,1)
75            END DO
76         END DO
77         ptrd2dx(jpi, : ) = 0._wp      ;      ptrd2dy(jpi, : ) = 0._wp
78         ptrd2dx( : ,jpj) = 0._wp      ;      ptrd2dy( : ,jpj) = 0._wp
79         !
80      CASE( 'TRA' )                    ! Tracers
81         ptrd2dx(:,:) = ptrd2dx(:,:) * tmask_i(:,:)
82         ptrd2dy(:,:) = ptrd2dy(:,:) * tmask_i(:,:)
83         !
84      END SELECT
85     
86      SELECT CASE( ctype )    !==  Basin averaged tracer/momentum trends  ==!
87      !
88      CASE( 'DYN' )                    ! Momentum
89         umo(ktrd) = 0._wp
90         vmo(ktrd) = 0._wp
91         !
92         SELECT CASE( ktrd )
93         CASE( jpdyn_trd_swf )         ! surface forcing
94            umo(ktrd) = SUM( ptrd2dx(:,:) * e1u(:,:) * e2u(:,:) * fse3u(:,:,1) )
95            vmo(ktrd) = SUM( ptrd2dy(:,:) * e1v(:,:) * e2v(:,:) * fse3v(:,:,1) )
96         END SELECT
97         !
98      CASE( 'TRA' )              ! Tracers
99         tmo(ktrd) = SUM( ptrd2dx(:,:) * e1e2t(:,:) * fse3t(:,:,1) )
100         smo(ktrd) = SUM( ptrd2dy(:,:) * e1e2t(:,:) * fse3t(:,:,1) )
101      END SELECT
102     
103      SELECT CASE( ctype )    !==  Basin averaged tracer/momentum square trends  ==!   (now field)
104      !
105      CASE( 'DYN' )              ! Momentum
106         hke(ktrd) = SUM(   un(:,:,1) * ptrd2dx(:,:) * e1u(:,:) * e2u(:,:) * fse3u(:,:,1)   &
107            &             + vn(:,:,1) * ptrd2dy(:,:) * e1v(:,:) * e2v(:,:) * fse3v(:,:,1)   )
108         !
109      CASE( 'TRA' )              ! Tracers
110         t2(ktrd) = SUM( ptrd2dx(:,:) * e1e2t(:,:) * fse3t(:,:,1) * tsn(:,:,1,jp_tem) )
111         s2(ktrd) = SUM( ptrd2dy(:,:) * e1e2t(:,:) * fse3t(:,:,1) * tsn(:,:,1,jp_sal) )
112         !     
113      END SELECT
114      !
115   END SUBROUTINE trd_2d
116
117
118   SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd, ctype )
119      !!---------------------------------------------------------------------
120      !!                  ***  ROUTINE trd_3d  ***
121      !!
122      !! ** Purpose : verify the basin averaged properties of tracers and/or
123      !!              momentum equations at every time step frequency nn_trd.
124      !!----------------------------------------------------------------------
125      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptrd3dx   ! Temperature or U trend
126      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptrd3dy   ! Salinity    or V trend
127      INTEGER,                          INTENT(in   ) ::   ktrd      ! momentum or tracer trend index
128      CHARACTER(len=3),                 INTENT(in   ) ::   ctype     ! momentum ('DYN') or tracers ('TRA') trends
129      !!
130      INTEGER  ::   ji, jj, jk   ! dummy loop indices
131      !!----------------------------------------------------------------------
132
133      SELECT CASE( ctype )    !==  Mask the trends  ==!
134      !
135      CASE( 'DYN' )              ! Momentum       
136         DO jk = 1, jpkm1
137            DO jj = 1, jpjm1
138               DO ji = 1, jpim1
139                  ptrd3dx(ji,jj,jk) = ptrd3dx(ji,jj,jk) * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)
140                  ptrd3dy(ji,jj,jk) = ptrd3dy(ji,jj,jk) * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
141               END DO
142            END DO
143         END DO
144         ptrd3dx(jpi, : ,:) = 0._wp      ;      ptrd3dy(jpi, : ,:) = 0._wp
145         ptrd3dx( : ,jpj,:) = 0._wp      ;      ptrd3dy( : ,jpj,:) = 0._wp
146         !
147      CASE( 'TRA' )              ! Tracers
148         DO jk = 1, jpkm1
149            ptrd3dx(:,:,jk) = ptrd3dx(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
150            ptrd3dy(:,:,jk) = ptrd3dy(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
151         END DO
152         !
153      END SELECT   
154
155      SELECT CASE( ctype )    !==  Basin averaged tracer/momentum trends  ==!
156      !
157      CASE( 'DYN' )              ! Momentum
158         umo(ktrd) = 0._wp
159         vmo(ktrd) = 0._wp
160         DO jk = 1, jpkm1
161            umo(ktrd) = umo(ktrd) + SUM( ptrd3dx(:,:,jk) * e1u(:,:) * e2u(:,:) * fse3u(:,:,jk) )
162            vmo(ktrd) = vmo(ktrd) + SUM( ptrd3dy(:,:,jk) * e1v(:,:) * e2v(:,:) * fse3v(:,:,jk) )
163         END DO
164         !
165      CASE( 'TRA' )              ! Tracers
166         tmo(ktrd) = 0._wp
167         smo(ktrd) = 0._wp
168         DO jk = 1, jpkm1
169            tmo(ktrd) = tmo(ktrd) + SUM( ptrd3dx(:,:,jk) * e1e2t(:,:) * fse3t(:,:,jk) )
170            smo(ktrd) = smo(ktrd) + SUM( ptrd3dy(:,:,jk) * e1e2t(:,:) * fse3t(:,:,jk) )
171         END DO
172         !
173      END SELECT
174
175      SELECT CASE( ctype )    !==  Basin averaged tracer/momentum square trends  ==!   (now field)
176      !
177      CASE( 'DYN' )              ! Momentum
178         hke(ktrd) = 0._wp
179         DO jk = 1, jpkm1
180            hke(ktrd) = hke(ktrd) + SUM(   un(:,:,jk) * ptrd3dx(:,:,jk) * e1u(:,:) * e2u(:,:) * fse3u(:,:,jk)   &
181               &                         + vn(:,:,jk) * ptrd3dy(:,:,jk) * e1v(:,:) * e2v(:,:) * fse3v(:,:,jk)   )
182         END DO
183         !
184      CASE( 'TRA' )              ! Tracers
185         t2(ktrd) = 0._wp
186         s2(ktrd) = 0._wp
187         DO jk = 1, jpkm1
188            t2(ktrd) = t2(ktrd) + SUM( ptrd3dx(:,:,jk) * tsn(:,:,jk,jp_tem) * e1e2t(:,:) * fse3t(:,:,jk) )
189            s2(ktrd) = s2(ktrd) + SUM( ptrd3dy(:,:,jk) * tsn(:,:,jk,jp_sal) * e1e2t(:,:) * fse3t(:,:,jk) )
190         END DO
191         !
192      END SELECT
193      !
194   END SUBROUTINE trd_3d
195
196
197   SUBROUTINE trd_icp_init
198      !!---------------------------------------------------------------------
199      !!                  ***  ROUTINE trd_icp_init  ***
200      !!
201      !! ** Purpose :   Read the namtrd namelist
202      !!----------------------------------------------------------------------
203      INTEGER  ::   ji, jj, jk   ! dummy loop indices
204      !!----------------------------------------------------------------------
205
206      IF(lwp) THEN
207         WRITE(numout,*)
208         WRITE(numout,*) 'trd_icp_init : integral constraints properties trends'
209         WRITE(numout,*) '~~~~~~~~~~~~~'
210      ENDIF
211
212      ! Total volume at t-points:
213      tvolt = 0._wp
214      DO jk = 1, jpkm1
215         tvolt = SUM( e1e2t(:,:) * fse3t(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:) )
216      END DO
217      IF( lk_mpp )   CALL mpp_sum( tvolt )   ! sum over the global domain
218
219      IF(lwp) WRITE(numout,*) '                total ocean volume at T-point   tvolt = ',tvolt
220
221#if  defined key_trddyn
222      ! Initialization of potential to kinetic energy conversion
223      rpktrd = 0._wp
224
225      ! Total volume at u-, v- points:
226      tvolu = 0._wp
227      tvolv = 0._wp
228
229      DO jk = 1, jpk
230         DO jj = 2, jpjm1
231            DO ji = fs_2, fs_jpim1   ! vector opt.
232               tvolu = tvolu + e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) * tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)
233               tvolv = tvolv + e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) * tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
234            END DO
235         END DO
236      END DO
237      IF( lk_mpp )   CALL mpp_sum( tvolu )   ! sums over the global domain
238      IF( lk_mpp )   CALL mpp_sum( tvolv )
239
240      IF(lwp) THEN
241         WRITE(numout,*) '                total ocean volume at U-point   tvolu = ',tvolu
242         WRITE(numout,*) '                total ocean volume at V-point   tvolv = ',tvolv
243      ENDIF
244#endif
245      !
246   END SUBROUTINE trd_icp_init
247
248
249   SUBROUTINE trd_dwr( kt )
250      !!---------------------------------------------------------------------
251      !!                  ***  ROUTINE trd_dwr  ***
252      !!
253      !! ** Purpose :  write dynamic trends in ocean.output
254      !!----------------------------------------------------------------------
255      !
256      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
257      !
258      INTEGER  ::   ji, jj, jk   ! dummy loop indices
259      REAL(wp) ::   zcof         ! local scalar
260      REAL(wp), POINTER, DIMENSION(:,:,:)  ::  zkx, zky, zkz, zkepe 
261      !!----------------------------------------------------------------------
262
263      CALL wrk_alloc( jpi, jpj, jpk, zkx, zky, zkz, zkepe )
264
265      ! I. Momentum trends
266      ! -------------------
267
268      IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
269
270         ! I.1 Conversion potential energy - kinetic energy
271         ! --------------------------------------------------
272         ! c a u t i o n here, trends are computed at kt+1 (now , but after the swap)
273         zkx  (:,:,:) = 0._wp
274         zky  (:,:,:) = 0._wp
275         zkz  (:,:,:) = 0._wp
276         zkepe(:,:,:) = 0._wp
277   
278         CALL eos( tsn, rhd, rhop )       ! now potential and in situ densities
279
280         zcof = 0.5_wp / rau0             ! Density flux at w-point
281         zkz(:,:,1) = 0._wp
282         DO jk = 2, jpk
283            zkz(:,:,jk) = e1e2t(:,:) * wn(:,:,jk) * ( rhop(:,:,jk) + rhop(:,:,jk-1) ) * tmask_i(:,:)
284         END DO
285         
286         zcof   = 0.5_wp / rau0           ! Density flux at u and v-points
287         DO jk = 1, jpkm1
288            DO jj = 1, jpjm1
289               DO ji = 1, jpim1
290                  zkx(ji,jj,jk) = zcof * e2u(ji,jj) * fse3u(ji,jj,jk) * un(ji,jj,jk) * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) )
291                  zky(ji,jj,jk) = zcof * e1v(ji,jj) * fse3v(ji,jj,jk) * vn(ji,jj,jk) * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) )
292               END DO
293            END DO
294         END DO
295         
296         DO jk = 1, jpkm1                 ! Density flux divergence at t-point
297            DO jj = 2, jpjm1
298               DO ji = 2, jpim1
299                  zkepe(ji,jj,jk) = - (  zkz(ji,jj,jk) - zkz(ji  ,jj  ,jk+1)               &
300                     &                 + zkx(ji,jj,jk) - zkx(ji-1,jj  ,jk  )               &
301                     &                 + zky(ji,jj,jk) - zky(ji  ,jj-1,jk  )   )           &
302                     &              / ( e1e2t(ji,jj) * fse3t(ji,jj,jk) ) * tmask(ji,jj,jk) * tmask_i(ji,jj)
303               END DO
304            END DO
305         END DO
306
307         ! I.2 Basin averaged kinetic energy trend
308         ! ----------------------------------------
309         peke = 0._wp
310         DO jk = 1, jpkm1
311            peke = peke + SUM( zkepe(:,:,jk) * fsdept(:,:,jk) * e1e2t(:,:) * fse3t(:,:,jk) )
312         END DO
313         peke = grav * peke
314
315         ! I.3 Sums over the global domain
316         ! ---------------------------------
317         IF( lk_mpp ) THEN
318            CALL mpp_sum( peke )
319            CALL mpp_sum( umo , jptot_dyn )
320            CALL mpp_sum( vmo , jptot_dyn )
321            CALL mpp_sum( hke , jptot_dyn )
322         ENDIF
323
324         ! I.2 Print dynamic trends in the ocean.output file
325         ! --------------------------------------------------
326
327         IF(lwp) THEN
328            WRITE (numout,*)
329            WRITE (numout,*)
330            WRITE (numout,9500) kt
331            WRITE (numout,9501) umo(jpicpd_hpg) / tvolu, vmo(jpicpd_hpg) / tvolv
332            WRITE (numout,9502) umo(jpicpd_keg) / tvolu, vmo(jpicpd_keg) / tvolv
333            WRITE (numout,9503) umo(jpicpd_rvo) / tvolu, vmo(jpicpd_rvo) / tvolv
334            WRITE (numout,9504) umo(jpicpd_pvo) / tvolu, vmo(jpicpd_pvo) / tvolv
335            WRITE (numout,9505) umo(jpicpd_ldf) / tvolu, vmo(jpicpd_ldf) / tvolv
336            WRITE (numout,9506) umo(jpicpd_had) / tvolu, vmo(jpicpd_had) / tvolv
337            WRITE (numout,9507) umo(jpicpd_zad) / tvolu, vmo(jpicpd_zad) / tvolv
338            WRITE (numout,9508) umo(jpicpd_zdf) / tvolu, vmo(jpicpd_zdf) / tvolv
339            WRITE (numout,9509) umo(jpicpd_spg) / tvolu, vmo(jpicpd_spg) / tvolv
340            WRITE (numout,9510) umo(jpicpd_swf) / tvolu, vmo(jpicpd_swf) / tvolv
341            WRITE (numout,9511) umo(jpicpd_dat) / tvolu, vmo(jpicpd_dat) / tvolv
342            WRITE (numout,9512) umo(jpicpd_bfr) / tvolu, vmo(jpicpd_bfr) / tvolv
343            WRITE (numout,9513)
344            WRITE (numout,9514)                                                 &
345            &     (  umo(jpicpd_hpg) + umo(jpicpd_keg) + umo(jpicpd_rvo) + umo(jpicpd_pvo) + umo(jpicpd_ldf)   &
346            &      + umo(jpicpd_had) + umo(jpicpd_zad) + umo(jpicpd_zdf) + umo(jpicpd_spg) + umo(jpicpd_dat)   &
347            &      + umo(jpicpd_swf) + umo(jpicpd_bfr) ) / tvolu,   &
348            &     (  vmo(jpicpd_hpg) + vmo(jpicpd_keg) + vmo(jpicpd_rvo) + vmo(jpicpd_pvo) + vmo(jpicpd_ldf)   &
349            &      + vmo(jpicpd_had) + vmo(jpicpd_zad) + vmo(jpicpd_zdf) + vmo(jpicpd_spg) + vmo(jpicpd_dat)   &
350            &      + vmo(jpicpd_swf) + vmo(jpicpd_bfr) ) / tvolv
351         ENDIF
352
353 9500    FORMAT(' momentum trend at it= ', i6, ' :', /' ==============================')
354 9501    FORMAT(' pressure gradient          u= ', e20.13, '    v= ', e20.13)
355 9502    FORMAT(' ke gradient                u= ', e20.13, '    v= ', e20.13)
356 9503    FORMAT(' relative vorticity term    u= ', e20.13, '    v= ', e20.13)
357 9504    FORMAT(' coriolis term              u= ', e20.13, '    v= ', e20.13)
358 9505    FORMAT(' horizontal diffusion       u= ', e20.13, '    v= ', e20.13)
359 9506    FORMAT(' horizontal advection       u= ', e20.13, '    v= ', e20.13)
360 9507    FORMAT(' vertical advection         u= ', e20.13, '    v= ', e20.13)
361 9508    FORMAT(' vertical diffusion         u= ', e20.13, '    v= ', e20.13)
362 9509    FORMAT(' surface pressure gradient  u= ', e20.13, '    v= ', e20.13)
363 9510    FORMAT(' surface wind forcing       u= ', e20.13, '    v= ', e20.13)
364 9511    FORMAT(' dampimg term               u= ', e20.13, '    v= ', e20.13)
365 9512    FORMAT(' bottom flux                u= ', e20.13, '    v= ', e20.13)
366 9513    FORMAT(' -----------------------------------------------------------------------------')
367 9514    FORMAT(' total trend                u= ', e20.13, '    v= ', e20.13)
368
369         IF(lwp) THEN
370            WRITE (numout,*)
371            WRITE (numout,*)
372            WRITE (numout,9520) kt
373            WRITE (numout,9521) hke(jpicpd_hpg) / tvolt
374            WRITE (numout,9522) hke(jpicpd_keg) / tvolt
375            WRITE (numout,9523) hke(jpicpd_rvo) / tvolt
376            WRITE (numout,9524) hke(jpicpd_pvo) / tvolt
377            WRITE (numout,9525) hke(jpicpd_ldf) / tvolt
378            WRITE (numout,9526) hke(jpicpd_had) / tvolt
379            WRITE (numout,9527) hke(jpicpd_zad) / tvolt
380            WRITE (numout,9528) hke(jpicpd_zdf) / tvolt
381            WRITE (numout,9529) hke(jpicpd_spg) / tvolt
382            WRITE (numout,9530) hke(jpicpd_swf) / tvolt
383            WRITE (numout,9531) hke(jpicpd_dat) / tvolt
384            WRITE (numout,9532) hke(jpicpd_bfr) / tvolt
385            WRITE (numout,9533)
386            WRITE (numout,9534)   &
387            &     (  hke(jpicpd_hpg) + hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_pvo) + hke(jpicpd_ldf)   &
388            &      + hke(jpicpd_had) + hke(jpicpd_zad) + hke(jpicpd_zdf) + hke(jpicpd_spg) + hke(jpicpd_dat)   &
389            &      + hke(jpicpd_swf) + hke(jpicpd_bfr) ) / tvolt
390         ENDIF
391
392 9520    FORMAT(' kinetic energy trend at it= ', i6, ' :', /' ====================================')
393 9521    FORMAT(' pressure gradient         u2= ', e20.13)
394 9522    FORMAT(' ke gradient               u2= ', e20.13)
395 9523    FORMAT(' relative vorticity term   u2= ', e20.13)
396 9524    FORMAT(' coriolis term             u2= ', e20.13)
397 9525    FORMAT(' horizontal diffusion      u2= ', e20.13)
398 9526    FORMAT(' horizontal advection      u2= ', e20.13)
399 9527    FORMAT(' vertical advection        u2= ', e20.13)
400 9528    FORMAT(' vertical diffusion        u2= ', e20.13)
401 9529    FORMAT(' surface pressure gradient u2= ', e20.13)
402 9530    FORMAT(' surface wind forcing      u2= ', e20.13)
403 9531    FORMAT(' dampimg term              u2= ', e20.13)
404 9532    FORMAT(' bottom flux               u2= ', e20.13)
405 9533    FORMAT(' --------------------------------------------------')
406 9534    FORMAT(' total trend               u2= ', e20.13)
407
408         IF(lwp) THEN
409            WRITE (numout,*)
410            WRITE (numout,*)
411            WRITE (numout,9540) kt
412            WRITE (numout,9541) ( hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_had) + hke(jpicpd_zad) ) / tvolt
413            WRITE (numout,9542) ( hke(jpicpd_keg) + hke(jpicpd_had) + hke(jpicpd_zad) ) / tvolt
414            WRITE (numout,9543) ( hke(jpicpd_pvo) ) / tvolt
415            WRITE (numout,9544) ( hke(jpicpd_rvo) ) / tvolt
416            WRITE (numout,9545) ( hke(jpicpd_spg) ) / tvolt
417            WRITE (numout,9546) ( hke(jpicpd_ldf) ) / tvolt
418            WRITE (numout,9547) ( hke(jpicpd_zdf) ) / tvolt
419            WRITE (numout,9548) ( hke(jpicpd_hpg) ) / tvolt, rpktrd / tvolt
420            WRITE (numout,*)
421            WRITE (numout,*)
422         ENDIF
423
424 9540    FORMAT(' energetic consistency at it= ', i6, ' :', /' =========================================')
425 9541    FORMAT(' 0 = non linear term(true if key_vorenergy or key_combined): ', e20.13)
426 9542    FORMAT(' 0 = ke gradient + horizontal + vertical advection         : ', e20.13)
427 9543    FORMAT(' 0 = coriolis term  (true if key_vorenergy or key_combined): ', e20.13)
428 9544    FORMAT(' 0 = uh.( rot(u) x uh ) (true if enstrophy conser.)        : ', e20.13)
429 9545    FORMAT(' 0 = surface pressure gradient                             : ', e20.13)
430 9546    FORMAT(' 0 > horizontal diffusion                                  : ', e20.13)
431 9547    FORMAT(' 0 > vertical diffusion                                    : ', e20.13)
432 9548    FORMAT(' pressure gradient u2 = - 1/rau0 u.dz(rhop)                : ', e20.13, '  u.dz(rhop) =', e20.13)
433         !
434         ! Save potential to kinetic energy conversion for next time step
435         rpktrd = peke
436         !
437      ENDIF
438      !
439      CALL wrk_dealloc( jpi, jpj, jpk, zkx, zky, zkz, zkepe )
440      !
441   END SUBROUTINE trd_dwr
442
443
444   SUBROUTINE trd_twr( kt )
445      !!---------------------------------------------------------------------
446      !!                  ***  ROUTINE trd_twr  ***
447      !!
448      !! ** Purpose :  write active tracers trends in ocean.output
449      !!----------------------------------------------------------------------
450      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
451      !!----------------------------------------------------------------------
452
453      ! I. Tracers trends
454      ! -----------------
455
456      IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
457
458         ! I.1 Sums over the global domain
459         ! -------------------------------
460         IF( lk_mpp ) THEN
461            CALL mpp_sum( tmo, jptot_tra )   
462            CALL mpp_sum( smo, jptot_tra )
463            CALL mpp_sum( t2 , jptot_tra )
464            CALL mpp_sum( s2 , jptot_tra )
465         ENDIF
466
467         ! I.2 Print tracers trends in the ocean.output file
468         ! --------------------------------------------------
469         
470         IF(lwp) THEN
471            WRITE (numout,*)
472            WRITE (numout,*)
473            WRITE (numout,9400) kt
474            WRITE (numout,9401)  tmo(jpicpt_xad) / tvolt, smo(jpicpt_xad) / tvolt
475            WRITE (numout,9411)  tmo(jpicpt_yad) / tvolt, smo(jpicpt_yad) / tvolt
476            WRITE (numout,9402)  tmo(jpicpt_zad) / tvolt, smo(jpicpt_zad) / tvolt
477            WRITE (numout,9403)  tmo(jpicpt_ldf) / tvolt, smo(jpicpt_ldf) / tvolt
478            WRITE (numout,9404)  tmo(jpicpt_zdf) / tvolt, smo(jpicpt_zdf) / tvolt
479            WRITE (numout,9405)  tmo(jpicpt_npc) / tvolt, smo(jpicpt_npc) / tvolt
480            WRITE (numout,9406)  tmo(jpicpt_dmp) / tvolt, smo(jpicpt_dmp) / tvolt
481            WRITE (numout,9407)  tmo(jpicpt_qsr) / tvolt
482            WRITE (numout,9408)  tmo(jpicpt_nsr) / tvolt, smo(jpicpt_nsr) / tvolt
483            WRITE (numout,9409) 
484            WRITE (numout,9410) (  tmo(jpicpt_xad) + tmo(jpicpt_yad) + tmo(jpicpt_zad) + tmo(jpicpt_ldf) + tmo(jpicpt_zdf)   &
485            &                    + tmo(jpicpt_npc) + tmo(jpicpt_dmp) + tmo(jpicpt_qsr) + tmo(jpicpt_nsr) ) / tvolt,   &
486            &                   (  smo(jpicpt_xad) + smo(jpicpt_yad) + smo(jpicpt_zad) + smo(jpicpt_ldf) + smo(jpicpt_zdf)   &
487            &                    + smo(jpicpt_npc) + smo(jpicpt_dmp)                   + smo(jpicpt_nsr) ) / tvolt
488         ENDIF
489
4909400     FORMAT(' tracer trend at it= ',i6,' :     temperature',   &
491              '              salinity',/' ============================')
4929401     FORMAT(' zonal      advection        ',e20.13,'     ',e20.13)
4939411     FORMAT(' meridional advection        ',e20.13,'     ',e20.13)
4949402     FORMAT(' vertical advection          ',e20.13,'     ',e20.13)
4959403     FORMAT(' horizontal diffusion        ',e20.13,'     ',e20.13)
4969404     FORMAT(' vertical diffusion          ',e20.13,'     ',e20.13)
4979405     FORMAT(' static instability mixing   ',e20.13,'     ',e20.13)
4989406     FORMAT(' damping term                ',e20.13,'     ',e20.13)
4999407     FORMAT(' penetrative qsr             ',e20.13)
5009408     FORMAT(' non solar radiation         ',e20.13,'     ',e20.13)
5019409     FORMAT(' -------------------------------------------------------------------------')
5029410     FORMAT(' total trend                 ',e20.13,'     ',e20.13)
503
504
505         IF(lwp) THEN
506            WRITE (numout,*)
507            WRITE (numout,*)
508            WRITE (numout,9420) kt
509            WRITE (numout,9421)   t2(jpicpt_xad) / tvolt, s2(jpicpt_xad) / tvolt
510            WRITE (numout,9431)   t2(jpicpt_yad) / tvolt, s2(jpicpt_yad) / tvolt
511            WRITE (numout,9422)   t2(jpicpt_zad) / tvolt, s2(jpicpt_zad) / tvolt
512            WRITE (numout,9423)   t2(jpicpt_ldf) / tvolt, s2(jpicpt_ldf) / tvolt
513            WRITE (numout,9424)   t2(jpicpt_zdf) / tvolt, s2(jpicpt_zdf) / tvolt
514            WRITE (numout,9425)   t2(jpicpt_npc) / tvolt, s2(jpicpt_npc) / tvolt
515            WRITE (numout,9426)   t2(jpicpt_dmp) / tvolt, s2(jpicpt_dmp) / tvolt
516            WRITE (numout,9427)   t2(jpicpt_qsr) / tvolt
517            WRITE (numout,9428)   t2(jpicpt_nsr) / tvolt, s2(jpicpt_nsr) / tvolt
518            WRITE (numout,9429)
519            WRITE (numout,9430) (  t2(jpicpt_xad) + t2(jpicpt_yad) + t2(jpicpt_zad) + t2(jpicpt_ldf) + t2(jpicpt_zdf)   &
520            &                    + t2(jpicpt_npc) + t2(jpicpt_dmp) + t2(jpicpt_qsr) + t2(jpicpt_nsr) ) / tvolt,   &
521            &                   (  s2(jpicpt_xad) + s2(jpicpt_yad) + s2(jpicpt_zad) + s2(jpicpt_ldf) + s2(jpicpt_zdf)   &
522            &                    + s2(jpicpt_npc) + s2(jpicpt_dmp)                  + s2(jpicpt_nsr) ) / tvolt
523         ENDIF
524
5259420     FORMAT(' tracer**2 trend at it= ', i6, ' :      temperature',   &
526            '               salinity', /, ' ===============================')
5279421     FORMAT(' zonal      advection      * t   ', e20.13, '     ', e20.13)
5289431     FORMAT(' meridional advection      * t   ', e20.13, '     ', e20.13)
5299422     FORMAT(' vertical advection        * t   ', e20.13, '     ', e20.13)
5309423     FORMAT(' horizontal diffusion      * t   ', e20.13, '     ', e20.13)
5319424     FORMAT(' vertical diffusion        * t   ', e20.13, '     ', e20.13)
5329425     FORMAT(' static instability mixing * t   ', e20.13, '     ', e20.13)
5339426     FORMAT(' damping term              * t   ', e20.13, '     ', e20.13)
5349427     FORMAT(' penetrative qsr           * t   ', e20.13)
5359428     FORMAT(' non solar radiation       * t   ', e20.13, '     ', e20.13)
5369429     FORMAT(' -----------------------------------------------------------------------------')
5379430     FORMAT(' total trend                *t = ', e20.13, '  *s = ', e20.13)
538
539
540         IF(lwp) THEN
541            WRITE (numout,*)
542            WRITE (numout,*)
543            WRITE (numout,9440) kt
544            WRITE (numout,9441) ( tmo(jpicpt_xad)+tmo(jpicpt_yad)+tmo(jpicpt_zad) )/tvolt,   &
545            &                   ( smo(jpicpt_xad)+smo(jpicpt_yad)+smo(jpicpt_zad) )/tvolt
546            WRITE (numout,9442)   tmo(jpicpt_zl1)/tvolt,  smo(jpicpt_zl1)/tvolt
547            WRITE (numout,9443)   tmo(jpicpt_ldf)/tvolt,  smo(jpicpt_ldf)/tvolt
548            WRITE (numout,9444)   tmo(jpicpt_zdf)/tvolt,  smo(jpicpt_zdf)/tvolt
549            WRITE (numout,9445)   tmo(jpicpt_npc)/tvolt,  smo(jpicpt_npc)/tvolt
550            WRITE (numout,9446) ( t2(jpicpt_xad)+t2(jpicpt_yad)+t2(jpicpt_zad) )/tvolt,    &
551            &                   ( s2(jpicpt_xad)+s2(jpicpt_yad)+s2(jpicpt_zad) )/tvolt
552            WRITE (numout,9447)   t2(jpicpt_ldf)/tvolt,   s2(jpicpt_ldf)/tvolt
553            WRITE (numout,9448)   t2(jpicpt_zdf)/tvolt,   s2(jpicpt_zdf)/tvolt
554            WRITE (numout,9449)   t2(jpicpt_npc)/tvolt,   s2(jpicpt_npc)/tvolt
555         ENDIF
556
5579440     FORMAT(' tracer consistency at it= ',i6,   &
558            ' :         temperature','                salinity',/,   &
559            ' ==================================')
5609441     FORMAT(' 0 = horizontal+vertical advection +    ',e20.13,'       ',e20.13)
5619442     FORMAT('     1st lev vertical advection         ',e20.13,'       ',e20.13)
5629443     FORMAT(' 0 = horizontal diffusion               ',e20.13,'       ',e20.13)
5639444     FORMAT(' 0 = vertical diffusion                 ',e20.13,'       ',e20.13)
5649445     FORMAT(' 0 = static instability mixing          ',e20.13,'       ',e20.13)
5659446     FORMAT(' 0 = horizontal+vertical advection * t  ',e20.13,'       ',e20.13)
5669447     FORMAT(' 0 > horizontal diffusion          * t  ',e20.13,'       ',e20.13)
5679448     FORMAT(' 0 > vertical diffusion            * t  ',e20.13,'       ',e20.13)
5689449     FORMAT(' 0 > static instability mixing     * t  ',e20.13,'       ',e20.13)
569         !
570      ENDIF
571      !
572   END SUBROUTINE trd_twr
573
574#   else
575   !!----------------------------------------------------------------------
576   !!   Default case :                                         Empty module
577   !!----------------------------------------------------------------------
578   INTERFACE trd_icp
579      MODULE PROCEDURE trd_2d, trd_3d
580   END INTERFACE
581
582CONTAINS
583   SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype )       ! Empty routine
584      REAL, DIMENSION(:,:) ::   ptrd2dx, ptrd2dy
585      INTEGER                     , INTENT(in   ) ::   ktrd         ! tracer trend index
586      CHARACTER(len=3)            , INTENT(in   ) ::   ctype        ! momentum ('DYN') or tracers ('TRA') trends
587      WRITE(*,*) 'trd_2d: You should not have seen this print! error ?', &
588          &       ptrd2dx(1,1), ptrd2dy(1,1), ktrd, ctype
589   END SUBROUTINE trd_2d
590   SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd , ctype )       ! Empty routine
591      REAL, DIMENSION(:,:,:) ::   ptrd3dx, ptrd3dy
592      INTEGER                     , INTENT(in   ) ::   ktrd         ! tracer trend index
593      CHARACTER(len=3)            , INTENT(in   ) ::   ctype        ! momentum ('DYN') or tracers ('TRA') trends
594      WRITE(*,*) 'trd_3d: You should not have seen this print! error ?', &
595          &       ptrd3dx(1,1,1), ptrd3dy(1,1,1), ktrd, ctype
596   END SUBROUTINE trd_3d
597   SUBROUTINE trd_icp_init               ! Empty routine
598   END SUBROUTINE trd_icp_init
599   SUBROUTINE trd_dwr( kt )          ! Empty routine
600      WRITE(*,*) 'trd_dwr: You should not have seen this print! error ?', kt
601   END SUBROUTINE trd_dwr
602   SUBROUTINE trd_twr( kt )          ! Empty routine
603      WRITE(*,*) 'trd_twr: You should not have seen this print! error ?', kt
604   END SUBROUTINE trd_twr
605#endif
606
607   !!======================================================================
608END MODULE trdicp
Note: See TracBrowser for help on using the repository browser.