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

source: vendors/AGRIF/dev_r12970_AGRIF_CMEMS/AGRIF_FILES/modupdate.F90 @ 13027

Last change on this file since 13027 was 13027, checked in by rblod, 4 years ago

New AGRIF library, see ticket #2129

  • Property svn:keywords set to Id
File size: 122.2 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!> Module Agrif_Update
24!>
25!> This module  contains procedures to update a parent grid from its child grids.
26!
27module Agrif_Update
28!
29    use Agrif_UpdateBasic
30    use Agrif_Arrays
31    use Agrif_CurgridFunctions
32    use Agrif_Mask
33#if defined AGRIF_MPI
34    use Agrif_Mpp
35#endif
36!
37    implicit none
38!
39    logical, private :: precomputedone(7) = .FALSE.
40!
41contains
42!
43!===================================================================================================
44!  subroutine Agrif_UpdateVariable
45!
46!> subroutine to set arguments for Agrif_UpdatenD
47!---------------------------------------------------------------------------------------------------
48subroutine Agrif_UpdateVariable ( parent, child, updateinf, updatesup, procname )
49!---------------------------------------------------------------------------------------------------
50    type(Agrif_Variable),  pointer      :: parent       !< Variable on the parent grid
51    type(Agrif_Variable),  pointer      :: child        !< Variable on the child grid
52    integer, dimension(6), intent(in)   :: updateinf    !< First positions where interpolations are calculated
53    integer, dimension(6), intent(in)   :: updatesup    !< Last  positions where interpolations are calculated
54    procedure()                         :: procname     !< Data recovery procedure
55!---------------------------------------------------------------------------------------------------
56    integer, dimension(6) :: nb_child           ! Number of cells on the child grid
57    integer, dimension(6) :: lb_child
58    integer, dimension(6) :: ub_child
59    integer, dimension(6) :: lb_parent
60    real   , dimension(6) ::  s_child           ! Child  grid position (s_root = 0)
61    real   , dimension(6) ::  s_parent          ! Parent grid position (s_root = 0)
62    real   , dimension(6) :: ds_child           ! Child  grid dx (ds_root = 1)
63    real   , dimension(6) :: ds_parent          ! Parent grid dx (ds_root = 1)
64    logical, dimension(6) :: do_update          ! Indicates if we perform update for each dimension
65    integer, dimension(6) :: posvar             ! Position of the variable on the cell (1 or 2)
66    integer, dimension(6) :: oldparentlbound, oldparentubound
67    integer               :: n, nbdim
68    logical               :: wholeupdate
69    type(Agrif_Variable), pointer :: root   ! Variable on the root grid
70!
71    root => child % root_var
72    nbdim = root % nbdim
73!
74    call PreProcessToInterpOrUpdate( parent,   child,       &
75                                     nb_child, ub_child,    &
76                                     lb_child, lb_parent,   &
77                                      s_child,  s_parent,   &
78                                     ds_child, ds_parent, nbdim, interp=.false. )
79!
80    do_update(:) = .true.
81    posvar(1:nbdim) = root % posvar(1:nbdim)
82!
83    do n = 1,nbdim
84!
85        if ( root % interptab(n) == 'N' ) then
86            posvar(n) = 1
87            do_update(n) = .false.
88            oldparentlbound(n) = parent % lb(n)
89            oldparentubound(n) = parent % ub(n)
90            parent % lb(n) = child % lb(n)
91            parent % ub(n) = child % ub(n)
92        end if
93!
94    enddo
95
96    wholeupdate = .FALSE.
97!
98    do n = 1,nbdim
99        if ( do_update(n) ) then
100            if ( (updateinf(n) > updatesup(n)) .OR. &
101                ((updateinf(n) == -99) .AND. (updatesup(n) == -99)) &
102               ) then
103                wholeupdate = .TRUE.
104            endif
105        endif
106    enddo
107!
108    IF (wholeupdate) THEN
109        call Agrif_UpdateWhole(parent, child,               &
110                updateinf(1:nbdim), updatesup(1:nbdim),     &
111                lb_child(1:nbdim), lb_parent(1:nbdim),      &
112                nb_child(1:nbdim), posvar(1:nbdim),         &
113                do_update(1:nbdim),                         &
114                s_child(1:nbdim),   s_parent(1:nbdim),      &
115                ds_child(1:nbdim), ds_parent(1:nbdim), nbdim, procname)
116    ELSE
117        call Agrif_UpdateBcnD(parent, child,                &
118                updateinf(1:nbdim), updatesup(1:nbdim),     &
119                lb_child(1:nbdim), lb_parent(1:nbdim),      &
120                nb_child(1:nbdim), posvar(1:nbdim),         &
121                do_update(1:nbdim),                         &
122                s_child(1:nbdim),   s_parent(1:nbdim),      &
123                ds_child(1:nbdim), ds_parent(1:nbdim), nbdim, procname)
124    ENDIF
125!
126    do n = 1,nbdim
127!
128        if ( root % interptab(n) == 'N' ) then  ! No space DIMENSION
129            parent % lb(n) = oldparentlbound(n)
130            parent % ub(n) = oldparentubound(n)
131        end if
132!
133    enddo
134!---------------------------------------------------------------------------------------------------
135end subroutine Agrif_UpdateVariable
136!===================================================================================================
137!
138!===================================================================================================
139!  subroutine Agrif_UpdateWhole
140!---------------------------------------------------------------------------------------------------
141subroutine Agrif_UpdateWhole ( parent, child, uinf, usup,       &
142                               lb_child, lb_parent,             &
143                               nb_child, posvar,                &
144                               do_update,                       &
145                               s_child,   s_parent,             &
146                               ds_child, ds_parent, nbdim, procname )
147!---------------------------------------------------------------------------------------------------
148#if defined AGRIF_MPI
149    include 'mpif.h'
150#endif
151!
152    type(Agrif_Variable),      pointer    :: parent         !< Variable on the parent grid
153    type(Agrif_Variable),      pointer    :: child          !< Variable on the child grid
154    integer, dimension(nbdim), intent(in) :: uinf           !< First positions where interpolations are calculated
155    integer, dimension(nbdim), intent(in) :: usup           !< Last  positions where interpolations are calculated
156    integer,                   intent(in) :: nbdim          !< Number of dimensions of the grid variable
157    integer, dimension(nbdim), intent(in) :: lb_child       !< Index of the first point inside the domain for the parent grid variable
158    integer, dimension(nbdim), intent(in) :: lb_parent      !< Index of the first point inside the domain for the child grid variable
159    integer, dimension(nbdim), intent(in) :: nb_child       !< Number of cells of the child grid
160    integer, dimension(nbdim), intent(in) :: posvar         !< Position of the variable on the cell (1 or 2)
161    logical, dimension(nbdim), intent(in) :: do_update      !< Indicates if we update for each dimension
162    real,    dimension(nbdim), intent(in) :: s_child        !< Positions of the child grid
163    real,    dimension(nbdim), intent(in) :: s_parent       !< Positions of the parent grid
164    real,    dimension(nbdim), intent(in) :: ds_child       !< Space steps of the child grid
165    real,    dimension(nbdim), intent(in) :: ds_parent      !< Space steps of the parent grid
166    procedure()                           :: procname       !< Data recovery procedure
167!
168    integer, dimension(nbdim)     :: type_update ! Type of update (copy or average)
169    integer, dimension(nbdim,2)   :: lubglob
170    integer, dimension(nbdim,2,2) :: indtab      ! limits of the child grid that will be used in the update scheme
171    integer, dimension(nbdim,2,2) :: indtruetab  ! grid variable where boundary conditions are
172    integer :: coeffraf, i
173    integer :: uinfloc, usuploc
174!
175    type_update = child % root_var % type_update(1:nbdim)
176!
177    do i = 1, nbdim
178!
179        if ( do_update(i) ) then
180!
181            coeffraf = nint(ds_parent(i)/ds_child(i))
182            uinfloc = 0
183            usuploc = nb_child(i)/coeffraf - 1
184
185            IF (posvar(i) == 1) THEN
186                usuploc = usuploc - 1
187            ENDIF
188
189            IF (uinf(i) > usup(i)) THEN
190                uinfloc = uinf(i)
191                usuploc = usuploc - uinf(i)
192            ENDIF
193
194            indtab(i,1,1) = lb_child(i) + (uinfloc + 1) * coeffraf
195            indtab(i,1,2) = lb_child(i) + (usuploc + 1) * coeffraf
196
197            IF ( posvar(i) == 1 ) THEN
198                IF ( type_update(i) == Agrif_Update_Full_Weighting ) THEN
199                    indtab(i,1,1) = indtab(i,1,1) - (coeffraf - 1)
200                    indtab(i,1,2) = indtab(i,1,2) + (coeffraf - 1)
201                ELSE IF ( type_update(i) /= Agrif_Update_Copy ) THEN
202                    indtab(i,1,1) = indtab(i,1,1) - coeffraf / 2
203                    indtab(i,1,2) = indtab(i,1,2) + coeffraf / 2
204                ENDIF
205            ELSE
206                indtab(i,1,1) = indtab(i,1,1) - coeffraf
207                indtab(i,1,2) = indtab(i,1,2) - 1
208    ! at this point, indices are OK for an average
209                IF ( type_update(i) == Agrif_Update_Full_Weighting ) THEN
210                    indtab(i,1,1) = indtab(i,1,1) - coeffraf / 2
211                    indtab(i,1,2) = indtab(i,1,2) + coeffraf / 2
212                ENDIF
213            ENDIF
214!
215        else    ! IF ( .not.do_update(i) ) THEN
216!
217            if ( posvar(i) == 1 ) then
218                indtab(i,1,1) = lb_child(i)
219                indtab(i,1,2) = lb_child(i) + nb_child(i)
220            else
221                indtab(i,1,1) = lb_child(i)
222                indtab(i,1,2) = lb_child(i) + nb_child(i) - 1
223            endif
224!
225        endif
226    enddo
227
228! lubglob contains the global lbound and ubound of the child array
229! lubglob(:,1) : global lbound for each dimension
230! lubglob(:,2) : global lbound for each dimension
231!
232    call Agrif_get_var_global_bounds(child, lubglob, nbdim)
233!
234    indtruetab(1:nbdim,1,1) = max(indtab(1:nbdim,1,1), lubglob(1:nbdim,1))
235    indtruetab(1:nbdim,1,2) = min(indtab(1:nbdim,1,2), lubglob(1:nbdim,2))
236!
237    call Agrif_UpdatenD(type_update, parent, child,                         &
238                        indtruetab(1:nbdim,1,1), indtruetab(1:nbdim,1,2),   &
239                        lb_child(1:nbdim), lb_parent(1:nbdim),              &
240                        s_child(1:nbdim), s_parent(1:nbdim),                &
241                        ds_child(1:nbdim), ds_parent(1:nbdim),              &
242#if defined AGRIF_MPI
243                        posvar, do_update,                                  &
244#endif
245                        nbdim, procname)
246!---------------------------------------------------------------------------------------------------
247end subroutine Agrif_UpdateWhole
248!===================================================================================================
249!
250!===================================================================================================
251!  subroutine Agrif_UpdateBcnd
252!---------------------------------------------------------------------------------------------------
253subroutine Agrif_UpdateBcnd ( parent, child, uinf, usup,        &
254                              lb_child, lb_parent,              &
255                              nb_child, posvar,                 &
256                              do_update,                        &
257                              s_child,  s_parent,               &
258                              ds_child, ds_parent, nbdim, procname )
259!---------------------------------------------------------------------------------------------------
260#if defined AGRIF_MPI
261    include 'mpif.h'
262#endif
263!
264    type(Agrif_Variable), pointer           :: parent       !< Variable on the parent grid
265    type(Agrif_Variable), pointer           :: child        !< Variable on the child grid
266    integer, dimension(nbdim), intent(in)   :: uinf         !< First positions where interpolations are calculated
267    integer, dimension(nbdim), intent(in)   :: usup         !< Last  positions where interpolations are calculated
268    integer                                 :: nbdim        !< Number of dimensions of the grid variable
269    integer, dimension(nbdim), intent(in)   :: lb_child     !< Index of the first point inside the domain for
270                                                            !!   the parent grid variable
271    integer, dimension(nbdim), intent(in)   :: lb_parent    !< Index of the first point inside the domain for
272                                                            !!   the child grid variable
273    integer, dimension(nbdim), intent(in)   :: nb_child     !< Number of cells of the child grid
274    integer, dimension(nbdim), intent(in)   :: posvar       !< Position of the variable on the cell (1 or 2)
275    logical, dimension(nbdim), intent(in)   :: do_update    !< Indicates if we update for each dimension
276    real,    dimension(nbdim), intent(in)   :: s_child      !< Positions of the child grid
277    real,    dimension(nbdim), intent(in)   :: s_parent     !< Positions of the parent grid
278    real,    dimension(nbdim), intent(in)   :: ds_child     !< Space steps of the child grid
279    real,    dimension(nbdim), intent(in)   :: ds_parent    !< Space steps of the parent grid
280    procedure()                             :: procname     !< Data recovery procedure
281    integer, dimension(6)  :: loctab_child      ! Indicates if the child grid has a common border
282                                                !    with the root grid
283    type(Agrif_Variable), pointer :: root_var   ! Variable on the root grid
284!
285    integer,dimension(nbdim)     :: type_update ! Type of update (copy or average)
286    integer,dimension(nbdim,2)   :: lubglob
287    integer                      :: i
288    integer,dimension(nbdim,2,2) :: indtab         ! Arrays indicating the limits of the child
289    integer,dimension(nbdim,2,2) :: indtruetab     ! grid variable where boundary conditions are
290    integer,dimension(nbdim,2,2,nbdim)   :: ptres  ! calculated
291    integer                      :: nb, ndir
292    integer :: coeffraf
293    integer :: n
294!
295    root_var => child % root_var
296    loctab_child(1:nbdim) = 0
297!
298    do n = 1,nbdim
299!
300        select case(root_var % interptab(n))
301!
302        case('x') ! x DIMENSION
303!
304            if (Agrif_Curgrid % NearRootBorder(1))      loctab_child(n) = -1
305            if (Agrif_Curgrid % DistantRootBorder(1))   loctab_child(n) = -2
306            if ((Agrif_Curgrid % NearRootBorder(1)) .AND. &
307                (Agrif_Curgrid % DistantRootBorder(1))) loctab_child(n) = -3
308!
309        case('y') ! y DIMENSION
310!
311            if (Agrif_Curgrid % NearRootBorder(2))      loctab_child(n) = -1
312            if (Agrif_Curgrid % DistantRootBorder(2))   loctab_child(n) = -2
313            if ((Agrif_Curgrid % NearRootBorder(2)) .AND. &
314                (Agrif_Curgrid % DistantRootBorder(2))) loctab_child(n) = -3
315!
316        case('z') ! z DIMENSION
317!
318            if (Agrif_Curgrid % NearRootBorder(3))      loctab_child(n) = -1
319            if (Agrif_Curgrid % DistantRootBorder(3))   loctab_child(n) = -2
320            if ((Agrif_Curgrid % NearRootBorder(3)) .AND. &
321                (Agrif_Curgrid % DistantRootBorder(3))) loctab_child(n) = -3
322!
323        case('N') ! No space DIMENSION
324!
325            loctab_child(n) = -3
326!
327        end select
328!
329    enddo
330
331    type_update = child % root_var % type_update(1:nbdim)
332!
333    DO i = 1, nbdim
334        coeffraf = nint(ds_parent(i)/ds_child(i))
335        indtab(i,1,1) = lb_child(i) + (uinf(i) + 1) * coeffraf
336        indtab(i,1,2) = lb_child(i) + (usup(i) + 1) * coeffraf
337
338        indtab(i,2,1) = lb_child(i) + nb_child(i) - (usup(i)+1) *  coeffraf
339        indtab(i,2,2) = lb_child(i) + nb_child(i) - (uinf(i)+1) *  coeffraf
340
341        IF (posvar(i) == 1) THEN
342            IF (type_update(i) == Agrif_Update_Full_Weighting) THEN
343                indtab(i,:,1) = indtab(i,:,1) - (coeffraf - 1)
344                indtab(i,:,2) = indtab(i,:,2) + (coeffraf - 1)
345            ELSE IF (type_update(i) /= Agrif_Update_Copy) THEN
346                indtab(i,:,1) = indtab(i,:,1) - coeffraf / 2
347                indtab(i,:,2) = indtab(i,:,2) + coeffraf / 2
348            ENDIF
349        ELSE
350            indtab(i,1,1) = indtab(i,1,1) - coeffraf
351            indtab(i,1,2) = indtab(i,1,2) - 1
352            indtab(i,2,2) = indtab(i,2,2) + coeffraf - 1
353            IF (type_update(i) == Agrif_Update_Full_Weighting) THEN
354                indtab(i,1,1) = indtab(i,1,1) - coeffraf/2
355                indtab(i,1,2) = indtab(i,1,2) + coeffraf/2
356                indtab(i,2,1) = indtab(i,2,1) - coeffraf/2
357                indtab(i,2,2) = indtab(i,2,2) + coeffraf/2
358            ENDIF
359        ENDIF
360    ENDDO
361!
362    call Agrif_get_var_global_bounds(child,lubglob,nbdim)
363!
364    indtruetab(1:nbdim,1,1) = max(indtab(1:nbdim,1,1),lubglob(1:nbdim,1))
365    indtruetab(1:nbdim,1,2) = max(indtab(1:nbdim,1,2),lubglob(1:nbdim,1))
366    indtruetab(1:nbdim,2,1) = min(indtab(1:nbdim,2,1),lubglob(1:nbdim,2))
367    indtruetab(1:nbdim,2,2) = min(indtab(1:nbdim,2,2),lubglob(1:nbdim,2))
368!
369    do nb = 1,nbdim
370        if ( do_update(nb) ) then
371            do ndir = 1,2
372              if (loctab_child(nb) /= (-ndir) .AND. loctab_child(nb) /= -3) then
373                ptres(nb,1,ndir,nb) = indtruetab(nb,ndir,1)
374                ptres(nb,2,ndir,nb) = indtruetab(nb,ndir,2)
375                do i = 1,nbdim
376                    if ( i /= nb ) then
377                        if ( do_update(i) ) then
378                            ptres(i,1,ndir,nb) = indtruetab(i,1,1)
379                            ptres(i,2,ndir,nb) = indtruetab(i,2,2)
380                        else
381                            if (posvar(i) == 1) then
382                                ptres(i,1,ndir,nb) = lb_child(i)
383                                ptres(i,2,ndir,nb) = lb_child(i) + nb_child(i)
384                            else
385                                ptres(i,1,ndir,nb) = lb_child(i)
386                                ptres(i,2,ndir,nb) = lb_child(i) + nb_child(i) - 1
387                            endif
388                        endif
389                    endif
390                enddo
391              endif
392            enddo
393        endif
394    enddo
395!
396    do nb = 1,nbdim
397        if ( do_update(nb) ) then
398            do ndir = 1,2
399              if (loctab_child(nb) /= (-ndir) .AND. loctab_child(nb) /= -3) then
400                call Agrif_UpdatenD(type_update, parent, child,             &
401                        ptres(1:nbdim,1,ndir,nb),ptres(1:nbdim,2,ndir,nb),  &
402                        lb_child(1:nbdim),lb_parent(1:nbdim),               &
403                        s_child(1:nbdim),s_parent(1:nbdim),                 &
404                        ds_child(1:nbdim),ds_parent(1:nbdim),               &
405#if defined AGRIF_MPI
406                        posvar,do_update,                                   &
407#endif
408                        nbdim,procname,nb,ndir)
409              endif
410            enddo
411        endif
412    enddo
413!---------------------------------------------------------------------------------------------------
414end subroutine Agrif_UpdateBcnd
415!===================================================================================================
416!
417!===================================================================================================
418!  subroutine Agrif_UpdatenD
419!
420!> updates a 2D grid variable on the parent grid of the current grid
421!---------------------------------------------------------------------------------------------------
422subroutine Agrif_UpdatenD ( type_update, parent, child,     &
423                            pttab, petab,                   &
424                            lb_child, lb_parent,            &
425                            s_child,  s_parent,             &
426                            ds_child, ds_parent,            &
427#if defined AGRIF_MPI
428                            posvar, do_update,              &
429#endif
430                            nbdim, procname, nb, ndir )
431!---------------------------------------------------------------------------------------------------
432#if defined AGRIF_MPI
433    include 'mpif.h'
434#endif
435!
436    type(Agrif_Variable), pointer           :: parent       !< Variable of the parent grid
437    type(Agrif_Variable), pointer           :: child        !< Variable of the child grid
438    integer,                   intent(in)   :: nbdim
439    integer, dimension(nbdim), intent(in)   :: type_update  !< Type of update (copy or average)
440    integer, dimension(nbdim), intent(in)   :: pttab        !< Index of the first point inside the domain
441    integer, dimension(nbdim), intent(in)   :: petab        !< Index of the first point inside the domain
442    integer, dimension(nbdim), intent(in)   :: lb_child  !< Index of the first point inside the domain for the child
443                                                            !!    grid variable
444    integer, dimension(nbdim), intent(in)   :: lb_parent !< Index of the first point inside the domain for the parent
445                                                            !!    grid variable
446    real,    dimension(nbdim), intent(in)   :: s_child      !< Positions of the child grid
447    real,    dimension(nbdim), intent(in)   :: s_parent     !< Positions of the parent grid
448    real,    dimension(nbdim), intent(in)   :: ds_child     !< Space steps of the child grid
449    real,    dimension(nbdim), intent(in)   :: ds_parent    !< Space steps of the parent grid
450    procedure()                             :: procname     !< Data recovery procedure
451    integer, optional,         intent(in)   :: nb, ndir
452!---------------------------------------------------------------------------------------------------
453    integer, dimension(nbdim)       :: pttruetab, cetruetab
454#if defined AGRIF_MPI
455    integer, dimension(nbdim)       :: posvar  !< Position of the variable on the cell (1 or 2)
456    logical, dimension(nbdim)       :: do_update
457#endif
458    integer, dimension(nbdim)       :: coords
459    integer, dimension(nbdim)       :: indmin, indmax
460    integer, dimension(nbdim)       :: indminglob, indmaxglob
461    real   , dimension(nbdim)       :: s_Child_temp, s_Parent_temp
462    integer, dimension(nbdim)       :: lowerbound,upperbound
463    integer, dimension(nbdim)       :: pttruetabwhole, cetruetabwhole
464    integer, dimension(nbdim,2,2)   :: childarray
465    integer, dimension(nbdim,2,2)   :: parentarray
466    integer,dimension(nbdim)     :: type_update_temp
467    logical :: memberin, member
468    integer :: nbin, ndirin
469    integer :: i, j,k,l,m
470    LOGICAL,DIMENSION(:),ALLOCATABLE :: member_chuncks
471    INTEGER,DIMENSION(:,:),ALLOCATABLE :: decal_chunks
472    INTEGER :: agrif_external_switch_index
473    INTEGER, DIMENSION(2) :: test_orientation
474!
475#if defined AGRIF_MPI
476!
477    integer,dimension(nbdim)    :: indminglob2,indmaxglob2
478    INTEGER, DIMENSION(:,:),ALLOCATABLE :: indminglob_chunks, indmaxglob_chunks
479    INTEGER, DIMENSION(:,:),ALLOCATABLE :: indminglob2_chunks,indmaxglob2_chunks
480    INTEGER, DIMENSION(:,:),ALLOCATABLE :: indminglob3_chunks,indmaxglob3_chunks
481    logical, dimension(0:Agrif_Nbprocs-1) :: sendtoproc1,recvfromproc1
482    logical, dimension(0:Agrif_Nbprocs-1) :: sendtoproc2,recvfromproc2
483    integer                               :: code, local_proc
484    integer, dimension(nbdim,4)           :: tab3
485    integer, dimension(nbdim,4,0:Agrif_Nbprocs-1) :: tab4
486    integer, dimension(nbdim,0:Agrif_Nbprocs-1,8) :: tab4t
487    integer, dimension(nbdim,0:Agrif_Nbprocs-1,8) :: tab5t
488    logical :: find_list_update
489    logical, dimension(0:Agrif_Nbprocs-1) :: memberinall, memberinall2
490    logical, dimension(1) :: memberin1
491!
492#endif
493!
494    type(Agrif_Variable), pointer, save :: tempC => NULL()       ! Temporary child grid variable
495    type(Agrif_Variable), pointer, save :: tempP => NULL()       ! Temporary parent grid variable
496    type(Agrif_Variable), pointer, save :: tempCextend => NULL() ! Temporary child
497
498    type(Agrif_Variable), pointer, save :: tempPextend => NULL() ! Temporary parent
499    type(Agrif_Variable), pointer, save :: tempPextend_chunk => NULL() ! Temporary parent
500    type(Agrif_Variable), pointer :: tempP_indic, tempP_average
501    type(Agrif_Variable), pointer :: tempC_indic
502    logical :: compute_average
503    real :: coeff_multi
504    integer :: nb_dimensions
505
506! CHUNK (e.g. periodicity)
507
508    INTEGER :: nb_chunks
509    INTEGER, DIMENSION(:,:,:,:), ALLOCATABLE :: parentarray_chunk
510    INTEGER, DIMENSION(:,:,:,:), ALLOCATABLE :: parentarray_chunk_decal
511    INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: bounds_chunks
512    logical,dimension(:),allocatable :: correction_required
513!
514
515!
516!   Get local lower and upper bound of the child variable
517    call Agrif_get_var_bounds_array(child, lowerbound, upperbound, nbdim)
518
519! here pttab and petab corresponds to the (global) indices of the points needed in the update
520! pttruetab and cetruetab contains only indices that are present on the local processor
521!
522    coords = child % root_var % coords
523!
524
525    call Agrif_Childbounds( nbdim, lowerbound, upperbound, pttab, petab, Agrif_Procrank,    &
526                            coords, pttruetab, cetruetab, memberin )
527
528        if (agrif_debug_update) then
529        print *,'************CHILDBOUNDS*********************************'
530#ifdef AGRIF_MPI
531         print *,'Processeur ',Agrif_Procrank
532#endif
533        print *,'memberin ',memberin
534        do i = 1 , nbdim
535        print *,'Direction ',i,' indices debut: ',pttab(i),pttruetab(i)
536        print *,'Direction ',i,' indices fin  : ',petab(i),cetruetab(i)
537        enddo
538        print *,'*********************************************'
539        endif
540
541    call Agrif_Prtbounds( nbdim, indminglob, indmaxglob, s_Parent_temp, s_Child_temp,       &
542                         s_child, ds_child, s_parent, ds_parent,                            &
543                         pttab, petab, lb_child, lb_parent                                  &
544#if defined AGRIF_MPI
545                       , posvar, type_update, do_update, pttruetabwhole, cetruetabwhole     &
546#endif
547            )
548
549        if (agrif_debug_update) then
550        print *,'************PRTBOUNDS*********************************'
551#ifdef AGRIF_MPI
552         print *,'Processeur ',Agrif_Procrank
553#endif
554        do i = 1 , nbdim
555        print *,'Direction ',i,' indminglob : ',indminglob(i)
556        print *,'Direction ',i,' indmaxglob : ',indmaxglob(i)
557        enddo
558
559        do i = 1 , nbdim
560        print *,'Direction ',i,' s_Parent_temp : ',s_Parent_temp(i)
561        print *,'Direction ',i,' s_Child_temp  : ',s_Child_temp(i)
562        enddo
563        print *,'*********************************************'
564        endif
565
566#if defined AGRIF_MPI
567!
568    IF (memberin) THEN
569        call Agrif_GlobalToLocalBounds(childarray,lowerbound,upperbound,    &
570                                       pttruetab,cetruetab, coords,         &
571                                       nbdim, Agrif_Procrank, member)
572    ENDIF
573
574        if (agrif_debug_update) then
575        print *,'************GlobalToLocalBounds******************'
576#ifdef AGRIF_MPI
577         print *,'Processeur ',Agrif_Procrank
578#endif
579        do i = 1 , nbdim
580        print *,'Direction ',i,' childarray global : ',childarray(i,1,1),childarray(i,2,1)
581        print *,'Direction ',i,' childarray local  : ',childarray(i,1,2),childarray(i,2,2)
582        enddo
583
584        print *,'*********************************************'
585        endif
586
587    call Agrif_Prtbounds(nbdim, indmin, indmax,                     &
588                         s_Parent_temp, s_Child_temp,               &
589                         s_child, ds_child, s_parent, ds_parent,    &
590                         pttruetab, cetruetab, lb_child, lb_parent, &
591                         posvar, type_update, do_update,            &
592                         pttruetabwhole, cetruetabwhole)
593
594        if (agrif_debug_update) then
595        print *,'************PRTBOUNDS II *********************************'
596#ifdef AGRIF_MPI
597         print *,'Processeur ',Agrif_Procrank
598#endif
599        do i = 1 , nbdim
600        print *,'Direction ',i,' indmin : ',indmin(i)
601        print *,'Direction ',i,' indmax : ',indmax(i)
602        enddo
603
604        do i = 1 , nbdim
605        print *,'Direction ',i,' s_Parent_temp : ',s_Parent_temp(i)
606        print *,'Direction ',i,' s_Child_temp  : ',s_Child_temp(i)
607        enddo
608        print *,'*********************************************'
609        endif
610!
611#else
612    indmin = indminglob
613    indmax = indmaxglob
614    pttruetabwhole = pttruetab
615    cetruetabwhole = cetruetab
616    childarray(:,1,2) = pttruetab
617    childarray(:,2,2) = cetruetab
618#endif
619
620    IF (.not.present(nb)) THEN
621        nbin=0
622        ndirin=0
623    ELSE
624        nbin = nb
625        ndirin = ndir
626    ENDIF
627
628    IF (memberin) THEN
629!
630        IF ( .not.associated(tempC) )  allocate(tempC)
631!
632        call Agrif_array_allocate(tempC,pttruetab,cetruetab,nbdim)
633        call Agrif_var_set_array_tozero(tempC,nbdim)
634
635        SELECT CASE (nbdim)
636        CASE(1)
637            CALL procname(tempC%array1,                 &
638                      childarray(1,1,2),childarray(1,2,2),.TRUE.,nbin,ndirin)
639        CASE(2)
640            CALL procname(tempC%array2,                 &
641                      childarray(1,1,2),childarray(1,2,2),  &
642                      childarray(2,1,2),childarray(2,2,2),.TRUE.,nbin,ndirin)
643        CASE(3)
644            CALL procname(tempC%array3,                 &
645                      childarray(1,1,2),childarray(1,2,2),  &
646                      childarray(2,1,2),childarray(2,2,2),  &
647                      childarray(3,1,2),childarray(3,2,2),.TRUE.,nbin,ndirin)
648        CASE(4)
649            CALL procname(tempC%array4,                 &
650                      childarray(1,1,2),childarray(1,2,2),  &
651                      childarray(2,1,2),childarray(2,2,2),  &
652                      childarray(3,1,2),childarray(3,2,2),  &
653                      childarray(4,1,2),childarray(4,2,2),.TRUE.,nbin,ndirin)
654        CASE(5)
655            CALL procname(tempC%array5,                 &
656                      childarray(1,1,2),childarray(1,2,2),  &
657                      childarray(2,1,2),childarray(2,2,2),  &
658                      childarray(3,1,2),childarray(3,2,2),  &
659                      childarray(4,1,2),childarray(4,2,2),  &
660                      childarray(5,1,2),childarray(5,2,2),.TRUE.,nbin,ndirin)
661        CASE(6)
662            CALL procname(tempC%array6,                 &
663                      childarray(1,1,2),childarray(1,2,2),  &
664                      childarray(2,1,2),childarray(2,2,2),  &
665                      childarray(3,1,2),childarray(3,2,2),  &
666                      childarray(4,1,2),childarray(4,2,2),  &
667                      childarray(5,1,2),childarray(5,2,2),  &
668                      childarray(6,1,2),childarray(6,2,2),.TRUE.,nbin,ndirin)
669        END SELECT
670!
671    ENDIF
672!
673#if defined AGRIF_MPI
674!
675!     tab2 contains the necessary limits of the parent grid for each processor
676
677    if (Associated(child%list_update)) then
678        call Agrif_Find_list_update(child%list_update,pttab,petab,                      &
679                                    lb_child,lb_parent,nbdim,                         &
680                                    find_list_update,tab4t,tab5t,memberinall,memberinall2,  &
681                                    sendtoproc1,recvfromproc1,sendtoproc2,recvfromproc2)
682    else
683        find_list_update = .FALSE.
684    endif
685
686    if (.not.find_list_update) then
687        tab3(:,1) = pttruetab(:)
688        tab3(:,2) = cetruetab(:)
689        tab3(:,3) = pttruetabwhole(:)
690        tab3(:,4) = cetruetabwhole(:)
691!
692        call MPI_ALLGATHER(tab3,4*nbdim,MPI_INTEGER,tab4,4*nbdim,MPI_INTEGER,Agrif_mpi_comm,code)
693
694        if ( .not.associated(tempCextend) ) allocate(tempCextend)
695        do k=0,Agrif_Nbprocs-1
696            do j=1,4
697                do i=1,nbdim
698                    tab4t(i,k,j) = tab4(i,j,k)
699                enddo
700            enddo
701        enddo
702
703        memberin1(1) = memberin
704        call MPI_ALLGATHER(memberin1,1,MPI_LOGICAL,memberinall,1,MPI_LOGICAL,Agrif_mpi_comm,code)
705
706        call Get_External_Data_first(tab4t(:,:,1),tab4t(:,:,2),tab4t(:,:,3),tab4t(:,:,4),   &
707                                     nbdim, memberinall, coords,                            &
708                                     sendtoproc1,recvfromproc1,                             &
709                                     tab4t(:,:,5),tab4t(:,:,6),tab4t(:,:,7),tab4t(:,:,8),   &
710                                     tab4t(:,:,1),tab4t(:,:,2))
711    endif
712
713    call ExchangeSameLevel(sendtoproc1,recvfromproc1,nbdim,         &
714            tab4t(:,:,3),tab4t(:,:,4),tab4t(:,:,5),tab4t(:,:,6),    &
715            tab4t(:,:,7),tab4t(:,:,8),memberin,tempC,tempCextend)
716
717#else
718    tempCextend => tempC
719#endif
720!
721!     Update of the parent grid (tempP) from the child grid (tempC)
722!
723    IF (memberin) THEN
724!
725        IF ( .not.associated(tempP) ) allocate(tempP)
726!
727        call Agrif_array_allocate(tempP,indmin,indmax,nbdim)
728!
729        if ( nbdim == 1 ) then
730            tempP % array1 = 0.
731            call Agrif_Update_1D_Recursive( type_update(1),  &
732                                            tempP%array1,   &
733                                            tempCextend%array1, &
734                                            indmin(1), indmax(1),   &
735                                            pttruetabwhole(1), cetruetabwhole(1),   &
736                                            s_Child_temp(1), s_Parent_temp(1),      &
737                                            ds_child(1), ds_parent(1) )
738
739            IF (Agrif_UseSpecialValueInUpdate) THEN
740            allocate(tempC_indic)
741            allocate(tempP_indic)
742            call Agrif_array_allocate(tempC_indic,lbound(tempCextend%array1),ubound(tempCextend%array1),nbdim)
743            call Agrif_array_allocate(tempP_indic,lbound(tempP%array1),ubound(tempP%array1),nbdim)
744
745            compute_average = .FALSE.
746            type_update_temp(1:nbdim) = type_update(1:nbdim)
747            IF (ANY(type_update(1:nbdim) == Agrif_Update_Full_Weighting)) THEN
748              compute_average = .TRUE.
749              allocate(tempP_average)
750              call Agrif_array_allocate(tempP_average,lbound(tempP%array1),ubound(tempP%array1),nbdim)
751              WHERE (type_update(1:nbdim) == Agrif_Update_Full_Weighting)
752                type_update_temp(1:nbdim) = Agrif_Update_Average
753              END WHERE
754              call Agrif_Update_1D_Recursive( type_update_temp(1),   &
755                                            tempP_average%array1,       &
756                                            tempCextend%array1, &
757                                            indmin(1), indmax(1),   &
758                                            pttruetabwhole(1), cetruetabwhole(1),   &
759                                            s_Child_temp(1), s_Parent_temp(1),      &
760                                            ds_child(1), ds_parent(1) )
761              coeff_multi = 1.
762              do nb_dimensions=1,nbdim
763                coeff_multi = coeff_multi * nint(ds_parent(nb_dimensions)/ds_child(nb_dimensions))
764              enddo
765            ENDIF
766
767            WHERE (tempCextend%array1 == Agrif_SpecialValueFineGrid)
768              tempC_indic%array1 = 0.
769            ELSEWHERE
770              tempC_indic%array1 = 1.
771            END WHERE
772
773            Agrif_UseSpecialValueInUpdate = .FALSE.
774            Agrif_Update_Weights = .TRUE.
775
776             call Agrif_Update_1D_Recursive( type_update_temp(1),   &
777                                            tempP_indic%array1,       &
778                                            tempC_indic%array1, &
779                                            indmin(1), indmax(1),   &
780                                            pttruetabwhole(1), cetruetabwhole(1),   &
781                                            s_Child_temp(1), s_Parent_temp(1),      &
782                                            ds_child(1), ds_parent(1) )
783
784           Agrif_UseSpecialValueInUpdate = .TRUE.
785           Agrif_Update_Weights = .FALSE.
786
787           IF (compute_average) THEN
788               WHERE (tempP_indic%array1 == 0.)
789                  tempP%array1 = Agrif_SpecialValueFineGrid
790               ELSEWHERE ((tempP_indic%array1 == coeff_multi).AND.(tempP%array1 /= Agrif_SpecialValueFineGrid))
791                  tempP%array1 = tempP%array1 /tempP_indic%array1
792               ELSEWHERE
793                  tempP%array1 = tempP_average%array1 /tempP_indic%array1
794               END WHERE
795
796           ELSE
797               WHERE (tempP_indic%array1 == 0.)
798                  tempP%array1 = Agrif_SpecialValueFineGrid
799               ELSEWHERE
800                  tempP%array1 = tempP%array1 /tempP_indic%array1
801               END WHERE
802            ENDIF
803
804            deallocate(tempP_indic%array1)
805            deallocate(tempC_indic%array1)
806            deallocate(tempC_indic)
807            deallocate(tempP_indic)
808            IF (compute_average) THEN
809              deallocate(tempP_average%array1)
810              deallocate(tempP_average)
811            ENDIF
812            ENDIF
813
814        endif
815        if ( nbdim == 2 ) then
816            call Agrif_Update_2D_Recursive( type_update(1:2),   &
817                                            tempP%array2,       &
818                                            tempCextend%array2, &
819                                            indmin(1:2), indmax(1:2),   &
820                                            pttruetabwhole(1:2), cetruetabwhole(1:2),   &
821                                            s_Child_temp(1:2), s_Parent_temp(1:2),      &
822                                            ds_child(1:2), ds_parent(1:2) )
823
824            IF (Agrif_UseSpecialValueInUpdate) THEN
825            allocate(tempC_indic)
826            allocate(tempP_indic)
827            call Agrif_array_allocate(tempC_indic,lbound(tempCextend%array2),ubound(tempCextend%array2),nbdim)
828            call Agrif_array_allocate(tempP_indic,lbound(tempP%array2),ubound(tempP%array2),nbdim)
829
830            compute_average = .FALSE.
831            type_update_temp(1:nbdim) = type_update(1:nbdim)
832            IF (ANY(type_update == Agrif_Update_Full_Weighting)) THEN
833              compute_average = .TRUE.
834              allocate(tempP_average)
835              call Agrif_array_allocate(tempP_average,lbound(tempP%array2),ubound(tempP%array2),nbdim)
836              WHERE (type_update(1:nbdim) == Agrif_Update_Full_Weighting)
837                type_update_temp(1:nbdim) = Agrif_Update_Average
838              END WHERE
839              call Agrif_Update_2D_Recursive( type_update_temp(1:2),   &
840                                            tempP_average%array2,       &
841                                            tempCextend%array2, &
842                                            indmin(1:2), indmax(1:2),   &
843                                            pttruetabwhole(1:2), cetruetabwhole(1:2),   &
844                                            s_Child_temp(1:2), s_Parent_temp(1:2),      &
845                                            ds_child(1:2), ds_parent(1:2) )
846              coeff_multi = 1.
847              do nb_dimensions=1,nbdim
848                coeff_multi = coeff_multi * nint(ds_parent(nb_dimensions)/ds_child(nb_dimensions))
849              enddo
850            ENDIF
851
852            WHERE (tempCextend%array2 == Agrif_SpecialValueFineGrid)
853              tempC_indic%array2 = 0.
854            ELSEWHERE
855              tempC_indic%array2 = 1.
856            END WHERE
857
858            Agrif_UseSpecialValueInUpdate = .FALSE.
859            Agrif_Update_Weights = .TRUE.
860
861            call Agrif_Update_2D_Recursive( type_update_temp(1:2),   &
862                                            tempP_indic%array2,       &
863                                            tempC_indic%array2, &
864                                            indmin(1:2), indmax(1:2),   &
865                                            pttruetabwhole(1:2), cetruetabwhole(1:2),   &
866                                            s_Child_temp(1:2), s_Parent_temp(1:2),      &
867                                            ds_child(1:2), ds_parent(1:2) )
868
869           Agrif_UseSpecialValueInUpdate = .TRUE.
870           Agrif_Update_Weights = .FALSE.
871
872           IF (compute_average) THEN
873               WHERE (tempP_indic%array2 == 0.)
874                  tempP%array2 = Agrif_SpecialValueFineGrid
875               ELSEWHERE ((tempP_indic%array2 == coeff_multi).AND.(tempP%array2 /= Agrif_SpecialValueFineGrid))
876                  tempP%array2 = tempP%array2 /tempP_indic%array2
877               ELSEWHERE
878                  tempP%array2 = tempP_average%array2 /tempP_indic%array2
879               END WHERE
880
881           ELSE
882               WHERE (tempP_indic%array2 == 0.)
883                  tempP%array2 = Agrif_SpecialValueFineGrid
884               ELSEWHERE
885                  tempP%array2 = tempP%array2 /tempP_indic%array2
886               END WHERE
887            ENDIF
888
889            deallocate(tempP_indic%array2)
890            deallocate(tempC_indic%array2)
891            deallocate(tempC_indic)
892            deallocate(tempP_indic)
893            IF (compute_average) THEN
894              deallocate(tempP_average%array2)
895              deallocate(tempP_average)
896            ENDIF
897            ENDIF
898
899        endif
900        if ( nbdim == 3 ) then
901            call Agrif_Update_3D_Recursive( type_update(1:3),   &
902                                            tempP%array3,       &
903                                            tempCextend%array3, &
904                                            indmin(1:3), indmax(1:3),   &
905                                            pttruetabwhole(1:3), cetruetabwhole(1:3),   &
906                                            s_Child_temp(1:3), s_Parent_temp(1:3),      &
907                                            ds_child(1:3), ds_parent(1:3) )
908
909            IF (Agrif_UseSpecialValueInUpdate) THEN
910            allocate(tempC_indic)
911            allocate(tempP_indic)
912            call Agrif_array_allocate(tempC_indic,lbound(tempCextend%array3),ubound(tempCextend%array3),nbdim)
913            call Agrif_array_allocate(tempP_indic,lbound(tempP%array3),ubound(tempP%array3),nbdim)
914
915            compute_average = .FALSE.
916            type_update_temp(1:nbdim) = type_update(1:nbdim)
917            IF (ANY(type_update == Agrif_Update_Full_Weighting)) THEN
918              compute_average = .TRUE.
919              allocate(tempP_average)
920              call Agrif_array_allocate(tempP_average,lbound(tempP%array3),ubound(tempP%array3),nbdim)
921              WHERE (type_update(1:nbdim) == Agrif_Update_Full_Weighting)
922                type_update_temp(1:nbdim) = Agrif_Update_Average
923              END WHERE
924              call Agrif_Update_3D_Recursive( type_update_temp(1:3),   &
925                                            tempP_average%array3,       &
926                                            tempCextend%array3, &
927                                            indmin(1:3), indmax(1:3),   &
928                                            pttruetabwhole(1:3), cetruetabwhole(1:3),   &
929                                            s_Child_temp(1:3), s_Parent_temp(1:3),      &
930                                            ds_child(1:3), ds_parent(1:3) )
931              coeff_multi = 1.
932              do nb_dimensions=1,nbdim
933                coeff_multi = coeff_multi * nint(ds_parent(nb_dimensions)/ds_child(nb_dimensions))
934              enddo
935            ENDIF
936
937            WHERE (tempCextend%array3 == Agrif_SpecialValueFineGrid)
938              tempC_indic%array3 = 0.
939            ELSEWHERE
940              tempC_indic%array3 = 1.
941            END WHERE
942
943            Agrif_UseSpecialValueInUpdate = .FALSE.
944            Agrif_Update_Weights = .TRUE.
945
946             call Agrif_Update_3D_Recursive( type_update_temp(1:3),   &
947                                            tempP_indic%array3,       &
948                                            tempC_indic%array3, &
949                                            indmin(1:3), indmax(1:3),   &
950                                            pttruetabwhole(1:3), cetruetabwhole(1:3),   &
951                                            s_Child_temp(1:3), s_Parent_temp(1:3),      &
952                                            ds_child(1:3), ds_parent(1:3) )
953
954           Agrif_UseSpecialValueInUpdate = .TRUE.
955           Agrif_Update_Weights = .FALSE.
956
957           IF (compute_average) THEN
958               WHERE (tempP_indic%array3 == 0.)
959                  tempP%array3 = Agrif_SpecialValueFineGrid
960               ELSEWHERE ((tempP_indic%array3 == coeff_multi).AND.(tempP%array3 /= Agrif_SpecialValueFineGrid))
961                  tempP%array3 = tempP%array3 /tempP_indic%array3
962               ELSEWHERE
963                  tempP%array3 = tempP_average%array3 /tempP_indic%array3
964               END WHERE
965
966           ELSE
967               WHERE (tempP_indic%array3 == 0.)
968                  tempP%array3 = Agrif_SpecialValueFineGrid
969               ELSEWHERE
970                  tempP%array3 = tempP%array3 /tempP_indic%array3
971               END WHERE
972            ENDIF
973
974            deallocate(tempP_indic%array3)
975            deallocate(tempC_indic%array3)
976            deallocate(tempC_indic)
977            deallocate(tempP_indic)
978            IF (compute_average) THEN
979              deallocate(tempP_average%array3)
980              deallocate(tempP_average)
981            ENDIF
982            ENDIF
983
984        endif
985        if ( nbdim == 4 ) then
986            call Agrif_Update_4D_Recursive( type_update(1:4),   &
987                                            tempP%array4,       &
988                                            tempCextend%array4, &
989                                            indmin(1:4), indmax(1:4),   &
990                                            pttruetabwhole(1:4), cetruetabwhole(1:4),   &
991                                            s_Child_temp(1:4), s_Parent_temp(1:4),      &
992                                            ds_child(1:4), ds_parent(1:4) )
993
994            IF (Agrif_UseSpecialValueInUpdate) THEN
995
996            allocate(tempC_indic)
997            allocate(tempP_indic)
998            call Agrif_array_allocate(tempC_indic,lbound(tempCextend%array4),ubound(tempCextend%array4),nbdim)
999            call Agrif_array_allocate(tempP_indic,lbound(tempP%array4),ubound(tempP%array4),nbdim)
1000
1001            compute_average = .FALSE.
1002            type_update_temp(1:nbdim) = type_update(1:nbdim)
1003            IF (ANY(type_update == Agrif_Update_Full_Weighting)) THEN
1004              compute_average = .TRUE.
1005              allocate(tempP_average)
1006              call Agrif_array_allocate(tempP_average,lbound(tempP%array4),ubound(tempP%array4),nbdim)
1007              WHERE (type_update(1:nbdim) == Agrif_Update_Full_Weighting)
1008                type_update_temp(1:nbdim) = Agrif_Update_Average
1009              END WHERE
1010              call Agrif_Update_4D_Recursive( type_update_temp(1:4),   &
1011                                            tempP_average%array4,       &
1012                                            tempCextend%array4, &
1013                                            indmin(1:4), indmax(1:4),   &
1014                                            pttruetabwhole(1:4), cetruetabwhole(1:4),   &
1015                                            s_Child_temp(1:4), s_Parent_temp(1:4),      &
1016                                            ds_child(1:4), ds_parent(1:4) )
1017              coeff_multi = 1.
1018              do nb_dimensions=1,nbdim
1019                coeff_multi = coeff_multi * nint(ds_parent(nb_dimensions)/ds_child(nb_dimensions))
1020              enddo
1021            ENDIF
1022
1023            WHERE (tempCextend%array4 == Agrif_SpecialValueFineGrid)
1024              tempC_indic%array4 = 0.
1025            ELSEWHERE
1026              tempC_indic%array4 = 1.
1027            END WHERE
1028
1029            Agrif_UseSpecialValueInUpdate = .FALSE.
1030            Agrif_Update_Weights = .TRUE.
1031
1032             call Agrif_Update_4D_Recursive( type_update_temp(1:4),   &
1033                                            tempP_indic%array4,       &
1034                                            tempC_indic%array4, &
1035                                            indmin(1:4), indmax(1:4),   &
1036                                            pttruetabwhole(1:4), cetruetabwhole(1:4),   &
1037                                            s_Child_temp(1:4), s_Parent_temp(1:4),      &
1038                                            ds_child(1:4), ds_parent(1:4) )
1039
1040           Agrif_UseSpecialValueInUpdate = .TRUE.
1041           Agrif_Update_Weights = .FALSE.
1042
1043           IF (compute_average) THEN
1044               WHERE (tempP_indic%array4 == 0.)
1045                  tempP%array4 = Agrif_SpecialValueFineGrid
1046               ELSEWHERE ((tempP_indic%array4 == coeff_multi).AND.(tempP%array4 /= Agrif_SpecialValueFineGrid))
1047                  tempP%array4 = tempP%array4 /tempP_indic%array4
1048               ELSEWHERE
1049                  tempP%array4 = tempP_average%array4 /tempP_indic%array4
1050               END WHERE
1051
1052           ELSE
1053               WHERE (tempP_indic%array4 == 0.)
1054                  tempP%array4 = Agrif_SpecialValueFineGrid
1055               ELSEWHERE
1056                  tempP%array4 = tempP%array4 /tempP_indic%array4
1057               END WHERE
1058            ENDIF
1059            deallocate(tempP_indic%array4)
1060            deallocate(tempC_indic%array4)
1061            deallocate(tempC_indic)
1062            deallocate(tempP_indic)
1063            IF (compute_average) THEN
1064              deallocate(tempP_average%array4)
1065              deallocate(tempP_average)
1066            ENDIF
1067            ENDIF
1068
1069        endif
1070        if ( nbdim == 5 ) then
1071            call Agrif_Update_5D_Recursive( type_update(1:5),   &
1072                                            tempP%array5,       &
1073                                            tempCextend%array5, &
1074                                            indmin(1:5), indmax(1:5),   &
1075                                            pttruetabwhole(1:5), cetruetabwhole(1:5),   &
1076                                            s_Child_temp(1:5), s_Parent_temp(1:5),      &
1077                                            ds_child(1:5), ds_parent(1:5) )
1078
1079            IF (Agrif_UseSpecialValueInUpdate) THEN
1080            allocate(tempC_indic)
1081            allocate(tempP_indic)
1082            call Agrif_array_allocate(tempC_indic,lbound(tempCextend%array5),ubound(tempCextend%array5),nbdim)
1083            call Agrif_array_allocate(tempP_indic,lbound(tempP%array5),ubound(tempP%array5),nbdim)
1084
1085            compute_average = .FALSE.
1086            type_update_temp(1:nbdim) = type_update(1:nbdim)
1087            IF (ANY(type_update == Agrif_Update_Full_Weighting)) THEN
1088              compute_average = .TRUE.
1089              allocate(tempP_average)
1090              call Agrif_array_allocate(tempP_average,lbound(tempP%array5),ubound(tempP%array5),nbdim)
1091              WHERE (type_update(1:nbdim) == Agrif_Update_Full_Weighting)
1092                type_update_temp(1:nbdim) = Agrif_Update_Average
1093              END WHERE
1094              call Agrif_Update_5D_Recursive( type_update_temp(1:5),   &
1095                                            tempP_average%array5,       &
1096                                            tempCextend%array5, &
1097                                            indmin(1:5), indmax(1:5),   &
1098                                            pttruetabwhole(1:5), cetruetabwhole(1:5),   &
1099                                            s_Child_temp(1:5), s_Parent_temp(1:5),      &
1100                                            ds_child(1:5), ds_parent(1:5) )
1101              coeff_multi = 1.
1102              do nb_dimensions=1,nbdim
1103                coeff_multi = coeff_multi * nint(ds_parent(nb_dimensions)/ds_child(nb_dimensions))
1104              enddo
1105            ENDIF
1106
1107            WHERE (tempCextend%array5 == Agrif_SpecialValueFineGrid)
1108              tempC_indic%array5 = 0.
1109            ELSEWHERE
1110              tempC_indic%array5 = 1.
1111            END WHERE
1112
1113            Agrif_UseSpecialValueInUpdate = .FALSE.
1114            Agrif_Update_Weights = .TRUE.
1115
1116             call Agrif_Update_5D_Recursive( type_update_temp(1:5),   &
1117                                            tempP_indic%array5,       &
1118                                            tempC_indic%array5, &
1119                                            indmin(1:5), indmax(1:5),   &
1120                                            pttruetabwhole(1:5), cetruetabwhole(1:5),   &
1121                                            s_Child_temp(1:5), s_Parent_temp(1:5),      &
1122                                            ds_child(1:5), ds_parent(1:5) )
1123
1124           Agrif_UseSpecialValueInUpdate = .TRUE.
1125           Agrif_Update_Weights = .FALSE.
1126
1127           IF (compute_average) THEN
1128               WHERE (tempP_indic%array5 == 0.)
1129                  tempP%array5 = Agrif_SpecialValueFineGrid
1130               ELSEWHERE ((tempP_indic%array5 == coeff_multi).AND.(tempP%array5 /= Agrif_SpecialValueFineGrid))
1131                  tempP%array5 = tempP%array5 /tempP_indic%array5
1132               ELSEWHERE
1133                  tempP%array5 = tempP_average%array5 /tempP_indic%array5
1134               END WHERE
1135
1136           ELSE
1137               WHERE (tempP_indic%array5 == 0.)
1138                  tempP%array5 = Agrif_SpecialValueFineGrid
1139               ELSEWHERE
1140                  tempP%array5 = tempP%array5 /tempP_indic%array5
1141               END WHERE
1142            ENDIF
1143
1144            deallocate(tempP_indic%array5)
1145            deallocate(tempC_indic%array5)
1146            deallocate(tempC_indic)
1147            deallocate(tempP_indic)
1148            IF (compute_average) THEN
1149              deallocate(tempP_average%array5)
1150              deallocate(tempP_average)
1151            ENDIF
1152            ENDIF
1153
1154        endif
1155        if ( nbdim == 6 ) then
1156            call Agrif_Update_6D_Recursive( type_update(1:6),   &
1157                                            tempP%array6,       &
1158                                            tempCextend%array6, &
1159                                            indmin(1:6), indmax(1:6),   &
1160                                            pttruetabwhole(1:6), cetruetabwhole(1:6),   &
1161                                            s_Child_temp(1:6), s_Parent_temp(1:6),      &
1162                                            ds_child(1:6), ds_parent(1:6) )
1163            IF (Agrif_UseSpecialValueInUpdate) THEN
1164            allocate(tempC_indic)
1165            allocate(tempP_indic)
1166            call Agrif_array_allocate(tempC_indic,lbound(tempCextend%array6),ubound(tempCextend%array6),nbdim)
1167            call Agrif_array_allocate(tempP_indic,lbound(tempP%array6),ubound(tempP%array6),nbdim)
1168
1169            compute_average = .FALSE.
1170            type_update_temp(1:nbdim) = type_update(1:nbdim)
1171            IF (ANY(type_update == Agrif_Update_Full_Weighting)) THEN
1172              compute_average = .TRUE.
1173              allocate(tempP_average)
1174              call Agrif_array_allocate(tempP_average,lbound(tempP%array6),ubound(tempP%array6),nbdim)
1175              type_update_temp(1:nbdim) = type_update
1176              WHERE (type_update(1:nbdim) == Agrif_Update_Full_Weighting)
1177                type_update_temp(1:nbdim) = Agrif_Update_Average
1178              END WHERE
1179              call Agrif_Update_6D_Recursive( type_update_temp(1:6),   &
1180                                            tempP_average%array6,       &
1181                                            tempCextend%array6, &
1182                                            indmin(1:6), indmax(1:6),   &
1183                                            pttruetabwhole(1:6), cetruetabwhole(1:6),   &
1184                                            s_Child_temp(1:6), s_Parent_temp(1:6),      &
1185                                            ds_child(1:6), ds_parent(1:6) )
1186              coeff_multi = 1.
1187              do nb_dimensions=1,nbdim
1188                coeff_multi = coeff_multi * nint(ds_parent(nb_dimensions)/ds_child(nb_dimensions))
1189              enddo
1190            ENDIF
1191
1192           IF (compute_average) THEN
1193               WHERE (tempP_indic%array6 == 0.)
1194                  tempP%array6 = Agrif_SpecialValueFineGrid
1195               ELSEWHERE ((tempP_indic%array6 == coeff_multi).AND.(tempP%array6 /= Agrif_SpecialValueFineGrid))
1196                  tempP%array6 = tempP%array6 /tempP_indic%array6
1197               ELSEWHERE
1198                  tempP%array6 = tempP_average%array6 /tempP_indic%array6
1199               END WHERE
1200
1201           ELSE
1202               WHERE (tempP_indic%array6 == 0.)
1203                  tempP%array6 = Agrif_SpecialValueFineGrid
1204               ELSEWHERE
1205                  tempP%array6 = tempP%array6 /tempP_indic%array6
1206               END WHERE
1207            ENDIF
1208
1209            Agrif_UseSpecialValueInUpdate = .FALSE.
1210            Agrif_Update_Weights = .TRUE.
1211
1212             call Agrif_Update_6D_Recursive( type_update_temp(1:6),   &
1213                                            tempP_indic%array6,       &
1214                                            tempC_indic%array6, &
1215                                            indmin(1:6), indmax(1:6),   &
1216                                            pttruetabwhole(1:6), cetruetabwhole(1:6),   &
1217                                            s_Child_temp(1:6), s_Parent_temp(1:6),      &
1218                                            ds_child(1:6), ds_parent(1:6) )
1219
1220           Agrif_UseSpecialValueInUpdate = .TRUE.
1221           Agrif_Update_Weights = .FALSE.
1222
1223            WHERE (tempP_indic%array6 == 0.)
1224              tempP%array6 = Agrif_SpecialValueFineGrid
1225            ELSEWHERE
1226              tempP%array6 = tempP%array6 /tempP_indic%array6
1227            END WHERE
1228
1229            deallocate(tempP_indic%array6)
1230            deallocate(tempC_indic%array6)
1231            deallocate(tempC_indic)
1232            deallocate(tempP_indic)
1233            IF (compute_average) THEN
1234              deallocate(tempP_average%array6)
1235              deallocate(tempP_average)
1236            ENDIF
1237            ENDIF
1238        endif
1239!
1240        call Agrif_array_deallocate(tempCextend,nbdim)
1241!
1242    ENDIF
1243
1244    if (agrif_debug_update .and. nbdim==2) then
1245        print *,'MINMAXUPDATE = ',minval(tempP%array2),maxval(tempP%array2)
1246    endif
1247
1248#if defined AGRIF_MPI
1249    local_proc = Agrif_Procrank
1250    call Agrif_get_var_bounds_array(parent,lowerbound,upperbound,nbdim)
1251    call Agrif_ChildGrid_to_ParentGrid()
1252
1253          parentarray(:,1,1) = indminglob
1254          parentarray(:,2,1) = indmaxglob
1255          parentarray(:,1,2) = indminglob
1256          parentarray(:,2,2) = indmaxglob
1257        if (associated(agrif_external_mapping)) then
1258          call agrif_external_mapping(nbdim,child%root_var % posvar(1),child%root_var % posvar(2), &
1259                                      parentarray,parentarray_chunk,correction_required,nb_chunks)
1260          allocate(decal_chunks(nb_chunks,nbdim))
1261          do i=1,nb_chunks
1262            decal_chunks(i,:)=parentarray_chunk(i,:,1,1)-parentarray_chunk(i,:,1,2)
1263          enddo
1264        else
1265            nb_chunks=1
1266            allocate(correction_required(nb_chunks))
1267            correction_required=.FALSE.
1268            allocate(parentarray_chunk(nb_chunks,nbdim,2,2))
1269            parentarray_chunk(1,:,:,:)=parentarray
1270            allocate(decal_chunks(nb_chunks,nbdim))
1271            decal_chunks=0
1272        endif
1273        if (agrif_debug_update) then
1274        print *,'AVANT PARENTCHILDBOUNDS'
1275        print *,'nombre de chunks ',nb_chunks
1276        do i=1,nb_chunks
1277          print *,'CHUNK Number : ',i
1278          do j=1,nbdim
1279           print *,'Direction ',j
1280           print *,'MIN MAX (2) = ',parentarray_chunk(i,j,1,2),parentarray_chunk(i,j,2,2)
1281           print *,'MIN MAX (1) = ',parentarray_chunk(i,j,1,1),parentarray_chunk(i,j,2,1)
1282          enddo
1283        enddo
1284        print *,'APRES PARENTCHILDBOUNDS'
1285        endif
1286
1287        allocate(indminglob_chunks(nb_chunks,nbdim))
1288        allocate(indmaxglob_chunks(nb_chunks,nbdim))
1289        allocate(indminglob2_chunks(nb_chunks,nbdim))
1290        allocate(indmaxglob2_chunks(nb_chunks,nbdim))
1291        allocate(indminglob3_chunks(nb_chunks,nbdim))
1292        allocate(indmaxglob3_chunks(nb_chunks,nbdim))
1293        allocate(member_chuncks(nb_chunks))
1294
1295        do i=1,nb_chunks
1296          indminglob_chunks(i,:) = parentarray_chunk(i,:,1,2)
1297          indmaxglob_chunks(i,:) = parentarray_chunk(i,:,2,2)
1298        enddo
1299
1300        do i=1,nb_chunks
1301              call Agrif_Childbounds(nbdim,lowerbound,upperbound,         &
1302                   indminglob_chunks(i,:),indmaxglob_chunks(i,:), local_proc, coords,   &
1303                   indminglob2_chunks(i,:),indmaxglob2_chunks(i,:),member_chuncks(i))
1304        enddo
1305
1306        if (agrif_debug_update) then
1307        print *,'************CHILDBOUNDSPARENTMPI*********************************'
1308#ifdef AGRIF_MPI
1309         print *,'Processeur ',Agrif_Procrank
1310#endif
1311        do j=1,nb_chunks
1312        print *,'Chunk number ',j
1313
1314        do i = 1 , nbdim
1315        print *,'Direction ',i,' indices debut: ',indminglob_chunks(j,i),indminglob2_chunks(j,i)
1316        print *,'Direction ',i,' indices fin  : ',indmaxglob_chunks(j,i),indmaxglob2_chunks(j,i)
1317        enddo
1318        enddo
1319        print *,'*********************************************'
1320        endif
1321
1322        allocate(parentarray_chunk_decal(nb_chunks,nbdim,2,2))
1323        do j=1,nb_chunks
1324        if (agrif_debug_update) print *,'CHUNK = ',j
1325        if (member_chuncks(j)) then
1326            ! call Agrif_GlobalToLocalBounds(parentarray_chunk(j,:,:,:),                               &
1327                                           ! lowerbound,  upperbound,                                  &
1328                                           ! indminglob2_chunks(j,:), indmaxglob2_chunks(j,:), coords, &
1329                                           ! nbdim, local_proc, member_chuncks(j),check_perio=.TRUE.)
1330
1331            call Agrif_GlobalToLocalBounds(parentarray_chunk(j,:,:,:),                               &
1332                                           lowerbound,  upperbound,                                  &
1333                                           indminglob2_chunks(j,:), indmaxglob2_chunks(j,:), coords, &
1334                                           nbdim, local_proc, member_chuncks(j))
1335
1336            if (correction_required(j)) then
1337                do i=1,2
1338                    test_orientation(1)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
1339                    parentarray_chunk(j,i,1,1),i)
1340                    test_orientation(2)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
1341                    parentarray_chunk(j,i,2,1),i)
1342                    parentarray_chunk_decal(j,i,1,1)=minval(test_orientation)
1343                    parentarray_chunk_decal(j,i,2,1)=maxval(test_orientation)
1344                enddo
1345                do i=3,nbdim
1346                  parentarray_chunk_decal(j,i,:,1) = parentarray_chunk(j,i,:,1)+decal_chunks(j,i)
1347                enddo
1348            else
1349            do i=1,nbdim
1350            parentarray_chunk_decal(j,i,:,1) = parentarray_chunk(j,i,:,1)+decal_chunks(j,i)
1351            enddo
1352            endif
1353
1354            if (agrif_debug_update) then
1355            do i=1,nbdim
1356            print *,'parentarray = ',i,parentarray_chunk(j,i,1,1),parentarray_chunk(j,i,2,1), &
1357                parentarray_chunk(j,i,1,2),parentarray_chunk(j,i,2,2)
1358                print *,'parentarraydecal = ',i,parentarray_chunk_decal(j,i,1,1),parentarray_chunk_decal(j,i,2,1)
1359            enddo
1360            endif
1361        endif
1362        enddo
1363
1364!     call Agrif_Childbounds(nbdim, lowerbound, upperbound,                   &
1365!                            indminglob,  indmaxglob,  local_proc, coords,    &
1366!                            indminglob2, indmaxglob2, member)
1367
1368!         if (agrif_debug_update) then
1369!         print *,'************CHILDBOUNDS PARENT*********************************'
1370! #ifdef AGRIF_MPI
1371!          print *,'Processeur ',Agrif_Procrank
1372! #endif
1373!         print *,'member ',member
1374!         do i = 1 , nbdim
1375!         print *,'Direction ',i,' indminglob2 : ',indminglob2(i)
1376!         print *,'Direction ',i,' indmaxglob2 : ',indmaxglob2(i)
1377!         enddo
1378!         print *,'*********************************************'
1379!         endif
1380!
1381!     IF (member) THEN
1382!         call Agrif_GlobalToLocalBounds(parentarray, lowerbound, upperbound, &
1383!                                        indminglob2, indmaxglob2, coords,    &
1384!                                        nbdim, local_proc, member,check_perio=.TRUE.)
1385!     ENDIF
1386
1387!         if (agrif_debug_update) then
1388!         print *,'************GlobalToLocalBounds II******************'
1389! #ifdef AGRIF_MPI
1390!          print *,'Processeur ',Agrif_Procrank
1391! #endif
1392!         do i = 1 , nbdim
1393!         print *,'Direction ',i,' parentarray global : ',parentarray(i,1,1),parentarray(i,2,1)
1394!         print *,'Direction ',i,' parentarray local  : ',parentarray(i,1,2),parentarray(i,2,2)
1395!         enddo
1396
1397!         print *,'*********************************************'
1398!         endif
1399
1400
1401        parentarray(:,1,:)=Huge(1)
1402        parentarray(:,2,:)=-Huge(1)
1403        indminglob2=Huge(1)
1404        indmaxglob2=-Huge(1)
1405        member = .FALSE.
1406        do j=1,nb_chunks
1407          if (member_chuncks(j)) then
1408            do i=1,nbdim
1409             parentarray(i,1,1) = min(parentarray(i,1,1),parentarray_chunk_decal(j,i,1,1))
1410             parentarray(i,1,2) = min(parentarray(i,1,2),parentarray_chunk(j,i,1,2))
1411             parentarray(i,2,1) = max(parentarray(i,2,1),parentarray_chunk_decal(j,i,2,1))
1412             parentarray(i,2,2) = max(parentarray(i,2,2),parentarray_chunk(j,i,2,2))
1413            enddo
1414
1415            if (correction_required(j)) then
1416                do i=1,2
1417                    test_orientation(1)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
1418                    indminglob2_chunks(j,i),i)
1419                    test_orientation(2)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
1420                    indmaxglob2_chunks(j,i),i)
1421                    indminglob2(i)=min(indminglob2(i),minval(test_orientation))
1422                    indmaxglob2(i)=max(indmaxglob2(i),maxval(test_orientation))
1423                enddo
1424
1425                do i=3,nbdim
1426             indminglob2(i)=min(indminglob2(i),indminglob2_chunks(j,i)+decal_chunks(j,i))
1427             indmaxglob2(i)=max(indmaxglob2(i),indmaxglob2_chunks(j,i)+decal_chunks(j,i))
1428                enddo
1429            else
1430                           do i=1,nbdim
1431             indminglob2(i)=min(indminglob2(i),indminglob2_chunks(j,i)+decal_chunks(j,i))
1432             indmaxglob2(i)=max(indmaxglob2(i),indmaxglob2_chunks(j,i)+decal_chunks(j,i))
1433                         enddo
1434         endif
1435
1436            member = .TRUE.
1437          endif
1438        enddo
1439
1440        call Agrif_ParentGrid_to_ChildGrid()
1441
1442        if (agrif_debug_update) then
1443        print *,'************ FINAL PARENTARRAY *****************'
1444#ifdef AGRIF_MPI
1445        print *,'Processeur ',Agrif_Procrank,' MEMBER = ',member
1446        do i=1,nbdim
1447         print *,'Direction ',i,' indices debut = ',parentarray(i,1,1),parentarray(i,1,2)
1448         print *,'Direction ',i,' indices fin = ',parentarray(i,2,1),parentarray(i,2,2)
1449        enddo
1450#endif
1451        endif
1452
1453        if (agrif_debug_update) then
1454        print *,'************ FINAL INDMINGLOB *****************'
1455#ifdef AGRIF_MPI
1456        print *,'Processeur ',Agrif_Procrank,' MEMBER = ',member
1457        do i=1,nbdim
1458         print *,'Direction ',i,' indices debut = ',indminglob2(i)
1459         print *,'Direction ',i,' indices fin = ',indmaxglob2(i)
1460        enddo
1461#endif
1462        endif
1463
1464    if (.not.find_list_update) then
1465        tab3(:,1) = indmin(:)
1466        tab3(:,2) = indmax(:)
1467        tab3(:,3) = indminglob2(:)
1468        tab3(:,4) = indmaxglob2(:)
1469!
1470        call MPI_ALLGATHER(tab3,4*nbdim,MPI_INTEGER,tab4,4*nbdim,MPI_INTEGER,Agrif_mpi_comm,code)
1471
1472        IF ( .not.associated(tempPextend) ) allocate(tempPextend)
1473        DO k=0,Agrif_Nbprocs-1
1474            do j=1,4
1475                do i=1,nbdim
1476                    tab5t(i,k,j) = tab4(i,j,k)
1477                enddo
1478            enddo
1479        enddo
1480
1481        memberin1(1) = member
1482        call MPI_ALLGATHER(memberin1,1,MPI_LOGICAL,memberinall2,1,MPI_LOGICAL,Agrif_mpi_comm,code)
1483        call Get_External_Data_first(tab5t(:,:,1),tab5t(:,:,2),tab5t(:,:,3),tab5t(:,:,4),   &
1484                                     nbdim, memberinall2, coords,                           &
1485                                     sendtoproc2, recvfromproc2,                            &
1486                                     tab5t(:,:,5),tab5t(:,:,6),tab5t(:,:,7),tab5t(:,:,8),   &
1487                                     tab5t(:,:,1),tab5t(:,:,2))
1488
1489        call Agrif_Addto_list_update(child%list_update,pttab,petab,lb_child,lb_parent,      &
1490                                     nbdim,tab4t,tab5t,memberinall,memberinall2,            &
1491                                     sendtoproc1,recvfromproc1,sendtoproc2,recvfromproc2)
1492
1493    endif
1494
1495    call ExchangeSameLevel(sendtoproc2,recvfromproc2,nbdim,                     &
1496                            tab5t(:,:,3),tab5t(:,:,4),tab5t(:,:,5),tab5t(:,:,6),&
1497                            tab5t(:,:,7),tab5t(:,:,8),member,tempP,tempPextend)
1498#else
1499    tempPextend => tempP
1500    parentarray(:,1,1) = indmin
1501    parentarray(:,2,1) = indmax
1502    parentarray(:,1,2) = indmin
1503    parentarray(:,2,2) = indmax
1504
1505        if (associated(agrif_external_mapping)) then
1506          call Agrif_ChildGrid_to_ParentGrid()
1507          call agrif_external_mapping(nbdim,child%root_var % posvar(1),child%root_var % posvar(2), &
1508                                      parentarray,parentarray_chunk,correction_required,nb_chunks)
1509          call Agrif_ParentGrid_to_ChildGrid()
1510          allocate(decal_chunks(nb_chunks,nbdim))
1511          do i=1,nb_chunks
1512            decal_chunks(i,:)=parentarray_chunk(i,:,1,1)-parentarray_chunk(i,:,1,2)
1513          enddo
1514        else
1515            nb_chunks=1
1516            allocate(correction_required(nb_chunks))
1517            correction_required=.FALSE.
1518            allocate(parentarray_chunk(nb_chunks,nbdim,2,2))
1519            parentarray_chunk(1,:,:,:)=parentarray
1520        endif
1521        if (agrif_debug_update) then
1522        print *,'AVANT PARENTCHILDBOUNDS'
1523        print *,'nombre de chunks ',nb_chunks
1524        do i=1,nb_chunks
1525          print *,'CHUNK Number : ',i
1526          do j=1,nbdim
1527           print *,'Direction ',j
1528           print *,'MIN MAX (2) = ',parentarray_chunk(i,j,1,2),parentarray_chunk(i,j,2,2)
1529           print *,'MIN MAX (1) = ',parentarray_chunk(i,j,1,1),parentarray_chunk(i,j,2,1)
1530          enddo
1531        enddo
1532        print *,'APRES PARENTCHILDBOUNDS'
1533        endif
1534        allocate(member_chuncks(nb_chunks))
1535        allocate(parentarray_chunk_decal(nb_chunks,nbdim,2,2))
1536        member_chuncks = .TRUE.
1537        member = .TRUE.
1538        do j=1,nb_chunks
1539        if (agrif_debug_update) print *,'CHUNK = ',j
1540        if (member_chuncks(j)) then
1541            do i=1,nbdim
1542            parentarray_chunk_decal(j,i,:,1) = parentarray_chunk(j,i,:,1)   !+decal_chunks(j,i)
1543            enddo
1544            if (agrif_debug_update) then
1545            do i=1,nbdim
1546            print *,'parentarray = ',i,parentarray_chunk(j,i,1,1),parentarray_chunk(j,i,2,1), &
1547                parentarray_chunk(j,i,1,2),parentarray_chunk(j,i,2,2)
1548                print *,'parentarraydecal = ',i,parentarray_chunk_decal(j,i,1,1),parentarray_chunk_decal(j,i,2,1)
1549            enddo
1550            endif
1551        endif
1552        enddo
1553
1554    member = .TRUE.
1555#endif
1556
1557        if (agrif_debug_update .and. nbdim==2) then
1558        print *,'MINMAXUPDATEEXTND = ',minval(tempPextend%array2),maxval(tempPextend%array2)
1559    endif
1560!
1561!   Special values on the child grid
1562    if ( Agrif_UseSpecialValueFineGrid ) then
1563!
1564!cc         noraftab(1:nbdim) =
1565!cc     &    child % root_var % interptab(1:nbdim) == 'N'
1566!
1567#if defined AGRIF_MPI
1568!
1569!          allocate(childvalues% var)
1570!
1571!          Call Agrif_array_allocate(childvalues%var,
1572!     &                      pttruetab,cetruetab,nbdim)
1573!          Call Agrif_var_full_copy_array(childvalues% var,
1574!     &                                tempC,
1575!     &                                nbdim)
1576!          Call Agrif_CheckMasknD(tempC,childvalues,
1577!     &                           pttruetab(1:nbdim),cetruetab(1:nbdim),
1578!     &                           pttruetab(1:nbdim),cetruetab(1:nbdim),
1579!     &                           noraftab(1:nbdim),nbdim)
1580!          Call Agrif_array_deallocate(childvalues% var,nbdim)
1581!         Deallocate(childvalues % var)
1582!
1583#else
1584!
1585!          Call Agrif_get_var_bounds_array(child,
1586!     &                              lowerbound,upperbound,nbdim)
1587!          Call Agrif_CheckMasknD(tempC,child,
1588!     &                           pttruetab(1:nbdim),cetruetab(1:nbdim),
1589!     &                           lowerbound,
1590!     &                           upperbound,
1591!     &                           noraftab(1:nbdim),nbdim)
1592!
1593#endif
1594!
1595    endif
1596!
1597!   Special values on the parent grid
1598    if (Agrif_UseSpecialValue) then
1599!
1600#if defined AGRIF_MPI
1601!
1602!          Call GiveAgrif_SpecialValueToTab_mpi(parent,tempP,
1603!     &                 parentarray,
1604!     &                 Agrif_SpecialValue,nbdim)
1605!
1606!
1607#else
1608!
1609!          Call GiveAgrif_SpecialValueToTab(parent,tempP,
1610!     &                  indmin,indmax,
1611!     &                  Agrif_SpecialValue,nbdim)
1612!
1613#endif
1614!
1615    endif
1616!
1617
1618
1619        call Agrif_ChildGrid_to_ParentGrid()
1620
1621        if (nb_chunks > 1) then
1622            allocate(tempPextend_chunk)
1623            SELECT CASE(nbdim)
1624            CASE(1)
1625                call Agrif_array_allocate(tempPextend_chunk,lbound(tempPextend%array1),ubound(tempPextend%array1),nbdim)
1626            CASE(2)
1627                call Agrif_array_allocate(tempPextend_chunk,lbound(tempPextend%array2),ubound(tempPextend%array2),nbdim)
1628            CASE(3)
1629                call Agrif_array_allocate(tempPextend_chunk,lbound(tempPextend%array3),ubound(tempPextend%array3),nbdim)
1630            CASE(4)
1631                call Agrif_array_allocate(tempPextend_chunk,lbound(tempPextend%array4),ubound(tempPextend%array4),nbdim)
1632            CASE(5)
1633                call Agrif_array_allocate(tempPextend_chunk,lbound(tempPextend%array5),ubound(tempPextend%array5),nbdim)
1634            CASE(6)
1635                call Agrif_array_allocate(tempPextend_chunk,lbound(tempPextend%array6),ubound(tempPextend%array6),nbdim)
1636            END SELECT
1637        else
1638            tempPextend_chunk => tempPextend
1639        endif
1640
1641        do i=1,nb_chunks
1642
1643        if (member_chuncks(i)) then
1644
1645            if (nb_chunks > 1) then
1646                call Agrif_var_copy_array (tempPextend_chunk, parentarray_chunk_decal(i,:,1,1), parentarray_chunk_decal(i,:,2,1), &
1647                                           tempPextend, parentarray_chunk_decal(i,:,1,1), parentarray_chunk_decal(i,:,2,1), nbdim )
1648            endif
1649!
1650        SELECT CASE(nbdim)
1651        CASE(1)
1652            ! call procname( tempPextend % array1(            &
1653            !         parentarray(1,1,1):parentarray(1,2,1)), &
1654            !         parentarray(1,1,2),parentarray(1,2,2),.FALSE.,nbin,ndirin)
1655
1656            call procname(tempPextend_chunk%array1(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1)),         &
1657                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),.FALSE.,nbin,ndirin)
1658
1659        CASE(2)
1660
1661            if (correction_required(i)) then
1662              call correct_field(tempPextend_chunk%array2(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1663                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1)), &
1664                      parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1), &
1665                      parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1))
1666            endif
1667
1668            call procname(tempPextend_chunk%array2(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1669                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1)),         &
1670                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
1671                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),.FALSE.,nbin,ndirin)
1672
1673        CASE(3)
1674
1675            if (correction_required(i)) then
1676                do k=parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2)
1677              call correct_field(tempPextend_chunk%array3(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1678                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1),k), &
1679                      parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1), &
1680                      parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1))
1681                enddo
1682            endif
1683
1684            call procname(tempPextend_chunk%array3(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1685                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1), &
1686                      parentarray_chunk_decal(i,3,1,1):parentarray_chunk_decal(i,3,2,1)),         &
1687                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
1688                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
1689                      parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),.FALSE.,nbin,ndirin)
1690
1691
1692
1693
1694        CASE(4)
1695
1696            if (correction_required(i)) then
1697                do l=parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2)
1698                do k=parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2)
1699              call correct_field(tempPextend_chunk%array4(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1700                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1),k,l), &
1701                      parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1), &
1702                      parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1))
1703                enddo
1704            enddo
1705            endif
1706
1707            call procname(tempPextend_chunk%array4(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1708                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1), &
1709                      parentarray_chunk_decal(i,3,1,1):parentarray_chunk_decal(i,3,2,1), &
1710                      parentarray_chunk_decal(i,4,1,1):parentarray_chunk_decal(i,4,2,1)),         &
1711                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
1712                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
1713                      parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),    &
1714                      parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2),.FALSE.,nbin,ndirin)
1715
1716        CASE(5)
1717
1718            if (correction_required(i)) then
1719                do m=parentarray_chunk(i,5,1,2),parentarray_chunk(i,5,2,2)
1720                do l=parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2)
1721                do k=parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2)
1722              call correct_field(tempPextend_chunk%array5(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1723                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1),k,l,m), &
1724                      parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1), &
1725                      parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1))
1726                enddo
1727            enddo
1728            enddo
1729            endif
1730
1731            call procname(tempPextend_chunk%array5(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1732                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1), &
1733                      parentarray_chunk_decal(i,3,1,1):parentarray_chunk_decal(i,3,2,1), &
1734                      parentarray_chunk_decal(i,4,1,1):parentarray_chunk_decal(i,4,2,1), &
1735                      parentarray_chunk_decal(i,5,1,1):parentarray_chunk_decal(i,5,2,1)),         &
1736                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
1737                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
1738                      parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),    &
1739                      parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2),    &
1740                      parentarray_chunk(i,5,1,2),parentarray_chunk(i,5,2,2),.FALSE.,nbin,ndirin)
1741
1742        CASE(6)
1743
1744            call procname(tempPextend_chunk%array6(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
1745                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1), &
1746                      parentarray_chunk_decal(i,3,1,1):parentarray_chunk_decal(i,3,2,1), &
1747                      parentarray_chunk_decal(i,4,1,1):parentarray_chunk_decal(i,4,2,1), &
1748                      parentarray_chunk_decal(i,5,1,1):parentarray_chunk_decal(i,5,2,1), &
1749                      parentarray_chunk_decal(i,6,1,1):parentarray_chunk_decal(i,6,2,1)),         &
1750                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
1751                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
1752                      parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),    &
1753                      parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2),    &
1754                      parentarray_chunk(i,5,1,2),parentarray_chunk(i,5,2,2),&
1755                      parentarray_chunk(i,6,1,2),parentarray_chunk(i,6,2,2),.FALSE.,nbin,ndirin)
1756
1757        END SELECT
1758        ENDIF
1759        enddo ! enddo i=1,nb_chunks
1760
1761        if (nb_chunks > 1) then
1762            call Agrif_array_deallocate(tempPextend_chunk,nbdim)
1763            deallocate(tempPextend_chunk)
1764        endif
1765!
1766        call Agrif_ParentGrid_to_ChildGrid()
1767!
1768        if (ANY(member_chuncks)) call Agrif_array_deallocate(tempPextend,nbdim)
1769!
1770
1771!
1772#if defined AGRIF_MPI
1773    IF (memberin) THEN
1774        call Agrif_array_deallocate(tempP,nbdim)
1775        call Agrif_array_deallocate(tempC,nbdim)
1776    ENDIF
1777#endif
1778!---------------------------------------------------------------------------------------------------
1779end subroutine Agrif_UpdatenD
1780!===================================================================================================
1781!
1782!===================================================================================================
1783!  subroutine Agrif_Prtbounds
1784!
1785!> calculates the bounds of the parent grid to be updated by the child grid
1786!---------------------------------------------------------------------------------------------------
1787subroutine Agrif_Prtbounds ( nbdim, indmin, indmax, s_Parent_temp, s_Child_temp,        &
1788                             s_child, ds_child, s_parent, ds_parent,                    &
1789                             pttruetab, cetruetab, lb_child, lb_parent            &
1790#if defined AGRIF_MPI
1791                            ,posvar, type_update, do_update,                &
1792                             pttruetabwhole, cetruetabwhole                             &
1793#endif
1794                )
1795!---------------------------------------------------------------------------------------------------
1796    integer,                   intent(in)   :: nbdim
1797    integer, dimension(nbdim), intent(out)  :: indmin, indmax
1798    real,    dimension(nbdim), intent(out)  :: s_Parent_temp, s_Child_temp
1799    real,    dimension(nbdim), intent(in)   :: s_child,  ds_child
1800    real,    dimension(nbdim), intent(in)   :: s_parent, ds_parent
1801    integer, dimension(nbdim), intent(in)   :: pttruetab, cetruetab
1802    integer, dimension(nbdim), intent(in)   :: lb_child, lb_parent
1803#if defined AGRIF_MPI
1804    integer, dimension(nbdim), intent(in)   :: posvar   !< Position of the variable on the cell (1 or 2)
1805    integer, dimension(nbdim), intent(in)   :: type_update
1806    logical, dimension(nbdim), intent(in)   :: do_update
1807    integer,dimension(nbdim), intent(out)   :: pttruetabwhole, cetruetabwhole
1808#endif
1809!
1810    real,dimension(nbdim) :: dim_newmin,dim_newmax
1811    integer :: i
1812#if defined AGRIF_MPI
1813    real    :: positionmin, positionmax
1814    integer :: imin, imax
1815    integer :: coeffraf
1816#endif
1817!
1818    do i = 1,nbdim
1819!
1820        dim_newmin(i) = s_child(i) + (pttruetab(i) - lb_child(i)) * ds_child(i)
1821        dim_newmax(i) = s_child(i) + (cetruetab(i) - lb_child(i)) * ds_child(i)
1822!
1823        indmin(i) = lb_parent(i) + agrif_ceiling((dim_newmin(i)-s_parent(i))/ds_parent(i))
1824        indmax(i) = lb_parent(i) + agrif_int(    (dim_newmax(i)-s_parent(i))/ds_parent(i))
1825!
1826#if defined AGRIF_MPI
1827        positionmin = s_parent(i) + (indmin(i)-lb_parent(i))*ds_parent(i)
1828        IF ( do_update(i) ) THEN
1829            IF (posvar(i) == 1) THEN
1830                IF ((type_update(i) == Agrif_Update_Average).OR.(type_update(i) == Agrif_Update_Max)) THEN
1831                    positionmin = positionmin - ds_parent(i)/2.
1832                ELSE IF (type_update(i) == Agrif_Update_Full_Weighting) THEN
1833                    positionmin = positionmin - (ds_parent(i)-ds_child(i))
1834                ENDIF
1835            ELSE
1836                IF (type_update(i) /= Agrif_Update_Full_Weighting) THEN
1837                    positionmin = positionmin - ds_parent(i)/2.
1838                ELSE
1839                    coeffraf = nint(ds_parent(i)/ds_child(i))
1840                    if (mod(coeffraf,2) == 1) then
1841                        positionmin = positionmin - (ds_parent(i)-ds_child(i))
1842                    else
1843                        positionmin = positionmin - (ds_parent(i)-ds_child(i))-ds_child(i)/2.
1844                    endif
1845                ENDIF
1846            ENDIF
1847        ENDIF
1848!
1849        imin = lb_child(i) + agrif_ceiling((positionmin-s_child(i))/ds_child(i))
1850        positionmin = s_child(i)  + (imin - lb_child(i)) * ds_child(i)
1851        positionmax = s_parent(i) + (indmax(i)-lb_parent(i))*ds_parent(i)
1852        pttruetabwhole(i) = imin
1853
1854        IF ( do_update(i) ) THEN
1855            IF (posvar(i) == 1) THEN
1856                IF ((type_update(i) == Agrif_Update_Average).OR.(type_update(i) == Agrif_Update_Max)) THEN
1857                    positionmax = positionmax  + ds_parent(i)/2.
1858                ELSE IF (type_update(i) == Agrif_Update_Full_Weighting) THEN
1859                    positionmax = positionmax  + (ds_parent(i)-ds_child(i))
1860                ENDIF
1861            ELSE
1862                IF (type_update(i) /= Agrif_Update_Full_Weighting) THEN
1863                    positionmax = positionmax  + ds_parent(i)/2.
1864                ELSE
1865                    coeffraf = nint(ds_parent(i)/ds_child(i))
1866                    if (mod(coeffraf,2) == 1) then
1867                        positionmax = positionmax + (ds_parent(i)-ds_child(i))
1868                    else
1869                        positionmax = positionmax + (ds_parent(i)-ds_child(i)) + ds_child(i)/2.
1870                    endif
1871                ENDIF
1872            ENDIF
1873        ENDIF
1874
1875        imax = lb_child(i) +agrif_int((positionmax-s_child(i))/ds_child(i))
1876        positionmax = s_child(i) + (imax - lb_child(i)) * ds_child(i)
1877        cetruetabwhole(i) = imax
1878#endif
1879!
1880        s_Parent_temp(i) = s_parent(i) + (indmin(i) - lb_parent(i)) * ds_parent(i)
1881        s_Child_temp(i)  = dim_newmin(i)
1882
1883#if defined AGRIF_MPI
1884        s_Child_temp(i) = positionmin
1885#endif
1886!
1887    enddo
1888!---------------------------------------------------------------------------------------------------
1889end subroutine Agrif_Prtbounds
1890!===================================================================================================
1891!
1892!===================================================================================================
1893!  subroutine Agrif_Update_1D_Recursive
1894!
1895!> Updates a 1D grid variable on the parent grid
1896!---------------------------------------------------------------------------------------------------
1897subroutine Agrif_Update_1D_Recursive ( type_update,         &
1898                                       tempP, tempC,        &
1899                                       indmin, indmax,      &
1900                                       lb_child, ub_child,  &
1901                                       s_child,  s_parent,  &
1902                                       ds_child, ds_parent )
1903!---------------------------------------------------------------------------------------------------
1904    integer,                            intent(in)  :: type_update            !< Type of update (copy or average)
1905    integer,                            intent(in)  :: indmin, indmax
1906    integer,                            intent(in)  :: lb_child, ub_child
1907    real,                               intent(in)  ::  s_child,  s_parent
1908    real,                               intent(in)  :: ds_child, ds_parent
1909    real, dimension(indmin:indmax),     intent(out) :: tempP
1910    real, dimension(lb_child:ub_child), intent(in)  :: tempC
1911!---------------------------------------------------------------------------------------------------
1912    call Agrif_UpdateBase(type_update,              &
1913                          tempP(indmin:indmax),     &
1914                          tempC(lb_child:ub_child), &
1915                          indmin, indmax,           &
1916                          lb_child, ub_child,       &
1917                          s_parent, s_child,        &
1918                          ds_parent, ds_child)
1919!---------------------------------------------------------------------------------------------------
1920end subroutine Agrif_Update_1D_Recursive
1921!===================================================================================================
1922!
1923!===================================================================================================
1924!  subroutine Agrif_Update_2D_Recursive
1925!
1926!> updates a 2D grid variable on the parent grid.
1927!! Calls #Agrif_Update_1D_Recursive and #Agrif_UpdateBase
1928!---------------------------------------------------------------------------------------------------
1929subroutine Agrif_Update_2D_Recursive ( type_update,         &
1930                                       tempP, tempC,        &
1931                                       indmin, indmax,      &
1932                                       lb_child, ub_child,  &
1933                                        s_child,  s_parent, &
1934                                       ds_child, ds_parent )
1935!---------------------------------------------------------------------------------------------------
1936    integer, dimension(2),          intent(in)  :: type_update            !< Type of update (copy or average)
1937    integer, dimension(2),          intent(in)  :: indmin, indmax
1938    integer, dimension(2),          intent(in)  :: lb_child, ub_child
1939    real,    dimension(2),          intent(in)  ::  s_child,  s_parent
1940    real,    dimension(2),          intent(in)  :: ds_child, ds_parent
1941    real,    dimension(          &
1942        indmin(1):indmax(1),     &
1943        indmin(2):indmax(2)),       intent(out) :: tempP
1944    real,    dimension(:,:),        intent(in)  :: tempC
1945!---------------------------------------------------------------------------------------------------
1946    real, dimension(indmin(1):indmax(1), lb_child(2):ub_child(2)) :: tabtemp
1947    real, dimension(indmin(2):indmax(2), indmin(1):indmax(1))     :: tempP_trsp
1948    real, dimension(lb_child(2):ub_child(2), indmin(1):indmax(1)) :: tabtemp_trsp
1949    integer :: i, j
1950    integer :: coeffraf
1951!
1952    tabtemp = 0.
1953    coeffraf = nint ( ds_parent(1) / ds_child(1) )
1954!
1955    IF((type_update(1) == Agrif_Update_Average) .AND. (coeffraf /= 1 )) THEN
1956!---CDIR NEXPAND
1957        if ( .NOT. precomputedone(1) ) then
1958            call Average1dPrecompute( ub_child(2)-lb_child(2)+1,    &
1959                                      indmax(1)-indmin(1)+1,              &
1960                                      ub_child(1)-lb_child(1)+1,    &
1961                                      s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
1962!            precomputedone(1) = .TRUE.
1963        endif
1964!---CDIR NEXPAND
1965        call Average1dAfterCompute( tabtemp, tempC, size(tabtemp), size(tempC), &
1966                    s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
1967!
1968    ELSE IF ((type_update(1) == Agrif_Update_Copy) .AND. (coeffraf /= 1 ))THEN
1969!---CDIR NEXPAND
1970         if ( .NOT. precomputedone(1) ) then
1971            call Agrif_basicupdate_copy1d_before( ub_child(2)-lb_child(2)+1, &
1972                                   indmax(1)-indmin(1)+1,           &
1973                                   ub_child(1)-lb_child(1)+1, &
1974                                   s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
1975!            precomputedone(1) = .TRUE.
1976        endif
1977!---CDIR NEXPAND
1978        call Agrif_basicupdate_copy1d_after(tabtemp,tempC,size(tabtemp),size(tempC),1)
1979!
1980    ELSE
1981        do j = lb_child(2),ub_child(2)
1982!
1983!---CDIR NEXPAND
1984            call Agrif_Update_1D_Recursive( type_update(1),             &
1985                                            tabtemp(:,j),               &
1986                                            tempC(:,j-lb_child(2)+1),   &
1987                                            indmin(1), indmax(1),       &
1988                                            lb_child(1),ub_child(1),    &
1989                                            s_child(1),  s_parent(1),   &
1990                                            ds_child(1),ds_parent(1))
1991        enddo
1992    ENDIF
1993!
1994    tabtemp_trsp = TRANSPOSE(tabtemp)
1995    coeffraf = nint(ds_parent(2)/ds_child(2))
1996!
1997    tempP_trsp = 0.
1998!
1999    IF((type_update(2) == Agrif_Update_Average) .AND. (coeffraf /= 1 )) THEN
2000!---CDIR NEXPAND
2001        if ( .NOT. precomputedone(2) ) then
2002            call Average1dPrecompute( indmax(1)-indmin(1)+1,          &
2003                                      indmax(2)-indmin(2)+1,          &
2004                                      ub_child(2)-lb_child(2)+1,&
2005                                      s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
2006!            precomputedone(2) = .TRUE.
2007        endif
2008!---CDIR NEXPAND
2009        call Average1dAfterCompute( tempP_trsp, tabtemp_trsp, size(tempP_trsp), size(tabtemp_trsp),&
2010                s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
2011!
2012    ELSE IF ((type_update(2) == Agrif_Update_Copy) .AND. (coeffraf /= 1 )) THEN
2013!---CDIR NEXPAND
2014        if ( .NOT. precomputedone(2) ) then
2015            call Agrif_basicupdate_copy1d_before( indmax(1)-indmin(1)+1,             &
2016                                   indmax(2)-indmin(2)+1,             &
2017                                   ub_child(2)-lb_child(2)+1,   &
2018                                   s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
2019!            precomputedone(2) = .TRUE.
2020        endif
2021!---CDIR NEXPAND
2022        call Agrif_basicupdate_copy1d_after( tempP_trsp, tabtemp_trsp, size(tempP_trsp), size(tabtemp_trsp),2)
2023!
2024    ELSE
2025        do i = indmin(1),indmax(1)
2026!
2027!---CDIR NEXPAND
2028            call Agrif_UpdateBase(type_update(2),                                       &
2029                                  tempP_trsp(indmin(2):indmax(2),i),            &
2030                                  tabtemp_trsp(lb_child(2):ub_child(2),i),&
2031                                  indmin(2),indmax(2),                          &
2032                                  lb_child(2),ub_child(2),                &
2033                                  s_parent(2),s_child(2),                       &
2034                                  ds_parent(2),ds_child(2))
2035!
2036        enddo
2037    ENDIF
2038!
2039    tempP = TRANSPOSE(tempP_trsp)
2040!---------------------------------------------------------------------------------------------------
2041end subroutine Agrif_Update_2D_Recursive
2042!===================================================================================================
2043!
2044subroutine Agrif_Update_2D_Recursive_ok ( type_update, &
2045                                        tempP, tempC, &
2046                                        indmin, indmax,   &
2047                                       lb_child, ub_child,                    &
2048                                       s_child, s_parent, ds_child, ds_parent )
2049!---------------------------------------------------------------------------------------------------
2050    INTEGER, DIMENSION(2), intent(in)   :: type_update            !< Type of update (copy or average)
2051    INTEGER, DIMENSION(2), intent(in)   :: indmin, indmax
2052    INTEGER, DIMENSION(2), intent(in)   :: lb_child, ub_child
2053    REAL,    DIMENSION(2), intent(in)   :: s_child,  s_parent
2054    REAL,    DIMENSION(2), intent(in)   :: ds_child, ds_parent
2055    REAL,    DIMENSION(                 &
2056                indmin(1):indmax(1),    &
2057                indmin(2):indmax(2)),           intent(out) :: tempP
2058    REAL, DIMENSION(                            &
2059                lb_child(1):ub_child(1),  &
2060                lb_child(2):ub_child(2)), intent(in)  :: tempC
2061!
2062    REAL, DIMENSION(indmin(1):indmax(1), lb_child(2):ub_child(2)) :: tabtemp
2063    INTEGER :: i
2064!
2065    do i = lb_child(2),ub_child(2)
2066        call Agrif_Update_1D_Recursive(type_update(1),                              &
2067                                       tabtemp(:, i),          &
2068                                       tempC(:,i),  &
2069                                       indmin(1),indmax(1),                 &
2070                                       lb_child(1),ub_child(1),       &
2071                                       s_child(1), s_parent(1),             &
2072                                      ds_child(1),ds_parent(1))
2073    enddo
2074!
2075    tempP = 0.
2076!
2077    do i = indmin(1),indmax(1)
2078        call Agrif_UpdateBase(type_update(2),                                       &
2079                              tempP(i,:),             &
2080                              tabtemp(i,:), &
2081                              indmin(2),indmax(2),                          &
2082                              lb_child(2),ub_child(2),                &
2083                              s_parent(2),s_child(2),                       &
2084                             ds_parent(2),ds_child(2))
2085    enddo
2086!---------------------------------------------------------------------------------------------------
2087end subroutine Agrif_Update_2D_Recursive_ok
2088!===================================================================================================
2089
2090!
2091!===================================================================================================
2092!  subroutine Agrif_Update_3D_Recursive
2093!
2094!> Updates a 3D grid variable on the parent grid.
2095!! Calls #Agrif_Update_2D_Recursive and #Agrif_UpdateBase.
2096!---------------------------------------------------------------------------------------------------
2097subroutine Agrif_Update_3D_Recursive ( type_update,         &
2098                                       tempP, tempC,        &
2099                                       indmin, indmax,      &
2100                                       lb_child, ub_child,  &
2101                                        s_child,  s_parent, &
2102                                       ds_child, ds_parent )
2103!---------------------------------------------------------------------------------------------------
2104    integer, dimension(3),          intent(in)  :: type_update            !< Type of update (copy or average)
2105    integer, dimension(3),          intent(in)  :: indmin, indmax
2106    integer, dimension(3),          intent(in)  :: lb_child, ub_child
2107    real,    dimension(3),          intent(in)  ::  s_child,  s_parent
2108    real,    dimension(3),          intent(in)  :: ds_child, ds_parent
2109    real,    dimension(          &
2110        indmin(1):indmax(1),     &
2111        indmin(2):indmax(2),     &
2112        indmin(3):indmax(3)),       intent(out) :: tempP
2113    real, dimension(             &
2114        lb_child(1):ub_child(1), &
2115        lb_child(2):ub_child(2), &
2116        lb_child(3):ub_child(3)),   intent(in)  :: tempC
2117!---------------------------------------------------------------------------------------------------
2118    real, dimension(            &
2119        indmin(1):indmax(1),    &
2120        indmin(2):indmax(2),    &
2121        lb_child(3):ub_child(3))    :: tabtemp
2122    integer :: i,j,k
2123    integer :: coeffraf,locind_child_left
2124    integer :: kuinf
2125!
2126    coeffraf = nint ( ds_parent(1) / ds_child(1) )
2127!
2128    if ((type_update(1) == Agrif_Update_Average) .AND. (coeffraf /= 1 )) then
2129!---CDIR NEXPAND
2130        call Average1dPrecompute(ub_child(2)-lb_child(2)+1,&
2131                                 indmax(1)-indmin(1)+1,&
2132                                 ub_child(1)-lb_child(1)+1,&
2133                                 s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
2134        precomputedone(1) = .TRUE.
2135    else if ((type_update(1) == Agrif_Update_Copy) .AND. (coeffraf /= 1 )) then
2136!---CDIR NEXPAND
2137         call Agrif_basicupdate_copy1d_before(ub_child(2)-lb_child(2)+1, &
2138                               indmax(1)-indmin(1)+1,           &
2139                               ub_child(1)-lb_child(1)+1, &
2140                               s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
2141        precomputedone(1) = .TRUE.
2142    endif
2143!
2144    coeffraf = nint ( ds_parent(2) / ds_child(2) )
2145!
2146    if ((type_update(2) == Agrif_Update_Average) .AND. (coeffraf /= 1 )) then
2147!---CDIR NEXPAND
2148        call Average1dPrecompute(indmax(1)-indmin(1)+1,&
2149                                 indmax(2)-indmin(2)+1,&
2150                                 ub_child(2)-lb_child(2)+1,&
2151                                 s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
2152        precomputedone(2) = .TRUE.
2153    else if ((type_update(2) == Agrif_Update_Copy) .AND. (coeffraf /= 1 )) then
2154!---CDIR NEXPAND
2155        call Agrif_basicupdate_copy1d_before( indmax(1)-indmin(1)+1,           &
2156                               indmax(2)-indmin(2)+1,           &
2157                               ub_child(2)-lb_child(2)+1, &
2158                               s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
2159        precomputedone(2) = .TRUE.
2160    endif
2161!
2162    do k = lb_child(3),ub_child(3)
2163        call Agrif_Update_2D_Recursive( type_update(1:2),tabtemp(:,:,k),tempC(:,:,k), &
2164                                        indmin(1:2),indmax(1:2),                &
2165                                        lb_child(1:2),ub_child(1:2),      &
2166                                        s_child(1:2),s_parent(1:2),             &
2167                                        ds_child(1:2),ds_parent(1:2) )
2168    enddo
2169!
2170    precomputedone(1) = .FALSE.
2171    precomputedone(2) = .FALSE.
2172!
2173    coeffraf = nint ( ds_parent(3) / ds_child(3) )
2174    locind_child_left = 1 + agrif_int((s_parent(3)-s_child(3))/ds_child(3))
2175!
2176    if (coeffraf == 1) then
2177        kuinf = lb_child(3)+locind_child_left-2
2178        do k=indmin(3),indmax(3)
2179            kuinf = kuinf + 1
2180            do j = indmin(2),indmax(2)
2181            do i = indmin(1),indmax(1)
2182                tempP(i,j,k) = tabtemp(i,j,kuinf)
2183            enddo
2184            enddo
2185        enddo
2186    else
2187        tempP = 0.
2188        do j = indmin(2),indmax(2)
2189        do i = indmin(1),indmax(1)
2190            call Agrif_UpdateBase(type_update(3),tempP(i,j,:),tabtemp(i,j,:),   &
2191                                  indmin(3),indmax(3),                  &
2192                                  lb_child(3),ub_child(3),        &
2193                                  s_parent(3),s_child(3),               &
2194                                  ds_parent(3),ds_child(3))
2195!
2196        enddo
2197        enddo
2198    endif
2199!---------------------------------------------------------------------------------------------------
2200end subroutine Agrif_Update_3D_Recursive
2201!===================================================================================================
2202!
2203!===================================================================================================
2204!  subroutine Agrif_Update_4D_Recursive
2205!
2206!> Updates a 4D grid variable on the parent grid.
2207!! Calls #Agrif_Update_3D_Recursive and #Agrif_UpdateBase.
2208!---------------------------------------------------------------------------------------------------
2209subroutine Agrif_Update_4D_Recursive ( type_update,         &
2210                                       tempP, tempC,        &
2211                                       indmin, indmax,      &
2212                                       lb_child, ub_child,  &
2213                                        s_child,  s_parent, &
2214                                       ds_child, ds_parent )
2215!---------------------------------------------------------------------------------------------------
2216    integer, dimension(4),          intent(in)  :: type_update            !< Type of update (copy or average)
2217    integer, dimension(4),          intent(in)  :: indmin, indmax
2218    integer, dimension(4),          intent(in)  :: lb_child, ub_child
2219    real,    dimension(4),          intent(in)  ::  s_child,  s_parent
2220    real,    dimension(4),          intent(in)  :: ds_child, ds_parent
2221    real,    dimension(          &
2222        indmin(1):indmax(1),     &
2223        indmin(2):indmax(2),     &
2224        indmin(3):indmax(3),     &
2225        indmin(4):indmax(4)),       intent(out) :: tempP
2226    real, dimension(             &
2227        lb_child(1):ub_child(1), &
2228        lb_child(2):ub_child(2), &
2229        lb_child(3):ub_child(3), &
2230        lb_child(4):ub_child(4)),   intent(in)  :: tempC
2231!---------------------------------------------------------------------------------------------------
2232    real, dimension(:,:,:,:), allocatable       :: tabtemp
2233    integer :: i,j,k,l
2234!
2235    allocate(tabtemp(indmin(1):indmax(1), &
2236                     indmin(2):indmax(2), &
2237                     indmin(3):indmax(3), &
2238                     lb_child(4):ub_child(4)))
2239!
2240    do l = lb_child(4), ub_child(4)
2241        call Agrif_Update_3D_Recursive(type_update(1:3),                    &
2242                                       tabtemp(indmin(1):indmax(1),         &
2243                                               indmin(2):indmax(2),         &
2244                                               indmin(3):indmax(3), l),     &
2245                                       tempC(lb_child(1):ub_child(1),       &
2246                                             lb_child(2):ub_child(2),       &
2247                                             lb_child(3):ub_child(3), l),   &
2248                                       indmin(1:3), indmax(1:3),            &
2249                                       lb_child(1:3), ub_child(1:3),        &
2250                                        s_child(1:3),  s_parent(1:3),       &
2251                                       ds_child(1:3), ds_parent(1:3))
2252    enddo
2253!
2254    tempP = 0.
2255!
2256    do k = indmin(3), indmax(3)
2257    do j = indmin(2), indmax(2)
2258    do i = indmin(1), indmax(1)
2259        call Agrif_UpdateBase(type_update(4),                               &
2260                              tempP(i,j,k,indmin(4):indmax(4)),             &
2261                              tabtemp(i,j,k,lb_child(4):ub_child(4)),       &
2262                              indmin(4), indmax(4),                         &
2263                              lb_child(4), ub_child(4),                     &
2264                               s_parent(4), s_child(4),                     &
2265                              ds_parent(4),ds_child(4) )
2266    enddo
2267    enddo
2268    enddo
2269!
2270    deallocate(tabtemp)
2271!---------------------------------------------------------------------------------------------------
2272end subroutine Agrif_Update_4D_Recursive
2273!===================================================================================================
2274!
2275!===================================================================================================
2276!  subroutine Agrif_Update_5D_Recursive
2277!
2278!> Updates a 5D grid variable on the parent grid.
2279!! Calls #Agrif_Update_4D_Recursive and #Agrif_UpdateBase.
2280!---------------------------------------------------------------------------------------------------
2281subroutine Agrif_Update_5D_Recursive ( type_update,         &
2282                                       tempP, tempC,        &
2283                                       indmin, indmax,      &
2284                                       lb_child, ub_child,  &
2285                                        s_child,  s_parent, &
2286                                       ds_child, ds_parent )
2287!---------------------------------------------------------------------------------------------------
2288    integer, dimension(5),          intent(in)  :: type_update            !< Type of update (copy or average)
2289    integer, dimension(5),          intent(in)  :: indmin, indmax
2290    integer, dimension(5),          intent(in)  :: lb_child, ub_child
2291    real,    dimension(5),          intent(in)  ::  s_child,  s_parent
2292    real,    dimension(5),          intent(in)  :: ds_child, ds_parent
2293    real,    dimension(          &
2294        indmin(1):indmax(1),     &
2295        indmin(2):indmax(2),     &
2296        indmin(3):indmax(3),     &
2297        indmin(4):indmax(4),     &
2298        indmin(5):indmax(5)),       intent(out) :: tempP
2299    real, dimension(             &
2300        lb_child(1):ub_child(1), &
2301        lb_child(2):ub_child(2), &
2302        lb_child(3):ub_child(3), &
2303        lb_child(4):ub_child(4), &
2304        lb_child(5):ub_child(5)),   intent(in)  :: tempC
2305!---------------------------------------------------------------------------------------------------
2306    real, dimension(:,:,:,:,:), allocatable     :: tabtemp
2307    integer :: i,j,k,l,m
2308!
2309    allocate(tabtemp(indmin(1):indmax(1), &
2310                     indmin(2):indmax(2), &
2311                     indmin(3):indmax(3), &
2312                     indmin(4):indmax(4), &
2313                     lb_child(5):ub_child(5)))
2314!
2315    do m = lb_child(5), ub_child(5)
2316        call Agrif_Update_4D_Recursive(type_update(1:4),                    &
2317                                       tabtemp(indmin(1):indmax(1),         &
2318                                               indmin(2):indmax(2),         &
2319                                               indmin(3):indmax(3),         &
2320                                               indmin(4):indmax(4), m),     &
2321                                       tempC(lb_child(1):ub_child(1),       &
2322                                             lb_child(2):ub_child(2),       &
2323                                             lb_child(3):ub_child(3),       &
2324                                             lb_child(4):ub_child(4), m),   &
2325                                       indmin(1:4),indmax(1:4),             &
2326                                       lb_child(1:4), ub_child(1:4),        &
2327                                        s_child(1:4), s_parent(1:4),        &
2328                                       ds_child(1:4), ds_parent(1:4))
2329    enddo
2330!
2331    tempP = 0.
2332!
2333    do l = indmin(4), indmax(4)
2334    do k = indmin(3), indmax(3)
2335    do j = indmin(2), indmax(2)
2336    do i = indmin(1), indmax(1)
2337        call Agrif_UpdateBase( type_update(5),                              &
2338                               tempP(i,j,k,l,indmin(5):indmax(5)),          &
2339                               tabtemp(i,j,k,l,lb_child(5):ub_child(5)),    &
2340                               indmin(5), indmax(5),                        &
2341                               lb_child(5), ub_child(5),                    &
2342                                s_parent(5), s_child(5),                    &
2343                               ds_parent(5),ds_child(5) )
2344    enddo
2345    enddo
2346    enddo
2347    enddo
2348!
2349    deallocate(tabtemp)
2350!---------------------------------------------------------------------------------------------------
2351end subroutine Agrif_Update_5D_Recursive
2352!===================================================================================================
2353!
2354!===================================================================================================
2355!  subroutine Agrif_Update_6D_Recursive
2356!
2357!> Updates a 6D grid variable on the parent grid.
2358!! Calls #Agrif_Update_5D_Recursive and #Agrif_UpdateBase.
2359!---------------------------------------------------------------------------------------------------
2360subroutine Agrif_Update_6D_Recursive ( type_update,         &
2361                                       tempP, tempC,        &
2362                                       indmin, indmax,      &
2363                                       lb_child, ub_child,  &
2364                                        s_child,  s_parent, &
2365                                       ds_child, ds_parent )
2366!---------------------------------------------------------------------------------------------------
2367    integer, dimension(6),          intent(in)  :: type_update            !< Type of update (copy or average)
2368    integer, dimension(6),          intent(in)  :: indmin, indmax
2369    integer, dimension(6),          intent(in)  :: lb_child, ub_child
2370    real,    dimension(6),          intent(in)  ::  s_child,  s_parent
2371    real,    dimension(6),          intent(in)  :: ds_child, ds_parent
2372    real,    dimension(          &
2373        indmin(1):indmax(1),     &
2374        indmin(2):indmax(2),     &
2375        indmin(3):indmax(3),     &
2376        indmin(4):indmax(4),     &
2377        indmin(5):indmax(5),     &
2378        indmin(6):indmax(6)),       intent(out) :: tempP
2379    real, dimension(             &
2380        lb_child(1):ub_child(1), &
2381        lb_child(2):ub_child(2), &
2382        lb_child(3):ub_child(3), &
2383        lb_child(4):ub_child(4), &
2384        lb_child(5):ub_child(5), &
2385        lb_child(6):ub_child(6)),   intent(in)  :: tempC
2386!---------------------------------------------------------------------------------------------------
2387    real, dimension(:,:,:,:,:,:), allocatable   :: tabtemp
2388    integer :: i,j,k,l,m,n
2389!
2390    allocate(tabtemp(indmin(1):indmax(1), &
2391                     indmin(2):indmax(2), &
2392                     indmin(3):indmax(3), &
2393                     indmin(4):indmax(4), &
2394                     indmin(5):indmax(5), &
2395                     lb_child(6):ub_child(6)))
2396!
2397    do n = lb_child(6),ub_child(6)
2398        call Agrif_Update_5D_Recursive(type_update(1:5),                    &
2399                                       tabtemp(indmin(1):indmax(1),         &
2400                                               indmin(2):indmax(2),         &
2401                                               indmin(3):indmax(3),         &
2402                                               indmin(4):indmax(4),         &
2403                                               indmin(5):indmax(5), n),     &
2404                                       tempC(lb_child(1):ub_child(1),       &
2405                                             lb_child(2):ub_child(2),       &
2406                                             lb_child(3):ub_child(3),       &
2407                                             lb_child(4):ub_child(4),       &
2408                                             lb_child(5):ub_child(5), n),   &
2409                                       indmin(1:5), indmax(1:5),            &
2410                                       lb_child(1:5),ub_child(1:5),         &
2411                                       s_child(1:5), s_parent(1:5),         &
2412                                      ds_child(1:5),ds_parent(1:5))
2413    enddo
2414!
2415    tempP = 0.
2416!
2417    do m = indmin(5), indmax(5)
2418    do l = indmin(4), indmax(4)
2419    do k = indmin(3), indmax(3)
2420    do j = indmin(2), indmax(2)
2421    do i = indmin(1), indmax(1)
2422        call Agrif_UpdateBase( type_update(6),                              &
2423                               tempP(i,j,k,l,m,indmin(6):indmax(6)),        &
2424                               tabtemp(i,j,k,l,m,lb_child(6):ub_child(6)),  &
2425                               indmin(6), indmax(6),                        &
2426                               lb_child(6),  ub_child(6),                   &
2427                                s_parent(6),  s_child(6),                   &
2428                               ds_parent(6), ds_child(6) )
2429    enddo
2430    enddo
2431    enddo
2432    enddo
2433    enddo
2434!
2435    deallocate(tabtemp)
2436!---------------------------------------------------------------------------------------------------
2437end subroutine Agrif_Update_6D_Recursive
2438!===================================================================================================
2439!
2440!===================================================================================================
2441!  subroutine Agrif_UpdateBase
2442!
2443!> Calls the updating method chosen by the user (copy, average or full-weighting).
2444!---------------------------------------------------------------------------------------------------
2445subroutine Agrif_UpdateBase ( type_update,              &
2446                              parent_tab, child_tab,    &
2447                              indmin, indmax,           &
2448                              lb_child, ub_child,       &
2449                              s_parent, s_child,        &
2450                              ds_parent, ds_child )
2451!---------------------------------------------------------------------------------------------------
2452    integer,                            intent(in) :: type_update
2453    integer,                            intent(in) :: indmin, indmax
2454    integer,                            intent(in) :: lb_child, ub_child
2455    real, dimension(indmin:indmax),     intent(out):: parent_tab
2456    real, dimension(lb_child:ub_child), intent(in) :: child_tab
2457    real,                               intent(in) :: s_parent,  s_child
2458    real,                               intent(in) :: ds_parent, ds_child
2459!---------------------------------------------------------------------------------------------------
2460    integer :: np       ! Length of parent array
2461    integer :: nc       ! Length of child  array
2462!
2463    np = indmax - indmin + 1
2464    nc = ub_child - lb_child + 1
2465!
2466    if     ( type_update == Agrif_Update_Copy ) then
2467!
2468        call Agrif_basicupdate_copy1d(              &
2469                    parent_tab, child_tab,          &
2470                    np,         nc,                 &
2471                    s_parent,  s_child,             &
2472                    ds_parent, ds_child )
2473!
2474    elseif ( type_update == Agrif_Update_Average ) then
2475!
2476        call Agrif_basicupdate_average1d(           &
2477                    parent_tab, child_tab,          &
2478                    np,         nc,                 &
2479                     s_parent,  s_child,            &
2480                    ds_parent, ds_child )
2481!
2482    elseif ( type_update == Agrif_Update_Max ) then
2483!
2484        call Agrif_basicupdate_max1d(           &
2485                    parent_tab, child_tab,          &
2486                    np,         nc,                 &
2487                     s_parent,  s_child,            &
2488                    ds_parent, ds_child )
2489    elseif ( type_update == Agrif_Update_Full_Weighting ) then
2490!
2491        call Agrif_basicupdate_full_weighting1D(    &
2492                    parent_tab, child_tab,          &
2493                    np,         nc,                 &
2494                    s_parent,  s_child,             &
2495                    ds_parent, ds_child )
2496!
2497    endif
2498!---------------------------------------------------------------------------------------------------
2499end subroutine Agrif_UpdateBase
2500!===================================================================================================
2501!
2502#if defined AGRIF_MPI
2503!===================================================================================================
2504!  subroutine Agrif_Find_list_update
2505!---------------------------------------------------------------------------------------------------
2506subroutine Agrif_Find_list_update ( list_update, pttab, petab, lb_child, lb_parent, nbdim, &
2507                                    find_list_update, tab4t, tab5t, memberinall, memberinall2,   &
2508                                    sendtoproc1, recvfromproc1, sendtoproc2, recvfromproc2 )
2509!---------------------------------------------------------------------------------------------------
2510    TYPE(Agrif_List_Interp_Loc),                   pointer     :: list_update
2511    INTEGER,                                       intent(in)  :: nbdim
2512    INTEGER, DIMENSION(nbdim),                     intent(in)  :: pttab, petab
2513    INTEGER, DIMENSION(nbdim),                     intent(in)  :: lb_child, lb_parent
2514    LOGICAL,                                       intent(out) :: find_list_update
2515    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1,8), intent(out) :: tab4t
2516    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1,8), intent(out) :: tab5t
2517    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1), intent(out) :: memberinall,memberinall2
2518    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1), intent(out) :: sendtoproc1,recvfromproc1
2519    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1), intent(out) :: sendtoproc2,recvfromproc2
2520!
2521    Type(Agrif_List_Interp_Loc), Pointer :: parcours
2522    INTEGER :: i
2523!
2524    find_list_update = .FALSE.
2525!
2526    parcours => list_update
2527
2528    Find_loop :  do while ( associated(parcours) )
2529        do i = 1,nbdim
2530            IF ((pttab(i) /= parcours%interp_loc%pttab(i)) .OR. &
2531                (petab(i) /= parcours%interp_loc%petab(i)) .OR. &
2532                (lb_child(i)  /= parcours%interp_loc%pttab_child(i)) .OR. &
2533                (lb_parent(i) /= parcours%interp_loc%pttab_parent(i))) THEN
2534                parcours => parcours%suiv
2535                cycle Find_loop
2536            ENDIF
2537        enddo
2538!
2539        tab4t = parcours%interp_loc%tab4t(1:nbdim,0:Agrif_Nbprocs-1,1:8)
2540        tab5t = parcours%interp_loc%tab5t(1:nbdim,0:Agrif_Nbprocs-1,1:8)
2541        memberinall  =  parcours%interp_loc%memberinall(0:Agrif_Nbprocs-1)
2542        memberinall2 =  parcours%interp_loc%memberinall2(0:Agrif_Nbprocs-1)
2543        sendtoproc1 =   parcours%interp_loc%sendtoproc1(0:Agrif_Nbprocs-1)
2544        sendtoproc2 =   parcours%interp_loc%sendtoproc2(0:Agrif_Nbprocs-1)
2545        recvfromproc1 = parcours%interp_loc%recvfromproc1(0:Agrif_Nbprocs-1)
2546        recvfromproc2 = parcours%interp_loc%recvfromproc2(0:Agrif_Nbprocs-1)
2547!
2548        find_list_update = .TRUE.
2549        exit Find_loop
2550!
2551    enddo Find_loop
2552!---------------------------------------------------------------------------------------------------
2553end subroutine Agrif_Find_list_update
2554!===================================================================================================
2555!
2556!===================================================================================================
2557!  subroutine Agrif_AddTo_list_update
2558!---------------------------------------------------------------------------------------------------
2559subroutine Agrif_AddTo_list_update ( list_update, pttab, petab, lb_child, lb_parent,  &
2560                                     nbdim, tab4t, tab5t, memberinall, memberinall2,        &
2561                                     sendtoproc1, recvfromproc1, sendtoproc2, recvfromproc2 )
2562!---------------------------------------------------------------------------------------------------
2563    TYPE(Agrif_List_Interp_Loc), pointer :: list_update
2564    INTEGER,                                        intent(in) :: nbdim
2565    INTEGER, DIMENSION(nbdim),                      intent(in) :: pttab, petab
2566    INTEGER, DIMENSION(nbdim),                      intent(in) :: lb_child, lb_parent
2567    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1,8),  intent(in) :: tab4t
2568    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1,8),  intent(in) :: tab5t
2569    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1),          intent(in) :: memberinall, memberinall2
2570    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1),          intent(in) :: sendtoproc1, recvfromproc1
2571    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1),          intent(in) :: sendtoproc2, recvfromproc2
2572!
2573    Type(Agrif_List_Interp_Loc), pointer :: parcours
2574!
2575    allocate(parcours)
2576    allocate(parcours%interp_loc)
2577
2578    parcours%interp_loc%pttab(1:nbdim) = pttab(1:nbdim)
2579    parcours%interp_loc%petab(1:nbdim) = petab(1:nbdim)
2580    parcours%interp_loc%pttab_child(1:nbdim)  = lb_child(1:nbdim)
2581    parcours%interp_loc%pttab_parent(1:nbdim) = lb_parent(1:nbdim)
2582
2583    allocate(parcours%interp_loc%tab4t(nbdim,0:Agrif_Nbprocs-1,8))
2584    allocate(parcours%interp_loc%tab5t(nbdim,0:Agrif_Nbprocs-1,8))
2585
2586    allocate(parcours%interp_loc%memberinall (0:Agrif_Nbprocs-1))
2587    allocate(parcours%interp_loc%memberinall2(0:Agrif_Nbprocs-1))
2588
2589    allocate(parcours%interp_loc%recvfromproc1(0:Agrif_Nbprocs-1))
2590    allocate(parcours%interp_loc%recvfromproc2(0:Agrif_Nbprocs-1))
2591    allocate(parcours%interp_loc%sendtoproc1(0:Agrif_Nbprocs-1))
2592    allocate(parcours%interp_loc%sendtoproc2(0:Agrif_Nbprocs-1))
2593
2594    parcours%interp_loc%tab4t = tab4t
2595    parcours%interp_loc%tab5t = tab5t
2596    parcours%interp_loc%memberinall   = memberinall
2597    parcours%interp_loc%memberinall2  = memberinall2
2598    parcours%interp_loc%sendtoproc1   = sendtoproc1
2599    parcours%interp_loc%sendtoproc2   = sendtoproc2
2600    parcours%interp_loc%recvfromproc1 = recvfromproc1
2601    parcours%interp_loc%recvfromproc2 = recvfromproc2
2602
2603    parcours%suiv => list_update
2604    list_update => parcours
2605!---------------------------------------------------------------------------------------------------
2606end subroutine Agrif_Addto_list_update
2607!===================================================================================================
2608#endif
2609!
2610end module Agrif_Update
Note: See TracBrowser for help on using the repository browser.