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.
Changeset 10087 for vendors/AGRIF/CMEMS_2020/AGRIF_FILES/modinterp.F90 – NEMO

Ignore:
Timestamp:
2018-09-05T15:33:44+02:00 (6 years ago)
Author:
rblod
Message:

update AGRIF library

File:
1 edited

Legend:

Unmodified
Added
Removed

  • vendors/AGRIF/CMEMS_2020/AGRIF_FILES/modinterp.F90

    r7752 r10087  
    2626module Agrif_Interpolation 
    2727! 
    28     use Agrif_InterpBasic 
    29     use Agrif_Arrays 
    30     use Agrif_Mask 
    31     use Agrif_CurgridFunctions 
     28  use Agrif_Init 
     29  use Agrif_Arrays 
     30  use Agrif_InterpBasic 
     31  use Agrif_User_Functions 
     32 
    3233#if defined AGRIF_MPI 
    3334    use Agrif_Mpp 
    3435#endif 
     36   
     37    use Agrif_Mask 
    3538! 
    3639    implicit none 
    … …  
    6669    integer, dimension(6) :: ub_child 
    6770    integer, dimension(6) :: lb_parent 
    68     real   , dimension(6) :: s_child,   s_parent 
    69     real   , dimension(6) :: ds_child, ds_parent 
     71    real(kind=8)   , dimension(6) :: s_child,   s_parent 
     72    real(kind=8)   , dimension(6) :: ds_child, ds_parent 
    7073    integer, dimension(child % root_var % nbdim,2,2) :: childarray 
    7174! 
    … …  
    115118    INTEGER, DIMENSION(nbdim), INTENT(in)   :: pttab_Parent !< Index of the first point inside the domain 
    116119                                                            !<    for the parent grid variable 
    117     REAL,    DIMENSION(nbdim), INTENT(in)   :: s_Child,s_Parent   !< Positions of the parent and child grids 
    118     REAL,    DIMENSION(nbdim), INTENT(in)   :: ds_Child,ds_Parent !< Space steps of the parent and child grids 
     120    REAL(kind=8),    DIMENSION(nbdim), INTENT(in)   :: s_Child,s_Parent   !< Positions of the parent and child grids 
     121    REAL(kind=8),    DIMENSION(nbdim), INTENT(in)   :: ds_Child,ds_Parent !< Space steps of the parent and child grids 
    119122    TYPE(Agrif_Variable),      pointer      :: restore            !< Indicates points where interpolation 
    120123    LOGICAL,                   INTENT(in)   :: torestore          !< Indicates if the array restore is used 
    … …  
    128131    INTEGER                       :: i,j,k,l,m,n 
    129132    INTEGER, DIMENSION(nbdim)     :: pttruetab,cetruetab 
    130     INTEGER, DIMENSION(nbdim)     :: indmin,     indmax 
     133    INTEGER, DIMENSION(nbdim)     :: indmin,     indmax, indmin_required_p, indmax_required_p 
    131134    INTEGER, DIMENSION(nbdim)     :: indminglob, indmaxglob 
    132135#if defined AGRIF_MPI 
    133136    INTEGER, DIMENSION(nbdim)     :: indminglob2,indmaxglob2 
     137    INTEGER, DIMENSION(nbdim)     :: indminglob3,indmaxglob3 
    134138#endif 
    135139    LOGICAL, DIMENSION(nbdim)     :: noraftab 
    136     REAL   , DIMENSION(nbdim)     :: s_Child_temp,s_Parent_temp 
     140    REAL(kind=8)   , DIMENSION(nbdim)     :: s_Child_temp,s_Parent_temp,s_Parent_temp_p 
    137141    INTEGER, DIMENSION(nbdim)     :: lowerbound, upperbound, coords 
    138142    INTEGER, DIMENSION(nbdim,2,2), INTENT(OUT) :: childarray 
    … …  
    148152    INTEGER, DIMENSION(nbdim,4,0:Agrif_Nbprocs-1)   :: tab4 
    149153    INTEGER, DIMENSION(nbdim,0:Agrif_Nbprocs-1,8)   :: tab4t 
     154    INTEGER,DIMENSION(nbdim,2) :: tab5 
     155    INTEGER,DIMENSION(nbdim,2,0:Agrif_Nbprocs-1) :: tab6 
     156    INTEGER,DIMENSION(nbdim,0:Agrif_Nbprocs-1,2) :: tab5t 
    150157    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1)           :: memberinall 
    151158    LOGICAL, DIMENSION(0:Agrif_Nbprocs-1)           :: sendtoproc1, recvfromproc1 
    … …  
    167174            child % list_interp,                                    & 
    168175            pttab, petab, pttab_Child, pttab_Parent, nbdim,         & 
    169             indmin, indmax, indminglob, indmaxglob,                 & 
     176            indmin, indmax, indmin_required_p, indmax_required_p,   & 
     177            indminglob, indmaxglob,                                 & 
    170178            pttruetab, cetruetab, memberin                          & 
    171179#if defined AGRIF_MPI 
    … …  
    174182#endif 
    175183        ) 
     184 
    176185! 
    177186    if (.not.find_list_interp) then 
    178187! 
     188! output : lowerbound and upperbound are the (local) lower and upper bounds of the child arrays 
     189 
    179190        call Agrif_get_var_bounds_array(child, lowerbound, upperbound, nbdim) 
     191 
     192! input : pttab, petab : global indexes where the interpolation is required 
     193! output : pttruetab, cetruetab : global indexes restricted to the bounds of the current processor  
     194! output : memberin is false if the current processor is not involved in the interpolation 
     195 
    180196        call Agrif_Childbounds(nbdim, lowerbound, upperbound,               & 
    181197                               pttab, petab, Agrif_Procrank, coords,        & 
    182198                               pttruetab, cetruetab, memberin) 
     199          
     200 
     201 
     202! output : indminglob, indmaxglob : global indexes required on the parent grid for the interpolation 
     203! output : s_Parent_temp, s_Child_temp : local s_Parent, s_Child relatively to indmin ou pttab 
    183204        call Agrif_Parentbounds(type_interp,nbdim,indminglob,indmaxglob,    & 
     205                                indmin_required_p, indmax_required_p,           & 
    184206                                s_Parent_temp,s_Child_temp,                 & 
    185207                                s_Child,ds_Child,                           & 
    … …  
    190212#if defined AGRIF_MPI 
    191213        if (memberin) then 
     214 
     215! output : indmin, indmax : global indexes required on the parent grid for the interpolation on the current processor (i.e. on pttruetab, cetruetab) 
     216! output : s_Parent_temp, s_Child_temp : local s_Parent, s_Child relatively to indmin ou pttruetab 
    192217            call Agrif_Parentbounds(type_interp,nbdim,indmin,indmax,        & 
     218                                    indmin_required_p, indmax_required_p,       & 
    193219                                    s_Parent_temp,s_Child_temp,             & 
    194220                                    s_Child,ds_Child,                       & 
    … …  
    200226 
    201227        local_proc = Agrif_Procrank 
     228 
     229! output : lowerbound and upperbound are the (local) lower and upper bounds of the parent arrays 
    202230        call Agrif_get_var_bounds_array(parent,lowerbound,upperbound,nbdim) 
    203231        call Agrif_ChildGrid_to_ParentGrid() 
    204 ! 
     232 
     233! input : indminglob,indmaxglob : global indexes where data are required for the interpolation 
     234! output : indminglob2,indmaxglob2 : global indexes restricted to the bounds of the current processor  
     235! output : member is false if the current processor does not need to send data 
     236! output : indminglob3,indmaxglob3 : global bounds on the current processor for the parent array 
     237 
    205238        call Agrif_Childbounds(nbdim,lowerbound,upperbound,                 & 
    206239                               indminglob,indmaxglob, local_proc, coords,   & 
    207                                indminglob2,indmaxglob2,member) 
     240                               indminglob2,indmaxglob2,member,              & 
     241                               indminglob3,indmaxglob3) 
    208242! 
    209243        if (member) then 
     244 
     245! output : parentarray 
     246! output : parentarray (:,:,2) : indminglob2, indmaxglob2 in term of local indexes on current processor 
     247! output : parentarray (:,:,1) : indminglob2, indmaxglob2 restricted to the current processor (different from indminglob2 ???) 
     248! output : member is .false. is the current processor has not data to send 
     249 
    210250            call Agrif_GlobalToLocalBounds(parentarray,                     & 
    211251                                           lowerbound,  upperbound,         & 
    … …  
    216256        call Agrif_ParentGrid_to_ChildGrid() 
    217257#else 
     258 
     259! In the sequentiel case, the following lines ensure that the bounds needed on the parent grid in the interpolation  
     260! do not exceed lower and upper bounds of the parent array 
     261 
     262! output : lowerbound and upperbound are the (local) lower and upper bounds of the parent arrays 
     263        call Agrif_get_var_bounds_array(parent,lowerbound,upperbound,nbdim) 
     264        call Agrif_ChildGrid_to_ParentGrid() 
     265 
     266! input : indminglob,indmaxglob : global indexes where data are required for the interpolation 
     267! output : indminglob2,indmaxglob2 : global indexes restricted to the bounds of the current processor  
     268! output : member is false if the current processor does not need to send data 
     269 
     270        call Agrif_Childbounds(nbdim,lowerbound,upperbound,                 & 
     271                               indminglob,indmaxglob, Agrif_Procrank, coords,   & 
     272                               indmin,indmax,member) 
     273 
     274        call Agrif_ParentGrid_to_ChildGrid() 
     275 
     276        indminglob = indmin 
     277        indmaxglob = indmax 
     278 
    218279        parentarray(:,1,1) = indminglob 
    219280        parentarray(:,2,1) = indmaxglob 
    220281        parentarray(:,1,2) = indminglob 
    221282        parentarray(:,2,2) = indmaxglob 
    222         indmin = indminglob 
    223         indmax = indmaxglob 
     283  
     284!       indmin = indminglob 
     285!        indmax = indmaxglob 
     286 
    224287        member = .TRUE. 
     288        s_Parent_temp = s_Parent + (indminglob - pttab_Parent) * ds_Parent 
     289 
    225290#endif 
     291!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
     292! Correct for non refined directions 
     293!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
     294        do i=1,nbdim 
     295          if (coords(i) == 0) then 
     296             indmin(i) = indminglob(i) 
     297             indmax(i) = indmaxglob(i) 
     298             pttruetab(i) = indminglob(i) 
     299             cetruetab(i) = indmaxglob(i) 
     300          endif 
     301        enddo 
    226302 
    227303    else 
    … …  
    231307        s_Child_temp  = s_Child + (pttruetab - pttab_Child) * ds_Child 
    232308#else 
     309 
     310! In the sequentiel case, the following lines ensure that the bounds needed on the parent grid in the interpolation  
     311! do not exceed lower and upper bounds of the parent array 
     312 
     313! output : lowerbound and upperbound are the (local) lower and upper bounds of the parent arrays 
     314        call Agrif_get_var_bounds_array(parent,lowerbound,upperbound,nbdim) 
     315        call Agrif_ChildGrid_to_ParentGrid() 
     316 
     317! input : indminglob,indmaxglob : global indexes where data are required for the interpolation 
     318! output : indminglob2,indmaxglob2 : global indexes restricted to the bounds of the current processor  
     319! output : member is false if the current processor does not need to send data 
     320 
     321        call Agrif_Childbounds(nbdim,lowerbound,upperbound,                 & 
     322                               indminglob,indmaxglob, Agrif_Procrank, coords,   & 
     323                               indmin,indmax,member) 
     324 
     325        call Agrif_ParentGrid_to_ChildGrid() 
     326 
     327        indminglob = indmin 
     328        indmaxglob = indmax 
     329 
    233330        parentarray(:,1,1) = indminglob 
    234331        parentarray(:,2,1) = indmaxglob 
    235332        parentarray(:,1,2) = indminglob 
    236333        parentarray(:,2,2) = indmaxglob 
    237         indmin = indminglob 
    238         indmax = indmaxglob 
     334 !       indmin = indminglob 
     335 !       indmax = indmaxglob 
    239336        member = .TRUE. 
    240337        s_Parent_temp = s_Parent + (indminglob - pttab_Parent) * ds_Parent 
    … …  
    246343        if (.not.associated(tempP)) allocate(tempP) 
    247344! 
     345 
    248346        call Agrif_array_allocate(tempP,parentarray(:,1,1),parentarray(:,2,1),nbdim) 
    249347        call Agrif_var_set_array_tozero(tempP,nbdim) 
    … …  
    286384                      parentarray(6,1,2),parentarray(6,2,2),.TRUE.,nb,ndir) 
    287385        end select 
     386 
    288387! 
    289388        call Agrif_ParentGrid_to_ChildGrid() 
    … …  
    298397        tab3(:,3) = indmin(:) 
    299398        tab3(:,4) = indmax(:) 
     399        tab5(:,1) = indminglob3(:) 
     400        tab5(:,2) = indmaxglob3(:) 
    300401! 
    301402        call MPI_ALLGATHER(tab3,4*nbdim,MPI_INTEGER,tab4,4*nbdim,MPI_INTEGER,Agrif_mpi_comm,code) 
    302  
     403        call MPI_ALLGATHER(tab5,2*nbdim,MPI_INTEGER,tab6,2*nbdim,MPI_INTEGER,Agrif_mpi_comm,code) 
    303404        if (.not.associated(tempPextend))   allocate(tempPextend) 
    304405 
    … …  
    311412        enddo 
    312413 
     414        do k=0,Agrif_Nbprocs-1 
     415          do j=1,2 
     416            do i=1,nbdim 
     417               tab5t(i,k,j) = tab6(i,j,k) 
     418            enddo 
     419          enddo 
     420        enddo 
     421       
    313422        memberin1(1) = memberin 
    314423        call MPI_ALLGATHER(memberin1,1,MPI_LOGICAL,memberinall,1,MPI_LOGICAL,Agrif_mpi_comm,code) 
    … …  
    319428                                     sendtoproc1,recvfromproc1,         & 
    320429                                     tab4t(:,:,5),tab4t(:,:,6),         & 
    321                                      tab4t(:,:,7),tab4t(:,:,8) ) 
     430                                     tab4t(:,:,7),tab4t(:,:,8),         & 
     431                                     tab5t(:,:,1),tab5t(:,:,2)) 
    322432    endif 
    323433 
    … …  
    333443                child%list_interp,pttab,petab,                  & 
    334444                pttab_Child,pttab_Parent,indmin,indmax,         & 
     445                indmin_required_p, indmax_required_p,           & 
    335446                indminglob,indmaxglob,                          & 
    336447                pttruetab,cetruetab,                            & 
    … …  
    345456    endif 
    346457! 
     458 
    347459    if (memberin) then 
    348460! 
    349461        if (.not.associated(tempC)) allocate(tempC) 
    350462! 
     463 
    351464        call Agrif_array_allocate(tempC,pttruetab,cetruetab,nbdim) 
     465 
    352466! 
    353467!       Special values on the parent grid 
    … …  
    357471! 
    358472            if (.not.associated(parentvalues))  allocate(parentvalues) 
    359 ! 
     473!t 
     474 
    360475            call Agrif_array_allocate(parentvalues,indmin,indmax,nbdim) 
    361476            call Agrif_var_full_copy_array(parentvalues,tempPextend,nbdim) 
    362477! 
    363             call Agrif_CheckMasknD(tempPextend,parentvalues,    & 
    364                     indmin(1:nbdim),indmax(1:nbdim),            & 
    365                     indmin(1:nbdim),indmax(1:nbdim),            & 
     478            call Agrif_CheckMasknD(tempPextend,parentvalues,                & 
     479                    indmin(1:nbdim),indmax(1:nbdim),                        & 
     480                    indmin(1:nbdim),indmax(1:nbdim),                        & 
     481                    indmin_required_p(1:nbdim),indmax_required_p(1:nbdim),  & 
    366482                    noraftab(1:nbdim),nbdim) 
    367483! 
    … …  
    391507                                                ds_Child(1:2),    ds_Parent(1:2) ) 
    392508            case(3) 
    393                 call Agrif_Interp_3D_recursive( type_interp(1:3),                       & 
    394                                                 tempPextend % array3,                   & 
    395                                                 tempC       % array3,                   & 
    396                                                 indmin(1:3), indmax(1:3),               & 
    397                                                 pttruetab(1:3),    cetruetab(1:3),      & 
    398                                                 s_Child_temp(1:3), s_Parent_temp(1:3),  & 
     509                s_Parent_temp_p = s_Parent + (indmin_required_p - pttab_Parent) * ds_Parent 
     510                call Agrif_Interp_3D_recursive( type_interp(1:3),                                 & 
     511                                                tempPextend % array3(                             & 
     512                                                indmin_required_p(1):indmax_required_p(1),        & 
     513                                                indmin_required_p(2):indmax_required_p(2),        & 
     514                                                indmin_required_p(3):indmax_required_p(3)),       & 
     515                                                tempC       % array3,                             & 
     516                                                indmin_required_p(1:3), indmax_required_p(1:3),   & 
     517                                                pttruetab(1:3),    cetruetab(1:3),                & 
     518                                                s_Child_temp(1:3), s_Parent_temp_p(1:3),          & 
    399519                                                ds_Child(1:3),    ds_Parent(1:3) ) 
     520 
    400521            case(4) 
    401                 call Agrif_Interp_4D_recursive( type_interp(1:4),                       & 
    402                                                 tempPextend % array4,                   & 
    403                                                 tempC       % array4,                   & 
    404                                                 indmin(1:4), indmax(1:4),               & 
    405                                                 pttruetab(1:4),    cetruetab(1:4),      & 
    406                                                 s_Child_temp(1:4), s_Parent_temp(1:4),  & 
     522                s_Parent_temp_p = s_Parent + (indmin_required_p - pttab_Parent) * ds_Parent 
     523                call Agrif_Interp_4D_recursive( type_interp(1:4),                                 & 
     524                                                tempPextend % array4(                             & 
     525                                                indmin_required_p(1):indmax_required_p(1),        & 
     526                                                indmin_required_p(2):indmax_required_p(2),        & 
     527                                                indmin_required_p(3):indmax_required_p(3),        & 
     528                                                indmin_required_p(4):indmax_required_p(4)),       & 
     529                                                tempC       % array4,                             & 
     530                                                indmin_required_p(1:4), indmax_required_p(1:4),   & 
     531                                                pttruetab(1:4),    cetruetab(1:4),                & 
     532                                                s_Child_temp(1:4), s_Parent_temp_p(1:4),          & 
    407533                                                ds_Child(1:4),    ds_Parent(1:4) ) 
    408534            case(5) 
    … …  
    595721        else    ! .not.to_restore 
    596722! 
     723 
    597724            if (memberin) then 
    598725    ! 
    … …  
    715842        endif 
    716843 
     844 
    717845        call Agrif_array_deallocate(tempPextend,nbdim) 
    718846        call Agrif_array_deallocate(tempC,nbdim) 
    … …  
    736864!--------------------------------------------------------------------------------------------------- 
    737865subroutine Agrif_Parentbounds ( type_interp, nbdim, indmin, indmax, & 
     866                                indmin_required,indmax_required,    & 
    738867                                s_Parent_temp, s_Child_temp,        & 
    739868                                s_Child, ds_Child,                  & 
    … …  
    745874    INTEGER,                   intent(in)  :: nbdim 
    746875    INTEGER, DIMENSION(nbdim), intent(out) :: indmin, indmax 
    747     REAL,    DIMENSION(nbdim), intent(out) :: s_Parent_temp, s_child_temp 
    748     REAL,    DIMENSION(nbdim), intent(in)  :: s_Child, ds_child 
    749     REAL,    DIMENSION(nbdim), intent(in)  :: s_Parent,ds_Parent 
     876    INTEGER, DIMENSION(nbdim), intent(out) :: indmin_required, indmax_required 
     877    REAL(kind=8),    DIMENSION(nbdim), intent(out) :: s_Parent_temp, s_child_temp 
     878    REAL(kind=8),    DIMENSION(nbdim), intent(in)  :: s_Child, ds_child 
     879    REAL(kind=8),    DIMENSION(nbdim), intent(in)  :: s_Parent,ds_Parent 
    750880    INTEGER, DIMENSION(nbdim), intent(in)  :: pttruetab, cetruetab 
    751881    INTEGER, DIMENSION(nbdim), intent(in)  :: pttab_Child, pttab_Parent 
    … …  
    753883    INTEGER, DIMENSION(nbdim), intent(in)  :: coords 
    754884! 
     885    REAL(kind=8) :: xpmin, xpmax 
     886    INTEGER :: coeffraf 
    755887    INTEGER :: i 
    756     REAL,DIMENSION(nbdim) :: dim_newmin, dim_newmax 
     888    REAL(kind=8),DIMENSION(nbdim) :: dim_newmin, dim_newmax 
    757889! 
    758890    dim_newmin = s_Child + (pttruetab - pttab_Child) * ds_Child 
    … …  
    763895        indmin(i) = pttab_Parent(i) + agrif_int((dim_newmin(i)-s_Parent(i))/ds_Parent(i)) 
    764896        indmax(i) = pttab_Parent(i) + agrif_ceiling((dim_newmax(i)-s_Parent(i))/ds_Parent(i)) 
     897 
     898        coeffraf = nint(ds_Parent(i)/ds_Child(i)) 
     899         
     900        indmin_required(i) = indmin(i) 
     901        indmax_required(i) = indmax(i) 
    765902! 
    766903!       Necessary for the Quadratic interpolation 
    767904! 
     905 
    768906        if ( (pttruetab(i) == cetruetab(i)) .and. (posvar(i) == 1) ) then 
     907            if (Agrif_UseSpecialValue) then 
     908               indmin(i) = indmin(i)-MaxSearch 
     909               indmax(i) = indmax(i)+MaxSearch 
     910            endif 
    769911        elseif ( coords(i) == 0 ) then  ! (interptab == 'N') 
    770912        elseif ( (type_interp(i) == Agrif_ppm)     .or.     & 
    … …  
    772914                 (type_interp(i) == Agrif_ppm_lim) .or.     & 
    773915                 (type_interp(i) == Agrif_weno) ) then 
    774             indmin(i) = indmin(i) - 2 
    775             indmax(i) = indmax(i) + 2 
    776  
     916                  
     917            if ((mod(coeffraf,2) == 0).AND.(posvar(i)==2)) then 
     918             
     919              xpmax = s_Parent(i)+(indmax(i)-pttab_Parent(i))*ds_Parent(i) 
     920              if (xpmax > dim_newmax(i)+ds_Child(i)) then 
     921                indmax(i) = indmax(i) + 1 
     922              else 
     923                indmax(i) = indmax(i) + 2 
     924              endif 
     925               
     926              xpmin = s_Parent(i)+(indmin(i)-pttab_Parent(i))*ds_Parent(i) 
     927              if (xpmin < dim_newmin(i)-ds_Child(i)) then 
     928                indmin(i) = indmin(i) - 1 
     929              else 
     930                indmin(i) = indmin(i) - 2 
     931              endif 
     932               
     933            else 
     934              indmin(i) = indmin(i) - 2 
     935              indmax(i) = indmax(i) + 2 
     936            endif 
     937 
     938            indmin_required(i) = indmin(i) 
     939            indmax_required(i) = indmax(i) 
     940         
    777941            if (Agrif_UseSpecialValue) then 
    778942               indmin(i) = indmin(i)-MaxSearch 
    779943               indmax(i) = indmax(i)+MaxSearch 
    780944            endif 
    781  
     945        elseif (type_interp(i) == Agrif_linearconservlim) then 
     946         
     947            if ((mod(coeffraf,2) == 0).AND.(posvar(i)==2)) then 
     948             
     949              xpmax = s_Parent(i)+(indmax(i)-pttab_Parent(i))*ds_Parent(i) 
     950              if (xpmax > dim_newmax(i)+ds_Child(i)) then 
     951                indmax(i) = indmax(i) 
     952              else 
     953                indmax(i) = indmax(i) + 1 
     954              endif 
     955               
     956              xpmin = s_Parent(i)+(indmin(i)-pttab_Parent(i))*ds_Parent(i) 
     957              if (xpmin < dim_newmin(i)-ds_Child(i)) then 
     958                indmin(i) = indmin(i) 
     959              else 
     960                indmin(i) = indmin(i) - 1 
     961              endif 
     962               
     963            else 
     964              indmin(i) = indmin(i) - 1 
     965              indmax(i) = indmax(i) + 1 
     966            endif 
     967 
     968            indmin_required(i) = indmin(i) 
     969            indmax_required(i) = indmax(i) 
     970         
     971            if (Agrif_UseSpecialValue) then 
     972               indmin(i) = indmin(i)-MaxSearch 
     973               indmax(i) = indmax(i)+MaxSearch 
     974            endif 
     975             
    782976        elseif ( (type_interp(i) /= Agrif_constant) .and.   & 
    783977                 (type_interp(i) /= Agrif_linear) ) then 
    784978            indmin(i) = indmin(i) - 1 
    785979            indmax(i) = indmax(i) + 1 
     980             
     981            indmin_required(i) = indmin(i) 
     982            indmax_required(i) = indmax(i) 
    786983 
    787984            if (Agrif_UseSpecialValue) then 
    … …  
    789986               indmax(i) = indmax(i)+MaxSearch 
    790987            endif 
    791  
    792988        elseif ( (type_interp(i) == Agrif_constant) .or.   & 
    793989                 (type_interp(i) == Agrif_linear) ) then 
     990            indmin_required(i) = indmin(i) 
     991            indmax_required(i) = indmax(i) 
    794992            if (Agrif_UseSpecialValue) then 
    795993               indmin(i) = indmin(i)-MaxSearch 
    796994               indmax(i) = indmax(i)+MaxSearch 
    797995            endif 
    798  
    799996        endif 
     997 
    800998! 
    801999    enddo 
    … …  
    8221020    integer,            intent(in)  :: indmin, indmax 
    8231021    integer,            intent(in)  :: pttab_child, petab_child 
    824     real,               intent(in)  :: s_child, s_parent 
    825     real,               intent(in)  :: ds_child, ds_parent 
     1022    real(kind=8),               intent(in)  :: s_child, s_parent 
     1023    real(kind=8),               intent(in)  :: ds_child, ds_parent 
    8261024    real, dimension(            & 
    8271025        indmin:indmax           & 
    … …  
    8571055    integer, dimension(2),              intent(in)  :: indmin, indmax 
    8581056    integer, dimension(2),              intent(in)  :: pttab_child, petab_child 
    859     real,    dimension(2),              intent(in)  :: s_child, s_parent 
    860     real,    dimension(2),              intent(in)  :: ds_child, ds_parent 
     1057    real(kind=8),    dimension(2),              intent(in)  :: s_child, s_parent 
     1058    real(kind=8),    dimension(2),              intent(in)  :: ds_child, ds_parent 
    8611059    real,    dimension(                 & 
    8621060        indmin(1):indmax(1),            & 
    … …  
    8751073        indmin(2):indmax(2),            & 
    8761074        pttab_child(1):petab_child(1))  :: tabtemp_trsp 
    877     integer                             :: i, j, coeffraf 
     1075    integer                             :: i, j, coeffraf, locind_child_left, ideb 
    8781076!--------------------------------------------------------------------------------------------------- 
    8791077! 
    … …  
    9001098                    s_parent(1),s_child(1),ds_parent(1),ds_child(1),1) 
    9011099!---CDIR NEXPAND 
    902         call PPM1dAfterCompute(tabin,tabtemp,size(tabin),size(tabtemp),1) 
     1100        call PPM1dAfterCompute(tabin,tabtemp,size(tabin),size(tabtemp),1,indchildppm(:,1),tabppm(:,:,1)) 
     1101    else if (coeffraf == 1) then 
     1102        locind_child_left = indmin(1) + agrif_int((s_child(1)-s_parent(1))/ds_parent(1)) 
     1103         
     1104            do j = indmin(2), indmax(2) 
     1105            ideb = locind_child_left 
     1106            do i = pttab_child(1), petab_child(1) 
     1107                tabtemp(i,j) = tabin(ideb,j) 
     1108                ideb = ideb + 1 
     1109            enddo 
     1110            enddo 
     1111 
    9031112    else 
    9041113        do j = indmin(2),indmax(2) 
    … …  
    9411150!---CDIR NEXPAND 
    9421151        call PPM1dAfterCompute(tabtemp_trsp, tabout_trsp,    & 
    943                                size(tabtemp_trsp), size(tabout_trsp), 2) 
     1152                               size(tabtemp_trsp), size(tabout_trsp), 2, & 
     1153                              indchildppm(:,2),tabppm(:,:,2)) 
    9441154    else 
    9451155        do i = pttab_child(1), petab_child(1) 
    … …  
    9761186    integer, dimension(3),              intent(in)  :: indmin, indmax 
    9771187    integer, dimension(3),              intent(in)  :: pttab_child, petab_child 
    978     real,    dimension(3),              intent(in)  :: s_child, s_parent 
    979     real,    dimension(3),              intent(in)  :: ds_child, ds_parent 
     1188    real(kind=8),    dimension(3),              intent(in)  :: s_child, s_parent 
     1189    real(kind=8),    dimension(3),              intent(in)  :: ds_child, ds_parent 
    9801190    real,    dimension(                 & 
    9811191        indmin(1):indmax(1),            & 
    … …  
    9911201        pttab_child(2):petab_child(2),  & 
    9921202        indmin(3):indmax(3))            :: tabtemp 
    993     integer                             :: i, j, k, coeffraf 
     1203    integer                             :: i, j, k, coeffraf,kp,kp1,kp2,kp3,kp4,kref 
    9941204    integer                             :: locind_child_left, kdeb 
     1205    real(kind=8)    :: ypos,globind_parent_left 
     1206    real(kind=8)    :: deltax, invdsparent 
     1207    real    :: t2,t3,t4,t5,t6,t7,t8 
     1208    integer :: locind_parent_left 
     1209 
    9951210! 
    9961211    coeffraf = nint ( ds_parent(1) / ds_child(1) ) 
    … …  
    10511266            enddo 
    10521267        enddo 
     1268    else if (type_interp(3) == Agrif_LAGRANGE) then 
     1269      invdsparent = 1./ds_parent(3) 
     1270      ypos = s_child(3) 
     1271      do k=pttab_child(3), petab_child(3) 
     1272        locind_parent_left = indmin(3)+agrif_int((ypos - s_parent(3))/ds_parent(3)) 
     1273        globind_parent_left = s_parent(3) + (locind_parent_left - indmin(3))*ds_parent(3) 
     1274        deltax = invdsparent*(ypos-globind_parent_left) 
     1275        deltax = nint(coeffraf*deltax)/real(coeffraf) 
     1276        ypos = ypos + ds_child(3) 
     1277 
     1278        if (abs(deltax) <= 0.0001) then 
     1279          do j = pttab_child(2), petab_child(2) 
     1280          do i = pttab_child(1), petab_child(1) 
     1281             tabout(i,j,k) = tabtemp(i,j,locind_parent_left) 
     1282          enddo 
     1283          enddo 
     1284        else 
     1285         t2 = deltax - 2. 
     1286        t3 = deltax - 1. 
     1287        t4 = deltax + 1. 
     1288 
     1289        t5 = -(1./6.)*deltax*t2*t3 
     1290        t6 = 0.5*t2*t3*t4 
     1291        t7 = -0.5*deltax*t2*t4 
     1292        t8 = (1./6.)*deltax*t3*t4       
     1293          do j = pttab_child(2), petab_child(2) 
     1294          do i = pttab_child(1), petab_child(1) 
     1295             tabout(i,j,k) = t5*tabtemp(i,j,locind_parent_left-1) + t6*tabtemp(i,j,locind_parent_left)    & 
     1296              +t7*tabtemp(i,j,locind_parent_left+1) + t8*tabtemp(i,j,locind_parent_left+2) 
     1297          enddo 
     1298          enddo 
     1299 
     1300        endif 
     1301 
     1302      enddo 
     1303    else if (type_interp(3) == Agrif_PPM) then 
     1304     call PPM1dPrecompute2d(1, & 
     1305                               indmax(3)-indmin(3)+1,           & 
     1306                               petab_child(3)-pttab_child(3)+1, & 
     1307                               s_parent(3),s_child(3),ds_parent(3),ds_child(3),1) 
     1308 
     1309     do k=pttab_child(3),petab_child(3) 
     1310        kref = k-pttab_child(3)+1 
     1311        kp=indmin(3)+indparentppm(kref,1)-1 
     1312        kp1 = kp + 1 
     1313        kp2 = kp1 + 1 
     1314        kp3 = kp2 + 1 
     1315        kp4 = kp3 + 1 
     1316        do j = pttab_child(2), petab_child(2) 
     1317        do i = pttab_child(1), petab_child(1) 
     1318         tabout(i,j,k) = tabppm(1,indchildppm(kref,1),1)*tabtemp(i,j,kp)   + & 
     1319                         tabppm(2,indchildppm(kref,1),1)*tabtemp(i,j,kp1)  + & 
     1320                         tabppm(3,indchildppm(kref,1),1)*tabtemp(i,j,kp2)  + & 
     1321                         tabppm(4,indchildppm(kref,1),1)*tabtemp(i,j,kp3)  + & 
     1322                         tabppm(5,indchildppm(kref,1),1)*tabtemp(i,j,kp4) 
     1323        enddo 
     1324        enddo 
     1325     enddo 
     1326 
    10531327    else 
     1328 
    10541329        do j = pttab_child(2), petab_child(2) 
    10551330        do i = pttab_child(1), petab_child(1) 
    … …  
    10631338        enddo 
    10641339        enddo 
     1340 
    10651341    endif 
    10661342!--------------------------------------------------------------------------------------------------- 
    … …  
    10831359    integer, dimension(4),              intent(in)  :: indmin, indmax 
    10841360    integer, dimension(4),              intent(in)  :: pttab_child, petab_child 
    1085     real,    dimension(4),              intent(in)  :: s_child, s_parent 
    1086     real,    dimension(4),              intent(in)  :: ds_child, ds_parent 
     1361    real(kind=8),    dimension(4),              intent(in)  :: s_child, s_parent 
     1362    real(kind=8),    dimension(4),              intent(in)  :: ds_child, ds_parent 
    10871363    real,    dimension(                 & 
    10881364        indmin(1):indmax(1),            & 
    … …  
    11021378        indmin(4):indmax(4))            :: tabtemp 
    11031379    integer                             :: i, j, k, l 
     1380 
     1381    real(kind=8)    :: ypos,globind_parent_left 
     1382    real(kind=8)    :: deltax, invdsparent 
     1383    real    :: t2,t3,t4,t5,t6,t7,t8 
     1384    integer :: locind_parent_left, coeffraf 
    11041385! 
    11051386    do l = indmin(4), indmax(4) 
    … …  
    11171398    enddo 
    11181399! 
     1400     if (type_interp(4) == Agrif_LAGRANGE) then 
     1401      coeffraf = nint(ds_parent(4)/ds_child(4)) 
     1402      invdsparent = 1./ds_parent(4) 
     1403      ypos = s_child(4) 
     1404      do l=pttab_child(4), petab_child(4) 
     1405        locind_parent_left = indmin(4)+agrif_int((ypos - s_parent(4))/ds_parent(4)) 
     1406        globind_parent_left = s_parent(4) + (locind_parent_left - indmin(4))*ds_parent(4) 
     1407        deltax = invdsparent*(ypos-globind_parent_left) 
     1408        deltax = nint(coeffraf*deltax)/real(coeffraf) 
     1409        ypos = ypos + ds_child(4) 
     1410 
     1411        if (abs(deltax) <= 0.0001) then 
     1412          do k = pttab_child(3), petab_child(3) 
     1413          do j = pttab_child(2), petab_child(2) 
     1414          do i = pttab_child(1), petab_child(1) 
     1415             tabout(i,j,k,l) = tabtemp(i,j,k,locind_parent_left) 
     1416          enddo 
     1417          enddo 
     1418          enddo 
     1419        else 
     1420         t2 = deltax - 2. 
     1421        t3 = deltax - 1. 
     1422        t4 = deltax + 1. 
     1423 
     1424        t5 = -(1./6.)*deltax*t2*t3 
     1425        t6 = 0.5*t2*t3*t4 
     1426        t7 = -0.5*deltax*t2*t4 
     1427        t8 = (1./6.)*deltax*t3*t4       
     1428          do k = pttab_child(3), petab_child(3) 
     1429          do j = pttab_child(2), petab_child(2) 
     1430          do i = pttab_child(1), petab_child(1) 
     1431             tabout(i,j,k,l) = t5*tabtemp(i,j,k,locind_parent_left-1) + t6*tabtemp(i,j,k,locind_parent_left)    & 
     1432              +t7*tabtemp(i,j,k,locind_parent_left+1) + t8*tabtemp(i,j,k,locind_parent_left+2) 
     1433          enddo 
     1434          enddo 
     1435          enddo 
     1436        endif 
     1437 
     1438      enddo 
     1439    else 
    11191440    do k = pttab_child(3), petab_child(3) 
    11201441    do j = pttab_child(2), petab_child(2) 
    … …  
    11301451    enddo 
    11311452    enddo 
     1453    endif 
    11321454!--------------------------------------------------------------------------------------------------- 
    11331455end subroutine Agrif_Interp_4D_recursive 
    … …  
    11491471    integer, dimension(5),              intent(in)  :: indmin, indmax 
    11501472    integer, dimension(5),              intent(in)  :: pttab_child, petab_child 
    1151     real,    dimension(5),              intent(in)  :: s_child, s_parent 
    1152     real,    dimension(5),              intent(in)  :: ds_child, ds_parent 
     1473    real(kind=8),    dimension(5),              intent(in)  :: s_child, s_parent 
     1474    real(kind=8),    dimension(5),              intent(in)  :: ds_child, ds_parent 
    11531475    real,    dimension(                 & 
    11541476        indmin(1):indmax(1),            & 
    … …  
    12221544    integer, dimension(6),                  intent(in)  :: indmin, indmax 
    12231545    integer, dimension(6),                  intent(in)  :: pttab_child, petab_child 
    1224     real,    dimension(6),                  intent(in)  :: s_child, s_parent 
    1225     real,    dimension(6),                  intent(in)  :: ds_child, ds_parent 
     1546    real(kind=8),    dimension(6),                  intent(in)  :: s_child, s_parent 
     1547    real(kind=8),    dimension(6),                  intent(in)  :: ds_child, ds_parent 
    12261548    real,    dimension(                 & 
    12271549        indmin(1):indmax(1),            & 
    … …  
    13011623    REAL, DIMENSION(indmin:indmax),           INTENT(IN)    :: parenttab 
    13021624    REAL, DIMENSION(pttab_child:petab_child), INTENT(OUT)   :: childtab 
    1303     REAL                                                    :: s_parent, s_child 
    1304     REAL                                                    :: ds_parent,ds_child 
     1625    REAL(kind=8)                                            :: s_parent, s_child 
     1626    REAL(kind=8)                                            :: ds_parent,ds_child 
    13051627! 
    13061628    if ( (indmin == indmax) .and. (pttab_child == petab_child) ) then 
    … …  
    13711693!--------------------------------------------------------------------------------------------------- 
    13721694function Agrif_Find_list_interp ( list_interp, pttab, petab, pttab_Child, pttab_Parent,     & 
    1373                                     nbdim, indmin, indmax, indminglob,  indmaxglob,         & 
     1695                                    nbdim, indmin, indmax, indmin_required_p, indmax_required_p, & 
     1696                                    indminglob,  indmaxglob,         & 
    13741697                                    pttruetab, cetruetab, memberin                          & 
    13751698#if defined AGRIF_MPI 
    … …  
    13821705    integer,                       intent(in)  :: nbdim 
    13831706    integer, dimension(nbdim),     intent(in)  :: pttab, petab, pttab_Child, pttab_Parent 
    1384     integer, dimension(nbdim),     intent(out) :: indmin, indmax 
     1707    integer, dimension(nbdim),     intent(out) :: indmin, indmax, indmin_required_p, indmax_required_p 
    13851708    integer, dimension(nbdim),     intent(out) :: indminglob, indmaxglob 
    13861709    integer, dimension(nbdim),     intent(out) :: pttruetab, cetruetab 
    … …  
    14211744        indmin = pil % indmin(1:nbdim) 
    14221745        indmax = pil % indmax(1:nbdim) 
     1746        indmin_required_p = pil % indmin_required_p(1:nbdim) 
     1747        indmax_required_p = pil % indmax_required_p(1:nbdim) 
    14231748 
    14241749        pttruetab = pil % pttruetab(1:nbdim) 
    … …  
    14521777!--------------------------------------------------------------------------------------------------- 
    14531778subroutine Agrif_AddTo_list_interp ( list_interp, pttab, petab, pttab_Child, pttab_Parent,  & 
    1454                                      indmin, indmax, indminglob, indmaxglob,                & 
     1779                                     indmin, indmax, indmin_required_p, indmax_required_p,  & 
     1780                                     indminglob, indmaxglob,                & 
    14551781                                     pttruetab, cetruetab,                                  & 
    14561782                                     memberin, nbdim                                        & 
    … …  
    14661792    integer                                 :: nbdim 
    14671793    integer, dimension(nbdim)               :: pttab, petab, pttab_Child, pttab_Parent 
    1468     integer, dimension(nbdim)               :: indmin,indmax 
     1794    integer, dimension(nbdim)               :: indmin,indmax, indmin_required_p, indmax_required_p 
    14691795    integer, dimension(nbdim)               :: indminglob, indmaxglob 
    14701796    integer, dimension(nbdim)               :: pttruetab, cetruetab 
    … …  
    14951821    pil % indmin(1:nbdim) = indmin(1:nbdim) 
    14961822    pil % indmax(1:nbdim) = indmax(1:nbdim) 
     1823 
     1824    pil % indmin_required_p(1:nbdim) = indmin_required_p(1:nbdim) 
     1825    pil % indmax_required_p(1:nbdim) = indmax_required_p(1:nbdim) 
    14971826 
    14981827    pil % memberin = memberin 
Note: See TracChangeset for help on using the changeset viewer.