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

source: vendors/AGRIF/dev/AGRIF_FILES/modinterp.F90

Last change on this file was 14975, checked in by jchanut, 3 years ago

#2638, merge new AGRIF library into trunk

  • Property svn:keywords set to Id
File size: 105.5 KB
Line 
1!
2! $Id$
3!
4!     AGRIF (Adaptive Grid Refinement In Fortran)
5!
6!     Copyright (C) 2003 Laurent Debreu (Laurent.Debreu@imag.fr)
7!                        Christophe Vouland (Christophe.Vouland@imag.fr)
8!
9!     This program is free software; you can redistribute it and/or modify
10!     it under the terms of the GNU General Public License as published by
11!     the Free Software Foundation; either version 2 of the License, or
12!     (at your option) any later version.
13!
14!     This program is distributed in the hope that it will be useful,
15!     but WITHOUT ANY WARRANTY; without even the implied warranty of
16!     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17!     GNU General Public License for more details.
18!
19!     You should have received a copy of the GNU General Public License
20!     along with this program; if not, write to the Free Software
21!     Foundation, Inc., 59 Temple Place- Suite 330, Boston, MA 02111-1307, USA.
22!
23!
24!> Module to initialize a fine grid from its parent grid, by using a space interpolation
25!
26module Agrif_Interpolation
27!
28    use Agrif_InterpBasic
29    use Agrif_Arrays
30    use Agrif_Mask
31    use Agrif_CurgridFunctions
32#if defined AGRIF_MPI
33    use Agrif_Mpp
34#endif
35!
36    implicit none
37!
38    logical, private:: precomputedone(7) = .FALSE.
39!
40    private :: Agrif_Parentbounds
41    private :: Agrif_Interp_1D_recursive, Agrif_Interp_2D_recursive, Agrif_Interp_3D_recursive
42    private :: Agrif_Interp_4D_recursive, Agrif_Interp_5D_recursive, Agrif_Interp_6D_recursive
43    private :: Agrif_InterpBase
44    private :: Agrif_Find_list_interp, Agrif_AddTo_list_interp
45!
46contains
47!
48!===================================================================================================
49!  subroutine Agrif_InterpVariable
50!
51!> Sets some arguments of subroutine Agrif_InterpnD, n being the dimension of the grid variable
52!---------------------------------------------------------------------------------------------------
53subroutine Agrif_InterpVariable ( parent, child, torestore, procname )
54!---------------------------------------------------------------------------------------------------
55    type(Agrif_Variable), pointer       :: parent       !< Variable on the parent grid
56    type(Agrif_Variable), pointer       :: child        !< Variable on the child grid
57    logical,              intent(in)    :: torestore    !< .false. indicates that the results of the
58                                                        !! interpolation are applied on the whole current grid
59    procedure()                         :: procname     !< Data recovery procedure
60!---------------------------------------------------------------------------------------------------
61    logical               :: memberin
62    integer               :: nbdim        ! Number of dimensions of the current grid
63    integer, dimension(6) :: type_interp  ! Type of interpolation (linear,spline,...)
64    integer, dimension(6) :: nb_child
65    integer, dimension(6) :: lb_child
66    integer, dimension(6) :: ub_child
67    integer, dimension(6) :: lb_parent
68    real(kind=8)   , dimension(6) :: s_child,   s_parent
69    real(kind=8)   , dimension(6) :: ds_child, ds_parent
70    integer, dimension(child % root_var % nbdim,2,2) :: childarray
71!
72    nbdim       = child % root_var % nbdim
73    type_interp = child % root_var % type_interp
74!
75    call PreProcessToInterpOrUpdate( parent,   child,       &
76                                     nb_child, ub_child,    &
77                                     lb_child, lb_parent,   &
78                                      s_child,  s_parent,   &
79                                     ds_child, ds_parent, nbdim, interp=.true.)
80!
81!   Call to a procedure of interpolation against the number of dimensions of the grid variable
82!
83    call Agrif_InterpnD(type_interp, parent, child, &
84                        lb_child, ub_child,         &
85                        lb_child, lb_parent,        &
86                        s_child,   s_parent,        &
87                        ds_child, ds_parent,        &
88                        child, torestore, nbdim,    &
89                        childarray, memberin,       &
90                        .false., procname, 0, 0)
91!---------------------------------------------------------------------------------------------------
92end subroutine Agrif_InterpVariable
93!===================================================================================================
94!
95!===================================================================================================
96!  subroutine Agrif_InterpnD
97!
98!> Interpolates a nD grid variable from its parent grid, by using a space interpolation
99!---------------------------------------------------------------------------------------------------
100subroutine Agrif_InterpnD ( type_interp, parent, child, pttab, petab, pttab_Child, pttab_Parent, &
101                            s_Child, s_Parent, ds_Child, ds_Parent, restore, torestore,          &
102                            nbdim, childarray, memberin, in_bc, procname, nb, ndir )
103!---------------------------------------------------------------------------------------------------
104#if defined AGRIF_MPI
105    include 'mpif.h'
106#endif
107!
108    INTEGER, DIMENSION(6),     INTENT(in)   :: type_interp  !< Type of interpolation !    (linear,...)
109    TYPE(Agrif_Variable),      pointer      :: parent       !< Variable of the parent grid
110    TYPE(Agrif_Variable),      pointer      :: child        !< Variable of the child grid
111    INTEGER, DIMENSION(nbdim), INTENT(in)   :: pttab        !< Index of the first point inside the domain
112    INTEGER, DIMENSION(nbdim), INTENT(in)   :: petab        !< Index of the first point inside the domain
113    INTEGER, DIMENSION(nbdim), INTENT(in)   :: pttab_Child  !< Index of the first point inside the domain
114                                                            !<    for the child grid variable
115    INTEGER, DIMENSION(nbdim), INTENT(in)   :: pttab_Parent !< Index of the first point inside the domain
116                                                            !<    for the parent grid variable
117    REAL(kind=8),    DIMENSION(nbdim), INTENT(in)   :: s_Child,s_Parent   !< Positions of the parent and child grids
118    REAL(kind=8),    DIMENSION(nbdim), INTENT(in)   :: ds_Child,ds_Parent !< Space steps of the parent and child grids
119    TYPE(Agrif_Variable),      pointer      :: restore            !< Indicates points where interpolation
120    LOGICAL,                   INTENT(in)   :: torestore          !< Indicates if the array restore is used
121    INTEGER,                   INTENT(in)   :: nbdim
122    LOGICAL,                   INTENT(out)  :: memberin
123    LOGICAL,                   INTENT(in)   :: in_bc              !< .true. if called from Agrif_CorrectVariable \n
124                                                                  !! .false. if called from Agrif_InterpVariable
125    procedure()                             :: procname           !< Data recovery procedure
126    INTEGER,                   INTENT(in)   :: nb, ndir
127!
128    INTEGER                       :: i,j,k,l,m,n
129    integer :: i1,j1,k1
130    INTEGER, DIMENSION(nbdim)     :: pttruetab,cetruetab
131    INTEGER, DIMENSION(nbdim)     :: indmin,     indmax
132    INTEGER, DIMENSION(nbdim)     :: indminglob, indmaxglob
133#if defined AGRIF_MPI
134    INTEGER, DIMENSION(nbdim)     :: indminglob2,indmaxglob2
135    INTEGER, DIMENSION(nbdim)     :: indminglob3,indmaxglob3
136    INTEGER, DIMENSION(:,:),ALLOCATABLE :: indminglob_chunks, indmaxglob_chunks
137    INTEGER, DIMENSION(:,:),ALLOCATABLE :: indminglob2_chunks,indmaxglob2_chunks
138    INTEGER, DIMENSION(:,:),ALLOCATABLE :: indminglob3_chunks,indmaxglob3_chunks
139#endif
140    LOGICAL, DIMENSION(nbdim)     :: noraftab
141    REAL(kind=8)   , DIMENSION(nbdim)     :: s_Child_temp,s_Parent_temp
142    INTEGER, DIMENSION(nbdim)     :: lowerbound, upperbound, coords
143    INTEGER, DIMENSION(nbdim,2,2), INTENT(OUT) :: childarray
144    INTEGER, DIMENSION(nbdim,2,2)              :: parentarray
145    INTEGER, DIMENSION(nbdim,2,2)              :: parentarray_decal
146    LOGICAL :: member
147    LOGICAL,DIMENSION(:),ALLOCATABLE :: member_chuncks
148    INTEGER,DIMENSION(:,:),ALLOCATABLE :: decal_chunks
149    LOGICAL :: find_list_interp
150!
151#if defined AGRIF_MPI
152!
153    INTEGER, PARAMETER          :: etiquette = 100
154    INTEGER                     :: code, local_proc
155    INTEGER, DIMENSION(nbdim,4)                     :: tab3
156    INTEGER, DIMENSION(nbdim,4,0:Agrif_Nbprocs-1)   :: tab4
157    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1,8)   :: tab4t
158    INTEGER, DIMENSION(nbdim,2) :: tab5
159    INTEGER, DIMENSION(nbdim,2,0:Agrif_Nbprocs-1) :: tab6
160    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1,2) :: tab5t
161    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1)           :: memberinall
162    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1)           :: sendtoproc1, recvfromproc1
163    LOGICAL, DIMENSION(1)                           :: memberin1
164    LOGICAL                                         :: memberout
165!
166#endif
167! CHUNK (periodicity)
168    INTEGER :: nb_chunks
169    INTEGER :: agrif_external_switch_index
170    INTEGER, DIMENSION(2) :: test_orientation
171    !INTEGER, DIMENSION(2,nbdim,2,2) :: parentarray_chunk
172    INTEGER, DIMENSION(:,:,:,:), ALLOCATABLE :: parentarray_chunk
173    INTEGER, DIMENSION(:,:,:,:), ALLOCATABLE :: parentarray_chunk_decal
174    INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: bounds_chunks
175    logical,dimension(:),allocatable :: correction_required
176!
177    type(Agrif_Variable), pointer, save :: tempC => NULL()        ! Temporary child grid variable
178    type(Agrif_Variable), pointer, save :: tempP => NULL()        ! Temporary parent grid variable
179    type(Agrif_Variable), pointer, save :: tempPextend => NULL()  ! Temporary parent grid variable
180    type(Agrif_Variable), pointer, save :: parentvalues => NULL()
181!
182    coords = child % root_var % coords
183!
184!   Boundaries of the current grid where interpolation is done
185    find_list_interp =                                              &
186        Agrif_Find_list_interp(                                     &
187            child % list_interp,                                    &
188            pttab, petab, pttab_Child, pttab_Parent, nbdim,         &
189            indmin, indmax, indminglob, indmaxglob,                 &
190            pttruetab, cetruetab, memberin,                         &
191            parentarray_chunk,parentarray_chunk_decal,decal_chunks,  &
192            correction_required,member_chuncks,nb_chunks            &
193#if defined AGRIF_MPI
194           ,indminglob2, indmaxglob2, parentarray,                  &
195            member, tab4t,memberinall,  sendtoproc1, recvfromproc1  &
196#endif
197        )
198!
199    if (.not.find_list_interp) then
200!
201        call Agrif_get_var_bounds_array(child, lowerbound, upperbound, nbdim, parent)
202        call Agrif_Childbounds(nbdim, lowerbound, upperbound,               &
203                               pttab, petab, Agrif_Procrank, coords,        &
204                               pttruetab, cetruetab, memberin)
205
206        if (agrif_debug_interp) then
207        print *,'************CHILDBOUNDS*********************************'
208#ifdef AGRIF_MPI
209         print *,'Processeur ',Agrif_Procrank
210#endif
211        print *,'memberin ',memberin
212        do i = 1 , nbdim
213        print *,'Direction ',i,' indices debut: ',pttab(i),pttruetab(i)
214        print *,'Direction ',i,' indices fin  : ',petab(i),cetruetab(i)
215        enddo
216        print *,'*********************************************'
217        endif
218
219        call Agrif_Parentbounds(type_interp,nbdim,indminglob,indmaxglob,    &
220                                s_Parent_temp,s_Child_temp,                 &
221                                s_Child,ds_Child,                           &
222                                s_Parent,ds_Parent,                         &
223                                pttab,petab,                                &
224                                pttab_Child,pttab_Parent,                   &
225                                child%root_var % posvar, coords)
226
227        if (agrif_debug_interp) then
228        print *,'************PARENTBOUNDS*********************************'
229#ifdef AGRIF_MPI
230         print *,'Processeur ',Agrif_Procrank
231#endif
232        do i = 1 , nbdim
233        print *,'Direction ',i,' indices debut: ',pttab(i),indminglob(i)
234        print *,'Direction ',i,' indices fin  : ',petab(i),indmaxglob(i)
235        enddo
236
237        do i = 1 , nbdim
238        print *,'Direction ',i,' s_parent_temp: ',s_parent_temp(i)
239        print *,'Direction ',i,' s_Child_temp : ',s_Child_temp(i)
240        enddo
241        print *,'*********************************************'
242        endif
243
244#if defined AGRIF_MPI
245        if (memberin) then
246            call Agrif_Parentbounds(type_interp,nbdim,indmin,indmax,        &
247                                    s_Parent_temp,s_Child_temp,             &
248                                    s_Child,ds_Child,                       &
249                                    s_Parent,ds_Parent,                     &
250                                    pttruetab,cetruetab,                    &
251                                    pttab_Child,pttab_Parent,               &
252                                    child%root_var % posvar, coords)
253
254        endif
255        if (agrif_debug_interp) then
256        print *,'************PARENTBOUNDSMPI*********************************'
257#ifdef AGRIF_MPI
258         print *,'Processeur ',Agrif_Procrank
259#endif
260        do i = 1 , nbdim
261        print *,'Direction ',i,' indices debut: ',pttruetab(i),indmin(i)
262        print *,'Direction ',i,' indices fin  : ',cetruetab(i),indmax(i)
263        enddo
264
265        do i = 1 , nbdim
266        print *,'Direction ',i,' s_parent_temp: ',s_parent_temp(i)
267        print *,'Direction ',i,' s_Child_temp : ',s_Child_temp(i)
268        enddo
269        print *,'*********************************************'
270        endif
271
272        local_proc = Agrif_Procrank
273        call Agrif_get_var_bounds_array(parent,lowerbound,upperbound,nbdim)
274        call Agrif_ChildGrid_to_ParentGrid()
275
276          parentarray(:,1,1) = indminglob
277          parentarray(:,2,1) = indmaxglob
278          parentarray(:,1,2) = indminglob
279          parentarray(:,2,2) = indmaxglob
280
281        if (associated(agrif_external_mapping)) then
282
283          call agrif_external_mapping(nbdim,child%root_var % posvar(1),child%root_var % posvar(2), &
284                                      parentarray,parentarray_chunk,correction_required,nb_chunks)
285          allocate(decal_chunks(nb_chunks,nbdim))
286          do i=1,nb_chunks
287            decal_chunks(i,:)=parentarray_chunk(i,:,1,1)-parentarray_chunk(i,:,1,2)
288          enddo
289        else
290            nb_chunks=1
291            allocate(correction_required(nb_chunks))
292            correction_required=.FALSE.
293            allocate(parentarray_chunk(nb_chunks,nbdim,2,2))
294            parentarray_chunk(1,:,:,:)=parentarray
295            allocate(decal_chunks(nb_chunks,nbdim))
296            decal_chunks=0
297        endif
298        if (agrif_debug_interp) then
299        print *,'AVANT PARENTCHILDBOUNDS'
300        print *,'nombre de chunks ',nb_chunks
301        do i=1,nb_chunks
302          print *,'CHUNK Number : ',i
303          do j=1,nbdim
304           print *,'Direction ',j
305           print *,'MIN MAX (2) = ',parentarray_chunk(i,j,1,2),parentarray_chunk(i,j,2,2)
306           print *,'MIN MAX (1) = ',parentarray_chunk(i,j,1,1),parentarray_chunk(i,j,2,1)
307          enddo
308        enddo
309        print *,'APRES PARENTCHILDBOUNDS'
310        endif
311
312        allocate(indminglob_chunks(nb_chunks,nbdim))
313        allocate(indmaxglob_chunks(nb_chunks,nbdim))
314        allocate(indminglob2_chunks(nb_chunks,nbdim))
315        allocate(indmaxglob2_chunks(nb_chunks,nbdim))
316        allocate(indminglob3_chunks(nb_chunks,nbdim))
317        allocate(indmaxglob3_chunks(nb_chunks,nbdim))
318        allocate(member_chuncks(nb_chunks))
319
320        do i=1,nb_chunks
321          indminglob_chunks(i,:) = parentarray_chunk(i,:,1,2)
322          indmaxglob_chunks(i,:) = parentarray_chunk(i,:,2,2)
323        enddo
324
325        do i=1,nb_chunks
326
327              call Agrif_Childbounds(nbdim,lowerbound,upperbound,         &
328                   indminglob_chunks(i,:),indmaxglob_chunks(i,:), local_proc, coords,   &
329                   indminglob2_chunks(i,:),indmaxglob2_chunks(i,:),member_chuncks(i),   &
330                   indminglob3_chunks(i,:),indmaxglob3_chunks(i,:))
331        enddo
332!
333        ! call Agrif_Childbounds(nbdim,lowerbound,upperbound,                 &
334                               ! indminglob,indmaxglob, local_proc, coords,   &
335                               ! indminglob2,indmaxglob2,member,              &
336                               ! indminglob3,indmaxglob3,check_perio=.TRUE.)
337
338        if (agrif_debug_interp) then
339        print *,'************CHILDBOUNDSPARENTMPI*********************************'
340#ifdef AGRIF_MPI
341         print *,'Processeur ',Agrif_Procrank
342#endif
343        do j=1,nb_chunks
344        print *,'Chunk number ',j
345
346        do i = 1 , nbdim
347        print *,'Direction ',i,' indices debut: ',indminglob_chunks(j,i),indminglob2_chunks(j,i),indminglob3_chunks(j,i)
348        print *,'Direction ',i,' indices fin  : ',indmaxglob_chunks(j,i),indmaxglob2_chunks(j,i),indmaxglob3_chunks(j,i)
349        enddo
350        enddo
351        print *,'*********************************************'
352        endif
353!
354        ! if (member) then
355            ! call Agrif_GlobalToLocalBounds(parentarray,                     &
356                                           ! lowerbound,  upperbound,         &
357                                           ! indminglob2, indmaxglob2, coords,&
358                                           ! nbdim, local_proc, member,check_perio=.TRUE.)
359            ! if (agrif_debug_interp) then
360            ! do i=1,nbdim
361            ! print *,'parentarray = ',i,parentarray(i,1,1),parentarray(i,2,1), &
362                ! parentarray(i,1,2),parentarray(i,2,2)
363            ! enddo
364            ! endif
365        ! endif
366
367        allocate(parentarray_chunk_decal(nb_chunks,nbdim,2,2))
368        do j=1,nb_chunks
369        if (agrif_debug_interp) print *,'CHUNK = ',j
370        if (member_chuncks(j)) then
371            ! call Agrif_GlobalToLocalBounds(parentarray_chunk(j,:,:,:),                               &
372                                           ! lowerbound,  upperbound,                                  &
373                                           ! indminglob2_chunks(j,:), indmaxglob2_chunks(j,:), coords, &
374                                           ! nbdim, local_proc, member_chuncks(j),check_perio=.TRUE.)
375
376            call Agrif_GlobalToLocalBounds(parentarray_chunk(j,:,:,:),                               &
377                                           lowerbound,  upperbound,                                  &
378                                           indminglob2_chunks(j,:), indmaxglob2_chunks(j,:), coords, &
379                                           nbdim, local_proc, member_chuncks(j))
380
381            if (correction_required(j)) then
382                do i=1,2
383                    test_orientation(1)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
384                    parentarray_chunk(j,i,1,1),i)
385                    test_orientation(2)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
386                    parentarray_chunk(j,i,2,1),i)
387                    parentarray_chunk_decal(j,i,1,1)=minval(test_orientation)
388                    parentarray_chunk_decal(j,i,2,1)=maxval(test_orientation)
389                enddo
390                do i=3,nbdim
391                  parentarray_chunk_decal(j,i,:,1) = parentarray_chunk(j,i,:,1)+decal_chunks(j,i)
392                enddo
393            else
394            do i=1,nbdim
395            parentarray_chunk_decal(j,i,:,1) = parentarray_chunk(j,i,:,1)+decal_chunks(j,i)
396            enddo
397            endif
398
399            if (agrif_debug_interp) then
400            do i=1,nbdim
401            print *,'parentarray = ',i,parentarray_chunk(j,i,1,1),parentarray_chunk(j,i,2,1), &
402                parentarray_chunk(j,i,1,2),parentarray_chunk(j,i,2,2)
403                print *,'parentarraydecal = ',i,parentarray_chunk_decal(j,i,1,1),parentarray_chunk_decal(j,i,2,1)
404            enddo
405            endif
406        endif
407        enddo
408
409        parentarray(:,1,:)=Huge(1)
410        parentarray(:,2,:)=-Huge(1)
411        indminglob2=Huge(1)
412        indmaxglob2=-Huge(1)
413        indminglob3=Huge(1)
414        indmaxglob3=-Huge(1)
415        member = .FALSE.
416        do j=1,nb_chunks
417          if (member_chuncks(j)) then
418                           do i=1,nbdim
419             parentarray(i,1,1) = min(parentarray(i,1,1),parentarray_chunk_decal(j,i,1,1))
420             parentarray(i,1,2) = min(parentarray(i,1,2),parentarray_chunk(j,i,1,2))
421             parentarray(i,2,1) = max(parentarray(i,2,1),parentarray_chunk_decal(j,i,2,1))
422             parentarray(i,2,2) = max(parentarray(i,2,2),parentarray_chunk(j,i,2,2))
423                         enddo
424            if (correction_required(j)) then
425               if (agrif_debug_interp) then
426                  do i=1,nbdim
427                     print *,'direction ',i
428                  print *,'glob2_chuk = ',indminglob2_chunks(j,i),indmaxglob2_chunks(j,i)
429                  print *,'glob3_chuk = ',indminglob3_chunks(j,i),indmaxglob3_chunks(j,i)
430               enddo
431               endif
432                do i=1,2
433                    test_orientation(1)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
434                    indminglob2_chunks(j,i),i)
435                    test_orientation(2)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
436                    indmaxglob2_chunks(j,i),i)
437                    indminglob2(i)=min(indminglob2(i),minval(test_orientation))
438                    indmaxglob2(i)=max(indmaxglob2(i),maxval(test_orientation))
439                enddo
440
441                do i=1,2
442                    test_orientation(1)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
443                    indminglob3_chunks(j,i),i)
444                    test_orientation(2)=agrif_external_switch_index(child%root_var % posvar(1),child%root_var % posvar(2), &
445                    indmaxglob3_chunks(j,i),i)
446                    indminglob3(i)=min(indminglob3(i),minval(test_orientation))
447                    indmaxglob3(i)=max(indmaxglob3(i),maxval(test_orientation))
448                enddo
449
450                do i=3,nbdim
451             indminglob2(i)=min(indminglob2(i),indminglob2_chunks(j,i)+decal_chunks(j,i))
452             indmaxglob2(i)=max(indmaxglob2(i),indmaxglob2_chunks(j,i)+decal_chunks(j,i))
453             indminglob3(i)=min(indminglob3(i),indminglob3_chunks(j,i)+decal_chunks(j,i))
454             indmaxglob3(i)=max(indmaxglob3(i),indmaxglob3_chunks(j,i)+decal_chunks(j,i))
455                enddo
456            else
457                           do i=1,nbdim
458             indminglob2(i)=min(indminglob2(i),indminglob2_chunks(j,i)+decal_chunks(j,i))
459             indmaxglob2(i)=max(indmaxglob2(i),indmaxglob2_chunks(j,i)+decal_chunks(j,i))
460             indminglob3(i)=min(indminglob3(i),indminglob3_chunks(j,i)+decal_chunks(j,i))
461             indmaxglob3(i)=max(indmaxglob3(i),indmaxglob3_chunks(j,i)+decal_chunks(j,i))
462                         enddo
463         endif
464
465            member = .TRUE.
466          endif
467        enddo
468
469        call Agrif_ParentGrid_to_ChildGrid()
470
471        if (agrif_debug_interp) then
472        print *,'************ FINAL PARENTARRAY *****************'
473#ifdef AGRIF_MPI
474        print *,'Processeur ',Agrif_Procrank,' MEMBER = ',member
475        do i=1,nbdim
476         print *,'Direction ',i,' indices debut = ',parentarray(i,1,1),parentarray(i,1,2)
477         print *,'Direction ',i,' indices fin = ',parentarray(i,2,1),parentarray(i,2,2)
478        enddo
479#endif
480        endif
481
482        if (agrif_debug_interp) then
483        print *,'************ FINAL INDMINGLOB *****************'
484#ifdef AGRIF_MPI
485        print *,'Processeur ',Agrif_Procrank,' MEMBER = ',member
486        do i=1,nbdim
487         print *,'Direction ',i,' indices debut = ',indminglob2(i),indminglob3(i)
488         print *,'Direction ',i,' indices fin = ',indmaxglob2(i),indmaxglob3(i)
489        enddo
490#endif
491        endif
492
493#else
494        parentarray(:,1,1) = indminglob
495        parentarray(:,2,1) = indmaxglob
496        parentarray(:,1,2) = indminglob
497        parentarray(:,2,2) = indmaxglob
498
499
500        if (associated(agrif_external_mapping)) then
501          call Agrif_ChildGrid_to_ParentGrid()
502          call agrif_external_mapping(nbdim,child%root_var % posvar(1),child%root_var % posvar(2), &
503                                      parentarray,parentarray_chunk,correction_required,nb_chunks)
504          call Agrif_ParentGrid_to_ChildGrid()
505          allocate(decal_chunks(nb_chunks,nbdim))
506          do i=1,nb_chunks
507            decal_chunks(i,:)=parentarray_chunk(i,:,1,1)-parentarray_chunk(i,:,1,2)
508          enddo
509        else
510            nb_chunks=1
511            allocate(correction_required(nb_chunks))
512            correction_required=.FALSE.
513            allocate(parentarray_chunk(nb_chunks,nbdim,2,2))
514            parentarray_chunk(1,:,:,:)=parentarray
515        endif
516        if (agrif_debug_interp) then
517        print *,'AVANT PARENTCHILDBOUNDS'
518        print *,'nombre de chunks ',nb_chunks
519        do i=1,nb_chunks
520          print *,'CHUNK Number : ',i
521          do j=1,nbdim
522           print *,'Direction ',j
523           print *,'MIN MAX (2) = ',parentarray_chunk(i,j,1,2),parentarray_chunk(i,j,2,2)
524           print *,'MIN MAX (1) = ',parentarray_chunk(i,j,1,1),parentarray_chunk(i,j,2,1)
525          enddo
526        enddo
527        print *,'APRES PARENTCHILDBOUNDS'
528        endif
529        allocate(member_chuncks(nb_chunks))
530        allocate(parentarray_chunk_decal(nb_chunks,nbdim,2,2))
531        member_chuncks = .TRUE.
532        member = .TRUE.
533        do j=1,nb_chunks
534        if (agrif_debug_interp) print *,'CHUNK = ',j
535        if (member_chuncks(j)) then
536            do i=1,nbdim
537            parentarray_chunk_decal(j,i,:,1) = parentarray_chunk(j,i,:,1)   !+decal_chunks(j,i)
538            if (agrif_debug_interp) then
539                print *,'ENCORE = ',parentarray_chunk(j,i,:,1),parentarray_chunk_decal(j,i,:,1)
540            endif
541            enddo
542            if (agrif_debug_interp) then
543            do i=1,nbdim
544            print *,'parentarray = ',i,parentarray_chunk(j,i,1,1),parentarray_chunk(j,i,2,1), &
545                parentarray_chunk(j,i,1,2),parentarray_chunk(j,i,2,2)
546                print *,'parentarraydecal = ',i,parentarray_chunk_decal(j,i,1,1),parentarray_chunk_decal(j,i,2,1)
547            enddo
548            endif
549        endif
550        enddo
551
552
553        indmin = indminglob
554        indmax = indmaxglob
555        member = .TRUE.
556#endif
557!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
558! Correct for non refined directions
559!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
560        do i=1,nbdim
561          if (coords(i) == 0) then
562             indmin(i) = indminglob(i)
563             indmax(i) = indmaxglob(i)
564             pttruetab(i) = indminglob(i)
565             cetruetab(i) = indmaxglob(i)
566          endif
567        enddo
568
569    else
570
571#if defined AGRIF_MPI
572        s_Parent_temp = s_Parent + (indmin - pttab_Parent) * ds_Parent
573        s_Child_temp  = s_Child + (pttruetab - pttab_Child) * ds_Child
574#else
575        parentarray(:,1,1) = indminglob
576        parentarray(:,2,1) = indmaxglob
577        parentarray(:,1,2) = indminglob
578        parentarray(:,2,2) = indmaxglob
579        indmin = indminglob
580        indmax = indmaxglob
581        member = .TRUE.
582        s_Parent_temp = s_Parent + (indminglob - pttab_Parent) * ds_Parent
583        s_Child_temp  = s_Child + (pttab - pttab_Child) * ds_Child
584#endif
585
586    endif
587
588        if (agrif_debug_interp) then
589        print *,'************SPARENTCHILD*********************************'
590#ifdef AGRIF_MPI
591         print *,'Processeur ',Agrif_Procrank
592#endif
593        do i = 1 , nbdim
594        print *,'Direction ',i,' s_Parent_temp : ',s_Parent_temp(i),indmin(i)
595        print *,'Direction ',i,' s_Child_temp  : ',s_Child_temp(i),pttruetab(i)
596        enddo
597        print *,'*********************************************'
598        endif
599
600    call Agrif_get_var_bounds_array(parent,lowerbound,upperbound,nbdim)
601!
602    if (member) then
603        if (.not.associated(tempP)) allocate(tempP)
604!
605        call Agrif_array_allocate(tempP,parentarray(:,1,1),parentarray(:,2,1),nbdim)
606        call Agrif_var_set_array_tozero(tempP,nbdim)
607    endif
608    Agrif_CurChildgrid=>Agrif_Curgrid
609    call Agrif_ChildGrid_to_ParentGrid()
610
611
612    do i=1,nb_chunks
613    if (agrif_debug_interp) then
614    print *,'PROCNAME POUR CHUCNK ',i
615    endif
616   
617    if (member_chuncks(i)) then
618        select case (nbdim)
619        case(1)
620            ! call procname(tempP%array1,                         &
621            !           parentarray(1,1,2),parentarray(1,2,2),.TRUE.,nb,ndir)
622
623            call procname(tempP%array1(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1)),         &
624                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),.TRUE.,nb,ndir)
625
626        case(2)
627            ! call procname(tempP%array2,                         &
628                      ! parentarray(1,1,2),parentarray(1,2,2),    &
629                      ! parentarray(2,1,2),parentarray(2,2,2),.TRUE.,nb,ndir)
630
631            call procname(tempP%array2(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
632                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1)),         &
633                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
634                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),.TRUE.,nb,ndir)
635            if (agrif_debug_interp) print *,'SORTIE DE PROCNAME'
636            if (correction_required(i)) then
637             call correct_field(tempP%array2(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
638                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1)), &
639                      parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1), &
640                      parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1))
641            endif
642
643        case(3)
644            call procname(tempP%array3(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
645                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1), &
646                      parentarray_chunk_decal(i,3,1,1):parentarray_chunk_decal(i,3,2,1)),                 &
647                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
648                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
649                      parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),.TRUE.,nb,ndir)
650
651            if (agrif_debug_interp) then
652                print *,'CHUNK = ',i
653                print *,'NBNDIR = ',nb,ndir,correction_required(i)
654                print *,'TEMPARRAY3 INDEX LOCAUX PUIS GLOBAUX'
655                print *,parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1)
656                print *,parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1)
657                print *,parentarray_chunk_decal(i,3,1,1),parentarray_chunk_decal(i,3,2,1)
658                print *,parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2)
659                print *,parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2)
660                print *,parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2)
661                do j1=parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1)
662                    do i1=parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1)
663                        print *,'valprocname = ',i1,j1,tempP%array3(i1,j1,1)
664                    enddo
665                enddo
666            endif
667            if (correction_required(i)) then
668                do k=parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2)
669             call correct_field(tempP%array3(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
670                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1),k), &
671                      parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1), &
672                      parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1))
673                enddo
674            if (agrif_debug_interp) then
675                do j1=parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1)
676                    do i1=parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1)
677                        print *,'valprocname apres correction = ',i1,j1,tempP%array3(i1,j1,1)
678                    enddo
679                enddo
680            endif
681            endif
682
683            ! call procname(tempP%array3,                         &
684                      ! parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
685                      ! parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
686                      ! parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),.TRUE.,nb,ndir)
687        case(4)
688
689            call procname(tempP%array4(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
690                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1), &
691                      parentarray_chunk_decal(i,3,1,1):parentarray_chunk_decal(i,3,2,1), &
692                      parentarray_chunk_decal(i,4,1,1):parentarray_chunk_decal(i,4,2,1)),                 &
693                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
694                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
695                      parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),    &
696                      parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2),.TRUE.,nb,ndir)
697
698            if (correction_required(i)) then
699                do l=parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2)
700                do k=parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2)
701             call correct_field(tempP%array4(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
702                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1),k,l), &
703                      parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1), &
704                      parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1))
705                enddo
706            enddo
707            endif
708
709            ! call procname(tempP%array4,                         &
710            !           parentarray(1,1,2),parentarray(1,2,2),    &
711            !           parentarray(2,1,2),parentarray(2,2,2),    &
712            !           parentarray(3,1,2),parentarray(3,2,2),    &
713            !           parentarray(4,1,2),parentarray(4,2,2),.TRUE.,nb,ndir)
714        case(5)
715
716            call procname(tempP%array5(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
717                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1), &
718                      parentarray_chunk_decal(i,3,1,1):parentarray_chunk_decal(i,3,2,1), &
719                      parentarray_chunk_decal(i,4,1,1):parentarray_chunk_decal(i,4,2,1), &
720                      parentarray_chunk_decal(i,5,1,1):parentarray_chunk_decal(i,5,2,1)),                 &
721                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
722                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
723                      parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),    &
724                      parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2),    &
725                      parentarray_chunk(i,5,1,2),parentarray_chunk(i,5,2,2),.TRUE.,nb,ndir)
726
727            if (correction_required(i)) then
728                do m=parentarray_chunk(i,5,1,2),parentarray_chunk(i,5,2,2)
729                do l=parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2)
730                do k=parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2)
731             call correct_field(tempP%array5(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
732                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1),k,l,m), &
733                      parentarray_chunk_decal(i,1,1,1),parentarray_chunk_decal(i,1,2,1), &
734                      parentarray_chunk_decal(i,2,1,1),parentarray_chunk_decal(i,2,2,1))
735                enddo
736            enddo
737            enddo
738            endif
739
740            ! call procname(tempP%array5,                         &
741            !           parentarray(1,1,2),parentarray(1,2,2),    &
742            !           parentarray(2,1,2),parentarray(2,2,2),    &
743            !           parentarray(3,1,2),parentarray(3,2,2),    &
744            !           parentarray(4,1,2),parentarray(4,2,2),    &
745            !           parentarray(5,1,2),parentarray(5,2,2),.TRUE.,nb,ndir)
746        case(6)
747
748            call procname(tempP%array6(parentarray_chunk_decal(i,1,1,1):parentarray_chunk_decal(i,1,2,1), &
749                      parentarray_chunk_decal(i,2,1,1):parentarray_chunk_decal(i,2,2,1), &
750                      parentarray_chunk_decal(i,3,1,1):parentarray_chunk_decal(i,3,2,1), &
751                      parentarray_chunk_decal(i,4,1,1):parentarray_chunk_decal(i,4,2,1), &
752                      parentarray_chunk_decal(i,5,1,1):parentarray_chunk_decal(i,5,2,1), &
753                      parentarray_chunk_decal(i,6,1,1):parentarray_chunk_decal(i,6,2,1)),                 &
754                      parentarray_chunk(i,1,1,2),parentarray_chunk(i,1,2,2),    &
755                      parentarray_chunk(i,2,1,2),parentarray_chunk(i,2,2,2),    &
756                      parentarray_chunk(i,3,1,2),parentarray_chunk(i,3,2,2),    &
757                      parentarray_chunk(i,4,1,2),parentarray_chunk(i,4,2,2),    &
758                      parentarray_chunk(i,5,1,2),parentarray_chunk(i,5,2,2),    &
759                      parentarray_chunk(i,6,1,2),parentarray_chunk(i,6,2,2),.TRUE.,nb,ndir)
760
761            ! call procname(tempP%array6,                         &
762            !           parentarray(1,1,2),parentarray(1,2,2),    &
763            !           parentarray(2,1,2),parentarray(2,2,2),    &
764            !           parentarray(3,1,2),parentarray(3,2,2),    &
765            !           parentarray(4,1,2),parentarray(4,2,2),    &
766            !           parentarray(5,1,2),parentarray(5,2,2),    &
767            !           parentarray(6,1,2),parentarray(6,2,2),.TRUE.,nb,ndir)
768        end select
769!
770!
771    endif
772    enddo
773    call Agrif_ParentGrid_to_ChildGrid()
774    nullify(Agrif_CurChildgrid)
775
776#if defined AGRIF_MPI
777    if (.not.find_list_interp) then
778!
779        tab3(:,1) = indminglob2(:)
780        tab3(:,2) = indmaxglob2(:)
781        tab3(:,3) = indmin(:)
782        tab3(:,4) = indmax(:)
783        tab5(:,1) = indminglob3(:)
784        tab5(:,2) = indmaxglob3(:)
785        if (agrif_debug_interp) then
786         print *,'********************'
787         print *,'MPI VARIABLES'
788         print *,'INDMINGLOB2'
789         do i=1,nbdim
790            print *,'Direction ',i,indminglob2(i),indmaxglob2(i)
791         enddo
792         print *,'INDMIN'
793         do i=1,nbdim
794            print *,'Direction ',i,indmin(i),indmax(i)
795         enddo
796         print *,'INDMINGLOB3'
797         do i=1,nbdim
798            print *,'Direction ',i,indminglob3(i),indmaxglob3(i)
799         enddo
800        endif
801!
802        call MPI_ALLGATHER(tab3,4*nbdim,MPI_INTEGER,tab4,4*nbdim,MPI_INTEGER,Agrif_mpi_comm,code)
803        call MPI_ALLGATHER(tab5,2*nbdim,MPI_INTEGER,tab6,2*nbdim,MPI_INTEGER,Agrif_mpi_comm,code)
804        if (.not.associated(tempPextend))   allocate(tempPextend)
805
806        do k=0,Agrif_Nbprocs-1
807            do j=1,4
808                do i=1,nbdim
809                    tab4t(i,k,j) = tab4(i,j,k)
810                enddo
811            enddo
812        enddo
813
814        do k=0,Agrif_Nbprocs-1
815          do j=1,2
816            do i=1,nbdim
817               tab5t(i,k,j) = tab6(i,j,k)
818            enddo
819          enddo
820        enddo
821
822        memberin1(1) = memberin
823        call MPI_ALLGATHER(memberin1,1,MPI_LOGICAL,memberinall,1,MPI_LOGICAL,Agrif_mpi_comm,code)
824
825        call Get_External_Data_first(tab4t(:,:,1),tab4t(:,:,2),         &
826                                     tab4t(:,:,3),tab4t(:,:,4),         &
827                                     nbdim,memberinall, coords,         &
828                                     sendtoproc1,recvfromproc1,         &
829                                     tab4t(:,:,5),tab4t(:,:,6),         &
830                                     tab4t(:,:,7),tab4t(:,:,8),         &
831                                     tab5t(:,:,1),tab5t(:,:,2))
832    endif
833
834    call ExchangeSameLevel(sendtoproc1,recvfromproc1,nbdim,         &
835            tab4t(:,:,3),tab4t(:,:,4),tab4t(:,:,5),tab4t(:,:,6),    &
836            tab4t(:,:,7),tab4t(:,:,8),memberin,tempP,tempPextend)
837#else
838    tempPextend => tempP
839#endif
840
841    if (.not.find_list_interp) then
842        call Agrif_Addto_list_interp(                           &
843                child%list_interp,pttab,petab,                  &
844                pttab_Child,pttab_Parent,indmin,indmax,         &
845                indminglob,indmaxglob,                          &
846                pttruetab,cetruetab,                            &
847                memberin,nbdim,                                 &
848                parentarray_chunk,parentarray_chunk_decal,decal_chunks,&
849                correction_required,member_chuncks,nb_chunks    &
850#if defined AGRIF_MPI
851               ,indminglob2,indmaxglob2,                        &
852                parentarray,                                    &
853                member,                                         &
854                tab4t,memberinall,sendtoproc1,recvfromproc1     &
855#endif
856        )
857    endif
858!
859    if (memberin) then
860!
861        if (.not.associated(tempC)) allocate(tempC)
862!
863        call Agrif_array_allocate(tempC,pttruetab,cetruetab,nbdim)
864!
865!       Special values on the parent grid
866        if (Agrif_UseSpecialValue) then
867!
868            noraftab(1:nbdim) = child % root_var % interptab(1:nbdim) == 'N'
869!
870            if (.not.associated(parentvalues))  allocate(parentvalues)
871!
872            call Agrif_array_allocate(parentvalues,indmin,indmax,nbdim)
873            call Agrif_var_full_copy_array(parentvalues,tempPextend,nbdim)
874!
875            call Agrif_CheckMasknD(tempPextend,parentvalues,    &
876                    indmin(1:nbdim),indmax(1:nbdim),            &
877                    indmin(1:nbdim),indmax(1:nbdim),            &
878                    noraftab(1:nbdim),nbdim)
879!
880            call Agrif_array_deallocate(parentvalues,nbdim)
881!
882        endif
883!
884!       Interpolation of the current grid
885!
886        if ( memberin ) then
887            select case(nbdim)
888            case(1)
889                call Agrif_Interp_1D_recursive( type_interp(1),                         &
890                                                tempPextend%array1,                     &
891                                                tempC%array1,                           &
892                                                indmin(1), indmax(1),                   &
893                                                pttruetab(1),    cetruetab(1),          &
894                                                s_Child_temp(1), s_Parent_temp(1),      &
895                                                ds_Child(1),     ds_Parent(1) )
896            case(2)
897                call Agrif_Interp_2D_recursive( type_interp(1:2),                       &
898                                                tempPextend % array2,                   &
899                                                tempC       % array2,                   &
900                                                indmin(1:2), indmax(1:2),               &
901                                                pttruetab(1:2),    cetruetab(1:2),      &
902                                                s_Child_temp(1:2), s_Parent_temp(1:2),  &
903                                                ds_Child(1:2),    ds_Parent(1:2) )
904            case(3)
905                        if (agrif_debug_interp) then
906        print *,'APRES ECHANGE'
907        print *,'nombre de chunks ',nb_chunks
908        print *,'indmin = ',indmin
909        print *,'indmax = ',indmax
910        do i=1,nb_chunks
911          print *,'CHUNK Number : ',i
912          print *,'MEMBER = ',member_chuncks(i)
913          print *,'Correction = ',correction_required(i)
914        enddo
915        endif
916
917                if (agrif_debug_interp) then
918                !    if ((nb==1).AND.(ndir==1)) then
919                        print *,'valeur parent = '
920                        do j=indmin(2),indmax(2)
921                            do i=indmin(1),indmax(1)
922                                print *,'par = ',i,j,tempPextend%array3(i,j,1)
923                            enddo
924                        enddo
925
926                 !   endif
927                endif
928                call Agrif_Interp_3D_recursive( type_interp(1:3),                       &
929                                                tempPextend % array3,                   &
930                                                tempC       % array3,                   &
931                                                indmin(1:3), indmax(1:3),               &
932                                                pttruetab(1:3),    cetruetab(1:3),      &
933                                                s_Child_temp(1:3), s_Parent_temp(1:3),  &
934                                                ds_Child(1:3),    ds_Parent(1:3) )
935                if (agrif_debug_interp) then
936                  !  if ((nb==1).AND.(ndir==1)) then
937                        print *,'valeur enfnat = '
938                        do j=pttruetab(2),cetruetab(2)
939                            do i=pttruetab(1),cetruetab(1)
940                                print *,'par = ',i,j,tempC%array3(i,j,1)
941                            enddo
942                        enddo
943
944                  !  endif
945                endif
946            case(4)
947                call Agrif_Interp_4D_recursive( type_interp(1:4),                       &
948                                                tempPextend % array4,                   &
949                                                tempC       % array4,                   &
950                                                indmin(1:4), indmax(1:4),               &
951                                                pttruetab(1:4),    cetruetab(1:4),      &
952                                                s_Child_temp(1:4), s_Parent_temp(1:4),  &
953                                                ds_Child(1:4),    ds_Parent(1:4) )
954            case(5)
955                call Agrif_Interp_5D_recursive( type_interp(1:5),                       &
956                                                tempPextend % array5,                   &
957                                                tempC       % array5,                   &
958                                                indmin(1:5), indmax(1:5),               &
959                                                pttruetab(1:5),    cetruetab(1:5),      &
960                                                s_Child_temp(1:5), s_Parent_temp(1:5),  &
961                                                ds_Child(1:5),    ds_Parent(1:5) )
962            case(6)
963                call Agrif_Interp_6D_recursive( type_interp(1:6),                       &
964                                                tempPextend % array6,                   &
965                                                tempC       % array6,                   &
966                                                indmin(1:6), indmax(1:6),               &
967                                                pttruetab(1:6),    cetruetab(1:6),      &
968                                                s_Child_temp(1:6), s_Parent_temp(1:6),  &
969                                                ds_Child(1:6),    ds_Parent(1:6) )
970            end select
971!
972            call Agrif_get_var_bounds_array(child,lowerbound,upperbound,nbdim,parent)
973
974#if defined AGRIF_MPI
975            call Agrif_GlobalToLocalBounds(childarray, lowerbound, upperbound,  &
976                                            pttruetab, cetruetab, coords,       &
977                                            nbdim, Agrif_Procrank, memberout)
978#else
979            childarray(:,1,1) = pttruetab
980            childarray(:,2,1) = cetruetab
981            childarray(:,1,2) = pttruetab
982            childarray(:,2,2) = cetruetab
983!cccccccccccccc       memberout = .TRUE.
984#endif
985!
986!           Special values on the child grid
987            if (Agrif_UseSpecialValueFineGrid) then
988                call GiveAgrif_SpecialValueToTab_mpi( child, tempC, childarray, Agrif_SpecialValueFineGrid,nbdim )
989            endif
990!
991        endif   ! ( memberin )
992!
993        if (torestore) then
994!
995#if defined AGRIF_MPI
996!
997            SELECT CASE (nbdim)
998            CASE (1)
999                do i = pttruetab(1),cetruetab(1)
1000!hildarrayAModifier     if (restore%restore1D(i) == 0)      &
1001!hildarrayAModifier         child%array1(childarray(i,1,2)) = tempC%array1(i)
1002                enddo
1003            CASE (2)
1004                do i = pttruetab(1),cetruetab(1)
1005                do j = pttruetab(2),cetruetab(2)
1006!hildarrayAModifier     if (restore%restore2D(i,j) == 0)    &
1007!hildarrayAModifier         child%array2(childarray(i,1,2), &
1008!hildarrayAModifier                          childarray(j,2,2)) = tempC%array2(i,j)
1009                enddo
1010                enddo
1011            CASE (3)
1012                do i = pttruetab(1),cetruetab(1)
1013                do j = pttruetab(2),cetruetab(2)
1014                do k = pttruetab(3),cetruetab(3)
1015!hildarrayAModifier     if (restore%restore3D(i,j,k) == 0)  &
1016!hildarrayAModifier         child%array3(childarray(i,1,2), &
1017!hildarrayAModifier                          childarray(j,2,2), &
1018!hildarrayAModifier                          childarray(k,3,2)) = tempC%array3(i,j,k)
1019                enddo
1020                enddo
1021                enddo
1022            CASE (4)
1023                do i = pttruetab(1),cetruetab(1)
1024                do j = pttruetab(2),cetruetab(2)
1025                do k = pttruetab(3),cetruetab(3)
1026                do l = pttruetab(4),cetruetab(4)
1027!hildarrayAModifier     if (restore%restore4D(i,j,k,l) == 0)    &
1028!hildarrayAModifier         child%array4(childarray(i,1,2),     &
1029!hildarrayAModifier                          childarray(j,2,2),     &
1030!hildarrayAModifier                          childarray(k,3,2),     &
1031!hildarrayAModifier                          childarray(l,4,2)) = tempC%array4(i,j,k,l)
1032                enddo
1033                enddo
1034                enddo
1035                enddo
1036            CASE (5)
1037                do i = pttruetab(1),cetruetab(1)
1038                do j = pttruetab(2),cetruetab(2)
1039                do k = pttruetab(3),cetruetab(3)
1040                do l = pttruetab(4),cetruetab(4)
1041                do m = pttruetab(5),cetruetab(5)
1042!hildarrayAModifier     if (restore%restore5D(i,j,k,l,m) == 0)  &
1043!hildarrayAModifier         child%array5(childarray(i,1,2),     &
1044!hildarrayAModifier                          childarray(j,2,2),     &
1045!hildarrayAModifier                          childarray(k,3,2),     &
1046!hildarrayAModifier                          childarray(l,4,2),     &
1047!hildarrayAModifier                          childarray(m,5,2)) = tempC%array5(i,j,k,l,m)
1048                enddo
1049                enddo
1050                enddo
1051                enddo
1052                enddo
1053            CASE (6)
1054                do i = pttruetab(1),cetruetab(1)
1055                do j = pttruetab(2),cetruetab(2)
1056                do k = pttruetab(3),cetruetab(3)
1057                do l = pttruetab(4),cetruetab(4)
1058                do m = pttruetab(5),cetruetab(5)
1059                do n = pttruetab(6),cetruetab(6)
1060!hildarrayAModifier     if (restore%restore6D(i,j,k,l,m,n) == 0)    &
1061!hildarrayAModifier         child%array6(childarray(i,1,2),         &
1062!hildarrayAModifier                          childarray(j,2,2),         &
1063!hildarrayAModifier                          childarray(k,3,2),         &
1064!hildarrayAModifier                          childarray(l,4,2),         &
1065!hildarrayAModifier                          childarray(m,5,2),         &
1066!hildarrayAModifier                          childarray(n,6,2)) = tempC%array6(i,j,k,l,m,n)
1067                enddo
1068                enddo
1069                enddo
1070                enddo
1071                enddo
1072                enddo
1073            END SELECT
1074!
1075#else
1076            select case (nbdim)
1077            case (1)
1078                do i = pttruetab(1),cetruetab(1)
1079                    if (restore%restore1D(i) == 0)          &
1080                        parray1(i) = tempC % array1(i)
1081                enddo
1082            case (2)
1083                do j = pttruetab(2),cetruetab(2)
1084                do i = pttruetab(1),cetruetab(1)
1085                    if (restore%restore2D(i,j) == 0)        &
1086                        parray2(i,j) = tempC % array2(i,j)
1087                enddo
1088                enddo
1089            case (3)
1090                do k = pttruetab(3),cetruetab(3)
1091                do j = pttruetab(2),cetruetab(2)
1092                do i = pttruetab(1),cetruetab(1)
1093                    if (restore%restore3D(i,j,k) == 0)      &
1094                        parray3(i,j,k) = tempC % array3(i,j,k)
1095                enddo
1096                enddo
1097                enddo
1098            case (4)
1099                do l = pttruetab(4),cetruetab(4)
1100                do k = pttruetab(3),cetruetab(3)
1101                do j = pttruetab(2),cetruetab(2)
1102                do i = pttruetab(1),cetruetab(1)
1103                    if (restore%restore4D(i,j,k,l) == 0)    &
1104                        parray4(i,j,k,l) = tempC % array4(i,j,k,l)
1105                enddo
1106                enddo
1107                enddo
1108                enddo
1109            case (5)
1110                do m = pttruetab(5),cetruetab(5)
1111                do l = pttruetab(4),cetruetab(4)
1112                do k = pttruetab(3),cetruetab(3)
1113                do j = pttruetab(2),cetruetab(2)
1114                do i = pttruetab(1),cetruetab(1)
1115                    if (restore%restore5D(i,j,k,l,m) == 0)  &
1116                        parray5(i,j,k,l,m) = tempC % array5(i,j,k,l,m)
1117                enddo
1118                enddo
1119                enddo
1120                enddo
1121                enddo
1122            case (6)
1123                do n = pttruetab(6),cetruetab(6)
1124                do m = pttruetab(5),cetruetab(5)
1125                do l = pttruetab(4),cetruetab(4)
1126                do k = pttruetab(3),cetruetab(3)
1127                do j = pttruetab(2),cetruetab(2)
1128                do i = pttruetab(1),cetruetab(1)
1129                    if (restore%restore6D(i,j,k,l,m,n) == 0)    &
1130                        parray6(i,j,k,l,m,n) = tempC % array6(i,j,k,l,m,n)
1131                enddo
1132                enddo
1133                enddo
1134                enddo
1135                enddo
1136                enddo
1137            end select
1138!
1139#endif
1140!
1141        else    ! .not.to_restore
1142!
1143            if (memberin) then
1144    !
1145                if ( .not.in_bc ) then
1146                    select case(nbdim)
1147                    case(1)
1148                        call procname(tempC % array1(            &
1149                                childarray(1,1,1):childarray(1,2,1)), &
1150                                childarray(1,1,2),childarray(1,2,2),.FALSE.,nb,ndir)
1151                    case(2)
1152                        call procname(                                &
1153                                tempC % array2(            &
1154                                childarray(1,1,1):childarray(1,2,1),  &
1155                                childarray(2,1,1):childarray(2,2,1)), &
1156                                childarray(1,1,2),childarray(1,2,2),  &
1157                                childarray(2,1,2),childarray(2,2,2),.FALSE.,nb,ndir)
1158                    case(3)
1159                        call procname(                                &
1160                                tempC % array3(            &
1161                                childarray(1,1,1):childarray(1,2,1),  &
1162                                childarray(2,1,1):childarray(2,2,1),  &
1163                                childarray(3,1,1):childarray(3,2,1)), &
1164                                childarray(1,1,2),childarray(1,2,2),  &
1165                                childarray(2,1,2),childarray(2,2,2),  &
1166                                childarray(3,1,2),childarray(3,2,2),.FALSE.,nb,ndir)
1167                    case(4)
1168                        call procname(                                &
1169                                tempC % array4(            &
1170                                childarray(1,1,1):childarray(1,2,1),  &
1171                                childarray(2,1,1):childarray(2,2,1),  &
1172                                childarray(3,1,1):childarray(3,2,1),  &
1173                                childarray(4,1,1):childarray(4,2,1)), &
1174                                childarray(1,1,2),childarray(1,2,2),  &
1175                                childarray(2,1,2),childarray(2,2,2),  &
1176                                childarray(3,1,2),childarray(3,2,2),  &
1177                                childarray(4,1,2),childarray(4,2,2),.FALSE.,nb,ndir)
1178                    case(5)
1179                        call procname(                                &
1180                                tempC % array5(            &
1181                                childarray(1,1,1):childarray(1,2,1),  &
1182                                childarray(2,1,1):childarray(2,2,1),  &
1183                                childarray(3,1,1):childarray(3,2,1),  &
1184                                childarray(4,1,1):childarray(4,2,1),  &
1185                                childarray(5,1,1):childarray(5,2,1)), &
1186                                childarray(1,1,2),childarray(1,2,2),  &
1187                                childarray(2,1,2),childarray(2,2,2),  &
1188                                childarray(3,1,2),childarray(3,2,2),  &
1189                                childarray(4,1,2),childarray(4,2,2),  &
1190                                childarray(5,1,2),childarray(5,2,2),.FALSE.,nb,ndir)
1191                    case(6)
1192                        call procname(                                &
1193                                tempC % array6(            &
1194                                childarray(1,1,1):childarray(1,2,1),  &
1195                                childarray(2,1,1):childarray(2,2,1),  &
1196                                childarray(3,1,1):childarray(3,2,1),  &
1197                                childarray(4,1,1):childarray(4,2,1),  &
1198                                childarray(5,1,1):childarray(5,2,1),  &
1199                                childarray(6,1,1):childarray(6,2,1)), &
1200                                childarray(1,1,2),childarray(1,2,2),  &
1201                                childarray(2,1,2),childarray(2,2,2),  &
1202                                childarray(3,1,2),childarray(3,2,2),  &
1203                                childarray(4,1,2),childarray(4,2,2),  &
1204                                childarray(5,1,2),childarray(5,2,2),  &
1205                                childarray(6,1,2),childarray(6,2,2),.FALSE.,nb,ndir)
1206                    end select
1207                else    ! we are in_bc
1208                    select case (nbdim)
1209                    case (1)
1210                        parray1(childarray(1,1,2):childarray(1,2,2)) =    &
1211                         tempC%array1(childarray(1,1,1):childarray(1,2,1))
1212                    case (2)
1213                        parray2(childarray(1,1,2):childarray(1,2,2),      &
1214                                      childarray(2,1,2):childarray(2,2,2)) =    &
1215                         tempC%array2(childarray(1,1,1):childarray(1,2,1),      &
1216                                      childarray(2,1,1):childarray(2,2,1))
1217                    case (3)
1218                        parray3(childarray(1,1,2):childarray(1,2,2),      &
1219                                      childarray(2,1,2):childarray(2,2,2),      &
1220                                      childarray(3,1,2):childarray(3,2,2)) =    &
1221                         tempC%array3(childarray(1,1,1):childarray(1,2,1),      &
1222                                      childarray(2,1,1):childarray(2,2,1),      &
1223                                      childarray(3,1,1):childarray(3,2,1))
1224                if (agrif_debug_interp) then
1225                    if ((nb==1).AND.(ndir==1)) then
1226                        print *,'valeur enfnat2 = '
1227                        do j=childarray(2,1,2),childarray(2,2,2)
1228                            do i=childarray(1,1,2),childarray(1,2,2)
1229                                print *,'par = ',i,j,parray3(i,j,1)
1230                            enddo
1231                        enddo
1232
1233                    endif
1234                endif
1235                    case (4)
1236                        parray4(childarray(1,1,2):childarray(1,2,2),      &
1237                                      childarray(2,1,2):childarray(2,2,2),      &
1238                                      childarray(3,1,2):childarray(3,2,2),      &
1239                                      childarray(4,1,2):childarray(4,2,2)) =    &
1240                         tempC%array4(childarray(1,1,1):childarray(1,2,1),      &
1241                                      childarray(2,1,1):childarray(2,2,1),      &
1242                                      childarray(3,1,1):childarray(3,2,1),      &
1243                                      childarray(4,1,1):childarray(4,2,1))
1244                    case (5)
1245                        parray5(childarray(1,1,2):childarray(1,2,2),      &
1246                                      childarray(2,1,2):childarray(2,2,2),      &
1247                                      childarray(3,1,2):childarray(3,2,2),      &
1248                                      childarray(4,1,2):childarray(4,2,2),      &
1249                                      childarray(5,1,2):childarray(5,2,2)) =    &
1250                         tempC%array5(childarray(1,1,1):childarray(1,2,1),      &
1251                                      childarray(2,1,1):childarray(2,2,1),      &
1252                                      childarray(3,1,1):childarray(3,2,1),      &
1253                                      childarray(4,1,1):childarray(4,2,1),      &
1254                                      childarray(5,1,1):childarray(5,2,1))
1255                    case (6)
1256                        parray6(childarray(1,1,2):childarray(1,2,2),      &
1257                                      childarray(2,1,2):childarray(2,2,2),      &
1258                                      childarray(3,1,2):childarray(3,2,2),      &
1259                                      childarray(4,1,2):childarray(4,2,2),      &
1260                                      childarray(5,1,2):childarray(5,2,2),      &
1261                                      childarray(6,1,2):childarray(6,2,2)) =    &
1262                         tempC%array6(childarray(1,1,1):childarray(1,2,1),      &
1263                                      childarray(2,1,1):childarray(2,2,1),      &
1264                                      childarray(3,1,1):childarray(3,2,1),      &
1265                                      childarray(4,1,1):childarray(4,2,1),      &
1266                                      childarray(5,1,1):childarray(5,2,1),      &
1267                                      childarray(6,1,1):childarray(6,2,1))
1268                    end select
1269                endif  ! < (.not.in_bc)
1270            endif  ! < memberin
1271!
1272        endif
1273
1274        call Agrif_array_deallocate(tempPextend,nbdim)
1275        call Agrif_array_deallocate(tempC,nbdim)
1276
1277    endif
1278!
1279!   Deallocations
1280#if defined AGRIF_MPI
1281    if (member) then
1282    if (agrif_debug_interp) then
1283    print *,'ALLCOATED 0 = ',allocated(tempP%array3),size(tempP%array3)
1284    endif
1285        call Agrif_array_deallocate(tempP,nbdim)
1286    endif
1287#endif
1288!---------------------------------------------------------------------------------------------------
1289end subroutine Agrif_InterpnD
1290!===================================================================================================
1291!
1292!===================================================================================================
1293!  subroutine Agrif_Parentbounds
1294!
1295!> Calculates the bounds of the parent grid for the interpolation of the child grid
1296!---------------------------------------------------------------------------------------------------
1297subroutine Agrif_Parentbounds ( type_interp, nbdim, indmin, indmax, &
1298                                s_Parent_temp, s_Child_temp,        &
1299                                s_Child, ds_Child,                  &
1300                                s_Parent,ds_Parent,                 &
1301                                pttruetab, cetruetab,               &
1302                                pttab_Child, pttab_Parent, posvar, coords )
1303!---------------------------------------------------------------------------------------------------
1304    INTEGER, DIMENSION(6),     intent(in)  :: type_interp
1305    INTEGER,                   intent(in)  :: nbdim
1306    INTEGER, DIMENSION(nbdim), intent(out) :: indmin, indmax
1307    REAL(kind=8),    DIMENSION(nbdim), intent(out) :: s_Parent_temp, s_child_temp
1308    REAL(kind=8),    DIMENSION(nbdim), intent(in)  :: s_Child, ds_child
1309    REAL(kind=8),    DIMENSION(nbdim), intent(in)  :: s_Parent,ds_Parent
1310    INTEGER, DIMENSION(nbdim), intent(in)  :: pttruetab, cetruetab
1311    INTEGER, DIMENSION(nbdim), intent(in)  :: pttab_Child, pttab_Parent
1312    INTEGER, DIMENSION(nbdim), intent(in)  :: posvar
1313    INTEGER, DIMENSION(nbdim), intent(in)  :: coords
1314!
1315    INTEGER :: i
1316    REAL(kind=8),DIMENSION(nbdim) :: dim_newmin, dim_newmax
1317!
1318    dim_newmin = s_Child + (pttruetab - pttab_Child) * ds_Child
1319    dim_newmax = s_Child + (cetruetab - pttab_Child) * ds_Child
1320!
1321    do i = 1,nbdim
1322!
1323        indmin(i) = pttab_Parent(i) + agrif_int((dim_newmin(i)-s_Parent(i))/ds_Parent(i))
1324        indmax(i) = pttab_Parent(i) + agrif_ceiling((dim_newmax(i)-s_Parent(i))/ds_Parent(i))
1325!
1326!       Necessary for the Quadratic interpolation
1327!
1328        if ( (pttruetab(i) == cetruetab(i)) .and. (posvar(i) == 1) ) then
1329        elseif ( coords(i) == 0 ) then  ! (interptab == 'N')
1330        elseif ( (type_interp(i) == Agrif_ppm)     .or.     &
1331                 (type_interp(i) == Agrif_eno)     .or.     &
1332                 (type_interp(i) == Agrif_ppm_lim) .or.     &
1333                 (type_interp(i) == Agrif_weno) ) then
1334            indmin(i) = indmin(i) - 2
1335            indmax(i) = indmax(i) + 2
1336
1337            if (Agrif_UseSpecialValue) then
1338               indmin(i) = indmin(i)-MaxSearch
1339               indmax(i) = indmax(i)+MaxSearch
1340            endif
1341
1342        elseif ( (type_interp(i) /= Agrif_constant) .and.   &
1343                 (type_interp(i) /= Agrif_linear) ) then
1344            indmin(i) = indmin(i) - 1
1345            indmax(i) = indmax(i) + 1
1346
1347            if (Agrif_UseSpecialValue) then
1348               indmin(i) = indmin(i)-MaxSearch
1349               indmax(i) = indmax(i)+MaxSearch
1350            endif
1351
1352        elseif ( (type_interp(i) == Agrif_constant) .or.   &
1353                 (type_interp(i) == Agrif_linear) ) then
1354            if (Agrif_UseSpecialValue) then
1355               indmin(i) = indmin(i)-MaxSearch
1356               indmax(i) = indmax(i)+MaxSearch
1357            endif
1358
1359        endif
1360!
1361    enddo
1362!
1363    s_Parent_temp = s_Parent + (indmin - pttab_Parent) * ds_Parent
1364    s_Child_temp  = s_Child + (pttruetab - pttab_Child) * ds_Child
1365!---------------------------------------------------------------------------------------------------
1366end subroutine Agrif_Parentbounds
1367!===================================================================================================
1368!
1369!===================================================================================================
1370!  subroutine Agrif_Interp_1D_Recursive
1371!
1372!> Subroutine for the interpolation of a 1D grid variable.
1373!> It calls Agrif_InterpBase.
1374!---------------------------------------------------------------------------------------------------
1375subroutine Agrif_Interp_1D_recursive ( type_interp, tabin, tabout,  &
1376                                       indmin, indmax,              &
1377                                       pttab_child, petab_child,    &
1378                                       s_child,  s_parent,          &
1379                                       ds_child, ds_parent )
1380!---------------------------------------------------------------------------------------------------
1381    integer,            intent(in)  :: type_interp
1382    integer,            intent(in)  :: indmin, indmax
1383    integer,            intent(in)  :: pttab_child, petab_child
1384    real(kind=8),               intent(in)  :: s_child, s_parent
1385    real(kind=8),               intent(in)  :: ds_child, ds_parent
1386    real, dimension(            &
1387        indmin:indmax           &
1388    ),                  intent(in)  :: tabin
1389    real, dimension(            &
1390        pttab_child:petab_child &
1391    ),                  intent(out) :: tabout
1392!---------------------------------------------------------------------------------------------------
1393    call Agrif_InterpBase(type_interp,                      &
1394                          tabin(indmin:indmax),             &
1395                          tabout(pttab_child:petab_child),  &
1396                          indmin, indmax,                   &
1397                          pttab_child, petab_child,         &
1398                          s_parent,    s_child,             &
1399                         ds_parent,   ds_child)
1400!---------------------------------------------------------------------------------------------------
1401end subroutine Agrif_Interp_1D_recursive
1402!===================================================================================================
1403!
1404!===================================================================================================
1405!  subroutine Agrif_Interp_2D_Recursive
1406!
1407!> Subroutine for the interpolation of a 2D grid variable.
1408!> It calls Agrif_Interp_1D_recursive and Agrif_InterpBase.
1409!---------------------------------------------------------------------------------------------------
1410subroutine Agrif_Interp_2D_recursive ( type_interp, tabin, tabout,  &
1411                                       indmin, indmax,              &
1412                                       pttab_child, petab_child,    &
1413                                       s_child,  s_parent,          &
1414                                       ds_child, ds_parent )
1415!---------------------------------------------------------------------------------------------------
1416    integer, dimension(2),              intent(in)  :: type_interp
1417    integer, dimension(2),              intent(in)  :: indmin, indmax
1418    integer, dimension(2),              intent(in)  :: pttab_child, petab_child
1419    real(kind=8),    dimension(2),              intent(in)  :: s_child, s_parent
1420    real(kind=8),    dimension(2),              intent(in)  :: ds_child, ds_parent
1421    real,    dimension(                 &
1422        indmin(1):indmax(1),            &
1423        indmin(2):indmax(2)),           intent(in)  :: tabin
1424    real,    dimension(                 &
1425        pttab_child(1):petab_child(1),  &
1426        pttab_child(2):petab_child(2)), intent(out) :: tabout
1427!---------------------------------------------------------------------------------------------------
1428    real, dimension(                    &
1429        pttab_child(1):petab_child(1),  &
1430        indmin(2):indmax(2))            :: tabtemp
1431    real, dimension(                    &
1432        pttab_child(2):petab_child(2),  &
1433        pttab_child(1):petab_child(1))  :: tabout_trsp
1434    real, dimension(                    &
1435        indmin(2):indmax(2),            &
1436        pttab_child(1):petab_child(1))  :: tabtemp_trsp
1437    integer                             :: i, j, coeffraf
1438!---------------------------------------------------------------------------------------------------
1439!
1440    coeffraf = nint ( ds_parent(1) / ds_child(1) )
1441!
1442    if ((type_interp(1) == Agrif_Linear) .and. (coeffraf /= 1)) then
1443!---CDIR NEXPAND
1444        if(.NOT. precomputedone(1))     &
1445            call Linear1dPrecompute2d(                  &
1446                    indmax(2)-indmin(2)+1,              &
1447                    indmax(1)-indmin(1)+1,              &
1448                    petab_child(1)-pttab_child(1)+1,    &
1449                    s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
1450!---CDIR NEXPAND
1451        call Linear1dAfterCompute(tabin,tabtemp,size(tabin),size(tabtemp),1)
1452!
1453    elseif ((type_interp(1) == Agrif_PPM) .and. (coeffraf /= 1)) then
1454!---CDIR NEXPAND
1455        if(.NOT. precomputedone(1))     &
1456            call PPM1dPrecompute2d(                     &
1457                    indmax(2)-indmin(2)+1,              &
1458                    indmax(1)-indmin(1)+1,              &
1459                    petab_child(1)-pttab_child(1)+1,        &
1460                    s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
1461!---CDIR NEXPAND
1462        call PPM1dAfterCompute(tabin,tabtemp,size(tabin),size(tabtemp),1)
1463    else
1464        do j = indmin(2),indmax(2)
1465!
1466!---CDIR NEXPAND
1467            call Agrif_Interp_1D_recursive(type_interp(1),                  &
1468                    tabin(indmin(1):indmax(1),j),               &
1469                    tabtemp(pttab_child(1):petab_child(1),j),   &
1470                    indmin(1),indmax(1),                    &
1471                    pttab_child(1),petab_child(1),          &
1472                    s_child(1), s_parent(1),                &
1473                    ds_child(1),ds_parent(1))
1474!
1475        enddo
1476    endif
1477
1478    coeffraf = nint(ds_parent(2)/ds_child(2))
1479    tabtemp_trsp = TRANSPOSE(tabtemp)
1480
1481    if ((type_interp(2) == Agrif_Linear) .and. (coeffraf /= 1)) then
1482!---CDIR NEXPAND
1483        if(.NOT. precomputedone(2))     &
1484            call Linear1dPrecompute2d(                  &
1485                    petab_child(1)-pttab_child(1)+1,    &
1486                    indmax(2)-indmin(2)+1,              &
1487                    petab_child(2)-pttab_child(2)+1,    &
1488                    s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
1489!---CDIR NEXPAND
1490        call Linear1dAfterCompute(tabtemp_trsp,tabout_trsp, &
1491                size(tabtemp_trsp),size(tabout_trsp),2)
1492
1493    elseif ((type_interp(2) == Agrif_PPM) .and. (coeffraf /= 1)) then
1494!---CDIR NEXPAND
1495        if(.NOT. precomputedone(2))     &
1496            call PPM1dPrecompute2d(                     &
1497                    petab_child(1)-pttab_child(1)+1,    &
1498                    indmax(2)-indmin(2)+1,              &
1499                    petab_child(2)-pttab_child(2)+1,    &
1500                    s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
1501!---CDIR NEXPAND
1502        call PPM1dAfterCompute(tabtemp_trsp, tabout_trsp,    &
1503                               size(tabtemp_trsp), size(tabout_trsp), 2)
1504    else
1505        do i = pttab_child(1), petab_child(1)
1506!
1507!---CDIR NEXPAND
1508            call Agrif_InterpBase(type_interp(2),                                   &
1509                                  tabtemp_trsp(indmin(2):indmax(2), i),             &
1510                                  tabout_trsp(pttab_child(2):petab_child(2), i),    &
1511                                  indmin(2), indmax(2),                             &
1512                                  pttab_child(2), petab_child(2),                   &
1513                                  s_parent(2),    s_child(2),                       &
1514                                 ds_parent(2),   ds_child(2) )
1515        enddo
1516    endif
1517!
1518    tabout = TRANSPOSE(tabout_trsp)
1519!---------------------------------------------------------------------------------------------------
1520end subroutine Agrif_Interp_2D_recursive
1521!===================================================================================================
1522!
1523!===================================================================================================
1524!  subroutine Agrif_Interp_3D_Recursive
1525!
1526!> Subroutine for the interpolation of a 3D grid variable.
1527!> It calls #Agrif_Interp_2D_recursive and #Agrif_InterpBase.
1528!---------------------------------------------------------------------------------------------------
1529subroutine Agrif_Interp_3D_recursive ( type_interp, tabin, tabout,  &
1530                                       indmin, indmax,              &
1531                                       pttab_child, petab_child,    &
1532                                       s_child,  s_parent,          &
1533                                       ds_child, ds_parent )
1534!---------------------------------------------------------------------------------------------------
1535    integer, dimension(3),              intent(in)  :: type_interp
1536    integer, dimension(3),              intent(in)  :: indmin, indmax
1537    integer, dimension(3),              intent(in)  :: pttab_child, petab_child
1538    real(kind=8),    dimension(3),              intent(in)  :: s_child, s_parent
1539    real(kind=8),    dimension(3),              intent(in)  :: ds_child, ds_parent
1540    real,    dimension(                 &
1541        indmin(1):indmax(1),            &
1542        indmin(2):indmax(2),            &
1543        indmin(3):indmax(3)),           intent(in)  :: tabin
1544    real,    dimension(                 &
1545        pttab_child(1):petab_child(1),  &
1546        pttab_child(2):petab_child(2),  &
1547        pttab_child(3):petab_child(3)), intent(out) :: tabout
1548!---------------------------------------------------------------------------------------------------
1549    real, dimension(                    &
1550        pttab_child(1):petab_child(1),  &
1551        pttab_child(2):petab_child(2),  &
1552        indmin(3):indmax(3))            :: tabtemp
1553    integer                             :: i, j, k, coeffraf
1554    integer                             :: locind_child_left, kdeb
1555!
1556    coeffraf = nint ( ds_parent(1) / ds_child(1) )
1557    if ( (type_interp(1) == Agrif_Linear) .and. (coeffraf/=1) ) then
1558        call Linear1dPrecompute2d(indmax(2)-indmin(2)+1,            &
1559                                  indmax(1)-indmin(1)+1,            &
1560                                  petab_child(1)-pttab_child(1)+1,  &
1561                                  s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
1562        precomputedone(1) = .TRUE.
1563    elseif ( (type_interp(1) == Agrif_PPM) .and. (coeffraf/=1) ) then
1564        call PPM1dPrecompute2d(indmax(2)-indmin(2)+1,           &
1565                               indmax(1)-indmin(1)+1,           &
1566                               petab_child(1)-pttab_child(1)+1, &
1567                               s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
1568        precomputedone(1) = .TRUE.
1569    endif
1570
1571    coeffraf = nint ( ds_parent(2) / ds_child(2) )
1572    if ( (type_interp(2) == Agrif_Linear) .and. (coeffraf/=1) ) then
1573        call Linear1dPrecompute2d(petab_child(1)-pttab_child(1)+1,  &
1574                                  indmax(2)-indmin(2)+1,            &
1575                                  petab_child(2)-pttab_child(2)+1,  &
1576                                  s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
1577        precomputedone(2) = .TRUE.
1578    elseif ( (type_interp(2) == Agrif_PPM) .and. (coeffraf/=1) ) then
1579        call PPM1dPrecompute2d(petab_child(1)-pttab_child(1)+1, &
1580                               indmax(2)-indmin(2)+1,           &
1581                               petab_child(2)-pttab_child(2)+1, &
1582                               s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
1583        precomputedone(2) = .TRUE.
1584    endif
1585!
1586    do k = indmin(3), indmax(3)
1587        call Agrif_Interp_2D_recursive(type_interp(1:2),                            &
1588                                       tabin(indmin(1):indmax(1),                   &
1589                                             indmin(2):indmax(2), k),               &
1590                                       tabtemp(pttab_child(1):petab_child(1),       &
1591                                               pttab_child(2):petab_child(2), k),   &
1592                                       indmin(1:2), indmax(1:2),                    &
1593                                       pttab_child(1:2), petab_child(1:2),          &
1594                                       s_child(1:2),     s_parent(1:2),             &
1595                                       ds_child(1:2),   ds_parent(1:2) )
1596    enddo
1597!
1598    precomputedone(1) = .FALSE.
1599    precomputedone(2) = .FALSE.
1600    coeffraf = nint(ds_parent(3)/ds_child(3))
1601!
1602    if ( coeffraf == 1 ) then
1603        locind_child_left = 1 + agrif_int((s_child(3)-s_parent(3))/ds_parent(3))
1604        kdeb = indmin(3)+locind_child_left-2
1605        do k = pttab_child(3),petab_child(3)
1606            kdeb = kdeb + 1
1607            do j = pttab_child(2), petab_child(2)
1608            do i = pttab_child(1), petab_child(1)
1609                tabout(i,j,k) = tabtemp(i,j,kdeb)
1610            enddo
1611            enddo
1612        enddo
1613    else
1614        do j = pttab_child(2), petab_child(2)
1615        do i = pttab_child(1), petab_child(1)
1616            call Agrif_InterpBase(type_interp(3),                                   &
1617                                  tabtemp(i,j,indmin(3):indmax(3)),                 &
1618                                  tabout(i,j,pttab_child(3):petab_child(3)),        &
1619                                  indmin(3), indmax(3),                             &
1620                                  pttab_child(3), petab_child(3),                   &
1621                                  s_parent(3),    s_child(3),                       &
1622                                 ds_parent(3),   ds_child(3) )
1623        enddo
1624        enddo
1625    endif
1626!---------------------------------------------------------------------------------------------------
1627end subroutine Agrif_Interp_3D_recursive
1628!===================================================================================================
1629!
1630!===================================================================================================
1631!  subroutine Agrif_Interp_4D_Recursive
1632!
1633!> Subroutine for the interpolation of a 4D grid variable.
1634!> It calls #Agrif_Interp_3D_recursive and #Agrif_InterpBase.
1635!---------------------------------------------------------------------------------------------------
1636subroutine Agrif_Interp_4D_recursive ( type_interp, tabin, tabout,  &
1637                                       indmin, indmax,              &
1638                                       pttab_child, petab_child,    &
1639                                       s_child,  s_parent,          &
1640                                       ds_child, ds_parent )
1641!---------------------------------------------------------------------------------------------------
1642    integer, dimension(4),              intent(in)  :: type_interp
1643    integer, dimension(4),              intent(in)  :: indmin, indmax
1644    integer, dimension(4),              intent(in)  :: pttab_child, petab_child
1645    real(kind=8),    dimension(4),              intent(in)  :: s_child, s_parent
1646    real(kind=8),    dimension(4),              intent(in)  :: ds_child, ds_parent
1647    real,    dimension(                 &
1648        indmin(1):indmax(1),            &
1649        indmin(2):indmax(2),            &
1650        indmin(3):indmax(3),            &
1651        indmin(4):indmax(4)),           intent(in)  :: tabin
1652    real,    dimension(                 &
1653        pttab_child(1):petab_child(1),  &
1654        pttab_child(2):petab_child(2),  &
1655        pttab_child(3):petab_child(3),  &
1656        pttab_child(4):petab_child(4)), intent(out) :: tabout
1657!---------------------------------------------------------------------------------------------------
1658    real, dimension(                    &
1659        pttab_child(1):petab_child(1),  &
1660        pttab_child(2):petab_child(2),  &
1661        pttab_child(3):petab_child(3),  &
1662        indmin(4):indmax(4))            :: tabtemp
1663    integer                             :: i, j, k, l
1664!
1665    do l = indmin(4), indmax(4)
1666        call Agrif_Interp_3D_recursive(type_interp(1:3),                            &
1667                                       tabin(indmin(1):indmax(1),                   &
1668                                             indmin(2):indmax(2),                   &
1669                                             indmin(3):indmax(3), l),               &
1670                                       tabtemp(pttab_child(1):petab_child(1),       &
1671                                               pttab_child(2):petab_child(2),       &
1672                                               pttab_child(3):petab_child(3), l),   &
1673                                       indmin(1:3), indmax(1:3),                    &
1674                                       pttab_child(1:3), petab_child(1:3),          &
1675                                       s_child(1:3),    s_parent(1:3),              &
1676                                       ds_child(1:3),  ds_parent(1:3) )
1677    enddo
1678!
1679    do k = pttab_child(3), petab_child(3)
1680    do j = pttab_child(2), petab_child(2)
1681    do i = pttab_child(1), petab_child(1)
1682        call Agrif_InterpBase(type_interp(4),                                       &
1683                              tabtemp(i,j,k,indmin(4):indmax(4)),                   &
1684                              tabout(i,j,k,pttab_child(4):petab_child(4)),          &
1685                              indmin(4), indmax(4),                                 &
1686                              pttab_child(4), petab_child(4),                       &
1687                              s_parent(4),    s_child(4),                           &
1688                             ds_parent(4),   ds_child(4) )
1689    enddo
1690    enddo
1691    enddo
1692!---------------------------------------------------------------------------------------------------
1693end subroutine Agrif_Interp_4D_recursive
1694!===================================================================================================
1695!
1696!===================================================================================================
1697!  subroutine Agrif_Interp_5D_Recursive
1698!
1699!> Subroutine for the interpolation of a 5D grid variable.
1700!> It calls #Agrif_Interp_4D_recursive and #Agrif_InterpBase.
1701!---------------------------------------------------------------------------------------------------
1702subroutine Agrif_Interp_5D_recursive ( type_interp, tabin, tabout,  &
1703                                       indmin, indmax,              &
1704                                       pttab_child, petab_child,    &
1705                                       s_child,  s_parent,          &
1706                                       ds_child, ds_parent )
1707!---------------------------------------------------------------------------------------------------
1708    integer, dimension(5),              intent(in)  :: type_interp
1709    integer, dimension(5),              intent(in)  :: indmin, indmax
1710    integer, dimension(5),              intent(in)  :: pttab_child, petab_child
1711    real(kind=8),    dimension(5),              intent(in)  :: s_child, s_parent
1712    real(kind=8),    dimension(5),              intent(in)  :: ds_child, ds_parent
1713    real,    dimension(                 &
1714        indmin(1):indmax(1),            &
1715        indmin(2):indmax(2),            &
1716        indmin(3):indmax(3),            &
1717        indmin(4):indmax(4),            &
1718        indmin(5):indmax(5)),           intent(in)  :: tabin
1719    real,    dimension(                 &
1720        pttab_child(1):petab_child(1),  &
1721        pttab_child(2):petab_child(2),  &
1722        pttab_child(3):petab_child(3),  &
1723        pttab_child(4):petab_child(4),  &
1724        pttab_child(5):petab_child(5)), intent(out) :: tabout
1725!---------------------------------------------------------------------------------------------------
1726    real, dimension(                    &
1727        pttab_child(1):petab_child(1),  &
1728        pttab_child(2):petab_child(2),  &
1729        pttab_child(3):petab_child(3),  &
1730        pttab_child(4):petab_child(4),  &
1731        indmin(5):indmax(5))            :: tabtemp
1732    integer                             :: i, j, k, l, m
1733!
1734    do m = indmin(5), indmax(5)
1735        call Agrif_Interp_4D_recursive(type_interp(1:4),                            &
1736                                       tabin(indmin(1):indmax(1),                   &
1737                                             indmin(2):indmax(2),                   &
1738                                             indmin(3):indmax(3),                   &
1739                                             indmin(4):indmax(4),m),                &
1740                                       tabtemp(pttab_child(1):petab_child(1),       &
1741                                               pttab_child(2):petab_child(2),       &
1742                                               pttab_child(3):petab_child(3),       &
1743                                               pttab_child(4):petab_child(4), m),   &
1744                                       indmin(1:4),indmax(1:4),                     &
1745                                       pttab_child(1:4), petab_child(1:4),          &
1746                                       s_child(1:4),     s_parent(1:4),             &
1747                                       ds_child(1:4),   ds_parent(1:4) )
1748    enddo
1749!
1750    do l = pttab_child(4), petab_child(4)
1751    do k = pttab_child(3), petab_child(3)
1752    do j = pttab_child(2), petab_child(2)
1753    do i = pttab_child(1), petab_child(1)
1754        call Agrif_InterpBase(type_interp(5),                                       &
1755                              tabtemp(i,j,k,l,indmin(5):indmax(5)),                 &
1756                              tabout(i,j,k,l,pttab_child(5):petab_child(5)),        &
1757                              indmin(5), indmax(5),                                 &
1758                              pttab_child(5), petab_child(5),                       &
1759                              s_parent(5),   s_child(5),                            &
1760                              ds_parent(5), ds_child(5) )
1761    enddo
1762    enddo
1763    enddo
1764    enddo
1765!---------------------------------------------------------------------------------------------------
1766end subroutine Agrif_Interp_5D_recursive
1767!===================================================================================================
1768!
1769!===================================================================================================
1770!  subroutine Agrif_Interp_6D_Recursive
1771!
1772!> Subroutine for the interpolation of a 6D grid variable.
1773!> It calls #Agrif_Interp_5D_recursive and Agrif_InterpBase.
1774!---------------------------------------------------------------------------------------------------
1775subroutine Agrif_Interp_6D_recursive ( type_interp, tabin, tabout,  &
1776                                       indmin, indmax,              &
1777                                       pttab_child, petab_child,    &
1778                                       s_child,  s_parent,          &
1779                                       ds_child, ds_parent )
1780!---------------------------------------------------------------------------------------------------
1781    integer, dimension(6),                  intent(in)  :: type_interp
1782    integer, dimension(6),                  intent(in)  :: indmin, indmax
1783    integer, dimension(6),                  intent(in)  :: pttab_child, petab_child
1784    real(kind=8),    dimension(6),                  intent(in)  :: s_child, s_parent
1785    real(kind=8),    dimension(6),                  intent(in)  :: ds_child, ds_parent
1786    real,    dimension(                 &
1787        indmin(1):indmax(1),            &
1788        indmin(2):indmax(2),            &
1789        indmin(3):indmax(3),            &
1790        indmin(4):indmax(4),            &
1791        indmin(5):indmax(5),            &
1792        indmin(6):indmax(6)),               intent(in)  :: tabin
1793    real,    dimension(                 &
1794        pttab_child(1):petab_child(1),  &
1795        pttab_child(2):petab_child(2),  &
1796        pttab_child(3):petab_child(3),  &
1797        pttab_child(4):petab_child(4),  &
1798        pttab_child(5):petab_child(5),  &
1799        pttab_child(6):petab_child(6)),     intent(out) :: tabout
1800!---------------------------------------------------------------------------------------------------
1801    real, dimension(                    &
1802        pttab_child(1):petab_child(1),  &
1803        pttab_child(2):petab_child(2),  &
1804        pttab_child(3):petab_child(3),  &
1805        pttab_child(4):petab_child(4),  &
1806        pttab_child(5):petab_child(5),  &
1807        indmin(6):indmax(6))            :: tabtemp
1808    integer                             :: i, j, k, l, m, n
1809!
1810    do n = indmin(6), indmax(6)
1811        call Agrif_Interp_5D_recursive(type_interp(1:5),                            &
1812                                       tabin(indmin(1):indmax(1),                   &
1813                                             indmin(2):indmax(2),                   &
1814                                             indmin(3):indmax(3),                   &
1815                                             indmin(4):indmax(4),                   &
1816                                             indmin(5):indmax(5), n),               &
1817                                        tabtemp(pttab_child(1):petab_child(1),      &
1818                                                pttab_child(2):petab_child(2),      &
1819                                                pttab_child(3):petab_child(3),      &
1820                                                pttab_child(4):petab_child(4),      &
1821                                                pttab_child(5):petab_child(5), n),  &
1822                                        indmin(1:5),indmax(1:5),                    &
1823                                        pttab_child(1:5), petab_child(1:5),         &
1824                                        s_child(1:5), s_parent(1:5),                &
1825                                        ds_child(1:5),ds_parent(1:5) )
1826    enddo
1827!
1828    do m = pttab_child(5), petab_child(5)
1829    do l = pttab_child(4), petab_child(4)
1830    do k = pttab_child(3), petab_child(3)
1831    do j = pttab_child(2), petab_child(2)
1832    do i = pttab_child(1), petab_child(1)
1833        call Agrif_InterpBase(type_interp(6),                                       &
1834                              tabtemp(i,j,k,l,m,indmin(6):indmax(6)),               &
1835                              tabout(i,j,k,l,m,pttab_child(6):petab_child(6)),      &
1836                              indmin(6), indmax(6),                                 &
1837                              pttab_child(6), petab_child(6),                       &
1838                              s_parent(6),   s_child(6),                            &
1839                              ds_parent(6), ds_child(6) )
1840    enddo
1841    enddo
1842    enddo
1843    enddo
1844    enddo
1845!---------------------------------------------------------------------------------------------------
1846end subroutine Agrif_Interp_6D_recursive
1847!===================================================================================================
1848!
1849!===================================================================================================
1850!  subroutine Agrif_InterpBase
1851!
1852!> Calls the interpolation method chosen by the user (linear, lagrange, spline, etc.).
1853!---------------------------------------------------------------------------------------------------
1854subroutine Agrif_InterpBase ( type_interp, parenttab, childtab, indmin, indmax, &
1855                              pttab_child, petab_child,                         &
1856                              s_parent, s_child, ds_parent, ds_child )
1857!---------------------------------------------------------------------------------------------------
1858    INTEGER                                                 :: type_interp
1859    INTEGER                                                 :: indmin, indmax
1860    INTEGER                                                 :: pttab_child, petab_child
1861    REAL, DIMENSION(indmin:indmax),           INTENT(IN)    :: parenttab
1862    REAL, DIMENSION(pttab_child:petab_child), INTENT(OUT)   :: childtab
1863    REAL(kind=8)                                            :: s_parent, s_child
1864    REAL(kind=8)                                            :: ds_parent,ds_child
1865!
1866    if ( (indmin == indmax) .and. (pttab_child == petab_child) ) then
1867!
1868        childtab(pttab_child) = parenttab(indmin)
1869!
1870    elseif (type_interp == Agrif_LINEAR) then    !       Linear interpolation
1871!
1872        call Agrif_basicinterp_linear1D(parenttab,childtab,                   &
1873                indmax-indmin+1,petab_child-pttab_child+1,  &
1874                s_parent,s_child,ds_parent,ds_child)
1875!
1876    elseif ( type_interp == Agrif_PPM ) then     !       PPM interpolation
1877
1878        call PPM1d(parenttab,childtab,                      &
1879                indmax-indmin+1,petab_child-pttab_child+1,  &
1880                s_parent,s_child,ds_parent,ds_child)
1881!
1882    elseif ( type_interp == Agrif_PPM_LIM ) then     !       PPM interpolation
1883
1884        call PPM1d_lim(parenttab,childtab,                      &
1885                indmax-indmin+1,petab_child-pttab_child+1,  &
1886                s_parent,s_child,ds_parent,ds_child)
1887!
1888    elseif (type_interp == Agrif_LAGRANGE) then  !       Lagrange interpolation
1889!
1890        call lagrange1D(parenttab,childtab,                 &
1891                indmax-indmin+1,petab_child-pttab_child+1,  &
1892                s_parent,s_child,ds_parent,ds_child)
1893!
1894    elseif (type_interp == Agrif_ENO) then       !       Eno interpolation
1895!
1896        call ENO1d(parenttab,childtab,                      &
1897                indmax-indmin+1,petab_child-pttab_child+1,  &
1898                s_parent,s_child,ds_parent,ds_child)
1899!
1900    elseif (type_interp == Agrif_WENO) then      !       Weno interpolation
1901!
1902        call WENO1d(parenttab,childtab,                     &
1903                indmax-indmin+1,petab_child-pttab_child+1,  &
1904                s_parent,s_child,ds_parent,ds_child)
1905!
1906    elseif (type_interp == Agrif_LINEARCONSERV) then !   Linear conservative interpolation
1907!
1908        call Linear1dConserv(parenttab,childtab,            &
1909                indmax-indmin+1,petab_child-pttab_child+1,  &
1910                s_parent,s_child,ds_parent,ds_child)
1911!
1912    elseif (type_interp == Agrif_LINEARCONSERVLIM) then !Linear conservative interpolation
1913!
1914        call Linear1dConservLim(parenttab,childtab,         &
1915                indmax-indmin+1,petab_child-pttab_child+1,  &
1916                s_parent,s_child,ds_parent,ds_child)
1917!
1918    elseif (type_interp == Agrif_CONSTANT) then
1919!
1920        call Constant1d(parenttab,childtab,                 &
1921                indmax-indmin+1,petab_child-pttab_child+1,  &
1922                s_parent,s_child,ds_parent,ds_child)
1923!
1924    endif
1925!---------------------------------------------------------------------------------------------------
1926end subroutine Agrif_InterpBase
1927!===================================================================================================
1928!
1929!===================================================================================================
1930!  subroutine Agrif_Find_list_interp
1931!---------------------------------------------------------------------------------------------------
1932function Agrif_Find_list_interp ( list_interp, pttab, petab, pttab_Child, pttab_Parent,     &
1933                                    nbdim, indmin, indmax, indminglob,  indmaxglob,         &
1934                                    pttruetab, cetruetab, memberin,                         &
1935                                    parentarray_chunk,parentarray_chunk_decal,decal_chunks,   &
1936                                    correction_required,member_chuncks,nb_chunks            &
1937#if defined AGRIF_MPI
1938                                   ,indminglob2, indmaxglob2, parentarray,                  &
1939                                    member, tab4t, memberinall, sendtoproc1, recvfromproc1  &
1940#endif
1941    ) result(find_list_interp)
1942!---------------------------------------------------------------------------------------------------
1943    type(Agrif_List_Interp_Loc),   pointer     :: list_interp
1944    integer,                       intent(in)  :: nbdim
1945    integer, dimension(nbdim),     intent(in)  :: pttab, petab, pttab_Child, pttab_Parent
1946    integer, dimension(nbdim),     intent(out) :: indmin, indmax
1947    integer, dimension(nbdim),     intent(out) :: indminglob, indmaxglob
1948    integer, dimension(nbdim),     intent(out) :: pttruetab, cetruetab
1949    logical,                       intent(out) :: memberin
1950    integer :: nb_chunks
1951    integer, dimension(:,:,:,:), allocatable :: parentarray_chunk
1952    integer, dimension(:,:,:,:), allocatable :: parentarray_chunk_decal
1953    integer, dimension(:,:),allocatable :: decal_chunks
1954    logical, dimension(:),allocatable :: correction_required
1955    logical, dimension(:),allocatable :: member_chuncks
1956#if defined AGRIF_MPI
1957    integer, dimension(nbdim),     intent(out) :: indminglob2, indmaxglob2
1958    integer, dimension(nbdim,2,2), intent(out) :: parentarray
1959    logical,                       intent(out) :: member
1960    integer, dimension(nbdim,0:Agrif_Nbprocs-1,8), intent(out) :: tab4t
1961    logical, dimension(0:Agrif_Nbprocs-1),         intent(out) :: memberinall
1962    logical, dimension(0:Agrif_Nbprocs-1),         intent(out) :: sendtoproc1, recvfromproc1
1963#endif
1964    logical :: find_list_interp
1965!
1966    integer :: i
1967    type(Agrif_List_Interp_Loc), pointer :: parcours
1968    type(Agrif_Interp_Loc),      pointer :: pil
1969
1970    find_list_interp = .false.
1971
1972    if ( .not. associated(list_interp) )    return
1973
1974    parcours => list_interp
1975    find_loop : do while ( associated(parcours) )
1976
1977        pil => parcours % interp_loc
1978
1979        do i = 1,nbdim
1980            if ( (pttab(i) /= pil % pttab(i)) .or. &
1981                 (petab(i) /= pil % petab(i)) .or. &
1982                 (pttab_child(i)  /= pil % pttab_child(i)) .or. &
1983                 (pttab_parent(i) /= pil % pttab_parent(i)) ) then
1984                parcours => parcours % suiv
1985                cycle find_loop
1986            endif
1987        enddo
1988
1989        indmin = pil % indmin(1:nbdim)
1990        indmax = pil % indmax(1:nbdim)
1991
1992        pttruetab = pil % pttruetab(1:nbdim)
1993        cetruetab = pil % cetruetab(1:nbdim)
1994
1995#if !defined AGRIF_MPI
1996        indminglob  = pil % indminglob(1:nbdim)
1997        indmaxglob  = pil % indmaxglob(1:nbdim)
1998#else
1999        indminglob  = pil % indminglob2(1:nbdim)
2000        indmaxglob  = pil % indmaxglob2(1:nbdim)
2001        indminglob2 = pil % indminglob2(1:nbdim)
2002        indmaxglob2 = pil % indmaxglob2(1:nbdim)
2003        parentarray = pil % parentarray(1:nbdim,:,:)
2004        member      = pil % member
2005        tab4t         = pil % tab4t(1:nbdim, 0:Agrif_Nbprocs-1, 1:8)
2006        memberinall   = pil % memberinall(0:Agrif_Nbprocs-1)
2007        sendtoproc1   = pil % sendtoproc1(0:Agrif_Nbprocs-1)
2008        recvfromproc1 = pil % recvfromproc1(0:Agrif_Nbprocs-1)
2009#endif
2010        memberin = pil % memberin
2011
2012! chunks
2013        nb_chunks = pil % nb_chunks
2014        Allocate(parentarray_chunk(nb_chunks,nbdim,2,2))
2015        parentarray_chunk = pil % parentarray_chunk
2016        Allocate(parentarray_chunk_decal(nb_chunks,nbdim,2,2))
2017        parentarray_chunk_decal = pil % parentarray_chunk_decal
2018        Allocate(correction_required(nb_chunks))
2019        correction_required = pil % correction_required
2020        Allocate(decal_chunks(nb_chunks,nbdim))
2021        decal_chunks = pil % decal_chunks
2022        Allocate(member_chuncks(nb_chunks))
2023        member_chuncks = pil % member_chuncks
2024
2025        find_list_interp = .true.
2026        exit find_loop
2027    enddo find_loop
2028!---------------------------------------------------------------------------------------------------
2029end function Agrif_Find_list_interp
2030!===================================================================================================
2031!
2032!===================================================================================================
2033!  subroutine Agrif_AddTo_list_interp
2034!---------------------------------------------------------------------------------------------------
2035subroutine Agrif_AddTo_list_interp ( list_interp, pttab, petab, pttab_Child, pttab_Parent,  &
2036                                     indmin, indmax, indminglob, indmaxglob,                &
2037                                     pttruetab, cetruetab,                                  &
2038                                     memberin, nbdim,                                       &
2039                                    parentarray_chunk,parentarray_chunk_decal,decal_chunks,  &
2040                                    correction_required,member_chuncks,nb_chunks            &
2041#if defined AGRIF_MPI
2042                                    ,indminglob2, indmaxglob2,                              &
2043                                     parentarray,                                           &
2044                                     member,                                                &
2045                                     tab4t, memberinall, sendtoproc1, recvfromproc1         &
2046#endif
2047    )
2048!---------------------------------------------------------------------------------------------------
2049    type(Agrif_List_Interp_Loc), pointer    :: list_interp
2050    integer                                 :: nbdim
2051    integer, dimension(nbdim)               :: pttab, petab, pttab_Child, pttab_Parent
2052    integer, dimension(nbdim)               :: indmin,indmax
2053    integer, dimension(nbdim)               :: indminglob, indmaxglob
2054    integer, dimension(nbdim)               :: pttruetab, cetruetab
2055    logical                                 :: memberin
2056    integer :: nb_chunks
2057    integer, dimension(:,:,:,:), allocatable :: parentarray_chunk
2058    integer, dimension(:,:,:,:), allocatable :: parentarray_chunk_decal
2059    integer, dimension(:,:),allocatable :: decal_chunks
2060    logical, dimension(:),allocatable :: correction_required
2061    logical, dimension(:),allocatable :: member_chuncks
2062#if defined AGRIF_MPI
2063    integer, dimension(nbdim,2,2)           :: parentarray
2064    logical                                 :: member
2065    integer, dimension(nbdim)                       :: indminglob2,indmaxglob2
2066    integer, dimension(nbdim,0:Agrif_Nbprocs-1,8)   :: tab4t
2067    logical, dimension(0:Agrif_Nbprocs-1)           :: memberinall
2068    logical, dimension(0:Agrif_Nbprocs-1)           :: sendtoproc1
2069    logical, dimension(0:Agrif_Nbprocs-1)           :: recvfromproc1
2070#endif
2071!
2072    type(Agrif_List_Interp_Loc), pointer    :: parcours
2073    type(Agrif_Interp_Loc),      pointer    :: pil
2074!
2075    allocate(parcours)
2076    allocate(parcours % interp_loc)
2077
2078    pil => parcours % interp_loc
2079
2080    pil % pttab(1:nbdim) = pttab(1:nbdim)
2081    pil % petab(1:nbdim) = petab(1:nbdim)
2082    pil % pttab_child(1:nbdim) = pttab_child(1:nbdim)
2083    pil % pttab_parent(1:nbdim) = pttab_parent(1:nbdim)
2084
2085    pil % indmin(1:nbdim) = indmin(1:nbdim)
2086    pil % indmax(1:nbdim) = indmax(1:nbdim)
2087
2088    pil % memberin = memberin
2089#if !defined AGRIF_MPI
2090    pil % indminglob(1:nbdim) = indminglob(1:nbdim)
2091    pil % indmaxglob(1:nbdim) = indmaxglob(1:nbdim)
2092#else
2093    pil % indminglob2(1:nbdim) = indminglob2(1:nbdim)
2094    pil % indmaxglob2(1:nbdim) = indmaxglob2(1:nbdim)
2095    pil % parentarray(1:nbdim,:,:) = parentarray(1:nbdim,:,:)
2096    pil % member = member
2097    allocate(pil % tab4t(nbdim, 0:Agrif_Nbprocs-1, 8))
2098    allocate(pil % memberinall(0:Agrif_Nbprocs-1))
2099    allocate(pil % sendtoproc1(0:Agrif_Nbprocs-1))
2100    allocate(pil % recvfromproc1(0:Agrif_Nbprocs-1))
2101    pil % tab4t         = tab4t
2102    pil % memberinall   = memberinall
2103    pil % sendtoproc1   = sendtoproc1
2104    pil % recvfromproc1 = recvfromproc1
2105#endif
2106
2107    pil % pttruetab(1:nbdim) = pttruetab(1:nbdim)
2108    pil % cetruetab(1:nbdim) = cetruetab(1:nbdim)
2109
2110! chunks
2111    pil % nb_chunks = nb_chunks
2112    allocate(pil % parentarray_chunk(nb_chunks,nbdim,2,2))
2113    allocate(pil % parentarray_chunk_decal(nb_chunks,nbdim,2,2))
2114    allocate(pil % correction_required(nb_chunks))
2115    allocate(pil % decal_chunks(nb_chunks,nbdim))
2116    allocate(pil % member_chuncks(nb_chunks))
2117
2118    pil % parentarray_chunk   = parentarray_chunk
2119    pil % parentarray_chunk_decal = parentarray_chunk_decal
2120    pil % correction_required = correction_required
2121    pil % decal_chunks        = decal_chunks
2122    pil % member_chuncks      = member_chuncks
2123
2124
2125    parcours % suiv => list_interp
2126    list_interp => parcours
2127!---------------------------------------------------------------------------------------------------
2128end subroutine Agrif_Addto_list_interp
2129!===================================================================================================
2130!
2131end module Agrif_Interpolation
Note: See TracBrowser for help on using the repository browser.