|
Medusa
1.1
Coordinate Free Mehless Method implementation
|
|
Go to the documentation of this file. 1 #ifndef MEDUSA_BITS_DOMAINS_GRAINDROPFILL_HPP_
2 #define MEDUSA_BITS_DOMAINS_GRAINDROPFILL_HPP_
17 template <
typename vec_t>
19 bot(bot), top(top), seed_(
get_seed()) {
20 for (
int i = 0; i <
dim; ++i) {
21 assert_msg(
top[i] >
bot[i],
"Bottom must be less than top, got bot = %s, top = %s.",
26 template <
typename vec_t>
27 template <
typename func_t>
30 if (type == 0) type = 1;
31 std::vector<vec_t> nodes = fillBox(h);
33 for (
const vec_t& p : nodes) {
34 if (!domain.
contains(p, boundary_tree))
continue;
37 if (r2 >= allowed_r*allowed_r) {
43 template <
typename vec_t>
44 template <
typename func_t>
46 assert_msg(dx > 0,
"Initial spacing must be positive, got %g. Did you set it correctly?", dx);
48 std::vector<vec_t> result;
52 for (
int i = 0; i <
dim-1; ++i) {
53 grid_size[i] =
iceil((top[i] - bot[i]) / dx) + 1;
54 assert_msg(grid_size[i] >= 1,
"Grid size in dimension %d is not positive.");
58 Eigen::Map<Eigen::Matrix<scalar_t, Eigen::Dynamic, 1>>
59 eig_z(z_coor.data().data(), z_coor.size());
60 eig_z += 1e-4*dx*Eigen::Matrix<scalar_t, Eigen::Dynamic, 1>::Random(z_coor.size());
64 scalar_t min_z = eig_z.minCoeff(&min_lin_idx);
65 Vec<int,
dim-1> min_idx = z_coor.multiIndex(min_lin_idx);
68 while (dot_nr < max_points && min_z <= top[
dim-1] + excess_factor*dx) {
70 for (
int i = 0; i <
dim-1; ++i) p[i] = bot[i] + min_idx[i] * dx;
77 int idx_rad =
ifloor(r / dx);
83 scalar_t height_update = r*r-((cnt-min_idx).
template cast<scalar_t>()*dx).squaredNorm();
84 cnt_lin = z_coor.linearIndex(cnt);
85 if (height_update > 0) {
86 z_coor[cnt_lin] = std::max(z_coor[cnt_lin], min_z + std::sqrt(height_update));
89 min_z = z_coor[min_lin_idx];
93 while (dist > idx_rad) {
94 for (
int i = 0; i <
dim-1; ++i) {
95 low[i] = (((min_idx[i] - 2*idx_rad) % grid_size[i]) + grid_size[i]) % grid_size[i];
96 high[i] = (min_idx[i] + 2*idx_rad + 1) % grid_size[i];
98 auto old_min_idx = min_idx;
101 cnt_lin = z_coor.linearIndex(cnt);
102 if (z_coor[cnt_lin] < min_z) {
103 min_z = z_coor[cnt_lin];
104 min_lin_idx = cnt_lin;
108 dist = (min_idx - old_min_idx).cwiseAbs().maxCoeff();
116 template <
typename func_t>
119 result.resize(std::min(result.size(), max_points));
125 #endif // MEDUSA_BITS_DOMAINS_GRAINDROPFILL_HPP_
Root namespace for the whole library.
Scalar scalar_t
Type of the elements, alias of Scalar.
bool incrementCyclicCounter(Vec< int, dim > &counter, const Vec< int, dim > &low, const Vec< int, dim > &high, const Vec< int, dim > &size)
Increments a multi-dimensional counter with given upper and lower limits and global upper size,...
Class representing domain discretization along with an associated shape.
Eigen::Matrix< scalar_t, dim, 1, Eigen::ColMajor|Eigen::AutoAlign, dim, 1 > Vec
Fixed size vector type, representing a mathematical 1d/2d/3d vector.
Class representing a static k-d tree data structure.
vec_t bot
Bottom left corner of the box to fill.
@ dim
Number of elements of this matrix.
std::vector< vec_t > fillBox(const func_t &h) const
Actual fill algorithm which fills [bot, top] box with spacing h.
const Range< vec_t > & positions() const
Returns positions of all nodes.
#define assert_msg(cond,...)
Assert with better error reporting.
@ dim
Dimensionality of the domain.
std::pair< Range< int >, Range< scalar_t > > query(const vec_t &point, int k=1) const
Find k nearest neighbors to given point.
int ifloor(T x)
Floors a floating point to an integer.
vec_t::scalar_t scalar_t
Scalar type;.
int iceil(T x)
Ceils a floating point to an integer.
Class representing a simple n-dimensional grid structure, which supports indexing and storing values.
void operator()(DomainDiscretization< vec_t > &domain, const func_t &h, int type=0) const
Fills domain with a quality node distribution.
bool contains(const vec_t &point) const
Returns true if point is inside the domain.
unsigned int get_seed()
Return a random seed.
Range< vec_t > discretizeLineWithDensity(const vec_t &p, const vec_t &q, const func_t &delta_r)
Returns nodes lying on the line segment pq with approximate distances delta_r.
int addInternalNode(const vec_t &point, int type)
Adds a single interior node with specified type to this discretization.
vec_t top
Top right corner of the box to fill.
bool incrementCounter(Vec< int, dim > &counter, const Vec< int, dim > &limit)
Increments a multi-dimensional counter with given limits.
GrainDropFill(const vec_t &bot, const vec_t &top)
Prepare to fill domain within [bot, top] bounding box.