Medusa  1.1
Coordinate Free Mehless Method implementation
test/domains/BasicRelax_test.cpp
#include "gtest/gtest.h"
namespace mm {
TEST(DomainEngines, RelaxWithConstDistributionTestBoundaryProjection0) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
double step = 1.0 / 55;
DomainDiscretization<Vec2d> c = sh.discretizeWithStep(step);
// no of boundary nodes before relax
int N_1 = c.boundary().size();
double min_spacing = 0.5 * step;
BasicRelax relax;
relax.iterations(20).initialHeat(5).numNeighbours(4).projectionType(BasicRelax::DO_NOT_PROJECT);
relax(c, min_spacing);
// no of boundary nodes after relax
int N_2 = c.boundary().size();
// actual tests
c.findSupport(FindClosest(2));
std::vector<double> tmp;
for (int i = 0; i < c.size(); ++i) tmp.push_back((c.pos(i) - c.pos(c.support(i, 1))).norm());
// check for overlapping nodes
EXPECT_GT(*std::min_element(std::begin(tmp), std::end(tmp)), min_spacing / 10);
// check if nodes are added on boundary
EXPECT_EQ(N_2, N_1);
}
TEST(DomainEngines, RelaxWithConstDistributionTestBoundaryProjection1) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
double step = 1.0 / 55.0;
DomainDiscretization<Vec2d> c = sh.discretizeWithStep(step);
Range<int> to_remove;
int j = 0;
for (int i : c.boundary()) {
if (j++ % 3 != 0) to_remove.push_back(i);
}
c.removeNodes(to_remove); // keep every fifth node
// no of boundary nodes before relax
int N_1 = c.boundary().size();
double min_spacing = 0.5 * step;
BasicRelax relax;
relax.iterations(50).initialHeat(1).numNeighbours(4).finalHeat(0.1)
.projectionType(BasicRelax::PROJECT_IN_DIRECTION).boundaryProjectionThreshold(0.35);
relax(c, min_spacing);
// no of boundary nodes after relax
int N_2 = c.boundary().size();
// actual tests
c.findSupport(FindClosest(2));
std::vector<double> tmp;
for (int i = 0; i < c.size(); ++i) tmp.push_back((c.pos(i) - c.pos(c.support(i, 1))).norm());
// check for overlapping nodes
EXPECT_GT(*std::min_element(std::begin(tmp), std::end(tmp)), min_spacing / 10);
// check if nodes are added on boundary
EXPECT_GT(N_2, N_1);
}
TEST(DomainEngines, RelaxWithConstDistributionTestBoundaryProjection2) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
double step = 1.0 / 55.0;
DomainDiscretization<Vec2d> c = sh.discretizeWithStep(step);
Range<int> to_remove;
int j = 0;
for (int i : c.boundary()) {
if (j++ % 3 != 0) to_remove.push_back(i);
}
c.removeNodes(to_remove); // keep every fifth node
// no of boundary nodes before relax
int N_1 = c.boundary().size();
double min_spacing = 0.5 * step;
BasicRelax relax;
relax.iterations(50).initialHeat(1).numNeighbours(5)
relax(c, min_spacing);
// no of boundary nodes after relax
int N_2 = c.boundary().size();
// actual tests
c.findSupport(FindClosest(2));
std::vector<double> tmp;
for (int i = 0; i < c.size(); ++i) tmp.push_back((c.pos(i) - c.pos(c.support(i, 1))).norm());
// check for overlapping nodes
EXPECT_GT(*std::min_element(std::begin(tmp), std::end(tmp)), min_spacing / 10);
// check if nodes are added on boundary
EXPECT_GT(N_2, N_1);
}
TEST(DomainEngines, RelaxWithConstDistributionTestBoundaryProjectiontreshold) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
double step = 1.0 / 55.0;
DomainDiscretization<Vec2d> c = sh.discretizeWithStep(step);
// no of boundary nodes before relax
int N_1 = c.boundary().size();
double min_spacing = 0.5 * step;
BasicRelax relax;
relax.iterations(20).initialHeat(1).numNeighbours(4).boundaryProjectionThreshold(1.5)
relax(c, min_spacing);
// no of boundary nodes after relax
int N_2 = c.boundary().size();
// actual tests
c.findSupport(FindClosest(2));
std::vector<double> tmp;
for (int i = 0; i < c.size(); ++i) tmp.push_back((c.pos(i) - c.pos(c.support(i, 1))).norm());
// check for overlapping nodes
EXPECT_GT(*std::min_element(std::begin(tmp), std::end(tmp)), min_spacing / 10);
// check if nodes are added on boundary
EXPECT_EQ(N_2, N_1);
}
TEST(DomainEngines, RelaxWithVariableDistribution) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
auto fill_density = [](const Vec2d& p) {
return (0.005 + (p[0] - 0.5) * (p[0] - 0.5) / 2 + (p[1] - 0.5) * (p[1] - 0.5) / 2);
};
DomainDiscretization<Vec2d> domain = sh.discretizeBoundaryWithStep(fill_density({r, 0.0}));
GeneralFill<Vec2d> fill_engine; fill_engine.seed(1);
BasicRelax relax_engine;
relax_engine.iterations(50).initialHeat(1).finalHeat(0.05).boundaryProjectionThreshold(0.75)
.projectionType(BasicRelax::PROJECT_IN_DIRECTION).numNeighbours(5);
domain.fill(fill_engine, fill_density);
domain.relax(relax_engine, fill_density);
domain.findSupport(FindClosest(2));
std::vector<double> tmp;
for (int i = 0; i < domain.size(); ++i) tmp.push_back(domain.dr(i));
// check for overlapping nodes
EXPECT_GT(*std::min_element(std::begin(tmp), std::end(tmp)), 0.001);
}
TEST(DomainEngines, RelaxDeathTest) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
BasicRelax relax;
relax.iterations(50).initialHeat(20).finalHeat(0.05).boundaryProjectionThreshold(0.75)
.projectionType(BasicRelax::PROJECT_IN_DIRECTION).numNeighbours(5);
auto fill_density = [](Vec2d p) -> double {
return (0.05 + (p[0] - 0.5) * (p[0] - 0.5) / 2 + (p[1] - 0.5) * (p[1] - 0.5) / 2);
};
double step = 1.0 / 15.0;
DomainDiscretization<Vec2d> c = sh.discretizeWithStep(step);
// check for overlapping nodes
EXPECT_DEATH(relax(c, fill_density), "No nodes in relax pool anymore, perhaps use lower heat");
}
TEST(DomainEngines, RelaxWithoutDistribution) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
BasicRelax relax;
relax.iterations(100).initialHeat(1).finalHeat(0).numNeighbours(3)
.projectionType(BasicRelax::PROJECT_IN_DIRECTION).boundaryProjectionThreshold(0.55);
auto fill_density = [](Vec2d p) -> double {
return (0.005 + (p[0] - 0.5) * (p[0] - 0.5) / 2 + (p[1] - 0.5) * (p[1] - 0.5) / 2);
};
DomainDiscretization<Vec2d> c = sh.discretizeBoundaryWithStep(fill_density(Vec2d({r, 0.0})));
GeneralFill<Vec2d> fill_engine;
fill_engine.seed(1);
fill_engine(c, fill_density);
relax(c);
c.findSupport(FindClosest(2));
std::vector<double> tmp;
for (int i = 0; i < c.size(); ++i) tmp.push_back((c.pos(i) - c.pos(c.support(i, 1))).norm());
// check for overlapping nodes
EXPECT_GT(*std::min_element(std::begin(tmp), std::end(tmp)), 0.001);
}
TEST(DomainEngines, DISABLED_RelaxWikiImage1) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
auto fill_density = [](Vec2d p) -> double {
return (0.0025 + (p[0] - 0.5) * (p[0] - 0.5) / 2 + (p[1] - 0.5) * (p[1] - 0.5) / 4);
};
prn(fill_density(Vec2d({0.25, 0})))
DomainDiscretization<Vec2d> c = sh.discretizeBoundaryWithStep(
fill_density(Vec2d({0.25, 0})) / 2);
GeneralFill<Vec2d> fill_engine;
fill_engine(c, fill_density);
BasicRelax relax;
relax.iterations(200).initialHeat(2).numNeighbours(3)
relax.boundaryProjectionThreshold(0.45).finalHeat(0.5);
relax(c, fill_density);
c.findSupport(FindClosest(2));
std::vector<double> tmp;
for (int i = 0; i < c.size(); ++i) tmp.push_back((c.pos(i) - c.pos(c.support(i, 1))).norm());
// check for overlapping nodes
EXPECT_GT(*std::min_element(std::begin(tmp), std::end(tmp)), 0.0005);
}
TEST(DomainEngines, DISABLED_RelaxWikiImage2) {
double r = 0.25;
BallShape<Vec2d> sh({0.5, 0.5}, r);
auto fill_density = [](Vec2d p) -> double {
return (0.0025 + (p[0] - 0.5) * (p[0] - 0.5) / 2 + (p[1] - 0.5) * (p[1] - 0.5) / 4);
};
prn(fill_density(Vec2d({0.25, 0})))
DomainDiscretization<Vec2d> c = sh.discretizeBoundaryWithStep(
fill_density(Vec2d({0.25, 0})) / 2);
GeneralFill<Vec2d> fill_engine;
fill_engine(c, fill_density);
BasicRelax relax;
relax.iterations(200).initialHeat(2).numNeighbours(3)
relax.boundaryProjectionThreshold(0.45).finalHeat(0.5);
relax(c);
c.findSupport(FindClosest(2));
std::vector<double> tmp;
for (int i = 0; i < c.size(); ++i) tmp.push_back((c.pos(i) - c.pos(c.support(i, 1))).norm());
// check for overlapping nodes
EXPECT_GT(*std::min_element(std::begin(tmp), std::end(tmp)), 0.0005);
}
} // namespace mm
mm
Root namespace for the whole library.
Definition: Gaussian.hpp:14
mm::BasicRelax::DO_NOT_PROJECT
@ DO_NOT_PROJECT
Escaped nodes are frozen outside the domain till the end of relax and then removed.
Definition: BasicRelax_fwd.hpp:51
mm::BasicRelax::PROJECT_IN_DIRECTION
@ PROJECT_IN_DIRECTION
Project on boundary in relax movement direction.
Definition: BasicRelax_fwd.hpp:56
FindClosest.hpp
prn
#define prn(...)
Prints a variable name and value to standard output.
Definition: assert.hpp:34
BallShape.hpp
mm::BasicRelax::PROJECT_BETWEEN_CLOSEST
@ PROJECT_BETWEEN_CLOSEST
Project between two closest boundary nodes – this one might result in divergent behaviour,...
Definition: BasicRelax_fwd.hpp:63
Vec.hpp
BasicRelax.hpp
mm::Vec2d
Vec< double, 2 > Vec2d
Convenience typedef for 2d vector of doubles.
Definition: Vec_fwd.hpp:34
GeneralFill.hpp