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.
trdicp.F90 in branches/2011/dev_r2787_LOCEAN3_TRA_TRP/NEMOGCM/NEMO/OPA_SRC/TRD – NEMO

source: branches/2011/dev_r2787_LOCEAN3_TRA_TRP/NEMOGCM/NEMO/OPA_SRC/TRD/trdicp.F90 @ 2789

Last change on this file since 2789 was 2789, checked in by cetlod, 13 years ago

Implementation of the merge of TRA/TRP : first guess, see ticket #842

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