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(-)
(limited to 'apps')
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 <vector>
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<T,DESCRIPTOR>& grid)
@@ -88,7 +87,7 @@ void prepareGeometry(Grid2D<T,DESCRIPTOR>& grid)
{
const Vector<T,2> origin {2*cylinderD, 2*cylinderD};
IndicatorCircle2D<T> obstacle(origin, cylinderD/2);
- sGeometry.rename(1,2,obstacle);
+ sGeometry.rename(1,5,obstacle);
}
sGeometry.clean();
@@ -99,6 +98,16 @@ void prepareGeometry(Grid2D<T,DESCRIPTOR>& grid)
clout << "Prepare Geometry ... OK" << std::endl;
}
+void disableRefinedArea(Grid2D<T,DESCRIPTOR>& coarseGrid,
+ RefiningGrid2D<T,DESCRIPTOR>& 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<T,DESCRIPTOR>& grid,
Dynamics<T, DESCRIPTOR>& bulkDynamics,
sOnLatticeBoundaryCondition2D<T,DESCRIPTOR>& sBoundaryCondition)
@@ -114,9 +123,10 @@ void prepareLattice(Grid2D<T,DESCRIPTOR>& grid,
sLattice.defineDynamics(sGeometry, 0, &instances::getNoDynamics<T,DESCRIPTOR>());
sLattice.defineDynamics(sGeometry, 1, &bulkDynamics); // bulk
- sLattice.defineDynamics(sGeometry, 2, &instances::getNoDynamics<T,DESCRIPTOR>());
+ sLattice.defineDynamics(sGeometry, 2, &bulkDynamics); // walls
sLattice.defineDynamics(sGeometry, 3, &bulkDynamics); // inflow
sLattice.defineDynamics(sGeometry, 4, &bulkDynamics); // outflow
+ sLattice.defineDynamics(sGeometry, 5, &instances::getBounceBack<T,DESCRIPTOR>()); // cylinder
sBoundaryCondition.addVelocityBoundary(sGeometry, 2, omega);
sBoundaryCondition.addVelocityBoundary(sGeometry, 3, omega);
@@ -125,7 +135,7 @@ void prepareLattice(Grid2D<T,DESCRIPTOR>& grid,
AnalyticalConst2D<T,T> rho0(1.0);
AnalyticalConst2D<T,T> 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<T,DESCRIPTOR>& 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<T,T> StartScale(iTmaxStart, 1.0);
+ PolynomialStartScale<T,T> StartScale(iTmaxStart, 1);
T iTvec[1] { T(iT) };
T frac[1] { };
@@ -162,7 +172,7 @@ void setBoundaryValues(Grid2D<T,DESCRIPTOR>& grid, int iT)
void getResults(const std::string& prefix,
Grid2D<T,DESCRIPTOR>& grid,
- int iT, Timer<T>& 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<T,DESCRIPTOR> coarseGrid(coarseCuboid, N, LatticeVelocity<T>(uLat), Re);
prepareGeometry(coarseGrid);
- const Vector<T,2> fineExtend {10*cylinderD, 2.75*cylinderD};
- const Vector<T,2> 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<T, DESCRIPTOR> coarseBulkDynamics(
coarseGrid.getConverter().getLatticeRelaxationFrequency(),
instances::getBulkMomenta<T, DESCRIPTOR>());
@@ -228,7 +219,13 @@ int main(int argc, char* argv[])
sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> coarseSBoundaryCondition(coarseGrid.getSuperLattice());
createLocalBoundaryCondition2D<T, DESCRIPTOR>(coarseSBoundaryCondition);
- prepareLattice(coarseGrid, coarseBulkDynamics, coarseSBoundaryCondition);
+ const Vector<T,2> fineExtend {16*cylinderD, 3.75*cylinderD};
+ const Vector<T,2> fineOrigin {0.75*cylinderD, (ly-fineExtend[1])/2};
+
+ auto& fineGrid = coarseGrid.refine(fineOrigin, fineExtend);
+ prepareGeometry(fineGrid);
+
+ disableRefinedArea(coarseGrid, fineGrid);
BGKdynamics<T, DESCRIPTOR> fineBulkDynamics(
fineGrid.getConverter().getLatticeRelaxationFrequency(),
@@ -237,17 +234,15 @@ int main(int argc, char* argv[])
sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> fineSBoundaryCondition(fineGrid.getSuperLattice());
createLocalBoundaryCondition2D<T, DESCRIPTOR>(fineSBoundaryCondition);
- const Vector<T,2> fineExtend2 {1.75*cylinderD, 1.75*cylinderD};
- const Vector<T,2> fineOrigin2 {2*cylinderD-fineExtend2[0]/2, 2*cylinderD-fineExtend2[1]/2};
+ prepareLattice(coarseGrid, coarseBulkDynamics, coarseSBoundaryCondition);
+
+ const Vector<T,2> fineExtend2 {11*cylinderD, 3*cylinderD};
+ const Vector<T,2> 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<T, DESCRIPTOR> fineBulkDynamics2(
fineGrid2.getConverter().getLatticeRelaxationFrequency(),
@@ -256,34 +251,50 @@ int main(int argc, char* argv[])
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);
+
prepareLattice(fineGrid2, fineBulkDynamics2, fineSBoundaryCondition2);
+ prepareLattice(fineGrid3, fineBulkDynamics3, fineSBoundaryCondition3);
+
clout << "Starting simulation..." << endl;
Timer<T> 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<T,DESCRIPTOR>& grid,
sLattice.defineDynamics(sGeometry, 0, &instances::getNoDynamics<T, DESCRIPTOR>());
sLattice.defineDynamics(sGeometry, 1, &bulkDynamics); // bulk
- sLattice.defineDynamics(sGeometry, 2, &instances::getNoDynamics<T, DESCRIPTOR>());
+ 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<T, DESCRIPTOR> coarseBulkDynamics(
coarseGrid.getConverter().getLatticeRelaxationFrequency(),
instances::getBulkMomenta<T, DESCRIPTOR>());
--
cgit v1.2.3