diff options
| -rw-r--r-- | apps/adrian/cylinder2d/cylinder2d.cpp | 33 |
1 files changed, 15 insertions, 18 deletions
diff --git a/apps/adrian/cylinder2d/cylinder2d.cpp b/apps/adrian/cylinder2d/cylinder2d.cpp index 066bb6f..5b67408 100644 --- a/apps/adrian/cylinder2d/cylinder2d.cpp +++ b/apps/adrian/cylinder2d/cylinder2d.cpp @@ -38,21 +38,21 @@ typedef double T; /// Setup geometry relative to cylinder diameter as defined by [SchaeferTurek96] const T cylinderD = 0.1; -const int N = 5; // resolution of the cylinder +const int N = 10; // resolution of the cylinder const T deltaR = cylinderD / N; // coarse lattice spacing const T lx = 22*cylinderD + deltaR; // length of the channel const T ly = 4.1*cylinderD + deltaR; // height of the channel const T cylinderX = 2*cylinderD; const T cylinderY = 2*cylinderD + deltaR/2; -const T Re = 20.; // Reynolds number +const T Re = 100.; // Reynolds number const T tau = 0.51; // relaxation time const T maxPhysT = 16.; // max. simulation time in s, SI unit const Characteristics<T> PhysCharacteristics( 0.1, // char. phys. length - 0.2, // char. phys. velocity - 0.1*0.2/Re, // phsy. kinematic viscosity + 1.0, // char. phys. velocity + 0.1/Re, // phsy. kinematic viscosity 1.0); // char. phys. density void prepareGeometry(Grid2D<T,DESCRIPTOR>& grid, Vector<T,2> origin, Vector<T,2> extend) @@ -143,19 +143,19 @@ 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::getBounceBack<T,DESCRIPTOR>()); // walls + sLattice.defineDynamics(sGeometry, 2, &bulkDynamics); // walls sLattice.defineDynamics(sGeometry, 3, &bulkDynamics); // inflow sLattice.defineDynamics(sGeometry, 4, &bulkDynamics); // outflow - sLattice.defineDynamics(sGeometry, 5, &bulkDynamics); // cylinder + sLattice.defineDynamics(sGeometry, 5, &instances::getBounceBack<T,DESCRIPTOR>()); // cylinder + sBoundaryCondition.addVelocityBoundary(sGeometry, 2, omega); sBoundaryCondition.addVelocityBoundary(sGeometry, 3, omega); sBoundaryCondition.addPressureBoundary(sGeometry, 4, omega); - sBoundaryCondition.addVelocityBoundary(sGeometry, 5, omega); AnalyticalConst2D<T,T> rho0(1.0); AnalyticalConst2D<T,T> u0(0.0, 0.0); - auto materials = sGeometry.getMaterialIndicator({1, 3, 4, 5}); + auto materials = sGeometry.getMaterialIndicator({1, 2, 3, 4}); sLattice.defineRhoU(materials, rho0, u0); sLattice.iniEquilibrium(materials, rho0, u0); @@ -238,11 +238,8 @@ void takeMeasurements(Grid2D<T,DESCRIPTOR>& grid) intpolatePressure(&pressureInFrontOfCylinder, point1); intpolatePressure(&pressureBehindCylinder, point2); - clout << "pressure1=" << pressureInFrontOfCylinder; - clout << "; pressure2=" << pressureBehindCylinder; - T pressureDrop = pressureInFrontOfCylinder - pressureBehindCylinder; - clout << "; pressureDrop=" << pressureDrop; + clout << "pressureDrop=" << pressureDrop; const int input[3] {}; T drag[dragF.getTargetDim()] {}; @@ -275,22 +272,22 @@ int main(int argc, char* argv[]) const auto coarseDeltaX = coarseGrid.getConverter().getPhysDeltaX(); - const Vector<T,2> fineExtend {10*cylinderD, ly-2*coarseDeltaX}; - const Vector<T,2> fineOrigin {1*coarseDeltaX, coarseDeltaX}; + const Vector<T,2> fineExtend {10*cylinderD, domainExtend[1]-4*coarseDeltaX}; + const Vector<T,2> fineOrigin {0.5*cylinderD, (domainExtend[1]-fineExtend[1])/2}; auto& fineGrid = coarseGrid.refine(fineOrigin, fineExtend); prepareGeometry(fineGrid, domainOrigin, domainExtend); disableRefinedArea(coarseGrid, fineGrid); - const Vector<T,2> fineExtend2 {6*cylinderD, fineGrid.getExtend()[1]-2*coarseDeltaX}; - const Vector<T,2> fineOrigin2 {2*coarseDeltaX, (ly-fineExtend2[1])/2}; + const Vector<T,2> fineExtend2 {6*cylinderD, fineGrid.getExtend()[1]-4*coarseDeltaX}; + const Vector<T,2> fineOrigin2 {0.75*cylinderD, (domainExtend[1]-fineExtend2[1])/2}; auto& fineGrid2 = fineGrid.refine(fineOrigin2, fineExtend2); prepareGeometry(fineGrid2, domainOrigin, domainExtend); disableRefinedArea(fineGrid, fineGrid2); const Vector<T,2> fineExtend3 {4*cylinderD, 2*cylinderD}; - const Vector<T,2> fineOrigin3 {1.2*cylinderD, (ly-fineExtend3[1])/2}; + const Vector<T,2> fineOrigin3 {1*cylinderD, (domainExtend[1]-fineExtend3[1])/2}; auto& fineGrid3 = fineGrid2.refine(fineOrigin3, fineExtend3); prepareGeometry(fineGrid3, domainOrigin, domainExtend); @@ -312,7 +309,7 @@ int main(int argc, char* argv[]) clout << "Total number of active cells: " << coarseGrid.getActiveVoxelN() << endl; clout << "Starting simulation..." << endl; - const int statIter = coarseGrid.getConverter().getLatticeTime(0.1); + const int statIter = coarseGrid.getConverter().getLatticeTime(0.05); Timer<T> timer( coarseGrid.getConverter().getLatticeTime(maxPhysT), coarseGrid.getSuperGeometry().getStatistics().getNvoxel()); |