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 trunk/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES – NEMO

source: trunk/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES/modinterp.F90 @ 5656

Last change on this file since 5656 was 5656, checked in by timgraham, 9 years ago

Merge of AGRIF branch (branches/2014/dev_r4765_CNRS_agrif) onto the trunk

  • Property svn:keywords set to Id
File size: 76.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   , dimension(6) :: s_child,   s_parent
69    real   , 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,    DIMENSION(nbdim), INTENT(in)   :: s_Child,s_Parent   !< Positions of the parent and child grids
118    REAL,    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, DIMENSION(nbdim)     :: pttruetab,cetruetab
130    INTEGER, DIMENSION(nbdim)     :: indmin,     indmax
131    INTEGER, DIMENSION(nbdim)     :: indminglob, indmaxglob
132#if defined AGRIF_MPI
133    INTEGER, DIMENSION(nbdim)     :: indminglob2,indmaxglob2
134#endif
135    LOGICAL, DIMENSION(nbdim)     :: noraftab
136    REAL   , DIMENSION(nbdim)     :: s_Child_temp,s_Parent_temp
137    INTEGER, DIMENSION(nbdim)     :: lowerbound, upperbound, coords
138    INTEGER, DIMENSION(nbdim,2,2), INTENT(OUT) :: childarray
139    INTEGER, DIMENSION(nbdim,2,2)              :: parentarray
140    LOGICAL :: member
141    LOGICAL :: find_list_interp
142!
143#if defined AGRIF_MPI
144!
145    INTEGER, PARAMETER          :: etiquette = 100
146    INTEGER                     :: code, local_proc
147    INTEGER, DIMENSION(nbdim,4)                     :: tab3
148    INTEGER, DIMENSION(nbdim,4,0:Agrif_Nbprocs-1)   :: tab4
149    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1,8)   :: tab4t
150    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1)           :: memberinall
151    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1)           :: sendtoproc1, recvfromproc1
152    LOGICAL, DIMENSION(1)                           :: memberin1
153    LOGICAL                                         :: memberout
154!
155#endif
156!
157    type(Agrif_Variable), pointer, save :: tempC => NULL()        ! Temporary child grid variable
158    type(Agrif_Variable), pointer, save :: tempP => NULL()        ! Temporary parent grid variable
159    type(Agrif_Variable), pointer, save :: tempPextend => NULL()  ! Temporary parent grid variable
160    type(Agrif_Variable), pointer, save :: parentvalues => NULL()
161!
162    coords = child % root_var % coords
163!
164!   Boundaries of the current grid where interpolation is done
165    find_list_interp =                                              &
166        Agrif_Find_list_interp(                                     &
167            child % list_interp,                                    &
168            pttab, petab, pttab_Child, pttab_Parent, nbdim,         &
169            indmin, indmax, indminglob, indmaxglob,                 &
170            pttruetab, cetruetab, memberin                          &
171#if defined AGRIF_MPI
172           ,indminglob2, indmaxglob2, parentarray,                  &
173            member, tab4t,memberinall,  sendtoproc1, recvfromproc1  &
174#endif
175        )
176!
177    if (.not.find_list_interp) then
178!
179        call Agrif_get_var_bounds_array(child, lowerbound, upperbound, nbdim)
180        call Agrif_Childbounds(nbdim, lowerbound, upperbound,               &
181                               pttab, petab, Agrif_Procrank, coords,        &
182                               pttruetab, cetruetab, memberin)
183        call Agrif_Parentbounds(type_interp,nbdim,indminglob,indmaxglob,    &
184                                s_Parent_temp,s_Child_temp,                 &
185                                s_Child,ds_Child,                           &
186                                s_Parent,ds_Parent,                         &
187                                pttab,petab,                                &
188                                pttab_Child,pttab_Parent,                   &
189                                child%root_var % posvar, coords)
190#if defined AGRIF_MPI
191        if (memberin) then
192            call Agrif_Parentbounds(type_interp,nbdim,indmin,indmax,        &
193                                    s_Parent_temp,s_Child_temp,             &
194                                    s_Child,ds_Child,                       &
195                                    s_Parent,ds_Parent,                     &
196                                    pttruetab,cetruetab,                    &
197                                    pttab_Child,pttab_Parent,               &
198                                    child%root_var % posvar, coords)
199        endif
200
201        local_proc = Agrif_Procrank
202        call Agrif_get_var_bounds_array(parent,lowerbound,upperbound,nbdim)
203        call Agrif_ChildGrid_to_ParentGrid()
204!
205        call Agrif_Childbounds(nbdim,lowerbound,upperbound,                 &
206                               indminglob,indmaxglob, local_proc, coords,   &
207                               indminglob2,indmaxglob2,member)
208!
209        if (member) then
210            call Agrif_GlobalToLocalBounds(parentarray,                     &
211                                           lowerbound,  upperbound,         &
212                                           indminglob2, indmaxglob2, coords,&
213                                           nbdim, local_proc, member)
214        endif
215
216        call Agrif_ParentGrid_to_ChildGrid()
217#else
218        parentarray(:,1,1) = indminglob
219        parentarray(:,2,1) = indmaxglob
220        parentarray(:,1,2) = indminglob
221        parentarray(:,2,2) = indmaxglob
222        indmin = indminglob
223        indmax = indmaxglob
224        member = .TRUE.
225#endif
226
227    else
228
229#if defined AGRIF_MPI
230        s_Parent_temp = s_Parent + (indmin - pttab_Parent) * ds_Parent
231        s_Child_temp  = s_Child + (pttruetab - pttab_Child) * ds_Child
232#else
233        parentarray(:,1,1) = indminglob
234        parentarray(:,2,1) = indmaxglob
235        parentarray(:,1,2) = indminglob
236        parentarray(:,2,2) = indmaxglob
237        indmin = indminglob
238        indmax = indmaxglob
239        member = .TRUE.
240        s_Parent_temp = s_Parent + (indminglob - pttab_Parent) * ds_Parent
241        s_Child_temp  = s_Child + (pttab - pttab_Child) * ds_Child
242#endif
243    endif
244!
245    if (member) then
246        if (.not.associated(tempP)) allocate(tempP)
247!
248        call Agrif_array_allocate(tempP,parentarray(:,1,1),parentarray(:,2,1),nbdim)
249        call Agrif_var_set_array_tozero(tempP,nbdim)
250
251        call Agrif_ChildGrid_to_ParentGrid()
252!
253        select case (nbdim)
254        case(1)
255            call procname(tempP%array1,                         &
256                      parentarray(1,1,2),parentarray(1,2,2),.TRUE.,nb,ndir)
257        case(2)
258            call procname(tempP%array2,                         &
259                      parentarray(1,1,2),parentarray(1,2,2),    &
260                      parentarray(2,1,2),parentarray(2,2,2),.TRUE.,nb,ndir)
261        case(3)
262            call procname(tempP%array3,                         &
263                      parentarray(1,1,2),parentarray(1,2,2),    &
264                      parentarray(2,1,2),parentarray(2,2,2),    &
265                      parentarray(3,1,2),parentarray(3,2,2),.TRUE.,nb,ndir)
266        case(4)
267            call procname(tempP%array4,                         &
268                      parentarray(1,1,2),parentarray(1,2,2),    &
269                      parentarray(2,1,2),parentarray(2,2,2),    &
270                      parentarray(3,1,2),parentarray(3,2,2),    &
271                      parentarray(4,1,2),parentarray(4,2,2),.TRUE.,nb,ndir)
272        case(5)
273            call procname(tempP%array5,                         &
274                      parentarray(1,1,2),parentarray(1,2,2),    &
275                      parentarray(2,1,2),parentarray(2,2,2),    &
276                      parentarray(3,1,2),parentarray(3,2,2),    &
277                      parentarray(4,1,2),parentarray(4,2,2),    &
278                      parentarray(5,1,2),parentarray(5,2,2),.TRUE.,nb,ndir)
279        case(6)
280            call procname(tempP%array6,                         &
281                      parentarray(1,1,2),parentarray(1,2,2),    &
282                      parentarray(2,1,2),parentarray(2,2,2),    &
283                      parentarray(3,1,2),parentarray(3,2,2),    &
284                      parentarray(4,1,2),parentarray(4,2,2),    &
285                      parentarray(5,1,2),parentarray(5,2,2),    &
286                      parentarray(6,1,2),parentarray(6,2,2),.TRUE.,nb,ndir)
287        end select
288!
289        call Agrif_ParentGrid_to_ChildGrid()
290!
291    endif
292
293#if defined AGRIF_MPI
294    if (.not.find_list_interp) then
295!
296        tab3(:,1) = indminglob2(:)
297        tab3(:,2) = indmaxglob2(:)
298        tab3(:,3) = indmin(:)
299        tab3(:,4) = indmax(:)
300!
301        call MPI_ALLGATHER(tab3,4*nbdim,MPI_INTEGER,tab4,4*nbdim,MPI_INTEGER,Agrif_mpi_comm,code)
302
303        if (.not.associated(tempPextend))   allocate(tempPextend)
304
305        do k=0,Agrif_Nbprocs-1
306            do j=1,4
307                do i=1,nbdim
308                    tab4t(i,k,j) = tab4(i,j,k)
309                enddo
310            enddo
311        enddo
312
313        memberin1(1) = memberin
314        call MPI_ALLGATHER(memberin1,1,MPI_LOGICAL,memberinall,1,MPI_LOGICAL,Agrif_mpi_comm,code)
315
316        call Get_External_Data_first(tab4t(:,:,1),tab4t(:,:,2),         &
317                                     tab4t(:,:,3),tab4t(:,:,4),         &
318                                     nbdim,memberinall, coords,         &
319                                     sendtoproc1,recvfromproc1,         &
320                                     tab4t(:,:,5),tab4t(:,:,6),         &
321                                     tab4t(:,:,7),tab4t(:,:,8) )
322    endif
323
324    call ExchangeSameLevel(sendtoproc1,recvfromproc1,nbdim,         &
325            tab4t(:,:,3),tab4t(:,:,4),tab4t(:,:,5),tab4t(:,:,6),    &
326            tab4t(:,:,7),tab4t(:,:,8),memberin,tempP,tempPextend)
327#else
328    tempPextend => tempP
329#endif
330
331    if (.not.find_list_interp) then
332        call Agrif_Addto_list_interp(                           &
333                child%list_interp,pttab,petab,                  &
334                pttab_Child,pttab_Parent,indmin,indmax,         &
335                indminglob,indmaxglob,                          &
336                pttruetab,cetruetab,                            &
337                memberin,nbdim                                  &
338#if defined AGRIF_MPI
339               ,indminglob2,indmaxglob2,                        &
340                parentarray,                                    &
341                member,                                         &
342                tab4t,memberinall,sendtoproc1,recvfromproc1     &
343#endif
344        )
345    endif
346!
347    if (memberin) then
348!
349        if (.not.associated(tempC)) allocate(tempC)
350!
351        call Agrif_array_allocate(tempC,pttruetab,cetruetab,nbdim)
352!
353!       Special values on the parent grid
354        if (Agrif_UseSpecialValue) then
355!
356            noraftab(1:nbdim) = child % root_var % interptab(1:nbdim) == 'N'
357!
358            if (.not.associated(parentvalues))  allocate(parentvalues)
359!
360            call Agrif_array_allocate(parentvalues,indmin,indmax,nbdim)
361            call Agrif_var_full_copy_array(parentvalues,tempPextend,nbdim)
362!
363            call Agrif_CheckMasknD(tempPextend,parentvalues,    &
364                    indmin(1:nbdim),indmax(1:nbdim),            &
365                    indmin(1:nbdim),indmax(1:nbdim),            &
366                    noraftab(1:nbdim),nbdim)
367!
368            call Agrif_array_deallocate(parentvalues,nbdim)
369!
370        endif
371!
372!       Interpolation of the current grid
373!
374        if ( memberin ) then
375            select case(nbdim)
376            case(1)
377                call Agrif_Interp_1D_recursive( type_interp(1),                         &
378                                                tempPextend%array1,                     &
379                                                tempC%array1,                           &
380                                                indmin(1), indmax(1),                   &
381                                                pttruetab(1),    cetruetab(1),          &
382                                                s_Child_temp(1), s_Parent_temp(1),      &
383                                                ds_Child(1),     ds_Parent(1) )
384            case(2)
385                call Agrif_Interp_2D_recursive( type_interp(1:2),                       &
386                                                tempPextend % array2,                   &
387                                                tempC       % array2,                   &
388                                                indmin(1:2), indmax(1:2),               &
389                                                pttruetab(1:2),    cetruetab(1:2),      &
390                                                s_Child_temp(1:2), s_Parent_temp(1:2),  &
391                                                ds_Child(1:2),    ds_Parent(1:2) )
392            case(3)
393                call Agrif_Interp_3D_recursive( type_interp(1:3),                       &
394                                                tempPextend % array3,                   &
395                                                tempC       % array3,                   &
396                                                indmin(1:3), indmax(1:3),               &
397                                                pttruetab(1:3),    cetruetab(1:3),      &
398                                                s_Child_temp(1:3), s_Parent_temp(1:3),  &
399                                                ds_Child(1:3),    ds_Parent(1:3) )
400            case(4)
401                call Agrif_Interp_4D_recursive( type_interp(1:4),                       &
402                                                tempPextend % array4,                   &
403                                                tempC       % array4,                   &
404                                                indmin(1:4), indmax(1:4),               &
405                                                pttruetab(1:4),    cetruetab(1:4),      &
406                                                s_Child_temp(1:4), s_Parent_temp(1:4),  &
407                                                ds_Child(1:4),    ds_Parent(1:4) )
408            case(5)
409                call Agrif_Interp_5D_recursive( type_interp(1:5),                       &
410                                                tempPextend % array5,                   &
411                                                tempC       % array5,                   &
412                                                indmin(1:5), indmax(1:5),               &
413                                                pttruetab(1:5),    cetruetab(1:5),      &
414                                                s_Child_temp(1:5), s_Parent_temp(1:5),  &
415                                                ds_Child(1:5),    ds_Parent(1:5) )
416            case(6)
417                call Agrif_Interp_6D_recursive( type_interp(1:6),                       &
418                                                tempPextend % array6,                   &
419                                                tempC       % array6,                   &
420                                                indmin(1:6), indmax(1:6),               &
421                                                pttruetab(1:6),    cetruetab(1:6),      &
422                                                s_Child_temp(1:6), s_Parent_temp(1:6),  &
423                                                ds_Child(1:6),    ds_Parent(1:6) )
424            end select
425!
426            call Agrif_get_var_bounds_array(child,lowerbound,upperbound,nbdim)
427
428#if defined AGRIF_MPI
429            call Agrif_GlobalToLocalBounds(childarray, lowerbound, upperbound,  &
430                                            pttruetab, cetruetab, coords,       &
431                                            nbdim, Agrif_Procrank, memberout)
432#else
433            childarray(:,1,1) = pttruetab
434            childarray(:,2,1) = cetruetab
435            childarray(:,1,2) = pttruetab
436            childarray(:,2,2) = cetruetab
437!cccccccccccccc       memberout = .TRUE.
438#endif
439!
440!           Special values on the child grid
441            if (Agrif_UseSpecialValueFineGrid) then
442                call GiveAgrif_SpecialValueToTab_mpi( child, tempC, childarray, Agrif_SpecialValueFineGrid,nbdim )
443            endif
444!
445        endif   ! ( memberin )
446!
447        if (torestore) then
448!
449#if defined AGRIF_MPI
450!
451            SELECT CASE (nbdim)
452            CASE (1)
453                do i = pttruetab(1),cetruetab(1)
454!hildarrayAModifier     if (restore%restore1D(i) == 0)      &
455!hildarrayAModifier         child%array1(childarray(i,1,2)) = tempC%array1(i)
456                enddo
457            CASE (2)
458                do i = pttruetab(1),cetruetab(1)
459                do j = pttruetab(2),cetruetab(2)
460!hildarrayAModifier     if (restore%restore2D(i,j) == 0)    &
461!hildarrayAModifier         child%array2(childarray(i,1,2), &
462!hildarrayAModifier                          childarray(j,2,2)) = tempC%array2(i,j)
463                enddo
464                enddo
465            CASE (3)
466                do i = pttruetab(1),cetruetab(1)
467                do j = pttruetab(2),cetruetab(2)
468                do k = pttruetab(3),cetruetab(3)
469!hildarrayAModifier     if (restore%restore3D(i,j,k) == 0)  &
470!hildarrayAModifier         child%array3(childarray(i,1,2), &
471!hildarrayAModifier                          childarray(j,2,2), &
472!hildarrayAModifier                          childarray(k,3,2)) = tempC%array3(i,j,k)
473                enddo
474                enddo
475                enddo
476            CASE (4)
477                do i = pttruetab(1),cetruetab(1)
478                do j = pttruetab(2),cetruetab(2)
479                do k = pttruetab(3),cetruetab(3)
480                do l = pttruetab(4),cetruetab(4)
481!hildarrayAModifier     if (restore%restore4D(i,j,k,l) == 0)    &
482!hildarrayAModifier         child%array4(childarray(i,1,2),     &
483!hildarrayAModifier                          childarray(j,2,2),     &
484!hildarrayAModifier                          childarray(k,3,2),     &
485!hildarrayAModifier                          childarray(l,4,2)) = tempC%array4(i,j,k,l)
486                enddo
487                enddo
488                enddo
489                enddo
490            CASE (5)
491                do i = pttruetab(1),cetruetab(1)
492                do j = pttruetab(2),cetruetab(2)
493                do k = pttruetab(3),cetruetab(3)
494                do l = pttruetab(4),cetruetab(4)
495                do m = pttruetab(5),cetruetab(5)
496!hildarrayAModifier     if (restore%restore5D(i,j,k,l,m) == 0)  &
497!hildarrayAModifier         child%array5(childarray(i,1,2),     &
498!hildarrayAModifier                          childarray(j,2,2),     &
499!hildarrayAModifier                          childarray(k,3,2),     &
500!hildarrayAModifier                          childarray(l,4,2),     &
501!hildarrayAModifier                          childarray(m,5,2)) = tempC%array5(i,j,k,l,m)
502                enddo
503                enddo
504                enddo
505                enddo
506                enddo
507            CASE (6)
508                do i = pttruetab(1),cetruetab(1)
509                do j = pttruetab(2),cetruetab(2)
510                do k = pttruetab(3),cetruetab(3)
511                do l = pttruetab(4),cetruetab(4)
512                do m = pttruetab(5),cetruetab(5)
513                do n = pttruetab(6),cetruetab(6)
514!hildarrayAModifier     if (restore%restore6D(i,j,k,l,m,n) == 0)    &
515!hildarrayAModifier         child%array6(childarray(i,1,2),         &
516!hildarrayAModifier                          childarray(j,2,2),         &
517!hildarrayAModifier                          childarray(k,3,2),         &
518!hildarrayAModifier                          childarray(l,4,2),         &
519!hildarrayAModifier                          childarray(m,5,2),         &
520!hildarrayAModifier                          childarray(n,6,2)) = tempC%array6(i,j,k,l,m,n)
521                enddo
522                enddo
523                enddo
524                enddo
525                enddo
526                enddo
527            END SELECT
528!
529#else
530            select case (nbdim)
531            case (1)
532                do i = pttruetab(1),cetruetab(1)
533                    if (restore%restore1D(i) == 0)          &
534                        parray1(i) = tempC % array1(i)
535                enddo
536            case (2)
537                do j = pttruetab(2),cetruetab(2)
538                do i = pttruetab(1),cetruetab(1)
539                    if (restore%restore2D(i,j) == 0)        &
540                        parray2(i,j) = tempC % array2(i,j)
541                enddo
542                enddo
543            case (3)
544                do k = pttruetab(3),cetruetab(3)
545                do j = pttruetab(2),cetruetab(2)
546                do i = pttruetab(1),cetruetab(1)
547                    if (restore%restore3D(i,j,k) == 0)      &
548                        parray3(i,j,k) = tempC % array3(i,j,k)
549                enddo
550                enddo
551                enddo
552            case (4)
553                do l = pttruetab(4),cetruetab(4)
554                do k = pttruetab(3),cetruetab(3)
555                do j = pttruetab(2),cetruetab(2)
556                do i = pttruetab(1),cetruetab(1)
557                    if (restore%restore4D(i,j,k,l) == 0)    &
558                        parray4(i,j,k,l) = tempC % array4(i,j,k,l)
559                enddo
560                enddo
561                enddo
562                enddo
563            case (5)
564                do m = pttruetab(5),cetruetab(5)
565                do l = pttruetab(4),cetruetab(4)
566                do k = pttruetab(3),cetruetab(3)
567                do j = pttruetab(2),cetruetab(2)
568                do i = pttruetab(1),cetruetab(1)
569                    if (restore%restore5D(i,j,k,l,m) == 0)  &
570                        parray5(i,j,k,l,m) = tempC % array5(i,j,k,l,m)
571                enddo
572                enddo
573                enddo
574                enddo
575                enddo
576            case (6)
577                do n = pttruetab(6),cetruetab(6)
578                do m = pttruetab(5),cetruetab(5)
579                do l = pttruetab(4),cetruetab(4)
580                do k = pttruetab(3),cetruetab(3)
581                do j = pttruetab(2),cetruetab(2)
582                do i = pttruetab(1),cetruetab(1)
583                    if (restore%restore6D(i,j,k,l,m,n) == 0)    &
584                        parray6(i,j,k,l,m,n) = tempC % array6(i,j,k,l,m,n)
585                enddo
586                enddo
587                enddo
588                enddo
589                enddo
590                enddo
591            end select
592!
593#endif
594!
595        else    ! .not.to_restore
596!
597            if (memberin) then
598    !
599                if ( .not.in_bc ) then
600                    select case(nbdim)
601                    case(1)
602                        call procname(tempC % array1(            &
603                                childarray(1,1,1):childarray(1,2,1)), &
604                                childarray(1,1,2),childarray(1,2,2),.FALSE.,nb,ndir)
605                    case(2)
606                        call procname(                                &
607                                tempC % array2(            &
608                                childarray(1,1,1):childarray(1,2,1),  &
609                                childarray(2,1,1):childarray(2,2,1)), &
610                                childarray(1,1,2),childarray(1,2,2),  &
611                                childarray(2,1,2),childarray(2,2,2),.FALSE.,nb,ndir)
612                    case(3)
613                        call procname(                                &
614                                tempC % array3(            &
615                                childarray(1,1,1):childarray(1,2,1),  &
616                                childarray(2,1,1):childarray(2,2,1),  &
617                                childarray(3,1,1):childarray(3,2,1)), &
618                                childarray(1,1,2),childarray(1,2,2),  &
619                                childarray(2,1,2),childarray(2,2,2),  &
620                                childarray(3,1,2),childarray(3,2,2),.FALSE.,nb,ndir)
621                    case(4)
622                        call procname(                                &
623                                tempC % array4(            &
624                                childarray(1,1,1):childarray(1,2,1),  &
625                                childarray(2,1,1):childarray(2,2,1),  &
626                                childarray(3,1,1):childarray(3,2,1),  &
627                                childarray(4,1,1):childarray(4,2,1)), &
628                                childarray(1,1,2),childarray(1,2,2),  &
629                                childarray(2,1,2),childarray(2,2,2),  &
630                                childarray(3,1,2),childarray(3,2,2),  &
631                                childarray(4,1,2),childarray(4,2,2),.FALSE.,nb,ndir)
632                    case(5)
633                        call procname(                                &
634                                tempC % array5(            &
635                                childarray(1,1,1):childarray(1,2,1),  &
636                                childarray(2,1,1):childarray(2,2,1),  &
637                                childarray(3,1,1):childarray(3,2,1),  &
638                                childarray(4,1,1):childarray(4,2,1),  &
639                                childarray(5,1,1):childarray(5,2,1)), &
640                                childarray(1,1,2),childarray(1,2,2),  &
641                                childarray(2,1,2),childarray(2,2,2),  &
642                                childarray(3,1,2),childarray(3,2,2),  &
643                                childarray(4,1,2),childarray(4,2,2),  &
644                                childarray(5,1,2),childarray(5,2,2),.FALSE.,nb,ndir)
645                    case(6)
646                        call procname(                                &
647                                tempC % array6(            &
648                                childarray(1,1,1):childarray(1,2,1),  &
649                                childarray(2,1,1):childarray(2,2,1),  &
650                                childarray(3,1,1):childarray(3,2,1),  &
651                                childarray(4,1,1):childarray(4,2,1),  &
652                                childarray(5,1,1):childarray(5,2,1),  &
653                                childarray(6,1,1):childarray(6,2,1)), &
654                                childarray(1,1,2),childarray(1,2,2),  &
655                                childarray(2,1,2),childarray(2,2,2),  &
656                                childarray(3,1,2),childarray(3,2,2),  &
657                                childarray(4,1,2),childarray(4,2,2),  &
658                                childarray(5,1,2),childarray(5,2,2),  &
659                                childarray(6,1,2),childarray(6,2,2),.FALSE.,nb,ndir)
660                    end select
661                else    ! we are in_bc
662                    select case (nbdim)
663                    case (1)
664                        parray1(childarray(1,1,2):childarray(1,2,2)) =    &
665                         tempC%array1(childarray(1,1,1):childarray(1,2,1))
666                    case (2)
667                        parray2(childarray(1,1,2):childarray(1,2,2),      &
668                                      childarray(2,1,2):childarray(2,2,2)) =    &
669                         tempC%array2(childarray(1,1,1):childarray(1,2,1),      &
670                                      childarray(2,1,1):childarray(2,2,1))
671                    case (3)
672                        parray3(childarray(1,1,2):childarray(1,2,2),      &
673                                      childarray(2,1,2):childarray(2,2,2),      &
674                                      childarray(3,1,2):childarray(3,2,2)) =    &
675                         tempC%array3(childarray(1,1,1):childarray(1,2,1),      &
676                                      childarray(2,1,1):childarray(2,2,1),      &
677                                      childarray(3,1,1):childarray(3,2,1))
678                    case (4)
679                        parray4(childarray(1,1,2):childarray(1,2,2),      &
680                                      childarray(2,1,2):childarray(2,2,2),      &
681                                      childarray(3,1,2):childarray(3,2,2),      &
682                                      childarray(4,1,2):childarray(4,2,2)) =    &
683                         tempC%array4(childarray(1,1,1):childarray(1,2,1),      &
684                                      childarray(2,1,1):childarray(2,2,1),      &
685                                      childarray(3,1,1):childarray(3,2,1),      &
686                                      childarray(4,1,1):childarray(4,2,1))
687                    case (5)
688                        parray5(childarray(1,1,2):childarray(1,2,2),      &
689                                      childarray(2,1,2):childarray(2,2,2),      &
690                                      childarray(3,1,2):childarray(3,2,2),      &
691                                      childarray(4,1,2):childarray(4,2,2),      &
692                                      childarray(5,1,2):childarray(5,2,2)) =    &
693                         tempC%array5(childarray(1,1,1):childarray(1,2,1),      &
694                                      childarray(2,1,1):childarray(2,2,1),      &
695                                      childarray(3,1,1):childarray(3,2,1),      &
696                                      childarray(4,1,1):childarray(4,2,1),      &
697                                      childarray(5,1,1):childarray(5,2,1))
698                    case (6)
699                        parray6(childarray(1,1,2):childarray(1,2,2),      &
700                                      childarray(2,1,2):childarray(2,2,2),      &
701                                      childarray(3,1,2):childarray(3,2,2),      &
702                                      childarray(4,1,2):childarray(4,2,2),      &
703                                      childarray(5,1,2):childarray(5,2,2),      &
704                                      childarray(6,1,2):childarray(6,2,2)) =    &
705                         tempC%array6(childarray(1,1,1):childarray(1,2,1),      &
706                                      childarray(2,1,1):childarray(2,2,1),      &
707                                      childarray(3,1,1):childarray(3,2,1),      &
708                                      childarray(4,1,1):childarray(4,2,1),      &
709                                      childarray(5,1,1):childarray(5,2,1),      &
710                                      childarray(6,1,1):childarray(6,2,1))
711                    end select
712                endif  ! < (.not.in_bc)
713            endif  ! < memberin
714!
715        endif
716
717        call Agrif_array_deallocate(tempPextend,nbdim)
718        call Agrif_array_deallocate(tempC,nbdim)
719
720    endif
721!
722!   Deallocations
723#if defined AGRIF_MPI
724    if (member) then
725        call Agrif_array_deallocate(tempP,nbdim)
726    endif
727#endif
728!---------------------------------------------------------------------------------------------------
729end subroutine Agrif_InterpnD
730!===================================================================================================
731!
732!===================================================================================================
733!  subroutine Agrif_Parentbounds
734!
735!> Calculates the bounds of the parent grid for the interpolation of the child grid
736!---------------------------------------------------------------------------------------------------
737subroutine Agrif_Parentbounds ( type_interp, nbdim, indmin, indmax, &
738                                s_Parent_temp, s_Child_temp,        &
739                                s_Child, ds_Child,                  &
740                                s_Parent,ds_Parent,                 &
741                                pttruetab, cetruetab,               &
742                                pttab_Child, pttab_Parent, posvar, coords )
743!---------------------------------------------------------------------------------------------------
744    INTEGER, DIMENSION(6),     intent(in)  :: type_interp
745    INTEGER,                   intent(in)  :: nbdim
746    INTEGER, DIMENSION(nbdim), intent(out) :: indmin, indmax
747    REAL,    DIMENSION(nbdim), intent(out) :: s_Parent_temp, s_child_temp
748    REAL,    DIMENSION(nbdim), intent(in)  :: s_Child, ds_child
749    REAL,    DIMENSION(nbdim), intent(in)  :: s_Parent,ds_Parent
750    INTEGER, DIMENSION(nbdim), intent(in)  :: pttruetab, cetruetab
751    INTEGER, DIMENSION(nbdim), intent(in)  :: pttab_Child, pttab_Parent
752    INTEGER, DIMENSION(nbdim), intent(in)  :: posvar
753    INTEGER, DIMENSION(nbdim), intent(in)  :: coords
754!
755    INTEGER :: i
756    REAL,DIMENSION(nbdim) :: dim_newmin, dim_newmax
757!
758    dim_newmin = s_Child + (pttruetab - pttab_Child) * ds_Child
759    dim_newmax = s_Child + (cetruetab - pttab_Child) * ds_Child
760!
761    do i = 1,nbdim
762!
763        indmin(i) = pttab_Parent(i) + agrif_int((dim_newmin(i)-s_Parent(i))/ds_Parent(i))
764        indmax(i) = pttab_Parent(i) + agrif_ceiling((dim_newmax(i)-s_Parent(i))/ds_Parent(i))
765!
766!       Necessary for the Quadratic interpolation
767!
768        if ( (pttruetab(i) == cetruetab(i)) .and. (posvar(i) == 1) ) then
769        elseif ( coords(i) == 0 ) then  ! (interptab == 'N')
770        elseif ( (type_interp(i) == Agrif_ppm)     .or.     &
771                 (type_interp(i) == Agrif_eno)     .or.     &
772                 (type_interp(i) == Agrif_ppm_lim) .or.     &
773                 (type_interp(i) == Agrif_weno) ) then
774            indmin(i) = indmin(i) - 2
775            indmax(i) = indmax(i) + 2
776        elseif ( (type_interp(i) /= Agrif_constant) .and.   &
777                 (type_interp(i) /= Agrif_linear) ) then
778            indmin(i) = indmin(i) - 1
779            indmax(i) = indmax(i) + 1
780        endif
781!
782    enddo
783!
784    s_Parent_temp = s_Parent + (indmin - pttab_Parent) * ds_Parent
785    s_Child_temp  = s_Child + (pttruetab - pttab_Child) * ds_Child
786!---------------------------------------------------------------------------------------------------
787end subroutine Agrif_Parentbounds
788!===================================================================================================
789!
790!===================================================================================================
791!  subroutine Agrif_Interp_1D_Recursive
792!
793!> Subroutine for the interpolation of a 1D grid variable.
794!> It calls Agrif_InterpBase.
795!---------------------------------------------------------------------------------------------------
796subroutine Agrif_Interp_1D_recursive ( type_interp, tabin, tabout,  &
797                                       indmin, indmax,              &
798                                       pttab_child, petab_child,    &
799                                       s_child,  s_parent,          &
800                                       ds_child, ds_parent )
801!---------------------------------------------------------------------------------------------------
802    integer,            intent(in)  :: type_interp
803    integer,            intent(in)  :: indmin, indmax
804    integer,            intent(in)  :: pttab_child, petab_child
805    real,               intent(in)  :: s_child, s_parent
806    real,               intent(in)  :: ds_child, ds_parent
807    real, dimension(            &
808        indmin:indmax           &
809    ),                  intent(in)  :: tabin
810    real, dimension(            &
811        pttab_child:petab_child &
812    ),                  intent(out) :: tabout
813!---------------------------------------------------------------------------------------------------
814    call Agrif_InterpBase(type_interp,                      &
815                          tabin(indmin:indmax),             &
816                          tabout(pttab_child:petab_child),  &
817                          indmin, indmax,                   &
818                          pttab_child, petab_child,         &
819                          s_parent,    s_child,             &
820                         ds_parent,   ds_child)
821!---------------------------------------------------------------------------------------------------
822end subroutine Agrif_Interp_1D_recursive
823!===================================================================================================
824!
825!===================================================================================================
826!  subroutine Agrif_Interp_2D_Recursive
827!
828!> Subroutine for the interpolation of a 2D grid variable.
829!> It calls Agrif_Interp_1D_recursive and Agrif_InterpBase.
830!---------------------------------------------------------------------------------------------------
831subroutine Agrif_Interp_2D_recursive ( type_interp, tabin, tabout,  &
832                                       indmin, indmax,              &
833                                       pttab_child, petab_child,    &
834                                       s_child,  s_parent,          &
835                                       ds_child, ds_parent )
836!---------------------------------------------------------------------------------------------------
837    integer, dimension(2),              intent(in)  :: type_interp
838    integer, dimension(2),              intent(in)  :: indmin, indmax
839    integer, dimension(2),              intent(in)  :: pttab_child, petab_child
840    real,    dimension(2),              intent(in)  :: s_child, s_parent
841    real,    dimension(2),              intent(in)  :: ds_child, ds_parent
842    real,    dimension(                 &
843        indmin(1):indmax(1),            &
844        indmin(2):indmax(2)),           intent(in)  :: tabin
845    real,    dimension(                 &
846        pttab_child(1):petab_child(1),  &
847        pttab_child(2):petab_child(2)), intent(out) :: tabout
848!---------------------------------------------------------------------------------------------------
849    real, dimension(                    &
850        pttab_child(1):petab_child(1),  &
851        indmin(2):indmax(2))            :: tabtemp
852    real, dimension(                    &
853        pttab_child(2):petab_child(2),  &
854        pttab_child(1):petab_child(1))  :: tabout_trsp
855    real, dimension(                    &
856        indmin(2):indmax(2),            &
857        pttab_child(1):petab_child(1))  :: tabtemp_trsp
858    integer                             :: i, j, coeffraf
859!---------------------------------------------------------------------------------------------------
860!
861    coeffraf = nint ( ds_parent(1) / ds_child(1) )
862!
863    if ((type_interp(1) == Agrif_Linear) .and. (coeffraf /= 1)) then
864!---CDIR NEXPAND
865        if(.NOT. precomputedone(1))     &
866            call Linear1dPrecompute2d(                  &
867                    indmax(2)-indmin(2)+1,              &
868                    indmax(1)-indmin(1)+1,              &
869                    petab_child(1)-pttab_child(1)+1,    &
870                    s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
871!---CDIR NEXPAND
872        call Linear1dAfterCompute(tabin,tabtemp,size(tabin),size(tabtemp),1)
873!
874    elseif ((type_interp(1) == Agrif_PPM) .and. (coeffraf /= 1)) then
875!---CDIR NEXPAND
876        if(.NOT. precomputedone(1))     &
877            call PPM1dPrecompute2d(                     &
878                    indmax(2)-indmin(2)+1,              &
879                    indmax(1)-indmin(1)+1,              &
880                    petab_child(1)-pttab_child(1)+1,        &
881                    s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
882!---CDIR NEXPAND
883        call PPM1dAfterCompute(tabin,tabtemp,size(tabin),size(tabtemp),1)
884    else
885        do j = indmin(2),indmax(2)
886!
887!---CDIR NEXPAND
888            call Agrif_Interp_1D_recursive(type_interp(1),                  &
889                    tabin(indmin(1):indmax(1),j),               &
890                    tabtemp(pttab_child(1):petab_child(1),j),   &
891                    indmin(1),indmax(1),                    &
892                    pttab_child(1),petab_child(1),          &
893                    s_child(1), s_parent(1),                &
894                    ds_child(1),ds_parent(1))
895!
896        enddo
897    endif
898
899    coeffraf = nint(ds_parent(2)/ds_child(2))
900    tabtemp_trsp = TRANSPOSE(tabtemp)
901
902    if ((type_interp(2) == Agrif_Linear) .and. (coeffraf /= 1)) then
903!---CDIR NEXPAND
904        if(.NOT. precomputedone(2))     &
905            call Linear1dPrecompute2d(                  &
906                    petab_child(1)-pttab_child(1)+1,    &
907                    indmax(2)-indmin(2)+1,              &
908                    petab_child(2)-pttab_child(2)+1,    &
909                    s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
910!---CDIR NEXPAND
911        call Linear1dAfterCompute(tabtemp_trsp,tabout_trsp, &
912                size(tabtemp_trsp),size(tabout_trsp),2)
913
914    elseif ((type_interp(2) == Agrif_PPM) .and. (coeffraf /= 1)) then
915!---CDIR NEXPAND
916        if(.NOT. precomputedone(2))     &
917            call PPM1dPrecompute2d(                     &
918                    petab_child(1)-pttab_child(1)+1,    &
919                    indmax(2)-indmin(2)+1,              &
920                    petab_child(2)-pttab_child(2)+1,    &
921                    s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
922!---CDIR NEXPAND
923        call PPM1dAfterCompute(tabtemp_trsp, tabout_trsp,    &
924                               size(tabtemp_trsp), size(tabout_trsp), 2)
925    else
926        do i = pttab_child(1), petab_child(1)
927!
928!---CDIR NEXPAND
929            call Agrif_InterpBase(type_interp(2),                                   &
930                                  tabtemp_trsp(indmin(2):indmax(2), i),             &
931                                  tabout_trsp(pttab_child(2):petab_child(2), i),    &
932                                  indmin(2), indmax(2),                             &
933                                  pttab_child(2), petab_child(2),                   &
934                                  s_parent(2),    s_child(2),                       &
935                                 ds_parent(2),   ds_child(2) )
936        enddo
937    endif
938!
939    tabout = TRANSPOSE(tabout_trsp)
940!---------------------------------------------------------------------------------------------------
941end subroutine Agrif_Interp_2D_recursive
942!===================================================================================================
943!
944!===================================================================================================
945!  subroutine Agrif_Interp_3D_Recursive
946!
947!> Subroutine for the interpolation of a 3D grid variable.
948!> It calls #Agrif_Interp_2D_recursive and #Agrif_InterpBase.
949!---------------------------------------------------------------------------------------------------
950subroutine Agrif_Interp_3D_recursive ( type_interp, tabin, tabout,  &
951                                       indmin, indmax,              &
952                                       pttab_child, petab_child,    &
953                                       s_child,  s_parent,          &
954                                       ds_child, ds_parent )
955!---------------------------------------------------------------------------------------------------
956    integer, dimension(3),              intent(in)  :: type_interp
957    integer, dimension(3),              intent(in)  :: indmin, indmax
958    integer, dimension(3),              intent(in)  :: pttab_child, petab_child
959    real,    dimension(3),              intent(in)  :: s_child, s_parent
960    real,    dimension(3),              intent(in)  :: ds_child, ds_parent
961    real,    dimension(                 &
962        indmin(1):indmax(1),            &
963        indmin(2):indmax(2),            &
964        indmin(3):indmax(3)),           intent(in)  :: tabin
965    real,    dimension(                 &
966        pttab_child(1):petab_child(1),  &
967        pttab_child(2):petab_child(2),  &
968        pttab_child(3):petab_child(3)), intent(out) :: tabout
969!---------------------------------------------------------------------------------------------------
970    real, dimension(                    &
971        pttab_child(1):petab_child(1),  &
972        pttab_child(2):petab_child(2),  &
973        indmin(3):indmax(3))            :: tabtemp
974    integer                             :: i, j, k, coeffraf
975    integer                             :: locind_child_left, kdeb
976!
977    coeffraf = nint ( ds_parent(1) / ds_child(1) )
978    if ( (type_interp(1) == Agrif_Linear) .and. (coeffraf/=1) ) then
979        call Linear1dPrecompute2d(indmax(2)-indmin(2)+1,            &
980                                  indmax(1)-indmin(1)+1,            &
981                                  petab_child(1)-pttab_child(1)+1,  &
982                                  s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
983        precomputedone(1) = .TRUE.
984    elseif ( (type_interp(1) == Agrif_PPM) .and. (coeffraf/=1) ) then
985        call PPM1dPrecompute2d(indmax(2)-indmin(2)+1,           &
986                               indmax(1)-indmin(1)+1,           &
987                               petab_child(1)-pttab_child(1)+1, &
988                               s_parent(1),s_child(1),ds_parent(1),ds_child(1),1)
989        precomputedone(1) = .TRUE.
990    endif
991
992    coeffraf = nint ( ds_parent(2) / ds_child(2) )
993    if ( (type_interp(2) == Agrif_Linear) .and. (coeffraf/=1) ) then
994        call Linear1dPrecompute2d(petab_child(1)-pttab_child(1)+1,  &
995                                  indmax(2)-indmin(2)+1,            &
996                                  petab_child(2)-pttab_child(2)+1,  &
997                                  s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
998        precomputedone(2) = .TRUE.
999    elseif ( (type_interp(2) == Agrif_PPM) .and. (coeffraf/=1) ) then
1000        call PPM1dPrecompute2d(petab_child(1)-pttab_child(1)+1, &
1001                               indmax(2)-indmin(2)+1,           &
1002                               petab_child(2)-pttab_child(2)+1, &
1003                               s_parent(2),s_child(2),ds_parent(2),ds_child(2),2)
1004        precomputedone(2) = .TRUE.
1005    endif
1006!
1007    do k = indmin(3), indmax(3)
1008        call Agrif_Interp_2D_recursive(type_interp(1:2),                            &
1009                                       tabin(indmin(1):indmax(1),                   &
1010                                             indmin(2):indmax(2), k),               &
1011                                       tabtemp(pttab_child(1):petab_child(1),       &
1012                                               pttab_child(2):petab_child(2), k),   &
1013                                       indmin(1:2), indmax(1:2),                    &
1014                                       pttab_child(1:2), petab_child(1:2),          &
1015                                       s_child(1:2),     s_parent(1:2),             &
1016                                       ds_child(1:2),   ds_parent(1:2) )
1017    enddo
1018!
1019    precomputedone(1) = .FALSE.
1020    precomputedone(2) = .FALSE.
1021    coeffraf = nint(ds_parent(3)/ds_child(3))
1022!
1023    if ( coeffraf == 1 ) then
1024        locind_child_left = 1 + agrif_int((s_child(3)-s_parent(3))/ds_parent(3))
1025        kdeb = indmin(3)+locind_child_left-2
1026        do k = pttab_child(3),petab_child(3)
1027            kdeb = kdeb + 1
1028            do j = pttab_child(2), petab_child(2)
1029            do i = pttab_child(1), petab_child(1)
1030                tabout(i,j,k) = tabtemp(i,j,kdeb)
1031            enddo
1032            enddo
1033        enddo
1034    else
1035        do j = pttab_child(2), petab_child(2)
1036        do i = pttab_child(1), petab_child(1)
1037            call Agrif_InterpBase(type_interp(3),                                   &
1038                                  tabtemp(i,j,indmin(3):indmax(3)),                 &
1039                                  tabout(i,j,pttab_child(3):petab_child(3)),        &
1040                                  indmin(3), indmax(3),                             &
1041                                  pttab_child(3), petab_child(3),                   &
1042                                  s_parent(3),    s_child(3),                       &
1043                                 ds_parent(3),   ds_child(3) )
1044        enddo
1045        enddo
1046    endif
1047!---------------------------------------------------------------------------------------------------
1048end subroutine Agrif_Interp_3D_recursive
1049!===================================================================================================
1050!
1051!===================================================================================================
1052!  subroutine Agrif_Interp_4D_Recursive
1053!
1054!> Subroutine for the interpolation of a 4D grid variable.
1055!> It calls #Agrif_Interp_3D_recursive and #Agrif_InterpBase.
1056!---------------------------------------------------------------------------------------------------
1057subroutine Agrif_Interp_4D_recursive ( type_interp, tabin, tabout,  &
1058                                       indmin, indmax,              &
1059                                       pttab_child, petab_child,    &
1060                                       s_child,  s_parent,          &
1061                                       ds_child, ds_parent )
1062!---------------------------------------------------------------------------------------------------
1063    integer, dimension(4),              intent(in)  :: type_interp
1064    integer, dimension(4),              intent(in)  :: indmin, indmax
1065    integer, dimension(4),              intent(in)  :: pttab_child, petab_child
1066    real,    dimension(4),              intent(in)  :: s_child, s_parent
1067    real,    dimension(4),              intent(in)  :: ds_child, ds_parent
1068    real,    dimension(                 &
1069        indmin(1):indmax(1),            &
1070        indmin(2):indmax(2),            &
1071        indmin(3):indmax(3),            &
1072        indmin(4):indmax(4)),           intent(in)  :: tabin
1073    real,    dimension(                 &
1074        pttab_child(1):petab_child(1),  &
1075        pttab_child(2):petab_child(2),  &
1076        pttab_child(3):petab_child(3),  &
1077        pttab_child(4):petab_child(4)), intent(out) :: tabout
1078!---------------------------------------------------------------------------------------------------
1079    real, dimension(                    &
1080        pttab_child(1):petab_child(1),  &
1081        pttab_child(2):petab_child(2),  &
1082        pttab_child(3):petab_child(3),  &
1083        indmin(4):indmax(4))            :: tabtemp
1084    integer                             :: i, j, k, l
1085!
1086    do l = indmin(4), indmax(4)
1087        call Agrif_Interp_3D_recursive(type_interp(1:3),                            &
1088                                       tabin(indmin(1):indmax(1),                   &
1089                                             indmin(2):indmax(2),                   &
1090                                             indmin(3):indmax(3), l),               &
1091                                       tabtemp(pttab_child(1):petab_child(1),       &
1092                                               pttab_child(2):petab_child(2),       &
1093                                               pttab_child(3):petab_child(3), l),   &
1094                                       indmin(1:3), indmax(1:3),                    &
1095                                       pttab_child(1:3), petab_child(1:3),          &
1096                                       s_child(1:3),    s_parent(1:3),              &
1097                                       ds_child(1:3),  ds_parent(1:3) )
1098    enddo
1099!
1100    do k = pttab_child(3), petab_child(3)
1101    do j = pttab_child(2), petab_child(2)
1102    do i = pttab_child(1), petab_child(1)
1103        call Agrif_InterpBase(type_interp(4),                                       &
1104                              tabtemp(i,j,k,indmin(4):indmax(4)),                   &
1105                              tabout(i,j,k,pttab_child(4):petab_child(4)),          &
1106                              indmin(4), indmax(4),                                 &
1107                              pttab_child(4), petab_child(4),                       &
1108                              s_parent(4),    s_child(4),                           &
1109                             ds_parent(4),   ds_child(4) )
1110    enddo
1111    enddo
1112    enddo
1113!---------------------------------------------------------------------------------------------------
1114end subroutine Agrif_Interp_4D_recursive
1115!===================================================================================================
1116!
1117!===================================================================================================
1118!  subroutine Agrif_Interp_5D_Recursive
1119!
1120!> Subroutine for the interpolation of a 5D grid variable.
1121!> It calls #Agrif_Interp_4D_recursive and #Agrif_InterpBase.
1122!---------------------------------------------------------------------------------------------------
1123subroutine Agrif_Interp_5D_recursive ( type_interp, tabin, tabout,  &
1124                                       indmin, indmax,              &
1125                                       pttab_child, petab_child,    &
1126                                       s_child,  s_parent,          &
1127                                       ds_child, ds_parent )
1128!---------------------------------------------------------------------------------------------------
1129    integer, dimension(5),              intent(in)  :: type_interp
1130    integer, dimension(5),              intent(in)  :: indmin, indmax
1131    integer, dimension(5),              intent(in)  :: pttab_child, petab_child
1132    real,    dimension(5),              intent(in)  :: s_child, s_parent
1133    real,    dimension(5),              intent(in)  :: ds_child, ds_parent
1134    real,    dimension(                 &
1135        indmin(1):indmax(1),            &
1136        indmin(2):indmax(2),            &
1137        indmin(3):indmax(3),            &
1138        indmin(4):indmax(4),            &
1139        indmin(5):indmax(5)),           intent(in)  :: tabin
1140    real,    dimension(                 &
1141        pttab_child(1):petab_child(1),  &
1142        pttab_child(2):petab_child(2),  &
1143        pttab_child(3):petab_child(3),  &
1144        pttab_child(4):petab_child(4),  &
1145        pttab_child(5):petab_child(5)), intent(out) :: tabout
1146!---------------------------------------------------------------------------------------------------
1147    real, dimension(                    &
1148        pttab_child(1):petab_child(1),  &
1149        pttab_child(2):petab_child(2),  &
1150        pttab_child(3):petab_child(3),  &
1151        pttab_child(4):petab_child(4),  &
1152        indmin(5):indmax(5))            :: tabtemp
1153    integer                             :: i, j, k, l, m
1154!
1155    do m = indmin(5), indmax(5)
1156        call Agrif_Interp_4D_recursive(type_interp(1:4),                            &
1157                                       tabin(indmin(1):indmax(1),                   &
1158                                             indmin(2):indmax(2),                   &
1159                                             indmin(3):indmax(3),                   &
1160                                             indmin(4):indmax(4),m),                &
1161                                       tabtemp(pttab_child(1):petab_child(1),       &
1162                                               pttab_child(2):petab_child(2),       &
1163                                               pttab_child(3):petab_child(3),       &
1164                                               pttab_child(4):petab_child(4), m),   &
1165                                       indmin(1:4),indmax(1:4),                     &
1166                                       pttab_child(1:4), petab_child(1:4),          &
1167                                       s_child(1:4),     s_parent(1:4),             &
1168                                       ds_child(1:4),   ds_parent(1:4) )
1169    enddo
1170!
1171    do l = pttab_child(4), petab_child(4)
1172    do k = pttab_child(3), petab_child(3)
1173    do j = pttab_child(2), petab_child(2)
1174    do i = pttab_child(1), petab_child(1)
1175        call Agrif_InterpBase(type_interp(5),                                       &
1176                              tabtemp(i,j,k,l,indmin(5):indmax(5)),                 &
1177                              tabout(i,j,k,l,pttab_child(5):petab_child(5)),        &
1178                              indmin(5), indmax(5),                                 &
1179                              pttab_child(5), petab_child(5),                       &
1180                              s_parent(5),   s_child(5),                            &
1181                              ds_parent(5), ds_child(5) )
1182    enddo
1183    enddo
1184    enddo
1185    enddo
1186!---------------------------------------------------------------------------------------------------
1187end subroutine Agrif_Interp_5D_recursive
1188!===================================================================================================
1189!
1190!===================================================================================================
1191!  subroutine Agrif_Interp_6D_Recursive
1192!
1193!> Subroutine for the interpolation of a 6D grid variable.
1194!> It calls #Agrif_Interp_5D_recursive and Agrif_InterpBase.
1195!---------------------------------------------------------------------------------------------------
1196subroutine Agrif_Interp_6D_recursive ( type_interp, tabin, tabout,  &
1197                                       indmin, indmax,              &
1198                                       pttab_child, petab_child,    &
1199                                       s_child,  s_parent,          &
1200                                       ds_child, ds_parent )
1201!---------------------------------------------------------------------------------------------------
1202    integer, dimension(6),                  intent(in)  :: type_interp
1203    integer, dimension(6),                  intent(in)  :: indmin, indmax
1204    integer, dimension(6),                  intent(in)  :: pttab_child, petab_child
1205    real,    dimension(6),                  intent(in)  :: s_child, s_parent
1206    real,    dimension(6),                  intent(in)  :: ds_child, ds_parent
1207    real,    dimension(                 &
1208        indmin(1):indmax(1),            &
1209        indmin(2):indmax(2),            &
1210        indmin(3):indmax(3),            &
1211        indmin(4):indmax(4),            &
1212        indmin(5):indmax(5),            &
1213        indmin(6):indmax(6)),               intent(in)  :: tabin
1214    real,    dimension(                 &
1215        pttab_child(1):petab_child(1),  &
1216        pttab_child(2):petab_child(2),  &
1217        pttab_child(3):petab_child(3),  &
1218        pttab_child(4):petab_child(4),  &
1219        pttab_child(5):petab_child(5),  &
1220        pttab_child(6):petab_child(6)),     intent(out) :: tabout
1221!---------------------------------------------------------------------------------------------------
1222    real, dimension(                    &
1223        pttab_child(1):petab_child(1),  &
1224        pttab_child(2):petab_child(2),  &
1225        pttab_child(3):petab_child(3),  &
1226        pttab_child(4):petab_child(4),  &
1227        pttab_child(5):petab_child(5),  &
1228        indmin(6):indmax(6))            :: tabtemp
1229    integer                             :: i, j, k, l, m, n
1230!
1231    do n = indmin(6), indmax(6)
1232        call Agrif_Interp_5D_recursive(type_interp(1:5),                            &
1233                                       tabin(indmin(1):indmax(1),                   &
1234                                             indmin(2):indmax(2),                   &
1235                                             indmin(3):indmax(3),                   &
1236                                             indmin(4):indmax(4),                   &
1237                                             indmin(5):indmax(5), n),               &
1238                                        tabtemp(pttab_child(1):petab_child(1),      &
1239                                                pttab_child(2):petab_child(2),      &
1240                                                pttab_child(3):petab_child(3),      &
1241                                                pttab_child(4):petab_child(4),      &
1242                                                pttab_child(5):petab_child(5), n),  &
1243                                        indmin(1:5),indmax(1:5),                    &
1244                                        pttab_child(1:5), petab_child(1:5),         &
1245                                        s_child(1:5), s_parent(1:5),                &
1246                                        ds_child(1:5),ds_parent(1:5) )
1247    enddo
1248!
1249    do m = pttab_child(5), petab_child(5)
1250    do l = pttab_child(4), petab_child(4)
1251    do k = pttab_child(3), petab_child(3)
1252    do j = pttab_child(2), petab_child(2)
1253    do i = pttab_child(1), petab_child(1)
1254        call Agrif_InterpBase(type_interp(6),                                       &
1255                              tabtemp(i,j,k,l,m,indmin(6):indmax(6)),               &
1256                              tabout(i,j,k,l,m,pttab_child(6):petab_child(6)),      &
1257                              indmin(6), indmax(6),                                 &
1258                              pttab_child(6), petab_child(6),                       &
1259                              s_parent(6),   s_child(6),                            &
1260                              ds_parent(6), ds_child(6) )
1261    enddo
1262    enddo
1263    enddo
1264    enddo
1265    enddo
1266!---------------------------------------------------------------------------------------------------
1267end subroutine Agrif_Interp_6D_recursive
1268!===================================================================================================
1269!
1270!===================================================================================================
1271!  subroutine Agrif_InterpBase
1272!
1273!> Calls the interpolation method chosen by the user (linear, lagrange, spline, etc.).
1274!---------------------------------------------------------------------------------------------------
1275subroutine Agrif_InterpBase ( type_interp, parenttab, childtab, indmin, indmax, &
1276                              pttab_child, petab_child,                         &
1277                              s_parent, s_child, ds_parent, ds_child )
1278!---------------------------------------------------------------------------------------------------
1279    INTEGER                                                 :: type_interp
1280    INTEGER                                                 :: indmin, indmax
1281    INTEGER                                                 :: pttab_child, petab_child
1282    REAL, DIMENSION(indmin:indmax),           INTENT(IN)    :: parenttab
1283    REAL, DIMENSION(pttab_child:petab_child), INTENT(OUT)   :: childtab
1284    REAL                                                    :: s_parent, s_child
1285    REAL                                                    :: ds_parent,ds_child
1286!
1287    if ( (indmin == indmax) .and. (pttab_child == petab_child) ) then
1288!
1289        childtab(pttab_child) = parenttab(indmin)
1290!
1291    elseif (type_interp == Agrif_LINEAR) then    !       Linear interpolation
1292!
1293        call Agrif_basicinterp_linear1D(parenttab,childtab,                   &
1294                indmax-indmin+1,petab_child-pttab_child+1,  &
1295                s_parent,s_child,ds_parent,ds_child)
1296!
1297    elseif ( type_interp == Agrif_PPM ) then     !       PPM interpolation
1298
1299        call PPM1d(parenttab,childtab,                      &
1300                indmax-indmin+1,petab_child-pttab_child+1,  &
1301                s_parent,s_child,ds_parent,ds_child)
1302!
1303    elseif ( type_interp == Agrif_PPM_LIM ) then     !       PPM interpolation
1304
1305        call PPM1d_lim(parenttab,childtab,                      &
1306                indmax-indmin+1,petab_child-pttab_child+1,  &
1307                s_parent,s_child,ds_parent,ds_child)
1308!
1309    elseif (type_interp == Agrif_LAGRANGE) then  !       Lagrange interpolation
1310!
1311        call lagrange1D(parenttab,childtab,                 &
1312                indmax-indmin+1,petab_child-pttab_child+1,  &
1313                s_parent,s_child,ds_parent,ds_child)
1314!
1315    elseif (type_interp == Agrif_ENO) then       !       Eno interpolation
1316!
1317        call ENO1d(parenttab,childtab,                      &
1318                indmax-indmin+1,petab_child-pttab_child+1,  &
1319                s_parent,s_child,ds_parent,ds_child)
1320!
1321    elseif (type_interp == Agrif_WENO) then      !       Weno interpolation
1322!
1323        call WENO1d(parenttab,childtab,                     &
1324                indmax-indmin+1,petab_child-pttab_child+1,  &
1325                s_parent,s_child,ds_parent,ds_child)
1326!
1327    elseif (type_interp == Agrif_LINEARCONSERV) then !   Linear conservative interpolation
1328!
1329        call Linear1dConserv(parenttab,childtab,            &
1330                indmax-indmin+1,petab_child-pttab_child+1,  &
1331                s_parent,s_child,ds_parent,ds_child)
1332!
1333    elseif (type_interp == Agrif_LINEARCONSERVLIM) then !Linear conservative interpolation
1334!
1335        call Linear1dConservLim(parenttab,childtab,         &
1336                indmax-indmin+1,petab_child-pttab_child+1,  &
1337                s_parent,s_child,ds_parent,ds_child)
1338!
1339    elseif (type_interp == Agrif_CONSTANT) then
1340!
1341        call Constant1d(parenttab,childtab,                 &
1342                indmax-indmin+1,petab_child-pttab_child+1,  &
1343                s_parent,s_child,ds_parent,ds_child)
1344!
1345    endif
1346!---------------------------------------------------------------------------------------------------
1347end subroutine Agrif_InterpBase
1348!===================================================================================================
1349!
1350!===================================================================================================
1351!  subroutine Agrif_Find_list_interp
1352!---------------------------------------------------------------------------------------------------
1353function Agrif_Find_list_interp ( list_interp, pttab, petab, pttab_Child, pttab_Parent,     &
1354                                    nbdim, indmin, indmax, indminglob,  indmaxglob,         &
1355                                    pttruetab, cetruetab, memberin                          &
1356#if defined AGRIF_MPI
1357                                   ,indminglob2, indmaxglob2, parentarray,                  &
1358                                    member, tab4t, memberinall, sendtoproc1, recvfromproc1  &
1359#endif
1360    ) result(find_list_interp)
1361!---------------------------------------------------------------------------------------------------
1362    type(Agrif_List_Interp_Loc),   pointer     :: list_interp
1363    integer,                       intent(in)  :: nbdim
1364    integer, dimension(nbdim),     intent(in)  :: pttab, petab, pttab_Child, pttab_Parent
1365    integer, dimension(nbdim),     intent(out) :: indmin, indmax
1366    integer, dimension(nbdim),     intent(out) :: indminglob, indmaxglob
1367    integer, dimension(nbdim),     intent(out) :: pttruetab, cetruetab
1368    logical,                       intent(out) :: memberin
1369#if defined AGRIF_MPI
1370    integer, dimension(nbdim),     intent(out) :: indminglob2, indmaxglob2
1371    integer, dimension(nbdim,2,2), intent(out) :: parentarray
1372    logical,                       intent(out) :: member
1373    integer, dimension(nbdim,0:Agrif_Nbprocs-1,8), intent(out) :: tab4t
1374    logical, dimension(0:Agrif_Nbprocs-1),         intent(out) :: memberinall
1375    logical, dimension(0:Agrif_Nbprocs-1),         intent(out) :: sendtoproc1, recvfromproc1
1376#endif
1377    logical :: find_list_interp
1378!
1379    integer :: i
1380    type(Agrif_List_Interp_Loc), pointer :: parcours
1381    type(Agrif_Interp_Loc),      pointer :: pil
1382
1383    find_list_interp = .false.
1384
1385    if ( .not. associated(list_interp) )    return
1386
1387    parcours => list_interp
1388    find_loop : do while ( associated(parcours) )
1389
1390        pil => parcours % interp_loc
1391
1392        do i = 1,nbdim
1393            if ( (pttab(i) /= pil % pttab(i)) .or. &
1394                 (petab(i) /= pil % petab(i)) .or. &
1395                 (pttab_child(i)  /= pil % pttab_child(i)) .or. &
1396                 (pttab_parent(i) /= pil % pttab_parent(i)) ) then
1397                parcours => parcours % suiv
1398                cycle find_loop
1399            endif
1400        enddo
1401
1402        indmin = pil % indmin(1:nbdim)
1403        indmax = pil % indmax(1:nbdim)
1404
1405        pttruetab = pil % pttruetab(1:nbdim)
1406        cetruetab = pil % cetruetab(1:nbdim)
1407
1408#if !defined AGRIF_MPI
1409        indminglob  = pil % indminglob(1:nbdim)
1410        indmaxglob  = pil % indmaxglob(1:nbdim)
1411#else
1412        indminglob  = pil % indminglob2(1:nbdim)
1413        indmaxglob  = pil % indmaxglob2(1:nbdim)
1414        indminglob2 = pil % indminglob2(1:nbdim)
1415        indmaxglob2 = pil % indmaxglob2(1:nbdim)
1416        parentarray = pil % parentarray(1:nbdim,:,:)
1417        member      = pil % member
1418        tab4t         = pil % tab4t(1:nbdim, 0:Agrif_Nbprocs-1, 1:8)
1419        memberinall   = pil % memberinall(0:Agrif_Nbprocs-1)
1420        sendtoproc1   = pil % sendtoproc1(0:Agrif_Nbprocs-1)
1421        recvfromproc1 = pil % recvfromproc1(0:Agrif_Nbprocs-1)
1422#endif
1423        memberin = pil % memberin
1424        find_list_interp = .true.
1425        exit find_loop
1426    enddo find_loop
1427!---------------------------------------------------------------------------------------------------
1428end function Agrif_Find_list_interp
1429!===================================================================================================
1430!
1431!===================================================================================================
1432!  subroutine Agrif_AddTo_list_interp
1433!---------------------------------------------------------------------------------------------------
1434subroutine Agrif_AddTo_list_interp ( list_interp, pttab, petab, pttab_Child, pttab_Parent,  &
1435                                     indmin, indmax, indminglob, indmaxglob,                &
1436                                     pttruetab, cetruetab,                                  &
1437                                     memberin, nbdim                                        &
1438#if defined AGRIF_MPI
1439                                    ,indminglob2, indmaxglob2,                              &
1440                                     parentarray,                                           &
1441                                     member,                                                &
1442                                     tab4t, memberinall, sendtoproc1, recvfromproc1         &
1443#endif
1444    )
1445!---------------------------------------------------------------------------------------------------
1446    type(Agrif_List_Interp_Loc), pointer    :: list_interp
1447    integer                                 :: nbdim
1448    integer, dimension(nbdim)               :: pttab, petab, pttab_Child, pttab_Parent
1449    integer, dimension(nbdim)               :: indmin,indmax
1450    integer, dimension(nbdim)               :: indminglob, indmaxglob
1451    integer, dimension(nbdim)               :: pttruetab, cetruetab
1452    logical                                 :: memberin
1453#if defined AGRIF_MPI
1454    integer, dimension(nbdim,2,2)           :: parentarray
1455    logical                                 :: member
1456    integer, dimension(nbdim)                       :: indminglob2,indmaxglob2
1457    integer, dimension(nbdim,0:Agrif_Nbprocs-1,8)   :: tab4t
1458    logical, dimension(0:Agrif_Nbprocs-1)           :: memberinall
1459    logical, dimension(0:Agrif_Nbprocs-1)           :: sendtoproc1
1460    logical, dimension(0:Agrif_Nbprocs-1)           :: recvfromproc1
1461#endif
1462!
1463    type(Agrif_List_Interp_Loc), pointer    :: parcours
1464    type(Agrif_Interp_Loc),      pointer    :: pil
1465!
1466    allocate(parcours)
1467    allocate(parcours % interp_loc)
1468
1469    pil => parcours % interp_loc
1470
1471    pil % pttab(1:nbdim) = pttab(1:nbdim)
1472    pil % petab(1:nbdim) = petab(1:nbdim)
1473    pil % pttab_child(1:nbdim) = pttab_child(1:nbdim)
1474    pil % pttab_parent(1:nbdim) = pttab_parent(1:nbdim)
1475
1476    pil % indmin(1:nbdim) = indmin(1:nbdim)
1477    pil % indmax(1:nbdim) = indmax(1:nbdim)
1478
1479    pil % memberin = memberin
1480#if !defined AGRIF_MPI
1481    pil % indminglob(1:nbdim) = indminglob(1:nbdim)
1482    pil % indmaxglob(1:nbdim) = indmaxglob(1:nbdim)
1483#else
1484    pil % indminglob2(1:nbdim) = indminglob2(1:nbdim)
1485    pil % indmaxglob2(1:nbdim) = indmaxglob2(1:nbdim)
1486    pil % parentarray(1:nbdim,:,:) = parentarray(1:nbdim,:,:)
1487    pil % member = member
1488    allocate(pil % tab4t(nbdim, 0:Agrif_Nbprocs-1, 8))
1489    allocate(pil % memberinall(0:Agrif_Nbprocs-1))
1490    allocate(pil % sendtoproc1(0:Agrif_Nbprocs-1))
1491    allocate(pil % recvfromproc1(0:Agrif_Nbprocs-1))
1492    pil % tab4t         = tab4t
1493    pil % memberinall   = memberinall
1494    pil % sendtoproc1   = sendtoproc1
1495    pil % recvfromproc1 = recvfromproc1
1496#endif
1497
1498    pil % pttruetab(1:nbdim) = pttruetab(1:nbdim)
1499    pil % cetruetab(1:nbdim) = cetruetab(1:nbdim)
1500
1501    parcours % suiv => list_interp
1502    list_interp => parcours
1503!---------------------------------------------------------------------------------------------------
1504end subroutine Agrif_Addto_list_interp
1505!===================================================================================================
1506!
1507end module Agrif_Interpolation
Note: See TracBrowser for help on using the repository browser.