1 #ifndef MEDUSA_BITS_DOMAINS_BALLSHAPE_HPP_
2 #define MEDUSA_BITS_DOMAINS_BALLSHAPE_HPP_
18 template <
typename vec_t>
21 if (type == 0) type = -1;
26 for (
const auto& p : points) {
32 template <
typename vec_t>
34 scalar_t step,
int internal_type,
int boundary_type)
const {
35 auto d = discretizeBoundaryWithStep(step, boundary_type);
36 if (internal_type == 0) internal_type = 1;
38 int num_of_points =
iceil(2*radius_/step) - 1;
39 for (
const vec_t& v :
linspace(bbox().first, bbox().second, {num_of_points},
false)) {
40 d.addInternalNode(v, internal_type);
42 }
else if (
dim == 2) {
44 int num_of_circles =
iceil(radius_ / step);
45 scalar_t dr = radius_ / num_of_circles;
46 d.addInternalNode(center_, internal_type);
47 for (
int i = 1; i < num_of_circles; ++i) {
51 for (
int j = 0; j < num_of_points; ++j) {
52 d.addInternalNode(center_ + r * vec_t({std::cos(j*dfi), std::sin(j*dfi)}),
56 }
else if (
dim >= 3) {
57 assert_msg(
false,
"This domain does not support filling with density, "
58 "use a GeneralSurfaceFill engine instead.");
63 template <
typename vec_t>
65 return os <<
"BallShape(" << center_.transpose() <<
", " << radius_ <<
")";
68 template <
typename vec_t>
72 if (type == 0) type = -1;
77 }
else if (
dim == 2) {
81 while (cur_alpha < max_alpha) {
82 scalar_t r = dr(center_ + radius_ * vec_t({std::cos(cur_alpha), std::sin(cur_alpha)}));
83 cur_alpha += 2*std::asin(r/2.0/radius_);
84 alphas.push_back(cur_alpha);
88 for (
int i = 0; i < alphas.
size()-1; ++i) {
90 vec_t normal = {std::cos(a), std::sin(a)};
93 }
else if (
dim == 3) {
96 double d = (1 + v.squaredNorm());
97 vec_t p = {2*v[0]/d, 2*v[1]/d, (1-v.squaredNorm())/d};
98 return d/2 * dr(center_ + radius_ * p) / radius_;
106 for (
int i : domain2d.interior()) {
107 auto v = domain2d.pos(i);
108 double d = (1 + v.squaredNorm());
109 vec_t p = {2*v[0]/d, 2*v[1]/d, (1-v.squaredNorm())/d};
113 for (
int i : domain2d.boundary()) {
114 auto v = domain2d.pos(i);
115 vec_t p = {v[0], v[1], 0};
121 double d = (1 + v.squaredNorm());
122 vec_t p = {2*v[0]/d, 2*v[1]/d, -(1-v.squaredNorm())/d};
123 return d/2 * dr(center_ + radius_ * p) / radius_;
126 for (
int i : domain2d.interior()) {
127 auto v = domain2d.pos(i);
128 double d = (1 + v.squaredNorm());
129 vec_t p = {2*v[0]/d, 2*v[1]/d, -(1-v.squaredNorm())/d};
133 assert_msg(
false,
"This domain does not support filling with density, "
134 "use a GeneralSurfaceFill engine instead.");
141 #endif // MEDUSA_BITS_DOMAINS_BALLSHAPE_HPP_