Skip to content

Commit

Permalink
Correcting tensor cross terms computation for periodic bcs (AMReX-Cod…
Browse files Browse the repository at this point in the history 
…es#1254)

## Summary
Correcting tensor cross terms computation for periodic bcs

## Additional background
Fixing a bug introduced in PR AMReX-Codes#1235. Now it computes in the correct way the cross terms in the tensor solve for problems with periodic bcs. 

## Checklist

The proposed changes:
- [x] fix a bug or incorrect behavior in AMReX
- [ ] add new capabilities to AMReX
- [x] changes answers in the test suite to more than roundoff level
- [x] are likely to significantly affect the results of downstream AMReX users
- [ ] are described in the proposed changes to the AMReX documentation, if appropriate
  • Loading branch information
@OscarAntepara
OscarAntepara authored Aug 11, 2020
1 parent a4ab677 commit 1197adc
Show file tree
Hide file tree
Showing 6 changed files with 67 additions and 54 deletions.
13 changes: 7 additions & 6 deletions Src/LinearSolvers/MLMG/AMReX_MLEBTensorOp.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -420,6 +420,7 @@ MLEBTensorOp::compCrossTerms(int amrlev, int mglev, MultiFab const& mf) const
const Box& domain = geom.Domain();
const auto dlo = amrex::lbound(domain);
const auto dhi = amrex::ubound(domain);
const GpuArray<int,AMREX_SPACEDIM>& is_periodic = geom.isPeriodicArray();
const auto& bcondloc = *m_bcondloc[amrlev][mglev];

Array<MultiFab,AMREX_SPACEDIM> const& etamf = m_b_coeffs[amrlev][mglev];
Expand Down Expand Up @@ -494,15 +495,15 @@ MLEBTensorOp::compCrossTerms(int amrlev, int mglev, MultiFab const& mf) const
AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM
( xbx, txbx,
{
mltensor_cross_terms_fx(txbx,fxfab,vfab,etaxfab,kapxfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fx(txbx,fxfab,vfab,etaxfab,kapxfab,dxinv,dlo,dhi,bct,is_periodic[0]);
}
, ybx, tybx,
{
mltensor_cross_terms_fy(tybx,fyfab,vfab,etayfab,kapyfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fy(tybx,fyfab,vfab,etayfab,kapyfab,dxinv,dlo,dhi,bct,is_periodic[1]);
}
, zbx, tzbx,
{
mltensor_cross_terms_fz(tzbx,fzfab,vfab,etazfab,kapzfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fz(tzbx,fzfab,vfab,etazfab,kapzfab,dxinv,dlo,dhi,bct,is_periodic[2]);
}
);
}
Expand All @@ -516,15 +517,15 @@ MLEBTensorOp::compCrossTerms(int amrlev, int mglev, MultiFab const& mf) const
AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM
( xbx, txbx,
{
mlebtensor_cross_terms_fx(txbx,fxfab,vfab,etaxfab,kapxfab,apx,flag,dxinv,dlo,dhi,bct);
mlebtensor_cross_terms_fx(txbx,fxfab,vfab,etaxfab,kapxfab,apx,flag,dxinv,dlo,dhi,bct,is_periodic[0]);
}
, ybx, tybx,
{
mlebtensor_cross_terms_fy(tybx,fyfab,vfab,etayfab,kapyfab,apy,flag,dxinv,dlo,dhi,bct);
mlebtensor_cross_terms_fy(tybx,fyfab,vfab,etayfab,kapyfab,apy,flag,dxinv,dlo,dhi,bct,is_periodic[1]);
}
, zbx, tzbx,
{
mlebtensor_cross_terms_fz(tzbx,fzfab,vfab,etazfab,kapzfab,apz,flag,dxinv,dlo,dhi,bct);
mlebtensor_cross_terms_fz(tzbx,fzfab,vfab,etazfab,kapzfab,apz,flag,dxinv,dlo,dhi,bct,is_periodic[2]);
}
);
}
Expand Down
14 changes: 8 additions & 6 deletions Src/LinearSolvers/MLMG/AMReX_MLEBTensor_2D_K.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,8 @@ void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
Array4<EBCellFlag const> const& flag,
GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
const Dim3& dlo, const Dim3& dhi,
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct) noexcept
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct,
const int& is_periodic) noexcept
{
const Real dyi = dxinv[1];
const auto lo = amrex::lbound(box);
Expand Down Expand Up @@ -52,7 +53,7 @@ void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
Real dvdy = (0.5*dyi) * ((vel(i ,jhip,0,1)-vel(i ,jhim,0,1))*whi
+(vel(i-1,jlop,0,1)-vel(i-1,jlom,0,1))*wlo);

if (i==dhi.x+1){
if (i==dhi.x+1 and !is_periodic){
if (bct(Orientation(Direction::x,Orientation::high),0)==AMREX_LO_NEUMANN){
dudy = (dyi) * ((vel(i-1,jlop,0,0)-vel(i-1,jlom,0,0))*wlo);
}
Expand All @@ -66,7 +67,7 @@ void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
dvdy = (dyi) * ((vel(i ,jhip,0,1)-vel(i ,jhim,0,1))*whi);
}
}
if (i==dlo.x){
if (i==dlo.x and !is_periodic){
if (bct(Orientation(Direction::x,Orientation::low),0)==AMREX_LO_NEUMANN){
dudy = (dyi) * ((vel(i ,jhip,0,0)-vel(i ,jhim,0,0))*whi);
}
Expand Down Expand Up @@ -101,7 +102,8 @@ void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
Array4<EBCellFlag const> const& flag,
GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
const Dim3& dlo, const Dim3& dhi,
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct) noexcept
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct,
const int& is_periodic) noexcept
{
const Real dxi = dxinv[0];
const auto lo = amrex::lbound(box);
Expand Down Expand Up @@ -131,7 +133,7 @@ void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
Real dvdx = (0.5*dxi) * ((vel(ihip,j ,0,1)-vel(ihim,j ,0,1))*whi
+(vel(ilop,j-1,0,1)-vel(ilom,j-1,0,1))*wlo);

if (j==dhi.y+1){
if (j==dhi.y+1 and !is_periodic){
if (bct(Orientation(Direction::y,Orientation::high),0)==AMREX_LO_NEUMANN){
dudx = (dxi) * ((vel(ilop,j-1,0,0)-vel(ilom,j-1,0,0))*wlo);
}
Expand All @@ -145,7 +147,7 @@ void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
dvdx = (dxi) * ((vel(ihip,j ,0,1)-vel(ihim,j ,0,1))*whi);
}
}
if (j==dlo.y){
if (j==dlo.y and !is_periodic){
if (bct(Orientation(Direction::y,Orientation::low),0)==AMREX_LO_NEUMANN){
dudx = (dxi) * ((vel(ihip,j ,0,0)-vel(ihim,j ,0,0))*whi);
}
Expand Down
33 changes: 18 additions & 15 deletions Src/LinearSolvers/MLMG/AMReX_MLEBTensor_3D_K.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,8 @@ void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
Array4<EBCellFlag const> const& flag,
GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
const Dim3& dlo, const Dim3& dhi,
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct) noexcept
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct,
const int& is_periodic) noexcept
{
const Real dyi = dxinv[1];
const Real dzi = dxinv[2];
Expand Down Expand Up @@ -55,7 +56,7 @@ void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
Real dvdy = (0.5*dyi) * ((vel(i ,jhip,k,1)-vel(i ,jhim,k,1))*whi
+(vel(i-1,jlop,k,1)-vel(i-1,jlom,k,1))*wlo);

if (i==dhi.x+1){
if (i==dhi.x+1 and !is_periodic){
if (bct(Orientation(Direction::x,Orientation::high),0)==AMREX_LO_NEUMANN){
dudy = (dyi) * ((vel(i-1,jlop,k,0)-vel(i-1,jlom,k,0))*wlo);
}
Expand All @@ -69,7 +70,7 @@ void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
dvdy = (dyi) * ((vel(i ,jhip,k,1)-vel(i ,jhim,k,1))*whi);
}
}
if (i==dlo.x){
if (i==dlo.x and !is_periodic){
if (bct(Orientation(Direction::x,Orientation::low),0)==AMREX_LO_NEUMANN){
dudy = (dyi) * ((vel(i ,jhip,k,0)-vel(i ,jhim,k,0))*whi);
}
Expand Down Expand Up @@ -97,7 +98,7 @@ void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
Real dwdz = (0.5*dzi) * ((vel(i ,j,khip,2)-vel(i ,j,khim,2))*whi
+(vel(i-1,j,klop,2)-vel(i-1,j,klom,2))*wlo);

if (i==dhi.x+1){
if (i==dhi.x+1 and !is_periodic){
if (bct(Orientation(Direction::x,Orientation::high),0)==AMREX_LO_NEUMANN){
dudz = (dzi) * ((vel(i-1,j,klop,0)-vel(i-1,j,klom,0))*wlo);
}
Expand All @@ -111,7 +112,7 @@ void mlebtensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
dwdz = (dzi) * ((vel(i ,j,khip,2)-vel(i ,j,khim,2))*whi);
}
}
if (i==dlo.x){
if (i==dlo.x and !is_periodic){
if (bct(Orientation(Direction::x,Orientation::low),0)==AMREX_LO_NEUMANN){
dudz = (dzi) * ((vel(i ,j,khip,0)-vel(i ,j,khim,0))*whi);
}
Expand Down Expand Up @@ -148,7 +149,8 @@ void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
Array4<EBCellFlag const> const& flag,
GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
const Dim3& dlo, const Dim3& dhi,
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct) noexcept
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct,
const int& is_periodic) noexcept
{
const Real dxi = dxinv[0];
const Real dzi = dxinv[2];
Expand Down Expand Up @@ -181,7 +183,7 @@ void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
Real dvdx = (0.5*dxi) * ((vel(ihip,j ,k,1)-vel(ihim,j ,k,1))*whi
+(vel(ilop,j-1,k,1)-vel(ilom,j-1,k,1))*wlo);

if (j==dhi.y+1){
if (j==dhi.y+1 and !is_periodic){
if (bct(Orientation(Direction::y,Orientation::high),0)==AMREX_LO_NEUMANN){
dudx = (dxi) * ((vel(ilop,j-1,k,0)-vel(ilom,j-1,k,0))*wlo);
}
Expand All @@ -195,7 +197,7 @@ void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
dvdx = (dxi) * ((vel(ihip,j ,k,1)-vel(ihim,j ,k,1))*whi);
}
}
if (j==dlo.y){
if (j==dlo.y and !is_periodic){
if (bct(Orientation(Direction::y,Orientation::low),0)==AMREX_LO_NEUMANN){
dudx = (dxi) * ((vel(ihip,j ,k,0)-vel(ihim,j ,k,0))*whi);
}
Expand Down Expand Up @@ -223,7 +225,7 @@ void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
Real dwdz = (0.5*dzi) * ((vel(i,j ,khip,2)-vel(i,j ,khim,2))*whi
+(vel(i,j-1,klop,2)-vel(i,j-1,klom,2))*wlo);

if (j==dhi.y+1){
if (j==dhi.y+1 and !is_periodic){
if (bct(Orientation(Direction::y,Orientation::high),1)==AMREX_LO_NEUMANN){
dvdz = (dzi) * ((vel(i,j-1,klop,1)-vel(i,j-1,klom,1))*wlo);
}
Expand All @@ -237,7 +239,7 @@ void mlebtensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
dwdz = (dzi) * ((vel(i,j ,khip,2)-vel(i,j ,khim,2))*whi);
}
}
if (j==dlo.y){
if (j==dlo.y and !is_periodic){
if (bct(Orientation(Direction::y,Orientation::low),1)==AMREX_LO_NEUMANN){
dvdz = (dzi) * ((vel(i,j ,khip,1)-vel(i,j ,khim,1))*whi);
}
Expand Down Expand Up @@ -274,7 +276,8 @@ void mlebtensor_cross_terms_fz (Box const& box, Array4<Real> const& fz,
Array4<EBCellFlag const> const& flag,
GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
const Dim3& dlo, const Dim3& dhi,
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct) noexcept
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct,
const int& is_periodic) noexcept
{
const Real dxi = dxinv[0];
const Real dyi = dxinv[1];
Expand Down Expand Up @@ -307,7 +310,7 @@ void mlebtensor_cross_terms_fz (Box const& box, Array4<Real> const& fz,
Real dwdx = (0.5*dxi) * ((vel(ihip,j,k ,2)-vel(ihim,j,k ,2))*whi
+(vel(ilop,j,k-1,2)-vel(ilom,j,k-1,2))*wlo);

if (k==dhi.z+1){
if (k==dhi.z+1 and !is_periodic){
if (bct(Orientation(Direction::z,Orientation::high),0)==AMREX_LO_NEUMANN){
dudx = (dxi) * ((vel(ilop,j,k-1,0)-vel(ilom,j,k-1,0))*wlo);
}
Expand All @@ -321,7 +324,7 @@ void mlebtensor_cross_terms_fz (Box const& box, Array4<Real> const& fz,
dwdx = (dxi) * ((vel(ihip,j,k ,2)-vel(ihim,j,k ,2))*whi);
}
}
if (k==dlo.z){
if (k==dlo.z and !is_periodic){
if (bct(Orientation(Direction::z,Orientation::low),0)==AMREX_LO_NEUMANN){
dudx = (dxi) * ((vel(ihip,j,k ,0)-vel(ihim,j,k ,0))*whi);
}
Expand Down Expand Up @@ -349,7 +352,7 @@ void mlebtensor_cross_terms_fz (Box const& box, Array4<Real> const& fz,
Real dwdy = (0.5*dyi) * ((vel(i,jhip,k ,2)-vel(i,jhim,k ,2))*whi
+(vel(i,jlop,k-1,2)-vel(i,jlom,k-1,2))*wlo);

if (k==dhi.z+1){
if (k==dhi.z+1 and !is_periodic){
if (bct(Orientation(Direction::z,Orientation::high),1)==AMREX_LO_NEUMANN){
dvdy = (dyi) * ((vel(i,jlop,k-1,1)-vel(i,jlom,k-1,1))*wlo);
}
Expand All @@ -363,7 +366,7 @@ void mlebtensor_cross_terms_fz (Box const& box, Array4<Real> const& fz,
dwdy = (dyi) * ((vel(i,jhip,k ,2)-vel(i,jhim,k ,2))*whi);
}
}
if (k==dlo.z){
if (k==dlo.z and !is_periodic){
if (bct(Orientation(Direction::z,Orientation::low),1)==AMREX_LO_NEUMANN){
dvdy = (dyi) * ((vel(i,jhip,k ,1)-vel(i,jhim,k ,1))*whi);
}
Expand Down
14 changes: 8 additions & 6 deletions Src/LinearSolvers/MLMG/AMReX_MLTensorOp.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -215,6 +215,7 @@ MLTensorOp::apply (int amrlev, int mglev, MultiFab& out, MultiFab& in, BCMode bc
const Box& domain = geom.Domain();
const auto dlo = amrex::lbound(domain);
const auto dhi = amrex::ubound(domain);
const GpuArray<int,AMREX_SPACEDIM>& is_periodic = geom.isPeriodicArray();
const auto& bcondloc = *m_bcondloc[amrlev][mglev];

Array<MultiFab,AMREX_SPACEDIM> const& etamf = m_b_coeffs[amrlev][mglev];
Expand Down Expand Up @@ -262,15 +263,15 @@ MLTensorOp::apply (int amrlev, int mglev, MultiFab& out, MultiFab& in, BCMode bc
AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM
( xbx, txbx,
{
mltensor_cross_terms_fx(txbx,fxfab,vfab,etaxfab,kapxfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fx(txbx,fxfab,vfab,etaxfab,kapxfab,dxinv,dlo,dhi,bct,is_periodic[0]);
}
, ybx, tybx,
{
mltensor_cross_terms_fy(tybx,fyfab,vfab,etayfab,kapyfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fy(tybx,fyfab,vfab,etayfab,kapyfab,dxinv,dlo,dhi,bct,is_periodic[1]);
}
, zbx, tzbx,
{
mltensor_cross_terms_fz(tzbx,fzfab,vfab,etazfab,kapzfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fz(tzbx,fzfab,vfab,etazfab,kapzfab,dxinv,dlo,dhi,bct,is_periodic[2]);
}
);

Expand Down Expand Up @@ -421,6 +422,7 @@ MLTensorOp::compFlux (int amrlev, const Array<MultiFab*,AMREX_SPACEDIM>& fluxes,
const Box& domain = geom.Domain();
const auto dlo = amrex::lbound(domain);
const auto dhi = amrex::ubound(domain);
const GpuArray<int,AMREX_SPACEDIM>& is_periodic = geom.isPeriodicArray();
const auto& bcondloc = *m_bcondloc[amrlev][mglev];

Array<MultiFab,AMREX_SPACEDIM> const& etamf = m_b_coeffs[amrlev][mglev];
Expand Down Expand Up @@ -467,15 +469,15 @@ MLTensorOp::compFlux (int amrlev, const Array<MultiFab*,AMREX_SPACEDIM>& fluxes,
AMREX_LAUNCH_HOST_DEVICE_LAMBDA_DIM
( xbx, txbx,
{
mltensor_cross_terms_fx(txbx,fxfab,vfab,etaxfab,kapxfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fx(txbx,fxfab,vfab,etaxfab,kapxfab,dxinv,dlo,dhi,bct,is_periodic[0]);
}
, ybx, tybx,
{
mltensor_cross_terms_fy(tybx,fyfab,vfab,etayfab,kapyfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fy(tybx,fyfab,vfab,etayfab,kapyfab,dxinv,dlo,dhi,bct,is_periodic[1]);
}
, zbx, tzbx,
{
mltensor_cross_terms_fz(tzbx,fzfab,vfab,etazfab,kapzfab,dxinv,dlo,dhi,bct);
mltensor_cross_terms_fz(tzbx,fzfab,vfab,etazfab,kapzfab,dxinv,dlo,dhi,bct,is_periodic[2]);
}
);

Expand Down
14 changes: 8 additions & 6 deletions Src/LinearSolvers/MLMG/AMReX_MLTensor_2D_K.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -131,7 +131,8 @@ void mltensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
Array4<Real const> const& kapx,
GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
const Dim3& dlo, const Dim3& dhi,
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct) noexcept
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct,
const int& is_periodic) noexcept
{
const Real dyi = dxinv[1];
const auto lo = amrex::lbound(box);
Expand All @@ -143,7 +144,7 @@ void mltensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
for (int i = lo.x; i <= hi.x; ++i) {
Real dudy = (vel(i,j+1,0,0)+vel(i-1,j+1,0,0)-vel(i,j-1,0,0)-vel(i-1,j-1,0,0))*(0.25*dyi);
Real dvdy = (vel(i,j+1,0,1)+vel(i-1,j+1,0,1)-vel(i,j-1,0,1)-vel(i-1,j-1,0,1))*(0.25*dyi);
if (i==dhi.x+1){
if (i==dhi.x+1 and !is_periodic){
if (bct(Orientation(Direction::x,Orientation::high),0)==AMREX_LO_NEUMANN){
dudy = (0.5*dyi) * (vel(i-1,j+1,0,0)-vel(i-1,j-1,0,0));
}
Expand All @@ -157,7 +158,7 @@ void mltensor_cross_terms_fx (Box const& box, Array4<Real> const& fx,
dvdy = (0.5*dyi) * (vel(i ,j+1,0,1)-vel(i ,j-1,0,1));
}
}
if (i==dlo.x){
if (i==dlo.x and !is_periodic){
if (bct(Orientation(Direction::x,Orientation::low),0)==AMREX_LO_NEUMANN){
dudy = (0.5*dyi) * (vel(i ,j+1,0,0)-vel(i ,j-1,0,0));
}
Expand Down Expand Up @@ -188,7 +189,8 @@ void mltensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
Array4<Real const> const& kapy,
GpuArray<Real,AMREX_SPACEDIM> const& dxinv,
const Dim3& dlo, const Dim3& dhi,
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct) noexcept
Array2D<BoundCond,0,2*AMREX_SPACEDIM-1,0,AMREX_SPACEDIM-1> const& bct,
const int& is_periodic) noexcept
{
const Real dxi = dxinv[0];
const auto lo = amrex::lbound(box);
Expand All @@ -200,7 +202,7 @@ void mltensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
for (int i = lo.x; i <= hi.x; ++i) {
Real dudx = (vel(i+1,j,0,0)+vel(i+1,j-1,0,0)-vel(i-1,j,0,0)-vel(i-1,j-1,0,0))*(0.25*dxi);
Real dvdx = (vel(i+1,j,0,1)+vel(i+1,j-1,0,1)-vel(i-1,j,0,1)-vel(i-1,j-1,0,1))*(0.25*dxi);
if (j==dhi.y+1){
if (j==dhi.y+1 and !is_periodic){
if (bct(Orientation(Direction::y,Orientation::high),0)==AMREX_LO_NEUMANN){
dudx = (0.5*dxi) * (vel(i+1,j-1,0,0)-vel(i-1,j-1,0,0));
}
Expand All @@ -214,7 +216,7 @@ void mltensor_cross_terms_fy (Box const& box, Array4<Real> const& fy,
dvdx = (0.5*dxi) * (vel(i+1,j ,0,1)-vel(i-1,j ,0,1));
}
}
if (j==dlo.y){
if (j==dlo.y and !is_periodic){
if (bct(Orientation(Direction::y,Orientation::low),0)==AMREX_LO_NEUMANN){
dudx = (0.5*dxi) * (vel(i+1,j ,0,0)-vel(i-1,j ,0,0));
}
Expand Down
Loading

0 comments on commit 1197adc

Please sign in to comment.