source: branches/2011/DEV_r2739_STFC_dCSE/NEMOGCM/NEMO/OPA_SRC/DIA/diaptr.F90 @ 3849

Last change on this file since 3849 was 3849, checked in by trackstand2, 8 years ago

Merge branch 'partitioner'

  • Property svn:keywords set to Id
File size: 44.6 KB
Line 
1MODULE diaptr
2   !!======================================================================
3   !!                       ***  MODULE  diaptr  ***
4   !! Ocean physics:  Computes meridonal transports and zonal means
5   !!=====================================================================
6   !! History :  1.0  ! 2003-09  (C. Talandier, G. Madec)  Original code
7   !!            2.0  ! 2006-01  (A. Biastoch)  Allow sub-basins computation
8   !!            3.2  ! 2010-03  (O. Marti, S. Flavoni) Add fields
9   !!            3.3  ! 2010-10  (G. Madec)  dynamical allocation
10   !!----------------------------------------------------------------------
11
12   !!----------------------------------------------------------------------
13   !!   dia_ptr      : Poleward Transport Diagnostics module
14   !!   dia_ptr_init : Initialization, namelist read
15   !!   dia_ptr_wri  : Output of poleward fluxes
16   !!   ptr_vjk      : "zonal" sum computation of a "meridional" flux array
17   !!   ptr_tjk      : "zonal" mean computation of a tracer field
18   !!   ptr_vj       : "zonal" and vertical sum computation of a "meridional" flux array
19   !!                   (Generic interface to ptr_vj_3d, ptr_vj_2d)
20   !!----------------------------------------------------------------------
21   USE oce              ! ocean dynamics and active tracers
22   USE dom_oce          ! ocean space and time domain
23   USE phycst           ! physical constants
24   USE ldftra_oce       ! ocean active tracers: lateral physics
25   USE dianam           !
26   USE iom              ! IOM library
27   USE ioipsl           ! IO-IPSL library
28   USE in_out_manager   ! I/O manager
29   USE lib_mpp          ! MPP library
30   USE lbclnk           ! lateral boundary condition - processor exchanges
31
32   IMPLICIT NONE
33   PRIVATE
34
35   INTERFACE ptr_vj
36      MODULE PROCEDURE ptr_vj_3d, ptr_vj_2d
37   END INTERFACE
38
39   PUBLIC   dia_ptr_init   ! call in opa module
40   PUBLIC   dia_ptr        ! call in step module
41   PUBLIC   ptr_vj         ! call by tra_ldf & tra_adv routines
42   PUBLIC   ptr_vjk        ! call by tra_ldf & tra_adv routines
43
44   !                                           !!** namelist  namptr  **
45   LOGICAL , PUBLIC ::   ln_diaptr  = .FALSE.   !: Poleward transport flag (T) or not (F)
46   LOGICAL , PUBLIC ::   ln_subbas  = .FALSE.   !: Atlantic/Pacific/Indian basins calculation
47   LOGICAL , PUBLIC ::   ln_diaznl  = .FALSE.   !: Add zonal means and meridional stream functions
48   LOGICAL , PUBLIC ::   ln_ptrcomp = .FALSE.   !: Add decomposition : overturning (and gyre, soon ...)
49   INTEGER , PUBLIC ::   nn_fptr    = 15        !: frequency of ptr computation  [time step]
50   INTEGER , PUBLIC ::   nn_fwri    = 15        !: frequency of ptr outputs      [time step]
51
52   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:) ::   htr_adv, htr_ldf, htr_ove   !: Heat TRansports (adv, diff, overturn.)
53   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:) ::   str_adv, str_ldf, str_ove   !: Salt TRansports (adv, diff, overturn.)
54   
55   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   btmsk                  ! T-point basin interior masks
56   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   btm30                  ! mask out Southern Ocean (=0 south of 30°S)
57   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   htr  , str             ! adv heat and salt transports (approx)
58   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   tn_jk, sn_jk , v_msf   ! i-mean T and S, j-Stream-Function
59   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sjk  , r1_sjk          ! i-mean i-k-surface and its inverse       
60   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   htr_eiv, str_eiv       ! bolus adv heat ans salt transports ('key_diaeiv')
61   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   v_msf_eiv              ! bolus j-streamfuction              ('key_diaeiv')
62
63
64   INTEGER ::   niter       !
65   INTEGER ::   nidom_ptr   !
66   INTEGER ::   numptr      ! logical unit for Poleward TRansports
67   INTEGER ::   nptr        ! = 1 (ln_subbas=F) or = 5 (glo, atl, pac, ind, ipc) (ln_subbas=T)
68
69   REAL(wp) ::   rc_sv    = 1.e-6_wp   ! conversion from m3/s to Sverdrup
70   REAL(wp) ::   rc_pwatt = 1.e-15_wp  ! conversion from W    to PW (further x rau0 x Cp)
71   REAL(wp) ::   rc_ggram = 1.e-6_wp   ! conversion from g    to Pg
72
73   REAL(wp), TARGET, DIMENSION(:),   ALLOCATABLE, SAVE :: p_fval1d
74   REAL(wp), TARGET, DIMENSION(:,:), ALLOCATABLE, SAVE :: p_fval2d
75
76   !! Integer, 1D workspace arrays. Not common enough to be implemented in
77   !! wrk_nemo module.
78   INTEGER, ALLOCATABLE, SAVE, DIMENSION(:) :: ndex  , ndex_atl     , ndex_pac     , ndex_ind     , ndex_ipc
79   INTEGER, ALLOCATABLE, SAVE, DIMENSION(:) ::         ndex_atl_30  , ndex_pac_30  , ndex_ind_30  , ndex_ipc_30
80   INTEGER, ALLOCATABLE, SAVE, DIMENSION(:) :: ndex_h, ndex_h_atl_30, ndex_h_pac_30, ndex_h_ind_30, ndex_h_ipc_30
81
82   !! * Control permutation of array indices
83#  include "oce_ftrans.h90"
84#  include "dom_oce_ftrans.h90"
85#  include "ldftra_oce_ftrans.h90"
86
87   !! * Substitutions
88#  include "domzgr_substitute.h90"
89#  include "vectopt_loop_substitute.h90"
90   !!----------------------------------------------------------------------
91   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
92   !! $Id$
93   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
94   !!----------------------------------------------------------------------
95CONTAINS
96
97   FUNCTION dia_ptr_alloc()
98      !!----------------------------------------------------------------------
99      !!                    ***  ROUTINE dia_ptr_alloc  ***
100      !!----------------------------------------------------------------------
101      INTEGER               ::   dia_ptr_alloc   ! return value
102      INTEGER, DIMENSION(5) ::   ierr
103      !!----------------------------------------------------------------------
104      ierr(:) = 0
105      !
106      ALLOCATE( btmsk(jpi,jpj,nptr) ,           &
107         &      htr_adv(jpj) , str_adv(jpj) ,   &
108         &      htr_ldf(jpj) , str_ldf(jpj) ,   &
109         &      htr_ove(jpj) , str_ove(jpj),    &
110         &      htr(jpj,nptr) , str(jpj,nptr) , &
111         &      tn_jk(jpj,jpk,nptr) , sn_jk (jpj,jpk,nptr) , v_msf(jpj,jpk,nptr) , &
112         &      sjk  (jpj,jpk,nptr) , r1_sjk(jpj,jpk,nptr) , STAT=ierr(1)  )
113         !
114#if defined key_diaeiv
115      ALLOCATE( htr_eiv(jpj,nptr) , str_eiv(jpj,nptr) , &
116         &      v_msf_eiv(jpj,jpk,nptr) , STAT=ierr(2) )
117#endif
118      ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(3))
119      !
120      ALLOCATE(ndex(jpj*jpk),        ndex_atl(jpj*jpk), ndex_pac(jpj*jpk), &
121         &     ndex_ind(jpj*jpk),    ndex_ipc(jpj*jpk),                    &
122         &     ndex_atl_30(jpj*jpk), ndex_pac_30(jpj*jpk), Stat=ierr(4))
123
124      ALLOCATE(ndex_ind_30(jpj*jpk), ndex_ipc_30(jpj*jpk),                   &
125         &     ndex_h(jpj),          ndex_h_atl_30(jpj), ndex_h_pac_30(jpj), &
126         &     ndex_h_ind_30(jpj),   ndex_h_ipc_30(jpj), Stat=ierr(5) )
127         !
128      dia_ptr_alloc = MAXVAL( ierr )
129      IF(lk_mpp)   CALL mpp_sum( dia_ptr_alloc )
130      !
131   END FUNCTION dia_ptr_alloc
132
133
134   FUNCTION ptr_vj_3d( pva )   RESULT ( p_fval )
135      !!----------------------------------------------------------------------
136      !!                    ***  ROUTINE ptr_vj_3d  ***
137      !!
138      !! ** Purpose :   i-k sum computation of a j-flux array
139      !!
140      !! ** Method  : - i-k sum of pva using the interior 2D vmask (vmask_i).
141      !!              pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
142      !!
143      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
144      !!----------------------------------------------------------------------
145!FTRANS pva :I :I :z
146!! DCSE_NEMO: work around deficiency in ftrans
147!     REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pva   ! mask flux array at V-point
148      REAL(wp) , INTENT(in)    ::   pva(jpi,jpj,jpk)           ! mask flux array at V-point
149      !!
150      INTEGER                  ::   ji, jj, jk   ! dummy loop arguments
151      INTEGER                  ::   ijpj         ! ???
152      REAL(wp), POINTER, DIMENSION(:) :: p_fval  ! function value
153      !!--------------------------------------------------------------------
154      !
155      p_fval => p_fval1d
156
157      ijpj = jpj
158      p_fval(:) = 0._wp
159#if defined key_z_first
160      DO jj = 2, jpjm1
161         DO ji = 2, jpim1
162            DO jk = 1, jpkm1
163#else
164      DO jk = 1, jpkm1
165         DO jj = 2, jpjm1
166            DO ji = fs_2, fs_jpim1   ! Vector opt.
167#endif
168               p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj) 
169            END DO
170         END DO
171      END DO
172#if defined key_mpp_mpi
173      IF(lk_mpp)   CALL mpp_sum( p_fval, ijpj, ncomm_znl)
174#endif
175      !
176   END FUNCTION ptr_vj_3d
177
178!FTRANS CLEAR
179   !! * Re-instate directives to control permutation of array indices
180#  include "oce_ftrans.h90"
181#  include "dom_oce_ftrans.h90"
182#  include "ldftra_oce_ftrans.h90"
183
184   FUNCTION ptr_vj_2d( pva )   RESULT ( p_fval )
185      !!----------------------------------------------------------------------
186      !!                    ***  ROUTINE ptr_vj_2d  ***
187      !!
188      !! ** Purpose :   "zonal" and vertical sum computation of a i-flux array
189      !!
190      !! ** Method  : - i-k sum of pva using the interior 2D vmask (vmask_i).
191      !!      pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
192      !!
193      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
194      !!----------------------------------------------------------------------
195      IMPLICIT none
196      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pva   ! mask flux array at V-point
197      !!
198      INTEGER                  ::   ji,jj       ! dummy loop arguments
199      INTEGER                  ::   ijpj        ! ???
200      REAL(wp), POINTER, DIMENSION(:) :: p_fval ! function value
201      !!--------------------------------------------------------------------
202      !
203      p_fval => p_fval1d
204
205      ijpj = jpj
206      p_fval(:) = 0._wp
207      DO jj = 2, jpjm1
208         DO ji = nldi, nlei   ! No vector optimisation here. Better use a mask ?
209            p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj)
210         END DO
211      END DO
212#if defined key_mpp_mpi
213      CALL mpp_sum( p_fval, ijpj, ncomm_znl )
214#endif
215      !
216   END FUNCTION ptr_vj_2d
217
218
219   FUNCTION ptr_vjk( pva, pmsk )   RESULT ( p_fval )
220      !!----------------------------------------------------------------------
221      !!                    ***  ROUTINE ptr_vjk  ***
222      !!
223      !! ** Purpose :   i-sum computation of a j-velocity array
224      !!
225      !! ** Method  : - i-sum of pva using the interior 2D vmask (vmask_i).
226      !!              pva is supposed to be a masked flux (i.e. * vmask)
227      !!
228      !! ** Action  : - p_fval: i-mean poleward flux of pva
229      !!----------------------------------------------------------------------
230#if defined key_mpp_mpi
231      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
232      USE wrk_nemo, ONLY:   zwork => wrk_1d_1
233#endif
234      !!
235      IMPLICIT none
236!FTRANS pva :I :I :z
237!! DCSE_NEMO: work around a deficiency in ftrans
238!     REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk)           ::   pva    ! mask flux array at V-point
239      REAL(wp) , INTENT(in)             ::   pva(jpi,jpj,jpk)             ! mask flux array at V-point
240      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)    , OPTIONAL ::   pmsk   ! Optional 2D basin mask
241      !!
242      INTEGER                           :: ji, jj, jk ! dummy loop arguments
243      REAL(wp), POINTER, DIMENSION(:,:) :: p_fval     ! return function value
244#if defined key_mpp_mpi
245      INTEGER, DIMENSION(1) ::   ish
246      INTEGER, DIMENSION(2) ::   ish2
247      INTEGER               ::   ijpjjpk
248#endif
249      !!--------------------------------------------------------------------
250      !
251#if defined key_mpp_mpi
252      IF( wrk_in_use(1, 1) ) THEN
253         CALL ctl_stop('ptr_vjk: ERROR - requested workspace array is unavailable')   ;   RETURN
254      END IF
255#endif
256
257      p_fval => p_fval2d
258
259      p_fval(:,:) = 0._wp
260
261      !
262      IF( PRESENT( pmsk ) ) THEN
263         DO jk = 1, jpkm1
264            DO jj = 2, jpjm1
265!!gm here, use of tmask_i  ==> no need of loop over nldi, nlei....
266               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
267                  p_fval(jj,jk) = p_fval(jj,jk) + pva(ji,jj,jk) * e1v(ji,jj) * fse3v(ji,jj,jk) * pmsk(ji,jj)
268               END DO
269            END DO
270         END DO
271      ELSE
272         DO jk = 1, jpkm1
273            DO jj = 2, jpjm1
274               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
275                  p_fval(jj,jk) = p_fval(jj,jk) + pva(ji,jj,jk) * e1v(ji,jj) * fse3v(ji,jj,jk) * tmask_i(ji,jj)
276               END DO
277            END DO
278         END DO
279      END IF
280      !
281#if defined key_mpp_mpi
282      ijpjjpk = jpj*jpk
283      ish(1) = ijpjjpk  ;   ish2(1) = jpj   ;   ish2(2) = jpk
284      zwork(1:ijpjjpk) = RESHAPE( p_fval, ish )
285      CALL mpp_sum( zwork, ijpjjpk, ncomm_znl )
286      p_fval(:,:) = RESHAPE( zwork, ish2 )
287#endif
288      !
289#if defined key_mpp_mpi
290      IF( wrk_not_released(1, 1) )   CALL ctl_stop('ptr_vjk: ERROR - failed to release workspace array')
291#endif
292      !
293   END FUNCTION ptr_vjk
294
295!FTRANS CLEAR
296   !! * Re-instate directives to control permutation of array indices
297#  include "oce_ftrans.h90"
298#  include "dom_oce_ftrans.h90"
299#  include "ldftra_oce_ftrans.h90"
300
301   FUNCTION ptr_tjk( pta, pmsk )   RESULT ( p_fval )
302      !!----------------------------------------------------------------------
303      !!                    ***  ROUTINE ptr_tjk  ***
304      !!
305      !! ** Purpose :   i-sum computation of e1t*e3t * a tracer field
306      !!
307      !! ** Method  : - i-sum of mj(pta) using tmask
308      !!
309      !! ** Action  : - p_fval: i-sum of e1t*e3t*pta
310      !!----------------------------------------------------------------------
311#if defined key_mpp_mpi
312      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
313      USE wrk_nemo, ONLY:   zwork => wrk_1d_1
314#endif
315      !!
316!FTRANS pta :I :I :z
317!! DCSE_NEMO: work around a deficiency in ftrans
318!     REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pta    ! tracer flux array at T-point
319      REAL(wp) , INTENT(in)             :: pta(jpi,jpj,jpk)     ! tracer flux array at T-point
320      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)     ::   pmsk   ! Optional 2D basin mask
321      !!
322      INTEGER                           :: ji, jj, jk   ! dummy loop arguments
323      REAL(wp), POINTER, DIMENSION(:,:) :: p_fval       ! return function value
324#if defined key_mpp_mpi
325      INTEGER, DIMENSION(1) ::   ish
326      INTEGER, DIMENSION(2) ::   ish2
327      INTEGER               ::   ijpjjpk
328#endif
329      !!--------------------------------------------------------------------
330      !
331#if defined key_mpp_mpi
332      IF( wrk_in_use(1, 1) ) THEN
333         CALL ctl_stop('ptr_tjk: requested workspace array unavailable')   ;   RETURN
334      ENDIF
335#endif
336
337      p_fval => p_fval2d
338
339      p_fval(:,:) = 0._wp
340#if defined key_z_first
341      DO jj = 2, jpjm1
342         DO ji =  nldi, nlei
343            DO jk = 1, jpkm1
344#else
345      DO jk = 1, jpkm1
346         DO jj = 2, jpjm1
347            DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
348#endif
349               p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * e1t(ji,jj) * fse3t(ji,jj,jk) * pmsk(ji,jj)
350            END DO
351         END DO
352      END DO
353#if defined key_mpp_mpi
354      ijpjjpk = jpj*jpk
355      ish(1) = jpj*jpk   ;   ish2(1) = jpj   ;   ish2(2) = jpk
356      zwork(1:ijpjjpk)= RESHAPE( p_fval, ish )
357      CALL mpp_sum( zwork, ijpjjpk, ncomm_znl )
358      p_fval(:,:)= RESHAPE( zwork, ish2 )
359#endif
360      !
361#if defined key_mpp_mpi
362      IF( wrk_not_released(1, 1) )   CALL ctl_stop('ptr_tjk: failed to release workspace array')
363#endif
364      !   
365   END FUNCTION ptr_tjk
366
367!FTRANS CLEAR
368   !! * Re-instate directives to control permutation of array indices
369#  include "oce_ftrans.h90"
370#  include "dom_oce_ftrans.h90"
371#  include "ldftra_oce_ftrans.h90"
372
373   SUBROUTINE dia_ptr( kt )
374      !!----------------------------------------------------------------------
375      !!                  ***  ROUTINE dia_ptr  ***
376      !!----------------------------------------------------------------------
377      USE oce,     vt  =>   ua   ! use ua as workspace
378!! DCSE_NEMO: see ticket 873
379      USE oce,     vs  =>   va   ! use va as workspace
380!! DCSE_NEMO: ua, va are re-named, so need additional directives
381!FTRANS vt vs :I :I :z
382      IMPLICIT none
383      !!
384      INTEGER, INTENT(in) ::   kt   ! ocean time step index
385      !
386      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
387      REAL(wp) ::   zv               ! local scalar
388      !!----------------------------------------------------------------------
389      !
390      IF( kt == nit000 .OR. MOD( kt, nn_fptr ) == 0 )   THEN
391         !
392         IF( MOD( kt, nn_fptr ) == 0 ) THEN 
393            !
394            IF( ln_diaznl ) THEN               ! i-mean temperature and salinity
395               DO jn = 1, nptr
396                  tn_jk(:,:,jn) = ptr_tjk( tn(:,:,:), btmsk(:,:,jn) ) * r1_sjk(:,:,jn)
397               END DO
398            ENDIF
399            !
400            !                          ! horizontal integral and vertical dz
401            !                                ! eulerian velocity
402            v_msf(:,:,1) = ptr_vjk( vn(:,:,:) ) 
403            DO jn = 2, nptr
404               v_msf(:,:,jn) = ptr_vjk( vn(:,:,:), btmsk(:,:,jn)*btm30(:,:) ) 
405            END DO
406#if defined key_diaeiv
407            DO jn = 1, nptr                  ! bolus velocity
408               v_msf_eiv(:,:,jn) = ptr_vjk( v_eiv(:,:,:), btmsk(:,:,jn) )   ! here no btm30 for MSFeiv
409            END DO
410            !                                ! add bolus stream-function to the eulerian one
411            v_msf(:,:,:) = v_msf(:,:,:) + v_msf_eiv(:,:,:)
412#endif
413            !
414            !                          ! Transports
415            !                                ! local heat & salt transports at T-points  ( tn*mj[vn+v_eiv] )
416            vt(:,:,jpk) = 0._wp   ;   vs(:,:,jpk) = 0._wp
417
418            DO jk= 1, jpkm1
419               DO jj = 2, jpj
420                  DO ji = 1, jpi
421#if defined key_diaeiv 
422                     zv = ( vn(ji,jj,jk) + vn(ji,jj-1,jk) + v_eiv(ji,jj,jk) + v_eiv(ji,jj-1,jk) ) * 0.5_wp
423#else
424                     zv = ( vn(ji,jj,jk) + vn(ji,jj-1,jk) ) * 0.5_wp
425#endif
426                     vt(:,jj,jk) = zv * tn(:,jj,jk)
427                     vs(:,jj,jk) = zv * sn(:,jj,jk)
428                  END DO
429               END DO
430            END DO
431!!gm useless as overlap areas are not used in ptr_vjk
432            CALL lbc_lnk( vs, 'V', -1. )   ;   CALL lbc_lnk( vt, 'V', -1. )
433!!gm
434            !                                ! heat & salt advective transports (approximation)
435            htr(:,1) = SUM( ptr_vjk( vt(:,:,:) ) , 2 ) * rc_pwatt   ! SUM over jk + conversion
436            str(:,1) = SUM( ptr_vjk( vs(:,:,:) ) , 2 ) * rc_ggram
437            DO jn = 2, nptr 
438               htr(:,jn) = SUM( ptr_vjk( vt(:,:,:), btmsk(:,:,jn)*btm30(:,:) ) , 2 ) * rc_pwatt   ! mask Southern Ocean
439               str(:,jn) = SUM( ptr_vjk( vs(:,:,:), btmsk(:,:,jn)*btm30(:,:) ) , 2 ) * rc_ggram   ! mask Southern Ocean
440            END DO
441
442            IF( ln_ptrcomp ) THEN            ! overturning transport
443               htr_ove(:) = SUM( v_msf(:,:,1) * tn_jk(:,:,1), 2 ) * rc_pwatt   ! SUM over jk + conversion
444               str_ove(:) = SUM( v_msf(:,:,1) * sn_jk(:,:,1), 2 ) * rc_ggram
445            END IF
446            !                                ! Advective and diffusive transport
447            htr_adv(:) = htr_adv(:) * rc_pwatt        ! these are computed in tra_adv... and tra_ldf... routines
448            htr_ldf(:) = htr_ldf(:) * rc_pwatt        ! here just the conversion in PW and Gg
449            str_adv(:) = str_adv(:) * rc_ggram
450            str_ldf(:) = str_ldf(:) * rc_ggram
451
452#if defined key_diaeiv
453            DO jn = 1, nptr                  ! Bolus component
454               htr_eiv(:,jn) = SUM( v_msf_eiv(:,:,jn) * tn_jk(:,:,jn), 2 ) * rc_pwatt   ! SUM over jk
455               str_eiv(:,jn) = SUM( v_msf_eiv(:,:,jn) * sn_jk(:,:,jn), 2 ) * rc_ggram   ! SUM over jk
456            END DO
457#endif
458            !                                ! "Meridional" Stream-Function
459            DO jn = 1, nptr
460               DO jk = 2, jpk 
461                  v_msf    (:,jk,jn) = v_msf    (:,jk-1,jn) + v_msf    (:,jk,jn)       ! Eulerian j-Stream-Function
462#if defined key_diaeiv
463                  v_msf_eiv(:,jk,jn) = v_msf_eiv(:,jk-1,jn) + v_msf_eiv(:,jk,jn)       ! Bolus    j-Stream-Function
464
465#endif
466               END DO
467            END DO
468            v_msf    (:,:,:) = v_msf    (:,:,:) * rc_sv       ! converte in Sverdrups
469#if defined key_diaeiv
470            v_msf_eiv(:,:,:) = v_msf_eiv(:,:,:) * rc_sv
471#endif
472         ENDIF
473         !
474         CALL dia_ptr_wri( kt )                        ! outputs
475         !
476      ENDIF
477      !
478      IF( kt == nitend )   CALL histclo( numptr )      ! Close the file
479      !
480   END SUBROUTINE dia_ptr
481
482!FTRANS CLEAR
483   !! * Re-instate directives to control permutation of array indices
484#  include "oce_ftrans.h90"
485#  include "dom_oce_ftrans.h90"
486#  include "ldftra_oce_ftrans.h90"
487
488   SUBROUTINE dia_ptr_init
489      !!----------------------------------------------------------------------
490      !!                  ***  ROUTINE dia_ptr_init  ***
491      !!                   
492      !! ** Purpose :   Initialization, namelist read
493      !!----------------------------------------------------------------------
494      INTEGER ::   jn           ! dummy loop indices
495      INTEGER ::   inum, ierr   ! local integers
496#if defined key_mpp_mpi
497      INTEGER, DIMENSION(1) :: iglo, iloc, iabsf, iabsl, ihals, ihale, idid
498#endif
499      !!
500      NAMELIST/namptr/ ln_diaptr, ln_diaznl, ln_subbas, ln_ptrcomp, nn_fptr, nn_fwri
501      !!----------------------------------------------------------------------
502
503      !                                      ! allocate dia_ptr arrays
504      IF( dia_ptr_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'lim_sbc_init : unable to allocate arrays' )
505
506      REWIND( numnam )                 ! Read Namelist namptr : poleward transport parameters
507      READ  ( numnam, namptr )
508
509      IF(lwp) THEN                     ! Control print
510         WRITE(numout,*)
511         WRITE(numout,*) 'dia_ptr_init : poleward transport and msf initialization'
512         WRITE(numout,*) '~~~~~~~~~~~~'
513         WRITE(numout,*) '   Namelist namptr : set ptr parameters'
514         WRITE(numout,*) '      Poleward heat & salt transport (T) or not (F)      ln_diaptr  = ', ln_diaptr
515         WRITE(numout,*) '      Overturning heat & salt transport                  ln_ptrcomp = ', ln_ptrcomp
516         WRITE(numout,*) '      T & S zonal mean and meridional stream function    ln_diaznl  = ', ln_diaznl 
517         WRITE(numout,*) '      Global (F) or glo/Atl/Pac/Ind/Indo-Pac basins      ln_subbas  = ', ln_subbas
518         WRITE(numout,*) '      Frequency of computation                           nn_fptr    = ', nn_fptr
519         WRITE(numout,*) '      Frequency of outputs                               nn_fwri    = ', nn_fwri
520      ENDIF
521
522      IF( ln_subbas ) THEN   ;   nptr = 5       ! Global, Atlantic, Pacific, Indian, Indo-Pacific
523      ELSE                   ;   nptr = 1       ! Global only
524      ENDIF
525
526      rc_pwatt = rc_pwatt * rau0 * rcp          ! conversion from K.s-1 to PetaWatt
527
528      IF( .NOT. ln_diaptr ) THEN       ! diaptr not used
529        RETURN
530      ENDIF
531     
532      IF( lk_mpp )   CALL mpp_ini_znl( numout )     ! Define MPI communicator for zonal sum
533
534      IF( ln_subbas ) THEN                ! load sub-basin mask
535         CALL iom_open( 'subbasins', inum )
536         CALL iom_get( inum, jpdom_data, 'atlmsk', btmsk(:,:,2) )   ! Atlantic basin
537         CALL iom_get( inum, jpdom_data, 'pacmsk', btmsk(:,:,3) )   ! Pacific  basin
538         CALL iom_get( inum, jpdom_data, 'indmsk', btmsk(:,:,4) )   ! Indian   basin
539         CALL iom_close( inum )
540         btmsk(:,:,5) = MAX ( btmsk(:,:,3), btmsk(:,:,4) )          ! Indo-Pacific basin
541         WHERE( gphit(:,:) < -30._wp)   ;   btm30(:,:) = 0._wp      ! mask out Southern Ocean
542#if defined key_z_first
543         ELSEWHERE                     ;   btm30(:,:) = tmask_1(:,:)
544#else
545         ELSEWHERE                     ;   btm30(:,:) = tmask(:,:,1)
546#endif
547         END WHERE
548      ENDIF
549      btmsk(:,:,1) = tmask_i(:,:)                                   ! global ocean
550     
551      DO jn = 1, nptr
552         btmsk(:,:,jn) = btmsk(:,:,jn) * tmask_i(:,:)               ! interior domain only
553      END DO
554     
555      IF( lk_vvl )   CALL ctl_stop( 'diaptr: error in vvl case as constant i-mean surface is used' )
556
557      !                                   ! i-sum of e1v*e3v surface and its inverse
558      DO jn = 1, nptr
559         sjk(:,:,jn) = ptr_tjk( tmask(:,:,:), btmsk(:,:,jn) )
560         r1_sjk(:,:,jn) = 0._wp
561         WHERE( sjk(:,:,jn) /= 0._wp )   r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn)
562      END DO
563
564#if defined key_mpp_mpi 
565      iglo (1) = jpjglo                   ! MPP case using MPI  ('key_mpp_mpi')
566      iloc (1) = nlcj
567      iabsf(1) = njmppt(narea)
568      iabsl(:) = iabsf(:) + iloc(:) - 1
569      ihals(1) = nldj - 1
570      ihale(1) = nlcj - nlej
571      idid (1) = 2
572      CALL flio_dom_set( jpnj, nproc/jpni, idid, iglo, iloc, iabsf, iabsl, ihals, ihale, 'BOX', nidom_ptr )
573#else
574      nidom_ptr = FLIO_DOM_NONE
575#endif
576      !
577   END SUBROUTINE dia_ptr_init
578
579
580   SUBROUTINE dia_ptr_wri( kt )
581      !!---------------------------------------------------------------------
582      !!                ***  ROUTINE dia_ptr_wri  ***
583      !!
584      !! ** Purpose :   output of poleward fluxes
585      !!
586      !! ** Method  :   NetCDF file
587      !!----------------------------------------------------------------------
588      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
589      USE wrk_nemo, ONLY:   zphi => wrk_1d_1, zfoo => wrk_1d_2    ! 1D workspace
590      USE wrk_nemo, ONLY:   z_1  => wrk_2d_1                      ! 2D      -
591      !!
592      INTEGER, INTENT(in) ::   kt   ! ocean time-step index
593      !!
594      INTEGER, SAVE ::   nhoridz, ndepidzt, ndepidzw
595      INTEGER, SAVE ::   ndim  , ndim_atl     , ndim_pac     , ndim_ind     , ndim_ipc
596      INTEGER, SAVE ::           ndim_atl_30  , ndim_pac_30  , ndim_ind_30  , ndim_ipc_30
597      INTEGER, SAVE ::   ndim_h, ndim_h_atl_30, ndim_h_pac_30, ndim_h_ind_30, ndim_h_ipc_30
598      !!
599      CHARACTER (len=40) ::   clhstnam, clop, clop_once, cl_comment   ! temporary names
600      INTEGER            ::   iline, it, itmod, ji, jj, jk            !
601#if defined key_iomput
602      INTEGER            ::   inum                                    ! temporary logical unit
603#endif
604      REAL(wp)           ::   zsto, zout, zdt, zjulian                ! temporary scalars
605      !!----------------------------------------------------------------------
606
607      IF( wrk_in_use(1, 1,2) .OR. wrk_in_use(2, 1) ) THEN
608         CALL ctl_stop('dia_ptr_wri: requested workspace arrays unavailable')   ;   RETURN
609      ENDIF
610
611      ! define time axis
612      it    = kt / nn_fptr
613      itmod = kt - nit000 + 1
614     
615      ! Initialization
616      ! --------------
617      IF( kt == nit000 ) THEN
618         niter = ( nit000 - 1 ) / nn_fptr
619         zdt = rdt
620         IF( nacc == 1 )   zdt = rdtmin
621         !
622         IF(lwp) THEN
623            WRITE(numout,*)
624            WRITE(numout,*) 'dia_ptr_wri : poleward transport and msf writing: initialization , niter = ', niter
625            WRITE(numout,*) '~~~~~~~~~~~~'
626         ENDIF
627
628         ! Reference latitude (used in plots)
629         ! ------------------
630         !                                           ! =======================
631         IF( cp_cfg == "orca" ) THEN                 !   ORCA configurations
632            !                                        ! =======================
633            IF( jp_cfg == 05  )   iline = 192   ! i-line that passes near the North Pole
634            IF( jp_cfg == 025 )   iline = 384   ! i-line that passes near the North Pole
635            IF( jp_cfg == 1   )   iline =  96   ! i-line that passes near the North Pole
636            IF( jp_cfg == 2   )   iline =  48   ! i-line that passes near the North Pole
637            IF( jp_cfg == 4   )   iline =  24   ! i-line that passes near the North Pole
638            zphi(1:jpj) = 0._wp
639            DO ji = mi0(iline), mi1(iline) 
640               zphi(1:jpj) = gphiv(ji,:)         ! if iline is in the local domain
641               ! Correct highest latitude for some configurations - will work if domain is parallelized in J ?
642               IF( jp_cfg == 05 ) THEN
643                  DO jj = mj0(jpjdta), mj1(jpjdta) 
644                     zphi( jj ) = zphi(mj0(jpjdta-1)) + ( zphi(mj0(jpjdta-1))-zphi(mj0(jpjdta-2)) ) * 0.5_wp
645                     zphi( jj ) = MIN( zphi(jj), 90._wp )
646                  END DO
647               END IF
648               IF( jp_cfg == 1 .OR. jp_cfg == 2 .OR. jp_cfg == 4 ) THEN
649                  DO jj = mj0(jpjdta-1), mj1(jpjdta-1) 
650                     zphi( jj ) = 88.5_wp
651                  END DO
652                  DO jj = mj0(jpjdta  ), mj1(jpjdta  ) 
653                     zphi( jj ) = 89.5_wp
654                  END DO
655               END IF
656            END DO
657            ! provide the correct zphi to all local domains
658#if defined key_mpp_mpi
659            CALL mpp_sum( zphi, jpj, ncomm_znl )       
660#endif
661            !                                        ! =======================
662         ELSE                                        !   OTHER configurations
663            !                                        ! =======================
664            zphi(1:jpj) = gphiv(1,:)             ! assume lat/lon coordinate, select the first i-line
665            !
666         ENDIF
667         !
668         ! Work only on westmost processor (will not work if mppini2 is used)
669#if defined key_mpp_mpi
670         IF( l_znl_root ) THEN 
671#endif
672            !
673            ! OPEN netcdf file
674            ! ----------------
675            ! Define frequency of output and means
676            zsto = nn_fptr * zdt
677            IF( ln_mskland )   THEN    ! put 1.e+20 on land (very expensive!!)
678               clop      = "ave(only(x))"
679               clop_once = "once(only(x))"
680            ELSE                       ! no use of the mask value (require less cpu time)
681               clop      = "ave(x)"       
682               clop_once = "once"
683            ENDIF
684
685            zout = nn_fwri * zdt
686            zfoo(1:jpj) = 0._wp
687
688            CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )  ! Compute julian date from starting date of the run
689            zjulian = zjulian - adatrj                         ! set calendar origin to the beginning of the experiment
690
691#if defined key_iomput
692            ! Requested by IPSL people, use by their postpro...
693            IF(lwp) THEN
694               CALL dia_nam( clhstnam, nn_fwri,' ' )
695               CALL ctl_opn( inum, 'date.file', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
696               WRITE(inum,*) clhstnam
697               CLOSE(inum)
698            ENDIF
699#endif
700
701            CALL dia_nam( clhstnam, nn_fwri, 'diaptr' )
702            IF(lwp)WRITE( numout,*)" Name of diaptr NETCDF file : ", clhstnam
703
704            ! Horizontal grid : zphi()
705            CALL histbeg(clhstnam, 1, zfoo, jpj, zphi,   &
706               1, 1, 1, jpj, niter, zjulian, zdt*nn_fptr, nhoridz, numptr, domain_id=nidom_ptr)
707            ! Vertical grids : gdept_0, gdepw_0
708            CALL histvert( numptr, "deptht", "Vertical T levels",   &
709               &                   "m", jpk, gdept_0, ndepidzt, "down" )
710            CALL histvert( numptr, "depthw", "Vertical W levels",   &
711               &                   "m", jpk, gdepw_0, ndepidzw, "down" )
712            !
713            CALL wheneq ( jpj*jpk, MIN(sjk(:,:,1), 1._wp), 1, 1., ndex  , ndim  )      ! Lat-Depth
714            CALL wheneq ( jpj    , MIN(sjk(:,1,1), 1._wp), 1, 1., ndex_h, ndim_h )     ! Lat
715
716            IF( ln_subbas ) THEN
717               z_1(:,1) = 1._wp
718               WHERE ( gphit(jpi/2,:) < -30._wp )   z_1(:,1) = 0._wp
719               DO jk = 2, jpk
720                  z_1(:,jk) = z_1(:,1)
721               END DO
722               !                       ! Atlantic (jn=2)
723               CALL wheneq ( jpj*jpk, MIN(sjk(:,:,2)         , 1._wp), 1, 1., ndex_atl     , ndim_atl      ) ! Lat-Depth
724               CALL wheneq ( jpj*jpk, MIN(sjk(:,:,2)*z_1(:,:), 1._wp), 1, 1., ndex_atl_30  , ndim_atl_30   ) ! Lat-Depth
725               CALL wheneq ( jpj    , MIN(sjk(:,1,2)*z_1(:,1), 1._wp), 1, 1., ndex_h_atl_30, ndim_h_atl_30 ) ! Lat
726               !                       ! Pacific (jn=3)
727               CALL wheneq ( jpj*jpk, MIN(sjk(:,:,3)         , 1._wp), 1, 1., ndex_pac     , ndim_pac      ) ! Lat-Depth
728               CALL wheneq ( jpj*jpk, MIN(sjk(:,:,3)*z_1(:,:), 1._wp), 1, 1., ndex_pac_30  , ndim_pac_30   ) ! Lat-Depth
729               CALL wheneq ( jpj    , MIN(sjk(:,1,3)*z_1(:,1), 1._wp), 1, 1., ndex_h_pac_30, ndim_h_pac_30 ) ! Lat
730               !                       ! Indian (jn=4)
731               CALL wheneq ( jpj*jpk, MIN(sjk(:,:,4)         , 1._wp), 1, 1., ndex_ind     , ndim_ind      ) ! Lat-Depth
732               CALL wheneq ( jpj*jpk, MIN(sjk(:,:,4)*z_1(:,:), 1._wp), 1, 1., ndex_ind_30  , ndim_ind_30   ) ! Lat-Depth
733               CALL wheneq ( jpj    , MIN(sjk(:,1,4)*z_1(:,1), 1._wp), 1, 1., ndex_h_ind_30, ndim_h_ind_30 ) ! Lat
734               !                       ! Indo-Pacific (jn=5)
735               CALL wheneq ( jpj*jpk, MIN(sjk(:,:,5)         , 1._wp), 1, 1., ndex_ipc     , ndim_ipc      ) ! Lat-Depth
736               CALL wheneq ( jpj*jpk, MIN(sjk(:,:,5)*z_1(:,:), 1._wp), 1, 1., ndex_ipc_30  , ndim_ipc_30   ) ! Lat-Depth
737               CALL wheneq ( jpj    , MIN(sjk(:,1,5)*z_1(:,1), 1._wp), 1, 1., ndex_h_ipc_30, ndim_h_ipc_30 ) ! Lat
738            ENDIF
739            !
740#if defined key_diaeiv
741            cl_comment = ' (Bolus part included)'
742#else
743            cl_comment = '                      '
744#endif
745            IF( ln_diaznl ) THEN             !  Zonal mean T and S
746               CALL histdef( numptr, "zotemglo", "Zonal Mean Temperature","C" ,   &
747                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
748               CALL histdef( numptr, "zosalglo", "Zonal Mean Salinity","PSU"  ,   &
749                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
750
751               CALL histdef( numptr, "zosrfglo", "Zonal Mean Surface","m^2"   ,   &
752                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
753               !
754               IF (ln_subbas) THEN
755                  CALL histdef( numptr, "zotematl", "Zonal Mean Temperature: Atlantic","C" ,   &
756                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
757                  CALL histdef( numptr, "zosalatl", "Zonal Mean Salinity: Atlantic","PSU"  ,   &
758                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
759                  CALL histdef( numptr, "zosrfatl", "Zonal Mean Surface: Atlantic","m^2"   ,   &
760                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
761
762                  CALL histdef( numptr, "zotempac", "Zonal Mean Temperature: Pacific","C"  ,   &
763                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
764                  CALL histdef( numptr, "zosalpac", "Zonal Mean Salinity: Pacific","PSU"   ,   &
765                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
766                  CALL histdef( numptr, "zosrfpac", "Zonal Mean Surface: Pacific","m^2"    ,   &
767                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
768
769                  CALL histdef( numptr, "zotemind", "Zonal Mean Temperature: Indian","C"   ,   &
770                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
771                  CALL histdef( numptr, "zosalind", "Zonal Mean Salinity: Indian","PSU"    ,   &
772                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
773                  CALL histdef( numptr, "zosrfind", "Zonal Mean Surface: Indian","m^2"     ,   &
774                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
775
776                  CALL histdef( numptr, "zotemipc", "Zonal Mean Temperature: Pacific+Indian","C" ,   &
777                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
778                  CALL histdef( numptr, "zosalipc", "Zonal Mean Salinity: Pacific+Indian","PSU"  ,   &
779                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop, zsto, zout )
780                  CALL histdef( numptr, "zosrfipc", "Zonal Mean Surface: Pacific+Indian","m^2"   ,   &
781                     1, jpj, nhoridz, jpk, 1, jpk, ndepidzt, 32, clop_once, zsto, zout )
782               ENDIF
783            ENDIF
784            !
785            !  Meridional Stream-Function (Eulerian and Bolus)
786            CALL histdef( numptr, "zomsfglo", "Meridional Stream-Function: Global"//TRIM(cl_comment),"Sv" ,   &
787               1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
788            IF( ln_subbas .AND. ln_diaznl ) THEN
789               CALL histdef( numptr, "zomsfatl", "Meridional Stream-Function: Atlantic"//TRIM(cl_comment),"Sv" ,   &
790                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
791               CALL histdef( numptr, "zomsfpac", "Meridional Stream-Function: Pacific"//TRIM(cl_comment),"Sv"  ,   &
792                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
793               CALL histdef( numptr, "zomsfind", "Meridional Stream-Function: Indian"//TRIM(cl_comment),"Sv"   ,   &
794                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
795               CALL histdef( numptr, "zomsfipc", "Meridional Stream-Function: Indo-Pacific"//TRIM(cl_comment),"Sv" ,&
796                  1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
797            ENDIF
798            !
799            !  Heat transport
800            CALL histdef( numptr, "sophtadv", "Advective Heat Transport"      ,   &
801               "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
802            CALL histdef( numptr, "sophtldf", "Diffusive Heat Transport"      ,   &
803               "PW",1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
804            IF ( ln_ptrcomp ) THEN
805               CALL histdef( numptr, "sophtove", "Overturning Heat Transport"    ,   &
806                  "PW",1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
807            END IF
808            IF( ln_subbas ) THEN
809               CALL histdef( numptr, "sohtatl", "Heat Transport Atlantic"//TRIM(cl_comment),  &
810                  "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
811               CALL histdef( numptr, "sohtpac", "Heat Transport Pacific"//TRIM(cl_comment) ,  &
812                  "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
813               CALL histdef( numptr, "sohtind", "Heat Transport Indian"//TRIM(cl_comment)  ,  &
814                  "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
815               CALL histdef( numptr, "sohtipc", "Heat Transport Pacific+Indian"//TRIM(cl_comment), &
816                  "PW", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
817            ENDIF
818            !
819            !  Salt transport
820            CALL histdef( numptr, "sopstadv", "Advective Salt Transport"      ,   &
821               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
822            CALL histdef( numptr, "sopstldf", "Diffusive Salt Transport"      ,   &
823               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
824            IF ( ln_ptrcomp ) THEN
825               CALL histdef( numptr, "sopstove", "Overturning Salt Transport"    ,   &
826                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
827            END IF
828#if defined key_diaeiv
829            ! Eddy induced velocity
830            CALL histdef( numptr, "zomsfeiv", "Bolus Meridional Stream-Function: global",   &
831               "Sv"      , 1, jpj, nhoridz, jpk, 1, jpk, ndepidzw, 32, clop, zsto, zout )
832            CALL histdef( numptr, "sophteiv", "Bolus Advective Heat Transport",   &
833               "PW"      , 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
834            CALL histdef( numptr, "sopsteiv", "Bolus Advective Salt Transport",   &
835               "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
836#endif
837            IF( ln_subbas ) THEN
838               CALL histdef( numptr, "sostatl", "Salt Transport Atlantic"//TRIM(cl_comment)      ,  &
839                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
840               CALL histdef( numptr, "sostpac", "Salt Transport Pacific"//TRIM(cl_comment)      ,   &
841                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
842               CALL histdef( numptr, "sostind", "Salt Transport Indian"//TRIM(cl_comment)      ,    &
843                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
844               CALL histdef( numptr, "sostipc", "Salt Transport Pacific+Indian"//TRIM(cl_comment),  &
845                  "Giga g/s", 1, jpj, nhoridz, 1, 1, 1, -99, 32, clop, zsto, zout )
846            ENDIF
847            !
848            CALL histend( numptr )
849            !
850         END IF
851#if defined key_mpp_mpi
852      END IF
853#endif
854
855#if defined key_mpp_mpi
856      IF( MOD( itmod, nn_fptr ) == 0 .AND. l_znl_root ) THEN
857#else
858      IF( MOD( itmod, nn_fptr ) == 0  ) THEN
859#endif
860         niter = niter + 1
861
862         IF( ln_diaznl ) THEN
863            CALL histwrite( numptr, "zosrfglo", niter, sjk  (:,:,1) , ndim, ndex )
864            CALL histwrite( numptr, "zotemglo", niter, tn_jk(:,:,1)  , ndim, ndex )
865            CALL histwrite( numptr, "zosalglo", niter, sn_jk(:,:,1)  , ndim, ndex )
866
867            IF (ln_subbas) THEN
868               CALL histwrite( numptr, "zosrfatl", niter, sjk(:,:,2), ndim_atl, ndex_atl )
869               CALL histwrite( numptr, "zosrfpac", niter, sjk(:,:,3), ndim_pac, ndex_pac )
870               CALL histwrite( numptr, "zosrfind", niter, sjk(:,:,4), ndim_ind, ndex_ind )
871               CALL histwrite( numptr, "zosrfipc", niter, sjk(:,:,5), ndim_ipc, ndex_ipc )
872
873               CALL histwrite( numptr, "zotematl", niter, tn_jk(:,:,2)  , ndim_atl, ndex_atl )
874               CALL histwrite( numptr, "zosalatl", niter, sn_jk(:,:,2)  , ndim_atl, ndex_atl )
875               CALL histwrite( numptr, "zotempac", niter, tn_jk(:,:,3)  , ndim_pac, ndex_pac )
876               CALL histwrite( numptr, "zosalpac", niter, sn_jk(:,:,3)  , ndim_pac, ndex_pac )
877               CALL histwrite( numptr, "zotemind", niter, tn_jk(:,:,4)  , ndim_ind, ndex_ind )
878               CALL histwrite( numptr, "zosalind", niter, sn_jk(:,:,4)  , ndim_ind, ndex_ind )
879               CALL histwrite( numptr, "zotemipc", niter, tn_jk(:,:,5)  , ndim_ipc, ndex_ipc )
880               CALL histwrite( numptr, "zosalipc", niter, sn_jk(:,:,5)  , ndim_ipc, ndex_ipc )
881            END IF
882         ENDIF
883
884         ! overturning outputs:
885         CALL histwrite( numptr, "zomsfglo", niter, v_msf(:,:,1), ndim, ndex )
886         IF( ln_subbas .AND. ln_diaznl ) THEN
887            CALL histwrite( numptr, "zomsfatl", niter, v_msf(:,:,2) , ndim_atl_30, ndex_atl_30 )
888            CALL histwrite( numptr, "zomsfpac", niter, v_msf(:,:,3) , ndim_pac_30, ndex_pac_30 )
889            CALL histwrite( numptr, "zomsfind", niter, v_msf(:,:,4) , ndim_ind_30, ndex_ind_30 )
890            CALL histwrite( numptr, "zomsfipc", niter, v_msf(:,:,5) , ndim_ipc_30, ndex_ipc_30 )
891         ENDIF
892#if defined key_diaeiv
893         CALL histwrite( numptr, "zomsfeiv", niter, v_msf_eiv(:,:,1), ndim  , ndex   )
894#endif
895
896         ! heat transport outputs:
897         IF( ln_subbas ) THEN
898            CALL histwrite( numptr, "sohtatl", niter, htr(:,2)  , ndim_h_atl_30, ndex_h_atl_30 )
899            CALL histwrite( numptr, "sohtpac", niter, htr(:,3)  , ndim_h_pac_30, ndex_h_pac_30 )
900            CALL histwrite( numptr, "sohtind", niter, htr(:,4)  , ndim_h_ind_30, ndex_h_ind_30 )
901            CALL histwrite( numptr, "sohtipc", niter, htr(:,5)  , ndim_h_ipc_30, ndex_h_ipc_30 )
902            CALL histwrite( numptr, "sostatl", niter, str(:,2)  , ndim_h_atl_30, ndex_h_atl_30 )
903            CALL histwrite( numptr, "sostpac", niter, str(:,3)  , ndim_h_pac_30, ndex_h_pac_30 )
904            CALL histwrite( numptr, "sostind", niter, str(:,4)  , ndim_h_ind_30, ndex_h_ind_30 )
905            CALL histwrite( numptr, "sostipc", niter, str(:,5)  , ndim_h_ipc_30, ndex_h_ipc_30 )
906         ENDIF
907
908         CALL histwrite( numptr, "sophtadv", niter, htr_adv     , ndim_h, ndex_h )
909         CALL histwrite( numptr, "sophtldf", niter, htr_ldf     , ndim_h, ndex_h )
910         CALL histwrite( numptr, "sopstadv", niter, str_adv     , ndim_h, ndex_h )
911         CALL histwrite( numptr, "sopstldf", niter, str_ldf     , ndim_h, ndex_h )
912         IF( ln_ptrcomp ) THEN
913            CALL histwrite( numptr, "sopstove", niter, str_ove(:) , ndim_h, ndex_h )
914            CALL histwrite( numptr, "sophtove", niter, htr_ove(:) , ndim_h, ndex_h )
915         ENDIF
916#if defined key_diaeiv
917         CALL histwrite( numptr, "sophteiv", niter, htr_eiv(:,1)  , ndim_h, ndex_h )
918         CALL histwrite( numptr, "sopsteiv", niter, str_eiv(:,1)  , ndim_h, ndex_h )
919#endif
920         !
921      ENDIF
922      !
923      IF( wrk_not_released(1, 1,2) .OR.    &
924          wrk_not_released(2, 1)    )   CALL ctl_stop('dia_ptr_wri: failed to release workspace arrays')
925      !
926  END SUBROUTINE dia_ptr_wri
927
928   !!======================================================================
929END MODULE diaptr
Note: See TracBrowser for help on using the repository browser.