New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
modupdate.F90 in vendors/AGRIF/CMEMS_2020/AGRIF_FILES – NEMO

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

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

update AGRIF library

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