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/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/TRD – NEMO

source: branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/TRD/trdicp.F90 @ 2082

Last change on this file since 2082 was 2082, checked in by cetlod, 14 years ago

Improve the merge of TRA-TRC, see ticket #717
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 32.8 KB
Line 
1MODULE trdicp
2   !!======================================================================
3   !!                       ***  MODULE  trdicp  ***
4   !! Ocean diagnostics:  ocean tracers and dynamic trends
5   !!=====================================================================
6   !! History :       !  91-12 (G. Madec)
7   !!                 !  92-06 (M. Imbard) add time step frequency
8   !!                 !  96-01 (G. Madec)  terrain following coordinates
9   !!            8.5  !  02-06 (G. Madec)  F90: Free form and module
10   !!            9.0  !  04-08 (C. Talandier) New trends organization
11   !!----------------------------------------------------------------------
12#if  defined key_trdtra   ||   defined key_trddyn   ||   defined key_esopa
13   !!----------------------------------------------------------------------
14   !!   'key_trdtra'  or                  active tracers trends diagnostics
15   !!   'key_trddyn'                            momentum trends diagnostics
16   !!----------------------------------------------------------------------
17   !!----------------------------------------------------------------------
18   !!   trd_icp          : compute the basin averaged properties for tra/dyn
19   !!   trd_dwr          : print dynmaic trends in ocean.output file
20   !!   trd_twr          : print tracers trends in ocean.output file
21   !!   trd_icp_init     : initialization step
22   !!----------------------------------------------------------------------
23   USE oce             ! ocean dynamics and tracers variables
24   USE dom_oce         ! ocean space and time domain variables
25   USE trdmod_oce      ! ocean variables trends
26   USE ldftra_oce      ! ocean active tracers: lateral physics
27   USE ldfdyn_oce      ! ocean dynamics: lateral physics
28   USE zdf_oce         ! ocean vertical physics
29   USE in_out_manager  ! I/O manager
30   USE lib_mpp         ! distibuted memory computing library
31   USE eosbn2          ! equation of state
32   USE phycst          ! physical constants
33
34   IMPLICIT NONE
35   PRIVATE
36
37   INTERFACE trd_icp
38      MODULE PROCEDURE trd_2d, trd_3d
39   END INTERFACE
40
41   PUBLIC   trd_icp       ! called by trdmod.F90
42   PUBLIC   trd_dwr       ! called by step.F90
43   PUBLIC   trd_twr       ! called by step.F90
44   PUBLIC   trd_icp_init  ! called by opa.F90
45
46   !! * Substitutions
47#  include "domzgr_substitute.h90"
48#  include "vectopt_loop_substitute.h90"
49   !!----------------------------------------------------------------------
50   !!   OPA 9.0 , LOCEAN-IPSL (2005)
51   !! $Id$
52   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
53   !!----------------------------------------------------------------------
54
55CONTAINS
56
57   SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype )
58      !!---------------------------------------------------------------------
59      !!                  ***  ROUTINE trd_2d  ***
60      !!
61      !! ** Purpose : verify the basin averaged properties of tracers and/or
62      !!              momentum equations at every time step frequency nn_trd.
63      !!----------------------------------------------------------------------
64      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ptrd2dx             ! Temperature or U trend
65      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   ptrd2dy             ! Salinity    or V trend
66      INTEGER                     , INTENT(in   ) ::   ktrd                ! tracer trend index
67      CHARACTER(len=3)            , INTENT(in   ) ::   ctype               ! momentum ('DYN') or tracers ('TRA') trends
68      !!
69      INTEGER  ::   ji, jj                                                 ! loop indices
70      REAL(wp) ::   zmsku, zbtu, zbt                                       ! temporary scalars
71      REAL(wp) ::   zmskv, zbtv                                            !    "         "
72      !!----------------------------------------------------------------------
73
74
75      ! 1. Mask trends
76      ! --------------
77
78      SELECT CASE( ctype )
79      !
80      CASE( 'DYN' )              ! Momentum
81         DO jj = 1, jpjm1
82            DO ji = 1, jpim1
83               zmsku = tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,1)
84               zmskv = tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,1)
85               ptrd2dx(ji,jj) = ptrd2dx(ji,jj) * zmsku
86               ptrd2dy(ji,jj) = ptrd2dy(ji,jj) * zmskv
87            END DO
88         END DO
89         ptrd2dx(jpi, : ) = 0.e0      ;      ptrd2dy(jpi, : ) = 0.e0
90         ptrd2dx( : ,jpj) = 0.e0      ;      ptrd2dy( : ,jpj) = 0.e0
91         !
92      CASE( 'TRA' )              ! Tracers
93         ptrd2dx(:,:) = ptrd2dx(:,:) * tmask_i(:,:)
94         ptrd2dy(:,:) = ptrd2dy(:,:) * tmask_i(:,:)
95         !
96      END SELECT
97     
98      ! 2. Basin averaged tracer/momentum trends
99      ! ----------------------------------------
100
101      SELECT CASE( ctype )
102      !
103      CASE( 'DYN' )              ! Momentum
104         umo(ktrd) = 0.e0
105         vmo(ktrd) = 0.e0
106         !
107         SELECT CASE( ktrd )
108         !
109         CASE( jpdyn_trd_swf )         ! surface forcing
110            DO jj = 1, jpj
111               DO ji = 1, jpi
112                  umo(ktrd) = umo(ktrd) + ptrd2dx(ji,jj) * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1)
113                  vmo(ktrd) = vmo(ktrd) + ptrd2dy(ji,jj) * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1)
114               END DO
115            END DO
116            !
117         END SELECT
118         !
119      CASE( 'TRA' )              ! Tracers
120         tmo(ktrd) = 0.e0
121         smo(ktrd) = 0.e0
122         DO jj = 1, jpj
123            DO ji = 1, jpi
124               zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
125               tmo(ktrd) =  tmo(ktrd) + ptrd2dx(ji,jj) * zbt
126               smo(ktrd) =  smo(ktrd) + ptrd2dy(ji,jj) * zbt
127            END DO
128         END DO
129         !
130      END SELECT
131     
132      ! 3. Basin averaged tracer/momentum square trends
133      ! ----------------------------------------------
134      ! c a u t i o n: field now
135     
136      SELECT CASE( ctype )
137      !
138      CASE( 'DYN' )              ! Momentum
139         hke(ktrd) = 0.e0
140         DO jj = 1, jpj
141            DO ji = 1, jpi
142               zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1)
143               zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1)
144               hke(ktrd) = hke(ktrd)   &
145               &   + un(ji,jj,1) * ptrd2dx(ji,jj) * zbtu &
146               &   + vn(ji,jj,1) * ptrd2dy(ji,jj) * zbtv
147            END DO
148         END DO
149         !
150      CASE( 'TRA' )              ! Tracers
151         t2(ktrd) = 0.e0
152         s2(ktrd) = 0.e0
153         DO jj = 1, jpj
154            DO ji = 1, jpi
155               zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
156               t2(ktrd) = t2(ktrd) + ptrd2dx(ji,jj) * zbt * tn(ji,jj,1)
157               s2(ktrd) = s2(ktrd) + ptrd2dy(ji,jj) * zbt * sn(ji,jj,1)
158            END DO
159         END DO
160         !     
161      END SELECT
162      !
163   END SUBROUTINE trd_2d
164
165
166   SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd, ctype )
167      !!---------------------------------------------------------------------
168      !!                  ***  ROUTINE trd_3d  ***
169      !!
170      !! ** Purpose : verify the basin averaged properties of tracers and/or
171      !!              momentum equations at every time step frequency nn_trd.
172      !!----------------------------------------------------------------------
173      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptrd3dx            ! Temperature or U trend
174      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptrd3dy            ! Salinity    or V trend
175      INTEGER,                          INTENT(in   ) ::   ktrd               ! momentum or tracer trend index
176      CHARACTER(len=3),                 INTENT(in   ) ::   ctype              ! momentum ('DYN') or tracers ('TRA') trends
177      !!
178      INTEGER ::   ji, jj, jk
179      REAL(wp) ::   zbt, zbtu, zbtv, zmsku, zmskv                             ! temporary scalars
180      !!----------------------------------------------------------------------
181
182      ! 1. Mask the trends
183      ! ------------------
184
185      SELECT CASE( ctype )
186      !
187      CASE( 'DYN' )              ! Momentum       
188         DO jk = 1, jpk
189            DO jj = 1, jpjm1
190               DO ji = 1, jpim1
191                  zmsku = tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)
192                  zmskv = tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
193                  ptrd3dx(ji,jj,jk) = ptrd3dx(ji,jj,jk) * zmsku
194                  ptrd3dy(ji,jj,jk) = ptrd3dy(ji,jj,jk) * zmskv
195               END DO
196            END DO
197         END DO
198         ptrd3dx(jpi, : ,:) = 0.e0      ;      ptrd3dy(jpi, : ,:) = 0.e0
199         ptrd3dx( : ,jpj,:) = 0.e0      ;      ptrd3dy( : ,jpj,:) = 0.e0
200         !
201      CASE( 'TRA' )              ! Tracers
202         DO jk = 1, jpk
203            ptrd3dx(:,:,jk) = ptrd3dx(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
204            ptrd3dy(:,:,jk) = ptrd3dy(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
205         END DO
206         !
207      END SELECT   
208
209      ! 2. Basin averaged tracer/momentum trends
210      ! ----------------------------------------
211     
212      SELECT CASE( ctype )
213      !
214      CASE( 'DYN' )              ! Momentum
215         umo(ktrd) = 0.e0
216         vmo(ktrd) = 0.e0
217         DO jk = 1, jpk
218            DO jj = 1, jpj
219               DO ji = 1, jpi
220                  zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
221                  zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
222                  umo(ktrd) = umo(ktrd) + ptrd3dx(ji,jj,jk) * zbtu
223                  vmo(ktrd) = vmo(ktrd) + ptrd3dy(ji,jj,jk) * zbtv
224               END DO
225            END DO
226         END DO
227         !
228      CASE( 'TRA' )              ! Tracers
229         tmo(ktrd) = 0.e0
230         smo(ktrd) = 0.e0
231         DO jk = 1, jpkm1
232            DO jj = 1, jpj
233               DO ji = 1, jpi
234                  zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 
235                  tmo(ktrd) =  tmo(ktrd) + ptrd3dx(ji,jj,jk) * zbt
236                  smo(ktrd) =  smo(ktrd) + ptrd3dy(ji,jj,jk) * zbt
237               END DO
238            END DO
239         END DO
240         !
241      END SELECT
242
243      ! 3. Basin averaged tracer/momentum square trends
244      ! -----------------------------------------------
245      ! c a u t i o n: field now
246     
247      SELECT CASE( ctype )
248      !
249      CASE( 'DYN' )              ! Momentum
250         hke(ktrd) = 0.e0
251         DO jk = 1, jpk
252            DO jj = 1, jpj
253               DO ji = 1, jpi
254                  zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
255                  zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
256                  hke(ktrd) = hke(ktrd)   &
257                  &   + un(ji,jj,jk) * ptrd3dx(ji,jj,jk) * zbtu &
258                  &   + vn(ji,jj,jk) * ptrd3dy(ji,jj,jk) * zbtv
259               END DO
260            END DO
261         END DO
262         !
263      CASE( 'TRA' )              ! Tracers
264         t2(ktrd) = 0.e0
265         s2(ktrd) = 0.e0
266         DO jk = 1, jpk
267            DO jj = 1, jpj
268               DO ji = 1, jpi
269                  zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
270                  t2(ktrd) = t2(ktrd) + ptrd3dx(ji,jj,jk) * zbt * tn(ji,jj,jk)
271                  s2(ktrd) = s2(ktrd) + ptrd3dy(ji,jj,jk) * zbt * sn(ji,jj,jk)
272               END DO
273            END DO
274         END DO
275         !
276      END SELECT
277      !
278   END SUBROUTINE trd_3d
279
280
281
282   SUBROUTINE trd_icp_init
283      !!---------------------------------------------------------------------
284      !!                  ***  ROUTINE trd_icp_init  ***
285      !!
286      !! ** Purpose :   Read the namtrd namelist
287      !!----------------------------------------------------------------------
288      INTEGER  ::   ji, jj, jk
289      REAL(wp) ::   zmskt
290#if  defined key_trddyn
291      REAL(wp) ::   zmsku, zmskv
292#endif
293      !!----------------------------------------------------------------------
294
295      IF(lwp) THEN
296         WRITE(numout,*)
297         WRITE(numout,*) 'trd_icp_init : integral constraints properties trends'
298         WRITE(numout,*) '~~~~~~~~~~~~~'
299      ENDIF
300
301      ! Total volume at t-points:
302      tvolt = 0.e0
303      DO jk = 1, jpkm1
304         DO jj = 2, jpjm1
305            DO ji = fs_2, fs_jpim1   ! vector opt.
306               zmskt = tmask(ji,jj,jk) * tmask_i(ji,jj)
307               tvolt = tvolt + zmskt * e1t(ji,jj) *e2t(ji,jj) * fse3t(ji,jj,jk)
308            END DO
309         END DO
310      END DO
311      IF( lk_mpp )   CALL mpp_sum( tvolt )   ! sum over the global domain
312
313      IF(lwp) WRITE(numout,*) '                total ocean volume at T-point   tvolt = ',tvolt
314
315#if  defined key_trddyn
316      ! Initialization of potential to kinetic energy conversion
317      rpktrd = 0.e0
318
319      ! Total volume at u-, v- points:
320      tvolu = 0.e0
321      tvolv = 0.e0
322
323      DO jk = 1, jpk
324         DO jj = 2, jpjm1
325            DO ji = fs_2, fs_jpim1   ! vector opt.
326               zmsku = tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)
327               zmskv = tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
328               tvolu = tvolu + zmsku * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
329               tvolv = tvolv + zmskv * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
330            END DO
331         END DO
332      END DO
333      IF( lk_mpp )   CALL mpp_sum( tvolu )   ! sums over the global domain
334      IF( lk_mpp )   CALL mpp_sum( tvolv )
335
336      IF(lwp) THEN
337         WRITE(numout,*) '                total ocean volume at U-point   tvolu = ',tvolu
338         WRITE(numout,*) '                total ocean volume at V-point   tvolv = ',tvolv
339      ENDIF
340#endif
341      !
342   END SUBROUTINE trd_icp_init
343
344
345   SUBROUTINE trd_dwr( kt )
346      !!---------------------------------------------------------------------
347      !!                  ***  ROUTINE trd_dwr  ***
348      !!
349      !! ** Purpose :  write dynamic trends in ocean.output
350      !!----------------------------------------------------------------------
351      INTEGER, INTENT(in) ::   kt                                  ! ocean time-step index
352      !!
353      INTEGER  ::   ji, jj, jk
354      REAL(wp) ::   ze1e2w, zcof, zbe1ru, zbe2rv, zbtr, ztz, zth   !    "      scalars
355      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zkepe, zkx, zky, zkz   ! temporary arrays
356      !!----------------------------------------------------------------------
357
358      ! I. Momentum trends
359      ! -------------------
360
361      IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
362
363         ! I.1 Conversion potential energy - kinetic energy
364         ! --------------------------------------------------
365         ! c a u t i o n here, trends are computed at kt+1 (now , but after the swap)
366
367         zkx(:,:,:)   = 0.e0
368         zky(:,:,:)   = 0.e0
369         zkz(:,:,:)   = 0.e0
370         zkepe(:,:,:) = 0.e0
371   
372         CALL eos( tsn, rhd, rhop )       ! now potential and in situ densities
373
374         ! Density flux at w-point
375         DO jk = 2, jpk
376            DO jj = 1, jpj
377               DO ji = 1, jpi
378                  ze1e2w = 0.5 * e1t(ji,jj) * e2t(ji,jj) * wn(ji,jj,jk) * tmask_i(ji,jj)
379                  zkz(ji,jj,jk) = ze1e2w / rau0 * ( rhop(ji,jj,jk) + rhop(ji,jj,jk-1) )
380               END DO
381            END DO
382         END DO
383         zkz(:,:,1) = 0.e0
384         
385         ! Density flux at u and v-points
386         DO jk = 1, jpk
387            DO jj = 1, jpjm1
388               DO ji = 1, jpim1
389                  zcof   = 0.5 / rau0
390                  zbe1ru = zcof * e2u(ji,jj) * fse3u(ji,jj,jk) * un(ji,jj,jk)
391                  zbe2rv = zcof * e1v(ji,jj) * fse3v(ji,jj,jk) * vn(ji,jj,jk)
392                  zkx(ji,jj,jk) = zbe1ru * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) )
393                  zky(ji,jj,jk) = zbe2rv * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) )
394               END DO
395            END DO
396         END DO
397         
398         ! Density flux divergence at t-point
399         DO jk = 1, jpkm1
400            DO jj = 2, jpjm1
401               DO ji = 2, jpim1
402                  zbtr = 1. / ( e1t(ji,jj)*e2t(ji,jj)*fse3t(ji,jj,jk) )
403                  ztz = - zbtr * (    zkz(ji,jj,jk) - zkz(ji,jj,jk+1) )
404                  zth = - zbtr * (  ( zkx(ji,jj,jk) - zkx(ji-1,jj,jk) )   &
405                    &             + ( zky(ji,jj,jk) - zky(ji,jj-1,jk) )  )
406                  zkepe(ji,jj,jk) = (zth + ztz) * tmask(ji,jj,jk) * tmask_i(ji,jj)
407               END DO
408            END DO
409         END DO
410         zkepe( : , : ,jpk) = 0.e0
411         zkepe( : ,jpj, : ) = 0.e0
412         zkepe(jpi, : , : ) = 0.e0
413
414         ! I.2 Basin averaged kinetic energy trend
415         ! ----------------------------------------
416         peke = 0.e0
417         DO jk = 1,jpk
418            DO jj = 1, jpj
419               DO ji = 1, jpi
420                  peke = peke + zkepe(ji,jj,jk) * grav * fsdept(ji,jj,jk)   &
421                     &                     * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
422               END DO
423            END DO
424         END DO
425
426         ! I.3 Sums over the global domain
427         ! ---------------------------------
428         IF( lk_mpp ) THEN
429            CALL mpp_sum( peke )
430            CALL mpp_sum( umo , jptot_dyn )
431            CALL mpp_sum( vmo , jptot_dyn )
432            CALL mpp_sum( hke , jptot_dyn )
433         ENDIF
434
435         ! I.2 Print dynamic trends in the ocean.output file
436         ! --------------------------------------------------
437
438         IF(lwp) THEN
439            WRITE (numout,*)
440            WRITE (numout,*)
441            WRITE (numout,9500) kt
442            WRITE (numout,9501) umo(jpicpd_hpg) / tvolu, vmo(jpicpd_hpg) / tvolv
443            WRITE (numout,9502) umo(jpicpd_keg) / tvolu, vmo(jpicpd_keg) / tvolv
444            WRITE (numout,9503) umo(jpicpd_rvo) / tvolu, vmo(jpicpd_rvo) / tvolv
445            WRITE (numout,9504) umo(jpicpd_pvo) / tvolu, vmo(jpicpd_pvo) / tvolv
446            WRITE (numout,9505) umo(jpicpd_ldf) / tvolu, vmo(jpicpd_ldf) / tvolv
447            WRITE (numout,9506) umo(jpicpd_had) / tvolu, vmo(jpicpd_had) / tvolv
448            WRITE (numout,9507) umo(jpicpd_zad) / tvolu, vmo(jpicpd_zad) / tvolv
449            WRITE (numout,9508) umo(jpicpd_zdf) / tvolu, vmo(jpicpd_zdf) / tvolv
450            WRITE (numout,9509) umo(jpicpd_spg) / tvolu, vmo(jpicpd_spg) / tvolv
451            WRITE (numout,9510) umo(jpicpd_swf) / tvolu, vmo(jpicpd_swf) / tvolv
452            WRITE (numout,9511) umo(jpicpd_dat) / tvolu, vmo(jpicpd_dat) / tvolv
453            WRITE (numout,9512) umo(jpicpd_bfr) / tvolu, vmo(jpicpd_bfr) / tvolv
454            WRITE (numout,9513)
455            WRITE (numout,9514)                                                 &
456            &     (  umo(jpicpd_hpg) + umo(jpicpd_keg) + umo(jpicpd_rvo) + umo(jpicpd_pvo) + umo(jpicpd_ldf)   &
457            &      + umo(jpicpd_had) + umo(jpicpd_zad) + umo(jpicpd_zdf) + umo(jpicpd_spg) + umo(jpicpd_dat)   &
458            &      + umo(jpicpd_swf) + umo(jpicpd_bfr) ) / tvolu,   &
459            &     (  vmo(jpicpd_hpg) + vmo(jpicpd_keg) + vmo(jpicpd_rvo) + vmo(jpicpd_pvo) + vmo(jpicpd_ldf)   &
460            &      + vmo(jpicpd_had) + vmo(jpicpd_zad) + vmo(jpicpd_zdf) + vmo(jpicpd_spg) + vmo(jpicpd_dat)   &
461            &      + vmo(jpicpd_swf) + vmo(jpicpd_bfr) ) / tvolv
462         ENDIF
463
464 9500    FORMAT(' momentum trend at it= ', i6, ' :', /' ==============================')
465 9501    FORMAT(' pressure gradient          u= ', e20.13, '    v= ', e20.13)
466 9502    FORMAT(' ke gradient                u= ', e20.13, '    v= ', e20.13)
467 9503    FORMAT(' relative vorticity term    u= ', e20.13, '    v= ', e20.13)
468 9504    FORMAT(' coriolis term              u= ', e20.13, '    v= ', e20.13)
469 9505    FORMAT(' horizontal diffusion       u= ', e20.13, '    v= ', e20.13)
470 9506    FORMAT(' horizontal advection       u= ', e20.13, '    v= ', e20.13)
471 9507    FORMAT(' vertical advection         u= ', e20.13, '    v= ', e20.13)
472 9508    FORMAT(' vertical diffusion         u= ', e20.13, '    v= ', e20.13)
473 9509    FORMAT(' surface pressure gradient  u= ', e20.13, '    v= ', e20.13)
474 9510    FORMAT(' surface wind forcing       u= ', e20.13, '    v= ', e20.13)
475 9511    FORMAT(' dampimg term               u= ', e20.13, '    v= ', e20.13)
476 9512    FORMAT(' bottom flux                u= ', e20.13, '    v= ', e20.13)
477 9513    FORMAT(' -----------------------------------------------------------------------------')
478 9514    FORMAT(' total trend                u= ', e20.13, '    v= ', e20.13)
479
480         IF(lwp) THEN
481            WRITE (numout,*)
482            WRITE (numout,*)
483            WRITE (numout,9520) kt
484            WRITE (numout,9521) hke(jpicpd_hpg) / tvolt
485            WRITE (numout,9522) hke(jpicpd_keg) / tvolt
486            WRITE (numout,9523) hke(jpicpd_rvo) / tvolt
487            WRITE (numout,9524) hke(jpicpd_pvo) / tvolt
488            WRITE (numout,9525) hke(jpicpd_ldf) / tvolt
489            WRITE (numout,9526) hke(jpicpd_had) / tvolt
490            WRITE (numout,9527) hke(jpicpd_zad) / tvolt
491            WRITE (numout,9528) hke(jpicpd_zdf) / tvolt
492            WRITE (numout,9529) hke(jpicpd_spg) / tvolt
493            WRITE (numout,9530) hke(jpicpd_swf) / tvolt
494            WRITE (numout,9531) hke(jpicpd_dat) / tvolt
495            WRITE (numout,9532) hke(jpicpd_bfr) / tvolt
496            WRITE (numout,9533)
497            WRITE (numout,9534)   &
498            &     (  hke(jpicpd_hpg) + hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_pvo) + hke(jpicpd_ldf)   &
499            &      + hke(jpicpd_had) + hke(jpicpd_zad) + hke(jpicpd_zdf) + hke(jpicpd_spg) + hke(jpicpd_dat)   &
500            &      + hke(jpicpd_swf) + hke(jpicpd_bfr) ) / tvolt
501         ENDIF
502
503 9520    FORMAT(' kinetic energy trend at it= ', i6, ' :', /' ====================================')
504 9521    FORMAT(' pressure gradient         u2= ', e20.13)
505 9522    FORMAT(' ke gradient               u2= ', e20.13)
506 9523    FORMAT(' relative vorticity term   u2= ', e20.13)
507 9524    FORMAT(' coriolis term             u2= ', e20.13)
508 9525    FORMAT(' horizontal diffusion      u2= ', e20.13)
509 9526    FORMAT(' horizontal advection      u2= ', e20.13)
510 9527    FORMAT(' vertical advection        u2= ', e20.13)
511 9528    FORMAT(' vertical diffusion        u2= ', e20.13)
512 9529    FORMAT(' surface pressure gradient u2= ', e20.13)
513 9530    FORMAT(' surface wind forcing      u2= ', e20.13)
514 9531    FORMAT(' dampimg term              u2= ', e20.13)
515 9532    FORMAT(' bottom flux               u2= ', e20.13)
516 9533    FORMAT(' --------------------------------------------------')
517 9534    FORMAT(' total trend               u2= ', e20.13)
518
519         IF(lwp) THEN
520            WRITE (numout,*)
521            WRITE (numout,*)
522            WRITE (numout,9540) kt
523            WRITE (numout,9541) ( hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_had) + hke(jpicpd_zad) ) / tvolt
524            WRITE (numout,9542) ( hke(jpicpd_keg) + hke(jpicpd_had) + hke(jpicpd_zad) ) / tvolt
525            WRITE (numout,9543) ( hke(jpicpd_pvo) ) / tvolt
526            WRITE (numout,9544) ( hke(jpicpd_rvo) ) / tvolt
527            WRITE (numout,9545) ( hke(jpicpd_spg) ) / tvolt
528            WRITE (numout,9546) ( hke(jpicpd_ldf) ) / tvolt
529            WRITE (numout,9547) ( hke(jpicpd_zdf) ) / tvolt
530            WRITE (numout,9548) ( hke(jpicpd_hpg) ) / tvolt, rpktrd / tvolt
531            WRITE (numout,*)
532            WRITE (numout,*)
533         ENDIF
534
535 9540    FORMAT(' energetic consistency at it= ', i6, ' :', /' =========================================')
536 9541    FORMAT(' 0 = non linear term(true if key_vorenergy or key_combined): ', e20.13)
537 9542    FORMAT(' 0 = ke gradient + horizontal + vertical advection         : ', e20.13)
538 9543    FORMAT(' 0 = coriolis term  (true if key_vorenergy or key_combined): ', e20.13)
539 9544    FORMAT(' 0 = uh.( rot(u) x uh ) (true if enstrophy conser.)        : ', e20.13)
540 9545    FORMAT(' 0 = surface pressure gradient                             : ', e20.13)
541 9546    FORMAT(' 0 > horizontal diffusion                                  : ', e20.13)
542 9547    FORMAT(' 0 > vertical diffusion                                    : ', e20.13)
543 9548    FORMAT(' pressure gradient u2 = - 1/rau0 u.dz(rhop)                : ', e20.13, '  u.dz(rhop) =', e20.13)
544         !
545         ! Save potential to kinetic energy conversion for next time step
546         rpktrd = peke
547         !
548      ENDIF
549      !
550   END SUBROUTINE trd_dwr
551
552
553   SUBROUTINE trd_twr( kt )
554      !!---------------------------------------------------------------------
555      !!                  ***  ROUTINE trd_twr  ***
556      !!
557      !! ** Purpose :  write active tracers trends in ocean.output
558      !!----------------------------------------------------------------------
559      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
560      !!----------------------------------------------------------------------
561
562      ! I. Tracers trends
563      ! -----------------
564
565      IF( MOD(kt,nn_trd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
566
567         ! I.1 Sums over the global domain
568         ! -------------------------------
569         IF( lk_mpp ) THEN
570            CALL mpp_sum( tmo, jptot_tra )   
571            CALL mpp_sum( smo, jptot_tra )
572            CALL mpp_sum( t2 , jptot_tra )
573            CALL mpp_sum( s2 , jptot_tra )
574         ENDIF
575
576         ! I.2 Print tracers trends in the ocean.output file
577         ! --------------------------------------------------
578         
579         IF(lwp) THEN
580            WRITE (numout,*)
581            WRITE (numout,*)
582            WRITE (numout,9400) kt
583            WRITE (numout,9401)  tmo(jpicpt_xad) / tvolt, smo(jpicpt_xad) / tvolt
584            WRITE (numout,9411)  tmo(jpicpt_yad) / tvolt, smo(jpicpt_yad) / tvolt
585            WRITE (numout,9402)  tmo(jpicpt_zad) / tvolt, smo(jpicpt_zad) / tvolt
586            WRITE (numout,9403)  tmo(jpicpt_ldf) / tvolt, smo(jpicpt_ldf) / tvolt
587            WRITE (numout,9404)  tmo(jpicpt_zdf) / tvolt, smo(jpicpt_zdf) / tvolt
588            WRITE (numout,9405)  tmo(jpicpt_npc) / tvolt, smo(jpicpt_npc) / tvolt
589            WRITE (numout,9406)  tmo(jpicpt_dmp) / tvolt, smo(jpicpt_dmp) / tvolt
590            WRITE (numout,9407)  tmo(jpicpt_qsr) / tvolt
591            WRITE (numout,9408)  tmo(jpicpt_nsr) / tvolt, smo(jpicpt_nsr) / tvolt
592            WRITE (numout,9409) 
593            WRITE (numout,9410) (  tmo(jpicpt_xad) + tmo(jpicpt_yad) + tmo(jpicpt_zad) + tmo(jpicpt_ldf) + tmo(jpicpt_zdf)   &
594            &                    + tmo(jpicpt_npc) + tmo(jpicpt_dmp) + tmo(jpicpt_qsr) + tmo(jpicpt_nsr) ) / tvolt,   &
595            &                   (  smo(jpicpt_xad) + smo(jpicpt_yad) + smo(jpicpt_zad) + smo(jpicpt_ldf) + smo(jpicpt_zdf)   &
596            &                    + smo(jpicpt_npc) + smo(jpicpt_dmp)                   + smo(jpicpt_nsr) ) / tvolt
597         ENDIF
598
5999400     FORMAT(' tracer trend at it= ',i6,' :     temperature',   &
600              '              salinity',/' ============================')
6019401     FORMAT(' zonal      advection        ',e20.13,'     ',e20.13)
6029411     FORMAT(' meridional advection        ',e20.13,'     ',e20.13)
6039402     FORMAT(' vertical advection          ',e20.13,'     ',e20.13)
6049403     FORMAT(' horizontal diffusion        ',e20.13,'     ',e20.13)
6059404     FORMAT(' vertical diffusion          ',e20.13,'     ',e20.13)
6069405     FORMAT(' static instability mixing   ',e20.13,'     ',e20.13)
6079406     FORMAT(' damping term                ',e20.13,'     ',e20.13)
6089407     FORMAT(' penetrative qsr             ',e20.13)
6099408     FORMAT(' non solar radiation         ',e20.13,'     ',e20.13)
6109409     FORMAT(' -------------------------------------------------------------------------')
6119410     FORMAT(' total trend                 ',e20.13,'     ',e20.13)
612
613
614         IF(lwp) THEN
615            WRITE (numout,*)
616            WRITE (numout,*)
617            WRITE (numout,9420) kt
618            WRITE (numout,9421)   t2(jpicpt_xad) / tvolt, s2(jpicpt_xad) / tvolt
619            WRITE (numout,9431)   t2(jpicpt_yad) / tvolt, s2(jpicpt_yad) / tvolt
620            WRITE (numout,9422)   t2(jpicpt_zad) / tvolt, s2(jpicpt_zad) / tvolt
621            WRITE (numout,9423)   t2(jpicpt_ldf) / tvolt, s2(jpicpt_ldf) / tvolt
622            WRITE (numout,9424)   t2(jpicpt_zdf) / tvolt, s2(jpicpt_zdf) / tvolt
623            WRITE (numout,9425)   t2(jpicpt_npc) / tvolt, s2(jpicpt_npc) / tvolt
624            WRITE (numout,9426)   t2(jpicpt_dmp) / tvolt, s2(jpicpt_dmp) / tvolt
625            WRITE (numout,9427)   t2(jpicpt_qsr) / tvolt
626            WRITE (numout,9428)   t2(jpicpt_nsr) / tvolt, s2(jpicpt_nsr) / tvolt
627            WRITE (numout,9429)
628            WRITE (numout,9430) (  t2(jpicpt_xad) + t2(jpicpt_yad) + t2(jpicpt_zad) + t2(jpicpt_ldf) + t2(jpicpt_zdf)   &
629            &                    + t2(jpicpt_npc) + t2(jpicpt_dmp) + t2(jpicpt_qsr) + t2(jpicpt_nsr) ) / tvolt,   &
630            &                   (  s2(jpicpt_xad) + s2(jpicpt_yad) + s2(jpicpt_zad) + s2(jpicpt_ldf) + s2(jpicpt_zdf)   &
631            &                    + s2(jpicpt_npc) + s2(jpicpt_dmp)                  + s2(jpicpt_nsr) ) / tvolt
632         ENDIF
633
6349420     FORMAT(' tracer**2 trend at it= ', i6, ' :      temperature',   &
635            '               salinity', /, ' ===============================')
6369421     FORMAT(' zonal      advection      * t   ', e20.13, '     ', e20.13)
6379431     FORMAT(' meridional advection      * t   ', e20.13, '     ', e20.13)
6389422     FORMAT(' vertical advection        * t   ', e20.13, '     ', e20.13)
6399423     FORMAT(' horizontal diffusion      * t   ', e20.13, '     ', e20.13)
6409424     FORMAT(' vertical diffusion        * t   ', e20.13, '     ', e20.13)
6419425     FORMAT(' static instability mixing * t   ', e20.13, '     ', e20.13)
6429426     FORMAT(' damping term              * t   ', e20.13, '     ', e20.13)
6439427     FORMAT(' penetrative qsr           * t   ', e20.13)
6449428     FORMAT(' non solar radiation       * t   ', e20.13, '     ', e20.13)
6459429     FORMAT(' -----------------------------------------------------------------------------')
6469430     FORMAT(' total trend                *t = ', e20.13, '  *s = ', e20.13)
647
648
649         IF(lwp) THEN
650            WRITE (numout,*)
651            WRITE (numout,*)
652            WRITE (numout,9440) kt
653            WRITE (numout,9441) ( tmo(jpicpt_xad)+tmo(jpicpt_yad)+tmo(jpicpt_zad) )/tvolt,   &
654            &                   ( smo(jpicpt_xad)+smo(jpicpt_yad)+smo(jpicpt_zad) )/tvolt
655            WRITE (numout,9442)   tmo(jpicpt_zl1)/tvolt,  smo(jpicpt_zl1)/tvolt
656            WRITE (numout,9443)   tmo(jpicpt_ldf)/tvolt,  smo(jpicpt_ldf)/tvolt
657            WRITE (numout,9444)   tmo(jpicpt_zdf)/tvolt,  smo(jpicpt_zdf)/tvolt
658            WRITE (numout,9445)   tmo(jpicpt_npc)/tvolt,  smo(jpicpt_npc)/tvolt
659            WRITE (numout,9446) ( t2(jpicpt_xad)+t2(jpicpt_yad)+t2(jpicpt_zad) )/tvolt,    &
660            &                   ( s2(jpicpt_xad)+s2(jpicpt_yad)+s2(jpicpt_zad) )/tvolt
661            WRITE (numout,9447)   t2(jpicpt_ldf)/tvolt,   s2(jpicpt_ldf)/tvolt
662            WRITE (numout,9448)   t2(jpicpt_zdf)/tvolt,   s2(jpicpt_zdf)/tvolt
663            WRITE (numout,9449)   t2(jpicpt_npc)/tvolt,   s2(jpicpt_npc)/tvolt
664         ENDIF
665
6669440     FORMAT(' tracer consistency at it= ',i6,   &
667            ' :         temperature','                salinity',/,   &
668            ' ==================================')
6699441     FORMAT(' 0 = horizontal+vertical advection +    ',e20.13,'       ',e20.13)
6709442     FORMAT('     1st lev vertical advection         ',e20.13,'       ',e20.13)
6719443     FORMAT(' 0 = horizontal diffusion               ',e20.13,'       ',e20.13)
6729444     FORMAT(' 0 = vertical diffusion                 ',e20.13,'       ',e20.13)
6739445     FORMAT(' 0 = static instability mixing          ',e20.13,'       ',e20.13)
6749446     FORMAT(' 0 = horizontal+vertical advection * t  ',e20.13,'       ',e20.13)
6759447     FORMAT(' 0 > horizontal diffusion          * t  ',e20.13,'       ',e20.13)
6769448     FORMAT(' 0 > vertical diffusion            * t  ',e20.13,'       ',e20.13)
6779449     FORMAT(' 0 > static instability mixing     * t  ',e20.13,'       ',e20.13)
678         !
679      ENDIF
680      !
681   END SUBROUTINE trd_twr
682
683#   else
684   !!----------------------------------------------------------------------
685   !!   Default case :                                         Empty module
686   !!----------------------------------------------------------------------
687   INTERFACE trd_icp
688      MODULE PROCEDURE trd_2d, trd_3d
689   END INTERFACE
690
691CONTAINS
692   SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype )       ! Empty routine
693      REAL, DIMENSION(:,:) ::   ptrd2dx, ptrd2dy
694      INTEGER                     , INTENT(in   ) ::   ktrd         ! tracer trend index
695      CHARACTER(len=3)            , INTENT(in   ) ::   ctype        ! momentum ('DYN') or tracers ('TRA') trends
696      WRITE(*,*) 'trd_2d: You should not have seen this print! error ?', &
697          &       ptrd2dx(1,1), ptrd2dy(1,1), ktrd, ctype
698   END SUBROUTINE trd_2d
699   SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd , ctype )       ! Empty routine
700      REAL, DIMENSION(:,:,:) ::   ptrd3dx, ptrd3dy
701      INTEGER                     , INTENT(in   ) ::   ktrd         ! tracer trend index
702      CHARACTER(len=3)            , INTENT(in   ) ::   ctype        ! momentum ('DYN') or tracers ('TRA') trends
703      WRITE(*,*) 'trd_3d: You should not have seen this print! error ?', &
704          &       ptrd3dx(1,1,1), ptrd3dy(1,1,1), ktrd, ctype
705   END SUBROUTINE trd_3d
706   SUBROUTINE trd_icp_init               ! Empty routine
707   END SUBROUTINE trd_icp_init
708   SUBROUTINE trd_dwr( kt )          ! Empty routine
709      WRITE(*,*) 'trd_dwr: You should not have seen this print! error ?', kt
710   END SUBROUTINE trd_dwr
711   SUBROUTINE trd_twr( kt )          ! Empty routine
712      WRITE(*,*) 'trd_twr: You should not have seen this print! error ?', kt
713   END SUBROUTINE trd_twr
714#endif
715
716   !!======================================================================
717END MODULE trdicp
Note: See TracBrowser for help on using the repository browser.