finished base for curling stability tests

added comments
This commit is contained in:
PavelMikus 2022-05-11 13:17:01 +02:00
parent 4144b73ccd
commit e39d14bf98

View file

@ -48,7 +48,7 @@ CurledFilament::CurledFilament(const Vec3f &position) :
struct Cell {
float volume;
float curled_height;
int island_id = std::numeric_limits<int>::max();
int accumulator_id = std::numeric_limits<int>::max();
};
class CentroidAccumulator {
@ -262,7 +262,7 @@ public:
return local_coords.x() >= 0 && local_coords.y() >= 0 && local_coords.x() < this->global_cell_count.x()
&& local_coords.y() < this->global_cell_count.y();
};
CentroidAccumulators accumulators(issues.supports_nedded.size() + 4);
CentroidAccumulators accumulators(issues.supports_nedded.size() + 4); // just an estimation; one for each support point from prev step, and 4 for the base
auto custom_comparator = [](const Vec2i &left, const Vec2i &right) {
if (left.x() == right.x()) {
return left.y() < right.y();
@ -270,26 +270,28 @@ public:
return left.x() < right.x();
};
int next_island_id = -1;
std::set<Vec2i, decltype(custom_comparator)> coords_to_check(custom_comparator);
//initialization of the bed accumulators from the bed cells (first layer of grid)
int next_accumulator_id = -1; // accumulators from the bed have negative ids, starting with -1. Accumulators generated by support points have nonegative ids, starting with 0, and sorted by height
// The reason is, that during merging of accumulators (when they meet at the upper cells), algorithm always keeps the one with the lower id (so bed is preffered), and discards the other
std::set<Vec2i, decltype(custom_comparator)> coords_to_check(custom_comparator); // set of coords to check for the current accumulator (search based on connectivity)
for (int y = 0; y < global_cell_count.y(); ++y) {
for (int x = 0; x < global_cell_count.x(); ++x) {
Cell &origin_cell = this->access_cell(Vec3i(x, y, 0));
if (origin_cell.volume > 0 && origin_cell.island_id == std::numeric_limits<int>::max()) {
CentroidAccumulator &acc = accumulators.create_accumulator(next_island_id, 0);
if (origin_cell.volume > 0 && origin_cell.accumulator_id == std::numeric_limits<int>::max()) { // cell has volume and no accumulator assigned yet
CentroidAccumulator &acc = accumulators.create_accumulator(next_accumulator_id, 0); // create new accumulator, height is 0, because we are on the bed
coords_to_check.clear();
coords_to_check.insert(Vec2i(x, y));
while (!coords_to_check.empty()) {
while (!coords_to_check.empty()) { // insert the origin cell coords in to the set, and search all connected cells with volume and without assigned accumulator, assign them to acc
Vec2i current_coords = *coords_to_check.begin();
coords_to_check.erase(coords_to_check.begin());
if (!validate_xy_coords(current_coords)) {
if (!validate_xy_coords(current_coords)) { // invalid coords, drop
continue;
}
Cell &cell = this->access_cell(Vec3i(current_coords.x(), current_coords.y(), 0));
if (cell.volume <= 0 || cell.island_id != std::numeric_limits<int>::max()) {
if (cell.volume <= 0 || cell.accumulator_id != std::numeric_limits<int>::max()) { // empty cell or already assigned, drop
continue;
}
cell.island_id = next_island_id;
cell.accumulator_id = next_accumulator_id; // update cell accumulator id, update the accumulator with the new cell data, add neighbours to queue
Vec3crd cell_center = this->get_cell_center(
Vec3i(current_coords.x(), current_coords.y(), local_z_index_offset));
acc.add_base_points( { Point(cell_center.head<2>()) });
@ -305,12 +307,15 @@ public:
}
}
}
next_island_id--;
next_accumulator_id--;
//base area is separated from the base convex hull - bed accumulators are initialized with convex hull area (TODO compute from number of covered cells instead )
// but support points are initialized with constant, and during merging, the base_areas are added.
acc.base_area = unscale<float>(unscale<float>(acc.base_hull().area())); //apply unscale 2x, it has units of area
}
}
}
// sort support points by height, so that their accumulators ids are also sorted by height
std::sort(issues.supports_nedded.begin(), issues.supports_nedded.end(),
[](const SupportPoint &left, const SupportPoint &right) {
return left.position.z() < right.position.z();
@ -318,21 +323,27 @@ public:
for (int index = 0; index < int(issues.supports_nedded.size()); ++index) {
Vec3i local_coords = this->to_local_cell_coords(
this->to_global_cell_coords(issues.supports_nedded[index].position));
this->access_cell(local_coords).island_id = index;
this->access_cell(local_coords).accumulator_id = index; // assign accumulator id (in case that multiple support points fall into the same cell, they are just overriden)
CentroidAccumulator &acc = accumulators.create_accumulator(index,
issues.supports_nedded[index].position.z());
acc.add_base_points( { Point(scaled(Vec2f(issues.supports_nedded[index].position.head<2>()))) });
acc.base_area = params.support_points_interface_area;
}
// add curling data to each cell
for (const CurledFilament &curling : issues.curling_up) {
this->access_cell(curling.position).curled_height += curling.estimated_height;
}
// keep map of modified accumulator for each layer, so that discarded accumulators are not further checked for stability
// the value of the map is list of cells with curling, that should be further checked for pressure stability with repsect to the accumulator
std::unordered_map<int, std::vector<Vec2i>> modified_acc_ids;
// At the end of the layer check, accumulators are further filtered, since merging causes that single accumulator can have mutliple ids.
// But each accumulator should be checked only once.
std::unordered_map<CentroidAccumulator*, std::vector<Vec2i>> filtered_active_accumulators;
modified_acc_ids.reserve(issues.supports_nedded.size() + 1);
// For each grid layer, cells are added to the accumulators and all active accumulators are checked of stability.
for (int z = 1; z < local_z_cell_count; ++z) {
std::cout << "current z: " << z << std::endl;
@ -341,39 +352,41 @@ public:
for (int x = 0; x < global_cell_count.x(); ++x) {
for (int y = 0; y < global_cell_count.y(); ++y) {
Cell &current = this->access_cell(Vec3i(x, y, z));
// distribute islands info
if (current.volume > 0 && current.island_id == std::numeric_limits<int>::max()) {
int min_island_id_found = std::numeric_limits<int>::max();
// distribute islands info - look for neighbours under the cell, and pick the smallest accumulator id
// also gather all ids, they will be merged to the smallest id accumualtor
if (current.volume > 0 && current.accumulator_id == std::numeric_limits<int>::max()) {
int min_accumulator_id_found = std::numeric_limits<int>::max();
std::unordered_set<int> ids_to_merge { };
for (int y_offset = -2; y_offset <= 2; ++y_offset) {
for (int x_offset = -2; x_offset <= 2; ++x_offset) {
if (validate_xy_coords(Vec2i(x + x_offset, y + y_offset))) {
Cell &under = this->access_cell(Vec3i(x + x_offset, y + y_offset, z - 1));
if (under.island_id < min_island_id_found) {
min_island_id_found = under.island_id;
if (under.accumulator_id < min_accumulator_id_found) {
min_accumulator_id_found = under.accumulator_id;
}
ids_to_merge.insert(under.island_id);
ids_to_merge.insert(under.accumulator_id);
}
}
}
// assign island and update its info
if (min_island_id_found < std::numeric_limits<int>::max()) {
// assign accumulator and update its info
if (min_accumulator_id_found < std::numeric_limits<int>::max()) { // accumulator id found
ids_to_merge.erase(std::numeric_limits<int>::max());
ids_to_merge.erase(min_island_id_found);
current.island_id = min_island_id_found;
ids_to_merge.erase(min_accumulator_id_found);
current.accumulator_id = min_accumulator_id_found; //assign accumualtor id to the cell
for (auto id : ids_to_merge) {
accumulators.merge_to(id, min_island_id_found);
accumulators.merge_to(id, min_accumulator_id_found); //merge other ids to the smallest id
}
CentroidAccumulator &acc = accumulators.access(min_island_id_found);
//update the acc with new point, and store in the modified accumulators map
CentroidAccumulator &acc = accumulators.access(min_accumulator_id_found);
acc.accumulated_value += current.volume
* unscale(this->get_cell_center(this->to_global_cell_coords(Vec3i(x, y, z)))).cast<
float>();
acc.accumulated_volume += current.volume;
modified_acc_ids.emplace(min_island_id_found, std::vector<Vec2i> { });
modified_acc_ids.emplace(min_accumulator_id_found, std::vector<Vec2i> { });
}
}
// distribute curling
// distribute curling (add curling from neighbours under, but also decrease but some factor)
if (current.volume > 0) {
float curled_height = 0;
for (int y_offset = -2; y_offset <= 2; ++y_offset) {
@ -385,18 +398,18 @@ public:
}
}
current.curled_height += std::max(0.0f,
float(curled_height - unscaled(this->cell_size.z()) / 2.0f));
float(curled_height - unscaled(this->cell_size.z()) / 1.5f));
current.curled_height = std::min(current.curled_height,
float(unscaled(this->cell_size.z()) * 2.0f));
std::cout << "Curling: " << current.curled_height << std::endl;
if (current.curled_height / unscaled(this->cell_size.z()) > 0.3f) {
modified_acc_ids[current.island_id].push_back( { x, y });
if (current.curled_height / unscaled(this->cell_size.z()) > 1.5f) { // just a magic threshold number.
modified_acc_ids[current.accumulator_id].push_back( { x, y });
}
}
}
}
//all cells of the grid layer checked, now further filter the modified accumulators, because multiple ids can point to the same acc
filtered_active_accumulators.clear();
for (const auto &pair : modified_acc_ids) {
CentroidAccumulator *acc = &accumulators.access(pair.first);
@ -411,142 +424,89 @@ public:
private:
void check_accumulators_stability(int z, CentroidAccumulators &accumulators,
std::unordered_map<CentroidAccumulator*, std::vector<Vec2i>> filtered_active_accumulators,
Issues &issues, const Params &params) {
Issues &issues, const Params &params) const {
for (const auto &pair : filtered_active_accumulators) {
std::cout << "Z: " << z << std::endl;
CentroidAccumulator &acc = *pair.first;
Vec3f centroid = acc.accumulated_value / acc.accumulated_volume;
Point base_centroid = acc.base_hull().centroid();
Vec3f base_centroid_3d { };
base_centroid_3d << unscale(base_centroid).cast<float>(), acc.base_height;
std::cout << "acc.accumulated_value : " << acc.accumulated_value.x() << " "
<< acc.accumulated_value.y() << " " << acc.accumulated_value.z() << std::endl;
std::cout << "acc.accumulated_volume : " << acc.accumulated_volume << std::endl;
std::cout << "centroid: " << centroid.x() << " " << centroid.y() << " " << centroid.z()
<< std::endl;
Point centroid_base_projection = Point(scaled(Vec2f(centroid.head<2>())));
Point pivot;
double distance_scaled_sq = std::numeric_limits<double>::max();
bool inside = true;
if (acc.base_hull().points.size() == 1) {
pivot = acc.base_hull().points[0];
distance_scaled_sq = (pivot - centroid_base_projection).squaredNorm();
inside = distance_scaled_sq < params.support_points_interface_area;
} else {
for (Line line : acc.base_hull().lines()) {
Point closest_point;
double dist_sq = line.distance_to_squared(centroid_base_projection, &closest_point);
if (dist_sq < distance_scaled_sq) {
pivot = closest_point;
distance_scaled_sq = dist_sq;
}
if ((centroid_base_projection - closest_point).cast<double>().dot(
line.normal().cast<double>())
> 0) {
inside = false;
}
}
}
Vec3f pivot_3d;
pivot_3d << unscale(pivot).cast<float>(), acc.base_height;
float embedded_distance = unscaled(sqrt(distance_scaled_sq));
float centroid_pivot_distance = (centroid - pivot_3d).norm();
float base_center_pivot_distance = float(
unscale(Vec2crd(acc.base_hull().centroid() - pivot)).norm());
std::cout << "centroid inside ? " << inside << " and embedded distance is: " << embedded_distance
std::cout << "base_centroid_3d: " << base_centroid_3d.x() << " " << base_centroid_3d.y() << " "
<< base_centroid_3d.z()
<< std::endl;
bool additional_supports_needed = false;
float sticking_force = acc.base_area
* (acc.base_height == 0 ? params.base_adhesion : params.support_adhesion);
float sticking_torque = base_center_pivot_distance * sticking_force;
std::cout << "sticking force: " << sticking_force << " sticking torque: " << sticking_torque
<< std::endl;
float xy_movement_force = acc.accumulated_volume * params.filament_density
* params.max_acceleration;
float xy_movement_torque = xy_movement_force * centroid_pivot_distance;
std::cout << "xy_movement_force: " << xy_movement_force << " xy_movement_torque: "
<< xy_movement_torque << std::endl;
float weight = acc.accumulated_volume * params.filament_density * params.gravity_constant;
float weight_torque = embedded_distance * weight;
if (!inside) {
weight_torque *= -1;
}
std::cout << "weight: " << weight << " weight_torque: " << weight_torque << std::endl;
float extruder_conflict_torque = params.tolerable_extruder_conflict_force * 2.0f
* centroid_pivot_distance;
std::cout << "extruder_conflict_torque: " << extruder_conflict_torque << std::endl;
float total_momentum = sticking_torque + weight_torque - xy_movement_torque
- extruder_conflict_torque;
additional_supports_needed = total_momentum < 0;
std::cout << "total_momentum: " << total_momentum << std::endl;
std::cout << "additional supports needed: " << additional_supports_needed << std::endl;
if (additional_supports_needed) {
Vec2f attractor_dir =
unscale(Vec2crd(inside ?
pivot - centroid_base_projection :
centroid_base_projection - pivot)).cast<float>().normalized();
Vec2f attractor = unscale(centroid_base_projection).cast<float>() + 10000 * attractor_dir;
std::cout << " attractor: " << attractor.x() << " | " << attractor.y() << std::endl;
double min_dist = std::numeric_limits<double>::max();
Vec3f support_point = centroid;
Vec2i coords = Vec2i(0, 0);
for (int y = 0; y < global_cell_count.y(); ++y) {
for (int x = 0; x < global_cell_count.x(); ++x) {
Cell &cell = this->access_cell(Vec3i(x, y, z));
if (cell.island_id != std::numeric_limits<int>::max() &&
&accumulators.access(cell.island_id) == &acc) {
Vec3f cell_center =
unscale(this->get_cell_center(this->to_global_cell_coords(Vec3i(x, y, z)))).cast<
float>();
float dist_sq = (cell_center.head<2>() - attractor).squaredNorm();
if (dist_sq < min_dist) {
min_dist = dist_sq;
support_point = cell_center;
coords = Vec2i(x, y);
}
// find the cell that is furthest from the base centroid ( its a heurstic to find a possible problems with balance without checking all layer cells)
//TODO better result are if first pivot is chosen as the closest point of the convex hull to the base centroid, and then cell furthest in the direction defined by
// the vector from base centroid to this pivot is taken.
double max_dist_sqr = 0;
Vec3f suspicious_point = centroid;
Vec2i coords = Vec2i(0, 0);
for (int y = 0; y < global_cell_count.y(); ++y) {
for (int x = 0; x < global_cell_count.x(); ++x) {
const Cell &cell = this->access_cell(Vec3i(x, y, z));
if (cell.accumulator_id != std::numeric_limits<int>::max() &&
&accumulators.access(cell.accumulator_id) == &acc) {
Vec3f cell_center =
unscale(this->get_cell_center(this->to_global_cell_coords(Vec3i(x, y, z)))).cast<
float>();
float dist_sq = (cell_center - base_centroid_3d).squaredNorm();
if (dist_sq > max_dist_sqr) {
max_dist_sqr = dist_sq;
suspicious_point = cell_center;
coords = Vec2i(x, y);
}
}
}
}
int final_height_coords = z;
while (final_height_coords > 0
&& this->access_cell(Vec3i(coords.x(), coords.y(), final_height_coords)).volume > 0) {
final_height_coords--;
// for the suspicious point, add movement force in xy (bed sliding, it is assumed that the worst direction is taken, for simplicity)
float xy_movement_force = acc.accumulated_volume * params.filament_density
* params.max_acceleration;
std::cout << "xy_movement_force: " << xy_movement_force << std::endl;
// also add weight (weight is the small factor, because the materials are very light. The weight torque will be computed much higher then what is real,
//since it does not push in the suspicoius point, but in centroid. Its approximation)
float weight = acc.accumulated_volume * params.filament_density * params.gravity_constant;
std::cout << "weight: " << weight << std::endl;
float force = this->check_point_stability_under_pressure(acc, base_centroid, suspicious_point,
xy_movement_force + weight + params.tolerable_extruder_conflict_force,
params);
if (force > 0) {
this->add_suppport_point(Vec3i(coords.x(), coords.y(), z), force, acc, issues, params);
}
for (const Vec2i &cell : pair.second) {
Vec3f pressure_point = unscale(
this->get_cell_center(
this->to_global_cell_coords(Vec3i(cell.x(), cell.y(), z)))).cast<float>();
float force = this->check_point_stability_under_pressure(acc, base_centroid, pressure_point,
params.max_curled_conflict_extruder_force, //TODO add linear scaling of the extruder force based on the curled height (but current data about curled height are kind of unreliable in scale)
params);
if (force > 0) {
this->add_suppport_point(Vec3i(cell.x(), cell.y(), z), force, acc, issues, params);
}
support_point.z() = unscaled(
(final_height_coords + this->local_z_index_offset) * this->cell_size.z());
float expected_force = total_momentum / (support_point - pivot_3d).norm();
std::cout << " new support point: " << support_point.x() << " | " << support_point.y() << " | "
<< support_point.z() << std::endl;
std::cout << " expected_force: " << expected_force << std::endl;
issues.supports_nedded.emplace_back(support_point, expected_force);
acc.add_base_points( { Point::new_scale(Vec2f(support_point.head<2>())) });
acc.base_area += params.support_points_interface_area;
}
}
}
float check_point_stability_under_pressure(CentroidAccumulator &acc, const Point &base_centroid,
const Vec3f &pressure_point, const Params &params) {
const Vec3f &pressure_point, float pressure_force, const Params &params) const {
Point pressure_base_projection = Point(scaled(Vec2f(pressure_point.head<2>())));
Point pivot;
double distance_scaled_sq = std::numeric_limits<double>::max();
bool inside = true;
// find pivot, which is the closest point of the accumulator base hull to pressure point (if the object should fall, it would be over this point)
if (acc.base_hull().points.size() == 1) {
pivot = acc.base_hull().points[0];
distance_scaled_sq = (pivot - pressure_base_projection).squaredNorm();
@ -565,19 +525,25 @@ private:
}
}
}
float embedded_distance = unscaled(sqrt(distance_scaled_sq));
// float embedded_distance = unscaled(sqrt(distance_scaled_sq));
float base_center_pivot_distance = float(unscale(Vec2crd(base_centroid - pivot)).norm());
Vec3f pivot_3d;
pivot_3d << unscale(pivot).cast<float>(), acc.base_height;
//sticking force estimated from the base area and support points
float sticking_force = acc.base_area
* (acc.base_height == 0 ? params.base_adhesion : params.support_adhesion);
float sticking_torque = sticking_force * base_center_pivot_distance;
float pressure_arm = inside ? pressure_point.z() - pivot_3d.z() : (pressure_point - pivot_3d).norm();
float pressure_force = 1.0f;
std::cout << "sticking force: " << sticking_force << std::endl;
std::cout << "pressure force: " << pressure_force << std::endl;
float sticking_torque = sticking_force * base_center_pivot_distance; // sticking torque is computed from the distance to the centroid
float pressure_arm = inside ? pressure_point.z() - pivot_3d.z() : (pressure_point - pivot_3d).norm(); // pressure arm is again higher then in reality,
// since it assumes the worst direction of the pressure force (perpendicular to the vector between pivot and pressure point)
float pressure_torque = pressure_arm * pressure_force;
std::cout << "sticking_torque: " << sticking_torque << std::endl;
std::cout << "pressure_torque: " << pressure_torque << std::endl;
if (sticking_torque < pressure_torque) {
return pressure_force;
} else {
@ -585,6 +551,27 @@ private:
}
}
void add_suppport_point(const Vec3i &local_coords, float expected_force, CentroidAccumulator &acc, Issues &issues,
const Params &params) const {
//add support point - but first finding the lowest full cell is needed, since putting support point to the current cell may end up inside the model,
// essentially doing nothing.
int final_height_coords = local_coords.z();
while (final_height_coords > 0
&& this->access_cell(this->to_global_cell_coords(Vec3i(local_coords.x(), local_coords.y(), final_height_coords))).volume > 0) {
final_height_coords--;
}
Vec3f support_point = unscale(this->get_cell_center(this->to_global_cell_coords(Vec3i(local_coords.x(), local_coords.y(), final_height_coords)))).cast<
float>();
std::cout << " new support point: " << support_point.x() << " | " << support_point.y() << " | "
<< support_point.z() << std::endl;
std::cout << " expected_force: " << expected_force << std::endl;
issues.supports_nedded.emplace_back(support_point, expected_force);
acc.add_base_points( { Point::new_scale(Vec2f(support_point.head<2>())) });
acc.base_area += params.support_points_interface_area;
}
#ifdef DEBUG_FILES
public:
void debug_export() const {
@ -603,18 +590,18 @@ public:
float max_volume = 0;
int min_island_id = 0;
int max_island_id = 0;
float max_curling_height = 0;
float max_curling_height = 0.5f;
for (int x = 0; x < global_cell_count.x(); ++x) {
for (int y = 0; y < global_cell_count.y(); ++y) {
for (int z = 0; z < local_z_cell_count; ++z) {
const Cell &cell = access_cell(Vec3i(x, y, z));
max_volume = std::max(max_volume, cell.volume);
if (cell.island_id != std::numeric_limits<int>::max()) {
min_island_id = std::min(min_island_id, cell.island_id);
max_island_id = std::max(max_island_id, cell.island_id);
if (cell.accumulator_id != std::numeric_limits<int>::max()) {
min_island_id = std::min(min_island_id, cell.accumulator_id);
max_island_id = std::max(max_island_id, cell.accumulator_id);
}
max_curling_height = std::max(max_curling_height, cell.curled_height);
// max_curling_height = std::max(max_curling_height, cell.curled_height);
}
}
}
@ -630,8 +617,8 @@ public:
fprintf(volume_grid_file, "v %f %f %f %f %f %f\n", center(0), center(1), center(2),
volume_color.x(), volume_color.y(), volume_color.z());
}
if (cell.island_id != std::numeric_limits<int>::max()) {
auto island_color = value_to_rgbf(min_island_id, max_island_id + 1, cell.island_id);
if (cell.accumulator_id != std::numeric_limits<int>::max()) {
auto island_color = value_to_rgbf(min_island_id, max_island_id + 1, cell.accumulator_id);
fprintf(islands_grid_file, "v %f %f %f %f %f %f\n", center(0), center(1),
center(2),
island_color.x(), island_color.y(), island_color.z());