Commit 5f97d9d2 authored by Dennis Gläser's avatar Dennis Gläser
Browse files

[test][disknetwork] add test for generation on shape

TODO: Implement intersection algorithms with shapes such that
constraints on the shape boundaries can be enforced
parent 61cd016b
......@@ -17,4 +17,11 @@ frackit_add_test(NAME test_generate_disk_network_cylinder
SOURCES test_generate_disk_network_cylinder.cc
LABELS entitynetwork sampling)
# test disk network creation (in shape from file)
frackit_add_test(NAME test_generate_disk_network_shape
SOURCES test_generate_disk_network_shape.cc
LABELS entitynetwork sampling
COMPILE_DEFINITIONS BREPFILE="${CMAKE_SOURCE_DIR}/test/entitynetwork/box.brep"
COMMAND ./test_generate_disk_network_shape)
set(CMAKE_BUILD_TYPE Debug)
DBRep_DrawableShape
CASCADE Topology V1, (c) Matra-Datavision
Locations 0
Curve2ds 24
1 0 0 1 0
1 0 0 1 0
1 1 0 0 -1
1 0 0 0 1
1 0 -1 1 0
1 0 0 1 0
1 0 0 0 -1
1 0 0 0 1
1 0 0 1 0
1 0 1 1 0
1 1 0 0 -1
1 1 0 0 1
1 0 -1 1 0
1 0 1 1 0
1 0 0 0 -1
1 1 0 0 1
1 0 0 0 1
1 0 0 1 0
1 1 0 0 1
1 0 0 1 0
1 0 0 0 1
1 0 1 1 0
1 1 0 0 1
1 0 1 1 0
Curves 12
1 0 0 0 0 0 1
1 0 0 1 -0 1 0
1 0 1 0 0 0 1
1 0 0 0 -0 1 0
1 1 0 0 0 0 1
1 1 0 1 -0 1 0
1 1 1 0 0 0 1
1 1 0 0 -0 1 0
1 0 0 0 1 0 -0
1 0 0 1 1 0 -0
1 0 1 0 1 0 -0
1 0 1 1 1 0 -0
Polygon3D 0
PolygonOnTriangulations 0
Surfaces 6
1 0 0 0 1 0 -0 0 0 1 0 -1 0
1 0 0 0 -0 1 0 0 0 1 1 0 -0
1 0 0 1 0 0 1 1 0 -0 -0 1 0
1 0 1 0 -0 1 0 0 0 1 1 0 -0
1 0 0 0 0 0 1 1 0 -0 -0 1 0
1 1 0 0 1 0 -0 0 0 1 0 -1 0
Triangulations 0
TShapes 35
Ve
1e-07
0 0 1
0 0
0101101
*
Ve
1e-07
0 0 0
0 0
0101101
*
Ve
1e-07
0 1 1
0 0
0101101
*
Ve
1e-07
0 1 0
0 0
0101101
*
Ve
1e-07
1 0 1
0 0
0101101
*
Ve
1e-07
1 0 0
0 0
0101101
*
Ve
1e-07
1 1 1
0 0
0101101
*
Ve
1e-07
1 1 0
0 0
0101101
*
Ed
1e-07 1 1 0
1 1 0 0 1
2 1 1 0 0 1
2 2 2 0 0 1
0
0101000
-35 0 +34 0 *
Ed
1e-07 1 1 0
1 2 0 0 1
2 3 1 0 0 1
2 4 3 0 0 1
0
0101000
-33 0 +35 0 *
Ed
1e-07 1 1 0
1 3 0 0 1
2 5 1 0 0 1
2 6 4 0 0 1
0
0101000
-33 0 +32 0 *
Ed
1e-07 1 1 0
1 4 0 0 1
2 7 1 0 0 1
2 8 5 0 0 1
0
0101000
-32 0 +34 0 *
Ed
1e-07 1 1 0
1 5 0 0 1
2 9 6 0 0 1
2 10 2 0 0 1
0
0101000
-31 0 +30 0 *
Ed
1e-07 1 1 0
1 6 0 0 1
2 11 6 0 0 1
2 12 3 0 0 1
0
0101000
-29 0 +31 0 *
Ed
1e-07 1 1 0
1 7 0 0 1
2 13 6 0 0 1
2 14 4 0 0 1
0
0101000
-29 0 +28 0 *
Ed
1e-07 1 1 0
1 8 0 0 1
2 15 6 0 0 1
2 16 5 0 0 1
0
0101000
-28 0 +30 0 *
Ed
1e-07 1 1 0
1 9 0 0 1
2 17 2 0 0 1
2 18 5 0 0 1
0
0101000
-30 0 +34 0 *
Ed
1e-07 1 1 0
1 10 0 0 1
2 19 2 0 0 1
2 20 3 0 0 1
0
0101000
-31 0 +35 0 *
Ed
1e-07 1 1 0
1 11 0 0 1
2 21 4 0 0 1
2 22 5 0 0 1
0
0101000
-28 0 +32 0 *
Ed
1e-07 1 1 0
1 12 0 0 1
2 23 4 0 0 1
2 24 3 0 0 1
0
0101000
-29 0 +33 0 *
Wi
0101100
-27 0 -26 0 +25 0 +24 0 *
Wi
0101100
-23 0 -22 0 +21 0 +20 0 *
Wi
0101100
-19 0 -23 0 +18 0 +27 0 *
Wi
0101100
-17 0 -21 0 +16 0 +25 0 *
Wi
0101100
-24 0 -17 0 +20 0 +19 0 *
Wi
0101100
-26 0 -16 0 +22 0 +18 0 *
Fa
0 1e-07 1 0
0111000
+15 0 *
Fa
0 1e-07 6 0
0111000
+14 0 *
Fa
0 1e-07 2 0
0111000
+13 0 *
Fa
0 1e-07 4 0
0111000
+12 0 *
Fa
0 1e-07 5 0
0111000
+11 0 *
Fa
0 1e-07 3 0
0111000
+10 0 *
Sh
0101100
-9 0 +8 0 -7 0 +6 0 -5 0 +4 0 *
So
0100000
+3 0 *
Co
1100000
-35 0 +34 0 -33 0 -32 0 +31 0 -30 0 +29 0 +28 0 +27 0 +26 0
-25 0 -24 0 -23 0 -22 0 +21 0 +20 0 +19 0 -18 0 -17 0 +16 0
-15 0 +14 0 -13 0 +12 0 -11 0 +10 0 -9 0 +8 0 -7 0 +6 0
-5 0 +4 0 +3 0 +2 0 *
+1 0
\ No newline at end of file
#include <iostream>
#include <random>
#include <fstream>
#include <BRepTools.hxx>
#include <BRep_Builder.hxx>
#include <TopoDS_Shape.hxx>
#include <TopoDS_Compound.hxx>
#include <TopoDS_Solid.hxx>
#include <frackit/geometry/box.hh>
#include <frackit/geometry/disk.hh>
#include <frackit/precision/precision.hh>
#include <frackit/magnitude/magnitude.hh>
#include <frackit/magnitude/containedmagnitude.hh>
#include <frackit/occ/breputilities.hh>
#include <frackit/sampling/pointsampling.hh>
#include <frackit/sampling/geometrysampling.hh>
#include <frackit/entitynetwork/constraints.hh>
//! create a disk network embedded in a cylinder
int main(int argc, char** argv)
{
//! print welcome message
std::cout << "\n\n"
<< "##########################################\n"
<< "## Creating a confined network of disks ##\n"
<< "##########################################\n"
<< "\n\n";
// brep file name is specified in CMakeLists.txt
std::ifstream domainShapeFile(BREPFILE);
if (!domainShapeFile)
throw std::runtime_error(std::string("Could not open shape file"));
TopoDS_Shape domainShape;
BRepTools::Read(domainShape, domainShapeFile, BRep_Builder());
// obtain the (single) solid contained in the file & get its boundaries
using namespace Frackit;
const auto solids = OCCUtilities::getSolids(domainShape);
assert(solids.size() == 1);
const auto& domain = solids[0];
const auto domainFaces = OCCUtilities::getFaces(domain);
const auto domainVolume = computeMagnitude(domain);
// create the disk samplers
using ctype = double;
using Disk = Disk<ctype>;
using DiskSampler = GeometrySampler<Disk>;
using PointSampler = GeometryPointSampler< Box<ctype> >;
// sample points within bounding box of domain
PointSampler pointSampler(OCCUtilities::getBoundingBox(domain));
// sampler for disks of orientation 1
DiskSampler diskSampler_1(std::normal_distribution<ctype>(0.35, 0.1),
std::normal_distribution<ctype>(0.225, 0.05),
std::normal_distribution<ctype>(toRadians(25.0), toRadians(5.0)),
std::normal_distribution<ctype>(toRadians(0.0), toRadians(5.0)),
std::normal_distribution<ctype>(toRadians(45.0), toRadians(5.0)));
// sampler for disks of orientation 1
DiskSampler diskSampler_2(std::normal_distribution<ctype>(0.35, 0.1),
std::normal_distribution<ctype>(0.225, 0.05),
std::normal_distribution<ctype>(toRadians(-35.0), toRadians(5.0)),
std::normal_distribution<ctype>(toRadians(0.0), toRadians(5.0)),
std::normal_distribution<ctype>(toRadians(-45.0), toRadians(5.0)));
//! containers to store created entities
std::vector<Disk> diskSet1;
std::vector<Disk> diskSet2;
//! enforce some constraints on the network
Frackit::EntityNetworkConstraints<ctype> constraintsOnSelf;
Frackit::EntityNetworkConstraints<ctype> constraintsOnOther;
Frackit::EntityNetworkConstraints<ctype> constraintsOnDomain;
// constraints among disks of the same set
constraintsOnSelf.setMinDistance(0.1);
constraintsOnSelf.setMinIntersectingAngle(toRadians(45.0));
constraintsOnSelf.setMinIntersectionMagnitude(0.01);
constraintsOnSelf.setMinIntersectionDistance(0.05);
// constraints with the other disk set
constraintsOnOther.setMinDistance(0.01);
constraintsOnOther.setMinIntersectingAngle(toRadians(15.0));
constraintsOnOther.setMinIntersectionMagnitude(0.01);
constraintsOnOther.setMinIntersectionDistance(0.05);
// constraints with the other disk set
constraintsOnDomain.setMinDistance(0.01);
constraintsOnDomain.setMinIntersectingAngle(toRadians(15.0));
constraintsOnDomain.setMinIntersectionMagnitude(0.01);
constraintsOnDomain.setMinIntersectionDistance(0.05);
//! create pseudo-random number generator to select set during creation
std::default_random_engine randomNumber{std::random_device{}()};
//! keep track of number of created & accepted entities
std::size_t total = 0;
std::size_t accepted = 0;
std::size_t accepted_1 = 0;
std::size_t accepted_2 = 0;
//! print header of status output
std::cout << "##############################################################################################################################\n"
<< "# No. of entities | No. rejected entities | Acceptance Ratio | Current entity density [m²/m³] | Progress [%] #\n"
<< "#----------------------------------------------------------------------------------------------------------------------------#\n";
ctype currentDensity = 0.0;
ctype currentDiskArea = 0.0;
const std::size_t numTargetEntities_1 = 5;
const std::size_t numTargetEntities_2 = 5;
while (accepted_1 != numTargetEntities_1 || accepted_2 != numTargetEntities_2)
{
bool createSecondary = randomNumber()%2;
if (!createSecondary && accepted_1 == numTargetEntities_1)
createSecondary = true;
else if (createSecondary && accepted_2 == numTargetEntities_2)
createSecondary = false;
auto disk = createSecondary ? diskSampler_1(pointSampler)
: diskSampler_2(pointSampler);
total++;
// We don't want ellipses of too large aspect ratio
if (disk.majorAxisLength() > disk.minorAxisLength()*3.0) continue;
auto& diskSetSelf = createSecondary ? diskSet2 : diskSet1;
const auto& diskSetOther = createSecondary ? diskSet1 : diskSet2;
// enforce constraints w.r.t. to the own disk set
if (!constraintsOnSelf.evaluate(diskSetSelf, disk)) continue;
// enforce constraints w.r.t. to the other disk set
if (!constraintsOnOther.evaluate(diskSetOther, disk)) continue;
// // enforce constraints w.r.t. the domain boundaries
// TODO: INTERSECTION ALGOS WITH SHAPE CLASSES
// bool violates = false;
// for (const auto& face : domainFaces)
// if (!constraintsOnDomain.evaluate(face, disk))
// { violates = true; break; }
// if (violates) continue;
// reject if intersection with domain is too small (here: 0.2m²)
const auto containedArea = computeContainedMagnitude(disk, domain);
if (containedArea < 0.1) continue;
// if we get here, the disk is admissible
diskSetSelf.push_back(std::move(disk));
accepted++;
if (createSecondary) accepted_2++;
else accepted_1++;
// compute new density (use minimum w.r.t. domain/network volume)
currentDiskArea += containedArea;
currentDensity = currentDiskArea/domainVolume;
std::cout << " " << std::setw(18) << std::setfill(' ')
<< std::to_string(diskSet1.size() + diskSet2.size()) + std::string(9, ' ')
<< "| " << std::setprecision(6) << std::setw(24) << std::setfill(' ')
<< std::to_string(total-accepted) + std::string(12, ' ')
<< "| " << std::setprecision(6) << std::setw(19) << std::setfill(' ')
<< std::to_string( double(double(accepted)/double(total)) ) + std::string(7, ' ')
<< "| " << std::setprecision(6) << std::setw(33) << std::setfill(' ')
<< std::to_string(currentDensity) + std::string(17, ' ')
<< "| " << std::string(3, ' ') + std::to_string( double(accepted_1+accepted_2)/double(numTargetEntities_1+numTargetEntities_2) ) << std::endl;
}
// create a single compound and write to .brep file
BRep_Builder b;
TopoDS_Compound c;
b.MakeCompound(c);
for (const auto& disk : diskSet1) b.Add(c, OCCUtilities::getShape(disk));
for (const auto& disk : diskSet2) b.Add(c, OCCUtilities::getShape(disk));
b.Add(c, domain);
BRepTools::Write(c, "disknetwork.brep");
// fragment the network and write in separate file
std::vector<TopoDS_Shape> allDisks;
for (const auto& disk : diskSet1) allDisks.push_back(OCCUtilities::getShape(disk));
for (const auto& disk : diskSet2) allDisks.push_back(OCCUtilities::getShape(disk));
const auto fragmentedNetwork = OCCUtilities::fragment(allDisks, Frackit::Precision<ctype>::confusion());
BRepTools::Write(fragmentedNetwork, "disknetwork_fragmented.brep");
const auto confinedNetwork = OCCUtilities::intersect(fragmentedNetwork,
domain,
Frackit::Precision<ctype>::confusion());
BRepTools::Write(confinedNetwork, "disknetwork_confined.brep");
std::vector<TopoDS_Shape> allShapes({confinedNetwork, domain});
BRepTools::Write(OCCUtilities::fragment(allShapes, Frackit::Precision<ctype>::confusion()),
"disknetwork_medium.brep");
std::cout << "All tests passed" << std::endl;
return 0;
}
Markdown is supported
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