init_bc Subroutine

public subroutine init_bc()

Uses

  • proc~~init_bc~~UsesGraph proc~init_bc init_bc mod_physics mod_physics proc~init_bc->mod_physics module~mod_comm_lib mod_comm_lib proc~init_bc->module~mod_comm_lib module~mod_global_parameters mod_global_parameters proc~init_bc->module~mod_global_parameters module~mod_basic_types mod_basic_types module~mod_global_parameters->module~mod_basic_types module~mod_connectivity mod_connectivity module~mod_global_parameters->module~mod_connectivity module~mod_constants mod_constants module~mod_global_parameters->module~mod_constants module~mod_physicaldata mod_physicaldata module~mod_global_parameters->module~mod_physicaldata module~mod_variables mod_variables module~mod_global_parameters->module~mod_variables mpi mpi module~mod_global_parameters->mpi module~mod_variables->module~mod_basic_types

Arguments

None

Calls

proc~~init_bc~~CallsGraph proc~init_bc init_bc proc~mpistop mpistop proc~init_bc->proc~mpistop mpi_abort mpi_abort proc~mpistop->mpi_abort

Called by

proc~~init_bc~~CalledByGraph proc~init_bc init_bc proc~initialize_vars initialize_vars proc~initialize_vars->proc~init_bc proc~initialize_amrvac initialize_amrvac proc~initialize_amrvac->proc~initialize_vars proc~main main proc~main->proc~initialize_amrvac program~agile agile program~agile->proc~main

Source Code

  subroutine init_bc()
    use mod_global_parameters
    use mod_physics, only: phys_req_diagonal, physics_type
    use mod_comm_lib, only: mpistop

    integer :: nghostcellsCo, interpolation_order
    integer :: nx1,nx2,nx3, nxCo1,nxCo2,nxCo3, ixGmin1,ixGmin2,ixGmin3,ixGmax1,&
         ixGmax2,ixGmax3, i1,i2,i3, ic1,ic2,ic3, inc1,inc2,inc3, idir

    ixGmin1=ixGlo1;ixGmin2=ixGlo2;ixGmin3=ixGlo3;ixGmax1=ixGhi1
    ixGmax2=ixGhi2;ixGmax3=ixGhi3;
    ixMmin1=ixGmin1+nghostcells;ixMmin2=ixGmin2+nghostcells
    ixMmin3=ixGmin3+nghostcells;ixMmax1=ixGmax1-nghostcells
    ixMmax2=ixGmax2-nghostcells;ixMmax3=ixGmax3-nghostcells;
    ixCoGmin1=1;ixCoGmin2=1;ixCoGmin3=1;
    ixCoGmax1=(ixGhi1-2*nghostcells)/2+2*nghostcells
    ixCoGmax2=(ixGhi2-2*nghostcells)/2+2*nghostcells
    ixCoGmax3=(ixGhi3-2*nghostcells)/2+2*nghostcells;
    ixCoGsmin1=0;ixCoGsmin2=0;ixCoGsmin3=0;
    ixCoGsmax1=ixCoGmax1;ixCoGsmax2=ixCoGmax2;ixCoGsmax3=ixCoGmax3;

    ixCoMmin1=ixCoGmin1+nghostcells;ixCoMmin2=ixCoGmin2+nghostcells
    ixCoMmin3=ixCoGmin3+nghostcells;ixCoMmax1=ixCoGmax1-nghostcells
    ixCoMmax2=ixCoGmax2-nghostcells;ixCoMmax3=ixCoGmax3-nghostcells;

    nx1=ixMmax1-ixMmin1+1;nx2=ixMmax2-ixMmin2+1;nx3=ixMmax3-ixMmin3+1;
    nxCo1=nx1/2;nxCo2=nx2/2;nxCo3=nx3/2;

    if(ghost_copy) then
       interpolation_order=1
    else
       interpolation_order=2
    end if
    nghostcellsCo=int((nghostcells+1)/2)

    if (nghostcellsCo+interpolation_order-1>nghostcells) then
       call mpistop("interpolation order for prolongation in getbc too high")
    end if

    ! (iib,i) index has following meanings: iib = 0 means it is not at any physical boundary
    ! iib=-1 means it is at the minimum side of a physical boundary
    ! iib= 1 means it is at the maximum side of a physical boundary
    ! i=-1 means subregion prepared for the neighbor at its minimum side
    ! i= 1 means subregion prepared for the neighbor at its maximum side

    ixS_srl_min1(:,-1)=ixMmin1
    ixS_srl_min1(:, 0)=ixMmin1
    ixS_srl_min1(:, 1)=ixMmax1+1-nghostcells
    ixS_srl_max1(:,-1)=ixMmin1-1+nghostcells
    ixS_srl_max1(:, 0)=ixMmax1
    ixS_srl_max1(:, 1)=ixMmax1

    ixR_srl_min1(:,-1)=1
    ixR_srl_min1(:, 0)=ixMmin1
    ixR_srl_min1(:, 1)=ixMmax1+1
    ixR_srl_max1(:,-1)=nghostcells
    ixR_srl_max1(:, 0)=ixMmax1
    ixR_srl_max1(:, 1)=ixGmax1

    ixS_r_min1(:,-1)=ixCoMmin1
    ixS_r_min1(:, 0)=ixCoMmin1
    ixS_r_min1(:, 1)=ixCoMmax1+1-nghostcells
    ixS_r_max1(:,-1)=ixCoMmin1-1+nghostcells
    ixS_r_max1(:, 0)=ixCoMmax1
    ixS_r_max1(:, 1)=ixCoMmax1

    ixR_r_min1(:, 0)=1
    ixR_r_min1(:, 1)=ixMmin1
    ixR_r_min1(:, 2)=ixMmin1+nxCo1
    ixR_r_min1(:, 3)=ixMmax1+1
    ixR_r_max1(:, 0)=nghostcells
    ixR_r_max1(:, 1)=ixMmin1-1+nxCo1
    ixR_r_max1(:, 2)=ixMmax1
    ixR_r_max1(:, 3)=ixGmax1

    ixS_p_min1(:, 0)=ixMmin1-(interpolation_order-1)
    ixS_p_min1(:, 1)=ixMmin1-(interpolation_order-1)
    ixS_p_min1(:, 2)=ixMmin1+nxCo1-nghostcellsCo-(interpolation_order-1)
    ixS_p_min1(:, 3)=ixMmax1+1-nghostcellsCo-(interpolation_order-1)
    ixS_p_max1(:, 0)=ixMmin1-1+nghostcellsCo+(interpolation_order-1)
    ixS_p_max1(:, 1)=ixMmin1-1+nxCo1+nghostcellsCo+(interpolation_order-1)
    ixS_p_max1(:, 2)=ixMmax1+(interpolation_order-1)
    ixS_p_max1(:, 3)=ixMmax1+(interpolation_order-1)

    if(.not.phys_req_diagonal) then
       ! exclude ghost-cell region when diagonal cells are unknown
       ixS_p_min1(:, 0)=ixMmin1
       ixS_p_max1(:, 3)=ixMmax1
       ixS_p_max1(:, 1)=ixMmin1-1+nxCo1+(interpolation_order-1)
       ixS_p_min1(:, 2)=ixMmin1+nxCo1-(interpolation_order-1)
    end if

    ixR_p_min1(:, 0)=ixCoMmin1-nghostcellsCo-(interpolation_order-1)
    ixR_p_min1(:, 1)=ixCoMmin1-(interpolation_order-1)
    ixR_p_min1(:, 2)=ixCoMmin1-nghostcellsCo-(interpolation_order-1)
    ixR_p_min1(:, 3)=ixCoMmax1+1-(interpolation_order-1)
    ixR_p_max1(:, 0)=nghostcells+(interpolation_order-1)
    ixR_p_max1(:, 1)=ixCoMmax1+nghostcellsCo+(interpolation_order-1)
    ixR_p_max1(:, 2)=ixCoMmax1+(interpolation_order-1)
    ixR_p_max1(:, 3)=ixCoMmax1+nghostcellsCo+(interpolation_order-1)

    if(.not.phys_req_diagonal) then
       ixR_p_max1(:, 0)=nghostcells
       ixR_p_min1(:, 3)=ixCoMmax1+1
       ixR_p_max1(:, 1)=ixCoMmax1+(interpolation_order-1)
       ixR_p_min1(:, 2)=ixCoMmin1-(interpolation_order-1)
    end if



    ixS_srl_min2(:,-1)=ixMmin2
    ixS_srl_min2(:, 0)=ixMmin2
    ixS_srl_min2(:, 1)=ixMmax2+1-nghostcells
    ixS_srl_max2(:,-1)=ixMmin2-1+nghostcells
    ixS_srl_max2(:, 0)=ixMmax2
    ixS_srl_max2(:, 1)=ixMmax2

    ixR_srl_min2(:,-1)=1
    ixR_srl_min2(:, 0)=ixMmin2
    ixR_srl_min2(:, 1)=ixMmax2+1
    ixR_srl_max2(:,-1)=nghostcells
    ixR_srl_max2(:, 0)=ixMmax2
    ixR_srl_max2(:, 1)=ixGmax2

    ixS_r_min2(:,-1)=ixCoMmin2
    ixS_r_min2(:, 0)=ixCoMmin2
    ixS_r_min2(:, 1)=ixCoMmax2+1-nghostcells
    ixS_r_max2(:,-1)=ixCoMmin2-1+nghostcells
    ixS_r_max2(:, 0)=ixCoMmax2
    ixS_r_max2(:, 1)=ixCoMmax2

    ixR_r_min2(:, 0)=1
    ixR_r_min2(:, 1)=ixMmin2
    ixR_r_min2(:, 2)=ixMmin2+nxCo2
    ixR_r_min2(:, 3)=ixMmax2+1
    ixR_r_max2(:, 0)=nghostcells
    ixR_r_max2(:, 1)=ixMmin2-1+nxCo2
    ixR_r_max2(:, 2)=ixMmax2
    ixR_r_max2(:, 3)=ixGmax2

    ixS_p_min2(:, 0)=ixMmin2-(interpolation_order-1)
    ixS_p_min2(:, 1)=ixMmin2-(interpolation_order-1)
    ixS_p_min2(:, 2)=ixMmin2+nxCo2-nghostcellsCo-(interpolation_order-1)
    ixS_p_min2(:, 3)=ixMmax2+1-nghostcellsCo-(interpolation_order-1)
    ixS_p_max2(:, 0)=ixMmin2-1+nghostcellsCo+(interpolation_order-1)
    ixS_p_max2(:, 1)=ixMmin2-1+nxCo2+nghostcellsCo+(interpolation_order-1)
    ixS_p_max2(:, 2)=ixMmax2+(interpolation_order-1)
    ixS_p_max2(:, 3)=ixMmax2+(interpolation_order-1)

    if(.not.phys_req_diagonal) then
       ! exclude ghost-cell region when diagonal cells are unknown
       ixS_p_min2(:, 0)=ixMmin2
       ixS_p_max2(:, 3)=ixMmax2
       ixS_p_max2(:, 1)=ixMmin2-1+nxCo2+(interpolation_order-1)
       ixS_p_min2(:, 2)=ixMmin2+nxCo2-(interpolation_order-1)
    end if

    ixR_p_min2(:, 0)=ixCoMmin2-nghostcellsCo-(interpolation_order-1)
    ixR_p_min2(:, 1)=ixCoMmin2-(interpolation_order-1)
    ixR_p_min2(:, 2)=ixCoMmin2-nghostcellsCo-(interpolation_order-1)
    ixR_p_min2(:, 3)=ixCoMmax2+1-(interpolation_order-1)
    ixR_p_max2(:, 0)=nghostcells+(interpolation_order-1)
    ixR_p_max2(:, 1)=ixCoMmax2+nghostcellsCo+(interpolation_order-1)
    ixR_p_max2(:, 2)=ixCoMmax2+(interpolation_order-1)
    ixR_p_max2(:, 3)=ixCoMmax2+nghostcellsCo+(interpolation_order-1)

    if(.not.phys_req_diagonal) then
       ixR_p_max2(:, 0)=nghostcells
       ixR_p_min2(:, 3)=ixCoMmax2+1
       ixR_p_max2(:, 1)=ixCoMmax2+(interpolation_order-1)
       ixR_p_min2(:, 2)=ixCoMmin2-(interpolation_order-1)
    end if



    ixS_srl_min3(:,-1)=ixMmin3
    ixS_srl_min3(:, 0)=ixMmin3
    ixS_srl_min3(:, 1)=ixMmax3+1-nghostcells
    ixS_srl_max3(:,-1)=ixMmin3-1+nghostcells
    ixS_srl_max3(:, 0)=ixMmax3
    ixS_srl_max3(:, 1)=ixMmax3

    ixR_srl_min3(:,-1)=1
    ixR_srl_min3(:, 0)=ixMmin3
    ixR_srl_min3(:, 1)=ixMmax3+1
    ixR_srl_max3(:,-1)=nghostcells
    ixR_srl_max3(:, 0)=ixMmax3
    ixR_srl_max3(:, 1)=ixGmax3

    ixS_r_min3(:,-1)=ixCoMmin3
    ixS_r_min3(:, 0)=ixCoMmin3
    ixS_r_min3(:, 1)=ixCoMmax3+1-nghostcells
    ixS_r_max3(:,-1)=ixCoMmin3-1+nghostcells
    ixS_r_max3(:, 0)=ixCoMmax3
    ixS_r_max3(:, 1)=ixCoMmax3

    ixR_r_min3(:, 0)=1
    ixR_r_min3(:, 1)=ixMmin3
    ixR_r_min3(:, 2)=ixMmin3+nxCo3
    ixR_r_min3(:, 3)=ixMmax3+1
    ixR_r_max3(:, 0)=nghostcells
    ixR_r_max3(:, 1)=ixMmin3-1+nxCo3
    ixR_r_max3(:, 2)=ixMmax3
    ixR_r_max3(:, 3)=ixGmax3

    ixS_p_min3(:, 0)=ixMmin3-(interpolation_order-1)
    ixS_p_min3(:, 1)=ixMmin3-(interpolation_order-1)
    ixS_p_min3(:, 2)=ixMmin3+nxCo3-nghostcellsCo-(interpolation_order-1)
    ixS_p_min3(:, 3)=ixMmax3+1-nghostcellsCo-(interpolation_order-1)
    ixS_p_max3(:, 0)=ixMmin3-1+nghostcellsCo+(interpolation_order-1)
    ixS_p_max3(:, 1)=ixMmin3-1+nxCo3+nghostcellsCo+(interpolation_order-1)
    ixS_p_max3(:, 2)=ixMmax3+(interpolation_order-1)
    ixS_p_max3(:, 3)=ixMmax3+(interpolation_order-1)

    if(.not.phys_req_diagonal) then
       ! exclude ghost-cell region when diagonal cells are unknown
       ixS_p_min3(:, 0)=ixMmin3
       ixS_p_max3(:, 3)=ixMmax3
       ixS_p_max3(:, 1)=ixMmin3-1+nxCo3+(interpolation_order-1)
       ixS_p_min3(:, 2)=ixMmin3+nxCo3-(interpolation_order-1)
    end if

    ixR_p_min3(:, 0)=ixCoMmin3-nghostcellsCo-(interpolation_order-1)
    ixR_p_min3(:, 1)=ixCoMmin3-(interpolation_order-1)
    ixR_p_min3(:, 2)=ixCoMmin3-nghostcellsCo-(interpolation_order-1)
    ixR_p_min3(:, 3)=ixCoMmax3+1-(interpolation_order-1)
    ixR_p_max3(:, 0)=nghostcells+(interpolation_order-1)
    ixR_p_max3(:, 1)=ixCoMmax3+nghostcellsCo+(interpolation_order-1)
    ixR_p_max3(:, 2)=ixCoMmax3+(interpolation_order-1)
    ixR_p_max3(:, 3)=ixCoMmax3+nghostcellsCo+(interpolation_order-1)

    if(.not.phys_req_diagonal) then
       ixR_p_max3(:, 0)=nghostcells
       ixR_p_min3(:, 3)=ixCoMmax3+1
       ixR_p_max3(:, 1)=ixCoMmax3+(interpolation_order-1)
       ixR_p_min3(:, 2)=ixCoMmin3-(interpolation_order-1)
    end if

    if (stagger_grid) then
       ! Staggered (face-allocated) variables
       do idir=1,ndim
          ixS_srl_stg_min1(idir,-1)=ixMmin1-kr(idir,1)
          ixS_srl_stg_max1(idir,-1)=ixMmin1-1+nghostcells
          ixS_srl_stg_min1(idir,0) =ixMmin1-kr(idir,1)
          ixS_srl_stg_max1(idir,0) =ixMmax1
          ixS_srl_stg_min1(idir,1) =ixMmax1-nghostcells+1-kr(idir,1)
          ixS_srl_stg_max1(idir,1) =ixMmax1

          ixR_srl_stg_min1(idir,-1)=1-kr(idir,1)
          ixR_srl_stg_max1(idir,-1)=nghostcells
          ixR_srl_stg_min1(idir,0) =ixMmin1-kr(idir,1)
          ixR_srl_stg_max1(idir,0) =ixMmax1
          ixR_srl_stg_min1(idir,1) =ixMmax1+1-kr(idir,1)
          ixR_srl_stg_max1(idir,1) =ixGmax1

          ixS_r_stg_min1(idir,-1)=ixCoMmin1-kr(idir,1)
          ixS_r_stg_max1(idir,-1)=ixCoMmin1-1+nghostcells
          ixS_r_stg_min1(idir,0) =ixCoMmin1-kr(idir,1)
          ixS_r_stg_max1(idir,0) =ixCoMmax1
          ixS_r_stg_min1(idir,1) =ixCoMmax1+1-nghostcells-kr(idir,1)
          ixS_r_stg_max1(idir,1) =ixCoMmax1

          ixR_r_stg_min1(idir,0)=1-kr(idir,1)
          ixR_r_stg_max1(idir,0)=nghostcells
          ixR_r_stg_min1(idir,1)=ixMmin1-kr(idir,1)
          ixR_r_stg_max1(idir,1)=ixMmin1-1+nxCo1
          ixR_r_stg_min1(idir,2)=ixMmin1+nxCo1-kr(idir,1)
          ixR_r_stg_max1(idir,2)=ixMmax1
          ixR_r_stg_min1(idir,3)=ixMmax1+1-kr(idir,1)
          ixR_r_stg_max1(idir,3)=ixGmax1

          ixS_srl_stg_min2(idir,-1)=ixMmin2-kr(idir,2)
          ixS_srl_stg_max2(idir,-1)=ixMmin2-1+nghostcells
          ixS_srl_stg_min2(idir,0) =ixMmin2-kr(idir,2)
          ixS_srl_stg_max2(idir,0) =ixMmax2
          ixS_srl_stg_min2(idir,1) =ixMmax2-nghostcells+1-kr(idir,2)
          ixS_srl_stg_max2(idir,1) =ixMmax2

          ixR_srl_stg_min2(idir,-1)=1-kr(idir,2)
          ixR_srl_stg_max2(idir,-1)=nghostcells
          ixR_srl_stg_min2(idir,0) =ixMmin2-kr(idir,2)
          ixR_srl_stg_max2(idir,0) =ixMmax2
          ixR_srl_stg_min2(idir,1) =ixMmax2+1-kr(idir,2)
          ixR_srl_stg_max2(idir,1) =ixGmax2

          ixS_r_stg_min2(idir,-1)=ixCoMmin2-kr(idir,2)
          ixS_r_stg_max2(idir,-1)=ixCoMmin2-1+nghostcells
          ixS_r_stg_min2(idir,0) =ixCoMmin2-kr(idir,2)
          ixS_r_stg_max2(idir,0) =ixCoMmax2
          ixS_r_stg_min2(idir,1) =ixCoMmax2+1-nghostcells-kr(idir,2)
          ixS_r_stg_max2(idir,1) =ixCoMmax2

          ixR_r_stg_min2(idir,0)=1-kr(idir,2)
          ixR_r_stg_max2(idir,0)=nghostcells
          ixR_r_stg_min2(idir,1)=ixMmin2-kr(idir,2)
          ixR_r_stg_max2(idir,1)=ixMmin2-1+nxCo2
          ixR_r_stg_min2(idir,2)=ixMmin2+nxCo2-kr(idir,2)
          ixR_r_stg_max2(idir,2)=ixMmax2
          ixR_r_stg_min2(idir,3)=ixMmax2+1-kr(idir,2)
          ixR_r_stg_max2(idir,3)=ixGmax2

          ixS_srl_stg_min3(idir,-1)=ixMmin3-kr(idir,3)
          ixS_srl_stg_max3(idir,-1)=ixMmin3-1+nghostcells
          ixS_srl_stg_min3(idir,0) =ixMmin3-kr(idir,3)
          ixS_srl_stg_max3(idir,0) =ixMmax3
          ixS_srl_stg_min3(idir,1) =ixMmax3-nghostcells+1-kr(idir,3)
          ixS_srl_stg_max3(idir,1) =ixMmax3

          ixR_srl_stg_min3(idir,-1)=1-kr(idir,3)
          ixR_srl_stg_max3(idir,-1)=nghostcells
          ixR_srl_stg_min3(idir,0) =ixMmin3-kr(idir,3)
          ixR_srl_stg_max3(idir,0) =ixMmax3
          ixR_srl_stg_min3(idir,1) =ixMmax3+1-kr(idir,3)
          ixR_srl_stg_max3(idir,1) =ixGmax3

          ixS_r_stg_min3(idir,-1)=ixCoMmin3-kr(idir,3)
          ixS_r_stg_max3(idir,-1)=ixCoMmin3-1+nghostcells
          ixS_r_stg_min3(idir,0) =ixCoMmin3-kr(idir,3)
          ixS_r_stg_max3(idir,0) =ixCoMmax3
          ixS_r_stg_min3(idir,1) =ixCoMmax3+1-nghostcells-kr(idir,3)
          ixS_r_stg_max3(idir,1) =ixCoMmax3

          ixR_r_stg_min3(idir,0)=1-kr(idir,3)
          ixR_r_stg_max3(idir,0)=nghostcells
          ixR_r_stg_min3(idir,1)=ixMmin3-kr(idir,3)
          ixR_r_stg_max3(idir,1)=ixMmin3-1+nxCo3
          ixR_r_stg_min3(idir,2)=ixMmin3+nxCo3-kr(idir,3)
          ixR_r_stg_max3(idir,2)=ixMmax3
          ixR_r_stg_min3(idir,3)=ixMmax3+1-kr(idir,3)
          ixR_r_stg_max3(idir,3)=ixGmax3

          if (idir==1) then
             ! Parallel components

             ixS_p_stg_min1(idir,0)=ixMmin1-1 ! -1 to make redundant
             ixS_p_stg_max1(idir,0)=ixMmin1-1+nghostcellsCo
             ixS_p_stg_min1(idir,1)=ixMmin1-1 ! -1 to make redundant
             ixS_p_stg_max1(idir,1)=ixMmin1-1+nxCo1+nghostcellsCo
             ixS_p_stg_min1(idir,2)=ixMmax1-nxCo1-nghostcellsCo
             ixS_p_stg_max1(idir,2)=ixMmax1
             ixS_p_stg_min1(idir,3)=ixMmax1-nghostcellsCo
             ixS_p_stg_max1(idir,3)=ixMmax1

             ixR_p_stg_min1(idir,0)=ixCoMmin1-1-nghostcellsCo
             ixR_p_stg_max1(idir,0)=ixCoMmin1-1
             ixR_p_stg_min1(idir,1)=ixCoMmin1-1 ! -1 to make redundant
             ixR_p_stg_max1(idir,1)=ixCoMmax1+nghostcellsCo
             ixR_p_stg_min1(idir,2)=ixCoMmin1-1-nghostcellsCo
             ixR_p_stg_max1(idir,2)=ixCoMmax1
             ixR_p_stg_min1(idir,3)=ixCoMmax1+1-1 ! -1 to make redundant
             ixR_p_stg_max1(idir,3)=ixCoMmax1+nghostcellsCo


             ixS_p_stg_min1(idir,0)=ixMmin1-1 ! -1 to make redundant
             ixS_p_stg_max1(idir,0)=ixMmin1-1+nghostcellsCo
             ixS_p_stg_min1(idir,1)=ixMmin1-1 ! -1 to make redundant
             ixS_p_stg_max1(idir,1)=ixMmin1-1+nxCo1+nghostcellsCo
             ixS_p_stg_min1(idir,2)=ixMmax1-nxCo1-nghostcellsCo
             ixS_p_stg_max1(idir,2)=ixMmax1
             ixS_p_stg_min1(idir,3)=ixMmax1-nghostcellsCo
             ixS_p_stg_max1(idir,3)=ixMmax1

             ixR_p_stg_min1(idir,0)=ixCoMmin1-1-nghostcellsCo
             ixR_p_stg_max1(idir,0)=ixCoMmin1-1
             ixR_p_stg_min1(idir,1)=ixCoMmin1-1 ! -1 to make redundant
             ixR_p_stg_max1(idir,1)=ixCoMmax1+nghostcellsCo
             ixR_p_stg_min1(idir,2)=ixCoMmin1-1-nghostcellsCo
             ixR_p_stg_max1(idir,2)=ixCoMmax1
             ixR_p_stg_min1(idir,3)=ixCoMmax1+1-1 ! -1 to make redundant
             ixR_p_stg_max1(idir,3)=ixCoMmax1+nghostcellsCo


             ixS_p_stg_min1(idir,0)=ixMmin1-1 ! -1 to make redundant
             ixS_p_stg_max1(idir,0)=ixMmin1-1+nghostcellsCo
             ixS_p_stg_min1(idir,1)=ixMmin1-1 ! -1 to make redundant
             ixS_p_stg_max1(idir,1)=ixMmin1-1+nxCo1+nghostcellsCo
             ixS_p_stg_min1(idir,2)=ixMmax1-nxCo1-nghostcellsCo
             ixS_p_stg_max1(idir,2)=ixMmax1
             ixS_p_stg_min1(idir,3)=ixMmax1-nghostcellsCo
             ixS_p_stg_max1(idir,3)=ixMmax1

             ixR_p_stg_min1(idir,0)=ixCoMmin1-1-nghostcellsCo
             ixR_p_stg_max1(idir,0)=ixCoMmin1-1
             ixR_p_stg_min1(idir,1)=ixCoMmin1-1 ! -1 to make redundant
             ixR_p_stg_max1(idir,1)=ixCoMmax1+nghostcellsCo
             ixR_p_stg_min1(idir,2)=ixCoMmin1-1-nghostcellsCo
             ixR_p_stg_max1(idir,2)=ixCoMmax1
             ixR_p_stg_min1(idir,3)=ixCoMmax1+1-1 ! -1 to make redundant
             ixR_p_stg_max1(idir,3)=ixCoMmax1+nghostcellsCo

          else

             ! Perpendicular component
             ixS_p_stg_min1(idir,0)=ixMmin1
             ixS_p_stg_max1(idir,0)=ixMmin1-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min1(idir,1)=ixMmin1
             ixS_p_stg_max1(idir,1)=ixMmin1-1+nxCo1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min1(idir,2)=ixMmax1+&
                  1-nxCo1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max1(idir,2)=ixMmax1
             ixS_p_stg_min1(idir,3)=ixMmax1+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max1(idir,3)=ixMmax1

             ixR_p_stg_min1(idir,0)=ixCoMmin1-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max1(idir,0)=ixCoMmin1-1
             ixR_p_stg_min1(idir,1)=ixCoMmin1
             ixR_p_stg_max1(idir,1)=ixCoMmax1+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min1(idir,2)=ixCoMmin1-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max1(idir,2)=ixCoMmax1
             ixR_p_stg_min1(idir,3)=ixCoMmax1+1
             ixR_p_stg_max1(idir,3)=ixCoMmax1+nghostcellsCo+&
                  (interpolation_order-1)


             ! Perpendicular component
             ixS_p_stg_min1(idir,0)=ixMmin1
             ixS_p_stg_max1(idir,0)=ixMmin1-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min1(idir,1)=ixMmin1
             ixS_p_stg_max1(idir,1)=ixMmin1-1+nxCo1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min1(idir,2)=ixMmax1+&
                  1-nxCo1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max1(idir,2)=ixMmax1
             ixS_p_stg_min1(idir,3)=ixMmax1+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max1(idir,3)=ixMmax1

             ixR_p_stg_min1(idir,0)=ixCoMmin1-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max1(idir,0)=ixCoMmin1-1
             ixR_p_stg_min1(idir,1)=ixCoMmin1
             ixR_p_stg_max1(idir,1)=ixCoMmax1+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min1(idir,2)=ixCoMmin1-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max1(idir,2)=ixCoMmax1
             ixR_p_stg_min1(idir,3)=ixCoMmax1+1
             ixR_p_stg_max1(idir,3)=ixCoMmax1+nghostcellsCo+&
                  (interpolation_order-1)


             ! Perpendicular component
             ixS_p_stg_min1(idir,0)=ixMmin1
             ixS_p_stg_max1(idir,0)=ixMmin1-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min1(idir,1)=ixMmin1
             ixS_p_stg_max1(idir,1)=ixMmin1-1+nxCo1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min1(idir,2)=ixMmax1+&
                  1-nxCo1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max1(idir,2)=ixMmax1
             ixS_p_stg_min1(idir,3)=ixMmax1+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max1(idir,3)=ixMmax1

             ixR_p_stg_min1(idir,0)=ixCoMmin1-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max1(idir,0)=ixCoMmin1-1
             ixR_p_stg_min1(idir,1)=ixCoMmin1
             ixR_p_stg_max1(idir,1)=ixCoMmax1+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min1(idir,2)=ixCoMmin1-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max1(idir,2)=ixCoMmax1
             ixR_p_stg_min1(idir,3)=ixCoMmax1+1
             ixR_p_stg_max1(idir,3)=ixCoMmax1+nghostcellsCo+&
                  (interpolation_order-1)

          end if
          if (idir==2) then
             ! Parallel components

             ixS_p_stg_min2(idir,0)=ixMmin2-1 ! -1 to make redundant
             ixS_p_stg_max2(idir,0)=ixMmin2-1+nghostcellsCo
             ixS_p_stg_min2(idir,1)=ixMmin2-1 ! -1 to make redundant
             ixS_p_stg_max2(idir,1)=ixMmin2-1+nxCo2+nghostcellsCo
             ixS_p_stg_min2(idir,2)=ixMmax2-nxCo2-nghostcellsCo
             ixS_p_stg_max2(idir,2)=ixMmax2
             ixS_p_stg_min2(idir,3)=ixMmax2-nghostcellsCo
             ixS_p_stg_max2(idir,3)=ixMmax2

             ixR_p_stg_min2(idir,0)=ixCoMmin2-1-nghostcellsCo
             ixR_p_stg_max2(idir,0)=ixCoMmin2-1
             ixR_p_stg_min2(idir,1)=ixCoMmin2-1 ! -1 to make redundant
             ixR_p_stg_max2(idir,1)=ixCoMmax2+nghostcellsCo
             ixR_p_stg_min2(idir,2)=ixCoMmin2-1-nghostcellsCo
             ixR_p_stg_max2(idir,2)=ixCoMmax2
             ixR_p_stg_min2(idir,3)=ixCoMmax2+1-1 ! -1 to make redundant
             ixR_p_stg_max2(idir,3)=ixCoMmax2+nghostcellsCo


             ixS_p_stg_min2(idir,0)=ixMmin2-1 ! -1 to make redundant
             ixS_p_stg_max2(idir,0)=ixMmin2-1+nghostcellsCo
             ixS_p_stg_min2(idir,1)=ixMmin2-1 ! -1 to make redundant
             ixS_p_stg_max2(idir,1)=ixMmin2-1+nxCo2+nghostcellsCo
             ixS_p_stg_min2(idir,2)=ixMmax2-nxCo2-nghostcellsCo
             ixS_p_stg_max2(idir,2)=ixMmax2
             ixS_p_stg_min2(idir,3)=ixMmax2-nghostcellsCo
             ixS_p_stg_max2(idir,3)=ixMmax2

             ixR_p_stg_min2(idir,0)=ixCoMmin2-1-nghostcellsCo
             ixR_p_stg_max2(idir,0)=ixCoMmin2-1
             ixR_p_stg_min2(idir,1)=ixCoMmin2-1 ! -1 to make redundant
             ixR_p_stg_max2(idir,1)=ixCoMmax2+nghostcellsCo
             ixR_p_stg_min2(idir,2)=ixCoMmin2-1-nghostcellsCo
             ixR_p_stg_max2(idir,2)=ixCoMmax2
             ixR_p_stg_min2(idir,3)=ixCoMmax2+1-1 ! -1 to make redundant
             ixR_p_stg_max2(idir,3)=ixCoMmax2+nghostcellsCo


             ixS_p_stg_min2(idir,0)=ixMmin2-1 ! -1 to make redundant
             ixS_p_stg_max2(idir,0)=ixMmin2-1+nghostcellsCo
             ixS_p_stg_min2(idir,1)=ixMmin2-1 ! -1 to make redundant
             ixS_p_stg_max2(idir,1)=ixMmin2-1+nxCo2+nghostcellsCo
             ixS_p_stg_min2(idir,2)=ixMmax2-nxCo2-nghostcellsCo
             ixS_p_stg_max2(idir,2)=ixMmax2
             ixS_p_stg_min2(idir,3)=ixMmax2-nghostcellsCo
             ixS_p_stg_max2(idir,3)=ixMmax2

             ixR_p_stg_min2(idir,0)=ixCoMmin2-1-nghostcellsCo
             ixR_p_stg_max2(idir,0)=ixCoMmin2-1
             ixR_p_stg_min2(idir,1)=ixCoMmin2-1 ! -1 to make redundant
             ixR_p_stg_max2(idir,1)=ixCoMmax2+nghostcellsCo
             ixR_p_stg_min2(idir,2)=ixCoMmin2-1-nghostcellsCo
             ixR_p_stg_max2(idir,2)=ixCoMmax2
             ixR_p_stg_min2(idir,3)=ixCoMmax2+1-1 ! -1 to make redundant
             ixR_p_stg_max2(idir,3)=ixCoMmax2+nghostcellsCo

          else

             ! Perpendicular component
             ixS_p_stg_min2(idir,0)=ixMmin2
             ixS_p_stg_max2(idir,0)=ixMmin2-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min2(idir,1)=ixMmin2
             ixS_p_stg_max2(idir,1)=ixMmin2-1+nxCo2+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min2(idir,2)=ixMmax2+&
                  1-nxCo2-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max2(idir,2)=ixMmax2
             ixS_p_stg_min2(idir,3)=ixMmax2+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max2(idir,3)=ixMmax2

             ixR_p_stg_min2(idir,0)=ixCoMmin2-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max2(idir,0)=ixCoMmin2-1
             ixR_p_stg_min2(idir,1)=ixCoMmin2
             ixR_p_stg_max2(idir,1)=ixCoMmax2+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min2(idir,2)=ixCoMmin2-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max2(idir,2)=ixCoMmax2
             ixR_p_stg_min2(idir,3)=ixCoMmax2+1
             ixR_p_stg_max2(idir,3)=ixCoMmax2+nghostcellsCo+&
                  (interpolation_order-1)


             ! Perpendicular component
             ixS_p_stg_min2(idir,0)=ixMmin2
             ixS_p_stg_max2(idir,0)=ixMmin2-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min2(idir,1)=ixMmin2
             ixS_p_stg_max2(idir,1)=ixMmin2-1+nxCo2+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min2(idir,2)=ixMmax2+&
                  1-nxCo2-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max2(idir,2)=ixMmax2
             ixS_p_stg_min2(idir,3)=ixMmax2+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max2(idir,3)=ixMmax2

             ixR_p_stg_min2(idir,0)=ixCoMmin2-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max2(idir,0)=ixCoMmin2-1
             ixR_p_stg_min2(idir,1)=ixCoMmin2
             ixR_p_stg_max2(idir,1)=ixCoMmax2+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min2(idir,2)=ixCoMmin2-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max2(idir,2)=ixCoMmax2
             ixR_p_stg_min2(idir,3)=ixCoMmax2+1
             ixR_p_stg_max2(idir,3)=ixCoMmax2+nghostcellsCo+&
                  (interpolation_order-1)


             ! Perpendicular component
             ixS_p_stg_min2(idir,0)=ixMmin2
             ixS_p_stg_max2(idir,0)=ixMmin2-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min2(idir,1)=ixMmin2
             ixS_p_stg_max2(idir,1)=ixMmin2-1+nxCo2+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min2(idir,2)=ixMmax2+&
                  1-nxCo2-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max2(idir,2)=ixMmax2
             ixS_p_stg_min2(idir,3)=ixMmax2+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max2(idir,3)=ixMmax2

             ixR_p_stg_min2(idir,0)=ixCoMmin2-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max2(idir,0)=ixCoMmin2-1
             ixR_p_stg_min2(idir,1)=ixCoMmin2
             ixR_p_stg_max2(idir,1)=ixCoMmax2+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min2(idir,2)=ixCoMmin2-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max2(idir,2)=ixCoMmax2
             ixR_p_stg_min2(idir,3)=ixCoMmax2+1
             ixR_p_stg_max2(idir,3)=ixCoMmax2+nghostcellsCo+&
                  (interpolation_order-1)

          end if
          if (idir==3) then
             ! Parallel components

             ixS_p_stg_min3(idir,0)=ixMmin3-1 ! -1 to make redundant
             ixS_p_stg_max3(idir,0)=ixMmin3-1+nghostcellsCo
             ixS_p_stg_min3(idir,1)=ixMmin3-1 ! -1 to make redundant
             ixS_p_stg_max3(idir,1)=ixMmin3-1+nxCo3+nghostcellsCo
             ixS_p_stg_min3(idir,2)=ixMmax3-nxCo3-nghostcellsCo
             ixS_p_stg_max3(idir,2)=ixMmax3
             ixS_p_stg_min3(idir,3)=ixMmax3-nghostcellsCo
             ixS_p_stg_max3(idir,3)=ixMmax3

             ixR_p_stg_min3(idir,0)=ixCoMmin3-1-nghostcellsCo
             ixR_p_stg_max3(idir,0)=ixCoMmin3-1
             ixR_p_stg_min3(idir,1)=ixCoMmin3-1 ! -1 to make redundant
             ixR_p_stg_max3(idir,1)=ixCoMmax3+nghostcellsCo
             ixR_p_stg_min3(idir,2)=ixCoMmin3-1-nghostcellsCo
             ixR_p_stg_max3(idir,2)=ixCoMmax3
             ixR_p_stg_min3(idir,3)=ixCoMmax3+1-1 ! -1 to make redundant
             ixR_p_stg_max3(idir,3)=ixCoMmax3+nghostcellsCo


             ixS_p_stg_min3(idir,0)=ixMmin3-1 ! -1 to make redundant
             ixS_p_stg_max3(idir,0)=ixMmin3-1+nghostcellsCo
             ixS_p_stg_min3(idir,1)=ixMmin3-1 ! -1 to make redundant
             ixS_p_stg_max3(idir,1)=ixMmin3-1+nxCo3+nghostcellsCo
             ixS_p_stg_min3(idir,2)=ixMmax3-nxCo3-nghostcellsCo
             ixS_p_stg_max3(idir,2)=ixMmax3
             ixS_p_stg_min3(idir,3)=ixMmax3-nghostcellsCo
             ixS_p_stg_max3(idir,3)=ixMmax3

             ixR_p_stg_min3(idir,0)=ixCoMmin3-1-nghostcellsCo
             ixR_p_stg_max3(idir,0)=ixCoMmin3-1
             ixR_p_stg_min3(idir,1)=ixCoMmin3-1 ! -1 to make redundant
             ixR_p_stg_max3(idir,1)=ixCoMmax3+nghostcellsCo
             ixR_p_stg_min3(idir,2)=ixCoMmin3-1-nghostcellsCo
             ixR_p_stg_max3(idir,2)=ixCoMmax3
             ixR_p_stg_min3(idir,3)=ixCoMmax3+1-1 ! -1 to make redundant
             ixR_p_stg_max3(idir,3)=ixCoMmax3+nghostcellsCo


             ixS_p_stg_min3(idir,0)=ixMmin3-1 ! -1 to make redundant
             ixS_p_stg_max3(idir,0)=ixMmin3-1+nghostcellsCo
             ixS_p_stg_min3(idir,1)=ixMmin3-1 ! -1 to make redundant
             ixS_p_stg_max3(idir,1)=ixMmin3-1+nxCo3+nghostcellsCo
             ixS_p_stg_min3(idir,2)=ixMmax3-nxCo3-nghostcellsCo
             ixS_p_stg_max3(idir,2)=ixMmax3
             ixS_p_stg_min3(idir,3)=ixMmax3-nghostcellsCo
             ixS_p_stg_max3(idir,3)=ixMmax3

             ixR_p_stg_min3(idir,0)=ixCoMmin3-1-nghostcellsCo
             ixR_p_stg_max3(idir,0)=ixCoMmin3-1
             ixR_p_stg_min3(idir,1)=ixCoMmin3-1 ! -1 to make redundant
             ixR_p_stg_max3(idir,1)=ixCoMmax3+nghostcellsCo
             ixR_p_stg_min3(idir,2)=ixCoMmin3-1-nghostcellsCo
             ixR_p_stg_max3(idir,2)=ixCoMmax3
             ixR_p_stg_min3(idir,3)=ixCoMmax3+1-1 ! -1 to make redundant
             ixR_p_stg_max3(idir,3)=ixCoMmax3+nghostcellsCo

          else

             ! Perpendicular component
             ixS_p_stg_min3(idir,0)=ixMmin3
             ixS_p_stg_max3(idir,0)=ixMmin3-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min3(idir,1)=ixMmin3
             ixS_p_stg_max3(idir,1)=ixMmin3-1+nxCo3+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min3(idir,2)=ixMmax3+&
                  1-nxCo3-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max3(idir,2)=ixMmax3
             ixS_p_stg_min3(idir,3)=ixMmax3+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max3(idir,3)=ixMmax3

             ixR_p_stg_min3(idir,0)=ixCoMmin3-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max3(idir,0)=ixCoMmin3-1
             ixR_p_stg_min3(idir,1)=ixCoMmin3
             ixR_p_stg_max3(idir,1)=ixCoMmax3+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min3(idir,2)=ixCoMmin3-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max3(idir,2)=ixCoMmax3
             ixR_p_stg_min3(idir,3)=ixCoMmax3+1
             ixR_p_stg_max3(idir,3)=ixCoMmax3+nghostcellsCo+&
                  (interpolation_order-1)


             ! Perpendicular component
             ixS_p_stg_min3(idir,0)=ixMmin3
             ixS_p_stg_max3(idir,0)=ixMmin3-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min3(idir,1)=ixMmin3
             ixS_p_stg_max3(idir,1)=ixMmin3-1+nxCo3+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min3(idir,2)=ixMmax3+&
                  1-nxCo3-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max3(idir,2)=ixMmax3
             ixS_p_stg_min3(idir,3)=ixMmax3+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max3(idir,3)=ixMmax3

             ixR_p_stg_min3(idir,0)=ixCoMmin3-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max3(idir,0)=ixCoMmin3-1
             ixR_p_stg_min3(idir,1)=ixCoMmin3
             ixR_p_stg_max3(idir,1)=ixCoMmax3+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min3(idir,2)=ixCoMmin3-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max3(idir,2)=ixCoMmax3
             ixR_p_stg_min3(idir,3)=ixCoMmax3+1
             ixR_p_stg_max3(idir,3)=ixCoMmax3+nghostcellsCo+&
                  (interpolation_order-1)


             ! Perpendicular component
             ixS_p_stg_min3(idir,0)=ixMmin3
             ixS_p_stg_max3(idir,0)=ixMmin3-1+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min3(idir,1)=ixMmin3
             ixS_p_stg_max3(idir,1)=ixMmin3-1+nxCo3+nghostcellsCo+&
                  (interpolation_order-1)
             ixS_p_stg_min3(idir,2)=ixMmax3+&
                  1-nxCo3-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max3(idir,2)=ixMmax3
             ixS_p_stg_min3(idir,3)=ixMmax3+&
                  1-nghostcellsCo-(interpolation_order-1)
             ixS_p_stg_max3(idir,3)=ixMmax3

             ixR_p_stg_min3(idir,0)=ixCoMmin3-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max3(idir,0)=ixCoMmin3-1
             ixR_p_stg_min3(idir,1)=ixCoMmin3
             ixR_p_stg_max3(idir,1)=ixCoMmax3+nghostcellsCo+&
                  (interpolation_order-1)
             ixR_p_stg_min3(idir,2)=ixCoMmin3-nghostcellsCo-(interpolation_order-&
                  1)
             ixR_p_stg_max3(idir,2)=ixCoMmax3
             ixR_p_stg_min3(idir,3)=ixCoMmax3+1
             ixR_p_stg_max3(idir,3)=ixCoMmax3+nghostcellsCo+&
                  (interpolation_order-1)

          end if

       end do
       ! calculate sizes for buffer arrays for siblings
       do i3=-1,1
          do i2=-1,1
             do i1=-1,1
                ! Staggered (face-allocated) variables
                do idir=1,ndim
                   sizes_srl_send_stg(idir,i1,i2,i3)=(ixS_srl_stg_max1(idir,&
                        i1)-ixS_srl_stg_min1(idir,i1)+1)*(ixS_srl_stg_max2(idir,&
                        i2)-ixS_srl_stg_min2(idir,i2)+1)*(ixS_srl_stg_max3(idir,&
                        i3)-ixS_srl_stg_min3(idir,i3)+1)
                   sizes_srl_recv_stg(idir,i1,i2,i3)=(ixR_srl_stg_max1(idir,&
                        i1)-ixR_srl_stg_min1(idir,i1)+1)*(ixR_srl_stg_max2(idir,&
                        i2)-ixR_srl_stg_min2(idir,i2)+1)*(ixR_srl_stg_max3(idir,&
                        i3)-ixR_srl_stg_min3(idir,i3)+1)
                   sizes_r_send_stg(idir,i1,i2,i3)=(ixS_r_stg_max1(idir,&
                        i1)-ixS_r_stg_min1(idir,i1)+1)*(ixS_r_stg_max2(idir,&
                        i2)-ixS_r_stg_min2(idir,i2)+1)*(ixS_r_stg_max3(idir,&
                        i3)-ixS_r_stg_min3(idir,i3)+1)
                end do
                sizes_srl_send_total(i1,i2,i3)=sum(sizes_srl_send_stg(:,i1,i2,i3))
                sizes_srl_recv_total(i1,i2,i3)=sum(sizes_srl_recv_stg(:,i1,i2,i3))
                sizes_r_send_total(i1,i2,i3)=sum(sizes_r_send_stg(:,i1,i2,i3))
             end do
          end do
       end do

       do i3=0,3
          do i2=0,3
             do i1=0,3
                ! Staggered (face-allocated) variables
                do idir=1,ndim
                   sizes_r_recv_stg(idir,i1,i2,i3)=(ixR_r_stg_max1(idir,&
                        i1)-ixR_r_stg_min1(idir,i1)+1)*(ixR_r_stg_max2(idir,&
                        i2)-ixR_r_stg_min2(idir,i2)+1)*(ixR_r_stg_max3(idir,&
                        i3)-ixR_r_stg_min3(idir,i3)+1)
                   sizes_p_send_stg(idir,i1,i2,i3)=(ixS_p_stg_max1(idir,&
                        i1)-ixS_p_stg_min1(idir,i1)+1)*(ixS_p_stg_max2(idir,&
                        i2)-ixS_p_stg_min2(idir,i2)+1)*(ixS_p_stg_max3(idir,&
                        i3)-ixS_p_stg_min3(idir,i3)+1)
                   sizes_p_recv_stg(idir,i1,i2,i3)=(ixR_p_stg_max1(idir,&
                        i1)-ixR_p_stg_min1(idir,i1)+1)*(ixR_p_stg_max2(idir,&
                        i2)-ixR_p_stg_min2(idir,i2)+1)*(ixR_p_stg_max3(idir,&
                        i3)-ixR_p_stg_min3(idir,i3)+1)
                end do
                sizes_r_recv_total(i1,i2,i3)=sum(sizes_r_recv_stg(:,i1,i2,i3))
                sizes_p_send_total(i1,i2,i3)=sum(sizes_p_send_stg(:,i1,i2,i3))
                sizes_p_recv_total(i1,i2,i3)=sum(sizes_p_recv_stg(:,i1,i2,i3))
             end do
          end do
       end do
    end if
    if(.not.stagger_grid .or. physics_type=='mf') then
       ! extend index range to physical boundary

       ixS_srl_min1(-1,0)=1
       ixS_srl_min1( 1,0)=ixMmin1
       ixS_srl_min1( 2,0)=1
       ixS_srl_max1(-1,0)=ixMmax1
       ixS_srl_max1( 1,0)=ixGmax1
       ixS_srl_max1( 2,0)=ixGmax1

       ixR_srl_min1(-1,0)=1
       ixR_srl_min1( 1,0)=ixMmin1
       ixR_srl_min1( 2,0)=1
       ixR_srl_max1(-1,0)=ixMmax1
       ixR_srl_max1( 1,0)=ixGmax1
       ixR_srl_max1( 2,0)=ixGmax1

       ixS_r_min1(-1,0)=1
       ixS_r_min1( 1,0)=ixCoMmin1
       ixS_r_max1(-1,0)=ixCoMmax1
       ixS_r_max1( 1,0)=ixCoGmax1

       ixR_r_min1(-1,1)=1
       ixR_r_max1(-1,1)=ixMmin1-1+nxCo1
       ixR_r_min1( 1,2)=ixMmin1+nxCo1
       ixR_r_max1( 1,2)=ixGmax1

       ixS_p_min1(-1,1)=1
       ixS_p_max1( 1,2)=ixGmax1

       ixR_p_min1(-1,1)=1
       ixR_p_max1( 1,2)=ixCoGmax1


       ixS_srl_min2(-1,0)=1
       ixS_srl_min2( 1,0)=ixMmin2
       ixS_srl_min2( 2,0)=1
       ixS_srl_max2(-1,0)=ixMmax2
       ixS_srl_max2( 1,0)=ixGmax2
       ixS_srl_max2( 2,0)=ixGmax2

       ixR_srl_min2(-1,0)=1
       ixR_srl_min2( 1,0)=ixMmin2
       ixR_srl_min2( 2,0)=1
       ixR_srl_max2(-1,0)=ixMmax2
       ixR_srl_max2( 1,0)=ixGmax2
       ixR_srl_max2( 2,0)=ixGmax2

       ixS_r_min2(-1,0)=1
       ixS_r_min2( 1,0)=ixCoMmin2
       ixS_r_max2(-1,0)=ixCoMmax2
       ixS_r_max2( 1,0)=ixCoGmax2

       ixR_r_min2(-1,1)=1
       ixR_r_max2(-1,1)=ixMmin2-1+nxCo2
       ixR_r_min2( 1,2)=ixMmin2+nxCo2
       ixR_r_max2( 1,2)=ixGmax2

       ixS_p_min2(-1,1)=1
       ixS_p_max2( 1,2)=ixGmax2

       ixR_p_min2(-1,1)=1
       ixR_p_max2( 1,2)=ixCoGmax2


       ixS_srl_min3(-1,0)=1
       ixS_srl_min3( 1,0)=ixMmin3
       ixS_srl_min3( 2,0)=1
       ixS_srl_max3(-1,0)=ixMmax3
       ixS_srl_max3( 1,0)=ixGmax3
       ixS_srl_max3( 2,0)=ixGmax3

       ixR_srl_min3(-1,0)=1
       ixR_srl_min3( 1,0)=ixMmin3
       ixR_srl_min3( 2,0)=1
       ixR_srl_max3(-1,0)=ixMmax3
       ixR_srl_max3( 1,0)=ixGmax3
       ixR_srl_max3( 2,0)=ixGmax3

       ixS_r_min3(-1,0)=1
       ixS_r_min3( 1,0)=ixCoMmin3
       ixS_r_max3(-1,0)=ixCoMmax3
       ixS_r_max3( 1,0)=ixCoGmax3

       ixR_r_min3(-1,1)=1
       ixR_r_max3(-1,1)=ixMmin3-1+nxCo3
       ixR_r_min3( 1,2)=ixMmin3+nxCo3
       ixR_r_max3( 1,2)=ixGmax3

       ixS_p_min3(-1,1)=1
       ixS_p_max3( 1,2)=ixGmax3

       ixR_p_min3(-1,1)=1
       ixR_p_max3( 1,2)=ixCoGmax3

    end if

    
    
    
    
    
    
    
    
    
    
    

  end subroutine init_bc