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.
diaptr.F90 in branches/UKMO/v3_6_extra_CMIP6_diagnostics/NEMOGCM/NEMO/OPA_SRC/DIA – NEMO

source: branches/UKMO/v3_6_extra_CMIP6_diagnostics/NEMOGCM/NEMO/OPA_SRC/DIA/diaptr.F90 @ 6652

Last change on this file since 6652 was 6652, checked in by timgraham, 8 years ago

Added heat transport from Bolus advection to diaptr (global and by sub basin).

File size: 23.9 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   !!            3.6  ! 2014-12  (C. Ethe) use of IOM
11   !!----------------------------------------------------------------------
12
13   !!----------------------------------------------------------------------
14   !!   dia_ptr      : Poleward Transport Diagnostics module
15   !!   dia_ptr_init : Initialization, namelist read
16   !!   ptr_sjk      : "zonal" mean computation of a field - tracer or flux array
17   !!   ptr_sj       : "zonal" and vertical sum computation of a "meridional" flux array
18   !!                   (Generic interface to ptr_sj_3d, ptr_sj_2d)
19   !!----------------------------------------------------------------------
20   USE oce              ! ocean dynamics and active tracers
21   USE dom_oce          ! ocean space and time domain
22   USE phycst           ! physical constants
23   USE ldftra_oce 
24   !
25   USE iom              ! IOM library
26   USE in_out_manager   ! I/O manager
27   USE lib_mpp          ! MPP library
28   USE timing           ! preformance summary
29
30   IMPLICIT NONE
31   PRIVATE
32
33   INTERFACE ptr_sj
34      MODULE PROCEDURE ptr_sj_3d, ptr_sj_2d
35   END INTERFACE
36
37   PUBLIC   ptr_sj         ! call by tra_ldf & tra_adv routines
38   PUBLIC   ptr_sjk        !
39   PUBLIC   dia_ptr_init   ! call in step module
40   PUBLIC   dia_ptr        ! call in step module
41   PUBLIC   dia_ptr_ohst_components        ! called from tra_ldf/tra_adv routines
42
43   !                                  !!** namelist  namptr  **
44   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   htr_adv, htr_ldf, htr_eiv   !: Heat TRansports (adv, diff, Bolus.)
45   REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) ::   str_adv, str_ldf, str_eiv   !: Salt TRansports (adv, diff, Bolus.)
46
47   LOGICAL, PUBLIC ::   ln_diaptr   !  Poleward transport flag (T) or not (F)
48   LOGICAL, PUBLIC ::   ln_subbas   !  Atlantic/Pacific/Indian basins calculation
49   INTEGER, PUBLIC ::   nptr        ! = 1 (l_subbas=F) or = 5 (glo, atl, pac, ind, ipc) (l_subbas=T)
50
51   REAL(wp) ::   rc_sv    = 1.e-6_wp   ! conversion from m3/s to Sverdrup
52   REAL(wp) ::   rc_pwatt = 1.e-15_wp  ! conversion from W    to PW (further x rau0 x Cp)
53   REAL(wp) ::   rc_ggram = 1.e-6_wp   ! conversion from g    to Pg
54
55   CHARACTER(len=3), ALLOCATABLE, SAVE, DIMENSION(:)     :: clsubb
56   REAL(wp),         ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk   ! T-point basin interior masks
57   REAL(wp),         ALLOCATABLE, SAVE, DIMENSION(:,:)   :: btm30   ! mask out Southern Ocean (=0 south of 30°S)
58
59   REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:)     :: p_fval1d
60   REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:,:)   :: p_fval2d
61
62
63   !! * Substitutions
64#  include "domzgr_substitute.h90"
65#  include "vectopt_loop_substitute.h90"
66   !!----------------------------------------------------------------------
67   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
68   !! $Id$
69   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
70   !!----------------------------------------------------------------------
71CONTAINS
72
73   SUBROUTINE dia_ptr( pvtr )
74      !!----------------------------------------------------------------------
75      !!                  ***  ROUTINE dia_ptr  ***
76      !!----------------------------------------------------------------------
77      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL ::   pvtr   ! j-effective transport
78      !
79      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
80      REAL(wp) ::   zv, zsfc               ! local scalar
81      REAL(wp), DIMENSION(jpi,jpj)     ::  z2d   ! 2D workspace
82      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  z3d   ! 3D workspace
83      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zmask   ! 3D workspace
84      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) ::  zts   ! 3D workspace
85      CHARACTER( len = 12 )  :: cl1
86      !!----------------------------------------------------------------------
87      !
88      IF( nn_timing == 1 )   CALL timing_start('dia_ptr')
89
90      !
91      IF( PRESENT( pvtr ) ) THEN
92         IF( iom_use("zomsfglo") ) THEN    ! effective MSF
93            z3d(1,:,:) = ptr_sjk( pvtr(:,:,:) )  ! zonal cumulative effective transport
94            DO jk = 2, jpkm1 
95              z3d(1,:,jk) = z3d(1,:,jk-1) + z3d(1,:,jk)   ! effective j-Stream-Function (MSF)
96            END DO
97            DO ji = 1, jpi
98               z3d(ji,:,:) = z3d(1,:,:)
99            ENDDO
100            cl1 = TRIM('zomsf'//clsubb(1) )
101            CALL iom_put( cl1, z3d * rc_sv )
102            DO jn = 2, nptr                                    ! by sub-basins
103               z3d(1,:,:) =  ptr_sjk( pvtr(:,:,:), btmsk(:,:,jn)*btm30(:,:) ) 
104               DO jk = 2, jpkm1 
105                  z3d(1,:,jk) = z3d(1,:,jk-1) + z3d(1,:,jk)    ! effective j-Stream-Function (MSF)
106               END DO
107               DO ji = 1, jpi
108                  z3d(ji,:,:) = z3d(1,:,:)
109               ENDDO
110               cl1 = TRIM('zomsf'//clsubb(jn) )
111               CALL iom_put( cl1, z3d * rc_sv )
112            END DO
113         ENDIF
114         !
115      ELSE
116         !
117         IF( iom_use("zotemglo") ) THEN    ! i-mean i-k-surface
118            DO jk = 1, jpkm1
119               DO jj = 1, jpj
120                  DO ji = 1, jpi
121                     zsfc = e1t(ji,jj) * fse3t(ji,jj,jk)
122                     zmask(ji,jj,jk)      = tmask(ji,jj,jk)      * zsfc
123                     zts(ji,jj,jk,jp_tem) = tsn(ji,jj,jk,jp_tem) * zsfc
124                     zts(ji,jj,jk,jp_sal) = tsn(ji,jj,jk,jp_sal) * zsfc
125                  ENDDO
126               ENDDO
127            ENDDO
128            DO jn = 1, nptr
129               zmask(1,:,:) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) )
130               cl1 = TRIM('zosrf'//clsubb(jn) )
131               CALL iom_put( cl1, zmask )
132               !
133               z3d(1,:,:) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) &
134                  &            / MAX( zmask(1,:,:), 10.e-15 )
135               DO ji = 1, jpi
136                  z3d(ji,:,:) = z3d(1,:,:)
137               ENDDO
138               cl1 = TRIM('zotem'//clsubb(jn) )
139               CALL iom_put( cl1, z3d )
140               !
141               z3d(1,:,:) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) &
142                  &            / MAX( zmask(1,:,:), 10.e-15 )
143               DO ji = 1, jpi
144                  z3d(ji,:,:) = z3d(1,:,:)
145               ENDDO
146               cl1 = TRIM('zosal'//clsubb(jn) )
147               CALL iom_put( cl1, z3d )
148            END DO
149         ENDIF
150         !
151         !                                ! Advective and diffusive heat and salt transport
152         IF( iom_use("sophtadv") .OR. iom_use("sopstadv") ) THEN   
153            z2d(1,:) = htr_adv(:,1) * rc_pwatt        !  (conversion in PW)
154            DO ji = 1, jpi
155               z2d(ji,:) = z2d(1,:)
156            ENDDO
157            cl1 = 'sophtadv'                 
158            CALL iom_put( TRIM(cl1), z2d )
159            z2d(1,:) = str_adv(:,1) * rc_ggram        ! (conversion in Gg)
160            DO ji = 1, jpi
161               z2d(ji,:) = z2d(1,:)
162            ENDDO
163            cl1 = 'sopstadv'
164            CALL iom_put( TRIM(cl1), z2d )
165            IF( ln_subbas ) THEN
166              DO jn=2,nptr
167               z2d(1,:) = htr_adv(:,jn) * rc_pwatt        !  (conversion in PW)
168               DO ji = 1, jpi
169                 z2d(ji,:) = z2d(1,:)
170               ENDDO
171               cl1 = TRIM('sophtadv_'//clsubb(jn))                 
172               CALL iom_put( cl1, z2d )
173               z2d(1,:) = str_adv(:,jn) * rc_ggram        ! (conversion in Gg)
174               DO ji = 1, jpi
175                  z2d(ji,:) = z2d(1,:)
176               ENDDO
177               cl1 = TRIM('sopstadv_'//clsubb(jn))                 
178               CALL iom_put( cl1, z2d )             
179              ENDDO
180            ENDIF
181         ENDIF
182         !
183         IF( iom_use("sophtldf") .OR. iom_use("sopstldf") ) THEN   
184            z2d(1,:) = htr_ldf(:,1) * rc_pwatt        !  (conversion in PW)
185            DO ji = 1, jpi
186               z2d(ji,:) = z2d(1,:)
187            ENDDO
188            cl1 = 'sophtldf'
189            CALL iom_put( TRIM(cl1), z2d )
190            z2d(1,:) = str_ldf(:,1) * rc_ggram        !  (conversion in Gg)
191            DO ji = 1, jpi
192               z2d(ji,:) = z2d(1,:)
193            ENDDO
194            cl1 = 'sopstldf'
195            CALL iom_put( TRIM(cl1), z2d )
196            IF( ln_subbas ) THEN
197              DO jn=2,nptr
198               z2d(1,:) = htr_ldf(:,jn) * rc_pwatt        !  (conversion in PW)
199               DO ji = 1, jpi
200                 z2d(ji,:) = z2d(1,:)
201               ENDDO
202               cl1 = TRIM('sophtldf_'//clsubb(jn))                 
203               CALL iom_put( cl1, z2d )
204               z2d(1,:) = str_ldf(:,jn) * rc_ggram        ! (conversion in Gg)
205               DO ji = 1, jpi
206                  z2d(ji,:) = z2d(1,:)
207               ENDDO
208               cl1 = TRIM('sopstldf_'//clsubb(jn))                 
209               CALL iom_put( cl1, z2d )             
210              ENDDO
211            ENDIF
212         ENDIF
213
214#ifdef key_diaeiv
215         IF(lk_traldf_eiv) THEN
216            IF( iom_use("sophteiv") .OR. iom_use("sopsteiv") ) THEN
217               z2d(1,:) = htr_eiv(:,1) * rc_pwatt        !  (conversion in PW)
218               DO ji = 1, jpi
219                  z2d(ji,:) = z2d(1,:)
220               ENDDO
221               cl1 = 'sophteiv'
222               CALL iom_put( TRIM(cl1), z2d )
223               z2d(1,:) = str_eiv(:,1) * rc_ggram        !  (conversion in Gg)
224               DO ji = 1, jpi
225                  z2d(ji,:) = z2d(1,:)
226               ENDDO
227               cl1 = 'sopsteiv'
228               CALL iom_put( TRIM(cl1), z2d )
229               IF( ln_subbas ) THEN
230                  DO jn=2,nptr
231                     z2d(1,:) = htr_eiv(:,jn) * rc_pwatt        !  (conversion in PW)
232                     DO ji = 1, jpi
233                        z2d(ji,:) = z2d(1,:)
234                     ENDDO
235                     cl1 = TRIM('sophteiv_'//clsubb(jn))                 
236                     CALL iom_put( cl1, z2d )
237                     z2d(1,:) = str_eiv(:,jn) * rc_ggram        ! (conversion in Gg)
238                     DO ji = 1, jpi
239                        z2d(ji,:) = z2d(1,:)
240                     ENDDO
241                     cl1 = TRIM('sopsteiv_'//clsubb(jn)) 
242                     CALL iom_put( cl1, z2d )             
243                  ENDDO
244               ENDIF
245            ENDIF
246         ENDIF
247#endif
248         !
249      ENDIF
250      !
251      IF( nn_timing == 1 )   CALL timing_stop('dia_ptr')
252      !
253   END SUBROUTINE dia_ptr
254
255
256   SUBROUTINE dia_ptr_init
257      !!----------------------------------------------------------------------
258      !!                  ***  ROUTINE dia_ptr_init  ***
259      !!                   
260      !! ** Purpose :   Initialization, namelist read
261      !!----------------------------------------------------------------------
262      INTEGER ::  jn           ! local integers
263      INTEGER ::  inum, ierr   ! local integers
264      INTEGER ::  ios          ! Local integer output status for namelist read
265      !!
266      NAMELIST/namptr/ ln_diaptr, ln_subbas
267      !!----------------------------------------------------------------------
268
269      REWIND( numnam_ref )              ! Namelist namptr in reference namelist : Poleward transport
270      READ  ( numnam_ref, namptr, IOSTAT = ios, ERR = 901)
271901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namptr in reference namelist', lwp )
272
273      REWIND( numnam_cfg )              ! Namelist namptr in configuration namelist : Poleward transport
274      READ  ( numnam_cfg, namptr, IOSTAT = ios, ERR = 902 )
275902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namptr in configuration namelist', lwp )
276      IF(lwm) WRITE ( numond, namptr )
277
278      IF(lwp) THEN                     ! Control print
279         WRITE(numout,*)
280         WRITE(numout,*) 'dia_ptr_init : poleward transport and msf initialization'
281         WRITE(numout,*) '~~~~~~~~~~~~'
282         WRITE(numout,*) '   Namelist namptr : set ptr parameters'
283         WRITE(numout,*) '      Poleward heat & salt transport (T) or not (F)      ln_diaptr  = ', ln_diaptr
284         WRITE(numout,*) '      Global (F) or glo/Atl/Pac/Ind/Indo-Pac basins      ln_subbas  = ', ln_subbas
285      ENDIF
286
287      IF( ln_diaptr ) THEN 
288         !
289         IF( ln_subbas ) THEN
290            nptr = 5            ! Global, Atlantic, Pacific, Indian, Indo-Pacific
291            ALLOCATE( clsubb(nptr) )
292            clsubb(1) = 'glo' ;  clsubb(2) = 'atl'  ;  clsubb(3) = 'pac'  ;  clsubb(4) = 'ind'  ;  clsubb(5) = 'ipc'
293         ELSE               
294            nptr = 1       ! Global only
295            ALLOCATE( clsubb(nptr) )
296            clsubb(1) = 'glo' 
297         ENDIF
298
299         !                                      ! allocate dia_ptr arrays
300         IF( dia_ptr_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dia_ptr_init : unable to allocate arrays' )
301
302         rc_pwatt = rc_pwatt * rau0_rcp          ! conversion from K.s-1 to PetaWatt
303
304         IF( lk_mpp )   CALL mpp_ini_znl( numout )     ! Define MPI communicator for zonal sum
305
306         IF( ln_subbas ) THEN                ! load sub-basin mask
307            CALL iom_open( 'subbasins', inum,  ldstop = .FALSE.  )
308            CALL iom_get( inum, jpdom_data, 'atlmsk', btmsk(:,:,2) )   ! Atlantic basin
309            CALL iom_get( inum, jpdom_data, 'pacmsk', btmsk(:,:,3) )   ! Pacific  basin
310            CALL iom_get( inum, jpdom_data, 'indmsk', btmsk(:,:,4) )   ! Indian   basin
311            CALL iom_close( inum )
312            btmsk(:,:,5) = MAX ( btmsk(:,:,3), btmsk(:,:,4) )          ! Indo-Pacific basin
313            WHERE( gphit(:,:) < -30._wp)   ;   btm30(:,:) = 0._wp      ! mask out Southern Ocean
314            ELSE WHERE                     ;   btm30(:,:) = ssmask(:,:)
315            END WHERE
316         ENDIF
317   
318         btmsk(:,:,1) = tmask_i(:,:)                                   ! global ocean
319     
320         DO jn = 1, nptr
321            btmsk(:,:,jn) = btmsk(:,:,jn) * tmask_i(:,:)               ! interior domain only
322         END DO
323
324         ! Initialise arrays to zero because diatpr is called before they are first calculated
325         ! Note that this means diagnostics will not be exactly correct when model run is restarted.
326         htr_adv(:,:) = 0._wp  ;  str_adv(:,:) =  0._wp 
327         htr_ldf(:,:) = 0._wp  ;  str_ldf(:,:) =  0._wp 
328         htr_eiv(:,:) = 0._wp  ;  str_eiv(:,:) =  0._wp 
329         !
330      ENDIF 
331      !
332   END SUBROUTINE dia_ptr_init
333
334   SUBROUTINE dia_ptr_ohst_components( ktra, cptr, pva ) 
335      !!----------------------------------------------------------------------
336      !!                    ***  ROUTINE dia_ptr_ohst_components  ***
337      !!----------------------------------------------------------------------
338      !! Wrapper for heat and salt transport calculations to calculate them for each basin
339      !! Called from all advection and/or diffusion routines
340      !!----------------------------------------------------------------------
341      INTEGER                         , INTENT(in )  :: ktra  ! tracer index
342      CHARACTER(len=3)                , INTENT(in)   :: cptr  ! transport type  'adv'/'ldf'/'eiv'
343      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)   :: pva   ! 3D input array of advection/diffusion
344      INTEGER                                        :: jn    !
345
346     
347      IF( cptr == 'adv' ) THEN
348         IF( ktra == jp_tem )  htr_adv(:,1) = ptr_sj( pva(:,:,:) )
349         IF( ktra == jp_sal )  str_adv(:,1) = ptr_sj( pva(:,:,:) )
350      ENDIF
351      IF( cptr == 'ldf' ) THEN
352         IF( ktra == jp_tem )  htr_ldf(:,1) = ptr_sj( pva(:,:,:) )
353         IF( ktra == jp_sal )  str_ldf(:,1) = ptr_sj( pva(:,:,:) )
354      ENDIF
355      IF( cptr == 'eiv' ) THEN
356         IF( ktra == jp_tem )  htr_eiv(:,1) = ptr_sj( pva(:,:,:) )
357         IF( ktra == jp_sal )  str_eiv(:,1) = ptr_sj( pva(:,:,:) )
358      ENDIF
359      !
360      IF( ln_subbas ) THEN
361         !
362         IF( cptr == 'adv' ) THEN
363             IF( ktra == jp_tem ) THEN
364                DO jn = 2, nptr
365                   htr_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) )
366                END DO
367             ENDIF
368             IF( ktra == jp_sal ) THEN
369                DO jn = 2, nptr
370                   str_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) )
371                END DO
372             ENDIF
373         ENDIF
374         IF( cptr == 'ldf' ) THEN
375             IF( ktra == jp_tem ) THEN
376                DO jn = 2, nptr
377                    htr_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) )
378                 END DO
379             ENDIF
380             IF( ktra == jp_sal ) THEN
381                DO jn = 2, nptr
382                   str_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) )
383                END DO
384             ENDIF
385         ENDIF
386         IF( cptr == 'eiv' ) THEN
387             IF( ktra == jp_tem ) THEN
388                DO jn = 2, nptr
389                    htr_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) )
390                 END DO
391             ENDIF
392             IF( ktra == jp_sal ) THEN
393                DO jn = 2, nptr
394                   str_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) )
395                END DO
396             ENDIF
397         ENDIF
398         !
399      ENDIF
400
401   END SUBROUTINE
402
403
404   FUNCTION dia_ptr_alloc()
405      !!----------------------------------------------------------------------
406      !!                    ***  ROUTINE dia_ptr_alloc  ***
407      !!----------------------------------------------------------------------
408      INTEGER               ::   dia_ptr_alloc   ! return value
409      INTEGER, DIMENSION(3) ::   ierr
410      !!----------------------------------------------------------------------
411      ierr(:) = 0
412      !
413      ALLOCATE( btmsk(jpi,jpj,nptr) ,           &
414         &      htr_adv(jpj,nptr) , str_adv(jpj,nptr) ,   &
415         &      htr_eiv(jpj,nptr) , str_eiv(jpj,nptr) ,   &
416         &      htr_ldf(jpj,nptr) , str_ldf(jpj,nptr) , STAT=ierr(1)  )
417         !
418      ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2))
419      !
420      ALLOCATE( btm30(jpi,jpj), STAT=ierr(3)  )
421
422         !
423      dia_ptr_alloc = MAXVAL( ierr )
424      IF(lk_mpp)   CALL mpp_sum( dia_ptr_alloc )
425      !
426   END FUNCTION dia_ptr_alloc
427
428
429   FUNCTION ptr_sj_3d( pva, pmsk )   RESULT ( p_fval )
430      !!----------------------------------------------------------------------
431      !!                    ***  ROUTINE ptr_sj_3d  ***
432      !!
433      !! ** Purpose :   i-k sum computation of a j-flux array
434      !!
435      !! ** Method  : - i-k sum of pva using the interior 2D vmask (vmask_i).
436      !!              pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
437      !!
438      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
439      !!----------------------------------------------------------------------
440      REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk)       ::   pva   ! mask flux array at V-point
441      REAL(wp), INTENT(in), DIMENSION(jpi,jpj), OPTIONAL ::   pmsk   ! Optional 2D basin mask
442      !
443      INTEGER                  ::   ji, jj, jk   ! dummy loop arguments
444      INTEGER                  ::   ijpj         ! ???
445      REAL(wp), POINTER, DIMENSION(:) :: p_fval  ! function value
446      !!--------------------------------------------------------------------
447      !
448      p_fval => p_fval1d
449
450      ijpj = jpj
451      p_fval(:) = 0._wp
452      IF( PRESENT( pmsk ) ) THEN
453         DO jk = 1, jpkm1
454            DO jj = 2, jpjm1
455               DO ji = fs_2, fs_jpim1   ! Vector opt.
456                  p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj) * pmsk(ji,jj)
457               END DO
458            END DO
459         END DO
460      ELSE
461         DO jk = 1, jpkm1
462            DO jj = 2, jpjm1
463               DO ji = fs_2, fs_jpim1   ! Vector opt.
464                  p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj) 
465               END DO
466            END DO
467         END DO
468      ENDIF
469#if defined key_mpp_mpi
470      IF(lk_mpp)   CALL mpp_sum( p_fval, ijpj, ncomm_znl)
471#endif
472      !
473   END FUNCTION ptr_sj_3d
474
475
476   FUNCTION ptr_sj_2d( pva, pmsk )   RESULT ( p_fval )
477      !!----------------------------------------------------------------------
478      !!                    ***  ROUTINE ptr_sj_2d  ***
479      !!
480      !! ** Purpose :   "zonal" and vertical sum computation of a i-flux array
481      !!
482      !! ** Method  : - i-k sum of pva using the interior 2D vmask (vmask_i).
483      !!      pva is supposed to be a masked flux (i.e. * vmask*e1v*e3v)
484      !!
485      !! ** Action  : - p_fval: i-k-mean poleward flux of pva
486      !!----------------------------------------------------------------------
487      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)           ::   pva   ! mask flux array at V-point
488      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj), OPTIONAL ::   pmsk   ! Optional 2D basin mask
489      !
490      INTEGER                  ::   ji,jj       ! dummy loop arguments
491      INTEGER                  ::   ijpj        ! ???
492      REAL(wp), POINTER, DIMENSION(:) :: p_fval ! function value
493      !!--------------------------------------------------------------------
494      !
495      p_fval => p_fval1d
496
497      ijpj = jpj
498      p_fval(:) = 0._wp
499      IF( PRESENT( pmsk ) ) THEN
500         DO jj = 2, jpjm1
501            DO ji = nldi, nlei   ! No vector optimisation here. Better use a mask ?
502               p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj) * pmsk(ji,jj)
503            END DO
504         END DO
505      ELSE
506         DO jj = 2, jpjm1
507            DO ji = nldi, nlei   ! No vector optimisation here. Better use a mask ?
508               p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj)
509            END DO
510         END DO
511      ENDIF
512#if defined key_mpp_mpi
513      CALL mpp_sum( p_fval, ijpj, ncomm_znl )
514#endif
515      !
516   END FUNCTION ptr_sj_2d
517
518
519   FUNCTION ptr_sjk( pta, pmsk )   RESULT ( p_fval )
520      !!----------------------------------------------------------------------
521      !!                    ***  ROUTINE ptr_sjk  ***
522      !!
523      !! ** Purpose :   i-sum computation of an array
524      !!
525      !! ** Method  : - i-sum of pva using the interior 2D vmask (vmask_i).
526      !!
527      !! ** Action  : - p_fval: i-mean poleward flux of pva
528      !!----------------------------------------------------------------------
529      !!
530      IMPLICIT none
531      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk)           ::   pta    ! mask flux array at V-point
532      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)    , OPTIONAL ::   pmsk   ! Optional 2D basin mask
533      !!
534      INTEGER                           :: ji, jj, jk ! dummy loop arguments
535      REAL(wp), POINTER, DIMENSION(:,:) :: p_fval     ! return function value
536#if defined key_mpp_mpi
537      INTEGER, DIMENSION(1) ::   ish
538      INTEGER, DIMENSION(2) ::   ish2
539      INTEGER               ::   ijpjjpk
540      REAL(wp), DIMENSION(jpj*jpk) ::   zwork    ! mask flux array at V-point
541#endif
542      !!--------------------------------------------------------------------
543      !
544      p_fval => p_fval2d
545
546      p_fval(:,:) = 0._wp
547      !
548      IF( PRESENT( pmsk ) ) THEN
549         DO jk = 1, jpkm1
550            DO jj = 2, jpjm1
551!!gm here, use of tmask_i  ==> no need of loop over nldi, nlei....
552               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
553                  p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * pmsk(ji,jj)
554               END DO
555            END DO
556         END DO
557      ELSE
558         DO jk = 1, jpkm1
559            DO jj = 2, jpjm1
560               DO ji =  nldi, nlei   ! No vector optimisation here. Better use a mask ?
561                  p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * tmask_i(ji,jj)
562               END DO
563            END DO
564         END DO
565      END IF
566      !
567#if defined key_mpp_mpi
568      ijpjjpk = jpj*jpk
569      ish(1) = ijpjjpk  ;   ish2(1) = jpj   ;   ish2(2) = jpk
570      zwork(1:ijpjjpk) = RESHAPE( p_fval, ish )
571      CALL mpp_sum( zwork, ijpjjpk, ncomm_znl )
572      p_fval(:,:) = RESHAPE( zwork, ish2 )
573#endif
574      !
575   END FUNCTION ptr_sjk
576
577
578   !!======================================================================
579END MODULE diaptr
Note: See TracBrowser for help on using the repository browser.