From f544a36b51c5616670afc6c0700e8ca0c5425e79 Mon Sep 17 00:00:00 2001 From: Adrian Kummerlaender Date: Tue, 22 Jan 2019 13:28:45 +0100 Subject: Tidy up refined cylinder2d --- apps/adrian/cylinder2d/cylinder2d.cpp | 119 ++++++++++++++++-------------- apps/adrian/poiseuille2d/poiseuille2d.cpp | 12 +-- 2 files changed, 69 insertions(+), 62 deletions(-) diff --git a/apps/adrian/cylinder2d/cylinder2d.cpp b/apps/adrian/cylinder2d/cylinder2d.cpp index 2e18ebd..ffab6a0 100644 --- a/apps/adrian/cylinder2d/cylinder2d.cpp +++ b/apps/adrian/cylinder2d/cylinder2d.cpp @@ -31,21 +31,20 @@ #include using namespace olb; -using namespace olb::descriptors; typedef double T; -#define DESCRIPTOR D2Q9Descriptor +#define DESCRIPTOR descriptors::D2Q9Descriptor /// Setup geometry relative to cylinder diameter as defined by [SchaeferTurek96] const T cylinderD = 1.0; -const int N = 10; // resolution of the model +const int N = 20; // 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.00333; // lattice velocity (also used as inflow velocity) -const T maxPhysT = 120.; // max. simulation time in s, SI unit +const T uLat = 0.1; // lattice velocity +const T maxPhysT = 60.; // max. simulation time in s, SI unit void prepareGeometry(Grid2D& grid) @@ -88,7 +87,7 @@ void prepareGeometry(Grid2D& grid) { const Vector origin {2*cylinderD, 2*cylinderD}; IndicatorCircle2D obstacle(origin, cylinderD/2); - sGeometry.rename(1,2,obstacle); + sGeometry.rename(1,5,obstacle); } sGeometry.clean(); @@ -99,6 +98,16 @@ void prepareGeometry(Grid2D& grid) clout << "Prepare Geometry ... OK" << std::endl; } +void disableRefinedArea(Grid2D& coarseGrid, + RefiningGrid2D& fineGrid) +{ + auto& sGeometry = coarseGrid.getSuperGeometry(); + auto refinedOverlap = fineGrid.getRefinedOverlap(); + sGeometry.rename(1,0,*refinedOverlap); + sGeometry.rename(2,0,*refinedOverlap); + sGeometry.rename(5,0,*refinedOverlap); +} + void prepareLattice(Grid2D& grid, Dynamics& bulkDynamics, sOnLatticeBoundaryCondition2D& sBoundaryCondition) @@ -114,9 +123,10 @@ void prepareLattice(Grid2D& grid, sLattice.defineDynamics(sGeometry, 0, &instances::getNoDynamics()); sLattice.defineDynamics(sGeometry, 1, &bulkDynamics); // bulk - sLattice.defineDynamics(sGeometry, 2, &instances::getNoDynamics()); + sLattice.defineDynamics(sGeometry, 2, &bulkDynamics); // walls sLattice.defineDynamics(sGeometry, 3, &bulkDynamics); // inflow sLattice.defineDynamics(sGeometry, 4, &bulkDynamics); // outflow + sLattice.defineDynamics(sGeometry, 5, &instances::getBounceBack()); // cylinder sBoundaryCondition.addVelocityBoundary(sGeometry, 2, omega); sBoundaryCondition.addVelocityBoundary(sGeometry, 3, omega); @@ -125,7 +135,7 @@ void prepareLattice(Grid2D& grid, AnalyticalConst2D rho0(1.0); AnalyticalConst2D u0(0.0, 0.0); - auto materials = sGeometry.getMaterialIndicator({1, 2, 3, 4}); + auto materials = sGeometry.getMaterialIndicator({1, 2, 3, 4, 5}); sLattice.defineRhoU(materials, rho0, u0); sLattice.iniEquilibrium(materials, rho0, u0); @@ -140,11 +150,11 @@ void setBoundaryValues(Grid2D& grid, int iT) auto& sGeometry = grid.getSuperGeometry(); auto& sLattice = grid.getSuperLattice(); - const int iTmaxStart = converter.getLatticeTime(0.05 * maxPhysT); - const int iTupdate = 20; + const int iTmaxStart = converter.getLatticeTime(0.2*maxPhysT); + const int iTupdate = 10; if ( iT % iTupdate == 0 && iT <= iTmaxStart ) { - PolynomialStartScale StartScale(iTmaxStart, 1.0); + PolynomialStartScale StartScale(iTmaxStart, 1); T iTvec[1] { T(iT) }; T frac[1] { }; @@ -162,7 +172,7 @@ void setBoundaryValues(Grid2D& grid, int iT) void getResults(const std::string& prefix, Grid2D& grid, - int iT, Timer& timer) + int iT) { OstreamManager clout(std::cout,"getResults"); @@ -180,8 +190,6 @@ void getResults(const std::string& prefix, vtmWriter.addFunctor(pressure); vtmWriter.addFunctor(knudsen); - const int statIter = converter.getLatticeTime(maxPhysT/10.); - if (iT==0) { vtmWriter.createMasterFile(); } @@ -189,13 +197,6 @@ void getResults(const std::string& prefix, if (iT%100==0) { vtmWriter.write(iT); } - - if (iT%statIter==0 ) { - timer.update(iT); - timer.printStep(); - - sLattice.getStatistics().print(iT,converter.getPhysTime(iT)); - } } int main(int argc, char* argv[]) @@ -211,16 +212,6 @@ int main(int argc, char* argv[]) Grid2D coarseGrid(coarseCuboid, N, LatticeVelocity(uLat), Re); prepareGeometry(coarseGrid); - const Vector fineExtend {10*cylinderD, 2.75*cylinderD}; - const Vector fineOrigin {0.75*cylinderD, (ly-fineExtend[1])/2}; - - auto& fineGrid = coarseGrid.refine(fineOrigin, fineExtend); - prepareGeometry(fineGrid); - - auto refinedOverlap = fineGrid.getRefinedOverlap(); - coarseGrid.getSuperGeometry().rename(1,0,*refinedOverlap); - coarseGrid.getSuperGeometry().rename(2,0,*refinedOverlap); - BGKdynamics coarseBulkDynamics( coarseGrid.getConverter().getLatticeRelaxationFrequency(), instances::getBulkMomenta()); @@ -228,7 +219,13 @@ int main(int argc, char* argv[]) sOnLatticeBoundaryCondition2D coarseSBoundaryCondition(coarseGrid.getSuperLattice()); createLocalBoundaryCondition2D(coarseSBoundaryCondition); - prepareLattice(coarseGrid, coarseBulkDynamics, coarseSBoundaryCondition); + const Vector fineExtend {16*cylinderD, 3.75*cylinderD}; + const Vector fineOrigin {0.75*cylinderD, (ly-fineExtend[1])/2}; + + auto& fineGrid = coarseGrid.refine(fineOrigin, fineExtend); + prepareGeometry(fineGrid); + + disableRefinedArea(coarseGrid, fineGrid); BGKdynamics fineBulkDynamics( fineGrid.getConverter().getLatticeRelaxationFrequency(), @@ -237,17 +234,15 @@ int main(int argc, char* argv[]) sOnLatticeBoundaryCondition2D fineSBoundaryCondition(fineGrid.getSuperLattice()); createLocalBoundaryCondition2D(fineSBoundaryCondition); - const Vector fineExtend2 {1.75*cylinderD, 1.75*cylinderD}; - const Vector fineOrigin2 {2*cylinderD-fineExtend2[0]/2, 2*cylinderD-fineExtend2[1]/2}; + prepareLattice(coarseGrid, coarseBulkDynamics, coarseSBoundaryCondition); + + const Vector fineExtend2 {11*cylinderD, 3*cylinderD}; + const Vector fineOrigin2 {1*cylinderD, 2*cylinderD-fineExtend2[1]/2}; auto& fineGrid2 = fineGrid.refine(fineOrigin2, fineExtend2); prepareGeometry(fineGrid2); - auto refinedOverlap2 = fineGrid2.getRefinedOverlap(); - fineGrid.getSuperGeometry().rename(1,0,*refinedOverlap2); - fineGrid.getSuperGeometry().rename(2,0,*refinedOverlap2); - - prepareLattice(fineGrid, fineBulkDynamics, fineSBoundaryCondition); + disableRefinedArea(fineGrid, fineGrid2); BGKdynamics fineBulkDynamics2( fineGrid2.getConverter().getLatticeRelaxationFrequency(), @@ -256,34 +251,50 @@ int main(int argc, char* argv[]) sOnLatticeBoundaryCondition2D fineSBoundaryCondition2(fineGrid2.getSuperLattice()); createLocalBoundaryCondition2D(fineSBoundaryCondition2); + prepareLattice(fineGrid, fineBulkDynamics, fineSBoundaryCondition); + + const Vector fineExtend3 {1.5*cylinderD, 1.5*cylinderD}; + const Vector fineOrigin3 {2*cylinderD-fineExtend3[0]/2, 2*cylinderD-fineExtend3[1]/2}; + + auto& fineGrid3 = fineGrid2.refine(fineOrigin3, fineExtend3); + prepareGeometry(fineGrid3); + + disableRefinedArea(fineGrid2, fineGrid3); + + BGKdynamics fineBulkDynamics3( + fineGrid3.getConverter().getLatticeRelaxationFrequency(), + instances::getBulkMomenta()); + + sOnLatticeBoundaryCondition2D fineSBoundaryCondition3(fineGrid3.getSuperLattice()); + createLocalBoundaryCondition2D(fineSBoundaryCondition3); + prepareLattice(fineGrid2, fineBulkDynamics2, fineSBoundaryCondition2); + prepareLattice(fineGrid3, fineBulkDynamics3, fineSBoundaryCondition3); + clout << "Starting simulation..." << endl; Timer timer( coarseGrid.getConverter().getLatticeTime(maxPhysT), coarseGrid.getSuperGeometry().getStatistics().getNvoxel()); timer.start(); + const int statIter = coarseGrid.getConverter().getLatticeTime(0.5); for (int iT = 0; iT < coarseGrid.getConverter().getLatticeTime(maxPhysT); ++iT) { setBoundaryValues(coarseGrid, iT); coarseGrid.collideAndStream(); - getResults( - "coarse_", - coarseGrid, - iT, - timer); - getResults( - "fine_", - fineGrid, - iT, - timer); - getResults( - "fine2_", - fineGrid2, - iT, - timer); + getResults("level0_", coarseGrid, iT); + getResults("level1_", fineGrid, iT); + getResults("level2_", fineGrid2, iT); + getResults("level3_", fineGrid3, iT); + + if (iT%statIter == 0) { + timer.update(iT); + timer.printStep(); + + coarseGrid.getSuperLattice().getStatistics().print(iT, coarseGrid.getConverter().getPhysTime(iT)); + } } timer.stop(); diff --git a/apps/adrian/poiseuille2d/poiseuille2d.cpp b/apps/adrian/poiseuille2d/poiseuille2d.cpp index 6d0c86f..b2fb48b 100644 --- a/apps/adrian/poiseuille2d/poiseuille2d.cpp +++ b/apps/adrian/poiseuille2d/poiseuille2d.cpp @@ -105,13 +105,13 @@ void prepareLattice(Grid2D& grid, sLattice.defineDynamics(sGeometry, 0, &instances::getNoDynamics()); sLattice.defineDynamics(sGeometry, 1, &bulkDynamics); // bulk - sLattice.defineDynamics(sGeometry, 2, &instances::getNoDynamics()); + sLattice.defineDynamics(sGeometry, 2, &bulkDynamics); // walls sLattice.defineDynamics(sGeometry, 3, &bulkDynamics); // inflow sLattice.defineDynamics(sGeometry, 4, &bulkDynamics); // outflow - sBoundaryCondition.addVelocityBoundary(sGeometry, 2, omega); - sBoundaryCondition.addVelocityBoundary(sGeometry, 3, omega); - sBoundaryCondition.addPressureBoundary(sGeometry, 4, omega); + sBoundaryCondition.addVelocityBoundary(sGeometry, 2, omega); // 0-velocity walls + sBoundaryCondition.addVelocityBoundary(sGeometry, 3, omega); // velocity inflow + sBoundaryCondition.addPressureBoundary(sGeometry, 4, omega); // pressure outflow const T Lx = converter.getLatticeLength(lx); const T Ly = converter.getLatticeLength(ly); @@ -201,10 +201,6 @@ int main(int argc, char* argv[]) prepareGeometry(fineGrid); - auto refinedOverlap = fineGrid.getRefinedOverlap(); - coarseGrid.getSuperGeometry().rename(1,0,*refinedOverlap); - coarseGrid.getSuperGeometry().rename(2,0,*refinedOverlap); - BGKdynamics coarseBulkDynamics( coarseGrid.getConverter().getLatticeRelaxationFrequency(), instances::getBulkMomenta()); -- cgit v1.2.3