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

source: branches/dev_004_VVL/NEMO/OPA_SRC/TRD/trdicp.F90 @ 1388

Last change on this file since 1388 was 1152, checked in by rblod, 16 years ago

Convert cvs header to svn Id, step II
  • Property svn:eol-style set to native
  • Property svn:keywords set to Id
File size: 33.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, clpas )
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.
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
69      !!
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                                            !    "         "
74      !!----------------------------------------------------------------------
75
76      ! Control of optional arguments
77      cpas = 'fst'
78      IF( PRESENT(clpas) )  cpas = clpas
79
80      ! 1. Mask trends
81      ! --------------
82
83      SELECT CASE( ctype )
84      !
85      CASE( 'DYN' )              ! Momentum
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
96         !
97      CASE( 'TRA' )              ! Tracers
98         ptrd2dx(:,:) = ptrd2dx(:,:) * tmask_i(:,:)
99         ptrd2dy(:,:) = ptrd2dy(:,:) * tmask_i(:,:)
100         !
101      END SELECT
102     
103      ! 2. Basin averaged tracer/momentum trends
104      ! ----------------------------------------
105
106      SELECT CASE( ctype )
107      !
108      CASE( 'DYN' )              ! Momentum
109         umo(ktrd) = 0.e0
110         vmo(ktrd) = 0.e0
111         !
112         SELECT CASE( ktrd )
113         !
114         CASE( jpdyn_trd_swf )         ! surface forcing
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
121            !
122         CASE( jpdyn_trd_bfr )         ! bottom friction fluxes
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
129            !
130         END SELECT
131         !
132      CASE( 'TRA' )              ! Tracers
133         IF( cpas == 'fst' )   THEN
134            tmo(ktrd) = 0.e0
135            smo(ktrd) = 0.e0
136         ENDIF
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
144         !
145      END SELECT
146     
147      ! 3. Basin averaged tracer/momentum square trends
148      ! ----------------------------------------------
149      ! c a u t i o n: field now
150     
151      SELECT CASE( ctype )
152      !
153      CASE( 'DYN' )              ! Momentum
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
164         !
165      CASE( 'TRA' )              ! Tracers
166         IF( cpas == 'fst' )   THEN
167            t2(ktrd) = 0.e0
168            s2(ktrd) = 0.e0
169         ENDIF
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
177         !     
178      END SELECT
179      !
180   END SUBROUTINE trd_2d
181
182
183   SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd, ctype, clpas )
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.
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
195      !!
196      INTEGER ::   ji, jj, jk
197      CHARACTER(len=3) ::   cpas                                              ! number of passage
198      REAL(wp) ::   zbt, zbtu, zbtv, zmsku, zmskv                             ! temporary scalars
199      !!----------------------------------------------------------------------
200
201      ! Control of optional arguments
202      cpas = 'fst'
203      IF( PRESENT(clpas) )  cpas = clpas
204
205      ! 1. Mask the trends
206      ! ------------------
207
208      SELECT CASE( ctype )
209      !
210      CASE( 'DYN' )              ! Momentum       
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
218               END DO
219            END DO
220         END DO
221         ptrd3dx(jpi, : ,:) = 0.e0      ;      ptrd3dy(jpi, : ,:) = 0.e0
222         ptrd3dx( : ,jpj,:) = 0.e0      ;      ptrd3dy( : ,jpj,:) = 0.e0
223         !
224      CASE( 'TRA' )              ! Tracers
225         DO jk = 1, jpk
226            ptrd3dx(:,:,jk) = ptrd3dx(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
227            ptrd3dy(:,:,jk) = ptrd3dy(:,:,jk) * tmask(:,:,jk) * tmask_i(:,:)
228         END DO
229         !
230      END SELECT   
231
232      ! 2. Basin averaged tracer/momentum trends
233      ! ----------------------------------------
234     
235      SELECT CASE( ctype )
236      !
237      CASE( 'DYN' )              ! Momentum
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
250         !
251      CASE( 'TRA' )              ! Tracers
252         IF( cpas == 'fst' )   THEN
253            tmo(ktrd) = 0.e0
254            smo(ktrd) = 0.e0
255         ENDIF
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
265         !
266      END SELECT
267
268      ! 3. Basin averaged tracer/momentum square trends
269      ! -----------------------------------------------
270      ! c a u t i o n: field now
271     
272      SELECT CASE( ctype )
273      !
274      CASE( 'DYN' )              ! Momentum
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
287         !
288      CASE( 'TRA' )              ! Tracers
289         IF( cpas == 'fst' )   THEN
290            t2(ktrd) = 0.e0
291            s2(ktrd) = 0.e0
292         ENDIF
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
302         !
303      END SELECT
304      !
305   END SUBROUTINE trd_3d
306
307
308
309   SUBROUTINE trd_icp_init
310      !!---------------------------------------------------------------------
311      !!                  ***  ROUTINE trd_icp_init  ***
312      !!
313      !! ** Purpose :   Read the namtrd namelist
314      !!----------------------------------------------------------------------
315      INTEGER  ::   ji, jj, jk
316      REAL(wp) ::   zmskt
317#if  defined key_trddyn
318      REAL(wp) ::   zmsku, zmskv
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
340      IF(lwp) WRITE(numout,*) '                total ocean volume at T-point   tvolt = ',tvolt
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
364         WRITE(numout,*) '                total ocean volume at U-point   tvolu = ',tvolu
365         WRITE(numout,*) '                total ocean volume at V-point   tvolv = ',tvolv
366      ENDIF
367#endif
368      !
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
377      !!----------------------------------------------------------------------
378      INTEGER, INTENT(in) ::   kt                                  ! ocean time-step index
379      !!
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
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
394         zkx(:,:,:)   = 0.e0
395         zky(:,:,:)   = 0.e0
396         zkz(:,:,:)   = 0.e0
397         zkepe(:,:,:) = 0.e0
398   
399         CALL eos( tn, sn, rhd, rhop )       ! now potential and in situ densities
400
401         ! Density flux at w-point
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
410         zkz(:,:,1) = 0.e0
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
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
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
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_had) / tvolu, vmo(jpicpd_had) / tvolv
475            WRITE (numout,9507) umo(jpicpd_zad) / tvolu, vmo(jpicpd_zad) / tvolv
476            WRITE (numout,9508) umo(jpicpd_zdf) / tvolu, vmo(jpicpd_zdf) / tvolv
477            WRITE (numout,9509) umo(jpicpd_spg) / tvolu, vmo(jpicpd_spg) / tvolv
478            WRITE (numout,9510) umo(jpicpd_swf) / tvolu, vmo(jpicpd_swf) / tvolv
479            WRITE (numout,9511) umo(jpicpd_dat) / tvolu, vmo(jpicpd_dat) / tvolv
480            WRITE (numout,9512) umo(jpicpd_bfr) / tvolu, vmo(jpicpd_bfr) / tvolv
481            WRITE (numout,9513)
482            WRITE (numout,9514)                                                 &
483            &     (  umo(jpicpd_hpg) + umo(jpicpd_keg) + umo(jpicpd_rvo) + umo(jpicpd_pvo) + umo(jpicpd_ldf)   &
484            &      + umo(jpicpd_had) + umo(jpicpd_zad) + umo(jpicpd_zdf) + umo(jpicpd_spg) + umo(jpicpd_dat)   &
485            &      + umo(jpicpd_swf) + umo(jpicpd_bfr) ) / tvolu,   &
486            &     (  vmo(jpicpd_hpg) + vmo(jpicpd_keg) + vmo(jpicpd_rvo) + vmo(jpicpd_pvo) + vmo(jpicpd_ldf)   &
487            &      + vmo(jpicpd_had) + vmo(jpicpd_zad) + vmo(jpicpd_zdf) + vmo(jpicpd_spg) + vmo(jpicpd_dat)   &
488            &      + vmo(jpicpd_swf) + vmo(jpicpd_bfr) ) / tvolv
489         ENDIF
490
491 9500    FORMAT(' momentum trend at it= ', i6, ' :', /' ==============================')
492 9501    FORMAT(' pressure gradient          u= ', e20.13, '    v= ', e20.13)
493 9502    FORMAT(' ke gradient                u= ', e20.13, '    v= ', e20.13)
494 9503    FORMAT(' relative vorticity term    u= ', e20.13, '    v= ', e20.13)
495 9504    FORMAT(' coriolis term              u= ', e20.13, '    v= ', e20.13)
496 9505    FORMAT(' horizontal diffusion       u= ', e20.13, '    v= ', e20.13)
497 9506    FORMAT(' horizontal advection       u= ', e20.13, '    v= ', e20.13)
498 9507    FORMAT(' vertical advection         u= ', e20.13, '    v= ', e20.13)
499 9508    FORMAT(' vertical diffusion         u= ', e20.13, '    v= ', e20.13)
500 9509    FORMAT(' surface pressure gradient  u= ', e20.13, '    v= ', e20.13)
501 9510    FORMAT(' surface wind forcing       u= ', e20.13, '    v= ', e20.13)
502 9511    FORMAT(' dampimg term               u= ', e20.13, '    v= ', e20.13)
503 9512    FORMAT(' bottom flux                u= ', e20.13, '    v= ', e20.13)
504 9513    FORMAT(' -----------------------------------------------------------------------------')
505 9514    FORMAT(' total trend                u= ', e20.13, '    v= ', e20.13)
506
507         IF(lwp) THEN
508            WRITE (numout,*)
509            WRITE (numout,*)
510            WRITE (numout,9520) kt
511            WRITE (numout,9521) hke(jpicpd_hpg) / tvolt
512            WRITE (numout,9522) hke(jpicpd_keg) / tvolt
513            WRITE (numout,9523) hke(jpicpd_rvo) / tvolt
514            WRITE (numout,9524) hke(jpicpd_pvo) / tvolt
515            WRITE (numout,9525) hke(jpicpd_ldf) / tvolt
516            WRITE (numout,9526) hke(jpicpd_had) / tvolt
517            WRITE (numout,9527) hke(jpicpd_zad) / tvolt
518            WRITE (numout,9528) hke(jpicpd_zdf) / tvolt
519            WRITE (numout,9529) hke(jpicpd_spg) / tvolt
520            WRITE (numout,9530) hke(jpicpd_swf) / tvolt
521            WRITE (numout,9531) hke(jpicpd_dat) / tvolt
522            WRITE (numout,9532)
523            WRITE (numout,9533)   &
524            &     (  hke(jpicpd_hpg) + hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_pvo) + hke(jpicpd_ldf)   &
525            &      + hke(jpicpd_had) + hke(jpicpd_zad) + hke(jpicpd_zdf) + hke(jpicpd_spg) + hke(jpicpd_dat)   &
526            &      + hke(jpicpd_swf) ) / tvolt
527         ENDIF
528
529 9520    FORMAT(' kinetic energy trend at it= ', i6, ' :', /' ====================================')
530 9521    FORMAT(' pressure gradient         u2= ', e20.13)
531 9522    FORMAT(' ke gradient               u2= ', e20.13)
532 9523    FORMAT(' relative vorticity term   u2= ', e20.13)
533 9524    FORMAT(' coriolis term             u2= ', e20.13)
534 9525    FORMAT(' horizontal diffusion      u2= ', e20.13)
535 9526    FORMAT(' horizontal advection      u2= ', e20.13)
536 9527    FORMAT(' vertical advection        u2= ', e20.13)
537 9528    FORMAT(' vertical diffusion        u2= ', e20.13)
538 9529    FORMAT(' surface pressure gradient u2= ', e20.13)
539 9530    FORMAT(' surface wind forcing      u2= ', e20.13)
540 9531    FORMAT(' dampimg term              u2= ', e20.13)
541 9532    FORMAT(' --------------------------------------------------')
542 9533    FORMAT(' total trend               u2= ', e20.13)
543
544         IF(lwp) THEN
545            WRITE (numout,*)
546            WRITE (numout,*)
547            WRITE (numout,9540) kt
548            WRITE (numout,9541) ( hke(jpicpd_keg) + hke(jpicpd_rvo) + hke(jpicpd_had) + hke(jpicpd_zad) ) / tvolt
549            WRITE (numout,9542) ( hke(jpicpd_keg) + hke(jpicpd_had) + hke(jpicpd_zad) ) / tvolt
550            WRITE (numout,9543) ( hke(jpicpd_pvo) ) / tvolt
551            WRITE (numout,9544) ( hke(jpicpd_rvo) ) / tvolt
552            WRITE (numout,9545) ( hke(jpicpd_spg) ) / tvolt
553            WRITE (numout,9546) ( hke(jpicpd_ldf) ) / tvolt
554            WRITE (numout,9547) ( hke(jpicpd_zdf) ) / tvolt
555            WRITE (numout,9548) ( hke(jpicpd_hpg) ) / tvolt, rpktrd / tvolt
556            WRITE (numout,*)
557            WRITE (numout,*)
558         ENDIF
559
560 9540    FORMAT(' energetic consistency at it= ', i6, ' :', /' =========================================')
561 9541    FORMAT(' 0 = non linear term(true if key_vorenergy or key_combined): ', e20.13)
562 9542    FORMAT(' 0 = ke gradient + horizontal + vertical advection         : ', e20.13)
563 9543    FORMAT(' 0 = coriolis term  (true if key_vorenergy or key_combined): ', e20.13)
564 9544    FORMAT(' 0 = uh.( rot(u) x uh ) (true if enstrophy conser.)        : ', e20.13)
565 9545    FORMAT(' 0 = surface pressure gradient                             : ', e20.13)
566 9546    FORMAT(' 0 > horizontal diffusion                                  : ', e20.13)
567 9547    FORMAT(' 0 > vertical diffusion                                    : ', e20.13)
568 9548    FORMAT(' pressure gradient u2 = - 1/rau0 u.dz(rhop)                : ', e20.13, '  u.dz(rhop) =', e20.13)
569         !
570         ! Save potential to kinetic energy conversion for next time step
571         rpktrd = peke
572         !
573      ENDIF
574      !
575   END SUBROUTINE trd_dwr
576
577
578   SUBROUTINE trd_twr( kt )
579      !!---------------------------------------------------------------------
580      !!                  ***  ROUTINE trd_twr  ***
581      !!
582      !! ** Purpose :  write active tracers trends in ocean.output
583      !!----------------------------------------------------------------------
584      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
585      !!----------------------------------------------------------------------
586
587      ! I. Tracers trends
588      ! -----------------
589
590      IF( MOD(kt,ntrd) == 0 .OR. kt == nit000 .OR. kt == nitend ) THEN
591
592         ! I.1 Sums over the global domain
593         ! -------------------------------
594         IF( lk_mpp ) THEN
595            CALL mpp_sum( tmo, jptot_tra )   
596            CALL mpp_sum( smo, jptot_tra )
597            CALL mpp_sum( t2 , jptot_tra )
598            CALL mpp_sum( s2 , jptot_tra )
599         ENDIF
600
601         ! I.2 Print tracers trends in the ocean.output file
602         ! --------------------------------------------------
603         
604         IF(lwp) THEN
605            WRITE (numout,*)
606            WRITE (numout,*)
607            WRITE (numout,9400) kt
608            WRITE (numout,9401) (tmo(jpicpt_xad)+tmo(jpicpt_yad))/ tvolt, (smo(jpicpt_xad)+smo(jpicpt_yad))/ tvolt
609            WRITE (numout,9402)  tmo(jpicpt_zad) / tvolt, smo(jpicpt_zad) / tvolt
610            WRITE (numout,9403)  tmo(jpicpt_ldf) / tvolt, smo(jpicpt_ldf) / tvolt
611            WRITE (numout,9404)  tmo(jpicpt_zdf) / tvolt, smo(jpicpt_zdf) / tvolt
612            WRITE (numout,9405)  tmo(jpicpt_npc) / tvolt, smo(jpicpt_npc) / tvolt
613            WRITE (numout,9406)  tmo(jpicpt_dmp) / tvolt, smo(jpicpt_dmp) / tvolt
614            WRITE (numout,9407)  tmo(jpicpt_qsr) / tvolt
615            WRITE (numout,9408)  tmo(jpicpt_nsr) / tvolt, smo(jpicpt_nsr) / tvolt
616            WRITE (numout,9409) 
617            WRITE (numout,9410) (  tmo(jpicpt_xad) + tmo(jpicpt_yad) + tmo(jpicpt_zad) + tmo(jpicpt_ldf) + tmo(jpicpt_zdf)   &
618            &                    + tmo(jpicpt_npc) + tmo(jpicpt_dmp) + tmo(jpicpt_qsr) + tmo(jpicpt_nsr) ) / tvolt,   &
619            &                   (  smo(jpicpt_xad) + smo(jpicpt_yad) + smo(jpicpt_zad) + smo(jpicpt_ldf) + smo(jpicpt_zdf)   &
620            &                    + smo(jpicpt_npc) + smo(jpicpt_dmp)                   + smo(jpicpt_nsr) ) / tvolt
621         ENDIF
622
6239400     FORMAT(' tracer trend at it= ',i6,' :     temperature',   &
624              '              salinity',/' ============================')
6259401     FORMAT(' horizontal advection        ',e20.13,'     ',e20.13)
6269402     FORMAT(' vertical advection          ',e20.13,'     ',e20.13)
6279403     FORMAT(' horizontal diffusion        ',e20.13,'     ',e20.13)
6289404     FORMAT(' vertical diffusion          ',e20.13,'     ',e20.13)
6299405     FORMAT(' static instability mixing   ',e20.13,'     ',e20.13)
6309406     FORMAT(' damping term                ',e20.13,'     ',e20.13)
6319407     FORMAT(' penetrative qsr             ',e20.13)
6329408     FORMAT(' non solar radiation         ',e20.13,'     ',e20.13)
6339409     FORMAT(' -------------------------------------------------------------------------')
6349410     FORMAT(' total trend                 ',e20.13,'     ',e20.13)
635
636
637         IF(lwp) THEN
638            WRITE (numout,*)
639            WRITE (numout,*)
640            WRITE (numout,9420) kt
641            WRITE (numout,9421) ( t2(jpicpt_xad)+t2(jpicpt_yad) )/ tvolt, ( s2(jpicpt_xad)+s2(jpicpt_yad) )/ tvolt
642            WRITE (numout,9422)   t2(jpicpt_zad) / tvolt, s2(jpicpt_zad) / tvolt
643            WRITE (numout,9423)   t2(jpicpt_ldf) / tvolt, s2(jpicpt_ldf) / tvolt
644            WRITE (numout,9424)   t2(jpicpt_zdf) / tvolt, s2(jpicpt_zdf) / tvolt
645            WRITE (numout,9425)   t2(jpicpt_npc) / tvolt, s2(jpicpt_npc) / tvolt
646            WRITE (numout,9426)   t2(jpicpt_dmp) / tvolt, s2(jpicpt_dmp) / tvolt
647            WRITE (numout,9427)   t2(jpicpt_qsr) / tvolt
648            WRITE (numout,9428)   t2(jpicpt_nsr) / tvolt, s2(jpicpt_nsr) / tvolt
649            WRITE (numout,9429)
650            WRITE (numout,9430) (  t2(jpicpt_xad) + t2(jpicpt_yad) + t2(jpicpt_zad) + t2(jpicpt_ldf) + t2(jpicpt_zdf)   &
651            &                    + t2(jpicpt_npc) + t2(jpicpt_dmp) + t2(jpicpt_qsr) + t2(jpicpt_nsr) ) / tvolt,   &
652            &                   (  s2(jpicpt_xad) + s2(jpicpt_yad) + s2(jpicpt_zad) + s2(jpicpt_ldf) + s2(jpicpt_zdf)   &
653            &                    + s2(jpicpt_npc) + s2(jpicpt_dmp)                  + s2(jpicpt_nsr) ) / tvolt
654         ENDIF
655
6569420     FORMAT(' tracer**2 trend at it= ', i6, ' :      temperature',   &
657            '               salinity', /, ' ===============================')
6589421     FORMAT(' horizontal advection      * t   ', e20.13, '     ', e20.13)
6599422     FORMAT(' vertical advection        * t   ', e20.13, '     ', e20.13)
6609423     FORMAT(' horizontal diffusion      * t   ', e20.13, '     ', e20.13)
6619424     FORMAT(' vertical diffusion        * t   ', e20.13, '     ', e20.13)
6629425     FORMAT(' static instability mixing * t   ', e20.13, '     ', e20.13)
6639426     FORMAT(' damping term              * t   ', e20.13, '     ', e20.13)
6649427     FORMAT(' penetrative qsr           * t   ', e20.13)
6659428     FORMAT(' non solar radiation       * t   ', e20.13, '     ', e20.13)
6669429     FORMAT(' -----------------------------------------------------------------------------')
6679430     FORMAT(' total trend                *t = ', e20.13, '  *s = ', e20.13)
668
669
670         IF(lwp) THEN
671            WRITE (numout,*)
672            WRITE (numout,*)
673            WRITE (numout,9440) kt
674            WRITE (numout,9441) ( tmo(jpicpt_xad)+tmo(jpicpt_yad)+tmo(jpicpt_zad) )/tvolt,   &
675            &                   ( smo(jpicpt_xad)+smo(jpicpt_yad)+smo(jpicpt_zad) )/tvolt
676            WRITE (numout,9442)   tmo(jpicpt_zl1)/tvolt,  smo(jpicpt_zl1)/tvolt
677            WRITE (numout,9443)   tmo(jpicpt_ldf)/tvolt,  smo(jpicpt_ldf)/tvolt
678            WRITE (numout,9444)   tmo(jpicpt_zdf)/tvolt,  smo(jpicpt_zdf)/tvolt
679            WRITE (numout,9445)   tmo(jpicpt_npc)/tvolt,  smo(jpicpt_npc)/tvolt
680            WRITE (numout,9446) ( t2(jpicpt_xad)+t2(jpicpt_yad)+t2(jpicpt_zad) )/tvolt,    &
681            &                   ( s2(jpicpt_xad)+s2(jpicpt_yad)+s2(jpicpt_zad) )/tvolt
682            WRITE (numout,9447)   t2(jpicpt_ldf)/tvolt,   s2(jpicpt_ldf)/tvolt
683            WRITE (numout,9448)   t2(jpicpt_zdf)/tvolt,   s2(jpicpt_zdf)/tvolt
684            WRITE (numout,9449)   t2(jpicpt_npc)/tvolt,   s2(jpicpt_npc)/tvolt
685         ENDIF
686
6879440     FORMAT(' tracer consistency at it= ',i6,   &
688            ' :         temperature','                salinity',/,   &
689            ' ==================================')
6909441     FORMAT(' 0 = horizontal+vertical advection +    ',e20.13,'       ',e20.13)
6919442     FORMAT('     1st lev vertical advection         ',e20.13,'       ',e20.13)
6929443     FORMAT(' 0 = horizontal diffusion               ',e20.13,'       ',e20.13)
6939444     FORMAT(' 0 = vertical diffusion                 ',e20.13,'       ',e20.13)
6949445     FORMAT(' 0 = static instability mixing          ',e20.13,'       ',e20.13)
6959446     FORMAT(' 0 = horizontal+vertical advection * t  ',e20.13,'       ',e20.13)
6969447     FORMAT(' 0 > horizontal diffusion          * t  ',e20.13,'       ',e20.13)
6979448     FORMAT(' 0 > vertical diffusion            * t  ',e20.13,'       ',e20.13)
6989449     FORMAT(' 0 > static instability mixing     * t  ',e20.13,'       ',e20.13)
699         !
700      ENDIF
701      !
702   END SUBROUTINE trd_twr
703
704#   else
705   !!----------------------------------------------------------------------
706   !!   Default case :                                         Empty module
707   !!----------------------------------------------------------------------
708   INTERFACE trd_icp
709      MODULE PROCEDURE trd_2d, trd_3d
710   END INTERFACE
711
712CONTAINS
713   SUBROUTINE trd_2d( ptrd2dx, ptrd2dy, ktrd , ctype, clpas )       ! Empty routine
714      REAL, DIMENSION(:,:) ::   ptrd2dx, ptrd2dy
715      INTEGER                     , INTENT(in   ) ::   ktrd         ! tracer trend index
716      CHARACTER(len=3)            , INTENT(in   ) ::   ctype        ! momentum ('DYN') or tracers ('TRA') trends
717      CHARACTER(len=3), INTENT(in), OPTIONAL ::   clpas             ! number of passage
718      WRITE(*,*) 'trd_2d: You should not have seen this print! error ?', &
719          &       ptrd2dx(1,1), ptrd2dy(1,1), ktrd, ctype, clpas
720   END SUBROUTINE trd_2d
721   SUBROUTINE trd_3d( ptrd3dx, ptrd3dy, ktrd , ctype, clpas )       ! Empty routine
722      REAL, DIMENSION(:,:,:) ::   ptrd3dx, ptrd3dy
723      INTEGER                     , INTENT(in   ) ::   ktrd         ! tracer trend index
724      CHARACTER(len=3)            , INTENT(in   ) ::   ctype        ! momentum ('DYN') or tracers ('TRA') trends
725      CHARACTER(len=3), INTENT(in), OPTIONAL ::   clpas             ! number of passage
726      WRITE(*,*) 'trd_3d: You should not have seen this print! error ?', &
727          &       ptrd3dx(1,1,1), ptrd3dy(1,1,1), ktrd, ctype, clpas
728   END SUBROUTINE trd_3d
729   SUBROUTINE trd_icp_init               ! Empty routine
730   END SUBROUTINE trd_icp_init
731   SUBROUTINE trd_dwr( kt )          ! Empty routine
732      WRITE(*,*) 'trd_dwr: You should not have seen this print! error ?', kt
733   END SUBROUTINE trd_dwr
734   SUBROUTINE trd_twr( kt )          ! Empty routine
735      WRITE(*,*) 'trd_twr: You should not have seen this print! error ?', kt
736   END SUBROUTINE trd_twr
737#endif
738
739   !!======================================================================
740END MODULE trdicp
Note: See TracBrowser for help on using the repository browser.