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

source: trunk/NEMOGCM/EXTERNAL/AGRIF/AGRIF_FILES/modbc.F @ 2731

Last change on this file since 2731 was 2731, checked in by rblod, 13 years ago

Changes for Agrif in MPI

  • Property svn:keywords set to Id
File size: 46.5 KB
RevLine 
[1901]1!
[2281]2! $Id$
[1901]3!
4C     AGRIF (Adaptive Grid Refinement In Fortran)
5C
6C     Copyright (C) 2003 Laurent Debreu (Laurent.Debreu@imag.fr)
7C                        Christophe Vouland (Christophe.Vouland@imag.fr)   
8C
9C     This program is free software; you can redistribute it and/or modify
10C     it under the terms of the GNU General Public License as published by
11C     the Free Software Foundation; either version 2 of the License, or
12C     (at your option) any later version.
13C
14C     This program is distributed in the hope that it will be useful,
15C     but WITHOUT ANY WARRANTY; without even the implied warranty of
16C     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17C     GNU General Public License for more details.
18C
19C     You should have received a copy of the GNU General Public License
20C     along with this program; if not, write to the Free Software
21C     Foundation, Inc., 59 Temple Place- Suite 330, Boston, MA 02111-1307, USA.
22C
23C
24C
25CCC   Module Agrif_Boundary 
26C 
27      Module Agrif_Boundary 
28C
29CCC   Description:
30CCC   Module to calculate the boundary conditions on the child grids from their
31CCC   parent grids.
32C
33C     Modules used: 
34C   
35      Use Agrif_Interpolation       
36C
37      IMPLICIT NONE
38C       
39      CONTAINS   
40C     Define procedures contained in this module
41C
42C
43     
44C
45C     **************************************************************************
46CCC   Subroutine Agrif_Interp_bc_1d
47C     **************************************************************************
48C 
49      Subroutine Agrif_Interp_bc_1d(TypeInterp,parent,child,tab,deb,fin,
50     &                              weight,pweight,procname)
51C
52CCC   Description:
53CCC   Subroutine to calculate the boundary conditions on a fine grid for a 1D 
54CCC   grid variable.
55C
56C     Declarations:
57C     
58     
59C
60C     Arguments     
61      External :: procname
62      Optional :: procname
63      INTEGER,DIMENSION(6,6) :: TypeInterp    ! TYPE of interpolation (linear,
64                                            ! lagrange, spline, ... )
65      TYPE(AGRIF_PVariable) :: parent       ! Variable on the parent grid
66      TYPE(AGRIF_PVariable) :: child        ! Variable on the child grid
67      TYPE(AGRIF_PVariable) :: childtemp    ! Temporary variable on the child
68                                            ! grid
69      INTEGER :: deb,fin                    ! Positions where interpolations are
70                                            ! done on the fine grid
71      REAL, DIMENSION(
72     &         child%var%lb(1):child%var%ub(1) 
73     &    ), Target :: tab ! Values of the grid variable
74      LOGICAL :: pweight                    ! Indicates if weight is used for 
75                                            ! the temporal interpolation
76      REAL :: weight                        ! Coefficient for the time
77                                            ! interpolation
78C
79C
80C     Definition of a temporary AGRIF_PVariable data TYPE representing the grid
81C     variable.   
82C
83      allocate(childtemp % var) 
84C
85      childtemp % var % root_var => child % var % root_var
86C     
87C     Values of the grid variable
[2715]88      childtemp % var % parray1 => tab 
[1901]89C
90C     Temporary results for the time interpolation before and after the space 
91C     interpolation 
92      childtemp % var % oldvalues2D => child % var % oldvalues2D
93C 
94C     Index indicating if a space interpolation is necessary
95      childtemp % var % interpIndex => child % var % interpIndex       
96      childtemp % var % Interpolationshouldbemade = 
97     &                 child % var % Interpolationshouldbemade   
98      childtemp % var % list_interp => child % var% list_interp 
99     
100      childtemp % var% lb = child % var % lb
101      childtemp % var% ub = child % var % ub         
102C
103C     Call to the procedure for the calculations of the boundary conditions
104      IF (present(procname)) THEN
105      Call Agrif_CorrectVariable
106     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight,procname)
107      ELSE
108      Call Agrif_CorrectVariable
109     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight)
110      ENDIF
111
112C   
113      child % var % oldvalues2D => childtemp % var % oldvalues2D
114      child % var % list_interp => childtemp % var %list_interp
115C         
116      deallocate(childtemp % var)
117C
118C       
119      End Subroutine Agrif_Interp_bc_1D
120C
121     
122C
123C
124C
125C     **************************************************************************
126CCC   Subroutine Agrif_Interp_bc_2d
127C     **************************************************************************
128C 
129      Subroutine Agrif_Interp_bc_2d(TypeInterp,parent,child,tab,deb,fin,
130     &                              weight,pweight,procname)
131C
132CCC   Description:
133CCC   Subroutine to calculate the boundary conditions on a fine grid for a 2D 
134CCC   grid variable.
135C
136C     Declarations:
137C     
138     
139C
140C     Arguments     
141      External :: procname
142      Optional :: procname
143      INTEGER,DIMENSION(6,6) :: TypeInterp    ! TYPE of interpolation (linear,
144                                            ! lagrange, spline, ... )
145      TYPE(AGRIF_PVariable) :: parent       ! Variable on the parent grid
146      TYPE(AGRIF_PVariable) :: child        ! Variable on the child grid
147      TYPE(AGRIF_PVariable) :: childtemp    ! Temporary variable on the child
148                                            ! grid
149      INTEGER :: deb,fin                    ! Positions where interpolations are
150                                            ! done on the fine grid
151      REAL, DIMENSION(
152     &         child%var%lb(1):child%var%ub(1),
153     &         child%var%lb(2):child%var%ub(2) 
154     &    ), Target :: tab ! Values of the grid variable
155      LOGICAL :: pweight                    ! Indicates if weight is used for 
156                                            ! the temporal interpolation
157      REAL :: weight                        ! Coefficient for the time
158                                            ! interpolation
159C
160C
161C     Definition of a temporary AGRIF_PVariable data TYPE representing the grid
162C     variable.   
163C
164      allocate(childtemp % var) 
165C
166      childtemp % var % root_var => child % var % root_var
167C     
168C     Values of the grid variable
[2715]169      childtemp % var % parray2 => tab 
[1901]170C
171C     Temporary results for the time interpolation before and after the space 
172C     interpolation 
173      childtemp % var % oldvalues2D => child % var % oldvalues2D
174C 
175C     Index indicating if a space interpolation is necessary
176      childtemp % var % interpIndex => child % var % interpIndex       
177      childtemp % var % Interpolationshouldbemade = 
178     &                 child % var % Interpolationshouldbemade   
179      childtemp % var % list_interp => child % var% list_interp 
180     
181      childtemp % var% lb = child % var % lb
182      childtemp % var% ub = child % var % ub         
183C
184C     Call to the procedure for the calculations of the boundary conditions
185      IF (present(procname)) THEN
186      Call Agrif_CorrectVariable
187     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight,procname)
188      ELSE
189      Call Agrif_CorrectVariable
190     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight)
191      ENDIF
192
193C   
194      child % var % oldvalues2D => childtemp % var % oldvalues2D
195      child % var % list_interp => childtemp % var %list_interp
196C         
197      deallocate(childtemp % var)
198C
199C       
200      End Subroutine Agrif_Interp_bc_2D
201C
202C
203C
204C     **************************************************************************
205CCC   Subroutine Agrif_Interp_bc_3d
206C     **************************************************************************
207C 
208      Subroutine Agrif_Interp_bc_3d(TypeInterp,parent,child,tab,deb,fin,
209     &                              weight,pweight,procname)           
210C
211CCC   Description:
212CCC   Subroutine to calculate the boundary conditions on a fine grid for a 3D 
213CCC   variable.
214C
215C     Declarations:
216C     
217     
218C
219C     Arguments 
220      External :: procname
221      Optional :: procname     
222      INTEGER,DIMENSION(6,6) :: TypeInterp    ! TYPE of interpolation (linear,
223                                            ! lagrange, spline, ... )
224      TYPE(AGRIF_PVariable) :: parent       ! Variable on the parent grid
225      TYPE(AGRIF_PVariable) :: child        ! Variable on the child grid
226      TYPE(AGRIF_PVariable) :: childtemp    ! Temporary variable on the child
227                                            ! grid
228      INTEGER :: deb,fin                    ! Positions where interpolations
229                                            ! are done on the fine grid
230      REAL, DIMENSION(
231     &         child%var%lb(1):child%var%ub(1),
232     &         child%var%lb(2):child%var%ub(2),
233     &         child%var%lb(3):child%var%ub(3)
234     &         ), Target :: tab ! Values of the grid variable
235      LOGICAL :: pweight                    ! Indicates if weight is used for
236                                            ! the temporal interpolation
237      REAL :: weight                        ! Coefficient for the time
238                                            ! interpolation
239C
240C
241C     Definition of a temporary AGRIF_PVariable data TYPE representing the grid 
242C     variable. 
243C
244      allocate(childtemp % var) 
245C
246      childtemp % var % root_var => child % var % root_var       
247C     
248C     Values of the grid variable
[2715]249      childtemp % var % parray3 => tab
[1901]250C
251C     Temporary results for the time interpolation before and after the space 
252C     interpolation 
253      childtemp % var % oldvalues2D => child % var % oldvalues2D
254C 
255C     Index indicating if a space interpolation is necessary
256      childtemp % var % interpIndex => child % var % interpIndex
257      childtemp % var % Interpolationshouldbemade = 
258     &                 child % var % Interpolationshouldbemade 
259      childtemp % var % list_interp => child % var% list_interp 
260     
261      childtemp % var% lb = child % var % lb
262      childtemp % var% ub = child % var % ub
263C
264C     Call to the procedure for the calculations of the boundary conditions     
265      IF (present(procname)) THEN
266      Call Agrif_CorrectVariable
267     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight,procname)
268      ELSE
269      Call Agrif_CorrectVariable
270     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight)
271      ENDIF
272C
273      child % var % oldvalues2D => childtemp % var % oldvalues2D
274      child % var % list_interp => childtemp % var %list_interp     
275C     
276      deallocate(childtemp % var)
277C
278C       
279      End Subroutine Agrif_Interp_bc_3D               
280C
281C
282C
283C
284C     **************************************************************************
285CCC   Subroutine Agrif_Interp_bc_4d
286C     **************************************************************************
287C 
288      Subroutine Agrif_Interp_bc_4d(TypeInterp,parent,child,tab,deb,fin,
289     &                             weight,pweight,procname)           
290C
291CCC   Description:
292CCC   Subroutine to calculate the boundary conditions on a fine grid for a 4D 
293CCC   grid variable.
294C
295C     Declarations:
296C     
297     
298C
299C     Arguments     
300      External :: procname
301      Optional :: procname     
302      INTEGER,DIMENSION(6,6) :: TypeInterp      ! TYPE of interpolation (linear,
303                                              ! lagrange, spline, ... )
304      TYPE(AGRIF_PVariable) :: parent         ! Variable on the parent grid
305      TYPE(AGRIF_PVariable) :: child          ! Variable on the child grid
306      TYPE(AGRIF_PVariable) :: childtemp      ! Temporary varaiable on the child
307                                              ! grid
308      INTEGER :: deb,fin                      ! Positions where interpolations
309                                              ! are done on the fine grid
310      REAL, DIMENSION(
311     &         child%var%lb(1):child%var%ub(1),
312     &         child%var%lb(2):child%var%ub(2),
313     &         child%var%lb(3):child%var%ub(3),
314     &         child%var%lb(4):child%var%ub(4)
315     &        ), Target :: tab ! Values of the grid variable
316      LOGICAL :: pweight                      ! Indicates if weight is used for
317                                              ! the temporal interpolation
318      REAL :: weight                          ! Coefficient for the time
319                                              ! interpolation
320C
321C
322C     Definition of a temporary AGRIF_PVariable data TYPE representing the grid 
323C     variable. 
324C
325      allocate(childtemp % var) 
326C
327      childtemp % var % root_var => child % var % root_var 
328C     
329C     Values of the grid variable
[2715]330      childtemp % var % parray4 => tab 
[1901]331C
332C     Temporary results for the time interpolation before and after the space 
333C     interpolation       
334      childtemp % var % oldvalues2D => child % var % oldvalues2D
335C 
336C     Index indicating if a space interpolation is necessary
337      childtemp % var % interpIndex => child % var % interpIndex       
338      childtemp % var % Interpolationshouldbemade = 
339     &                 child % var % Interpolationshouldbemade 
340      childtemp % var % list_interp => child % var% list_interp 
341     
342      childtemp % var% lb = child % var % lb
343      childtemp % var% ub = child % var % ub               
344C
345C     Call to the procedure for the calculations of the boundary conditions
346      IF (present(procname)) THEN
347      Call Agrif_CorrectVariable
348     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight,procname)
349      ELSE
350      Call Agrif_CorrectVariable
351     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight)
352      ENDIF
353C
354      child % var % oldvalues2D => childtemp % var % oldvalues2D
355      child % var % list_interp => childtemp % var %list_interp     
356C     
357      deallocate(childtemp % var)
358C
359C       
360      End Subroutine Agrif_Interp_bc_4D
361C
362C
363C
364C     **************************************************************************
365CCC   Subroutine Agrif_Interp_bc_5d
366C     **************************************************************************
367C 
368      Subroutine Agrif_Interp_bc_5d(TypeInterp,parent,child,tab,deb,fin,
369     &                             weight,pweight,procname)           
370C
371CCC   Description:
372CCC   Subroutine to calculate the boundary conditions on a fine grid for a 5D 
373CCC   grid variable.
374C
375C     Declarations:
376C     
377     
378C
379C     Arguments     
380      External :: procname
381      Optional :: procname     
382      INTEGER,DIMENSION(6,6) :: TypeInterp      ! TYPE of interpolation (linear,
383                                              ! lagrange, spline, ... )
384      TYPE(AGRIF_PVariable) :: parent         ! Variable on the parent grid
385      TYPE(AGRIF_PVariable) :: child          ! Variable on the child grid
386      TYPE(AGRIF_PVariable) :: childtemp      ! Temporary varaiable on the child
387                                              ! grid
388      INTEGER :: deb,fin                      ! Positions where interpolations
389                                              ! are done on the fine grid
390      REAL, DIMENSION(
391     &         child%var%lb(1):child%var%ub(1),
392     &         child%var%lb(2):child%var%ub(2),
393     &         child%var%lb(3):child%var%ub(3),
394     &         child%var%lb(4):child%var%ub(4),
395     &         child%var%lb(5):child%var%ub(5)
396     &         ), Target :: tab ! Values of the grid variable
397      LOGICAL :: pweight                      ! Indicates if weight is used for
398                                              ! the temporal interpolation
399      REAL :: weight                          ! Coefficient for the time
400                                              ! interpolation
401C
402C
403C     Definition of a temporary AGRIF_PVariable data TYPE representing the grid 
404C     variable. 
405C
406      allocate(childtemp % var) 
407C
408      childtemp % var % root_var => child % var % root_var 
409C     
410C     Values of the grid variable
[2715]411      childtemp % var % parray5 => tab 
[1901]412C
413C     Temporary results for the time interpolation before and after the space 
414C     interpolation       
415      childtemp % var % oldvalues2D => child % var % oldvalues2D
416C 
417C     Index indicating if a space interpolation is necessary
418      childtemp % var % interpIndex => child % var % interpIndex       
419      childtemp % var % Interpolationshouldbemade = 
420     &                 child % var % Interpolationshouldbemade 
421      childtemp % var % list_interp => child % var% list_interp 
422     
423      childtemp % var% lb = child % var % lb
424      childtemp % var% ub = child % var % ub
425                     
426C
427C     Call to the procedure for the calculations of the boundary conditions 
428      IF (present(procname)) THEN
429      Call Agrif_CorrectVariable
430     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight,procname)
431      ELSE
432      Call Agrif_CorrectVariable
433     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight)
434      ENDIF
435C
436      child % var % oldvalues2D => childtemp % var % oldvalues2D
437      child % var % list_interp => childtemp % var %list_interp     
438C     
439      deallocate(childtemp % var)
440C
441C       
442      End Subroutine Agrif_Interp_bc_5D
443C
444C
445C
446C
447C     **************************************************************************
448CCC   Subroutine Agrif_Interp_bc_6d
449C     **************************************************************************
450C 
451      Subroutine Agrif_Interp_bc_6d(TypeInterp,parent,child,tab,deb,fin,
452     &                             weight,pweight)           
453C
454CCC   Description:
455CCC   Subroutine to calculate the boundary conditions on a fine grid for a 6D 
456CCC   grid variable.
457C
458C     Declarations:
459C     
460     
461C
462C     Arguments     
463      INTEGER,DIMENSION(6,6) :: TypeInterp      ! TYPE of interpolation (linear,
464                                              ! lagrange, spline, ... )
465      TYPE(AGRIF_PVariable) :: parent         ! Variable on the parent grid
466      TYPE(AGRIF_PVariable) :: child          ! Variable on the child grid
467      TYPE(AGRIF_PVariable) :: childtemp      ! Temporary varaiable on the child
468                                              ! grid
469      INTEGER :: deb,fin                      ! Positions where interpolations
470                                              ! are done on the fine grid
471      REAL, DIMENSION(
472     &         child%var%lb(1):child%var%ub(1),
473     &         child%var%lb(2):child%var%ub(2),
474     &         child%var%lb(3):child%var%ub(3),
475     &         child%var%lb(4):child%var%ub(4),
476     &         child%var%lb(5):child%var%ub(5),
477     &         child%var%lb(6):child%var%ub(6)
478     &         ), Target :: tab ! Values of the grid variable
479      LOGICAL :: pweight                      ! Indicates if weight is used for
480                                              ! the temporal interpolation
481      REAL :: weight                          ! Coefficient for the time
482                                              ! interpolation
483C
484C
485C     Definition of a temporary AGRIF_PVariable data TYPE representing the grid 
486C     variable. 
487C
488      allocate(childtemp % var) 
489C
490      childtemp % var % root_var => child % var % root_var 
491C     
492C     Values of the grid variable
[2715]493      childtemp % var % parray6 => tab 
[1901]494C
495C     Temporary results for the time interpolation before and after the space 
496C     interpolation       
497      childtemp % var % oldvalues2D => child % var % oldvalues2D
498C 
499C     Index indicating if a space interpolation is necessary
500      childtemp % var % interpIndex => child % var % interpIndex       
501      childtemp % var % Interpolationshouldbemade = 
502     &                 child % var % Interpolationshouldbemade 
503      childtemp % var % list_interp => child % var% list_interp 
504     
505      childtemp % var% lb = child % var % lb
506      childtemp % var% ub = child % var % ub               
507C
508C     Call to the procedure for the calculations of the boundary conditions
509      Call Agrif_CorrectVariable
510     &     (TypeInterp,parent,childtemp,deb,fin,pweight,weight)
511C
512      child % var % oldvalues2D => childtemp % var % oldvalues2D
513      child % var % list_interp => childtemp % var %list_interp     
514C     
515      deallocate(childtemp % var)
516C
517C       
518      End Subroutine Agrif_Interp_bc_6D
519C
520C
521C     **************************************************************************
522CCC   Subroutine Agrif_CorrectVariable
523C     **************************************************************************
524C
525      Subroutine AGRIF_CorrectVariable(TypeInterp,parent,child,deb,fin,
526     &                                 pweight,weight,procname)
527C
528CCC   Description:
529CCC   Subroutine to calculate the boundary conditions on a fine grid.
530C
531C     Declarations:
532C     
533     
534C
535C     Arguments
536      External :: procname
537      Optional :: procname
538      TYPE(AGRIF_PVariable) :: parent         ! Variable on the parent grid
539      TYPE(AGRIF_PVariable) :: child          ! Variable on the child grid
540      INTEGER,DIMENSION(6,6)  :: TypeInterp     ! TYPE of interpolation
541                                              !    (linear,lagrange,...)
542      INTEGER               :: deb,fin        ! Positions where boundary
543                                              !    conditions are calculated
544      LOGICAL               :: pweight        ! Indicates if weight is used
545                                              !    for the time interpolation
546      REAL                  :: weight         ! Coefficient for the time
547                                              !    interpolation
548C
549C     Local scalars
550      TYPE(Agrif_Grid)    , Pointer :: Agrif_Child_Gr,Agrif_Parent_Gr
551      TYPE(AGRIF_Variable), Pointer :: root   ! Variable on the root grid
552      INTEGER                       :: nbdim  ! Number of dimensions of
553                                              !    the grid variable
554      INTEGER                       :: n
555      INTEGER,DIMENSION(6)          :: pttab_child  ! Index of the first point
556                                                    !    inside the domain for
557                                                    !    the child grid variable
558      INTEGER,DIMENSION(6)          :: pttab_parent ! Index of the first point
559                                                    !    inside the domain for
560                                                    !    the parent grid
561                                                    !    variable
562      INTEGER,DIMENSION(6)          :: nbtab_Child  ! Number of the cells
563      INTEGER,DIMENSION(6)          :: posvartab_Child    ! Position of the
564                                                    !    variable on the cell
565      INTEGER,DIMENSION(6)          :: loctab_Child ! Indicates if the child
566                                                    !    grid has a common
567                                                    !    border with the root
568                                                    !    grid     
569      REAL, DIMENSION(6)            :: s_child,s_parent   ! Positions of the
570                                                    !    parent and child grids
571      REAL, DIMENSION(6)            :: ds_child,ds_parent ! Space steps of the
572                                                    !    parent and child grids
573C
574C
575      loctab_child(:) = 0
576C
577      Agrif_Child_Gr => Agrif_Curgrid
578      Agrif_Parent_Gr => Agrif_Curgrid % parent
579      root => child % var % root_var 
580      nbdim = root % nbdim
581C
582      do n = 1,nbdim
583        posvartab_child(n) = root % posvar(n)     
584      enddo
585C
586C
587      do n = 1,nbdim
588C
589        Select case(root % interptab(n))
590C
591          case('x') ! x DIMENSION
592C
593            nbtab_Child(n) = Agrif_Child_Gr % nb(1)
594            pttab_Child(n) = root % point(1)
595            pttab_Parent(n) = root % point(1)
596            s_Child(n) = Agrif_Child_Gr % Agrif_x(1)
597            s_Parent(n) = Agrif_Parent_Gr % Agrif_x(1)
598            ds_Child(n) = Agrif_Child_Gr % Agrif_d(1)
599            ds_Parent(n) = Agrif_Parent_Gr % Agrif_d(1)
600            if (root % posvar(n) == 2) then
601                s_Child(n) = s_Child(n) + ds_Child(n)/2.
602                s_Parent(n) = s_Parent(n) + ds_Parent(n)/2.
603            endif
604C
605            if (Agrif_CURGRID % NearRootBorder(1)) 
606     &         loctab_child(n) = -1
607            if (Agrif_CURGRID % DistantRootBorder(1))
608     &         loctab_child(n) = -2
609            if ((Agrif_CURGRID % NearRootBorder(1)) .AND. 
610     &          (Agrif_CURGRID % DistantRootBorder(1))) 
611     &         loctab_child(n) = -3
612C
613          case('y') ! y DIMENSION     
614C
615            nbtab_Child(n) = Agrif_Child_Gr % nb(2)
616            pttab_Child(n) = root % point(2)
617            pttab_Parent(n) = root % point(2)
618            s_Child(n) = Agrif_Child_Gr % Agrif_x(2)
619            s_Parent(n) = Agrif_Parent_Gr % Agrif_x(2) 
620            ds_Child(n) = Agrif_Child_Gr % Agrif_d(2)
621            ds_Parent(n) = Agrif_Parent_Gr % Agrif_d(2)
622            if (root % posvar(n) == 2) then       
623                s_Child(n) = s_Child(n) + ds_Child(n)/2.
624                s_Parent(n) = s_Parent(n) + ds_Parent(n)/2.
625            endif       
626C
627            if (Agrif_CURGRID % NearRootBorder(2)) 
628     &         loctab_child(n) = -1
629            if (Agrif_CURGRID % DistantRootBorder(2)) 
630     &         loctab_child(n) = -2
631            if ((Agrif_CURGRID % NearRootBorder(2)) .AND. 
632     &          (Agrif_CURGRID % DistantRootBorder(2))) 
633     &         loctab_child(n) = -3
634C
635          case('z') ! z DIMENSION
636C
637            nbtab_Child(n) = Agrif_Child_Gr % nb(3)
638            pttab_Child(n) = root % point(3)
639            pttab_Parent(n) = root % point(3)
640            s_Child(n) = Agrif_Child_Gr % Agrif_x(3)
641            s_Parent(n) = Agrif_Parent_Gr % Agrif_x(3)
642            ds_Child(n) = Agrif_Child_Gr % Agrif_d(3)
643            ds_Parent(n) = Agrif_Parent_Gr % Agrif_d(3)
644            if (root % posvar(n) == 2) then       
645                s_Child(n) = s_Child(n) + ds_Child(n)/2.
646                s_Parent(n) = s_Parent(n) + ds_Parent(n)/2.
647            endif       
648C
649            if (Agrif_CURGRID % NearRootBorder(3)) 
650     &         loctab_child(n) = -1
651            if (Agrif_CURGRID % DistantRootBorder(3)) 
652     &         loctab_child(n) = -2
653            if ((Agrif_CURGRID % NearRootBorder(3)) .AND. 
654     &          (Agrif_CURGRID % DistantRootBorder(3))) 
655     &         loctab_child(n) = -3
656C
657          case('N') ! No space DIMENSION     
658C
659            nbtab_Child(n) = child % var % ub(n) - child % var % lb(n)
660            pttab_Child(n) = child % var % lb(n)
661C
662C           No interpolation but only a copy of the values of the grid variable
663C     
664            posvartab_child(n) = 1
665            pttab_Parent(n)= pttab_Child(n)
666            s_Child(n)=0.
667            s_Parent(n)=0. 
668            ds_Child(n)=1.
669            ds_Parent(n)=1.
670            loctab_child(n) = -3
671C
672        End select
673C
674      enddo
675C
676         IF (present(procname)) THEN
677          Call AGRIF_CorrectnD
678     &         (TypeInterp,parent,child,deb,fin,pweight,weight,
679     &          pttab_Child(1:nbdim),pttab_Parent(1:nbdim),
680     &          nbtab_Child(1:nbdim),posvartab_Child(1:nbdim),
681     &          loctab_Child(1:nbdim),
682     &          s_Child(1:nbdim),s_Parent(1:nbdim),
683     &          ds_Child(1:nbdim),ds_Parent(1:nbdim),nbdim,procname)
684        ELSE
685          Call AGRIF_CorrectnD
686     &         (TypeInterp,parent,child,deb,fin,pweight,weight,
687     &          pttab_Child(1:nbdim),pttab_Parent(1:nbdim),
688     &          nbtab_Child(1:nbdim),posvartab_Child(1:nbdim),
689     &          loctab_Child(1:nbdim),
690     &          s_Child(1:nbdim),s_Parent(1:nbdim),
691     &          ds_Child(1:nbdim),ds_Parent(1:nbdim),nbdim)
692         ENDIF
693C
694C
695      End subroutine AGRIF_CorrectVariable
696C
697C     **************************************************************************
698CCC   Subroutine Agrif_Correctnd
699C     **************************************************************************
700C
701      Subroutine AGRIF_Correctnd(TypeInterp,parent,child,deb,fin,
702     &                           pweight,weight,
703     &                           pttab_child,pttab_Parent,
704     &                           nbtab_Child,posvartab_Child,
705     &                           loctab_Child,
706     &                           s_Child,s_Parent,
707     &                           ds_Child,ds_Parent,nbdim,procname)
708C
709CCC   Description:
710CCC   Subroutine to calculate the boundary conditions for a nD grid variable on 
711CCC   a fine grid by using a space and time interpolations; it is called by the 
712CCC   Agrif_CorrectVariable procedure.
713C
714C
715C     Declarations:
716C
717     
718C
[2019]719#ifdef key_mpp_mpi
[1901]720C
[2019]721      INCLUDE 'mpif.h'
[1901]722C
723#endif
724C
725C     Arguments
726      External :: procname
727      Optional :: procname
728      INTEGER,DIMENSION(6,6) :: TypeInterp ! TYPE of interpolation (linear,
729                                         !   spline,...) 
730      TYPE(AGRIF_PVariable)    :: parent ! Variable on the parent grid
731      TYPE(AGRIF_PVariable)    :: child  ! Variable on the child grid
732      INTEGER                  :: deb,fin ! Positions where interpolations
733                                         !    are done
734      LOGICAL                  :: pweight ! Indicates if weight is used for
735                                         !    the temporal interpolation
736      REAL                     :: weight ! Coefficient for the temporal
737                                         !    interpolation
738      INTEGER                  :: nbdim  ! Number of dimensions of the grid
739                                         !    variable
740      INTEGER,DIMENSION(nbdim) :: pttab_child ! Index of the first point inside
741                                         !    the domain for the parent
742                                         !    grid variable
743      INTEGER,DIMENSION(nbdim) :: pttab_Parent ! Index of the first point
744                                         !   inside the domain for the
745                                         !   child grid variable
746      INTEGER,DIMENSION(nbdim) :: nbtab_Child ! Number of cells of the child
747                                         !    grid
748      INTEGER,DIMENSION(nbdim) :: posvartab_Child ! Position of the grid
749                                         !    variable (1 or 2)
750      INTEGER,DIMENSION(nbdim) :: loctab_Child ! Indicates if the child
751                                        !    grid has a common border with
752                                        !    the root grid
753      REAL   ,DIMENSION(nbdim) :: s_Child,s_Parent ! Positions of the parent
754                                        !   and child grids
755      REAL   ,DIMENSION(nbdim) :: ds_Child,ds_Parent ! Space steps of the
756                                        !    parent and child grids
757C
758C     Local variables
759      TYPE(AGRIF_PVariable)        :: restore ! Variable on the parent     
760      INTEGER,DIMENSION(nbdim,2)   :: lubglob
761      INTEGER                      :: i         
762      INTEGER                      :: kindex  ! Index used for safeguard
763                                       !    and time interpolation
764      INTEGER,DIMENSION(nbdim,2,2) :: indtab ! Arrays indicating the limits
765                                       !    of the child     
766      INTEGER,DIMENSION(nbdim,2,2) :: indtruetab ! grid variable where
767                                       !   boundary conditions are
768      INTEGER,DIMENSION(nbdim,2,2,nbdim)   :: ptres,ptres2 ! calculated
769      INTEGER                      :: nb,ndir,n,sizetab(1)
770      REAL, DIMENSION(:), Allocatable :: tab ! Array used for the interpolation
771      REAL    :: c1t,c2t               ! Coefficients for the time interpolation
772                                       !    (c2t=1-c1t)
773C
[2019]774#ifdef key_mpp_mpi
[1901]775C
776      INTEGER,DIMENSION(nbdim)   :: lower,upper
777      INTEGER,DIMENSION(nbdim)   :: ltab,utab
778      INTEGER,DIMENSION(nbdim)   :: lb,ub
779      INTEGER,DIMENSION(nbdim,2) :: iminmaxg
780      INTEGER                    :: code
781C
782#endif           
783C     
784C
785      indtab(1:nbdim,2,1) = pttab_child(1:nbdim) + nbtab_child(1:nbdim)
786     &          + deb
787      indtab(1:nbdim,2,2) = indtab(1:nbdim,2,1) + ( fin - deb )
788       
789      indtab(1:nbdim,1,1) = pttab_child(1:nbdim) - fin
790      indtab(1:nbdim,1,2) = pttab_child(1:nbdim) - deb
791                 
792      WHERE (posvartab_child(1:nbdim) == 2)
793        indtab(1:nbdim,1,1) = indtab(1:nbdim,1,1) - 1
794        indtab(1:nbdim,1,2) = indtab(1:nbdim,1,2) - 1
795      END WHERE
796     
797
[2019]798#if !defined key_mpp_mpi
[1901]799      Call Agrif_nbdim_Get_bound_dimension(child%var,lubglob(:,1),
800     &              lubglob(:,2),nbdim)
801C
802#else
803C     
804      Call Agrif_nbdim_Get_bound_dimension(child%var,lb,ub,nbdim)
805     
806      DO i = 1,nbdim
807C
808        Call Agrif_Invloc(lb(i),Agrif_Procrank,i,iminmaxg(i,1))
809        Call Agrif_Invloc(ub(i),Agrif_Procrank,i,iminmaxg(i,2))
810C
811      ENDDO
812C
813      iminmaxg(1:nbdim,2) = - iminmaxg(1:nbdim,2)
814     
815      CALL MPI_ALLREDUCE(iminmaxg,lubglob,2*nbdim,MPI_INTEGER,MPI_MIN,
[2731]816     &                     MPI_COMM_AGRIF,code) 
[1901]817     
818      lubglob(1:nbdim,2) = - lubglob(1:nbdim,2)     
819C
820#endif
821C     
822      indtruetab(1:nbdim,1,1) = max(indtab(1:nbdim,1,1),
823     &     lubglob(1:nbdim,1))
824      indtruetab(1:nbdim,1,2) = max(indtab(1:nbdim,1,2),
825     &     lubglob(1:nbdim,1))
826      indtruetab(1:nbdim,2,1) = min(indtab(1:nbdim,2,1),
827     &     lubglob(1:nbdim,2))
828      indtruetab(1:nbdim,2,2) = min(indtab(1:nbdim,2,2),
829     &     lubglob(1:nbdim,2))
830
831                             
832C 
833C
834      do nb = 1,nbdim
835C
836        do ndir = 1,2     
837C
838          if (loctab_child(nb) /= (-ndir) 
839     &        .AND. loctab_child(nb) /= -3) then
840C           
841              do n = 1,2
842C
843                ptres(nb,n,ndir,nb) = indtruetab(nb,ndir,n) 
844C
845              enddo             
846C
847              do i = 1,nbdim
848C     
849                if (i .NE. nb) then     
850C
851                    if (loctab_child(i) == -1 
852     &                                 .OR. loctab_child(i) == -3) then
853C
854                        ptres(i,1,ndir,nb) = pttab_child(i)
855C
856                      else
857C
858                        ptres(i,1,ndir,nb) = indtruetab(i,1,1)
859C
860                    endif
861C
862                    if (loctab_child(i) == -2 
863     &                                 .OR. loctab_child(i) == -3) then
864C
865                        if (posvartab_child(i) == 1) then
866C
867                            ptres(i,2,ndir,nb) = pttab_child(i)
868     &                                                 + nbtab_child(i)
869C
870                          else
871C
872                            ptres(i,2,ndir,nb) = pttab_child(i)
873     &                                             + nbtab_child(i) - 1
874C
875                        endif                             
876C
877                      else
878C
879                        ptres(i,2,ndir,nb) = indtruetab(i,2,2)
880C
881                    endif                       
882C     
883                endif
884C     
885              enddo
886     
887C
[2019]888#if defined key_mpp_mpi
[1901]889                 Call Agrif_nbdim_Get_bound_dimension
890     &                    (child%var,lower,upper,nbdim)
891
892              do i = 1,nbdim
893C
894                 Call GetLocalBoundaries(ptres(i,1,ndir,nb),
895     &                                  ptres(i,2,ndir,nb),i,
896     &                                  lower(i),upper(i),
897     &                                  ltab(i),utab(i))
898                 ptres2(i,1,ndir,nb) = max(ltab(i),lower(i)) 
899                 ptres2(i,2,ndir,nb) = min(utab(i),upper(i))
900                 if ((i == nb) .AND. (ndir == 1)) then
901                     ptres2(i,2,ndir,nb) = max(utab(i),lower(i))
902                 elseif ((i == nb) .AND. (ndir == 2)) then
903                     ptres2(i,1,ndir,nb) = min(ltab(i),upper(i)) 
904                 endif
905C
906              enddo
907#else
908              ptres2(:,:,ndir,nb) = ptres(:,:,ndir,nb)
909#endif     
910           
911        endif
912     
913        enddo
914       enddo
915C
916      if (child % var % interpIndex 
917     &        /= Agrif_Curgrid % parent % ngridstep .OR.
918     &    child%var%Interpolationshouldbemade ) then
919C
920C     Space interpolation 
921C     
922      kindex = 1   
923     
924
925C
926      do nb = 1,nbdim
927C
928        do ndir = 1,2
929C
930          if (loctab_child(nb) /= (-ndir) 
931     &        .AND. loctab_child(nb) /= -3) then 
932C
933              IF (present(procname)) THEN
934              Call Agrif_InterpnD
935     &             (TYPEInterp(nb,:),parent,child,
936     &              ptres(1:nbdim,1,ndir,nb),ptres(1:nbdim,2,ndir,nb),
937     &              pttab_child(1:nbdim),pttab_Parent(1:nbdim),
938     &              s_Child(1:nbdim),s_Parent(1:nbdim),
939     &              ds_Child(1:nbdim),ds_Parent(1:nbdim),
940     &              restore,.FALSE.,nbdim,procname)
941              ELSE
942              Call Agrif_InterpnD             
943     &             (TYPEInterp(nb,:),parent,child,
944     &              ptres(1:nbdim,1,ndir,nb),ptres(1:nbdim,2,ndir,nb),
945     &              pttab_child(1:nbdim),pttab_Parent(1:nbdim),
946     &              s_Child(1:nbdim),s_Parent(1:nbdim),
947     &              ds_Child(1:nbdim),ds_Parent(1:nbdim),
948     &              restore,.FALSE.,nbdim)
949              ENDIF
950     
951              IF (.NOT. child%var%interpolationshouldbemade) THEN
952C     
953C             Safeguard of the values of the grid variable (at times n and n+1
954C                on the parent grid)
955C     
956              sizetab(1) = 1
957C
958              do i = 1,nbdim
959C         
960                sizetab(1) = sizetab(1)
961     &              * (ptres2(i,2,ndir,nb)-ptres2(i,1,ndir,nb)+1)
962C     
963              enddo
964
965              Call saveAfterInterp(child,
966     &           ptres2(:,:,ndir,nb),kindex,sizetab(1),nbdim)
967C
968           ENDIF
969C
970          endif
971C
972        enddo       
973C
974      enddo     
975C
976C
977      child % var % interpIndex = Agrif_Curgrid % parent % ngridstep
978C
979C
980      endif 
981     
982              IF (.NOT. child%var%interpolationshouldbemade) THEN
983C
984C
985C     Calculation of the coefficients c1t and c2t for the temporary
986C        interpolation
987      if (pweight) then
988C
989          c1t = weight
990C
991        else
992C
993          c1t = (REAL(AGRIF_Nbstepint()) + 1.) / Agrif_Rhot()
994C
995      endif
996C                                 
997      c2t = 1. - c1t
998C           
999C     Time interpolation
1000C
1001      kindex = 1
1002C
1003      do nb = 1,nbdim
1004C
1005        do ndir = 1,2     
1006C
1007          if (loctab_child(nb) /= (-ndir) 
1008     &        .AND. loctab_child(nb) /= -3) then
1009
1010              Call timeInterpolation
1011     &             (child,ptres2(:,:,ndir,nb),kindex,c1t,c2t,nbdim)
1012          endif
1013C
1014        enddo
1015C     
1016      enddo
1017C     
1018
1019       ENDIF
1020C 
1021      End Subroutine Agrif_Correctnd
1022C
1023C
1024C     **************************************************************************
1025CCC   Subroutine saveAfterInterp
1026C     **************************************************************************
1027C
1028      Subroutine saveAfterInterp(child,bounds,kindex,newsize,nbdim)
1029C
1030CCC   Descritpion:
1031CCC   Subroutine used to save the values of the grid variable on the fine grid 
1032CCC   after the space interpolation. 
1033C
1034C     Declarations:
1035C
1036     
1037C
1038C     argument
1039      TYPE (AGRIF_PVariable) :: child   ! The fine grid variable
1040      INTEGER                :: kindex  ! Index indicating where this safeguard
1041                                        ! is done on the fine grid
1042      INTEGER :: nbdim, newsize
1043      INTEGER,DIMENSION(nbdim,2) :: bounds
1044C
1045C     Local scalars
1046      INTEGER                      :: ir,jr,kr,lr,mr,nr
1047C
1048C
1049C     Allocation of the array oldvalues2d
1050
1051C
1052      if (newsize .LE. 0) return
1053C
1054      Call Agrif_Checksize
1055     &     (child,kindex+newsize) 
1056
1057      if (child % var % interpIndex 
1058     &        /= Agrif_Curgrid % parent % ngridstep ) then
1059       child%var%oldvalues2d(1,kindex:kindex+newsize-1)=
1060     &        child%var%oldvalues2d(2,kindex:kindex+newsize-1)
1061      endif
1062
1063      SELECT CASE (nbdim)
1064      CASE (1)
1065
1066!CDIR ALTCODE
1067         do ir=bounds(1,1),bounds(1,2)
1068            child%var%oldvalues2d(2,kindex) =
[2715]1069     &           child%var%parray1(ir)
[1901]1070            kindex = kindex + 1
1071         enddo       
1072C
1073      CASE (2)
1074
1075        do jr=bounds(2,1),bounds(2,2)
1076!CDIR ALTCODE
1077           do ir=bounds(1,1),bounds(1,2)
1078            child%var%oldvalues2d(2,kindex) =
[2715]1079     &           child%var%parray2(ir,jr)
[1901]1080            kindex = kindex + 1
1081           enddo
1082         enddo     
1083C
1084      CASE (3)
1085        do kr=bounds(3,1),bounds(3,2)
1086          do jr=bounds(2,1),bounds(2,2)
1087!CDIR ALTCODE
1088             do ir=bounds(1,1),bounds(1,2)
1089            child%var%oldvalues2d(2,kindex) =
[2715]1090     &           child%var%parray3(ir,jr,kr)
[1901]1091            kindex = kindex + 1
1092             enddo
1093           enddo
1094         enddo       
1095C
1096      CASE (4)
1097          do lr=bounds(4,1),bounds(4,2)
1098           do kr=bounds(3,1),bounds(3,2)
1099             do jr=bounds(2,1),bounds(2,2)
1100!CDIR ALTCODE
1101               do ir=bounds(1,1),bounds(1,2)
1102            child%var%oldvalues2d(2,kindex) =
[2715]1103     &           child%var%parray4(ir,jr,kr,lr)
[1901]1104            kindex = kindex + 1
1105               enddo
1106             enddo
1107           enddo
1108         enddo         
1109C
1110      CASE (5)
1111         do mr=bounds(5,1),bounds(5,2)
1112           do lr=bounds(4,1),bounds(4,2)
1113            do kr=bounds(3,1),bounds(3,2) 
1114              do jr=bounds(2,1),bounds(2,2)
1115!CDIR ALTCODE
1116                 do ir=bounds(1,1),bounds(1,2)
1117            child%var%oldvalues2d(2,kindex) =
[2715]1118     &           child%var%parray5(ir,jr,kr,lr,mr)
[1901]1119            kindex = kindex + 1
1120                 enddo
1121               enddo
1122             enddo
1123           enddo
1124         enddo   
1125C
1126      CASE (6)
1127        do nr=bounds(6,1),bounds(6,2)
1128           do mr=bounds(5,1),bounds(5,2)
1129             do lr=bounds(4,1),bounds(4,2)
1130               do kr=bounds(3,1),bounds(3,2)
1131                do jr=bounds(2,1),bounds(2,2)
1132!CDIR ALTCODE
1133                   do ir=bounds(1,1),bounds(1,2)
1134            child%var%oldvalues2d(2,kindex) =
[2715]1135     &           child%var%parray6(ir,jr,kr,lr,mr,nr)
[1901]1136            kindex = kindex + 1
1137                   enddo
1138                 enddo
1139               enddo
1140             enddo
1141           enddo
1142         enddo       
1143      END SELECT
1144C
1145C                                                 
1146      End subroutine saveAfterInterp
1147C
1148C
1149C
1150C     **************************************************************************
1151CCC   Subroutine timeInterpolation
1152C     **************************************************************************
1153C
1154      Subroutine timeInterpolation(child,bounds,kindex,c1t,c2t,nbdim) 
1155C
1156CCC   Descritpion:
1157CCC   Subroutine for a linear time interpolation on the child grid. 
1158C
1159C     Declarations:
1160C
1161     
1162C
1163C     argument
1164      TYPE (AGRIF_PVariable) :: child  ! The fine grid variable
1165      INTEGER :: nbdim
1166      INTEGER,DIMENSION(nbdim,2) :: bounds
1167      INTEGER                :: kindex ! Index indicating the values of the fine
1168                                       ! grid got before and after the space
1169                                       ! interpolation and used for the time
1170                                       ! interpolation
1171      REAL                   :: c1t,c2t! coefficients for the time interpolation
1172                                       ! (c2t=1-c1t) 
1173C
1174C     Local aruments     
1175      INTEGER :: i
1176C     Local scalars
1177      INTEGER                      :: ir,jr,kr,lr,mr,nr     
1178C
1179C
1180     
1181      SELECT CASE (nbdim)
1182      CASE (1)
1183
1184!CDIR ALTCODE
1185         do ir=bounds(1,1),bounds(1,2)
[2715]1186                child%var%parray1(ir) =
[1901]1187     &           c2t*child % var % oldvalues2d(1,kindex)   
1188     &         + c1t*child % var % oldvalues2d(2,kindex)   
1189            kindex = kindex + 1
1190         enddo       
1191C
1192      CASE (2)
1193
1194        do jr=bounds(2,1),bounds(2,2)
1195!CDIR ALTCODE
1196           do ir=bounds(1,1),bounds(1,2)
[2715]1197                child%var%parray2(ir,jr) =
[1901]1198     &           c2t*child % var % oldvalues2d(1,kindex)   
1199     &         + c1t*child % var % oldvalues2d(2,kindex) 
1200            kindex = kindex + 1
1201           enddo
1202         enddo     
1203C
1204      CASE (3)
1205        do kr=bounds(3,1),bounds(3,2)
1206          do jr=bounds(2,1),bounds(2,2)
1207!CDIR ALTCODE
1208             do ir=bounds(1,1),bounds(1,2)
[2715]1209                child%var%parray3(ir,jr,kr) =
[1901]1210     &           c2t*child % var % oldvalues2d(1,kindex)   
1211     &         + c1t*child % var % oldvalues2d(2,kindex) 
1212            kindex = kindex + 1
1213             enddo
1214           enddo
1215         enddo       
1216C
1217      CASE (4)
1218          do lr=bounds(4,1),bounds(4,2)
1219           do kr=bounds(3,1),bounds(3,2)
1220             do jr=bounds(2,1),bounds(2,2)
1221!CDIR ALTCODE
1222               do ir=bounds(1,1),bounds(1,2)
[2715]1223                child%var%parray4(ir,jr,kr,lr) =
[1901]1224     &           c2t*child % var % oldvalues2d(1,kindex)   
1225     &         + c1t*child % var % oldvalues2d(2,kindex) 
1226            kindex = kindex + 1
1227               enddo
1228             enddo
1229           enddo
1230         enddo         
1231C
1232      CASE (5)
1233         do mr=bounds(5,1),bounds(5,2)
1234           do lr=bounds(4,1),bounds(4,2)
1235            do kr=bounds(3,1),bounds(3,2) 
1236              do jr=bounds(2,1),bounds(2,2)
1237!CDIR ALTCODE
1238                 do ir=bounds(1,1),bounds(1,2)
[2715]1239                child%var%parray5(ir,jr,kr,lr,mr) =
[1901]1240     &           c2t*child % var % oldvalues2d(1,kindex)   
1241     &         + c1t*child % var % oldvalues2d(2,kindex) 
1242            kindex = kindex + 1
1243                 enddo
1244               enddo
1245             enddo
1246           enddo
1247         enddo   
1248C
1249      CASE (6)
1250        do nr=bounds(6,1),bounds(6,2)
1251           do mr=bounds(5,1),bounds(5,2)
1252             do lr=bounds(4,1),bounds(4,2)
1253               do kr=bounds(3,1),bounds(3,2)
1254                do jr=bounds(2,1),bounds(2,2)
1255!CDIR ALTCODE
1256                   do ir=bounds(1,1),bounds(1,2)
[2715]1257                child%var%parray6(ir,jr,kr,lr,mr,nr) =
[1901]1258     &           c2t*child % var % oldvalues2d(1,kindex)   
1259     &         + c1t*child % var % oldvalues2d(2,kindex) 
1260            kindex = kindex + 1
1261                   enddo
1262                 enddo
1263               enddo
1264             enddo
1265           enddo
1266         enddo       
1267      END SELECT
1268                                                   
1269C
1270C
1271      End subroutine timeInterpolation     
1272C
1273C
1274C
1275C     **************************************************************************
1276CCC   Subroutine Agrif_Checksize
1277C     **************************************************************************
1278C
1279      Subroutine Agrif_Checksize(child,newsize)
1280C
1281CCC   Descritpion:
1282CCC   Subroutine used in the saveAfterInterp procedure to allocate the 
1283CCC   oldvalues2d array. 
1284C
1285C     Declarations:
1286C
1287     
1288C     
1289C     TYPE argument
1290      TYPE (AGRIF_PVariable) :: child  ! The fine grid variable
1291C 
1292C     Scalar arguments
1293      INTEGER :: newsize               ! Size of the domains where the boundary
1294                                       ! conditions are calculated 
1295C
1296C     Local arrays
1297      REAL, DIMENSION(:,:), Allocatable :: tempoldvalues ! Temporary array
1298C               
1299C
1300      if (.NOT. associated(child % var % oldvalues2d)) then
1301C
1302          allocate(child % var % oldvalues2d(2,newsize))
1303C 
1304          child % var % oldvalues2d=0.
1305C
1306        else
1307C
1308          if (SIZE(child % var % oldvalues2d,2) < newsize) then   
1309C
1310              allocate(tempoldvalues(2,SIZE(child % var %
1311     &                                    oldvalues2d,2)))
1312C
1313              tempoldvalues = child % var % oldvalues2d
1314C
1315              deallocate(child % var % oldvalues2d)
1316C
1317              allocate(child % var % oldvalues2d(2,newsize))
1318C           
1319              child%var%oldvalues2d=0.
1320C
1321              child % var % oldvalues2d(:,1:SIZE(tempoldvalues,2)) = 
1322     &        tempoldvalues(:,:)
1323C
1324              deallocate(tempoldvalues)
1325C
1326          endif
1327C
1328      endif                                                   
1329C
1330C
1331      End  Subroutine Agrif_Checksize
1332C
1333C
1334C
1335C     
1336      End Module AGRIF_boundary
1337
Note: See TracBrowser for help on using the repository browser.