Commit 2e6512fe authored by Dennis Gläser's avatar Dennis Gläser
Browse files

Merge branch 'feature/examples-readmes' into 'master'

Feature/examples readmes

See merge request DennisGlaeser/frackit!47
parents 549dc65b bf1c2330
......@@ -61,8 +61,8 @@ blocks enables users to fully customize each step to their needs. Three exemplar
applications are contained in this repository:
* [Example 1][5] Generation of a simple network consisting of quadrilaterals with two main orientations.
* [Example 2][6] Generation of a network of quadrilaterals embedded in a cylindrical domain, mimicking a core sample of a fractured porous medium.
* [Example 3][4] Generation of a network consisting of both disks and quadrilaterals, embedded in one layer of a domain that composed of three layers.
* [Example 2][6] Generation of a network of quadrilaterals embedded in a cylindrical domain (mimicking a core sample).
* [Example 3][4] Generation of a network consisting of both disks and quadrilaterals, confined and contained in one layer of a domain that is composed of three layers.
Documentation
......@@ -81,14 +81,17 @@ Installation
Please note that the following requirements need to be installed:
* OpenCascade (>= 7.3.0)
* CMake (>2.8.12)
* CMake (>3.0)
* C, C++ compiler (C++17 required)
* Optional: Doxygen (>= 1.8)
### Installation of OpenCascade
Frackit requires the [OpenCascade][2] library to be installed on your system.
You can download the source code at https://www.opencascade.com/content/download-center,
and for details on the installation we refer to TODO:PUT_LINK.
You can download the source code [HERE][2],
and details on the installation can be found [HERE][10].
Please note that [OpenCascade][2] requires further 3rd party libraries, of which
the mandatory ones are Tcl/Tk and FreeType (see this [link][11]).
On Ubuntu, both of these can be installed from the command line.
### Building Frackit under Linux
After [OpenCascade][2] and the other requirements listed above have been installed,
......@@ -173,7 +176,6 @@ Contributing
=============
Contributions are highly appreciated.
Before you start coding, please have a look at the [styleguide][doc/styleguide.md].
For bug reports, please file an [issue](https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/issues).
If you want to contribute with new features of improvements of existing code, please
......@@ -186,6 +188,8 @@ of the Frackit repository as the target branch.
before merging.
* If you have developer status you don't need to do a fork and you can create branches directly.
In this project, we follow the [styleguide][2] of the [DuMuX][0] project.
Please have a look at these before you start coding your contributions.
[0]: https://dumux.org
[1]: http://gmsh.info/
......@@ -197,3 +201,7 @@ before merging.
[7]: https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/tree/master/appl/
[8]: http://www.doxygen.org/index.html
[9]: https://www.gnu.org/licenses/gpl-3.0.en.html
[10]: https://www.opencascade.com/doc/occt-6.9.1/overview/html/occt_dev_guides__building_cmake.html
[11]: https://www.opencascade.com/doc/occt-6.9.1/overview/html/occt_dev_guides__building_3rdparty_linux.html
[12]: https://git.iws.uni-stuttgart.de/dumux-repositories/dumux/blob/master/doc/styleguide.md
<!--- Example picture --->
<p align="center">
<img src="../../doc/img/example1_network.png" alt="frackit example 1" width="800"/>
</p>
Example 1
=========
In this exemplary application, a network of quadrilateral fractures is generated
within the unit cube (see image above). The main file containing the source code
to this example is the file `example1.cc` which is located in this folder.
Two main orientations are considered for the quadrilaterals, for both of which
a corresponding instance of the `QuadrilateralSampler` class is created.
For example, we instantiate the sampler class for the second orientation with:
```cpp
using QuadSampler = QuadrilateralSampler<worldDimension>;
using Distro = std::normal_distribution<ctype>;
QuadSampler quadSampler2(makeUniformPointSampler(domain), // use a new point sampler instance!
Distro(toRadians(0.0), toRadians(5.0)), // strike angle: mean value & standard deviation
Distro(toRadians(0.0), toRadians(5.0)), // dip angle: mean value & standard deviation
Distro(0.5, 0.1), // edge length: mean value & standard deviation
0.05); // threshold for minimum edge length
```
The first constructor argument is a point sampler with which the center points of
the quadrilaterals are sampled. Here we use uniformly sampled points in the unit
cube, which is represented by an instance of the `Box` class, stored in the
variable `domain`. The second and third arguments define the distributions for
the strike and dip angle, where in this case we use uniform distributions with
a mean value of 0° and a standard deviation of 5°. The fourth argument is the
distribution to be used for sampling the edge lengths, while the last argument
defines a minimum value below which the edge length must not fall.
The quadrilaterals are then sampled from the two samplers `quadSampler1` and
`quadSampler2`, using the `()` operator:
```cpp
auto quad = sampleIntoSet1 ? quadSampler1() : quadSampler2();
```
In this example we use the boolean variable `sampleIntoSet1` to determine from
which sampler we should sample the next quadrilateral (more details follow below).
The variable `quad` holds a new candidate for an entity of the network, however,
However, we want to enforce certain constraints such as a minimum distance between
entities. For this we use instances of the `EntityNetworkConstraints` class and
configure it as desired. For example, the constraints on entities of the same
orientation are defined in this example as follows:
```cpp
EntityNetworkConstraints constraintsOnSelf;
constraintsOnSelf.setMinDistance(0.05);
constraintsOnSelf.setMinIntersectingAngle(toRadians(30.0));
constraintsOnSelf.setMinIntersectionMagnitude(0.05);
constraintsOnSelf.setMinIntersectionDistance(0.05);
```
In the main loop of quadrilateral generation, the fulfilment of these constraints is
evaluated against the other quadrilaterals with:
```cpp
auto& entitySet = sampleIntoSet1 ? entitySet1 : entitySet2;
if (!constraintsOnSelf.evaluate(entitySet, quad))
{ status.increaseRejectedCounter(); continue; }
```
where `entityset1` and `entitySet2` are of type `std::vector<Quadrilateral>` and
store all quadrilaterals that are accepted. The function `evaluate` of the
`EntityNetworkConstraints` class evaluates the constraints for `quad` against all
entities contained in `entitySet` and returns `true` only if there no violation of
any of the defined constraints has been found. After an admissible quadrilateral
has been generated, the line
```cpp
// sample into the other set the next time
sampleIntoSet1 = !sampleIntoSet1;
```
at the end of the loop makes sure that a quadrilateral of the other orientation
is sampled next. In [Example 3][0] we will get to know how to use helper classes
that store different entity sets and automatically sample from various sampler
classes such that this does not have to be done manually.
After the desired number of entities has been generated, the entities are cast
into an entity network using the builder class:
```cpp
EntityNetworkBuilder builder;
builder.addEntities(entitySet1);
builder.addEntities(entitySet2);
const auto network = builder.build();
```
This network can then be written to disk, for example in [Gmsh][1] (.geo) file format:
```cpp
GmshWriter writer(network);
writer.write("network", // filename of the .geo files (will add extension .geo automatically)
0.1); // element size to be used
```
[0]: https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/tree/master/appl/example3
[1]: http://gmsh.info/
<!--- Example picture --->
<p align="center">
<img src="../../doc/img/example2_network.png" alt="frackit example 1" width="800"/>
</p>
Example 2
=========
In this exemplary application, a network of quadrilaterals, following the same
distributions as the network of [example 1][0], is created. However, in this
example we want to confined the network to a cylindrical domain (see image above).
We represent the domain via an instance of the internal cylinder class:
```cpp
// In this example we consider a cylindrical domain
Cylinder<ctype> domain(/*radius*/0.5, /*height*/1.0);
```
This creates an axis-aligned cylinder where the x- and y- axis form the basis of
the bottom face of the cylinder, while the cylinder axis is aligned with the
z-axis. Note that there are other ways to construct arbitrarily-oriented cylinders,
and for details on this we refer to the [class documentation][2].
In addition to constraints among entities of the different orientations, we also
want to enforce constraints with respect to the boundary of the cylindrical domain.
In this example, we simply reuse the constraints of entities of the same orientation
and write
```cpp
auto constraintsOnBoundary = constraintsOnSelf;
```
to create a new constraints object. Within the loop in which the entities are created
these constraints are enforced:
```cpp
// we want to enforce constraints also w.r.t. to the cylinder boundary
if (!constraintsOnBoundary.evaluate(domain.topFace(), quad))
{ status.increaseRejectedCounter(); continue; }
if (!constraintsOnBoundary.evaluate(domain.bottomFace(), quad))
{ status.increaseRejectedCounter(); continue; }
if (!constraintsOnBoundary.evaluate(domain.lateralFace(), quad))
{ status.increaseRejectedCounter(); continue; }
```
As you can see, the `Cylinder` class provides functions for obtaining the representations
of the top, bottom and lateral boundaries. The top and bottom boundaries are represented
by instances of the `Disk` class, while the lateral surface is described by the class
`CylinderSurface`. Moreover, we want to reject all those quadrilaterals of which only
a very small portion is inside the cylinder. This is done in the lines
```cpp
const auto containedArea = computeContainedMagnitude(quad, domain);
// reject if this is too small (< 0.01 m^2)
if (containedArea < 0.01)
{ status.increaseRejectedCounter(); continue; }
```
where `computeContainedMagnitude()` is a free function that returns the length/area/volume
of the part of a geometry that is contained inside of another geometry.
Finally, after the desired number of entities has been created, we cast the entities
into an instance of a `ContainedEntityNetwork`, using the corresponding builder class:
```cpp
ContainedEntityNetworkBuilder builder;
// define the domain (single sub-domain) and give it a unique id
builder.addConfiningSubDomain(domain, Id(1));
```
The class `ContainedEntityNetworkBuilder` returns an instance of a `ContainedEntityNetwork`
when the `build()` function is called (see below). This network implementation contains information
on (sub-)domains and which entities are embedded in which (sub-)domain. Each (sub-)domain
receives a unique identifier by means of the `Id` class. By calling
`builder.addConfiningSubDomain(domain, Id(1));`, we define the (sub-)domain to be
confining, i.e. entities that are added to this (sub-)domain will be confined to it
by cutting away all parts that lie outside the (sub-)domain. In contrast to that,
one could call `builder.addSubDomain(domain, Id(1));`, in which case embedded networks
will not be confined (see [example 3][2]). Entities are associated with the (sub-)
domain they are embedded in, and are added to the builder class by writing
```cpp
// define entities to be embedded in this domain
builder.addSubDomainEntities(entitySet1, Id(1));
builder.addSubDomainEntities(entitySet2, Id(1));
const auto network = builder.build();
```
The variable `network` now holds an instance of the `ContainedEntityNetwork`,
which we can also pass to the `GmshWriter`:
```cpp
GmshWriter writer(network);
writer.write("network", // filename of the .geo files (will add extension .geo automatically)
0.1, // element size to be used on the quadrilaterals
0.2); // element size to be used in the domain away from the quadrilaterals
```
As you can see, one can specify different mesh sizes to be used on the fracture
entities and in the rest of the domain.
[0]: https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/tree/master/appl/example1
[1]: https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/tree/master/geometry/cylinder.hh
[2]: https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/tree/master/appl/example3
<!--- Example picture --->
<p align="center">
<img src="../../doc/img/example3_network.png" alt="frackit example 1" width="800"/>
</p>
Example 3
=========
In contrast to the [previous example][0], we now want to create a network within
a complex domain that cannot be represented by the internal geometry classes.
The domain has been created with [Gmsh][1] and has been saved in .brep file format.
In a first step, we want to read the file and parse it into an instance of the
`TopoDS_Shape` class of [OpenCascade][2]. To this end, we use the utility function
provided by Frackit in the namespace OCCUtilities:
```cpp
const auto domainShape = OCCUtilities::readShape(BREPFILE);
```
BREPFILE is a preprocessor variable that is substituted by cmake with the actual
path to the .brep file used for this test (see the file `CMakeLists.txt` in this folder).
The domain specified in the .brep file consists of three layers (see above image),
and we want to construct an entity network only in the center one. The variable
`domainShape` contains the information of all layers, and we can use the utility
function
```cpp
const auto solids = OCCUtilities::getSolids(domainShape);
```
to extract the three layers, represented by instances of the `TopoDS_Solid` class
of [OpenCascade][2]. Note that the function `getSolids()` returns an instance of
`std::vector<TopoDS_Solid>`, and thus, we obtain the middle layer with
```cpp
const auto& networkDomain = solids[1];
```
Note that this requires knowledge about the ordering of the solids in the .brep file.
As in the previous examples, we need a point sampler that samples the points to be
used as the center points of the fracture entities. In this example, we again want to
uniformly sample the points within a box, for which we use the bounding box of
`networkDomain`. This can be obtained using the utility function:
```cpp
const auto domainBBox = OCCUtilities::getBoundingBox(networkDomain);
```
With this, we create a quadrilateral sampler in a similar way as in the previous examples
```cpp
QuadrilateralSampler<3> quadSampler(makeUniformPointSampler(domainBBox), // sampler for quadrilateral center points
Distro(toRadians(45.0), toRadians(5.0)), // strike angle: mean value & standard deviation
Distro(toRadians(90.0), toRadians(5.0)), // dip angle: mean value & standard deviation
Distro(45.0, 5.0), // edge length: mean value & standard deviation
5.0); // threshold for minimum edge length
```
For the entities of the other orientation, we want to use `Disk` objects in this case.
The instantiation of the corrensponding sampler class looks like this:
```cpp
DiskSampler diskSampler(makeUniformPointSampler(domainBBox), // sampler for disk center points
Distro(30.0, 6.5), // major axis length: mean value & standard deviation
Distro(24.0, 4.5), // minor axis length: mean value & standard deviation
Distro(toRadians(0.0), toRadians(7.5)), // rotation around x-axis: mean value & standard deviation
Distro(toRadians(0.0), toRadians(7.5)), // rotation around y-axis: mean value & standard deviation
Distro(toRadians(0.0), toRadians(7.5))); // rotation around z-axis: mean value & standard deviation
```
Here, the second and third constructor arguments specify the distributions used
for the major and minor axis lengths of the elliptical disks. The last three
constructor arguments specify the distributions used to determine the orientation
of the disks (for details we refer to the [class documentation][3]).
We use the `MultiGeometrySampler` class to facilitate sampling from several sampler
class with possibly varying geometry types. Arbitrary many samplers can be added to
this class, with the following syntax:
```cpp
// Define ids for the two entity sets
const Id diskSetId(1); // we give the set of orientation one, consisting of disks, the id 1
const Id quadSetId(2); // we give the set of orientation two, consisting of quadrilaterals, the id 2
using Disk = Disk<ctype>;
using Quad = Quadrilateral<ctype, 3>;
// Sampler that samples both geometry types (disks & quads)
// In this one can define an arbitrary number of samplers, each
// of which is associated with an entity set with a unique id
MultiGeometrySampler<Disk, Quad> multiSampler;
multiSampler.addGeometrySampler(diskSampler, diskSetId);
multiSampler.addGeometrySampler(quadSampler, quadSetId);
```
Each sampler is associated with a unique identifier and an error is thrown if it
is tried to add a sampler with an identifier that is already taken. Sampling
occurs again by using the `()` operator, but, in this case it receives an instance
of `Id`, in which it will store the identifier of the sampler from which the
geometry has been generated. That is, the code
```cpp
Id id;
auto geom = multiSampler(id);
```
stores in `id` the identifier of the sampler from which `geom` was sampled.
The variable `geom` holds an instance of an abstract geometry class as the
return type of the `()` operator of the `MultiGeometrySampler` must be uniquely
defined. Hoever, instances of the abstract geometry class can be cast back into
the actual geometry (in this cas `Disk` or `Quad`).
In this context, another useful class is the `MultiGeometryEntitySet`. It can
store arbitrarily many sets of entities of different types. These
types have to be provided as template arguments to the class:
```cpp
MultiGeometryEntitySet<Disk, Quad> entitySets;
```
Entities are added passing a unique identifier, that is, the code
```cpp
Disk disk = diskSampler();
entitySets.addEntity(disk, diskSetId);
```
adds a disk to an entity set with the identifier `diskSetId`. A respective overload
for this function is also implemented for the abstract geometry class mentioned above.
Thus, the geometries sampled from the `MultiGeometrySampler` can directly be inserted
into the instance of `MultiGeometryEntitySet`:
```cpp
Id id;
auto geom = multiSampler(id);
entitySets.addEntity(geom, id);
```
This allows the generation and storage of multiple geometry types of multiple orientations
in a compact way.
In the previous examples, the constraints, for a new entity candidate, were evaluated
manually against the previously admitted entities obtained from the different samplers.
This becomes increasingly cumbersome the more directions are to be considered, especially
if individual constraints are chosen between entities obtained from different samplers.
In such cases, the `EntityNetworkConstraintsMatrix` class can be used to facilitate
the constraints evaluation. It can be used in conjunction with the `MultiGeometryEntitySet`
class and allows for the definition of an arbitrary number of constraints,
which are to be evaluated among entities of specific orientations. In this example, we
create three different instances of the `EntityNetworkConstraints` class, and add them
to the matrix with:
```cpp
EntityNetworkConstraintsMatrix<Constraints> constraintsMatrix;
constraintsMatrix.addConstraints(constraints1, // constraint instance
IdPair(diskSetId, diskSetId)); // set between which to use these constraints
constraintsMatrix.addConstraints(constraints2, // constraint instance
IdPair(quadSetId, quadSetId)); // set between which to use these constraints
constraintsMatrix.addConstraints(constraintsOnOther, // constraint instance
{IdPair(diskSetId, quadSetId), // sets between which to use these constraints
IdPair(quadSetId, diskSetId)}); // sets between which to use these constraints
```
The `IdPair` objects are used to tell the constraints matrix between entities of
which sets, identified by ids, the individual constraints must hold. The evaluation
of all constraints can then be written in one line:
```cpp
// enforce constraints w.r.t. to the other entities
if (!constraintsMatrix.evaluate(entitySets, geom, id))
{ status.increaseRejectedCounter(); continue; }
```
In this function call, `id` holds the identifier of the set for which the candidate
`geom` was sampled. Internally, the `EntityNetworkConstraintsMatrix` instance will
check `geom` against the other entity sets using all constraints that were defined for this identifier. For instance, in this example a disk-shaped geometry sampled from
the sampler with id `diskSetId`, will be checked against all entities of the set with
id `diskSetId` using `constraints1` and against all entities of the set with id `quadSetId`
using `constraintsOnOther`.
After successful generation of the desired number of entities, we again want to
construct an entity network from the raw entities and write it out in [Gmsh][1]
file format. In a first step, we again create an instance of the network builder
class for contained networks and define the sub-domains:
```cpp
ContainedEntityNetworkBuilder builder;
// add sub-domains
builder.addConfiningSubDomain(solids[0], Id(1));
builder.addConfiningSubDomain(networkDomain, Id(2));
builder.addConfiningSubDomain(solids[2], Id(3));
```
We have given the identifier with index 2 to the (sub-)domain for which we have created
the entity network. Using the convenience function available in the `MultiGeometryEntitySet`
class, we can add all entities to the builder by writing
```cpp
ContainedEntityNetworkBuilder builder;
entitySets.exportEntitySets(builder, Id(2));
```
where the second argument states the (sub-)domain for which the entities should
be added to the builder (see [example 2][4]). If `entitySets` holds entity sets
associated with multiple (sub-)domains, or, if only some of the entity sets should
be added to the builder, one can use:
```cpp
ContainedEntityNetworkBuilder builder;
entitySets.exportEntitySets({diskSetId, quadSetId}, builder, Id(2));
```
With this call only the entities of those sets that are given in the list of ids
inside the brackets (`{}`) are exported. An entity network can then be built and
written out in the same manner as shown in the other examples. At the end of the
main file `example3.cc`, it is shown how different types of networks
(contained-confined, contained-unconfined, etc.) can be created from the raw entities.
[0]: https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/tree/master/appl/example2
[1]: http://gmsh.info/
[2]: https://www.opencascade.com/content/download-center
[3]: https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/tree/master/sampling/disksampler.hh
[4]: https://git.iws.uni-stuttgart.de/DennisGlaeser/frackit/tree/master/appl/example2
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment