summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--apps/adrian/cylinder2d/cylinder2d.cpp89
-rw-r--r--src/refinement/coupler2D.hh6
2 files changed, 48 insertions, 47 deletions
diff --git a/apps/adrian/cylinder2d/cylinder2d.cpp b/apps/adrian/cylinder2d/cylinder2d.cpp
index 739c02c..33028c0 100644
--- a/apps/adrian/cylinder2d/cylinder2d.cpp
+++ b/apps/adrian/cylinder2d/cylinder2d.cpp
@@ -42,10 +42,7 @@ const T ly = 2.0; // height of the channel
const int N = 100; // resolution of the model
const T Re = 1000.; // Reynolds number
const T baseTau = 0.52; // Relaxation time of coarsest grid
-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
-
+const T maxPhysT = 120.; // max. simulation time in s, SI unit
void prepareGeometry(Grid2D<T,DESCRIPTOR>& grid)
{
@@ -83,10 +80,10 @@ void prepareGeometry(Grid2D<T,DESCRIPTOR>& grid)
sGeometry.rename(1,4,outflow);
}
- // Set material number for vertically centered obstacle
+ // Set material number for vertically centered cylinder
{
const Vector<T,2> origin {1.25, ly/2};
- IndicatorCircle2D<T> obstacle(origin, 0.15);
+ IndicatorCircle2D<T> obstacle(origin, 0.2);
sGeometry.rename(1,2,obstacle);
}
@@ -120,35 +117,50 @@ void prepareLattice(Grid2D<T,DESCRIPTOR>& grid,
sBoundaryCondition.addVelocityBoundary(sGeometry, 3, omega);
sBoundaryCondition.addPressureBoundary(sGeometry, 4, omega);
- const T Lx = converter.getLatticeLength(lx);
- const T Ly = converter.getLatticeLength(ly);
- const T p0 = 8.*converter.getLatticeViscosity()*converter.getCharLatticeVelocity()*Lx/(Ly*Ly);
- AnalyticalLinear2D<T,T> rho(-p0/lx*DESCRIPTOR<T>::invCs2, 0, p0*DESCRIPTOR<T>::invCs2+1);
-
- const T maxVelocity = converter.getCharLatticeVelocity();
- const T radius = ly/2;
- std::vector<T> axisPoint{0, ly/2};
- std::vector<T> axisDirection{1, 0};
- Poiseuille2D<T> u(axisPoint, axisDirection, maxVelocity, radius);
-
- sLattice.defineRhoU(sGeometry, 1, rho, u);
- sLattice.iniEquilibrium(sGeometry, 1, rho, u);
- sLattice.defineRhoU(sGeometry, 2, rho, u);
- sLattice.iniEquilibrium(sGeometry, 2, rho, u);
+ AnalyticalConst2D<T,T> rho0(1.0);
+ AnalyticalConst2D<T,T> u0(0.0, 0.0);
- sLattice.defineRhoU(sGeometry, 3, rho, u);
- sLattice.iniEquilibrium(sGeometry, 3, rho, u);
- sLattice.defineRhoU(sGeometry, 4, rho, u);
- sLattice.iniEquilibrium(sGeometry, 4, rho, u);
+ auto materials = sGeometry.getMaterialIndicator({1, 3, 4});
+ sLattice.defineRhoU(materials, rho0, u0);
+ sLattice.iniEquilibrium(materials, rho0, u0);
sLattice.initialize();
clout << "Prepare lattice ... OK" << std::endl;
}
+void setBoundaryValues(Grid2D<T,DESCRIPTOR>& grid, int iT)
+{
+ auto& converter = grid.getConverter();
+ auto& sGeometry = grid.getSuperGeometry();
+ auto& sLattice = grid.getSuperLattice();
+
+ const int iTmaxStart = converter.getLatticeTime(0.1 * maxPhysT);
+ const int iTupdate = 50;
+
+ if ( iT % iTupdate == 0 && iT <= iTmaxStart ) {
+ OstreamManager clout(std::cout,"setBoundaryValues");
+ clout << "Update inflow." << std::endl;
+
+ PolynomialStartScale<T,T> StartScale(iTmaxStart, 1.0);
+
+ T iTvec[1] { T(iT) };
+ T frac[1] { };
+ StartScale(frac, iTvec);
+
+ const T maxVelocity = converter.getCharLatticeVelocity() * frac[0];
+ const T radius = ly/2;
+ std::vector<T> axisPoint{0, ly/2};
+ std::vector<T> axisDirection{1, 0};
+ Poiseuille2D<T> u(axisPoint, axisDirection, maxVelocity, radius);
+
+ sLattice.defineU(sGeometry, 3, u);
+ }
+}
+
void getResults(const std::string& prefix,
Grid2D<T,DESCRIPTOR>& grid,
- int iT, Timer<T>& timer, bool hasConverged)
+ int iT, Timer<T>& timer)
{
OstreamManager clout(std::cout,"getResults");
@@ -174,7 +186,7 @@ void getResults(const std::string& prefix,
vtmWriter.write(iT);
}
- if (iT%statIter==0 || hasConverged) {
+ if (iT%statIter==0 ) {
timer.update(iT);
timer.printStep();
@@ -221,7 +233,7 @@ int main(int argc, char* argv[])
sOnLatticeBoundaryCondition2D<T, DESCRIPTOR> fineSBoundaryCondition(fineGrid.getSuperLattice());
createLocalBoundaryCondition2D<T, DESCRIPTOR>(fineSBoundaryCondition);
- const Vector<T,2> fineExtend2 {0.8, 0.6};
+ const Vector<T,2> fineExtend2 {1.2, 0.8};
const Vector<T,2> fineOrigin2 {0.95, (ly-fineExtend2[1])/2};
auto& fineGrid2 = fineGrid.refine(fineOrigin2, fineExtend2);
@@ -242,20 +254,14 @@ int main(int argc, char* argv[])
prepareLattice(fineGrid2, fineBulkDynamics2, fineSBoundaryCondition2);
- clout << "starting simulation..." << endl;
+ clout << "Starting simulation..." << endl;
Timer<T> timer(
coarseGrid.getConverter().getLatticeTime(maxPhysT),
coarseGrid.getSuperGeometry().getStatistics().getNvoxel());
- util::ValueTracer<T> converge(
- fineGrid2.getConverter().getLatticeTime(physInterval),
- residuum);
timer.start();
for (int iT = 0; iT < coarseGrid.getConverter().getLatticeTime(maxPhysT); ++iT) {
- if (converge.hasConverged()) {
- clout << "Simulation converged." << endl;
- break;
- }
+ setBoundaryValues(coarseGrid, iT);
coarseGrid.collideAndStream();
@@ -263,22 +269,17 @@ int main(int argc, char* argv[])
"coarse_",
coarseGrid,
iT,
- timer,
- converge.hasConverged());
+ timer);
getResults(
"fine_",
fineGrid,
iT,
- timer,
- converge.hasConverged());
+ timer);
getResults(
"fine2_",
fineGrid2,
iT,
- timer,
- converge.hasConverged());
-
- converge.takeValue(fineGrid2.getSuperLattice().getStatistics().getAverageEnergy(), true);
+ timer);
}
timer.stop();
diff --git a/src/refinement/coupler2D.hh b/src/refinement/coupler2D.hh
index 78aad46..92c8d26 100644
--- a/src/refinement/coupler2D.hh
+++ b/src/refinement/coupler2D.hh
@@ -201,9 +201,9 @@ void FineCoupler2D<T,DESCRIPTOR>::couple()
}
for (int y=1; y < this->_coarseSize-2; ++y) {
- const auto rho = order2interpolation(_c2f_rho, y);
- const auto u = order2interpolation(_c2f_u, y);
- const auto fneq = order2interpolation(_c2f_fneq, y);
+ const auto rho = order4interpolation(_c2f_rho, y);
+ const auto u = order4interpolation(_c2f_u, y);
+ const auto fneq = order4interpolation(_c2f_fneq, y);
const T uSqr = u*u;