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.
prtctl.F90 in NEMO/trunk/NEMOGCM/NEMO/OCE_SRC/IOM – NEMO

source: NEMO/trunk/NEMOGCM/NEMO/OCE_SRC/IOM/prtctl.F90 @ 9594

Last change on this file since 9594 was 9570, checked in by nicolasmartin, 6 years ago

Global renaming for core routines (./NEMO)

  • Folders
    • LIM_SRC_3 -> ICE_SRC
    • OPA_SRC -> OCE_SRC
  • CPP key: key_lim3 -> key_si3
  • Modules, (sub)routines and variables names
    • MPI: mpi_comm_opa -> mpi_comm_oce, MPI_COMM_OPA -> MPI_COMM_OCE, mpi_init_opa -> mpi_init_oce
    • AGRIF: agrif_opa_* -> agrif_oce_*, agrif_lim3_* -> agrif_si3_* and few more
    • TOP-PISCES: p.zlim -> p.zice, namp.zlim -> namp.zice
  • Comments
    • NEMO/OPA -> NEMO/OCE
    • ESIM|LIM3 -> SI3
  • Property svn:keywords set to Id
File size: 24.8 KB
Line 
1MODULE prtctl
2   !!======================================================================
3   !!                       ***  MODULE prtctl   ***
4   !! Ocean system : print all SUM trends for each processor domain
5   !!======================================================================
6   !! History :  9.0  !  05-07  (C. Talandier) original code
7   !!            3.4  !  11-11  (C. Harris) decomposition changes for running with CICE
8   !!----------------------------------------------------------------------
9   USE dom_oce          ! ocean space and time domain variables
10#if defined key_nemocice_decomp
11   USE ice_domain_size, only: nx_global, ny_global
12#endif
13   USE in_out_manager   ! I/O manager
14   USE lib_mpp          ! distributed memory computing
15
16   IMPLICIT NONE
17   PRIVATE
18
19   INTEGER , DIMENSION(:), ALLOCATABLE, SAVE ::   numid
20   INTEGER , DIMENSION(:), ALLOCATABLE, SAVE ::   nlditl , nldjtl    ! first, last indoor index for each i-domain
21   INTEGER , DIMENSION(:), ALLOCATABLE, SAVE ::   nleitl , nlejtl    ! first, last indoor index for each j-domain
22   INTEGER , DIMENSION(:), ALLOCATABLE, SAVE ::   nimpptl, njmpptl   ! i-, j-indexes for each processor
23   INTEGER , DIMENSION(:), ALLOCATABLE, SAVE ::   nlcitl , nlcjtl    ! dimensions of every subdomain
24   INTEGER , DIMENSION(:), ALLOCATABLE, SAVE ::   ibonitl, ibonjtl   !
25
26   REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE ::   t_ctll , s_ctll    ! previous tracer trend values
27   REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE ::   u_ctll , v_ctll    ! previous velocity trend values
28
29   INTEGER ::   ktime   ! time step
30
31   PUBLIC prt_ctl         ! called by all subroutines
32   PUBLIC prt_ctl_info    ! called by all subroutines
33   PUBLIC prt_ctl_init    ! called by opa.F90
34   PUBLIC sub_dom         ! called by opa.F90
35
36   !!----------------------------------------------------------------------
37   !! NEMO/OCE 3.3 , NEMO Consortium (2010)
38   !! $Id$
39   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
40   !!----------------------------------------------------------------------
41CONTAINS
42
43   SUBROUTINE prt_ctl (tab2d_1, tab3d_1, mask1, clinfo1, tab2d_2, tab3d_2,   &
44      &                                  mask2, clinfo2, kdim, clinfo3 )
45      !!----------------------------------------------------------------------
46      !!                     ***  ROUTINE prt_ctl  ***
47      !!
48      !! ** Purpose : - print sum control of 2D or 3D arrays over the same area
49      !!                in mono and mpp case. This way can be usefull when
50      !!                debugging a new parametrization in mono or mpp.
51      !!
52      !! ** Method  : 2 possibilities exist when setting the ln_ctl parameter to
53      !!                .true. in the ocean namelist:
54      !!              - to debug a MPI run .vs. a mono-processor one;
55      !!                the control print will be done over each sub-domain.
56      !!                The nictl[se] and njctl[se] parameters in the namelist must
57      !!                be set to zero and [ij]splt to the corresponding splitted
58      !!                domain in MPI along respectively i-, j- directions.
59      !!              - to debug a mono-processor run over the whole domain/a specific area;
60      !!                in the first case the nictl[se] and njctl[se] parameters must be set
61      !!                to zero else to the indices of the area to be controled. In both cases
62      !!                isplt and jsplt must be set to 1.
63      !!              - All arguments of the above calling sequence are optional so their
64      !!                name must be explicitly typed if used. For instance if the 3D
65      !!                array tn(:,:,:) must be passed through the prt_ctl subroutine,
66      !!                it must looks like: CALL prt_ctl(tab3d_1=tn).
67      !!
68      !!                    tab2d_1 : first 2D array
69      !!                    tab3d_1 : first 3D array
70      !!                    mask1   : mask (3D) to apply to the tab[23]d_1 array
71      !!                    clinfo1 : information about the tab[23]d_1 array
72      !!                    tab2d_2 : second 2D array
73      !!                    tab3d_2 : second 3D array
74      !!                    mask2   : mask (3D) to apply to the tab[23]d_2 array
75      !!                    clinfo2 : information about the tab[23]d_2 array
76      !!                    kdim    : k- direction for 3D arrays
77      !!                    clinfo3 : additional information
78      !!----------------------------------------------------------------------
79      REAL(wp), DIMENSION(:,:)  , INTENT(in), OPTIONAL ::   tab2d_1
80      REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL ::   tab3d_1
81      REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL ::   mask1
82      CHARACTER (len=*)         , INTENT(in), OPTIONAL ::   clinfo1
83      REAL(wp), DIMENSION(:,:)  , INTENT(in), OPTIONAL ::   tab2d_2
84      REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL ::   tab3d_2
85      REAL(wp), DIMENSION(:,:,:), INTENT(in), OPTIONAL ::   mask2
86      CHARACTER (len=*)         , INTENT(in), OPTIONAL ::   clinfo2
87      INTEGER                   , INTENT(in), OPTIONAL ::   kdim
88      CHARACTER (len=*)         , INTENT(in), OPTIONAL ::   clinfo3
89      !
90      CHARACTER (len=15) :: cl2
91      INTEGER ::  jn, sind, eind, kdir,j_id
92      REAL(wp) :: zsum1, zsum2, zvctl1, zvctl2
93      REAL(wp), DIMENSION(jpi,jpj)     :: ztab2d_1, ztab2d_2
94      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmask1, zmask2, ztab3d_1, ztab3d_2
95      !!----------------------------------------------------------------------
96
97      ! Arrays, scalars initialization
98      kdir      = jpkm1
99      cl2       = ''
100      zsum1     = 0.e0
101      zsum2     = 0.e0
102      zvctl1    = 0.e0
103      zvctl2    = 0.e0
104      ztab2d_1(:,:)   = 0.e0
105      ztab2d_2(:,:)   = 0.e0
106      ztab3d_1(:,:,:) = 0.e0
107      ztab3d_2(:,:,:) = 0.e0
108      zmask1  (:,:,:) = 1.e0
109      zmask2  (:,:,:) = 1.e0
110
111      ! Control of optional arguments
112      IF( PRESENT(clinfo2) )   cl2                  = clinfo2
113      IF( PRESENT(kdim)    )   kdir                 = kdim
114      IF( PRESENT(tab2d_1) )   ztab2d_1(:,:)        = tab2d_1(:,:)
115      IF( PRESENT(tab2d_2) )   ztab2d_2(:,:)        = tab2d_2(:,:)
116      IF( PRESENT(tab3d_1) )   ztab3d_1(:,:,1:kdir) = tab3d_1(:,:,1:kdir)
117      IF( PRESENT(tab3d_2) )   ztab3d_2(:,:,1:kdir) = tab3d_2(:,:,1:kdir)
118      IF( PRESENT(mask1)   )   zmask1  (:,:,:)      = mask1  (:,:,:)
119      IF( PRESENT(mask2)   )   zmask2  (:,:,:)      = mask2  (:,:,:)
120
121      IF( lk_mpp .AND. jpnij > 1 ) THEN       ! processor number
122         sind = narea
123         eind = narea
124      ELSE                                    ! processors total number
125         sind = 1
126         eind = ijsplt
127      ENDIF
128
129      ! Loop over each sub-domain, i.e. the total number of processors ijsplt
130      DO jn = sind, eind
131         ! Set logical unit
132         j_id = numid(jn - narea + 1)
133         ! Set indices for the SUM control
134         IF( .NOT. lsp_area ) THEN
135            IF (lk_mpp .AND. jpnij > 1)   THEN
136               nictls = MAX(  1, nlditl(jn) )
137               nictle = MIN(jpi, nleitl(jn) )
138               njctls = MAX(  1, nldjtl(jn) )
139               njctle = MIN(jpj, nlejtl(jn) )
140               ! Do not take into account the bound of the domain
141               IF( ibonitl(jn) == -1 .OR. ibonitl(jn) == 2 ) nictls = MAX(2, nictls)
142               IF( ibonjtl(jn) == -1 .OR. ibonjtl(jn) == 2 ) njctls = MAX(2, njctls)
143               IF( ibonitl(jn) ==  1 .OR. ibonitl(jn) == 2 ) nictle = MIN(nictle, nleitl(jn) - 1)
144               IF( ibonjtl(jn) ==  1 .OR. ibonjtl(jn) == 2 ) njctle = MIN(njctle, nlejtl(jn) - 1)
145            ELSE
146               nictls = MAX(  1, nimpptl(jn) - 1 + nlditl(jn) )
147               nictle = MIN(jpi, nimpptl(jn) - 1 + nleitl(jn) )
148               njctls = MAX(  1, njmpptl(jn) - 1 + nldjtl(jn) )
149               njctle = MIN(jpj, njmpptl(jn) - 1 + nlejtl(jn) )
150               ! Do not take into account the bound of the domain
151               IF( ibonitl(jn) == -1 .OR. ibonitl(jn) == 2 ) nictls = MAX(2, nictls)
152               IF( ibonjtl(jn) == -1 .OR. ibonjtl(jn) == 2 ) njctls = MAX(2, njctls)
153               IF( ibonitl(jn) ==  1 .OR. ibonitl(jn) == 2 ) nictle = MIN(nictle, nimpptl(jn) + nleitl(jn) - 2)
154               IF( ibonjtl(jn) ==  1 .OR. ibonjtl(jn) == 2 ) njctle = MIN(njctle, njmpptl(jn) + nlejtl(jn) - 2)
155            ENDIF
156         ENDIF
157
158         IF( PRESENT(clinfo3)) THEN
159            IF ( clinfo3 == 'tra' )  THEN
160               zvctl1 = t_ctll(jn)
161               zvctl2 = s_ctll(jn)
162            ELSEIF ( clinfo3 == 'dyn' )   THEN
163               zvctl1 = u_ctll(jn)
164               zvctl2 = v_ctll(jn)
165            ENDIF
166         ENDIF
167
168         ! Compute the sum control
169         ! 2D arrays
170         IF( PRESENT(tab2d_1) )   THEN
171            zsum1 = SUM( ztab2d_1(nictls:nictle,njctls:njctle)*zmask1(nictls:nictle,njctls:njctle,1) )
172            zsum2 = SUM( ztab2d_2(nictls:nictle,njctls:njctle)*zmask2(nictls:nictle,njctls:njctle,1) )
173         ENDIF
174
175         ! 3D arrays
176         IF( PRESENT(tab3d_1) )   THEN
177            zsum1 = SUM( ztab3d_1(nictls:nictle,njctls:njctle,1:kdir)*zmask1(nictls:nictle,njctls:njctle,1:kdir) )
178            zsum2 = SUM( ztab3d_2(nictls:nictle,njctls:njctle,1:kdir)*zmask2(nictls:nictle,njctls:njctle,1:kdir) )
179         ENDIF
180
181         ! Print the result
182         IF( PRESENT(clinfo3) )   THEN
183            WRITE(j_id,FMT='(a,D23.16,3x,a,D23.16)')clinfo1, zsum1-zvctl1, cl2, zsum2-zvctl2
184            SELECT CASE( clinfo3 )
185            CASE ( 'tra-ta' ) 
186               t_ctll(jn) = zsum1
187            CASE ( 'tra' ) 
188                t_ctll(jn) = zsum1
189                s_ctll(jn) = zsum2
190            CASE ( 'dyn' ) 
191                u_ctll(jn) = zsum1
192                v_ctll(jn) = zsum2 
193            END SELECT
194         ELSEIF ( PRESENT(clinfo2) .OR. PRESENT(tab2d_2) .OR. PRESENT(tab3d_2) )   THEN
195            WRITE(j_id,FMT='(a,D23.16,3x,a,D23.16)')clinfo1, zsum1, cl2, zsum2
196         ELSE
197            WRITE(j_id,FMT='(a,D23.16)')clinfo1, zsum1
198         ENDIF
199
200      ENDDO
201      !
202   END SUBROUTINE prt_ctl
203
204
205   SUBROUTINE prt_ctl_info (clinfo1, ivar1, clinfo2, ivar2, itime)
206      !!----------------------------------------------------------------------
207      !!                     ***  ROUTINE prt_ctl_info  ***
208      !!
209      !! ** Purpose : - print information without any computation
210      !!
211      !! ** Action  : - input arguments
212      !!                    clinfo1 : information about the ivar1
213      !!                    ivar1   : value to print
214      !!                    clinfo2 : information about the ivar2
215      !!                    ivar2   : value to print
216      !!----------------------------------------------------------------------
217      CHARACTER (len=*), INTENT(in)           ::   clinfo1
218      INTEGER          , INTENT(in), OPTIONAL ::   ivar1
219      CHARACTER (len=*), INTENT(in), OPTIONAL ::   clinfo2
220      INTEGER          , INTENT(in), OPTIONAL ::   ivar2
221      INTEGER          , INTENT(in), OPTIONAL ::   itime
222      !
223      INTEGER :: jn, sind, eind, iltime, j_id
224      !!----------------------------------------------------------------------
225
226      IF( lk_mpp .AND. jpnij > 1 ) THEN       ! processor number
227         sind = narea
228         eind = narea
229      ELSE                                    ! total number of processors
230         sind = 1
231         eind = ijsplt
232      ENDIF
233
234      ! Set to zero arrays at each new time step
235      IF( PRESENT(itime) )   THEN
236         iltime = itime
237         IF( iltime > ktime )   THEN
238            t_ctll(:) = 0.e0   ;   s_ctll(:) = 0.e0
239            u_ctll(:) = 0.e0   ;   v_ctll(:) = 0.e0
240            ktime = iltime
241         ENDIF
242      ENDIF
243
244      ! Loop over each sub-domain, i.e. number of processors ijsplt
245      DO jn = sind, eind
246         !
247         j_id = numid(jn - narea + 1)         ! Set logical unit
248         !
249         IF( PRESENT(ivar1) .AND. PRESENT(clinfo2) .AND. PRESENT(ivar2) )   THEN
250            WRITE(j_id,*)clinfo1, ivar1, clinfo2, ivar2
251         ELSEIF ( PRESENT(ivar1) .AND. PRESENT(clinfo2) .AND. .NOT. PRESENT(ivar2) )   THEN
252            WRITE(j_id,*)clinfo1, ivar1, clinfo2
253         ELSEIF ( PRESENT(ivar1) .AND. .NOT. PRESENT(clinfo2) .AND. PRESENT(ivar2) )   THEN
254            WRITE(j_id,*)clinfo1, ivar1, ivar2
255         ELSEIF ( PRESENT(ivar1) .AND. .NOT. PRESENT(clinfo2) .AND. .NOT. PRESENT(ivar2) )   THEN
256            WRITE(j_id,*)clinfo1, ivar1
257         ELSE
258            WRITE(j_id,*)clinfo1
259         ENDIF
260         !
261      END DO
262      !
263   END SUBROUTINE prt_ctl_info
264
265
266   SUBROUTINE prt_ctl_init
267      !!----------------------------------------------------------------------
268      !!                     ***  ROUTINE prt_ctl_init  ***
269      !!
270      !! ** Purpose :   open ASCII files & compute indices
271      !!----------------------------------------------------------------------
272      INTEGER ::   jn, sind, eind, j_id
273      CHARACTER (len=28) :: clfile_out
274      CHARACTER (len=23) :: clb_name
275      CHARACTER (len=19) :: cl_run
276      !!----------------------------------------------------------------------
277
278      ! Allocate arrays
279      ALLOCATE( nlditl(ijsplt) , nleitl(ijsplt) , nimpptl(ijsplt) , ibonitl(ijsplt) ,   &
280         &      nldjtl(ijsplt) , nlejtl(ijsplt) , njmpptl(ijsplt) , ibonjtl(ijsplt) ,   &
281         &      nlcitl(ijsplt) , t_ctll(ijsplt) , u_ctll (ijsplt) ,                     &
282         &      nlcjtl(ijsplt) , s_ctll(ijsplt) , v_ctll (ijsplt)                       )
283
284      ! Initialization
285      t_ctll(:) = 0.e0
286      s_ctll(:) = 0.e0
287      u_ctll(:) = 0.e0
288      v_ctll(:) = 0.e0
289      ktime = 1
290
291      IF( lk_mpp .AND. jpnij > 1 ) THEN
292         sind = narea
293         eind = narea
294         clb_name = "('mpp.output_',I4.4)"
295         cl_run = 'MULTI processor run'
296         ! use indices for each area computed by mpp_init subroutine
297         nlditl(1:jpnij) = nldit(:) 
298         nleitl(1:jpnij) = nleit(:) 
299         nldjtl(1:jpnij) = nldjt(:) 
300         nlejtl(1:jpnij) = nlejt(:) 
301         !
302         nimpptl(1:jpnij) = nimppt(:)
303         njmpptl(1:jpnij) = njmppt(:)
304         !
305         nlcitl(1:jpnij) = nlcit(:)
306         nlcjtl(1:jpnij) = nlcjt(:)
307         !
308         ibonitl(1:jpnij) = ibonit(:)
309         ibonjtl(1:jpnij) = ibonjt(:)
310      ELSE
311         sind = 1
312         eind = ijsplt
313         clb_name = "('mono.output_',I4.4)"
314         cl_run = 'MONO processor run '
315         ! compute indices for each area as done in mpp_init subroutine
316         CALL sub_dom
317      ENDIF
318
319      ALLOCATE( numid(eind-sind+1) )
320
321      DO jn = sind, eind
322         WRITE(clfile_out,FMT=clb_name) jn-1
323         CALL ctl_opn( numid(jn -narea + 1), clfile_out, 'REPLACE', 'FORMATTED', 'SEQUENTIAL', 1, numout, .FALSE. )
324         j_id = numid(jn -narea + 1)
325         WRITE(j_id,*)
326         WRITE(j_id,*) '                 L O D Y C - I P S L'
327         WRITE(j_id,*) '                     O P A model'
328         WRITE(j_id,*) '            Ocean General Circulation Model'
329         WRITE(j_id,*) '               version OPA 9.0  (2005) '
330         WRITE(j_id,*)
331         WRITE(j_id,*) '                   PROC number: ', jn
332         WRITE(j_id,*)
333         WRITE(j_id,FMT="(19x,a20)")cl_run
334
335         ! Print the SUM control indices
336         IF( .NOT. lsp_area )   THEN
337            nictls = nimpptl(jn) + nlditl(jn) - 1
338            nictle = nimpptl(jn) + nleitl(jn) - 1
339            njctls = njmpptl(jn) + nldjtl(jn) - 1
340            njctle = njmpptl(jn) + nlejtl(jn) - 1
341         ENDIF
342         WRITE(j_id,*) 
343         WRITE(j_id,*) 'prt_ctl :  Sum control indices'
344         WRITE(j_id,*) '~~~~~~~'
345         WRITE(j_id,*)
346         WRITE(j_id,9000)'                                nlej   = ', nlejtl(jn), '              '
347         WRITE(j_id,9000)'                  ------------- njctle = ', njctle, ' -------------'
348         WRITE(j_id,9001)'                  |                                       |'
349         WRITE(j_id,9001)'                  |                                       |'
350         WRITE(j_id,9001)'                  |                                       |'
351         WRITE(j_id,9002)'           nictls = ', nictls,  '                           nictle = ', nictle
352         WRITE(j_id,9002)'           nldi   = ', nlditl(jn),  '                           nlei   = ', nleitl(jn)
353         WRITE(j_id,9001)'                  |                                       |'
354         WRITE(j_id,9001)'                  |                                       |'
355         WRITE(j_id,9001)'                  |                                       |'
356         WRITE(j_id,9004)'  njmpp  = ',njmpptl(jn),'   ------------- njctls = ', njctls, ' -------------'
357         WRITE(j_id,9003)'           nimpp  = ', nimpptl(jn), '        nldj   = ', nldjtl(jn), '              '
358         WRITE(j_id,*)
359         WRITE(j_id,*)
360
3619000     FORMAT(a41,i4.4,a14)
3629001     FORMAT(a59)
3639002     FORMAT(a20,i4.4,a36,i3.3)
3649003     FORMAT(a20,i4.4,a17,i4.4)
3659004     FORMAT(a11,i4.4,a26,i4.4,a14)
366      END DO
367      !
368   END SUBROUTINE prt_ctl_init
369
370
371   SUBROUTINE sub_dom
372      !!----------------------------------------------------------------------
373      !!                  ***  ROUTINE sub_dom  ***
374      !!                   
375      !! ** Purpose :   Lay out the global domain over processors.
376      !!                CAUTION:
377      !!                This part has been extracted from the mpp_init
378      !!                subroutine and names of variables/arrays have been
379      !!                slightly changed to avoid confusion but the computation
380      !!                is exactly the same. Any modification about indices of
381      !!                each sub-domain in the mppini.F90 module should be reported
382      !!                here.
383      !!
384      !! ** Method  :   Global domain is distributed in smaller local domains.
385      !!                Periodic condition is a function of the local domain position
386      !!                (global boundary or neighbouring domain) and of the global
387      !!                periodic
388      !!                Type :         jperio global periodic condition
389      !!                               nperio local  periodic condition
390      !!
391      !! ** Action  : - set domain parameters
392      !!                    nimpp     : longitudinal index
393      !!                    njmpp     : latitudinal  index
394      !!                    nperio    : lateral condition type
395      !!                    narea     : number for local area
396      !!                    nlcil      : first dimension
397      !!                    nlcjl      : second dimension
398      !!                    nbondil    : mark for "east-west local boundary"
399      !!                    nbondjl    : mark for "north-south local boundary"
400      !!
401      !! History :
402      !!        !  94-11  (M. Guyon)  Original code
403      !!        !  95-04  (J. Escobar, M. Imbard)
404      !!        !  98-02  (M. Guyon)  FETI method
405      !!        !  98-05  (M. Imbard, J. Escobar, L. Colombet )  SHMEM and MPI versions
406      !!   8.5  !  02-08  (G. Madec)  F90 : free form
407      !!----------------------------------------------------------------------
408      INTEGER ::   ji, jj, jn               ! dummy loop indices
409      INTEGER ::   &
410         ii, ij,                         &  ! temporary integers
411         irestil, irestjl,               &  !    "          "
412         ijpi  , ijpj, nlcil,            &  ! temporary logical unit
413         nlcjl , nbondil, nbondjl,       &
414         nrecil, nrecjl, nldil, nleil, nldjl, nlejl
415
416      INTEGER, DIMENSION(jpi,jpj) ::   iimpptl, ijmpptl, ilcitl, ilcjtl   ! workspace
417      REAL(wp) ::   zidom, zjdom            ! temporary scalars
418      INTEGER ::   inum                     ! local logical unit
419      !!----------------------------------------------------------------------
420
421      !
422      !
423      !  1. Dimension arrays for subdomains
424      ! -----------------------------------
425      !  Computation of local domain sizes ilcitl() ilcjtl()
426      !  These dimensions depend on global sizes isplt,jsplt and jpiglo,jpjglo
427      !  The subdomains are squares leeser than or equal to the global
428      !  dimensions divided by the number of processors minus the overlap
429      !  array (cf. par_oce.F90).
430
431#if defined key_nemocice_decomp
432      ijpi = ( nx_global+2-2*nn_hls + (isplt-1) ) / isplt + 2*nn_hls
433      ijpj = ( ny_global+2-2*nn_hls + (jsplt-1) ) / jsplt + 2*nn_hls 
434#else
435      ijpi = ( jpiglo-2*nn_hls + (isplt-1) ) / isplt + 2*nn_hls
436      ijpj = ( jpjglo-2*nn_hls + (jsplt-1) ) / jsplt + 2*nn_hls
437#endif
438
439
440      nrecil  = 2 * nn_hls
441      nrecjl  = 2 * nn_hls
442      irestil = MOD( jpiglo - nrecil , isplt )
443      irestjl = MOD( jpjglo - nrecjl , jsplt )
444
445      IF(  irestil == 0 )   irestil = isplt
446#if defined key_nemocice_decomp
447
448      ! In order to match CICE the size of domains in NEMO has to be changed
449      ! The last line of blocks (west) will have fewer points
450      DO jj = 1, jsplt 
451         DO ji=1, isplt-1 
452            ilcitl(ji,jj) = ijpi 
453         END DO
454         ilcitl(isplt,jj) = jpiglo - (isplt - 1) * (ijpi - nrecil)
455      END DO 
456
457#else
458
459      DO jj = 1, jsplt
460         DO ji = 1, irestil
461            ilcitl(ji,jj) = ijpi
462         END DO
463         DO ji = irestil+1, isplt
464            ilcitl(ji,jj) = ijpi -1
465         END DO
466      END DO
467
468#endif
469     
470      IF( irestjl == 0 )   irestjl = jsplt
471#if defined key_nemocice_decomp 
472
473      ! Same change to domains in North-South direction as in East-West.
474      DO ji = 1, isplt 
475         DO jj=1, jsplt-1 
476            ilcjtl(ji,jj) = ijpj 
477         END DO
478         ilcjtl(ji,jsplt) = jpjglo - (jsplt - 1) * (ijpj - nrecjl)
479      END DO 
480
481#else
482
483      DO ji = 1, isplt
484         DO jj = 1, irestjl
485            ilcjtl(ji,jj) = ijpj
486         END DO
487         DO jj = irestjl+1, jsplt
488            ilcjtl(ji,jj) = ijpj -1
489         END DO
490      END DO
491
492#endif
493      zidom = nrecil
494      DO ji = 1, isplt
495         zidom = zidom + ilcitl(ji,1) - nrecil
496      END DO
497      IF(lwp) WRITE(numout,*)
498      IF(lwp) WRITE(numout,*)' sum ilcitl(i,1) = ', zidom, ' jpiglo = ', jpiglo
499     
500      zjdom = nrecjl
501      DO jj = 1, jsplt
502         zjdom = zjdom + ilcjtl(1,jj) - nrecjl
503      END DO
504      IF(lwp) WRITE(numout,*)' sum ilcitl(1,j) = ', zjdom, ' jpjglo = ', jpjglo
505      IF(lwp) WRITE(numout,*)
506     
507
508      !  2. Index arrays for subdomains
509      ! -------------------------------
510
511      iimpptl(:,:) = 1
512      ijmpptl(:,:) = 1
513     
514      IF( isplt > 1 ) THEN
515         DO jj = 1, jsplt
516            DO ji = 2, isplt
517               iimpptl(ji,jj) = iimpptl(ji-1,jj) + ilcitl(ji-1,jj) - nrecil
518            END DO
519         END DO
520      ENDIF
521
522      IF( jsplt > 1 ) THEN
523         DO jj = 2, jsplt
524            DO ji = 1, isplt
525               ijmpptl(ji,jj) = ijmpptl(ji,jj-1)+ilcjtl(ji,jj-1)-nrecjl
526            END DO
527         END DO
528      ENDIF
529     
530      ! 3. Subdomain description
531      ! ------------------------
532
533      DO jn = 1, ijsplt
534         ii = 1 + MOD( jn-1, isplt )
535         ij = 1 + (jn-1) / isplt
536         nimpptl(jn) = iimpptl(ii,ij)
537         njmpptl(jn) = ijmpptl(ii,ij)
538         nlcitl (jn) = ilcitl (ii,ij)     
539         nlcil       = nlcitl (jn)     
540         nlcjtl (jn) = ilcjtl (ii,ij)     
541         nlcjl       = nlcjtl (jn)
542         nbondjl = -1                                    ! general case
543         IF( jn   >  isplt          )   nbondjl = 0      ! first row of processor
544         IF( jn   >  (jsplt-1)*isplt )  nbondjl = 1     ! last  row of processor
545         IF( jsplt == 1             )   nbondjl = 2      ! one processor only in j-direction
546         ibonjtl(jn) = nbondjl
547         
548         nbondil = 0                                     !
549         IF( MOD( jn, isplt ) == 1 )   nbondil = -1      !
550         IF( MOD( jn, isplt ) == 0 )   nbondil =  1      !
551         IF( isplt            == 1 )   nbondil =  2      ! one processor only in i-direction
552         ibonitl(jn) = nbondil
553         
554         nldil =  1   + nn_hls
555         nleil = nlcil - nn_hls
556         IF( nbondil == -1 .OR. nbondil == 2 )   nldil = 1
557         IF( nbondil ==  1 .OR. nbondil == 2 )   nleil = nlcil
558         nldjl =  1   + nn_hls
559         nlejl = nlcjl - nn_hls
560         IF( nbondjl == -1 .OR. nbondjl == 2 )   nldjl = 1
561         IF( nbondjl ==  1 .OR. nbondjl == 2 )   nlejl = nlcjl
562         nlditl(jn) = nldil
563         nleitl(jn) = nleil
564         nldjtl(jn) = nldjl
565         nlejtl(jn) = nlejl
566      END DO
567      !
568      ! Save processor layout in layout_prtctl.dat file
569      IF(lwp) THEN
570         CALL ctl_opn( inum, 'layout_prtctl.dat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE., narea )
571         WRITE(inum,'(a)') 'nproc nlcil nlcjl nldil nldjl nleil nlejl nimpptl njmpptl ibonitl ibonjtl'
572         !
573         DO jn = 1, ijsplt
574            WRITE(inum,'(i5,6i6,4i8)') jn-1,nlcitl(jn),  nlcjtl(jn), &
575               &                            nlditl(jn),  nldjtl(jn), &
576               &                            nleitl(jn),  nlejtl(jn), &
577               &                           nimpptl(jn), njmpptl(jn), &
578               &                           ibonitl(jn), ibonjtl(jn)
579         END DO
580         CLOSE(inum)   
581      END IF
582      !
583      !
584   END SUBROUTINE sub_dom
585
586   !!======================================================================
587END MODULE prtctl
Note: See TracBrowser for help on using the repository browser.