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.
modcluster.F in branches/UKMO/smagorinsky_diagnostics/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES – NEMO

source: branches/UKMO/smagorinsky_diagnostics/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES/modcluster.F @ 5699

Last change on this file since 5699 was 5699, checked in by davestorkey, 9 years ago

Remove SVN keyword updating from UKMO/smagorinsky_diagnostics branch.

File size: 47.9 KB
Line 
1!
2! $Id: modcluster.F 2715 2011-03-30 15:58:35Z rblod $
3!
4C     AGRIF (Adaptive Grid Refinement In Fortran)
5C
6C     Copyright (C) 2003 Laurent Debreu (Laurent.Debreu@imag.fr)
7C                        Christophe Vouland (Christophe.Vouland@imag.fr)   
8C
9C     This program is free software; you can redistribute it and/or modify
10C     it under the terms of the GNU General Public License as published by
11C     the Free Software Foundation; either version 2 of the License, or
12C     (at your option) any later version.
13C
14C     This program is distributed in the hope that it will be useful,
15C     but WITHOUT ANY WARRANTY; without even the implied warranty of
16C     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17C     GNU General Public License for more details.
18C
19C     You should have received a copy of the GNU General Public License
20C     along with this program; if not, write to the Free Software
21C     Foundation, Inc., 59 Temple Place- Suite 330, Boston, MA 02111-1307, USA.
22C
23C
24C
25CCC   Module AGRIF_Clustering 
26C
27C
28      Module Agrif_Clustering 
29C
30CCC   Description:
31CCC   Module to create and initialize the grid hierarchy from the 
32CCC   AGRIF_FixedGrids.in file.
33C
34C     Modules used:
35C 
36      Use Agrif_Curgridfunctions 
37      Use Agrif_Init_Vars
38      Use Agrif_Save
39C             
40      IMPLICIT NONE
41C       
42      Contains
43C     Define procedures contained in this module
44C
45C
46C
47C     **************************************************************************
48CCC   Subroutine Agrif_Cluster_All 
49C     **************************************************************************
50C
51      Recursive Subroutine Agrif_Cluster_All(g,coarsegrid)
52C
53CCC   Description:
54CCC   Subroutine for the clustering. A temporary grid hierarchy, pointed by 
55CCC   coarsegrid, is created.
56C
57CC    Method:
58CC           
59C
60C     Declarations:
61C
62C     Pointer arguments   
63      TYPE(AGRIF_grid)     ,pointer   :: g        ! Pointer on the current grid
64      TYPE(AGRIF_rectangle),pointer   :: coarsegrid
65C     
66C     Local pointer         
67      TYPE(AGRIF_lrectangle),pointer  :: parcours 
68      TYPE(AGRIF_grid)      ,pointer  :: newgrid 
69      REAL                            :: g_eps
70      INTEGER                         :: iii
71C     
72      g_eps = huge(1.)
73      do iii = 1 , Agrif_Probdim
74         g_eps = min(g_eps,g%Agrif_d(iii))
75      enddo
76C
77      g_eps = g_eps/100.                 
78C
79C     Necessary condition for clustering 
80      do iii = 1 , Agrif_Probdim
81         if (g%Agrif_d(iii)/Agrif_coeffref(iii).LT.
82     &                   (Agrif_mind(iii)-g_eps)) Return
83      enddo
84C
85      Nullify(coarsegrid%childgrids)
86C
87      Call Agrif_ClusterGridnD(g,coarsegrid) 
88C
89      parcours => coarsegrid % childgrids                   
90C     
91      do while (associated(parcours))
92C     
93C       Newgrid is created. It is a copy of a fine grid created previously by
94C          clustering.     
95        Allocate(newgrid)
96C     
97        Nullify(newgrid%child_grids)
98C
99        do iii = 1 , Agrif_Probdim
100           newgrid % nb(iii) = (parcours % r % imax(iii) - 
101     &                          parcours % r % imin(iii)) *
102     &                          Agrif_Coeffref(iii)
103C       
104           newgrid % Agrif_x(iii) = g%Agrif_x(iii) + 
105     &                         (parcours %r % imin(iii) -1) 
106     &                         *g%Agrif_d(iii)
107C 
108           newgrid % Agrif_d(iii) = g%Agrif_d(iii) / Agrif_Coeffref(iii)
109C
110        enddo
111C
112        if ( Agrif_Probdim .EQ. 1 ) then
113           allocate(newgrid%tabpoint1D(newgrid%nb(1)+1))
114           newgrid%tabpoint1D = 0
115        endif
116C
117        if ( Agrif_Probdim .EQ. 2 ) then
118           allocate(newgrid%tabpoint2D(newgrid%nb(1)+1,
119     &                                 newgrid%nb(2)+1))
120           newgrid%tabpoint2D = 0     
121        endif
122C
123        if ( Agrif_Probdim .EQ. 3 ) then
124           allocate(newgrid%tabpoint3D(newgrid%nb(1)+1,
125     &                                 newgrid%nb(2)+1,
126     &                                 newgrid%nb(3)+1))
127           newgrid%tabpoint3D = 0
128        endif
129C       Points detection on newgrid
130        Call Agrif_TabpointsnD(Agrif_mygrid,newgrid) 
131C
132C       Recursive call to Agrif_Cluster_All           
133        Call Agrif_Cluster_All (newgrid, parcours % r)     
134C
135        parcours => parcours % next 
136C     
137        if ( Agrif_Probdim .EQ. 1 ) Deallocate(newgrid%tabpoint1D)
138        if ( Agrif_Probdim .EQ. 2 ) Deallocate(newgrid%tabpoint2D)
139        if ( Agrif_Probdim .EQ. 3 ) Deallocate(newgrid%tabpoint3D)
140C
141        Deallocate(newgrid)   
142C
143      enddo
144C
145C     
146      End Subroutine Agrif_Cluster_All
147C
148C     **************************************************************************
149CCC   Subroutine Agrif_TabpointsnD
150C     **************************************************************************
151C
152      Recursive Subroutine Agrif_TabpointsnD(g,newgrid)
153C
154CCC   Description:
155CCC   Copy on newgrid of points detected on the grid hierarchy pointed by g.   
156C
157CC    Method:
158CC           
159C
160C     Declarations:
161C
162C     Arguments
163      TYPE(Agrif_Grid),pointer   :: g,newgrid
164C
165C     Local variables     
166      TYPE(Agrif_Pgrid),pointer  :: parcours
167C
168      REAL                  :: g_eps,newgrid_eps,eps
169      REAL   , DIMENSION(3) :: newmin,newmax
170      REAL   , DIMENSION(3) :: gmin,gmax
171      REAL   , DIMENSION(3) :: xmin
172      INTEGER, DIMENSION(3) :: igmin,inewmin
173      INTEGER, DIMENSION(3) :: inewmax
174      INTEGER               :: iii
175      INTEGER               :: i,j,k
176      INTEGER               :: i0,j0,k0
177C
178C     
179      parcours => g % child_grids
180C
181      do While(associated(parcours))
182        Call Agrif_TabpointsnD(parcours%gr,newgrid)
183        parcours => parcours % next
184      enddo     
185C
186      g_eps = 10.
187      newgrid_eps = 10.
188C     
189      do iii = 1 , Agrif_probdim
190         g_eps = min( g_eps , g % Agrif_d(iii) )
191         newgrid_eps = min(newgrid_eps,newgrid%Agrif_d(iii))
192      enddo
193C
194      eps = min(g_eps,newgrid_eps)/100.                 
195C                   
196      do iii = 1 , Agrif_probdim
197         if (newgrid%Agrif_d(iii) .LT. (g%Agrif_d(iii)-eps)) Return
198C                   
199         gmin(iii) = g%Agrif_x(iii)
200         gmax(iii) = g%Agrif_x(iii) + g%nb(iii) * g%Agrif_d(iii)
201C 
202         newmin(iii) = newgrid%Agrif_x(iii)
203         newmax(iii) = newgrid%Agrif_x(iii) + newgrid%nb(iii) *
204     &                 newgrid%Agrif_d(iii)
205C     
206         if (gmax(iii) .LT. newmin(iii)) Return
207C
208         if (gmin(iii) .GT. newmax(iii)) Return
209C           
210         inewmin(iii) = 1 - floor(-(max(gmin(iii),newmin(iii))-
211     &                    newmin(iii))
212     &                   /newgrid%Agrif_d(iii)) 
213C         
214         xmin(iii) = newgrid%Agrif_x(iii) + (inewmin(iii)-1)*
215     &              newgrid%Agrif_d(iii)         
216C       
217         igmin(iii) = 1 + nint((xmin(iii)-gmin(iii))/g%Agrif_d(iii))
218C
219         inewmax(iii) = 1 + int((min(gmax(iii),newmax(iii))-
220     &                  newmin(iii))/newgrid%Agrif_d(iii))
221      enddo
222C               
223      if ( Agrif_probdim .EQ. 1 ) then
224         i0 = igmin(1)
225         do i = inewmin(1),inewmax(1)
226            newgrid%tabpoint1D(i) = max(
227     &                           newgrid%tabpoint1D(i),
228     &                           g%tabpoint1D(i0))
229         enddo
230         i0 = i0 + int(newgrid%Agrif_d(1)/g%Agrif_d(1))
231      endif
232C
233      if ( Agrif_probdim .EQ. 2 ) then
234         i0 = igmin(1)
235         do i = inewmin(1),inewmax(1)
236            j0 = igmin(2)
237            do j = inewmin(2),inewmax(2)
238               newgrid%tabpoint2D(i,j) = max(
239     &                           newgrid%tabpoint2D(i,j),
240     &                           g%tabpoint2D(i0,j0))
241               j0 = j0 + int(newgrid%Agrif_d(2)/g%Agrif_d(2))
242            enddo
243            i0 = i0 + int(newgrid%Agrif_d(1)/g%Agrif_d(1))
244         enddo
245      endif
246C
247      if ( Agrif_probdim .EQ. 3 ) then
248         i0 = igmin(1)
249         do i = inewmin(1),inewmax(1)
250            j0 = igmin(2)
251            do j = inewmin(2),inewmax(2)
252               k0 = igmin(3)
253               do k = inewmin(3),inewmax(3)
254                  newgrid%tabpoint3D(i,j,k) = max(
255     &                           newgrid%tabpoint3D(i,j,k),
256     &                           g%tabpoint3D(i0,j0,k0))
257                  k0 = k0 + int(newgrid%Agrif_d(3)/g%Agrif_d(3))
258               enddo                                       
259               j0 = j0 + int(newgrid%Agrif_d(2)/g%Agrif_d(2))
260            enddo
261            i0 = i0 + int(newgrid%Agrif_d(1)/g%Agrif_d(1))
262         enddo
263      endif
264C
265      Return
266C
267C
268      End Subroutine Agrif_TabpointsnD
269C
270C     **************************************************************************
271CCC   Subroutine Agrif_ClusterGridnD
272C     **************************************************************************
273C     
274      Subroutine Agrif_ClusterGridnD(g,coarsegrid)
275C
276CCC   Description:
277CCC   Clustering on the grid pointed by g.
278C
279CC    Method:
280CC           
281C
282C     Declarations:
283C
284C     Pointer arguments
285      TYPE(AGRIF_grid)     ,pointer  :: g          ! Pointer on the current grid
286      TYPE(AGRIF_rectangle),pointer  :: coarsegrid
287C
288C     Local variables     
289      TYPE(Agrif_rectangle) :: newrect
290      TYPE(Agrif_Variable)  :: newflag
291C
292      INTEGER               :: i,j,k
293      INTEGER ,DIMENSION(3) :: sx
294      INTEGER               :: bufferwidth,flagpoints
295      INTEGER               :: n1,n2,m1,m2,o1,o2
296      INTEGER               :: iii
297C
298C     
299      bufferwidth = int(Agrif_Minwidth/5.)
300C     
301      do iii = 1 , Agrif_probdim
302         sx(iii) = g % nb(iii) + 1
303      enddo
304C                 
305      if ( Agrif_probdim .EQ. 1 ) then
306         allocate(newflag%iarray1(sx(1)))
307         newflag%iarray1 = 0     
308      endif
309      if ( Agrif_probdim .EQ. 2 ) then
310         allocate(newflag%iarray2(sx(1),sx(2)))
311         newflag%iarray2 = 0     
312      endif
313      if ( Agrif_probdim .EQ. 3 ) then
314         allocate(newflag%iarray3(sx(1),sx(2),sx(3)))
315         newflag%iarray3 = 0     
316      endif
317C
318      flagpoints = 0 
319C     
320      if (bufferwidth>0) then
321C     
322         if ( Agrif_probdim .EQ. 1 ) then
323            do i = bufferwidth,sx(1)-bufferwidth+1     
324               if (g % tabpoint1D(i) .EQ. 1) then
325                  m1 = i - bufferwidth + 1
326                  m2 = i + bufferwidth - 1                 
327                  flagpoints = flagpoints + 1
328                  newflag%iarray1(m1:m2)=1 
329               endif
330            enddo   
331         endif
332C
333        if ( Agrif_probdim .EQ. 2 ) then
334           do i = bufferwidth,sx(1)-bufferwidth+1     
335               do j = bufferwidth,sx(2)-bufferwidth+1
336                  if (g % tabpoint2D(i,j) .EQ. 1) then
337                     n1 = j - bufferwidth + 1
338                     n2 = j + bufferwidth - 1
339                     m1 = i - bufferwidth + 1
340                     m2 = i + bufferwidth - 1                 
341                     flagpoints = flagpoints + 1
342                     newflag%iarray2(m1:m2,n1:n2)=1 
343                  endif 
344               enddo   
345            enddo
346         endif
347C
348        if ( Agrif_probdim .EQ. 3 ) then
349            do i = bufferwidth,sx(1)-bufferwidth+1     
350               do j = bufferwidth,sx(2)-bufferwidth+1
351                  do k = bufferwidth,sx(3)-bufferwidth+1
352                     if (g % tabpoint3D(i,j,k) .EQ. 1) then
353                        o1 = k - bufferwidth + 1
354                        o2 = k + bufferwidth - 1
355                        n1 = j - bufferwidth + 1
356                        n2 = j + bufferwidth - 1
357                        m1 = i - bufferwidth + 1
358                        m2 = i + bufferwidth - 1                 
359                        flagpoints = flagpoints + 1
360                        newflag%iarray3(m1:m2,n1:n2,o1:o2)=1 
361                     endif
362                  enddo 
363               enddo   
364            enddo
365        endif
366      else
367          flagpoints = 1       
368C
369          if ( Agrif_probdim .EQ. 1 ) then
370          newflag%iarray1 = g % tabpoint1D
371          endif
372          if ( Agrif_probdim .EQ. 2 ) then
373          newflag%iarray2 = g % tabpoint2D
374          endif
375          if ( Agrif_probdim .EQ. 3 ) then
376          newflag%iarray3 = g % tabpoint3D
377          endif
378      endif
379C       
380      if (flagpoints .EQ. 0) then       
381          if ( Agrif_probdim .EQ. 1 ) deallocate(newflag%iarray1) 
382          if ( Agrif_probdim .EQ. 2 ) deallocate(newflag%iarray2) 
383          if ( Agrif_probdim .EQ. 3 ) deallocate(newflag%iarray3) 
384          Return
385      endif       
386C     
387      do iii = 1 , Agrif_probdim
388         newrect % imin(iii) = 1
389         newrect % imax(iii) = sx(iii)
390      enddo
391C     
392      Call Agrif_Clusternd(newflag,
393     &                     coarsegrid%childgrids,newrect) 
394C     
395      if ( Agrif_probdim .EQ. 1 ) deallocate(newflag%iarray1)
396      if ( Agrif_probdim .EQ. 2 ) deallocate(newflag%iarray2)
397      if ( Agrif_probdim .EQ. 3 ) deallocate(newflag%iarray3)
398C 
399C   
400      End Subroutine Agrif_ClusterGridnD
401C
402C     **************************************************************************
403CCC   Subroutine Agrif_ClusternD
404C     **************************************************************************
405C
406      Recursive subroutine Agrif_Clusternd(flag,boxlib,oldB) 
407C
408CCC   Description:
409CCC   Clustering on the grid pointed by oldB.
410C
411CC    Method:
412CC           
413C
414C     Declarations:
415C
416C     Arguments 
417      TYPE(Agrif_rectangle) ::  oldB   
418      TYPE(Agrif_Variable)  :: flag
419c      INTEGER,DIMENSION(oldB%imin(1):oldB%imax(1),
420c     &                  oldB%imin(2):oldB%imax(2)) :: flag
421      TYPE(Agrif_lrectangle),pointer :: boxlib     
422C
423C     Local variables     
424      TYPE(Agrif_lrectangle),pointer :: tempbox,parcbox,parcbox2       
425      TYPE(Agrif_rectangle) :: newB,newB2
426      INTEGER :: i,j,k,iii
427      INTEGER,DIMENSION(:),allocatable :: i_sig 
428      INTEGER,DIMENSION(:),allocatable :: j_sig       
429      INTEGER,DIMENSION(:),allocatable :: k_sig       
430      INTEGER,DIMENSION(3) :: ipu,ipl
431      INTEGER,DIMENSION(3) :: istr,islice
432      REAL :: cureff
433      REAL :: neweff
434      INTEGER :: ValMax,ValSum,TailleTab
435      INTEGER :: nbpoints,nbpointsflag 
436      LOGICAL :: test
437C         
438      allocate(i_sig(oldB%imin(1):oldB%imax(1)))
439      if ( Agrif_probdim .GE. 2 ) 
440     & allocate(j_sig(oldB%imin(2):oldB%imax(2)))
441      if ( Agrif_probdim .EQ. 3 ) 
442     & allocate(k_sig(oldB%imin(3):oldB%imax(3)))
443C
444      test = .FALSE.
445      do iii = 1 , Agrif_probdim
446         test = test .OR. ( (oldB%imax(iii)-oldB%imin(iii)+1) 
447     &                    .LT. Agrif_Minwidth)
448      enddo
449      if ( test ) Return
450C     
451      if ( Agrif_probdim .EQ. 1 ) i_sig = flag%iarray1
452      if ( Agrif_probdim .EQ. 2 ) then
453         do i = oldB%imin(1),oldB%imax(1)
454           i_sig(i) = SUM(flag%iarray2(i,
455     &                            oldB%imin(2):oldB%imax(2)))
456         enddo
457         do j = oldB%imin(2),oldB%imax(2)
458           j_sig(j) = SUM(flag%iarray2(
459     &                            oldB%imin(1):oldB%imax(1),j))
460         enddo     
461      endif
462      if ( Agrif_probdim .EQ. 3 ) then
463         do i = oldB%imin(1),oldB%imax(1)
464            i_sig(i) = SUM(flag%iarray3(i,
465     &                            oldB%imin(2):oldB%imax(2),
466     &                            oldB%imin(3):oldB%imax(3)))
467         enddo
468         do j = oldB%imin(2),oldB%imax(2)
469            j_sig(j) = SUM(flag%iarray3(
470     &                          oldB%imin(1):oldB%imax(1),j,
471     &                          oldB%imin(3):oldB%imax(3)))
472         enddo
473         do k = oldB%imin(3),oldB%imax(3)
474            k_sig(k) = SUM(flag%iarray3(
475     &                          oldB%imin(1):oldB%imax(1),
476     &                          oldB%imin(2):oldB%imax(2),k))
477         enddo     
478      endif
479C                   
480      do iii = 1 , Agrif_probdim
481         ipl(iii) = oldB%imin(iii)
482         ipu(iii) = oldB%imax(iii)
483      enddo
484C           
485      Call Agrif_Clusterprune(i_sig,ipl(1),ipu(1)) 
486      if ( Agrif_probdim .GE. 2 ) 
487     &   Call Agrif_Clusterprune(j_sig,ipl(2),ipu(2))     
488      if ( Agrif_probdim .EQ. 3 ) 
489     &   Call Agrif_Clusterprune(k_sig,ipl(3),ipu(3))     
490C     
491      test = .TRUE.
492      do iii = 1 , Agrif_probdim
493         test = test .AND. (ipl(iii).EQ.oldB%imin(iii))
494         test = test .AND. (ipu(iii).EQ.oldB%imax(iii))
495      enddo
496     
497      if (.NOT. test) then
498          do iii = 1 , Agrif_probdim
499             newB%imin(iii) = ipl(iii)
500             newB%imax(iii) = ipu(iii)
501          enddo
502C       
503          if ( Agrif_probdim .EQ. 1 ) 
504     &       nbpoints = SUM(flag%iarray1(newB%imin(1):newB%imax(1)))   
505          if ( Agrif_probdim .EQ. 2 ) 
506     &       nbpoints = SUM(flag%iarray2(newB%imin(1):newB%imax(1),
507     &                           newB%imin(2):newB%imax(2)))   
508          if ( Agrif_probdim .EQ. 3 ) 
509     &       nbpoints = SUM(flag%iarray3(newB%imin(1):newB%imax(1),
510     &                           newB%imin(2):newB%imax(2),   
511     &                           newB%imin(3):newB%imax(3)))
512C   
513          if ( Agrif_probdim .EQ. 1 ) 
514     &       TailleTab = newB%imax(1)-newB%imin(1)+1
515          if ( Agrif_probdim .EQ. 2 ) 
516     &       TailleTab = (newB%imax(1)-newB%imin(1)+1)*
517     &                   (newB%imax(2)-newB%imin(2)+1)
518          if ( Agrif_probdim .EQ. 3 ) 
519     &       TailleTab = (newB%imax(1)-newB%imin(1)+1)*
520     &                   (newB%imax(2)-newB%imin(2)+1)*
521     &                   (newB%imax(3)-newB%imin(3)+1)
522C
523          neweff = REAL(nbpoints)/TailleTab
524C     
525          if (nbpoints.GT.0) then
526C       
527              if ((neweff .GT .Agrif_efficiency)) then
528                  Call Agrif_Add_Rectangle(newB,boxlib)
529                  Return
530              else
531C           
532                  tempbox => boxlib
533                  newB2 = newB
534                  Call Agrif_Clusternd(flag,
535     &                      boxlib,newB)
536C     
537C                 Compute new efficiency
538C
539                  cureff = neweff
540                  parcbox2 => boxlib
541                  nbpoints = 0
542                  nbpointsflag = 0
543C
544                  do While (associated(parcbox2))
545                    if (associated(parcbox2,tempbox)) Exit
546                    newB = parcbox2%r
547C
548                    if ( Agrif_probdim .EQ. 1 ) Valsum = 
549     &                                 SUM(flag%iarray1(
550     &                                 newB%imin(1):newB%imax(1)))
551                    if ( Agrif_probdim .EQ. 2 ) Valsum = 
552     &                                 SUM(flag%iarray2(
553     &                                 newB%imin(1):newB%imax(1),
554     &                                 newB%imin(2):newB%imax(2)))
555                    if ( Agrif_probdim .EQ. 3 ) Valsum = 
556     &                                 SUM(flag%iarray3(
557     &                                 newB%imin(1):newB%imax(1),
558     &                                 newB%imin(2):newB%imax(2),
559     &                                 newB%imin(3):newB%imax(3)))
560C
561                    nbpointsflag = nbpointsflag + ValSum
562                    if ( Agrif_probdim .EQ. 1 ) 
563     &                 TailleTab = newB%imax(1)-newB%imin(1)+1
564                    if ( Agrif_probdim .EQ. 2 ) 
565     &                 TailleTab = (newB%imax(1)-newB%imin(1)+1)*
566     &                             (newB%imax(2)-newB%imin(2)+1)
567                    if ( Agrif_probdim .EQ. 3 ) 
568     &                 TailleTab = (newB%imax(1)-newB%imin(1)+1)*
569     &                             (newB%imax(2)-newB%imin(2)+1)*
570     &                             (newB%imax(3)-newB%imin(3)+1)
571                    nbpoints = nbpoints + TailleTab
572                    parcbox2 => parcbox2%next 
573                  enddo 
574C coefficient 1.05 avant 1.15 possibilité de laisser choix à l utilisateur
575                  if  (REAL(nbpointsflag)/REAL(nbpoints)
576     &                 .LT.(1.0001*cureff)) then
577                      parcbox2 => boxlib           
578                      do While (associated(parcbox2))
579                        if (associated(parcbox2,tempbox)) Exit
580                        deallocate(parcbox2%r)
581                        parcbox2 => parcbox2%next             
582                      enddo     
583                      boxlib => tempbox
584                      Call Agrif_Add_Rectangle(newB2,boxlib)
585                      Return
586                  endif
587              endif         
588          endif
589          Return
590      endif
591C       
592      do iii = 1 , Agrif_Probdim
593         istr(iii) = oldB%imax(iii)
594         islice(iii) = oldB%imin(iii)
595      enddo     
596C     
597      Call Agrif_Clusterslice(i_sig,islice(1),istr(1)) 
598      if ( Agrif_probdim .GE. 2 ) 
599     & Call Agrif_Clusterslice(j_sig,islice(2),istr(2)) 
600      if ( Agrif_probdim .EQ. 3 ) 
601     & Call Agrif_Clusterslice(k_sig,islice(3),istr(3)) 
602C     
603      ValSum = 0
604      do iii = 1 , Agrif_Probdim
605         Valsum = valSum + islice(iii)
606      enddo
607C     
608      if ( Valsum .EQ. -Agrif_Probdim ) then
609          Call Agrif_Add_Rectangle(oldB,boxlib) 
610          Return
611      endif
612C       
613      nullify(tempbox)
614      tempbox => boxlib
615      if ( Agrif_probdim .EQ. 1 ) 
616     &      cureff  = oldB%imax(1)-oldB%imin(1)+1
617      if ( Agrif_probdim .EQ. 2 ) 
618     &      cureff  = (oldB%imax(1)-oldB%imin(1)+1)*
619     &                (oldB%imax(2)-oldB%imin(2)+1)
620      if ( Agrif_probdim .EQ. 3 ) 
621     &       cureff = (oldB%imax(1)-oldB%imin(1)+1)*
622     &                (oldB%imax(2)-oldB%imin(2)+1)*
623     &                (oldB%imax(3)-oldB%imin(3)+1)
624      Nullify(parcbox)
625C 
626      do iii = 1 , Agrif_Probdim
627          newB%imax(iii) = oldB%imax(iii)           
628          newB%imin(iii) = oldB%imin(iii)
629      enddo     
630C
631      ValMax = 0
632      do iii = 1 , Agrif_Probdim
633         ValMax = Max(ValMax,istr(iii))
634      enddo
635C
636      if (istr(1) .EQ. ValMax ) then
637          newB%imax(1) = islice(1)
638          Call Agrif_Add_Rectangle(newB,parcbox)                   
639          newB%imin(1) = islice(1)+1
640          newB%imax(1) = oldB%imax(1)
641          Call Agrif_Add_Rectangle(newB,parcbox)
642      elseif ( Agrif_probdim .GE. 2 ) then
643         if (istr(2) .EQ. ValMax ) then
644            newB%imax(2) = islice(2)
645            Call Agrif_Add_Rectangle(newB,parcbox)
646            newB%imin(2) = islice(2)+1
647            newB%imax(2) = oldB%imax(2)
648            Call Agrif_Add_Rectangle(newB,parcbox)
649         elseif ( Agrif_probdim .EQ. 3 ) then
650            newB%imax(3) = islice(3)
651            Call Agrif_Add_Rectangle(newB,parcbox)
652            newB%imin(3) = islice(3)+1
653            newB%imax(3) = oldB%imax(3)
654            Call Agrif_Add_Rectangle(newB,parcbox)
655         endif
656      endif
657C     
658      do While (associated(parcbox))
659        newB = parcbox%r           
660C
661        if ( Agrif_probdim .EQ. 1 ) nbpoints = 
662     &                                 SUM(flag%iarray1(
663     &                                 newB%imin(1):newB%imax(1)))
664        if ( Agrif_probdim .EQ. 2 ) nbpoints = 
665     &                                 SUM(flag%iarray2(
666     &                                 newB%imin(1):newB%imax(1),
667     &                                 newB%imin(2):newB%imax(2)))
668        if ( Agrif_probdim .EQ. 3 ) nbpoints = 
669     &                                 SUM(flag%iarray3(
670     &                                 newB%imin(1):newB%imax(1),
671     &                                 newB%imin(2):newB%imax(2),
672     &                                 newB%imin(3):newB%imax(3)))
673C     
674       if ( Agrif_probdim .EQ. 1 ) 
675     &                 TailleTab = newB%imax(1)-newB%imin(1)+1
676       if ( Agrif_probdim .EQ. 2 ) 
677     &                 TailleTab = (newB%imax(1)-newB%imin(1)+1)*
678     &                             (newB%imax(2)-newB%imin(2)+1)
679       if ( Agrif_probdim .EQ. 3 ) 
680     &                 TailleTab = (newB%imax(1)-newB%imin(1)+1)*
681     &                             (newB%imax(2)-newB%imin(2)+1)*
682     &                             (newB%imax(3)-newB%imin(3)+1)
683
684        neweff = REAL(nbpoints) / TailleTab
685C     
686        if (nbpoints .GT. 0) then
687C     
688            if ((neweff .GT .Agrif_efficiency)) then
689                Call Agrif_Add_Rectangle(newB,boxlib)
690            else
691                tempbox => boxlib
692                newB2 = newB
693                Call Agrif_Clusternd(flag,
694     &                      boxlib,newB)
695C
696C               Compute new efficiency
697C
698                cureff = neweff
699                parcbox2 => boxlib
700                nbpoints = 0
701                nbpointsflag = 0
702C
703                do While (associated(parcbox2))
704                  if (associated(parcbox2,tempbox)) Exit
705                  newB = parcbox2%r
706C
707                  if ( Agrif_probdim .EQ. 1 ) ValSum = 
708     &                                 SUM(flag%iarray1(
709     &                                 newB%imin(1):newB%imax(1)))
710                  if ( Agrif_probdim .EQ. 2 ) ValSum = 
711     &                                 SUM(flag%iarray2(
712     &                                 newB%imin(1):newB%imax(1),
713     &                                 newB%imin(2):newB%imax(2)))
714                  if ( Agrif_probdim .EQ. 3 ) ValSum = 
715     &                                 SUM(flag%iarray3(
716     &                                 newB%imin(1):newB%imax(1),
717     &                                 newB%imin(2):newB%imax(2),
718     &                                 newB%imin(3):newB%imax(3)))
719C
720                  nbpointsflag = nbpointsflag + ValSum
721C
722                    if ( Agrif_probdim .EQ. 1 ) 
723     &                 TailleTab = newB%imax(1)-newB%imin(1)+1
724                    if ( Agrif_probdim .EQ. 2 ) 
725     &                 TailleTab = (newB%imax(1)-newB%imin(1)+1)*
726     &                             (newB%imax(2)-newB%imin(2)+1)
727                    if ( Agrif_probdim .EQ. 3 ) 
728     &                 TailleTab = (newB%imax(1)-newB%imin(1)+1)*
729     &                             (newB%imax(2)-newB%imin(2)+1)*
730     &                             (newB%imax(3)-newB%imin(3)+1)
731
732                  nbpoints = nbpoints + TailleTab
733C             
734                  parcbox2 => parcbox2%next 
735                enddo   
736C           
737                if  (REAL(nbpointsflag)/REAL(nbpoints)
738     &               .LT.(1.15*cureff)) then
739                    parcbox2 => boxlib           
740                    do While (associated(parcbox2))
741                      if (associated(parcbox2,tempbox)) Exit
742                      deallocate(parcbox2%r)
743                      parcbox2 => parcbox2%next             
744                    enddo
745                    boxlib => tempbox
746                    Call Agrif_Add_Rectangle(newB2,boxlib)
747                endif
748            endif
749        endif
750        parcbox => parcbox%next
751      enddo
752C     
753C     
754      Return
755C     
756      End Subroutine Agrif_Clusternd 
757C
758C     **************************************************************************
759CCC   Subroutine Agrif_Clusterslice
760C     **************************************************************************
761C                           
762      Subroutine Agrif_Clusterslice(sig,slice,str) 
763C
764C
765CCC   Description:
766CCC   
767C
768CC    Method:
769CC           
770C
771C     Declarations:
772C
773C     Arguments 
774      INTEGER                      :: slice,str
775      INTEGER,DIMENSION(slice:str) :: sig 
776C
777C     Local variables     
778      INTEGER                      :: ideb,ifin
779      INTEGER                      :: i,t,a,di,ts 
780      INTEGER,DIMENSION(slice:str) :: lap 
781C
782C   
783      ideb = slice
784      ifin = str
785C   
786      if (SIZE(sig) <= 2*Agrif_Minwidth) then
787          str = -1 
788          slice = -1 
789          Return
790      endif
791C     
792      t = -1 
793      a = -1 
794C
795      do i = ideb+Agrif_Minwidth-1,ifin-Agrif_Minwidth
796          if (sig(i) .EQ. 0) then
797              if ((i-ideb) < (ifin-i)) then
798                  di = i - ideb
799                else
800                  di = ifin - i 
801              endif
802C
803              if (di > t) then
804                a = i 
805                t = di 
806              endif
807         endif
808      enddo
809C     
810      if (a .NE. (-1)) then
811          slice = a 
812          str = t 
813          Return
814      endif
815C     
816      t = -1 
817      a = -1 
818C     
819      do i = ideb+1,ifin-1 
820        lap(i) = sig(i+1) + sig(i-1) - 2*sig(i) 
821      enddo
822C     
823      do i = ideb + Agrif_Minwidth-1,ifin-Agrif_Minwidth
824        if ((lap(i+1)*lap(i)) .LE. 0) then
825            ts = ABS(lap(i+1) - lap(i)) 
826            if (ts > t) then
827                t = ts 
828                a = i
829            endif
830       endif
831      enddo
832C     
833      if (a .EQ. (ideb + Agrif_Minwidth - 1)) then
834          a = -1
835          t = -1
836      endif
837C     
838      slice = a 
839      str = t 
840C
841C     
842      End Subroutine Agrif_Clusterslice       
843C
844C
845C   
846C     **************************************************************************
847CCC   Subroutine Agrif_Clusterprune
848C     **************************************************************************
849C     
850      Subroutine Agrif_Clusterprune(sig,pl,pu) 
851C
852C
853CCC   Description:
854CCC   
855C
856CC    Method:
857CC           
858C
859C     Declarations:
860C
861C     Arguments
862      INTEGER                  :: pl,pu 
863      INTEGER,DIMENSION(pl:pu) :: sig
864C
865C     Local variables     
866      INTEGER :: ideb,ifin       
867      INTEGER :: diff,addl,addu,udist,ldist 
868C
869C     
870      ideb = pl
871      ifin = pu
872C     
873      if (SIZE(sig) <= Agrif_Minwidth) then
874          return
875      endif
876C       
877      do While ((sig(pl) .EQ. 0) .AND. (pl < ifin)) 
878        pl = pl + 1 
879      enddo 
880C       
881      do While ((sig(pu) .EQ. 0) .AND. (pu > ideb)) 
882        pu = pu - 1 
883      enddo 
884C       
885      if ((pu-pl) < Agrif_Minwidth) then
886          diff = Agrif_Minwidth - (pu - pl + 1) 
887          udist = ifin - pu 
888          ldist = pl - ideb
889          addl = diff / 2 
890          addu = diff - addl 
891          if (addu > udist) then
892              addu = udist 
893              addl = diff - addu 
894          endif
895C         
896          if (addl > ldist) then
897              addl = ldist 
898              addu = diff - addl
899          endif
900C         
901          pu = pu + addu 
902          pl = pl - addl         
903C     
904      endif
905C       
906C
907      End Subroutine Agrif_Clusterprune
908C
909C
910C             
911C     **************************************************************************
912CCC   Subroutine Agrif_Add_Rectangle 
913C     **************************************************************************
914C       
915      Subroutine Agrif_Add_Rectangle(R,LR)
916C
917CCC   Description:
918CCC   Subroutine to add the Agrif_Rectangle R in a list managed by LR.
919C
920C     Declarations:
921C
922C     Arguments               
923      TYPE(AGRIF_rectangle)           :: R 
924      TYPE(AGRIF_lrectangle), Pointer :: LR
925C
926C     Local variable       
927      TYPE(AGRIF_lrectangle), Pointer :: newrect 
928C
929      INTEGER                         :: iii 
930C
931C       
932      allocate(newrect)
933      allocate(newrect % r)
934C
935      newrect % r = R
936C     
937      do iii = 1 , Agrif_Probdim
938         newrect % r % spaceref(iii) = Agrif_Coeffref(iii)
939         newrect % r % timeref(iii) = Agrif_Coeffreft(iii)
940      enddo
941C             
942      newrect % r % number = -1           
943      Nullify(newrect % r % childgrids)
944      newrect % next => LR 
945      LR => newrect 
946C
947C     
948      End Subroutine Agrif_Add_Rectangle   
949C 
950C
951C 
952C     **************************************************************************
953CCC   Subroutine Agrif_Read_Fix_Grd 
954C     **************************************************************************
955C     
956      Recursive Subroutine Agrif_Read_Fix_Grd(coarsegrid,j,nunit)
957C
958CCC   Description:
959CCC   Subroutine to create the grid hierarchy from the reading of the 
960CCC   AGRIF_FixedGrids.in file.
961C
962CC    Method:
963CC    Recursive subroutine and creation of a first grid hierarchy from the
964CC    reading of the AGRIF_FixedGrids.in file.       
965C
966C     Declarations:
967C
968C     Pointer argument     
969      TYPE(AGRIF_rectangle), Pointer   :: coarsegrid ! Pointer on the first grid
970                                                     ! of the grid hierarchy
971C                                   
972C     Scalar arguments                                   
973      INTEGER                          :: j          ! Number of the new grid
974      INTEGER                          :: nunit      ! unit associated with file
975C     
976C     Local variables     
977      TYPE(AGRIF_rectangle)            :: newrect    ! Pointer on a new grid   
978      TYPE(AGRIF_lrectangle), Pointer  :: parcours   ! Pointer for the recursive
979                                                     !    procedure
980      TYPE(AGRIF_lrectangle), Pointer  :: newlrect
981      TYPE(AGRIF_lrectangle), Pointer  :: end_list       
982      INTEGER                          :: i          ! for each child grid
983      INTEGER                          :: nb_grids   ! Number of child grids
984      INTEGER                          :: iii
985C         
986C
987      Nullify(newrect%childgrids)
988C     
989C     Reading of the number of child grids       
990      read(nunit,*) nb_grids
991C     
992C     coarsegrid%nbgridchild = nb_grids
993C     
994      allocate(end_list)
995C
996      nullify(end_list % r)
997      nullify(end_list % next)
998C     
999      coarsegrid % childgrids => end_list
1000C
1001C     Reading of imin(1),imax(1),imin(2),imax(2),imin(3),imax(3), and space and
1002C        time refinement factors for each child grid.
1003C     Creation and addition of the new grid to the grid hierarchy. 
1004C
1005      do i = 1,nb_grids
1006        allocate(newlrect)     
1007        newrect % number = j   ! Number of the grid
1008C
1009        if ( Agrif_USE_ONLY_FIXED_GRIDS .EQ. 0 ) then
1010           if (Agrif_Probdim == 3) then 
1011            read(nunit,*) newrect % imin(1), newrect % imax(1),
1012     &                    newrect % imin(2), newrect % imax(2),
1013     &                    newrect % imin(3), newrect % imax(3),
1014     &                    newrect % spaceref(1),newrect % spaceref(2),
1015     &                    newrect % spaceref(3),
1016     &                    newrect % timeref(1),newrect % timeref(2),
1017     &                    newrect % timeref(3)
1018           elseif (Agrif_Probdim == 2) then           
1019            read(nunit,*) newrect % imin(1),newrect % imax(1),
1020     &                  newrect % imin(2),newrect % imax(2),
1021     &                    newrect % spaceref(1),newrect % spaceref(2),
1022     &                    newrect % timeref(1),newrect % timeref(2)
1023           elseif (Agrif_Probdim == 1) then
1024            read(nunit,*) newrect % imin(1), newrect % imax(1),
1025     &                    newrect % spaceref(1),
1026     &                    newrect % timeref(1)
1027           endif
1028        else
1029           if (Agrif_Probdim == 3) then 
1030            read(nunit,*) newrect % imin(1), newrect % imax(1),
1031     &                    newrect % imin(2), newrect % imax(2),
1032     &                    newrect % imin(3), newrect % imax(3),
1033     &                    newrect % spaceref(1),newrect % spaceref(2),
1034     &                    newrect % spaceref(3),
1035     &                    newrect % timeref(1)
1036           elseif (Agrif_Probdim == 2) then           
1037            read(nunit,*) newrect % imin(1),newrect % imax(1),
1038     &                  newrect % imin(2),newrect % imax(2),
1039     &                    newrect % spaceref(1),newrect % spaceref(2),
1040     &                    newrect % timeref(1)
1041           elseif (Agrif_Probdim == 1) then
1042            read(nunit,*) newrect % imin(1), newrect % imax(1),
1043     &                    newrect % spaceref(1),
1044     &                    newrect % timeref(1)
1045           endif
1046C
1047           if ( Agrif_probdim .GE. 2 ) then
1048              do iii = 2 , Agrif_probdim
1049                 newrect % timeref(iii) = newrect % timeref(1) 
1050              enddo
1051           endif
1052C
1053        endif
1054C
1055C       Addition to the grid hierarchy
1056C
1057        nullify(newrect % childgrids)
1058        j = j + 1
1059        Allocate(newlrect%r)
1060        newlrect % r = newrect 
1061        nullify(newlrect % next)
1062        end_list % next => newlrect
1063        end_list => end_list % next
1064      enddo
1065C     
1066      coarsegrid % childgrids => coarsegrid % childgrids % next
1067      parcours => coarsegrid % childgrids
1068C
1069C     Recursive operation to create the grid hierarchy branch by branch
1070C
1071      do while (associated(parcours))
1072        call Agrif_Read_Fix_Grd (parcours % r,j,nunit)
1073        parcours => parcours % next
1074      enddo
1075C     
1076C
1077      End Subroutine Agrif_Read_Fix_Grd       
1078C       
1079C 
1080C
1081C     **************************************************************************
1082CCC   Subroutine Agrif_Create_Grids 
1083C     **************************************************************************
1084C
1085      Recursive Subroutine Agrif_Create_Grids(g,coarsegrid)
1086C
1087CCC   Description:
1088CCC   Subroutine to create the grid hierarchy (g) from the one created with the
1089CCC   Agrif_Read_Fix_Grd or Agrif_Cluster_All procedures (coarsegrid).
1090C
1091CC    Method:
1092CC    Recursive subroutine.       
1093C
1094C     Declarations:
1095C
1096C     Pointer arguments       
1097      TYPE(AGRIF_grid)     , Pointer  :: g          ! Pointer on the root coarse
1098                                                    ! grid
1099      TYPE(AGRIF_rectangle), Pointer  :: coarsegrid ! Pointer on the root coarse
1100                                                    ! grid of the grid hierarchy
1101                                                    ! created with the
1102                                                    ! Agrif_Read_Fix_Grd
1103                                                    ! subroutine
1104C
1105C     Local pointers
1106      TYPE(Agrif_grid)      , Pointer :: newgrid    ! New grid
1107      TYPE(Agrif_pgrid)     , Pointer :: newpgrid
1108      TYPE(Agrif_pgrid)     , Pointer :: parcours2
1109      TYPE(Agrif_lrectangle), Pointer :: parcours
1110      TYPE(Agrif_pgrid)     , Pointer :: end_list
1111      TYPE(Agrif_pgrid)     , Pointer :: parcours3
1112C
1113C     Local scalars     
1114      LOGICAL                         :: nullliste
1115      INTEGER                         :: iii
1116      INTEGER                         :: moving_grid_id = 0     
1117
1118C
1119      parcours3 => g % child_grids
1120C     
1121      if (associated(parcours3)) then
1122          do While (associated(parcours3 % next))
1123            parcours3 => parcours3 % next
1124          enddo
1125          end_list => parcours3
1126          nullliste=.FALSE.
1127        else
1128          allocate(end_list)
1129          nullify(end_list % gr)
1130          nullify(end_list % next)     
1131          g % child_grids => end_list 
1132          parcours3 => end_list     
1133          nullliste=.TRUE.
1134      endif
1135C     
1136      parcours => coarsegrid % childgrids
1137C     
1138C     Creation of the grid hierarchy from the one created by using the 
1139C     Agrif_Read_Fix_Grd subroutine
1140C
1141      do while (associated(parcours))
1142        allocate(newgrid)       
1143        moving_grid_id=moving_grid_id+1
1144        newgrid % grid_id = moving_grid_id
1145        do iii = 1 , Agrif_Probdim
1146           newgrid % spaceref(iii) = parcours % r % spaceref(iii)
1147           newgrid % timeref(iii) = parcours % r % timeref(iii)
1148        enddo
1149C
1150        do iii = 1 , Agrif_Probdim
1151          newgrid % nb(iii) = (parcours % r % imax(iii) 
1152     &                       - parcours % r % imin(iii)) 
1153     &                       * parcours % r % spaceref(iii)
1154C     
1155          newgrid % ix(iii) =  parcours % r % imin(iii)
1156C
1157          newgrid % Agrif_d(iii) = g % Agrif_d(iii) 
1158     &                  / REAL(newgrid % spaceref(iii))
1159C
1160          newgrid % Agrif_x(iii) = g % Agrif_x(iii) +
1161     &      (parcours % r % imin(iii) - 1)* g % Agrif_d(iii)
1162C     
1163        enddo
1164C
1165C       Pointer on the parent grid                       
1166C
1167        newgrid % parent => g     
1168
1169C       Level of the current grid
1170        newgrid % level = newgrid % parent % level + 1
1171        if (newgrid % level > Agrif_MaxLevelLoc) then
1172          Agrif_MaxLevelLoc = newgrid%level
1173        endif
1174
1175C     
1176C       Grid pointed by newgrid is a fixed grid     
1177C
1178        if (parcours % r % number .GT. 0) then     
1179            newgrid % fixed = .true. 
1180          else
1181            newgrid % fixed = .false.
1182        endif
1183C
1184C       Number of the grid pointed by newgrid
1185        newgrid % fixedrank = parcours % r % number
1186C           
1187C       No time calculation on this grid         
1188        newgrid % ngridstep = 0     
1189C
1190C       Test indicating if the current grid has a common border with the root 
1191C       coarse grid in the x direction
1192        do iii = 1 , Agrif_Probdim
1193           newgrid % NearRootBorder(iii) = .FALSE.
1194C       
1195           if ((newgrid % parent % NearRootBorder(iii)) .AND. 
1196     &         (newgrid % ix(iii) == 1)) then
1197               newgrid % NearRootBorder(iii) = .TRUE.
1198           endif
1199C
1200           newgrid % DistantRootBorder(iii) = .FALSE.
1201C 
1202           if ((newgrid % parent % DistantRootBorder(iii)) .AND. 
1203     &         (newgrid % ix(iii) + 
1204     &         (newgrid % nb(iii)/newgrid % spaceref(iii)) 
1205     &          - 1  == newgrid % parent % nb(iii))) then
1206               newgrid % DistantRootBorder(iii) = .TRUE.
1207           endif
1208        enddo
1209C
1210C       Writing in output files
1211C
1212        newgrid % oldgrid = .FALSE.     
1213C
1214C     
1215C       Definition of the CHARACTERistics of the variable of the grid pointed by
1216C       newgrid 
1217        Call Agrif_Create_Var (newgrid)
1218C     
1219C       Instanciation of the grid pointed by newgrid and its variables
1220        Call Agrif_Instance (newgrid)
1221C       
1222C       Nullify the variable of the grid pointed by newgrid
1223C
1224C
1225C       Addition of this grid to the grid hierarchy 
1226C     
1227        nullify(newgrid % child_grids)
1228        allocate(newpgrid) 
1229        newpgrid % gr => newgrid
1230        nullify(newpgrid % next)
1231        end_list % next => newpgrid
1232        end_list => end_list % next               
1233        parcours => parcours % next
1234C
1235C       Updating the total number of fixed grids
1236        if (newgrid % fixed) then
1237            AGRIF_nbfixedgrids = AGRIF_nbfixedgrids + 1
1238        endif
1239C                 
1240      enddo
1241C
1242C
1243      if (nullliste) then
1244          g % child_grids => g % child_grids % next
1245          parcours2 => g % child_grids
1246          deallocate(parcours3)
1247        else
1248          parcours2 => parcours3 % next
1249      endif
1250C     
1251      parcours => coarsegrid % childgrids
1252C
1253C     Recursive call to the subroutine Agrif_Create_Fixed_Grids to create the
1254C     grid hierarchy
1255C     
1256      do while (associated(parcours))
1257        Call Agrif_Create_Grids (parcours2 % gr,parcours % r)
1258        parcours => parcours % next
1259        parcours2 => parcours2 % next
1260      enddo 
1261C     
1262      Return     
1263C                 
1264      End Subroutine Agrif_Create_Grids
1265C
1266C     
1267C
1268C     **************************************************************************
1269CCC   Subroutine Agrif_Init_Hierarchy 
1270C     **************************************************************************
1271C
1272      Recursive Subroutine Agrif_Init_Hierarchy(g) 
1273C
1274CCC   Description:
1275CCC   Subroutine to initialize all the grids except the root coarse grid (this 
1276CCC   one, pointed by AGRIF_mygrid, is initialized by the subroutine
1277CCC   Agrif_Init_Grids defi ned in the module Agrif_Util and called in the main 
1278CCC   program ). 
1279C
1280CC    Method:
1281CC    Recursive subroutine.       
1282C
1283C     Declarations:
1284C
1285C     Pointer argument     
1286      TYPE(AGRIF_grid), Pointer  :: g         ! Pointer on the root coarse grid
1287C     
1288C     Local variables     
1289      TYPE(AGRIF_pgrid), Pointer :: parcours  ! Pointer for the recursive call
1290      LOGICAL                    :: Init_Hierarchy
1291C     
1292C
1293      parcours=>g%child_grids
1294C     
1295      do while (associated(parcours))
1296        Init_Hierarchy = .false.
1297        if ( AGRIF_USE_FIXED_GRIDS .EQ. 1 .OR. 
1298     &       AGRIF_USE_ONLY_FIXED_GRIDS .EQ. 1 ) then
1299           if ((parcours%gr%fixed) 
1300     &         .AND. (Agrif_mygrid%ngridstep == 0)) then
1301              Init_Hierarchy = .true.
1302           endif     
1303       endif
1304C
1305       if (.NOT. parcours%gr%fixed) Init_Hierarchy = .true.
1306       if (parcours % gr % oldgrid) Init_Hierarchy = .false.
1307C
1308       if (Init_Hierarchy) then
1309C
1310C           Instanciation of the grid pointed by parcours%gr and its variables 
1311            Call Agrif_Instance (parcours % gr)
1312C   
1313C           Allocation of the arrays containing values of the variables of the
1314C           grid pointed by parcours%gr     
1315C     
1316            Call Agrif_Allocation (parcours % gr)   
1317C           
1318            Call Agrif_initialisations(parcours % gr)             
1319C       
1320            Call Agrif_Instance(parcours % gr)
1321C
1322            if ( Agrif_USE_ONLY_FIXED_GRIDS .EQ. 0 ) then
1323              Call Agrif_Allocate_Restore (parcours % gr)
1324            endif
1325C
1326            if ( Agrif_USE_ONLY_FIXED_GRIDS .EQ. 0 ) then
1327C              Initialization by copy of the grids created by clustering 
1328               Call AGRIF_CopyFromold_All (parcours%gr,
1329     &                                     Agrif_oldmygrid)
1330            endif
1331C
1332C           Initialization by interpolation 
1333C           (this routine is written by the user) 
1334            Call Agrif_InitValues()
1335C
1336            if ( Agrif_USE_ONLY_FIXED_GRIDS .EQ. 0 ) then
1337               Call Agrif_Free_Restore (parcours % gr)
1338            endif
1339C
1340       endif                                 
1341       parcours => parcours % next 
1342C         
1343      enddo
1344C
1345      parcours => g % child_grids
1346C     
1347C     Recursive operation to initialize all the grids
1348      do while (associated(parcours))
1349        Call Agrif_Init_Hierarchy (parcours%gr)
1350        parcours => parcours%next
1351      enddo
1352C     
1353      End Subroutine Agrif_Init_Hierarchy 
1354C
1355C     **************************************************************************
1356CCC   Subroutine Agrif_Allocate_Restore 
1357C     **************************************************************************
1358C
1359      Subroutine Agrif_Allocate_Restore(Agrif_Gr)
1360C     
1361C     
1362C     Modules used:
1363C
1364      TYPE(AGRIF_grid), Pointer  :: Agrif_Gr   ! Pointer on the root coarse grid
1365C     
1366      INTEGER                    :: i
1367C
1368        do i = 1 , Agrif_NbVariables
1369            if ( Agrif_Mygrid%tabvars(i)%var % restaure ) then
1370            if ( Agrif_Gr%tabvars(i)%var % nbdim .EQ. 1 ) then   
1371            Allocate( Agrif_Gr%tabvars(i)%var% 
1372     &    Restore1D(lbound(Agrif_Gr%tabvars(i)%var%array1,1) 
1373     &    :ubound(Agrif_Gr%tabvars(i)%var%array1,1))) 
1374            Agrif_Gr%tabvars(i)%var%Restore1D = 0
1375C
1376            endif
1377            if ( Agrif_Gr%tabvars(i)%var % nbdim .EQ. 2 ) then
1378            Allocate( Agrif_Gr%tabvars(i)%var%Restore2D( 
1379     &      lbound(Agrif_Gr%tabvars(i)%var%array2,1):         
1380     &      ubound(Agrif_Gr%tabvars(i)%var%array2,1),         
1381     &      lbound(Agrif_Gr%tabvars(i)%var%array2,2)         
1382     &      :ubound(Agrif_Gr%tabvars(i)%var%array2,2)))         
1383            Agrif_Gr%tabvars(i)%var%Restore2D = 0
1384C
1385            endif
1386            if ( Agrif_Mygrid%tabvars(i)%var % nbdim .EQ. 3 ) then
1387C     
1388            Allocate( Agrif_Gr%tabvars(i)%var%Restore3D( 
1389     &      lbound(Agrif_Gr%tabvars(i)%var%array3,1):       
1390     &      ubound(Agrif_Gr%tabvars(i)%var%array3,1),       
1391     &      lbound(Agrif_Gr%tabvars(i)%var%array3,2):       
1392     &      ubound(Agrif_Gr%tabvars(i)%var%array3,2),       
1393     &      lbound(Agrif_Gr%tabvars(i)%var%array3,3):       
1394     &      ubound(Agrif_Gr%tabvars(i)%var%array3,3))) 
1395            Agrif_Gr%tabvars(i)%var%Restore3D = 0
1396            endif
1397C 
1398            endif
1399          enddo
1400C
1401      Return
1402C
1403C     
1404      End Subroutine Agrif_Allocate_Restore     
1405C
1406C
1407C
1408C
1409C     **************************************************************************
1410CCC   Subroutine Agrif_Free_Restore 
1411C     **************************************************************************
1412C
1413      Subroutine Agrif_Free_Restore(Agrif_Gr)
1414C
1415C     
1416C     Pointer argument     
1417      TYPE(AGRIF_grid), Pointer  :: Agrif_Gr   ! Pointer on the root coarse grid
1418      INTEGER :: i     
1419C
1420      do i = 1 , Agrif_NbVariables
1421         if ( Agrif_Mygrid % tabvars(i) % var % restaure) then
1422C 
1423            if (associated(Agrif_Gr%tabvars(i)%var%Restore1D)) then
1424               Deallocate(Agrif_Gr%tabvars(i)%var%Restore1D)
1425            endif
1426            if (associated(Agrif_Gr%tabvars(i)%var%Restore2D)) then
1427               Deallocate(Agrif_Gr%tabvars(i)%var%Restore2D)
1428            endif
1429            if (associated(Agrif_Gr%tabvars(i)%var%Restore3D)) then
1430               Deallocate(Agrif_Gr%tabvars(i)%var%Restore3D)
1431            endif
1432            if (associated(Agrif_Gr%tabvars(i)%var%Restore4D)) then
1433               Deallocate(Agrif_Gr%tabvars(i)%var%Restore4D)
1434            endif
1435            if (associated(Agrif_Gr%tabvars(i)%var%Restore5D)) then
1436               Deallocate(Agrif_Gr%tabvars(i)%var%Restore5D)
1437            endif
1438            if (associated(Agrif_Gr%tabvars(i)%var%Restore6D)) then
1439               Deallocate(Agrif_Gr%tabvars(i)%var%Restore6D)
1440            endif
1441C 
1442        endif
1443      enddo
1444C 
1445      Return
1446C
1447C     
1448      End Subroutine Agrif_Free_Restore
1449C
1450C
1451C
1452      End Module Agrif_Clustering
Note: See TracBrowser for help on using the repository browser.