diff --git a/apps/tutorials/lbm/04_LBComplexGeometry.cpp b/apps/tutorials/lbm/04_LBComplexGeometry.cpp
index e33fb6ee07b7970c3283a592616d161069fd9677..fdab4509e2256b1d09030ed9c08376d10d0af674 100644
--- a/apps/tutorials/lbm/04_LBComplexGeometry.cpp
+++ b/apps/tutorials/lbm/04_LBComplexGeometry.cpp
@@ -57,15 +57,18 @@ namespace walberla
{
uint_t numGhostLayers = uint_t(1);
+//! [typedefs]
using LatticeModel_T = lbm::D3Q27< lbm::collision_model::SRT >;
using Stencil_T = LatticeModel_T::Stencil;
using CommunicationStencil_T = LatticeModel_T::CommunicationStencil;
+//! [typedefs]
using PdfField_T = lbm::PdfField< LatticeModel_T >;
using flag_t = walberla::uint8_t;
using FlagField_T = FlagField< flag_t >;
+//! [vertexToFaceColor]
template< typename MeshType >
void vertexToFaceColor(MeshType& mesh, const typename MeshType::Color& defaultColor)
{
@@ -93,6 +96,7 @@ void vertexToFaceColor(MeshType& mesh, const typename MeshType::Color& defaultCo
mesh.set_color(*faceIt, useVertexColor ? vertexColor : defaultColor);
}
}
+//! [vertexToFaceColor]
int main(int argc, char** argv)
{
@@ -115,10 +119,12 @@ int main(int argc, char** argv)
const double remainingTimeLoggerFrequency =
parameters.getParameter< double >("remainingTimeLoggerFrequency", 3.0); // in seconds
+ //! [parseDomainParameters]
// read domain parameters
auto domainParameters = walberlaEnv.config()->getOneBlock("DomainSetup");
std::string meshFile = domainParameters.getParameter< std::string >("meshFile");
+ //! [parseDomainParameters]
Vector3< uint_t > domainScaling =
domainParameters.getParameter< Vector3< uint_t > >("domainScaling", Vector3< uint_t >(1));
@@ -134,19 +140,25 @@ int main(int argc, char** argv)
WALBERLA_LOG_INFO_ON_ROOT("Using mesh from " << meshFile << ".")
+ //! [readMesh]
// read in mesh with vertex colors on a single process and broadcast it
auto mesh = make_shared< mesh::TriangleMesh >();
mesh->request_vertex_colors();
mesh::readAndBroadcast(meshFile, *mesh);
+ //! [readMesh]
// color faces according to vertices
vertexToFaceColor(*mesh, mesh::TriangleMesh::Color(255, 255, 255));
+ //! [triDist]
// add information to mesh that is required for computing signed distances from a point to a triangle
auto triDist = make_shared< mesh::TriangleDistance< mesh::TriangleMesh > >(mesh);
+ //! [triDist]
+ //! [octree]
// building distance octree
auto distanceOctree = make_shared< mesh::DistanceOctree< mesh::TriangleMesh > >(triDist);
+ //! [octree]
WALBERLA_LOG_INFO_ON_ROOT("Octree has height " << distanceOctree->height())
@@ -157,17 +169,25 @@ int main(int argc, char** argv)
/// CREATE BLOCK FOREST ///
///////////////////////////
+ //! [aabb]
auto aabb = computeAABB(*mesh);
aabb.scale(domainScaling);
aabb.setCenter(aabb.center() + 0.2 * Vector3< real_t >(aabb.xSize(), 0, 0));
+ //! [aabb]
+ //! [bfc]
// create and configure block forest creator
mesh::ComplexGeometryStructuredBlockforestCreator bfc(aabb, Vector3< real_t >(dx),
mesh::makeExcludeMeshInterior(distanceOctree, dx));
+
bfc.setPeriodicity(periodicity);
+ //! [bfc]
+
+ //! [blockForest]
// create block forest
auto blocks = bfc.createStructuredBlockForest(cellsPerBlock);
+ //! [blockForest]
////////////////////////////////////
/// CREATE AND INITIALIZE FIELDS ///
@@ -183,6 +203,7 @@ int main(int argc, char** argv)
/// BOUNDARY HANDLING ///
/////////////////////////
+ //! [DefaultBoundaryHandling]
// create and initialize boundary handling
const FlagUID fluidFlagUID("Fluid");
@@ -197,9 +218,9 @@ int main(int argc, char** argv)
boundariesConfig.getParameter< real_t >("pressure0", real_c(1.0)),
boundariesConfig.getParameter< real_t >("pressure1", real_c(1.0)));
- geometry::initBoundaryHandling< BHFactory::BoundaryHandling >(*blocks, boundaryHandlingId, boundariesConfig);
- geometry::setNonBoundaryCellsToDomain< BHFactory::BoundaryHandling >(*blocks, boundaryHandlingId);
+ //! [DefaultBoundaryHandling]
+ //! [colorToBoundary]
// set NoSlip UID to boundaries that we colored
mesh::ColorToBoundaryMapper< mesh::TriangleMesh > colorToBoundaryMapper(
(mesh::BoundaryInfo(BHFactory::getNoSlipBoundaryUID())));
@@ -208,23 +229,30 @@ int main(int argc, char** argv)
// mark boundaries
auto boundaryLocations = colorToBoundaryMapper.addBoundaryInfoToMesh(*mesh);
+ //! [colorToBoundary]
+ //! [VTKMesh]
// write mesh info to file
mesh::VTKMeshWriter< mesh::TriangleMesh > meshWriter(mesh, "meshBoundaries", 1);
meshWriter.addDataSource(make_shared< mesh::BoundaryUIDFaceDataSource< mesh::TriangleMesh > >(boundaryLocations));
meshWriter.addDataSource(make_shared< mesh::ColorFaceDataSource< mesh::TriangleMesh > >());
meshWriter.addDataSource(make_shared< mesh::ColorVertexDataSource< mesh::TriangleMesh > >());
meshWriter();
+ //! [VTKMesh]
+ //! [boundarySetup]
// voxelize mesh
mesh::BoundarySetup boundarySetup(blocks, makeMeshDistanceFunction(distanceOctree), numGhostLayers);
// set fluid cells
boundarySetup.setDomainCells< BHFactory::BoundaryHandling >(boundaryHandlingId, mesh::BoundarySetup::OUTSIDE);
- // set up boundaries
+ // set up inflow/outflow/wall boundaries from DefaultBoundaryHandlingFactory
+ geometry::initBoundaryHandling< BHFactory::BoundaryHandling >(*blocks, boundaryHandlingId, boundariesConfig);
+ // set up obstacle boundaries from file
boundarySetup.setBoundaries< BHFactory::BoundaryHandling >(
boundaryHandlingId, makeBoundaryLocationFunction(distanceOctree, boundaryLocations), mesh::BoundarySetup::INSIDE);
+ //! [boundarySetup]
//////////////////////////////////
/// SET UP SWEEPS AND TIMELOOP ///
@@ -258,8 +286,9 @@ int main(int argc, char** argv)
//////////////////
// add VTK output to time loop
- uint_t vtkWriteFrequency = parameters.getParameter("writeFrequency", uint_c(0));
- auto vtkOutput = vtk::createVTKOutput_BlockData(*blocks, "fluid_field", vtkWriteFrequency, numGhostLayers, false,
+ auto VTKParams = walberlaEnv.config()->getBlock("VTK");
+ uint_t vtkWriteFrequency = VTKParams.getBlock("fluid_field").getParameter("writeFrequency", uint_t(0));
+ auto vtkOutput = vtk::createVTKOutput_BlockData(*blocks, "fluid_field", vtkWriteFrequency, uint_t(0), false,
"vtk_out", "simulation_step", false, true, true, false, 0);
field::FlagFieldCellFilter< FlagField_T > fluidFilter(flagFieldId);
diff --git a/apps/tutorials/lbm/04_LBComplexGeometry.dox b/apps/tutorials/lbm/04_LBComplexGeometry.dox
index 02ee5fcba38e813c5041b83f6a6a5a65aa21f128..dd1d7e74ea03aa1fb39565c727493508bf346d49 100644
--- a/apps/tutorials/lbm/04_LBComplexGeometry.dox
+++ b/apps/tutorials/lbm/04_LBComplexGeometry.dox
@@ -4,7 +4,8 @@ namespace walberla {
\page tutorial_lbm04 Tutorial - LBM 4: Complex Geometry
\brief A configurable application for LBM simulations with complex geometries
-DISCLAIMER: for this tutorial, you must have built `waLBerla` with `WALBERLA_BUILD_WITH_OPENMESH` enabled!
+
+> _NOTE:_ for this tutorial, you must have built `waLBerla` with `WALBERLA_BUILD_WITH_OPENMESH` enabled!
\section tutorial04_overview Overview
@@ -13,7 +14,7 @@ It has following features:
- make use of default LBM sweeps from the lbm module
- can handle some basic boundary conditions: no slip, free slip, pressure, and velocity boundaries
-- initialize geometry of the domain importing a mesh file
+- initialize geometry of the domain importing a mesh file (must be exported as triangle mesh in this tutorial)
- the boundary and geometry can be configured in a parameter file
The LBM routine mostly follows what has been demonstrated in LBM Tutorial 1. However, as our aim to load arbitrary complex geometry from a mesh file, we cannot specify the geometry via the Boundaries Block in the configuration file.
@@ -29,11 +30,7 @@ In this application, we will again simulate a channel flow around an obstacle. A
\subsection tutorial04_latticemodel Lattice Model
In this tutorial, we will use the standard D3Q27 stencil with a SRT collision model.
-\code
-using LatticeModel_T = lbm::D3Q27<lbm::collision_model::SRT>;
-using Stencil_T = LatticeModel_T::Stencil;
-using CommunicationStencil_T = LatticeModel_T::CommunicationStencil;
-\endcode
+\snippet 04_LBComplexGeometry.cpp typedefs
As this is only the basic LBM setup, other stencils and collision methods can be included to simulate more accurate, stable and complex flows.
@@ -42,142 +39,51 @@ As this is only the basic LBM setup, other stencils and collision methods can be
When building waLBerla with the CMAKE configuration `WALBERLA_BUILD_WITH_OPENMESH` enabled, one can use its capabilities for loading, storing and manipulating meshes.
OpenMesh supports many common mesh data formats. We will use an `.obj` file.
+> _*NOTE:*_ with the functions used in this tutorial, only triangle meshes are processed correctly. For more information, please refer to the walberla::mesh module.
Reading meshes into waLBerla is relatively straightforward, even for large scale simulations.
Firstly, we define a block in our parameter file responsible for the mesh and the domain:
-\code
-DomainSetup
-{
- meshFile bunny.obj;
- domainScaling < 10, 3, 1 >;
- dx 3;
- cellsPerBlock < 16, 16, 16 >;
- periodic < 0, 0, 1 >;
-}
-\endcode
+\snippet 04_LBComplexGeometry.prm domainSetup
We will comment on the meaning behind these parameters when they are first used in the code snippets.
The name of the file in which our mesh is stored (here: `bunny.obj` as we want to simulate a flow around the Stanford Bunny) is parsed as usual:
-\code
- auto domainParameters = walberlaEnv.config()->getOneBlock( "DomainSetup" );
- std::string meshFile = domainParameters.getParameter< std::string >( "meshFile" );
-\endcode
+\snippet 04_LBComplexGeometry.cpp parseDomainParameters
Afterwards, the mesh is read in on a single process and broadcasted to all other processes with
-\code
- auto mesh = make_shared< mesh::TriangleMesh >();
- mesh->request_vertex_colors();
- mesh::readAndBroadcast( meshFile, *mesh );
-\endcode
+\snippet 04_LBComplexGeometry.cpp readMesh
Note that we have requested vertex colors from the mesh for the purpose of coloring the faces of the object according to its vertices in the next step.
Therefore we define a helper function `vertexToFaceColor` that iterates over all faces and colors them in their vertex color. If no uniform coloring of the vertices is given, a default color is taken.
-\code
-template< typename MeshType >
-void vertexToFaceColor( MeshType & mesh, const typename MeshType::Color & defaultColor )
-{
- WALBERLA_CHECK( mesh.has_vertex_colors() );
- mesh.request_face_colors();
-
- for( auto faceIt = mesh.faces_begin(); faceIt != mesh.faces_end(); ++faceIt )
- {
- typename MeshType::Color vertexColor;
-
- bool useVertexColor = true;
-
- auto vertexIt = mesh.fv_iter( *faceIt );
- WALBERLA_ASSERT( vertexIt.is_valid() );
-
- vertexColor = mesh.color( *vertexIt );
-
- ++vertexIt;
- while( vertexIt.is_valid() && useVertexColor )
- {
- if( vertexColor != mesh.color( *vertexIt ) )
- useVertexColor = false;
- ++vertexIt;
- }
-
- mesh.set_color( *faceIt, useVertexColor ? vertexColor : defaultColor );
- }
-}
-\endcode
+\snippet 04_LBComplexGeometry.cpp vertexToFaceColor
+
After calling this function, we prepare for building the distance octree by precalculating information (e.g. face normals) of the mesh that is required for computing the singed distances from a point to a triangle:
-\code
-auto triDist = make_shared< mesh::TriangleDistance< mesh::TriangleMesh > >( mesh );
-\endcode
+\snippet 04_LBComplexGeometry.cpp triDist
From this information we can finally build the distance octree. It stores information about how close or far boundaries are to each other. Later, this information could be used for e.g. adaptive mesh refinement (note that this will not be covered in this tutorial).
-\code
-auto distanceOctree = make_shared< mesh::DistanceOctree< mesh::TriangleMesh > >( triDist );
-\endcode
+\snippet 04_LBComplexGeometry.cpp octree
Even though we have successfully loaded the complex geometry and set up the corresponding distance octree, we have not defined our computational LB domain yet.
In this tutorial, the LB domain is defined relatively to the loaded geometry. Henceforth, we calculate the axis-aligned bounding box of the geometry and scale it to our needs.
Here, we chose our channel to be 10x3x1 times the size of the Stanford Bunny. This scaling is defined in the parameter file (parameter: domainScaling).
As the bunny will be placed in the center of the bounding box, we shift the center to the left such that the bunny will be nearer to the inflow.
-\code
-auto aabb = computeAABB( *mesh );
-aabb.scale( domainScaling );
-aabb.setCenter( aabb.center() + 0.2 * Vector3<real_t>(aabb.xSize(), 0, 0) );
-\endcode
-
-Finally, we use a convenient built-in data structure that will be responsible for the creation of the structured block forest.
-\code
-mesh::ComplexGeometryStructuredBlockforestCreator bfc( aabb, Vector3< real_t >( dx ), mesh::makeExcludeMeshInterior( distanceOctree, dx ) );
-\endcode
+\snippet 04_LBComplexGeometry.cpp aabb
+
+Finally, we use a convenient built-in data structure that will be responsible for the creation of the structured block forest and set the periodicity according to the parameter file.
+\snippet 04_LBComplexGeometry.cpp bfc
The `ComplexGeometryStructuredBlockforestCreator` takes as arguments the axis-aligned bounding box of the domain, the cell sizes `dx` and an exclusion function.
In this tutorial, we want to exclude the interior of the Stanford Bunny with a maximum error of `dx`.
-After setting the periodicity obtained from the parameter file
-\code
-bfc.setPeriodicity( periodicity );
-\endcode
-we create the structured block forest on which the sweeps will be performed
-\code
-auto blocks = bfc.createStructuredBlockForest( cellsPerBlock );
-\endcode
+Afterwards, we create the structured block forest on which the sweeps will be performed
+\snippet 04_LBComplexGeometry.cpp blockForest
\subsection tutorial04_boundaryhandling Boundary Handling
The rest of the tutorial mainly follows Tutorial 01 and should therefore be self-explanatory. The only adaption that has to be made is the treatment of the boundaries.
Whereas the boundary conditions of the basic domain (inflow, outflow, no-slip) can be again treated by the lbm::DefaultBoundaryHandlingFactory, the no-slip boundary conditions of the bunny need special attention.
After the standard routine of the lbm::DefaultBoundaryHandlingFactory
-\code
-const FlagUID fluidFlagUID( "Fluid" );
-
-auto boundariesConfig = walberlaEnv.config()->getOneBlock( "Boundaries" );
-
-typedef lbm::DefaultBoundaryHandlingFactory< LatticeModel_T, FlagField_T > BHFactory;
-
-BlockDataID boundaryHandlingId = BHFactory::addBoundaryHandlingToStorage( blocks, "boundary handling", flagFieldId, pdfFieldId, fluidFlagUID,
- boundariesConfig.getParameter< Vector3<real_t> >( "velocity0", Vector3<real_t>() ),
- boundariesConfig.getParameter< Vector3<real_t> >( "velocity1", Vector3<real_t>() ),
- boundariesConfig.getParameter< real_t > ( "pressure0", real_c( 1.0 ) ),
- boundariesConfig.getParameter< real_t > ( "pressure1", real_c( 1.0 ) ) );
-
-geometry::initBoundaryHandling<BHFactory::BoundaryHandling>( *blocks, boundaryHandlingId, boundariesConfig );
-geometry::setNonBoundaryCellsToDomain<BHFactory::BoundaryHandling> ( *blocks, boundaryHandlingId );
-\endcode
+\snippet 04_LBComplexGeometry.cpp DefaultBoundaryHandling
we address the complex geometry. In a first step, the boundaries that we colored in `vertexToFaceColor` are equipped with the no-slip boundary UID and the boundary information is added to the mesh.
-\code
-mesh::ColorToBoundaryMapper<mesh::TriangleMesh> colorToBoundaryMapper(( mesh::BoundaryInfo( BHFactory::getNoSlipBoundaryUID() ) ));
-colorToBoundaryMapper.set( mesh::TriangleMesh::Color(255,255,255), mesh::BoundaryInfo( BHFactory::getNoSlipBoundaryUID() ) );
-
-auto boundaryLocations = colorToBoundaryMapper.addBoundaryInfoToMesh( *mesh );
-\endcode
+\snippet 04_LBComplexGeometry.cpp colorToBoundary
In order to have the mesh information available at postprocessing, we write this information to a VTK file as
-\code
-mesh::VTKMeshWriter< mesh::TriangleMesh > meshWriter( mesh, "meshBoundaries", 1 );
-meshWriter.addDataSource( make_shared< mesh::BoundaryUIDFaceDataSource< mesh::TriangleMesh > >( boundaryLocations ) );
-meshWriter.addDataSource( make_shared< mesh::ColorFaceDataSource< mesh::TriangleMesh > >() );
-meshWriter.addDataSource( make_shared< mesh::ColorVertexDataSource< mesh::TriangleMesh > >() );
-meshWriter();
-\endcode
+\snippet 04_LBComplexGeometry.cpp VTKMesh
Lastly, the boundary information of the complex geometry is added to the structured block forest and the fluid cells and the boundaries are set accordingly.
-\code
-mesh::BoundarySetup boundarySetup( blocks, makeMeshDistanceFunction( distanceOctree ), numGhostLayers );
-
-boundarySetup.setDomainCells< BHFactory::BoundaryHandling > ( boundaryHandlingId, mesh::BoundarySetup::OUTSIDE );
-
-boundarySetup.setBoundaries<BHFactory::BoundaryHandling>( boundaryHandlingId, makeBoundaryLocationFunction( distanceOctree, boundaryLocations ), mesh::BoundarySetup::INSIDE );
-\endcode
+\snippet 04_LBComplexGeometry.cpp boundarySetup
\section tutorial04_outlook Outlook
Now that you have seen how to load and process complex geometries from a mesh file, you can play around with this tutorial and adapt it to your needs.
diff --git a/apps/tutorials/lbm/04_LBComplexGeometry.prm b/apps/tutorials/lbm/04_LBComplexGeometry.prm
index bb88c074fccb1df34f1f293487a1ce09b92958b4..5509aa3a81a44a4594b3894804639486aa47c979 100644
--- a/apps/tutorials/lbm/04_LBComplexGeometry.prm
+++ b/apps/tutorials/lbm/04_LBComplexGeometry.prm
@@ -9,6 +9,7 @@ Parameters
remainingTimeLoggerFrequency 3; // in seconds
}
+//! [domainSetup]
DomainSetup
{
meshFile bunny.obj;
@@ -18,6 +19,7 @@ DomainSetup
cellsPerBlock < 16, 16, 16 >;
periodic < 0, 0, 1 >;
}
+//! [domainSetup]
StabilityChecker
{