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.F90 in vendors/AGRIF/CMEMS_2020/AGRIF_FILES – NEMO

source: vendors/AGRIF/CMEMS_2020/AGRIF_FILES/modcluster.F90 @ 10087

Last change on this file since 10087 was 10087, checked in by rblod, 6 years ago

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