Test refinement around obstacle

Coupler construction is to be wrapped by `Grid::refine` for refinement
areas that are fully inside the coarse domain.

1 files changed, 34 insertions, 16 deletions

diff --git a/apps/adrian/poiseuille2d/poiseuille2d.cpp b/apps/adrian/poiseuille2d/poiseuille2d.cpp
index 6c50339..f160135 100644
--- a/apps/adrian/poiseuille2d/poiseuille2d.cpp
+++ b/apps/adrian/poiseuille2d/poiseuille2d.cpp
@@ -38,7 +38,7 @@ typedef double T;
 const T lx  = 4.;             // length of the channel
 const T ly  = 1.;             // height of the channel
 const int N = 30;             // resolution of the model
-const T Re  = 100.;            // Reynolds number
+const T Re  = 100.;           // Reynolds number
 const T maxPhysT = 60.;       // max. simulation time in s, SI unit
 const T physInterval = 0.25;  // interval for the convergence check in s
 const T residuum = 1e-5;      // residuum for the convergence check
@@ -68,8 +68,8 @@ void prepareGeometry(UnitConverter<T,DESCRIPTOR> const& converter,
   IndicatorCuboid2D<T> outflow(extend, origin);
   superGeometry.rename(1,4,1,outflow);
 
-  //IndicatorCuboid2D<T> obstacle(Vector<T,2>{0.5,0.2}, Vector<T,2>{1.75,0.4});
-  //superGeometry.rename(1,2,obstacle);
+  IndicatorCuboid2D<T> obstacle(Vector<T,2>{0.2,0.1}, Vector<T,2>{1.25,0.45});
+  superGeometry.rename(1,5,obstacle);
 
   // Removes all not needed boundary voxels outside the surface
   superGeometry.clean();
@@ -365,8 +365,11 @@ public:
   {
     OLB_ASSERT(util::nearZero(extend[0]) || util::nearZero(extend[1]), "Coupling is only implemented alongside unit vectors");
 
-    _coarse.getCuboidGeometry().getFloorLatticeR(_origin, _coarseOrigin);
-    _fine.getCuboidGeometry().getFloorLatticeR(_origin, _fineOrigin);
+    std::vector<int> tmp(3,0);
+    _coarse.getCuboidGeometry().getLatticeR(_origin.toStdVector(), tmp);
+    _coarseOrigin = tmp;
+    _fine.getCuboidGeometry().getLatticeR(_origin.toStdVector(), tmp);
+    _fineOrigin = tmp;
   }
 };
 
@@ -611,8 +614,8 @@ int main(int argc, char* argv[])
   Vector<T,2> coarseExtend { lx, ly  };
   IndicatorCuboid2D<T> coarseCuboid(coarseExtend, coarseOrigin);
 
-  Vector<T,2> fineOrigin   { 0.25 * lx, 0.0 };
-  Vector<T,2> fineExtend   { 0.5  * lx, ly };
+  Vector<T,2> fineOrigin   { 0.25 * lx, 0.3 };
+  Vector<T,2> fineExtend   { 0.5  * lx, 0.4 };
   IndicatorCuboid2D<T> fineCuboid(fineExtend, fineOrigin);
 
   auto coarseGrid = Grid<T,DESCRIPTOR>::make(coarseCuboid, N, 0.8, Re);
@@ -621,9 +624,12 @@ int main(int argc, char* argv[])
   prepareGeometry(coarseGrid->getConverter(), coarseGrid->getSuperGeometry());
   prepareGeometry(fineGrid->getConverter(), fineGrid->getSuperGeometry());
 
-  IndicatorCuboid2D<T> overlapCuboid(Vector<T,2>{fineExtend[0]-2*coarseDeltaX, fineExtend[1]}, Vector<T,2>{fineOrigin[0]+coarseDeltaX, fineOrigin[1]});
+  fineGrid->getSuperGeometry().rename(2,1);
+  fineGrid->getSuperGeometry().rename(5,2);
+  IndicatorCuboid2D<T> overlapCuboid(fineExtend - 4*coarseDeltaX, fineOrigin + 2*coarseDeltaX);
   coarseGrid->getSuperGeometry().rename(1,0,overlapCuboid);
   coarseGrid->getSuperGeometry().rename(2,0,overlapCuboid);
+  coarseGrid->getSuperGeometry().rename(5,0,overlapCuboid);
 
   Dynamics<T, DESCRIPTOR>* coarseBulkDynamics;
   coarseBulkDynamics = new BGKdynamics<T, DESCRIPTOR>(
@@ -664,15 +670,15 @@ int main(int argc, char* argv[])
     residuum);
   timer.start();
 
-  const Vector<T,2> originC2F = fineOrigin;
-  const Vector<T,2> extendC2F { 0, fineExtend[1] };
-  FineCoupler<T,DESCRIPTOR> fineCoupler0(*coarseGrid, *fineGrid, fineOrigin, extendC2F);
-  FineCoupler<T,DESCRIPTOR> fineCoupler1(*coarseGrid, *fineGrid, fineOrigin + Vector<T,2>(fineExtend[0],0), extendC2F);
+  FineCoupler<T,DESCRIPTOR> fineCoupler0(*coarseGrid, *fineGrid, fineOrigin, Vector<T,2>{0,fineExtend[1]+coarseDeltaX});
+  FineCoupler<T,DESCRIPTOR> fineCoupler1(*coarseGrid, *fineGrid, fineOrigin + Vector<T,2>{fineExtend[0],0}, Vector<T,2>{0,fineExtend[1]+coarseDeltaX});
+  FineCoupler<T,DESCRIPTOR> fineCoupler2(*coarseGrid, *fineGrid, fineOrigin + Vector<T,2>(0,fineExtend[1]), Vector<T,2>{fineExtend[0]+coarseDeltaX,0});
+  FineCoupler<T,DESCRIPTOR> fineCoupler3(*coarseGrid, *fineGrid, fineOrigin, Vector<T,2>{fineExtend[0]+coarseDeltaX,0});
 
-  const Vector<T,2> originF2C { fineOrigin[0] + coarseDeltaX, fineOrigin[1] };
-  const Vector<T,2> extendF2C { 0, fineExtend[1] };
-  CoarseCoupler<T,DESCRIPTOR> coarseCoupler0(*coarseGrid, *fineGrid, originF2C, extendF2C);
-  CoarseCoupler<T,DESCRIPTOR> coarseCoupler1(*coarseGrid, *fineGrid, fineOrigin + Vector<T,2>(fineExtend[0]-coarseDeltaX,0), extendF2C);
+  CoarseCoupler<T,DESCRIPTOR> coarseCoupler0(*coarseGrid, *fineGrid, fineOrigin + coarseDeltaX, Vector<T,2>{0,fineExtend[1]-2*coarseDeltaX});
+  CoarseCoupler<T,DESCRIPTOR> coarseCoupler1(*coarseGrid, *fineGrid, fineOrigin + Vector<T,2>(fineExtend[0]-coarseDeltaX,coarseDeltaX), Vector<T,2>{0,fineExtend[1]-2*coarseDeltaX});
+  CoarseCoupler<T,DESCRIPTOR> coarseCoupler2(*coarseGrid, *fineGrid, fineOrigin + Vector<T,2>{coarseDeltaX,fineExtend[1]-coarseDeltaX}, Vector<T,2>{fineExtend[0]-2*coarseDeltaX,0});
+  CoarseCoupler<T,DESCRIPTOR> coarseCoupler3(*coarseGrid, *fineGrid, fineOrigin + coarseDeltaX, Vector<T,2>{fineExtend[0]-2*coarseDeltaX,0});
 
   for (int iT = 0; iT < coarseGrid->getConverter().getLatticeTime(maxPhysT); ++iT) {
     if (converge.hasConverged()) {
@@ -682,6 +688,8 @@ int main(int argc, char* argv[])
 
     fineCoupler0.store();
     fineCoupler1.store();
+    fineCoupler2.store();
+    fineCoupler3.store();
     coarseGrid->getSuperLattice().collideAndStream();
 
     fineGrid->getSuperLattice().collideAndStream();
@@ -689,15 +697,25 @@ int main(int argc, char* argv[])
     fineCoupler0.couple();
     fineCoupler1.interpolate();
     fineCoupler1.couple();
+    fineCoupler2.interpolate();
+    fineCoupler2.couple();
+    fineCoupler3.interpolate();
+    fineCoupler3.couple();
 
     fineGrid->getSuperLattice().collideAndStream();
     fineCoupler0.store();
     fineCoupler0.couple();
     fineCoupler1.store();
     fineCoupler1.couple();
+    fineCoupler2.store();
+    fineCoupler2.couple();
+    fineCoupler3.store();
+    fineCoupler3.couple();
 
     coarseCoupler0.couple();
     coarseCoupler1.couple();
+    coarseCoupler2.couple();
+    coarseCoupler3.couple();
 
     getResults(
       "coarse_",