diff options
-rw-r--r-- | apps/adrian/cylinder2d/optimized_grid/cylinder2d.cpp | 47 | ||||
-rw-r--r-- | src/refinement/coupler2D.hh | 12 | ||||
-rw-r--r-- | src/refinement/grid2D.h | 6 | ||||
-rw-r--r-- | src/refinement/grid2D.hh | 12 |
4 files changed, 52 insertions, 25 deletions
diff --git a/apps/adrian/cylinder2d/optimized_grid/cylinder2d.cpp b/apps/adrian/cylinder2d/optimized_grid/cylinder2d.cpp index e5a934c..d516394 100644 --- a/apps/adrian/cylinder2d/optimized_grid/cylinder2d.cpp +++ b/apps/adrian/cylinder2d/optimized_grid/cylinder2d.cpp @@ -137,7 +137,6 @@ void prepareLattice(Grid2D<T,DESCRIPTOR>& grid) instances::getBulkMomenta<T,DESCRIPTOR>()))); sOnLatticeBoundaryCondition2D<T,DESCRIPTOR>& sBoundaryCondition = grid.getOnLatticeBoundaryCondition(); - //createInterpBoundaryCondition2D<T,DESCRIPTOR>(sBoundaryCondition); createLocalBoundaryCondition2D<T,DESCRIPTOR>(sBoundaryCondition); const T omega = converter.getLatticeRelaxationFrequency(); @@ -251,26 +250,8 @@ void takeMeasurements(Grid2D<T,DESCRIPTOR>& grid, int iT, bool print=true) } } -int main(int argc, char* argv[]) +void setupRefinement(Grid2D<T,DESCRIPTOR>& coarseGrid, Vector<T,2> domainOrigin, Vector<T,2> domainExtend) { - olbInit(&argc, &argv); - singleton::directories().setOutputDir("./tmp/"); - OstreamManager clout(std::cout,"main"); - - const Vector<T,2> coarseOrigin {0.0, 0.0}; - const Vector<T,2> coarseExtend {lx, ly}; - IndicatorCuboid2D<T> coarseCuboid(coarseExtend, coarseOrigin); - - Grid2D<T,DESCRIPTOR> coarseGrid( - coarseCuboid, - RelaxationTime<T>(tau), - N, - PhysCharacteristics); - const Vector<T,2> domainOrigin = coarseGrid.getSuperGeometry().getStatistics().getMinPhysR(0); - const Vector<T,2> domainExtend = coarseGrid.getSuperGeometry().getStatistics().getPhysExtend(0); - - prepareGeometry(coarseGrid, domainOrigin, domainExtend); - const auto coarseDeltaX = coarseGrid.getConverter().getPhysDeltaX(); const Vector<T,2> fineOutflowExtend {1*cylinderD, domainExtend[1]}; @@ -331,6 +312,28 @@ int main(int argc, char* argv[]) auto& fineGrid3 = fineGrid2.refine(fineOrigin3, fineExtend3); prepareGeometry(fineGrid3, domainOrigin, domainExtend); disableRefinedArea(fineGrid2, fineGrid3); +} + +int main(int argc, char* argv[]) +{ + olbInit(&argc, &argv); + singleton::directories().setOutputDir("./tmp/"); + OstreamManager clout(std::cout,"main"); + + const Vector<T,2> coarseOrigin {0.0, 0.0}; + const Vector<T,2> coarseExtend {lx, ly}; + IndicatorCuboid2D<T> coarseCuboid(coarseExtend, coarseOrigin); + + Grid2D<T,DESCRIPTOR> coarseGrid( + coarseCuboid, + RelaxationTime<T>(tau), + N, + PhysCharacteristics); + const Vector<T,2> domainOrigin = coarseGrid.getSuperGeometry().getStatistics().getMinPhysR(0); + const Vector<T,2> domainExtend = coarseGrid.getSuperGeometry().getStatistics().getPhysExtend(0); + + prepareGeometry(coarseGrid, domainOrigin, domainExtend); + setupRefinement(coarseGrid, domainOrigin, domainExtend); coarseGrid.forEachGrid(prepareLattice); @@ -343,6 +346,8 @@ int main(int argc, char* argv[]) coarseGrid.getSuperGeometry().getStatistics().getNvoxel()); timer.start(); + Grid2D<T,DESCRIPTOR>& cylinderGrid = coarseGrid.locate({cylinderX, cylinderY}); + for (int iT = 0; iT <= coarseGrid.getConverter().getLatticeTime(maxPhysT); ++iT) { setBoundaryValues(coarseGrid, iT); @@ -356,7 +361,7 @@ int main(int argc, char* argv[]) getResults(grid, id, iT); }); - takeMeasurements(fineGrid3, iT); + takeMeasurements(cylinderGrid, iT); } } diff --git a/src/refinement/coupler2D.hh b/src/refinement/coupler2D.hh index 92c8d26..373a6de 100644 --- a/src/refinement/coupler2D.hh +++ b/src/refinement/coupler2D.hh @@ -73,7 +73,8 @@ Coupler2D<T,DESCRIPTOR>::Coupler2D(Grid2D<T,DESCRIPTOR>& coarse, Grid2D<T,DESCRI const auto& fineGeometry = _fine.getCuboidGeometry(); const T deltaX = _fine.getConverter().getPhysDeltaX(); - const Vector<T,2> stepPhysR = _vertical ? Vector<T,2> {0, deltaX} : Vector<T,2> {deltaX, 0}; + const Vector<T,2> stepPhysR = _vertical ? Vector<T,2> {0, deltaX} : + Vector<T,2> {deltaX, 0}; for (int i=0; i < _fineSize; ++i) { if (i % 2 == 0) { @@ -212,7 +213,8 @@ void FineCoupler2D<T,DESCRIPTOR>::couple() fineLattice.get(finePos, fineCell); for (int iPop=0; iPop < DESCRIPTOR<T>::q; ++iPop) { - fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho[0], u.data, uSqr) + this->getScalingFactor() * fneq[iPop]; + fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho[0], u.data, uSqr) + + this->getScalingFactor() * fneq[iPop]; } fineLattice.set(finePos, fineCell); @@ -230,7 +232,8 @@ void FineCoupler2D<T,DESCRIPTOR>::couple() fineLattice.get(finePos, fineCell); for (int iPop=0; iPop < DESCRIPTOR<T>::q; ++iPop) { - fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho[0], u.data, uSqr) + this->getScalingFactor() * fneq[iPop]; + fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho[0], u.data, uSqr) + + this->getScalingFactor() * fneq[iPop]; } fineLattice.set(finePos, fineCell); @@ -248,7 +251,8 @@ void FineCoupler2D<T,DESCRIPTOR>::couple() fineLattice.get(finePos, fineCell); for (int iPop=0; iPop < DESCRIPTOR<T>::q; ++iPop) { - fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho[0], u.data, uSqr) + this->getScalingFactor() * fneq[iPop]; + fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho[0], u.data, uSqr) + + this->getScalingFactor() * fneq[iPop]; } fineLattice.set(finePos, fineCell); diff --git a/src/refinement/grid2D.h b/src/refinement/grid2D.h index 6cbbf6f..f538850 100644 --- a/src/refinement/grid2D.h +++ b/src/refinement/grid2D.h @@ -147,6 +147,12 @@ public: void forEachGrid(std::function<void(Grid2D<T,DESCRIPTOR>&)>&& f); void forEachGrid(const std::string& id, std::function<void(Grid2D<T,DESCRIPTOR>&,const std::string&)>&& f); + /// Returns the finest grid representing a physical position + /** + * Only works if pos is actually contained in a node of the refinement tree. + **/ + Grid2D<T,DESCRIPTOR>& locate(Vector<T,2> pos); + std::size_t getActiveVoxelN() const; }; diff --git a/src/refinement/grid2D.hh b/src/refinement/grid2D.hh index fdf6461..9ca8aaf 100644 --- a/src/refinement/grid2D.hh +++ b/src/refinement/grid2D.hh @@ -281,6 +281,18 @@ void Grid2D<T,DESCRIPTOR>::forEachGrid( } template <typename T, template<typename> class DESCRIPTOR> +Grid2D<T,DESCRIPTOR>& Grid2D<T,DESCRIPTOR>::locate(Vector<T,2> pos) +{ + int iC; + for (auto& grid : _fineGrids) { + if (grid->getCuboidGeometry().getC(pos, iC)) { + return grid->locate(pos); + } + } + return *this; +} + +template <typename T, template<typename> class DESCRIPTOR> std::size_t Grid2D<T,DESCRIPTOR>::getActiveVoxelN() const { std::size_t n = _geometry->getStatistics().getNvoxel(); |