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 NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/TRD – NEMO

source: NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/TRD/trdtra.F90 @ 12340

Last change on this file since 12340 was 12340, checked in by acc, 4 years ago

Branch 2019/dev_r11943_MERGE_2019. This commit introduces basic do loop macro
substitution to the 2019 option 1, merge branch. These changes have been SETTE
tested. The only addition is the do_loop_substitute.h90 file in the OCE directory but
the macros defined therein are used throughout the code to replace identifiable, 2D-
and 3D- nested loop opening and closing statements with single-line alternatives. Code
indents are also adjusted accordingly.

The following explanation is taken from comments in the new header file:

This header file contains preprocessor definitions and macros used in the do-loop
substitutions introduced between version 4.0 and 4.2. The primary aim of these macros
is to assist in future applications of tiling to improve performance. This is expected
to be achieved by alternative versions of these macros in selected locations. The
initial introduction of these macros simply replaces all identifiable nested 2D- and
3D-loops with single line statements (and adjusts indenting accordingly). Do loops
are identifiable if they comform to either:

DO jk = ....

DO jj = .... DO jj = ...

DO ji = .... DO ji = ...
. OR .
. .

END DO END DO

END DO END DO

END DO

and white-space variants thereof.

Additionally, only loops with recognised jj and ji loops limits are treated; these are:
Lower limits of 1, 2 or fs_2
Upper limits of jpi, jpim1 or fs_jpim1 (for ji) or jpj, jpjm1 or fs_jpjm1 (for jj)

The macro naming convention takes the form: DO_2D_BT_LR where:

B is the Bottom offset from the PE's inner domain;
T is the Top offset from the PE's inner domain;
L is the Left offset from the PE's inner domain;
R is the Right offset from the PE's inner domain

So, given an inner domain of 2,jpim1 and 2,jpjm1, a typical example would replace:

DO jj = 2, jpj

DO ji = 1, jpim1
.
.

END DO

END DO

with:

DO_2D_01_10
.
.
END_2D

similar conventions apply to the 3D loops macros. jk loop limits are retained
through macro arguments and are not restricted. This includes the possibility of
strides for which an extra set of DO_3DS macros are defined.

In the example definition below the inner PE domain is defined by start indices of
(kIs, kJs) and end indices of (kIe, KJe)

#define DO_2D_00_00 DO jj = kJs, kJe ; DO ji = kIs, kIe
#define END_2D END DO ; END DO

TO DO:


Only conventional nested loops have been identified and replaced by this step. There are constructs such as:

DO jk = 2, jpkm1

z2d(:,:) = z2d(:,:) + e3w(:,:,jk,Kmm) * z3d(:,:,jk) * wmask(:,:,jk)

END DO

which may need to be considered.

  • Property svn:keywords set to Id
File size: 21.9 KB
Line 
1MODULE trdtra
2   !!======================================================================
3   !!                       ***  MODULE  trdtra  ***
4   !! Ocean diagnostics:  ocean tracers trends pre-processing
5   !!=====================================================================
6   !! History :  3.3  !  2010-06  (C. Ethe) creation for the TRA/TRC merge
7   !!            3.5  !  2012-02  (G. Madec) update the comments
8   !!----------------------------------------------------------------------
9
10   !!----------------------------------------------------------------------
11   !!   trd_tra       : pre-process the tracer trends
12   !!   trd_tra_adv   : transform a div(U.T) trend into a U.grad(T) trend
13   !!   trd_tra_mng   : tracer trend manager: dispatch to the diagnostic modules
14   !!   trd_tra_iom   : output 3D tracer trends using IOM
15   !!----------------------------------------------------------------------
16   USE oce            ! ocean dynamics and tracers variables
17   USE dom_oce        ! ocean domain
18   USE sbc_oce        ! surface boundary condition: ocean
19   USE zdf_oce        ! ocean vertical physics
20   USE trd_oce        ! trends: ocean variables
21   USE trdtrc         ! ocean passive mixed layer tracers trends
22   USE trdglo         ! trends: global domain averaged
23   USE trdpen         ! trends: Potential ENergy
24   USE trdmxl         ! ocean active mixed layer tracers trends
25   USE ldftra         ! ocean active tracers lateral physics
26   USE ldfslp
27   USE zdfddm         ! vertical physics: double diffusion
28   USE phycst         ! physical constants
29   !
30   USE in_out_manager ! I/O manager
31   USE iom            ! I/O manager library
32   USE lib_mpp        ! MPP library
33
34   IMPLICIT NONE
35   PRIVATE
36
37   PUBLIC   trd_tra   ! called by all tra_... modules
38
39   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   trdtx, trdty, trdt   ! use to store the temperature trends
40   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_evd  ! store avt_evd to calculate EVD trend
41
42   !! * Substitutions
43#  include "vectopt_loop_substitute.h90"
44#  include "do_loop_substitute.h90"
45   !!----------------------------------------------------------------------
46   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
47   !! $Id$
48   !! Software governed by the CeCILL license (see ./LICENSE)
49   !!----------------------------------------------------------------------
50CONTAINS
51
52   INTEGER FUNCTION trd_tra_alloc()
53      !!---------------------------------------------------------------------
54      !!                  ***  FUNCTION trd_tra_alloc  ***
55      !!---------------------------------------------------------------------
56      ALLOCATE( trdtx(jpi,jpj,jpk) , trdty(jpi,jpj,jpk) , trdt(jpi,jpj,jpk) , avt_evd(jpi,jpj,jpk), STAT= trd_tra_alloc )
57      !
58      CALL mpp_sum ( 'trdtra', trd_tra_alloc )
59      IF( trd_tra_alloc /= 0 )   CALL ctl_stop( 'STOP', 'trd_tra_alloc: failed to allocate arrays' )
60   END FUNCTION trd_tra_alloc
61
62
63   SUBROUTINE trd_tra( kt, Kmm, Krhs, ctype, ktra, ktrd, ptrd, pu, ptra )
64      !!---------------------------------------------------------------------
65      !!                  ***  ROUTINE trd_tra  ***
66      !!
67      !! ** Purpose : pre-process tracer trends
68      !!
69      !! ** Method  : - mask the trend
70      !!              - advection (ptra present) converte the incoming flux (U.T)
71      !!              into trend (U.T => -U.grat(T)=div(U.T)-T.div(U)) through a
72      !!              call to trd_tra_adv
73      !!              - 'TRA' case : regroup T & S trends
74      !!              - send the trends to trd_tra_mng (trdtrc) for further processing
75      !!----------------------------------------------------------------------
76      INTEGER                         , INTENT(in)           ::   kt      ! time step
77      CHARACTER(len=3)                , INTENT(in)           ::   ctype   ! tracers trends type 'TRA'/'TRC'
78      INTEGER                         , INTENT(in)           ::   ktra    ! tracer index
79      INTEGER                         , INTENT(in)           ::   ktrd    ! tracer trend index
80      INTEGER                         , INTENT(in)           ::   Kmm, Krhs ! time level indices
81      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)           ::   ptrd    ! tracer trend  or flux
82      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   pu      ! now velocity
83      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   ptra    ! now tracer variable
84      !
85      INTEGER ::   jk   ! loop indices
86      REAL(wp),        DIMENSION(jpi,jpj,jpk) ::   ztrds             ! 3D workspace
87      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zwt, zws, ztrdt   ! 3D workspace
88      !!----------------------------------------------------------------------
89      !     
90      IF( .NOT. ALLOCATED( trdtx ) ) THEN      ! allocate trdtra arrays
91         IF( trd_tra_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'trd_tra : unable to allocate arrays' )
92      ENDIF
93
94      IF( ctype == 'TRA' .AND. ktra == jp_tem ) THEN   !==  Temperature trend  ==!
95         !
96         SELECT CASE( ktrd )
97         !                            ! advection: transform the advective flux into a trend
98         CASE( jptra_xad )   ;   CALL trd_tra_adv( ptrd, pu, ptra, 'X', trdtx, Kmm ) 
99         CASE( jptra_yad )   ;   CALL trd_tra_adv( ptrd, pu, ptra, 'Y', trdty, Kmm ) 
100         CASE( jptra_zad )   ;   CALL trd_tra_adv( ptrd, pu, ptra, 'Z', trdt, Kmm )
101         CASE( jptra_bbc,    &        ! qsr, bbc: on temperature only, send to trd_tra_mng
102            &  jptra_qsr )   ;   trdt(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
103                                 ztrds(:,:,:) = 0._wp
104                                 CALL trd_tra_mng( trdt, ztrds, ktrd, kt, Kmm )
105 !!gm Gurvan, verify the jptra_evd trend please !
106         CASE( jptra_evd )   ;   avt_evd(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
107         CASE DEFAULT                 ! other trends: masked trends
108            trdt(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)              ! mask & store
109         END SELECT
110         !
111      ENDIF
112
113      IF( ctype == 'TRA' .AND. ktra == jp_sal ) THEN      !==  Salinity trends  ==!
114         !
115         SELECT CASE( ktrd )
116         !                            ! advection: transform the advective flux into a trend
117         !                            !            and send T & S trends to trd_tra_mng
118         CASE( jptra_xad  )   ;   CALL trd_tra_adv( ptrd , pu  , ptra, 'X'  , ztrds, Kmm ) 
119                                  CALL trd_tra_mng( trdtx, ztrds, ktrd, kt, Kmm   )
120         CASE( jptra_yad  )   ;   CALL trd_tra_adv( ptrd , pu  , ptra, 'Y'  , ztrds, Kmm ) 
121                                  CALL trd_tra_mng( trdty, ztrds, ktrd, kt, Kmm   )
122         CASE( jptra_zad  )   ;   CALL trd_tra_adv( ptrd , pu  , ptra, 'Z'  , ztrds, Kmm ) 
123                                  CALL trd_tra_mng( trdt , ztrds, ktrd, kt, Kmm   )
124         CASE( jptra_zdfp )           ! diagnose the "PURE" Kz trend (here: just before the swap)
125            !                         ! iso-neutral diffusion case otherwise jptra_zdf is "PURE"
126            ALLOCATE( zwt(jpi,jpj,jpk), zws(jpi,jpj,jpk), ztrdt(jpi,jpj,jpk) )
127            !
128            zwt(:,:, 1 ) = 0._wp   ;   zws(:,:, 1 ) = 0._wp            ! vertical diffusive fluxes
129            zwt(:,:,jpk) = 0._wp   ;   zws(:,:,jpk) = 0._wp
130            DO jk = 2, jpk
131               zwt(:,:,jk) = avt(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk)
132               zws(:,:,jk) = avs(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk)
133            END DO
134            !
135            ztrdt(:,:,jpk) = 0._wp   ;   ztrds(:,:,jpk) = 0._wp
136            DO jk = 1, jpkm1
137               ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / e3t(:,:,jk,Kmm)
138               ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / e3t(:,:,jk,Kmm) 
139            END DO
140            CALL trd_tra_mng( ztrdt, ztrds, jptra_zdfp, kt, Kmm ) 
141            !
142            !                         ! Also calculate EVD trend at this point.
143            zwt(:,:,:) = 0._wp   ;   zws(:,:,:) = 0._wp            ! vertical diffusive fluxes
144            DO jk = 2, jpk
145               zwt(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_tem,Krhs) - ts(:,:,jk,jp_tem,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk)
146               zws(:,:,jk) = avt_evd(:,:,jk) * ( ts(:,:,jk-1,jp_sal,Krhs) - ts(:,:,jk,jp_sal,Krhs) ) / e3w(:,:,jk,Kmm) * tmask(:,:,jk)
147            END DO
148            !
149            ztrdt(:,:,jpk) = 0._wp   ;   ztrds(:,:,jpk) = 0._wp
150            DO jk = 1, jpkm1
151               ztrdt(:,:,jk) = ( zwt(:,:,jk) - zwt(:,:,jk+1) ) / e3t(:,:,jk,Kmm)
152               ztrds(:,:,jk) = ( zws(:,:,jk) - zws(:,:,jk+1) ) / e3t(:,:,jk,Kmm) 
153            END DO
154            CALL trd_tra_mng( ztrdt, ztrds, jptra_evd, kt, Kmm ) 
155            !
156            DEALLOCATE( zwt, zws, ztrdt )
157            !
158         CASE DEFAULT                 ! other trends: mask and send T & S trends to trd_tra_mng
159            ztrds(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
160            CALL trd_tra_mng( trdt, ztrds, ktrd, kt, Kmm ) 
161         END SELECT
162      ENDIF
163
164      IF( ctype == 'TRC' ) THEN                           !==  passive tracer trend  ==!
165         !
166         SELECT CASE( ktrd )
167         !                            ! advection: transform the advective flux into a masked trend
168         CASE( jptra_xad )   ;   CALL trd_tra_adv( ptrd , pu , ptra, 'X', ztrds, Kmm ) 
169         CASE( jptra_yad )   ;   CALL trd_tra_adv( ptrd , pu , ptra, 'Y', ztrds, Kmm ) 
170         CASE( jptra_zad )   ;   CALL trd_tra_adv( ptrd , pu , ptra, 'Z', ztrds, Kmm ) 
171         CASE DEFAULT                 ! other trends: just masked
172                                 ztrds(:,:,:) = ptrd(:,:,:) * tmask(:,:,:)
173         END SELECT
174         !                            ! send trend to trd_trc
175         CALL trd_trc( ztrds, ktra, ktrd, kt, Kmm ) 
176         !
177      ENDIF
178      !
179   END SUBROUTINE trd_tra
180
181
182   SUBROUTINE trd_tra_adv( pf, pu, pt, cdir, ptrd, Kmm )
183      !!---------------------------------------------------------------------
184      !!                  ***  ROUTINE trd_tra_adv  ***
185      !!
186      !! ** Purpose :   transformed a advective flux into a masked advective trends
187      !!
188      !! ** Method  :   use the following transformation: -div(U.T) = - U grad(T) + T.div(U)
189      !!       i-advective trends = -un. di-1[T] = -( di-1[fi] - tn di-1[un] )
190      !!       j-advective trends = -un. di-1[T] = -( dj-1[fi] - tn dj-1[un] )
191      !!       k-advective trends = -un. di+1[T] = -( dk+1[fi] - tn dk+1[un] )
192      !!                where fi is the incoming advective flux.
193      !!----------------------------------------------------------------------
194      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pf      ! advective flux in one direction
195      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pu      ! now velocity   in one direction
196      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) ::   pt      ! now or before tracer
197      CHARACTER(len=1)                , INTENT(in   ) ::   cdir    ! X/Y/Z direction
198      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(  out) ::   ptrd    ! advective trend in one direction
199      INTEGER,  INTENT(in)                            ::   Kmm     ! time level index
200      !
201      INTEGER  ::   ji, jj, jk   ! dummy loop indices
202      INTEGER  ::   ii, ij, ik   ! index shift as function of the direction
203      !!----------------------------------------------------------------------
204      !
205      SELECT CASE( cdir )      ! shift depending on the direction
206      CASE( 'X' )   ;   ii = 1   ;   ij = 0   ;   ik = 0      ! i-trend
207      CASE( 'Y' )   ;   ii = 0   ;   ij = 1   ;   ik = 0      ! j-trend
208      CASE( 'Z' )   ;   ii = 0   ;   ij = 0   ;   ik =-1      ! k-trend
209      END SELECT
210      !
211      !                        ! set to zero uncomputed values
212      ptrd(jpi,:,:) = 0._wp   ;   ptrd(1,:,:) = 0._wp
213      ptrd(:,jpj,:) = 0._wp   ;   ptrd(:,1,:) = 0._wp
214      ptrd(:,:,jpk) = 0._wp
215      !
216      DO_3D_00_00( 1, jpkm1 )
217         ptrd(ji,jj,jk) = - (     pf (ji,jj,jk) - pf (ji-ii,jj-ij,jk-ik)                        &
218           &                  - ( pu(ji,jj,jk) - pu(ji-ii,jj-ij,jk-ik) ) * pt(ji,jj,jk)  )   &
219           &              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) * tmask(ji,jj,jk)
220      END_3D
221      !
222   END SUBROUTINE trd_tra_adv
223
224
225   SUBROUTINE trd_tra_mng( ptrdx, ptrdy, ktrd, kt, Kmm )
226      !!---------------------------------------------------------------------
227      !!                  ***  ROUTINE trd_tra_mng  ***
228      !!
229      !! ** Purpose :   Dispatch all tracer trends computation, e.g. 3D output,
230      !!                integral constraints, potential energy, and/or
231      !!                mixed layer budget.
232      !!----------------------------------------------------------------------
233      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdx   ! Temperature or U trend
234      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdy   ! Salinity    or V trend
235      INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
236      INTEGER                   , INTENT(in   ) ::   kt      ! time step
237      INTEGER                   , INTENT(in   ) ::   Kmm     ! time level index
238      !!----------------------------------------------------------------------
239
240      IF( neuler == 0 .AND. kt == nit000    ) THEN   ;   r2dt =      rdt      ! = rdt (restart with Euler time stepping)
241      ELSEIF(               kt <= nit000 + 1) THEN   ;   r2dt = 2. * rdt      ! = 2 rdt (leapfrog)
242      ENDIF
243
244      !                   ! 3D output of tracers trends using IOM interface
245      IF( ln_tra_trd )   CALL trd_tra_iom ( ptrdx, ptrdy, ktrd, kt, Kmm )
246
247      !                   ! Integral Constraints Properties for tracers trends                                       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
248      IF( ln_glo_trd )   CALL trd_glo( ptrdx, ptrdy, ktrd, 'TRA', kt, Kmm )
249
250      !                   ! Potential ENergy trends
251      IF( ln_PE_trd  )   CALL trd_pen( ptrdx, ptrdy, ktrd, kt, r2dt, Kmm )
252
253      !                   ! Mixed layer trends for active tracers
254      IF( ln_tra_mxl )   THEN   
255         !-----------------------------------------------------------------------------------------------
256         ! W.A.R.N.I.N.G :
257         ! jptra_ldf : called by traldf.F90
258         !                 at this stage we store:
259         !                  - the lateral geopotential diffusion (here, lateral = horizontal)
260         !                  - and the iso-neutral diffusion if activated
261         ! jptra_zdf : called by trazdf.F90
262         !                 * in case of iso-neutral diffusion we store the vertical diffusion component in the
263         !                   lateral trend including the K_z contrib, which will be removed later (see trd_mxl)
264         !-----------------------------------------------------------------------------------------------
265
266         SELECT CASE ( ktrd )
267         CASE ( jptra_xad )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_xad, '3D' )   ! zonal    advection
268         CASE ( jptra_yad )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_yad, '3D' )   ! merid.   advection
269         CASE ( jptra_zad )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_zad, '3D' )   ! vertical advection
270         CASE ( jptra_ldf )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_ldf, '3D' )   ! lateral  diffusion
271         CASE ( jptra_bbl )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_bbl, '3D' )   ! bottom boundary layer
272         CASE ( jptra_zdf )
273            IF( ln_traldf_iso ) THEN ; CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_ldf, '3D' )   ! lateral  diffusion (K_z)
274            ELSE                   ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_zdf, '3D' )   ! vertical diffusion (K_z)
275            ENDIF
276         CASE ( jptra_dmp )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_dmp, '3D' )   ! internal 3D restoring (tradmp)
277         CASE ( jptra_qsr )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_for, '3D' )   ! air-sea : penetrative sol radiat
278         CASE ( jptra_nsr )        ;   ptrdx(:,:,2:jpk) = 0._wp   ;   ptrdy(:,:,2:jpk) = 0._wp
279                                       CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_for, '2D' )   ! air-sea : non penetr sol radiation
280         CASE ( jptra_bbc )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_bbc, '3D' )   ! bottom bound cond (geoth flux)
281         CASE ( jptra_npc )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_npc, '3D' )   ! non penetr convect adjustment
282         CASE ( jptra_atf )        ;   CALL trd_mxl_zint( ptrdx, ptrdy, jpmxl_atf, '3D' )   ! asselin time filter (last trend)
283                                   !
284                                       CALL trd_mxl( kt, r2dt )                             ! trends: Mixed-layer (output)
285         END SELECT
286         !
287      ENDIF
288      !
289   END SUBROUTINE trd_tra_mng
290
291
292   SUBROUTINE trd_tra_iom( ptrdx, ptrdy, ktrd, kt, Kmm )
293      !!---------------------------------------------------------------------
294      !!                  ***  ROUTINE trd_tra_iom  ***
295      !!
296      !! ** Purpose :   output 3D tracer trends using IOM
297      !!----------------------------------------------------------------------
298      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdx   ! Temperature or U trend
299      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   ptrdy   ! Salinity    or V trend
300      INTEGER                   , INTENT(in   ) ::   ktrd    ! tracer trend index
301      INTEGER                   , INTENT(in   ) ::   kt      ! time step
302      INTEGER                   , INTENT(in   ) ::   Kmm     ! time level index
303      !!
304      INTEGER ::   ji, jj, jk   ! dummy loop indices
305      INTEGER ::   ikbu, ikbv   ! local integers
306      REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   z2dx, z2dy   ! 2D workspace
307      !!----------------------------------------------------------------------
308      !
309!!gm Rq: mask the trends already masked in trd_tra, but lbc_lnk should probably be added
310      !
311      ! Trends evaluated every time step that could go to the standard T file and can be output every ts into a 1ts file if 1ts output is selected
312      SELECT CASE( ktrd )
313      ! This total trend is done every time step
314      CASE( jptra_tot  )   ;   CALL iom_put( "ttrd_tot" , ptrdx )           ! model total trend
315                               CALL iom_put( "strd_tot" , ptrdy )
316      END SELECT
317      !
318      ! These trends are done every second time step. When 1ts output is selected must go different (2ts) file from standard T-file
319      IF( MOD( kt, 2 ) == 0 ) THEN
320         SELECT CASE( ktrd )
321         CASE( jptra_xad  )   ;   CALL iom_put( "ttrd_xad"  , ptrdx )        ! x- horizontal advection
322                                  CALL iom_put( "strd_xad"  , ptrdy )
323         CASE( jptra_yad  )   ;   CALL iom_put( "ttrd_yad"  , ptrdx )        ! y- horizontal advection
324                                  CALL iom_put( "strd_yad"  , ptrdy )
325         CASE( jptra_zad  )   ;   CALL iom_put( "ttrd_zad"  , ptrdx )        ! z- vertical   advection
326                                  CALL iom_put( "strd_zad"  , ptrdy )
327                                  IF( ln_linssh ) THEN                   ! cst volume : adv flux through z=0 surface
328                                     ALLOCATE( z2dx(jpi,jpj), z2dy(jpi,jpj) )
329                                     z2dx(:,:) = ww(:,:,1) * ts(:,:,1,jp_tem,Kmm) / e3t(:,:,1,Kmm)
330                                     z2dy(:,:) = ww(:,:,1) * ts(:,:,1,jp_sal,Kmm) / e3t(:,:,1,Kmm)
331                                     CALL iom_put( "ttrd_sad", z2dx )
332                                     CALL iom_put( "strd_sad", z2dy )
333                                     DEALLOCATE( z2dx, z2dy )
334                                  ENDIF
335         CASE( jptra_totad  ) ;   CALL iom_put( "ttrd_totad", ptrdx )        ! total   advection
336                                  CALL iom_put( "strd_totad", ptrdy )
337         CASE( jptra_ldf  )   ;   CALL iom_put( "ttrd_ldf"  , ptrdx )        ! lateral diffusion
338                                  CALL iom_put( "strd_ldf"  , ptrdy )
339         CASE( jptra_zdf  )   ;   CALL iom_put( "ttrd_zdf"  , ptrdx )        ! vertical diffusion (including Kz contribution)
340                                  CALL iom_put( "strd_zdf"  , ptrdy )
341         CASE( jptra_zdfp )   ;   CALL iom_put( "ttrd_zdfp" , ptrdx )        ! PURE vertical diffusion (no isoneutral contribution)
342                                  CALL iom_put( "strd_zdfp" , ptrdy )
343         CASE( jptra_evd )    ;   CALL iom_put( "ttrd_evd"  , ptrdx )        ! EVD trend (convection)
344                                  CALL iom_put( "strd_evd"  , ptrdy )
345         CASE( jptra_dmp  )   ;   CALL iom_put( "ttrd_dmp"  , ptrdx )        ! internal restoring (damping)
346                                  CALL iom_put( "strd_dmp"  , ptrdy )
347         CASE( jptra_bbl  )   ;   CALL iom_put( "ttrd_bbl"  , ptrdx )        ! bottom boundary layer
348                                  CALL iom_put( "strd_bbl"  , ptrdy )
349         CASE( jptra_npc  )   ;   CALL iom_put( "ttrd_npc"  , ptrdx )        ! static instability mixing
350                                  CALL iom_put( "strd_npc"  , ptrdy )
351         CASE( jptra_bbc  )   ;   CALL iom_put( "ttrd_bbc"  , ptrdx )        ! geothermal heating   (only on temperature)
352         CASE( jptra_nsr  )   ;   CALL iom_put( "ttrd_qns"  , ptrdx(:,:,1) ) ! surface forcing + runoff (ln_rnf=T)
353                                  CALL iom_put( "strd_cdt"  , ptrdy(:,:,1) )        ! output as 2D surface fields
354         CASE( jptra_qsr  )   ;   CALL iom_put( "ttrd_qsr"  , ptrdx )        ! penetrative solar radiat. (only on temperature)
355         END SELECT
356         ! the Asselin filter trend  is also every other time step but needs to be lagged one time step
357         ! Even when 1ts output is selected can go to the same (2ts) file as the trends plotted every even time step.
358      ELSE IF( MOD( kt, 2 ) == 1 ) THEN
359         SELECT CASE( ktrd )
360         CASE( jptra_atf  )   ;   CALL iom_put( "ttrd_atf" , ptrdx )        ! asselin time Filter
361                                  CALL iom_put( "strd_atf" , ptrdy )
362         END SELECT
363      END IF
364      !
365   END SUBROUTINE trd_tra_iom
366
367   !!======================================================================
368END MODULE trdtra
Note: See TracBrowser for help on using the repository browser.