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

source: trunk/NEMO/OPA_SRC/TRD/trdicp.F90 @ 881

Last change on this file since 881 was 719, checked in by ctlod, 17 years ago

get back to the nemo_v2_3 version for trunk

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 33.4 KB
RevLine 
[215]1MODULE trdicp
2   !!======================================================================
3   !!                       ***  MODULE  trdicp  ***
4   !! Ocean diagnostics:  ocean tracers and dynamic trends
5   !!=====================================================================
[503]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   !!----------------------------------------------------------------------
[215]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   !!----------------------------------------------------------------------
[503]18   !!   trd_icp          : compute the basin averaged properties for tra/dyn
[215]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
[503]37   INTERFACE trd_icp
[215]38      MODULE PROCEDURE trd_2d, trd_3d
39   END INTERFACE
40
[503]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
[215]45
46   !! * Substitutions
47#  include "domzgr_substitute.h90"
48#  include "vectopt_loop_substitute.h90"
49   !!----------------------------------------------------------------------
[247]50   !!   OPA 9.0 , LOCEAN-IPSL (2005)
[719]51   !! $Header$
[503]52   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
[215]53   !!----------------------------------------------------------------------
54
55CONTAINS
56
[503]57   SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype, clpas )
[215]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 ntrd.
[503]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      CHARACTER(len=3)            , INTENT(in   ), OPTIONAL ::   clpas     ! number of passage
[215]69      !!
[503]70      INTEGER  ::   ji, jj                                                 ! loop indices
71      CHARACTER(len=3) ::   cpas                                           ! number of passage
72      REAL(wp) ::   zmsku, zbtu, zbt                                       ! temporary scalars
73      REAL(wp) ::   zmskv, zbtv                                            !    "         "
[215]74      !!----------------------------------------------------------------------
75
[503]76      ! Control of optional arguments
77      cpas = 'fst'
78      IF( PRESENT(clpas) )  cpas = clpas
[215]79
[503]80      ! 1. Mask trends
81      ! --------------
[215]82
[503]83      SELECT CASE( ctype )
84      !
85      CASE( 'DYN' )              ! Momentum
[215]86         DO jj = 1, jpjm1
87            DO ji = 1, jpim1
88               zmsku = tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,1)
89               zmskv = tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,1)
90               ptrd2dx(ji,jj) = ptrd2dx(ji,jj) * zmsku
91               ptrd2dy(ji,jj) = ptrd2dy(ji,jj) * zmskv
92            END DO
93         END DO
94         ptrd2dx(jpi, : ) = 0.e0      ;      ptrd2dy(jpi, : ) = 0.e0
95         ptrd2dx( : ,jpj) = 0.e0      ;      ptrd2dy( : ,jpj) = 0.e0
[503]96         !
97      CASE( 'TRA' )              ! Tracers
[215]98         ptrd2dx(:,:) = ptrd2dx(:,:) * tmask_i(:,:)
99         ptrd2dy(:,:) = ptrd2dy(:,:) * tmask_i(:,:)
[503]100         !
[215]101      END SELECT
102     
[503]103      ! 2. Basin averaged tracer/momentum trends
104      ! ----------------------------------------
[215]105
[503]106      SELECT CASE( ctype )
107      !
108      CASE( 'DYN' )              ! Momentum
[215]109         umo(ktrd) = 0.e0
110         vmo(ktrd) = 0.e0
[503]111         !
112         SELECT CASE( ktrd )
113         !
114         CASE( jpdyn_trd_swf )         ! surface forcing
[215]115            DO jj = 1, jpj
116               DO ji = 1, jpi
117                  umo(ktrd) = umo(ktrd) + ptrd2dx(ji,jj) * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1)
118                  vmo(ktrd) = vmo(ktrd) + ptrd2dy(ji,jj) * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1)
119               END DO
120            END DO
[503]121            !
122         CASE( jpdyn_trd_bfr )         ! bottom friction fluxes
[215]123            DO jj = 1, jpj
124               DO ji = 1, jpi
125                  umo(ktrd) = umo(ktrd) + ptrd2dx(ji,jj)
126                  vmo(ktrd) = vmo(ktrd) + ptrd2dy(ji,jj)
127               END DO
128            END DO
[503]129            !
[215]130         END SELECT
[503]131         !
132      CASE( 'TRA' )              ! Tracers
133         IF( cpas == 'fst' )   THEN
134            tmo(ktrd) = 0.e0
135            smo(ktrd) = 0.e0
136         ENDIF
[215]137         DO jj = 1, jpj
138            DO ji = 1, jpi
139               zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
140               tmo(ktrd) =  tmo(ktrd) + ptrd2dx(ji,jj) * zbt
141               smo(ktrd) =  smo(ktrd) + ptrd2dy(ji,jj) * zbt
142            END DO
143         END DO
[503]144         !
[215]145      END SELECT
146     
[503]147      ! 3. Basin averaged tracer/momentum square trends
148      ! ----------------------------------------------
[215]149      ! c a u t i o n: field now
150     
[503]151      SELECT CASE( ctype )
152      !
153      CASE( 'DYN' )              ! Momentum
[215]154         hke(ktrd) = 0.e0
155         DO jj = 1, jpj
156            DO ji = 1, jpi
157               zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1)
158               zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1)
159               hke(ktrd) = hke(ktrd)   &
160               &   + un(ji,jj,1) * ptrd2dx(ji,jj) * zbtu &
161               &   + vn(ji,jj,1) * ptrd2dy(ji,jj) * zbtv
162            END DO
163         END DO
[503]164         !
165      CASE( 'TRA' )              ! Tracers
166         IF( cpas == 'fst' )   THEN
167            t2(ktrd) = 0.e0
168            s2(ktrd) = 0.e0
169         ENDIF
[215]170         DO jj = 1, jpj
171            DO ji = 1, jpi
172               zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
173               t2(ktrd) = t2(ktrd) + ptrd2dx(ji,jj) * zbt * tn(ji,jj,1)
174               s2(ktrd) = s2(ktrd) + ptrd2dy(ji,jj) * zbt * sn(ji,jj,1)
175            END DO
176         END DO
[503]177         !     
[215]178      END SELECT
[503]179      !
[215]180   END SUBROUTINE trd_2d
181
182
[503]183   SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd, ctype, clpas )
[215]184      !!---------------------------------------------------------------------
185      !!                  ***  ROUTINE trd_3d  ***
186      !!
187      !! ** Purpose : verify the basin averaged properties of tracers and/or
188      !!              momentum equations at every time step frequency ntrd.
[503]189      !!----------------------------------------------------------------------
190      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptrd3dx            ! Temperature or U trend
191      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptrd3dy            ! Salinity    or V trend
192      INTEGER,                          INTENT(in   ) ::   ktrd               ! momentum or tracer trend index
193      CHARACTER(len=3),                 INTENT(in   ) ::   ctype              ! momentum ('DYN') or tracers ('TRA') trends
194      CHARACTER(len=3),                 INTENT(in   ), OPTIONAL ::   clpas    ! number of passage
[215]195      !!
196      INTEGER ::   ji, jj, jk
[503]197      CHARACTER(len=3) ::   cpas                                              ! number of passage
198      REAL(wp) ::   zbt, zbtu, zbtv, zmsku, zmskv                             ! temporary scalars
[215]199      !!----------------------------------------------------------------------
200
[503]201      ! Control of optional arguments
202      cpas = 'fst'
203      IF( PRESENT(clpas) )  cpas = clpas
[215]204
[503]205      ! 1. Mask the trends
206      ! ------------------
[215]207
[503]208      SELECT CASE( ctype )
209      !
210      CASE( 'DYN' )              ! Momentum       
[215]211         DO jk = 1, jpk
212            DO jj = 1, jpjm1
213               DO ji = 1, jpim1
214                  zmsku = tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)
215                  zmskv = tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
216                  ptrd3dx(ji,jj,jk) = ptrd3dx(ji,jj,jk) * zmsku
217                  ptrd3dy(ji,jj,jk) = ptrd3dy(ji,jj,jk) * zmskv
[503]218               END DO
219            END DO
220         END DO
[215]221         ptrd3dx(jpi, : ,:) = 0.e0      ;      ptrd3dy(jpi, : ,:) = 0.e0
222         ptrd3dx( : ,jpj,:) = 0.e0      ;      ptrd3dy( : ,jpj,:) = 0.e0
[503]223         !
224      CASE( 'TRA' )              ! Tracers
[215]225         DO jk = 1, jpk
226            ptrd3dx(:,:,jk) = ptrd3dx(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
227            ptrd3dy(:,:,jk) = ptrd3dy(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
[503]228         END DO
229         !
[215]230      END SELECT   
231
[503]232      ! 2. Basin averaged tracer/momentum trends
233      ! ----------------------------------------
[215]234     
[503]235      SELECT CASE( ctype )
236      !
237      CASE( 'DYN' )              ! Momentum
[215]238         umo(ktrd) = 0.e0
239         vmo(ktrd) = 0.e0
240         DO jk = 1, jpk
241            DO jj = 1, jpj
242               DO ji = 1, jpi
243                  zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
244                  zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
245                  umo(ktrd) = umo(ktrd) + ptrd3dx(ji,jj,jk) * zbtu
246                  vmo(ktrd) = vmo(ktrd) + ptrd3dy(ji,jj,jk) * zbtv
247               END DO
248            END DO
249         END DO
[503]250         !
251      CASE( 'TRA' )              ! Tracers
252         IF( cpas == 'fst' )   THEN
253            tmo(ktrd) = 0.e0
254            smo(ktrd) = 0.e0
255         ENDIF
[215]256         DO jk = 1, jpkm1
257            DO jj = 1, jpj
258               DO ji = 1, jpi
259                  zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 
260                  tmo(ktrd) =  tmo(ktrd) + ptrd3dx(ji,jj,jk) * zbt
261                  smo(ktrd) =  smo(ktrd) + ptrd3dy(ji,jj,jk) * zbt
262               END DO
263            END DO
264         END DO
[503]265         !
[215]266      END SELECT
267
[503]268      ! 3. Basin averaged tracer/momentum square trends
269      ! -----------------------------------------------
[215]270      ! c a u t i o n: field now
271     
[503]272      SELECT CASE( ctype )
273      !
274      CASE( 'DYN' )              ! Momentum
[215]275         hke(ktrd) = 0.e0
276         DO jk = 1, jpk
277            DO jj = 1, jpj
278               DO ji = 1, jpi
279                  zbtu = e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
280                  zbtv = e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
281                  hke(ktrd) = hke(ktrd)   &
282                  &   + un(ji,jj,jk) * ptrd3dx(ji,jj,jk) * zbtu &
283                  &   + vn(ji,jj,jk) * ptrd3dy(ji,jj,jk) * zbtv
284               END DO
285            END DO
286         END DO
[503]287         !
288      CASE( 'TRA' )              ! Tracers
289         IF( cpas == 'fst' )   THEN
290            t2(ktrd) = 0.e0
291            s2(ktrd) = 0.e0
292         ENDIF
[215]293         DO jk = 1, jpk
294            DO jj = 1, jpj
295               DO ji = 1, jpi
296                  zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
297                  t2(ktrd) = t2(ktrd) + ptrd3dx(ji,jj,jk) * zbt * tn(ji,jj,jk)
298                  s2(ktrd) = s2(ktrd) + ptrd3dy(ji,jj,jk) * zbt * sn(ji,jj,jk)
299               END DO
300            END DO
301         END DO
[503]302         !
[215]303      END SELECT
[503]304      !
[215]305   END SUBROUTINE trd_3d
306
307
308
309   SUBROUTINE trd_icp_init
310      !!---------------------------------------------------------------------
311      !!                  ***  ROUTINE trd_icp_init  ***
312      !!
[503]313      !! ** Purpose :   Read the namtrd namelist
[215]314      !!----------------------------------------------------------------------
[503]315      INTEGER  ::   ji, jj, jk
[215]316      REAL(wp) ::   zmskt
317#if  defined key_trddyn
[503]318      REAL(wp) ::   zmsku, zmskv
[215]319#endif
320      !!----------------------------------------------------------------------
321
322      IF(lwp) THEN
323         WRITE(numout,*)
324         WRITE(numout,*) 'trd_icp_init : integral constraints properties trends'
325         WRITE(numout,*) '~~~~~~~~~~~~~'
326      ENDIF
327
328      ! Total volume at t-points:
329      tvolt = 0.e0
330      DO jk = 1, jpkm1
331         DO jj = 2, jpjm1
332            DO ji = fs_2, fs_jpim1   ! vector opt.
333               zmskt = tmask(ji,jj,jk) * tmask_i(ji,jj)
334               tvolt = tvolt + zmskt * e1t(ji,jj) *e2t(ji,jj) * fse3t(ji,jj,jk)
335            END DO
336         END DO
337      END DO
338      IF( lk_mpp )   CALL mpp_sum( tvolt )   ! sum over the global domain
339
[503]340      IF(lwp) WRITE(numout,*) '                total ocean volume at T-point   tvolt = ',tvolt
[215]341
342#if  defined key_trddyn
343      ! Initialization of potential to kinetic energy conversion
344      rpktrd = 0.e0
345
346      ! Total volume at u-, v- points:
347      tvolu = 0.e0
348      tvolv = 0.e0
349
350      DO jk = 1, jpk
351         DO jj = 2, jpjm1
352            DO ji = fs_2, fs_jpim1   ! vector opt.
353               zmsku = tmask_i(ji+1,jj  ) * tmask_i(ji,jj) * umask(ji,jj,jk)
354               zmskv = tmask_i(ji  ,jj+1) * tmask_i(ji,jj) * vmask(ji,jj,jk)
355               tvolu = tvolu + zmsku * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk)
356               tvolv = tvolv + zmskv * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk)
357            END DO
358         END DO
359      END DO
360      IF( lk_mpp )   CALL mpp_sum( tvolu )   ! sums over the global domain
361      IF( lk_mpp )   CALL mpp_sum( tvolv )
362
363      IF(lwp) THEN
[503]364         WRITE(numout,*) '                total ocean volume at U-point   tvolu = ',tvolu
365         WRITE(numout,*) '                total ocean volume at V-point   tvolv = ',tvolv
[215]366      ENDIF
367#endif
[503]368      !
[215]369   END SUBROUTINE trd_icp_init
370
371
372   SUBROUTINE trd_dwr( kt )
373      !!---------------------------------------------------------------------
374      !!                  ***  ROUTINE trd_dwr  ***
375      !!
376      !! ** Purpose :  write dynamic trends in ocean.output
[503]377      !!----------------------------------------------------------------------
378      INTEGER, INTENT(in) ::   kt                                  ! ocean time-step index
[215]379      !!
[503]380      INTEGER  ::   ji, jj, jk
381      REAL(wp) ::   ze1e2w, zcof, zbe1ru, zbe2rv, zbtr, ztz, zth   !    "      scalars
382      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zkepe, zkx, zky, zkz   ! temporary arrays
[215]383      !!----------------------------------------------------------------------
384
385      ! I. Momentum trends
386      ! -------------------
387
388      IF( MOD(kt,ntrd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
389
390         ! I.1 Conversion potential energy - kinetic energy
391         ! --------------------------------------------------
392         ! c a u t i o n here, trends are computed at kt+1 (now , but after the swap)
393
[503]394         zkx(:,:,:)   = 0.e0
395         zky(:,:,:)   = 0.e0
396         zkz(:,:,:)   = 0.e0
[215]397         zkepe(:,:,:) = 0.e0
398   
399         CALL eos( tn, sn, rhd, rhop )       ! now potential and in situ densities
400
[503]401         ! Density flux at w-point
[215]402         DO jk = 2, jpk
403            DO jj = 1, jpj
404               DO ji = 1, jpi
405                  ze1e2w = 0.5 * e1t(ji,jj) * e2t(ji,jj) * wn(ji,jj,jk) * tmask_i(ji,jj)
406                  zkz(ji,jj,jk) = ze1e2w / rau0 * ( rhop(ji,jj,jk) + rhop(ji,jj,jk-1) )
407               END DO
408            END DO
409         END DO
[503]410         zkz(:,:,1) = 0.e0
[215]411         
412         ! Density flux at u and v-points
413         DO jk = 1, jpk
414            DO jj = 1, jpjm1
415               DO ji = 1, jpim1
416                  zcof   = 0.5 / rau0
417                  zbe1ru = zcof * e2u(ji,jj) * fse3u(ji,jj,jk) * un(ji,jj,jk)
418                  zbe2rv = zcof * e1v(ji,jj) * fse3v(ji,jj,jk) * vn(ji,jj,jk)
419                  zkx(ji,jj,jk) = zbe1ru * ( rhop(ji,jj,jk) + rhop(ji+1,jj,jk) )
420                  zky(ji,jj,jk) = zbe2rv * ( rhop(ji,jj,jk) + rhop(ji,jj+1,jk) )
421               END DO
422            END DO
423         END DO
424         
425         ! Density flux divergence at t-point
426         DO jk = 1, jpkm1
427            DO jj = 2, jpjm1
428               DO ji = 2, jpim1
429                  zbtr = 1. / ( e1t(ji,jj)*e2t(ji,jj)*fse3t(ji,jj,jk) )
430                  ztz = - zbtr * (    zkz(ji,jj,jk) - zkz(ji,jj,jk+1) )
431                  zth = - zbtr * (  ( zkx(ji,jj,jk) - zkx(ji-1,jj,jk) )   &
432                    &             + ( zky(ji,jj,jk) - zky(ji,jj-1,jk) )  )
433                  zkepe(ji,jj,jk) = (zth + ztz) * tmask(ji,jj,jk) * tmask_i(ji,jj)
434               END DO
435            END DO
436         END DO
437         zkepe( : , : ,jpk) = 0.e0
438         zkepe( : ,jpj, : ) = 0.e0
439         zkepe(jpi, : , : ) = 0.e0
440
441         ! I.2 Basin averaged kinetic energy trend
442         ! ----------------------------------------
443         peke = 0.e0
444         DO jk = 1,jpk
445            DO jj = 1, jpj
446               DO ji = 1, jpi
447                  peke = peke + zkepe(ji,jj,jk) * grav * fsdept(ji,jj,jk)   &
448                     &                     * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
449               END DO
450            END DO
451         END DO
452
453         ! I.3 Sums over the global domain
454         ! ---------------------------------
455         IF( lk_mpp ) THEN
[503]456            CALL mpp_sum( peke )
457            CALL mpp_sum( umo , jptot_dyn )
458            CALL mpp_sum( vmo , jptot_dyn )
459            CALL mpp_sum( hke , jptot_dyn )
460         ENDIF
[215]461
462         ! I.2 Print dynamic trends in the ocean.output file
463         ! --------------------------------------------------
464
465         IF(lwp) THEN
466            WRITE (numout,*)
467            WRITE (numout,*)
468            WRITE (numout,9500) kt
[503]469            WRITE (numout,9501) umo(jpicpd_hpg) / tvolu, vmo(jpicpd_hpg) / tvolv
470            WRITE (numout,9502) umo(jpicpd_keg) / tvolu, vmo(jpicpd_keg) / tvolv
471            WRITE (numout,9503) umo(jpicpd_rvo) / tvolu, vmo(jpicpd_rvo) / tvolv
472            WRITE (numout,9504) umo(jpicpd_pvo) / tvolu, vmo(jpicpd_pvo) / tvolv
473            WRITE (numout,9505) umo(jpicpd_ldf) / tvolu, vmo(jpicpd_ldf) / tvolv
474            WRITE (numout,9506) umo(jpicpd_zad) / tvolu, vmo(jpicpd_zad) / tvolv
475            WRITE (numout,9507) umo(jpicpd_zdf) / tvolu, vmo(jpicpd_zdf) / tvolv
476            WRITE (numout,9508) umo(jpicpd_spg) / tvolu, vmo(jpicpd_spg) / tvolv
477            WRITE (numout,9509) umo(jpicpd_swf) / tvolu, vmo(jpicpd_swf) / tvolv
478            WRITE (numout,9510) umo(jpicpd_dat) / tvolu, vmo(jpicpd_dat) / tvolv
479            WRITE (numout,9511) umo(jpicpd_bfr) / tvolu, vmo(jpicpd_bfr) / tvolv
[215]480            WRITE (numout,9512)
481            WRITE (numout,9513)                                                 &
[503]482            &     (  umo(jpicpd_hpg) + umo(jpicpd_keg) + umo(jpicpd_rvo) + umo(jpicpd_pvo) + umo(jpicpd_ldf)   &
483            &      + umo(jpicpd_zad) + umo(jpicpd_zdf) + umo(jpicpd_spg) + umo(jpicpd_dat) + umo(jpicpd_swf)   &
484            &      + umo(jpicpd_bfr) ) / tvolu,   &
485            &     (  vmo(jpicpd_hpg) + vmo(jpicpd_keg) + vmo(jpicpd_rvo) + vmo(jpicpd_pvo) + vmo(jpicpd_ldf)   &
486            &      + vmo(jpicpd_zad) + vmo(jpicpd_zdf) + vmo(jpicpd_spg) + vmo(jpicpd_dat) + vmo(jpicpd_swf)   &
487            &      + vmo(jpicpd_bfr) ) / tvolv
[215]488         ENDIF
489
490 9500    FORMAT(' momentum trend at it= ', i6, ' :', /' ==============================')
491 9501    FORMAT(' pressure gradient          u= ', e20.13, '    v= ', e20.13)
492 9502    FORMAT(' ke gradient                u= ', e20.13, '    v= ', e20.13)
493 9503    FORMAT(' relative vorticity term    u= ', e20.13, '    v= ', e20.13)
494 9504    FORMAT(' coriolis term              u= ', e20.13, '    v= ', e20.13)
495 9505    FORMAT(' horizontal diffusion       u= ', e20.13, '    v= ', e20.13)
496 9506    FORMAT(' vertical advection         u= ', e20.13, '    v= ', e20.13)
497 9507    FORMAT(' vertical diffusion         u= ', e20.13, '    v= ', e20.13)
498 9508    FORMAT(' surface pressure gradient  u= ', e20.13, '    v= ', e20.13)
[503]499 9509    FORMAT(' surface wind forcing       u= ', e20.13, '    v= ', e20.13)
[215]500 9510    FORMAT(' dampimg term               u= ', e20.13, '    v= ', e20.13)
501 9511    FORMAT(' bottom flux                u= ', e20.13, '    v= ', e20.13)
502 9512    FORMAT(' -----------------------------------------------------------------------------')
503 9513    FORMAT(' total trend                u= ', e20.13, '    v= ', e20.13)
504
505         IF(lwp) THEN
506            WRITE (numout,*)
507            WRITE (numout,*)
508            WRITE (numout,9520) kt
[503]509            WRITE (numout,9521) hke(jpicpd_hpg) / tvolt
510            WRITE (numout,9522) hke(jpicpd_keg) / tvolt
511            WRITE (numout,9523) hke(jpicpd_rvo) / tvolt
512            WRITE (numout,9524) hke(jpicpd_pvo) / tvolt
513            WRITE (numout,9525) hke(jpicpd_ldf) / tvolt
514            WRITE (numout,9526) hke(jpicpd_zad) / tvolt
515            WRITE (numout,9527) hke(jpicpd_zdf) / tvolt
516            WRITE (numout,9528) hke(jpicpd_spg) / tvolt
517            WRITE (numout,9529) hke(jpicpd_swf) / tvolt
518            WRITE (numout,9530) hke(jpicpd_dat) / tvolt
[215]519            WRITE (numout,9531)
520            WRITE (numout,9532)   &
[503]521            &     (  hke(jpicpd_hpg) + hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_pvo) + hke(jpicpd_ldf)   &
522            &      + hke(jpicpd_zad) + hke(jpicpd_zdf) + hke(jpicpd_spg) + hke(jpicpd_dat) + hke(jpicpd_swf) ) / tvolt
[215]523         ENDIF
524
525 9520    FORMAT(' kinetic energy trend at it= ', i6, ' :', /' ====================================')
526 9521    FORMAT(' pressure gradient         u2= ', e20.13)
527 9522    FORMAT(' ke gradient               u2= ', e20.13)
528 9523    FORMAT(' relative vorticity term   u2= ', e20.13)
529 9524    FORMAT(' coriolis term             u2= ', e20.13)
530 9525    FORMAT(' horizontal diffusion      u2= ', e20.13)
531 9526    FORMAT(' vertical advection        u2= ', e20.13)
532 9527    FORMAT(' vertical diffusion        u2= ', e20.13)
533 9528    FORMAT(' surface pressure gradient u2= ', e20.13)
[503]534 9529    FORMAT(' surface wind forcing      u2= ', e20.13)
[215]535 9530    FORMAT(' dampimg term              u2= ', e20.13)
536 9531    FORMAT(' --------------------------------------------------')
537 9532    FORMAT(' total trend               u2= ', e20.13)
538
539         IF(lwp) THEN
540            WRITE (numout,*)
541            WRITE (numout,*)
542            WRITE (numout,9540) kt
[503]543            WRITE (numout,9541) ( hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_zad) ) / tvolt
544            WRITE (numout,9542) ( hke(jpicpd_keg) + hke(jpicpd_zad) ) / tvolt
545            WRITE (numout,9543) ( hke(jpicpd_pvo) ) / tvolt
546            WRITE (numout,9544) ( hke(jpicpd_rvo) ) / tvolt
547            WRITE (numout,9545) ( hke(jpicpd_spg) ) / tvolt
548            WRITE (numout,9546) ( hke(jpicpd_ldf) ) / tvolt
549            WRITE (numout,9547) ( hke(jpicpd_zdf) ) / tvolt
550            WRITE (numout,9548) ( hke(jpicpd_hpg) ) / tvolt, rpktrd / tvolt
551            WRITE (numout,*)
552            WRITE (numout,*)
[215]553         ENDIF
554
555 9540    FORMAT(' energetic consistency at it= ', i6, ' :', /' =========================================')
556 9541    FORMAT(' 0 = non linear term(true if key_vorenergy or key_combined): ', e20.13)
[503]557 9542    FORMAT(' 0 = ke gradient + vertical advection                      : ', e20.13)
[215]558 9543    FORMAT(' 0 = coriolis term  (true if key_vorenergy or key_combined): ', e20.13)
[503]559 9544    FORMAT(' 0 = uh.( rot(u) x uh ) (true if enstrophy conser.)        : ', e20.13)
560 9545    FORMAT(' 0 = surface pressure gradient                             : ', e20.13)
561 9546    FORMAT(' 0 > horizontal diffusion                                  : ', e20.13)
562 9547    FORMAT(' 0 > vertical diffusion                                    : ', e20.13)
563 9548    FORMAT(' pressure gradient u2 = - 1/rau0 u.dz(rhop)                : ', e20.13, '  u.dz(rhop) =', e20.13)
564         !
[215]565         ! Save potential to kinetic energy conversion for next time step
566         rpktrd = peke
[503]567         !
[215]568      ENDIF
[503]569      !
[215]570   END SUBROUTINE trd_dwr
571
572
573   SUBROUTINE trd_twr( kt )
574      !!---------------------------------------------------------------------
575      !!                  ***  ROUTINE trd_twr  ***
576      !!
577      !! ** Purpose :  write active tracers trends in ocean.output
578      !!----------------------------------------------------------------------
[503]579      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
[215]580      !!----------------------------------------------------------------------
581
582      ! I. Tracers trends
583      ! -----------------
584
585      IF( MOD(kt,ntrd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
586
587         ! I.1 Sums over the global domain
588         ! -------------------------------
589         IF( lk_mpp ) THEN
[503]590            CALL mpp_sum( tmo, jptot_tra )   
591            CALL mpp_sum( smo, jptot_tra )
592            CALL mpp_sum( t2 , jptot_tra )
593            CALL mpp_sum( s2 , jptot_tra )
[215]594         ENDIF
595
596         ! I.2 Print tracers trends in the ocean.output file
597         ! --------------------------------------------------
598         
599         IF(lwp) THEN
600            WRITE (numout,*)
601            WRITE (numout,*)
602            WRITE (numout,9400) kt
[503]603            WRITE (numout,9401) (tmo(jpicpt_xad)+tmo(jpicpt_yad))/ tvolt, (smo(jpicpt_xad)+smo(jpicpt_yad))/ tvolt
604            WRITE (numout,9402)  tmo(jpicpt_zad) / tvolt, smo(jpicpt_zad) / tvolt
605            WRITE (numout,9403)  tmo(jpicpt_ldf) / tvolt, smo(jpicpt_ldf) / tvolt
606            WRITE (numout,9404)  tmo(jpicpt_zdf) / tvolt, smo(jpicpt_zdf) / tvolt
607            WRITE (numout,9405)  tmo(jpicpt_npc) / tvolt, smo(jpicpt_npc) / tvolt
608            WRITE (numout,9406)  tmo(jpicpt_dmp) / tvolt, smo(jpicpt_dmp) / tvolt
609            WRITE (numout,9407)  tmo(jpicpt_qsr) / tvolt
610            WRITE (numout,9408)  tmo(jpicpt_nsr) / tvolt, smo(jpicpt_nsr) / tvolt
611            WRITE (numout,9409) 
612            WRITE (numout,9410) (  tmo(jpicpt_xad) + tmo(jpicpt_yad) + tmo(jpicpt_zad) + tmo(jpicpt_ldf) + tmo(jpicpt_zdf)   &
613            &                    + tmo(jpicpt_npc) + tmo(jpicpt_dmp) + tmo(jpicpt_qsr) + tmo(jpicpt_nsr) ) / tvolt,   &
614            &                   (  smo(jpicpt_xad) + smo(jpicpt_yad) + smo(jpicpt_zad) + smo(jpicpt_ldf) + smo(jpicpt_zdf)   &
615            &                    + smo(jpicpt_npc) + smo(jpicpt_dmp)                   + smo(jpicpt_nsr) ) / tvolt
[215]616         ENDIF
617
6189400     FORMAT(' tracer trend at it= ',i6,' :     temperature',   &
619              '              salinity',/' ============================')
6209401     FORMAT(' horizontal advection        ',e20.13,'     ',e20.13)
6219402     FORMAT(' vertical advection          ',e20.13,'     ',e20.13)
6229403     FORMAT(' horizontal diffusion        ',e20.13,'     ',e20.13)
6239404     FORMAT(' vertical diffusion          ',e20.13,'     ',e20.13)
[503]6249405     FORMAT(' static instability mixing   ',e20.13,'     ',e20.13)
[215]6259406     FORMAT(' damping term                ',e20.13,'     ',e20.13)
[503]6269407     FORMAT(' penetrative qsr             ',e20.13)
6279408     FORMAT(' non solar radiation         ',e20.13,'     ',e20.13)
[215]6289409     FORMAT(' -------------------------------------------------------------------------')
6299410     FORMAT(' total trend                 ',e20.13,'     ',e20.13)
630
631
632         IF(lwp) THEN
633            WRITE (numout,*)
634            WRITE (numout,*)
635            WRITE (numout,9420) kt
[503]636            WRITE (numout,9421) ( t2(jpicpt_xad)+t2(jpicpt_yad) )/ tvolt, ( s2(jpicpt_xad)+s2(jpicpt_yad) )/ tvolt
637            WRITE (numout,9422)   t2(jpicpt_zad) / tvolt, s2(jpicpt_zad) / tvolt
638            WRITE (numout,9423)   t2(jpicpt_ldf) / tvolt, s2(jpicpt_ldf) / tvolt
639            WRITE (numout,9424)   t2(jpicpt_zdf) / tvolt, s2(jpicpt_zdf) / tvolt
640            WRITE (numout,9425)   t2(jpicpt_npc) / tvolt, s2(jpicpt_npc) / tvolt
641            WRITE (numout,9426)   t2(jpicpt_dmp) / tvolt, s2(jpicpt_dmp) / tvolt
642            WRITE (numout,9427)   t2(jpicpt_qsr) / tvolt
643            WRITE (numout,9428)   t2(jpicpt_nsr) / tvolt, s2(jpicpt_nsr) / tvolt
[215]644            WRITE (numout,9429)
[503]645            WRITE (numout,9430) (  t2(jpicpt_xad) + t2(jpicpt_yad) + t2(jpicpt_zad) + t2(jpicpt_ldf) + t2(jpicpt_zdf)   &
646            &                    + t2(jpicpt_npc) + t2(jpicpt_dmp) + t2(jpicpt_qsr) + t2(jpicpt_nsr) ) / tvolt,   &
647            &                   (  s2(jpicpt_xad) + s2(jpicpt_yad) + s2(jpicpt_zad) + s2(jpicpt_ldf) + s2(jpicpt_zdf)   &
648            &                    + s2(jpicpt_npc) + s2(jpicpt_dmp)                  + s2(jpicpt_nsr) ) / tvolt
[215]649         ENDIF
650
6519420     FORMAT(' tracer**2 trend at it= ', i6, ' :      temperature',   &
652            '               salinity', /, ' ===============================')
6539421     FORMAT(' horizontal advection      * t   ', e20.13, '     ', e20.13)
6549422     FORMAT(' vertical advection        * t   ', e20.13, '     ', e20.13)
6559423     FORMAT(' horizontal diffusion      * t   ', e20.13, '     ', e20.13)
6569424     FORMAT(' vertical diffusion        * t   ', e20.13, '     ', e20.13)
[503]6579425     FORMAT(' static instability mixing * t   ', e20.13, '     ', e20.13)
[215]6589426     FORMAT(' damping term              * t   ', e20.13, '     ', e20.13)
[503]6599427     FORMAT(' penetrative qsr           * t   ', e20.13)
6609428     FORMAT(' non solar radiation       * t   ', e20.13, '     ', e20.13)
[215]6619429     FORMAT(' -----------------------------------------------------------------------------')
6629430     FORMAT(' total trend                *t = ', e20.13, '  *s = ', e20.13)
663
664
665         IF(lwp) THEN
666            WRITE (numout,*)
667            WRITE (numout,*)
668            WRITE (numout,9440) kt
[503]669            WRITE (numout,9441) ( tmo(jpicpt_xad)+tmo(jpicpt_yad)+tmo(jpicpt_zad) )/tvolt,   &
670            &                   ( smo(jpicpt_xad)+smo(jpicpt_yad)+smo(jpicpt_zad) )/tvolt
671            WRITE (numout,9442)   tmo(jpicpt_zl1)/tvolt,  smo(jpicpt_zl1)/tvolt
672            WRITE (numout,9443)   tmo(jpicpt_ldf)/tvolt,  smo(jpicpt_ldf)/tvolt
673            WRITE (numout,9444)   tmo(jpicpt_zdf)/tvolt,  smo(jpicpt_zdf)/tvolt
674            WRITE (numout,9445)   tmo(jpicpt_npc)/tvolt,  smo(jpicpt_npc)/tvolt
675            WRITE (numout,9446) ( t2(jpicpt_xad)+t2(jpicpt_yad)+t2(jpicpt_zad) )/tvolt,    &
676            &                   ( s2(jpicpt_xad)+s2(jpicpt_yad)+s2(jpicpt_zad) )/tvolt
677            WRITE (numout,9447)   t2(jpicpt_ldf)/tvolt,   s2(jpicpt_ldf)/tvolt
678            WRITE (numout,9448)   t2(jpicpt_zdf)/tvolt,   s2(jpicpt_zdf)/tvolt
679            WRITE (numout,9449)   t2(jpicpt_npc)/tvolt,   s2(jpicpt_npc)/tvolt
[215]680         ENDIF
681
6829440     FORMAT(' tracer consistency at it= ',i6,   &
683            ' :         temperature','                salinity',/,   &
684            ' ==================================')
[503]6859441     FORMAT(' 0 = horizontal+vertical advection +    ',e20.13,'       ',e20.13)
6869442     FORMAT('     1st lev vertical advection         ',e20.13,'       ',e20.13)
6879443     FORMAT(' 0 = horizontal diffusion               ',e20.13,'       ',e20.13)
6889444     FORMAT(' 0 = vertical diffusion                 ',e20.13,'       ',e20.13)
6899445     FORMAT(' 0 = static instability mixing          ',e20.13,'       ',e20.13)
6909446     FORMAT(' 0 = horizontal+vertical advection * t  ',e20.13,'       ',e20.13)
6919447     FORMAT(' 0 > horizontal diffusion          * t  ',e20.13,'       ',e20.13)
6929448     FORMAT(' 0 > vertical diffusion            * t  ',e20.13,'       ',e20.13)
6939449     FORMAT(' 0 > static instability mixing     * t  ',e20.13,'       ',e20.13)
694         !
[215]695      ENDIF
[503]696      !
[215]697   END SUBROUTINE trd_twr
698
699#   else
700   !!----------------------------------------------------------------------
701   !!   Default case :                                         Empty module
702   !!----------------------------------------------------------------------
[601]703   INTERFACE trd_icp
704      MODULE PROCEDURE trd_2d, trd_3d
705   END INTERFACE
706
[215]707CONTAINS
[521]708   SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype, clpas )       ! Empty routine
[503]709      REAL, DIMENSION(:,:) ::   ptrd2dx, ptrd2dy
[601]710      INTEGER                     , INTENT(in   ) ::   ktrd         ! tracer trend index
711      CHARACTER(len=3)            , INTENT(in   ) ::   ctype        ! momentum ('DYN') or tracers ('TRA') trends
712      CHARACTER(len=3), INTENT(in), OPTIONAL ::   clpas             ! number of passage
[521]713      WRITE(*,*) 'trd_2d: You should not have seen this print! error ?', &
714          &       ptrd2dx(1,1), ptrd2dy(1,1), ktrd, ctype, clpas
[215]715   END SUBROUTINE trd_2d
[521]716   SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd , ctype, clpas )       ! Empty routine
[503]717      REAL, DIMENSION(:,:,:) ::   ptrd3dx, ptrd3dy
[601]718      INTEGER                     , INTENT(in   ) ::   ktrd         ! tracer trend index
719      CHARACTER(len=3)            , INTENT(in   ) ::   ctype        ! momentum ('DYN') or tracers ('TRA') trends
720      CHARACTER(len=3), INTENT(in), OPTIONAL ::   clpas             ! number of passage
[521]721      WRITE(*,*) 'trd_3d: You should not have seen this print! error ?', &
722          &       ptrd3dx(1,1,1), ptrd3dy(1,1,1), ktrd, ctype, clpas
[215]723   END SUBROUTINE trd_3d
724   SUBROUTINE trd_icp_init               ! Empty routine
725   END SUBROUTINE trd_icp_init
726   SUBROUTINE trd_dwr( kt )          ! Empty routine
727      WRITE(*,*) 'trd_dwr: You should not have seen this print! error ?', kt
728   END SUBROUTINE trd_dwr
729   SUBROUTINE trd_twr( kt )          ! Empty routine
730      WRITE(*,*) 'trd_twr: You should not have seen this print! error ?', kt
731   END SUBROUTINE trd_twr
732#endif
733
734   !!======================================================================
735END MODULE trdicp
Note: See TracBrowser for help on using the repository browser.