1 files changed, 49 insertions, 34 deletions

diff --git a/apps/adrian/cylinder2d/cylinder2d.cpp b/apps/adrian/cylinder2d/cylinder2d.cpp
index ffab6a0..bf25417 100644
--- a/apps/adrian/cylinder2d/cylinder2d.cpp
+++ b/apps/adrian/cylinder2d/cylinder2d.cpp
@@ -39,11 +39,11 @@ typedef double T;
 /// Setup geometry relative to cylinder diameter as defined by [SchaeferTurek96]
 const T cylinderD = 1.0;
 
-const int N = 20;              // resolution of the model
+const int N = 5;               // resolution of the model
 const T lx  = 22  * cylinderD; // length of the channel
 const T ly  = 4.1 * cylinderD; // height of the channel
 const T Re  = 100.;            // Reynolds number
-const T uLat = 0.1;            // lattice velocity
+const T uLat = 0.05;           // lattice velocity
 const T maxPhysT = 60.;        // max. simulation time in s, SI unit
 
 
@@ -67,13 +67,13 @@ void prepareGeometry(Grid2D<T,DESCRIPTOR>& grid)
     IndicatorCuboid2D<T> lowerWall(wallExtend, wallOrigin);
     sGeometry.rename(1,2,lowerWall);
 
-    IndicatorCuboid2D<T> upperWall(wallExtend, wallOrigin + Vector<T,2> {0,ly});
+    IndicatorCuboid2D<T> upperWall(wallExtend, wallOrigin + Vector<T,2> {0,ly-physSpacing/2});
     sGeometry.rename(1,2,upperWall);
   }
 
   // Set material number for inflow and outflow
   {
-    const Vector<T,2> extend { physSpacing/2,  ly};
+    const Vector<T,2> extend { physSpacing/2,  ly-physSpacing/2};
     const Vector<T,2> origin {-physSpacing/4, -physSpacing/4};
 
     IndicatorCuboid2D<T> inflow(extend, origin);
@@ -216,27 +216,58 @@ int main(int argc, char* argv[])
     coarseGrid.getConverter().getLatticeRelaxationFrequency(),
     instances::getBulkMomenta<T, DESCRIPTOR>());
 
-  sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> coarseSBoundaryCondition(coarseGrid.getSuperLattice());
-  createLocalBoundaryCondition2D<T, DESCRIPTOR>(coarseSBoundaryCondition);
+  sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> coarseBoundaryCondition(coarseGrid.getSuperLattice());
+  createLocalBoundaryCondition2D<T, DESCRIPTOR>(coarseBoundaryCondition);
 
-  const Vector<T,2> fineExtend {16*cylinderD, 3.75*cylinderD};
-  const Vector<T,2> fineOrigin {0.75*cylinderD, (ly-fineExtend[1])/2};
+  const auto coarseDeltaX = coarseGrid.getConverter().getPhysDeltaX();
+
+  const Vector<T,2> fineExtend {8*cylinderD, ly-2*coarseDeltaX};
+  const Vector<T,2> fineOrigin {0.75*cylinderD, coarseDeltaX};
 
   auto& fineGrid = coarseGrid.refine(fineOrigin, fineExtend);
   prepareGeometry(fineGrid);
 
+  const Vector<T,2> fineExtendB {10*coarseDeltaX, ly};
+  const Vector<T,2> fineOriginB {lx-10*coarseDeltaX, 0};
+
+  auto& fineOutflowGrid = coarseGrid.refine(fineOriginB, fineExtendB, false);
+  prepareGeometry(fineOutflowGrid);
+
+  {
+    const Vector<T,2> origin = fineOutflowGrid.getOrigin();
+    const Vector<T,2> extend = fineOutflowGrid.getExtend();
+
+    const Vector<T,2> extendY = {0,extend[1]};
+
+    coarseGrid.addFineCoupling(fineOutflowGrid, origin, extendY);
+    coarseGrid.addCoarseCoupling(fineOutflowGrid, origin + Vector<T,2> {coarseDeltaX,0}, extendY);
+
+    IndicatorCuboid2D<T> refined(extend, origin + Vector<T,2> {2*coarseDeltaX,0});
+
+    coarseGrid.getSuperGeometry().rename(1,0,refined);
+    coarseGrid.getSuperGeometry().rename(2,0,refined);
+    coarseGrid.getSuperGeometry().rename(4,0,refined);
+  }
+
   disableRefinedArea(coarseGrid, fineGrid);
 
   BGKdynamics<T, DESCRIPTOR> fineBulkDynamics(
     fineGrid.getConverter().getLatticeRelaxationFrequency(),
     instances::getBulkMomenta<T, DESCRIPTOR>());
 
-  sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> fineSBoundaryCondition(fineGrid.getSuperLattice());
-  createLocalBoundaryCondition2D<T, DESCRIPTOR>(fineSBoundaryCondition);
+  sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> fineBoundaryCondition(fineGrid.getSuperLattice());
+  createLocalBoundaryCondition2D<T, DESCRIPTOR>(fineBoundaryCondition);
 
-  prepareLattice(coarseGrid, coarseBulkDynamics, coarseSBoundaryCondition);
+  BGKdynamics<T, DESCRIPTOR> fineOutflowBulkDynamics(
+    fineOutflowGrid.getConverter().getLatticeRelaxationFrequency(),
+    instances::getBulkMomenta<T, DESCRIPTOR>());
+
+  sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> fineOutflowBoundaryCondition(fineOutflowGrid.getSuperLattice());
+  createLocalBoundaryCondition2D<T, DESCRIPTOR>(fineOutflowBoundaryCondition);
+
+  prepareLattice(coarseGrid, coarseBulkDynamics, coarseBoundaryCondition);
 
-  const Vector<T,2> fineExtend2 {11*cylinderD, 3*cylinderD};
+  const Vector<T,2> fineExtend2 {4*cylinderD, 2*cylinderD};
   const Vector<T,2> fineOrigin2 {1*cylinderD, 2*cylinderD-fineExtend2[1]/2};
 
   auto& fineGrid2 = fineGrid.refine(fineOrigin2, fineExtend2);
@@ -248,29 +279,13 @@ int main(int argc, char* argv[])
     fineGrid2.getConverter().getLatticeRelaxationFrequency(),
     instances::getBulkMomenta<T, DESCRIPTOR>());
 
-  sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> fineSBoundaryCondition2(fineGrid2.getSuperLattice());
-  createLocalBoundaryCondition2D<T, DESCRIPTOR>(fineSBoundaryCondition2);
-
-  prepareLattice(fineGrid, fineBulkDynamics, fineSBoundaryCondition);
-
-  const Vector<T,2> fineExtend3 {1.5*cylinderD, 1.5*cylinderD};
-  const Vector<T,2> fineOrigin3 {2*cylinderD-fineExtend3[0]/2, 2*cylinderD-fineExtend3[1]/2};
-
-  auto& fineGrid3 = fineGrid2.refine(fineOrigin3, fineExtend3);
-  prepareGeometry(fineGrid3);
-
-  disableRefinedArea(fineGrid2, fineGrid3);
-
-  BGKdynamics<T, DESCRIPTOR> fineBulkDynamics3(
-    fineGrid3.getConverter().getLatticeRelaxationFrequency(),
-    instances::getBulkMomenta<T, DESCRIPTOR>());
-
-  sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> fineSBoundaryCondition3(fineGrid3.getSuperLattice());
-  createLocalBoundaryCondition2D<T, DESCRIPTOR>(fineSBoundaryCondition3);
+  sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> fineBoundaryCondition2(fineGrid2.getSuperLattice());
+  createLocalBoundaryCondition2D<T, DESCRIPTOR>(fineBoundaryCondition2);
 
-  prepareLattice(fineGrid2, fineBulkDynamics2, fineSBoundaryCondition2);
+  prepareLattice(fineGrid, fineBulkDynamics, fineBoundaryCondition);
+  prepareLattice(fineOutflowGrid, fineOutflowBulkDynamics, fineOutflowBoundaryCondition);
 
-  prepareLattice(fineGrid3, fineBulkDynamics3, fineSBoundaryCondition3);
+  prepareLattice(fineGrid2, fineBulkDynamics2, fineBoundaryCondition2);
 
   clout << "Starting simulation..." << endl;
   Timer<T> timer(
@@ -286,8 +301,8 @@ int main(int argc, char* argv[])
 
     getResults("level0_", coarseGrid, iT);
     getResults("level1_", fineGrid,   iT);
+    getResults("level1_outflow_", fineOutflowGrid, iT);
     getResults("level2_", fineGrid2,  iT);
-    getResults("level3_", fineGrid3,  iT);
 
     if (iT%statIter == 0) {
       timer.update(iT);