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.
trdtra.F90 in trunk/NEMO/OPA_SRC/TRD – NEMO

source: trunk/NEMO/OPA_SRC/TRD/trdtra.F90 @ 3

Last change on this file since 3 was 3, checked in by opalod, 20 years ago

Initial revision

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 14.8 KB
Line 
1MODULE trdtra
2   !!======================================================================
3   !!                       ***  MODULE  trdtra  ***
4   !! Ocean diagnostics:  ocean tracer trends
5   !!=====================================================================
6#if  defined key_trdtra   ||   defined key_esopa
7   !!----------------------------------------------------------------------
8   !!   'key_trdtra'  or                           tracer trend diagnostics
9   !!   'key_trdmld'                          mixed layer trend diagnostics
10   !!----------------------------------------------------------------------
11
12   !!----------------------------------------------------------------------
13   !!   trd_tra      : verify the basin averaged properties for tracers
14   !!   trd_tra_init : ???
15   !!----------------------------------------------------------------------
16   !! * Modules used
17   USE oce             ! ocean dynamics and tracers variables
18   USE dom_oce         ! ocean space and time domain variables
19   USE trdtra_oce      ! ocean active tracer trend variables
20   USE trddyn_oce      ! ocean dynamics trend variables
21   USE ldftra_oce      ! ocean active tracers: lateral physics
22   USE ldfdyn_oce      ! ocean dynamics: lateral physics
23   USE zdf_oce         ! ocean vertical physics
24   USE in_out_manager  ! I/O manager
25   USE lib_mpp         ! distibuted memory computing library
26
27   IMPLICIT NONE
28   PRIVATE
29
30   !! * Routine accessibility
31   PUBLIC trd_tra      ! called by step.F90
32   PUBLIC trd_tra_init ! called by opa.F90
33
34   !! * Shared module variables
35   LOGICAL, PUBLIC, PARAMETER ::   lk_trdtra = .TRUE.   ! momentum trend flag
36
37   !! * Substitutions
38#  include "domzgr_substitute.h90"
39#  include "vectopt_loop_substitute.h90"
40   !!----------------------------------------------------------------------
41   !!   OPA 9.0 , LODYC-IPSL  (2003)
42   !!----------------------------------------------------------------------
43
44CONTAINS
45
46   SUBROUTINE trd_tra( kt )
47      !!---------------------------------------------------------------------
48      !!                  ***  ROUTINE trd_tra  ***
49      !!
50      !! ** Purpose :   verify the basin averaged properties of the tracers
51      !!   equations at every time step frequency ntrd.
52      !!
53      !!   Method :
54      !!
55      !! History :
56      !!        !  91-12 (G. Madec)
57      !!        !  92-06 (M. Imbard) add time step frequency
58      !!        !  96-01 (G. Madec) terrain following coordinates
59      !!   8.5  !  02-06  (G. Madec)  F90: Free form and module
60      !!----------------------------------------------------------------------
61      !! * Arguments
62      INTEGER, INTENT( in ) ::   kt      ! ocean time-step index
63
64      !! * Local declarations
65      INTEGER :: ji, jj, jk, jn
66      REAL(wp)                ::   zbt 
67      REAL(wp), DIMENSION(10) ::   ztmo, zsmo, zt2, zs2
68
69      NAMELIST/namtrd/ ntrd, nctls
70      !!----------------------------------------------------------------------
71
72      ! 0. Initialization
73      ! -----------------
74
75      IF( kt == nit000 ) THEN
76
77         ! namelist namtrd : trend diagnostic
78         REWIND( numnam )
79         READ  ( numnam, namtrd )
80
81         IF(lwp) THEN
82            WRITE(numout,*) 'trd_tra : read namelist namtrd'
83            WRITE(numout,*) '~~~~~~~'
84            WRITE(numout,*) ' time step frequency trend       ntrd  = ', ntrd
85            WRITE(numout,*) ' '
86         ENDIF
87
88         ! Total volume at t-points:
89         tvolt = 0.
90         DO jk = 1, jpkm1
91            DO jj = 2, jpjm1
92               DO ji = fs_2, fs_jpim1   ! vector opt.
93                  tvolt = tvolt + e1t(ji,jj) *e2t(ji,jj) * fse3t(ji,jj,jk)* tmask(ji,jj,jk) * tmask_i(ji,jj)
94               END DO
95            END DO
96         END DO
97#if defined key_mpp
98         CALL mpp_sum( tvolt )
99#endif
100         IF(lwp) THEN
101            WRITE(numout,*)
102            WRITE(numout,*) '***trd_tra'
103            WRITE(numout,*) ' frequency ntrd = ',ntrd
104            WRITE(numout,*) ' '
105            WRITE(numout,*) ' tvolt = ',tvolt
106            WRITE(numout,*) ' '
107         ENDIF
108
109      ENDIF
110
111
112      ! 1. Advective trends and forcing trend
113      ! -------------------------------------
114
115      IF( MOD(kt,ntrd) == 0 .OR. kt  == nit000 .OR. kt == nitend ) THEN
116         
117         ! 1.1 Mask the forcing trend and substract it from the vertical diffusion trend
118         flxtrd(:,:,1) = flxtrd(:,:,1) * tmask_i(:,:)
119         flxtrd(:,:,2) = flxtrd(:,:,2) * tmask_i(:,:)
120         
121         ! 1.2 Mask the trends
122         DO jn = 1, 6
123            DO jk = 1, jpk
124               ttrd(:,:,jk,jn) = ttrd(:,:,jk,jn)* tmask(:,:,jk) * tmask_i(:,:)
125               strd(:,:,jk,jn) = strd(:,:,jk,jn)* tmask(:,:,jk) * tmask_i(:,:)
126            END DO
127         END DO
128         
129         DO jk = 1, jpk
130            ttrd(:,:,jk,7) = ttrd(:,:,jk,7) * tmask(:,:,jk) * tmask_i(:,:)
131         END DO
132
133
134         ! 2. Basin averaged tracer trend
135         ! ------------------------------
136         
137         DO jn = 1, 6
138            ztmo(jn) = 0.
139            zsmo(jn) = 0.
140            DO jk = 1, jpkm1
141               DO jj = 1, jpj
142                  DO ji = 1, jpi
143                     zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) 
144                     ztmo(jn) =  ztmo(jn) + ttrd(ji,jj,jk,jn) * zbt
145                     zsmo(jn) =  zsmo(jn) + strd(ji,jj,jk,jn) * zbt
146                  END DO
147               END DO
148            END DO
149         END DO
150         
151         ztmo(7) = 0.
152         DO jk = 1, jpk
153            DO jj = 1, jpj
154               DO ji = 1, jpi
155                  zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
156                  ztmo(7) =  ztmo(7) + ttrd(ji,jj,jk,7) * zbt
157               END DO
158            END DO
159         END DO
160         
161         ztmo(8) = 0.
162         zsmo(8) = 0.
163         DO jj = 1, jpj
164            DO ji = 1, jpi
165               zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
166               ztmo(8) =  ztmo(8) + flxtrd(ji,jj,1) * zbt
167               zsmo(8) =  zsmo(8) + flxtrd(ji,jj,2) * zbt
168            END DO
169         END DO
170         
171         
172         ! 3. Basin averaged tracer square trend
173         ! -------------------------------------
174         ! c a u t i o n: field before, because after the array swap
175         
176         DO jn = 1, 6
177            zt2(jn) = 0.e0
178            zs2(jn) = 0.e0
179            DO jk = 1, jpk
180               DO jj = 1, jpj
181                  DO ji = 1, jpi
182                     zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
183                     zt2(jn) = zt2(jn) + ttrd(ji,jj,jk,jn) * zbt * tb(ji,jj,jk)
184                     zs2(jn) = zs2(jn) + strd(ji,jj,jk,jn) * zbt * sb(ji,jj,jk)
185                  END DO
186               END DO
187            END DO
188         END DO
189         
190         zt2(7) = 0.e0
191         DO jk = 1, jpk
192            DO jj = 1, jpj
193               DO ji = 1, jpi
194                  zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)
195                  zt2(7) = zt2(7) + ttrd(ji,jj,jk,7) * zbt * tb(ji,jj,jk)
196               END DO
197            END DO
198         END DO
199         
200         zt2(8) = 0.e0
201         zs2(8) = 0.e0
202         DO jj = 1, jpj
203            DO ji = 1, jpi
204               zbt = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
205               zt2(8) = zt2(8) + flxtrd(ji,jj,1) * zbt * tb(ji,jj,1)
206               zs2(8) = zs2(8) + flxtrd(ji,jj,2) * zbt * sb(ji,jj,1)
207            END DO
208         END DO
209         
210#if defined key_mpp
211         CALL mpp_sum( ztmo, 10 )
212         CALL mpp_sum( zsmo, 10 )
213         CALL mpp_sum( zt2 , 10 )
214         CALL mpp_sum( zs2 , 10 )
215#endif
216         
217         ! 4. Print
218         ! --------
219         
220         IF(lwp) THEN
221            WRITE (numout,*)
222            WRITE (numout,*)
223            WRITE (numout,9400) kt
224            WRITE (numout,9401) ztmo(1) / tvolt, zsmo(1) / tvolt
225            WRITE (numout,9402) ztmo(2) / tvolt, zsmo(2) / tvolt
226            WRITE (numout,9403) ztmo(3) / tvolt, zsmo(3) / tvolt
227            WRITE (numout,9404) ztmo(4) / tvolt, zsmo(4) / tvolt
228            WRITE (numout,9405) ztmo(5) / tvolt, zsmo(5) / tvolt
229            WRITE (numout,9406) ztmo(6) / tvolt, zsmo(6) / tvolt
230            WRITE (numout,9407) ztmo(7) / tvolt
231            WRITE (numout,9408) ztmo(8) / tvolt, zsmo(8) / tvolt
232            WRITE (numout,9409)
233            WRITE (numout,9410) (  ztmo(1) + ztmo(2) + ztmo(3) + ztmo(4)              &
234            &                    + ztmo(5) + ztmo(6) + ztmo(7) + ztmo(8) ) / tvolt,   &
235            &                   (  zsmo(1) + zsmo(2) + zsmo(3) + zsmo(4)              &
236            &                    + zsmo(5) + zsmo(6)           + zsmo(8) ) / tvolt
237         ENDIF
238
2399400     FORMAT(' tracer trend at it= ',i6,' :     temperature',   &
240              '              salinity',/' ============================')
2419401     FORMAT(' horizontal advection        ',e20.13,'     ',e20.13)
2429402     FORMAT(' vertical advection          ',e20.13,'     ',e20.13)
2439403     FORMAT(' horizontal diffusion        ',e20.13,'     ',e20.13)
2449404     FORMAT(' vertical diffusion          ',e20.13,'     ',e20.13)
2459405     FORMAT(' STATIC instability mixing   ',e20.13,'     ',e20.13)
2469406     FORMAT(' damping term                ',e20.13,'     ',e20.13)
2479407     FORMAT(' penetrative qsr             ',e20.13,'     ',e20.13)
2489408     FORMAT(' forcing term                ',e20.13,'     ',e20.13)
2499409     FORMAT(' -------------------------------------------------------------------------')
2509410     FORMAT(' total trend                 ',e20.13,'     ',e20.13)
251
252
253         IF(lwp) THEN
254            WRITE (numout,*)
255            WRITE (numout,*)
256            WRITE (numout,9420) kt
257            WRITE (numout,9421) zt2(1) / tvolt, zs2(1) / tvolt
258            WRITE (numout,9422) zt2(2) / tvolt, zs2(2) / tvolt
259            WRITE (numout,9423) zt2(3) / tvolt, zs2(3) / tvolt
260            WRITE (numout,9424) zt2(4) / tvolt, zs2(4) / tvolt
261            WRITE (numout,9425) zt2(5) / tvolt, zs2(5) / tvolt
262            WRITE (numout,9426) zt2(6) / tvolt, zs2(6) / tvolt
263            WRITE (numout,9427) zt2(7) / tvolt
264            WRITE (numout,9428) zt2(8) / tvolt, zs2(8) / tvolt
265            WRITE (numout,9429)
266            WRITE (numout,9430) (  zt2(1) + zt2(2) + zt2(3) + zt2(4)              &
267            &                    + zt2(5) + zt2(6) + zt2(7) + zt2(8) ) / tvolt,   &
268            &                   (  zs2(1) + zs2(2) + zs2(3) + zs2(4)              &
269            &                    + zs2(5) + zs2(6)          + zs2(8) ) / tvolt
270         ENDIF
271
2729420     FORMAT(' tracer**2 trend at it= ', i6, ' :      temperature',   &
273            '               salinity', /, ' ===============================')
2749421     FORMAT(' horizontal advection      * t   ', e20.13, '     ', e20.13)
2759422     FORMAT(' vertical advection        * t   ', e20.13, '     ', e20.13)
2769423     FORMAT(' horizontal diffusion      * t   ', e20.13, '     ', e20.13)
2779424     FORMAT(' vertical diffusion        * t   ', e20.13, '     ', e20.13)
2789425     FORMAT(' STATIC instability mixing * t   ', e20.13, '     ', e20.13)
2799426     FORMAT(' damping term              * t   ', e20.13, '     ', e20.13)
2809427     FORMAT(' penetrative qsr           * t   ', e20.13, '     ', e20.13)
2819428     FORMAT(' forcing term              * t   ', e20.13, '     ', e20.13)
2829429     FORMAT(' -----------------------------------------------------------------------------')
2839430     FORMAT(' total trend                *t = ', e20.13, '  *s = ', e20.13)
284
285
286         IF(lwp) THEN
287            WRITE (numout,*)
288            WRITE (numout,*)
289            WRITE (numout,9440) kt
290            WRITE (numout,9441) ( ztmo(1)+ztmo(2) )/tvolt, ( zsmo(1)+zsmo(2) )/tvolt
291            WRITE (numout,9442)   ztmo(3)/tvolt,  zsmo(3)/tvolt
292            WRITE (numout,9443)   ztmo(4)/tvolt,  zsmo(4)/tvolt
293            WRITE (numout,9444)   ztmo(5)/tvolt,  zsmo(5)/tvolt
294            WRITE (numout,9445) ( zt2(1)+zt2(2) )/tvolt, ( zs2(1)+zs2(2) )/tvolt
295            WRITE (numout,9446)   zt2(3)/tvolt,   zs2(3)/tvolt
296            WRITE (numout,9447)   zt2(4)/tvolt,   zs2(4)/tvolt
297            WRITE (numout,9448)   zt2(5)/tvolt,   zs2(5)/tvolt
298         ENDIF
299
3009440     FORMAT(' tracer consistency at it= ',i6,   &
301            ' :         temperature','                salinity',/,   &
302            ' ==================================')
3039441     FORMAT(' 0 = horizontal+vertical advection      ',e20.13,'       ',e20.13)
3049442     FORMAT(' 0 = horizontal diffusion               ',e20.13,'       ',e20.13)
3059443     FORMAT(' 0 = vertical diffusion                 ',e20.13,'       ',e20.13)
3069444     FORMAT(' 0 = static instability mixing          ',e20.13,'       ',e20.13)
3079445     FORMAT(' 0 = horizontal+vertical advection * t  ',e20.13,'       ',e20.13)
3089446     FORMAT(' 0 > horizontal diffusion          * t  ',e20.13,'       ',e20.13)
3099447     FORMAT(' 0 > vertical diffusion            * t  ',e20.13,'       ',e20.13)
3109448     FORMAT(' 0 > static instability mixing     * t  ',e20.13,'       ',e20.13)
311      ENDIF
312
313   END SUBROUTINE trd_tra
314
315
316   SUBROUTINE trd_tra_init
317      !!---------------------------------------------------------------------
318      !!                  ***  ROUTINE trd_tra_init  ***
319      !!
320      !! ** Purpose :   
321      !!
322      !! ** Method  :
323      !!
324      !! History :
325      !!   9.0  !  03-09  (G. Madec)  Original code
326      !!----------------------------------------------------------------------
327      !! * Local declarations
328      INTEGER :: ji, jj, jk
329
330      NAMELIST/namtrd/ ntrd, nctls
331      !!----------------------------------------------------------------------
332
333      ! set to zero the tracers trends
334      ttrd  (:,:,:,:) = 0.e0
335      strd  (:,:,:,:) = 0.e0
336      ttrdh (:,:,:,:) = 0.e0
337      strdh (:,:,:,:) = 0.e0
338      flxtrd(:,:,  :) = 0.e0
339
340      ! namelist namtrd : trend diagnostic
341      REWIND( numnam )
342      READ  ( numnam, namtrd )
343
344      IF(lwp) THEN
345         WRITE(numout,*)
346         WRITE(numout,*) 'trd_tra : read namelist namtrd'
347         WRITE(numout,*) '~~~~~~~'
348         WRITE(numout,*) '          time step frequency trend       ntrd  = ', ntrd
349      ENDIF
350
351      ! Total volume at t-points:
352      tvolt = 0.
353      DO jk = 1, jpkm1
354         DO jj = 2, jpjm1
355            DO ji = fs_2, fs_jpim1   ! vector opt.
356               tvolt = tvolt + e1t(ji,jj) *e2t(ji,jj) * fse3t(ji,jj,jk)* tmask(ji,jj,jk) * tmask_i(ji,jj)
357            END DO
358         END DO
359      END DO
360#if defined key_mpp
361      CALL mpp_sum( tvolt )
362#endif
363      IF(lwp) THEN
364         WRITE(numout,*) '          total ocean volume at T-point   tvolt = ',tvolt
365      ENDIF
366
367   END SUBROUTINE trd_tra_init
368
369#   else
370   !!----------------------------------------------------------------------
371   !!   Default case :                                         Empty module
372   !!----------------------------------------------------------------------
373   LOGICAL, PUBLIC, PARAMETER ::   lk_trdtra = .FALSE.   ! momentum trend flag
374CONTAINS
375   SUBROUTINE trd_tra( kt )        ! Empty routine
376      WRITE(*,*) kt
377   END SUBROUTINE trd_tra
378   SUBROUTINE trd_tra_init         ! Empty routine
379   END SUBROUTINE trd_tra_init
380#endif
381
382   !!======================================================================
383END MODULE trdtra
Note: See TracBrowser for help on using the repository browser.