Time estimator fixes (#15304)

* Plan corners with junction deviation where the firmware uses it

The time estimator only ever had the classic per-axis jerk model, which limits a
corner by the largest single-axis component of the velocity change. That is
anisotropic: the same corner is allowed sqrt(2) more speed on a diagonal than on
an axis, which paints a four-lobed ripple around every circular wall in the
actual speed and actual flow views, worst on small parts whose walls are made of
short segments.

Klipper has no classic jerk at all and Marlin 2 has none while M205 J is in use;
both plan corners with junction deviation, which sees only the corner angle. Add
that model and use it for those machines:

  - Klipper: derived from the square corner velocity, as the firmware does
    (jd = scv^2 * (sqrt(2) - 1) / max_accel), reading the scv from
    machine_max_jerk_x, where process_SET_VELOCITY_LIMIT() already stores
    SQUARE_CORNER_VELOCITY.
  - Marlin 2: machine_max_junction_deviation, which was already loaded into the
    machine limits but never reached the planner.
  - Every other flavor keeps the classic jerk path unchanged.

The model has no per-axis jerk floor, so this also drops the hard slow spot the
estimator drew at the start of every loop from machine_max_jerk_e.

Toolpaths are unaffected: on a full export the only lines that change are M73.

The junction deviation maths, including Marlin's JD_HANDLE_SMALL_SEGMENTS arc
approximation, is ported from PrusaSlicer's src/libslic3r/GCode/GCodeProcessor.cpp.
The Klipper mapping is not in PrusaSlicer, which ignores SET_VELOCITY_LIMIT.

* Add tests for junction deviation corner planning

Cover the three properties the change rests on:

  - a right angle on Klipper is planned at exactly the square corner velocity,
    the identity that makes the scv to junction deviation mapping correct, and a
    shallow corner is planned far faster than per-axis jerk allows;
  - junction deviation gives the same speed whatever the corner's orientation,
    while classic jerk keeps its sqrt(2) spread, which is the four-lobed ripple;
  - machines that do not plan with junction deviation are provably untouched,
    including a Marlin 2 printer that has it disabled.
This commit is contained in:
Ian Bassi
2026-08-20 09:16:02 -03:00
committed by GitHub
parent dc75ce6811
commit aaa8e98bb0
3 changed files with 307 additions and 20 deletions

View File

@@ -298,6 +298,7 @@ void GCodeProcessor::TimeMachine::State::reset()
//BBS
enter_direction = { 0.0f, 0.0f, 0.0f };
exit_direction = { 0.0f, 0.0f, 0.0f };
jd_unit_vec = { 0.0f, 0.0f, 0.0f, 0.0f };
}
void GCodeProcessor::TimeMachine::CustomGCodeTime::reset()
@@ -5036,6 +5037,10 @@ void GCodeProcessor::process_G1(const std::array<std::optional<double>, 4>& axes
if (!is_extrusion_only_move(delta_pos))
curr.enter_direction = curr.enter_direction / norm;
curr.exit_direction = curr.enter_direction;
curr.jd_unit_vec = Vec4f(static_cast<float>(delta_pos[X]) * inv_distance,
static_cast<float>(delta_pos[Y]) * inv_distance,
static_cast<float>(delta_pos[Z]) * inv_distance,
static_cast<float>(delta_pos[E]) * inv_distance);
TimeBlock block;
block.move_type = type;
@@ -5118,22 +5123,32 @@ void GCodeProcessor::process_G1(const std::array<std::optional<double>, 4>& axes
block.acceleration = acceleration;
// calculates block exit feedrate
curr.safe_feedrate = block.feedrate_profile.cruise;
static const float PREVIOUS_FEEDRATE_THRESHOLD = 0.0001f;
const bool has_prev_move = !blocks.empty() && prev.feedrate > PREVIOUS_FEEDRATE_THRESHOLD;
for (unsigned char a = X; a <= E; ++a) {
float axis_max_jerk = get_axis_max_jerk(static_cast<PrintEstimatedStatistics::ETimeMode>(i), static_cast<Axis>(a));
if (curr.abs_axis_feedrate[a] > axis_max_jerk)
curr.safe_feedrate = std::min(curr.safe_feedrate, axis_max_jerk);
// Orca: junction deviation where the firmware uses it (Klipper always, Marlin 2 with M205 J).
// Negative leaves the classic jerk path below unchanged.
const float vmax_junction_jd = calc_vmax_junction_deviation(block, prev, curr, has_prev_move,
static_cast<PrintEstimatedStatistics::ETimeMode>(i));
const bool use_junction_deviation = vmax_junction_jd >= 0.0f;
// calculates block exit feedrate. Junction deviation has no per axis jerk floor, so a move is
// free to start from rest.
curr.safe_feedrate = use_junction_deviation ? 0.0f : block.feedrate_profile.cruise;
if (!use_junction_deviation) {
for (unsigned char a = X; a <= E; ++a) {
float axis_max_jerk = get_axis_max_jerk(static_cast<PrintEstimatedStatistics::ETimeMode>(i), static_cast<Axis>(a));
if (curr.abs_axis_feedrate[a] > axis_max_jerk)
curr.safe_feedrate = std::min(curr.safe_feedrate, axis_max_jerk);
}
}
block.feedrate_profile.exit = curr.safe_feedrate;
static const float PREVIOUS_FEEDRATE_THRESHOLD = 0.0001f;
// calculates block entry feedrate
float vmax_junction = curr.safe_feedrate;
if (!blocks.empty() && prev.feedrate > PREVIOUS_FEEDRATE_THRESHOLD) {
float vmax_junction = use_junction_deviation ? vmax_junction_jd : curr.safe_feedrate;
if (!use_junction_deviation && has_prev_move) {
bool prev_speed_larger = prev.feedrate > block.feedrate_profile.cruise;
float smaller_speed_factor = prev_speed_larger ? (block.feedrate_profile.cruise / prev.feedrate) : (prev.feedrate / block.feedrate_profile.cruise);
// Pick the smaller of the nominal speeds. Higher speed shall not be achieved at the junction during coasting.
@@ -5400,6 +5415,10 @@ void GCodeProcessor::process_VG1(const GCodeReader::GCodeLine& line)
if (!is_extrusion_only_move(delta_pos))
curr.enter_direction = curr.enter_direction / norm;
curr.exit_direction = curr.enter_direction;
curr.jd_unit_vec = Vec4f(static_cast<float>(delta_pos[X]) * inv_distance,
static_cast<float>(delta_pos[Y]) * inv_distance,
static_cast<float>(delta_pos[Z]) * inv_distance,
static_cast<float>(delta_pos[E]) * inv_distance);
TimeBlock block;
block.move_type = type;
@@ -5480,22 +5499,32 @@ void GCodeProcessor::process_VG1(const GCodeReader::GCodeLine& line)
block.acceleration = acceleration;
// calculates block exit feedrate
curr.safe_feedrate = block.feedrate_profile.cruise;
static const float PREVIOUS_FEEDRATE_THRESHOLD = 0.0001f;
const bool has_prev_move = !blocks.empty() && prev.feedrate > PREVIOUS_FEEDRATE_THRESHOLD;
for (unsigned char a = X; a <= E; ++a) {
float axis_max_jerk = get_axis_max_jerk(static_cast<PrintEstimatedStatistics::ETimeMode>(i), static_cast<Axis>(a));
if (curr.abs_axis_feedrate[a] > axis_max_jerk)
curr.safe_feedrate = std::min(curr.safe_feedrate, axis_max_jerk);
// Orca: junction deviation where the firmware uses it (Klipper always, Marlin 2 with M205 J).
// Negative leaves the classic jerk path below unchanged.
const float vmax_junction_jd = calc_vmax_junction_deviation(block, prev, curr, has_prev_move,
static_cast<PrintEstimatedStatistics::ETimeMode>(i));
const bool use_junction_deviation = vmax_junction_jd >= 0.0f;
// calculates block exit feedrate. Junction deviation has no per axis jerk floor, so a move is
// free to start from rest.
curr.safe_feedrate = use_junction_deviation ? 0.0f : block.feedrate_profile.cruise;
if (!use_junction_deviation) {
for (unsigned char a = X; a <= E; ++a) {
float axis_max_jerk = get_axis_max_jerk(static_cast<PrintEstimatedStatistics::ETimeMode>(i), static_cast<Axis>(a));
if (curr.abs_axis_feedrate[a] > axis_max_jerk)
curr.safe_feedrate = std::min(curr.safe_feedrate, axis_max_jerk);
}
}
block.feedrate_profile.exit = curr.safe_feedrate;
static const float PREVIOUS_FEEDRATE_THRESHOLD = 0.0001f;
// calculates block entry feedrate
float vmax_junction = curr.safe_feedrate;
if (!blocks.empty() && prev.feedrate > PREVIOUS_FEEDRATE_THRESHOLD) {
float vmax_junction = use_junction_deviation ? vmax_junction_jd : curr.safe_feedrate;
if (!use_junction_deviation && has_prev_move) {
bool prev_speed_larger = prev.feedrate > block.feedrate_profile.cruise;
float smaller_speed_factor = prev_speed_larger ? (block.feedrate_profile.cruise / prev.feedrate) : (prev.feedrate / block.feedrate_profile.cruise);
// Pick the smaller of the nominal speeds. Higher speed shall not be achieved at the junction during coasting.
@@ -7168,6 +7197,86 @@ float GCodeProcessor::get_axis_max_jerk_with_jd(PrintEstimatedStatistics::ETimeM
return get_axis_max_jerk_with_jd(mode, axis, get_acceleration(mode));
}
float GCodeProcessor::get_junction_deviation(PrintEstimatedStatistics::ETimeMode mode, float acceleration) const
{
const size_t id = static_cast<size_t>(mode);
// Klipper has no classic jerk: jd = scv^2 * (sqrt(2) - 1) / max_accel
// (toolhead.py::_calc_junction_deviation). Passing the block acceleration back in makes it cancel
// in calc_vmax_junction_deviation(), leaving the identity v == scv at a 90 degree corner.
if (m_flavor == gcfKlipper) {
// machine_max_jerk_x holds the square corner velocity; process_SET_VELOCITY_LIMIT() writes it.
const float scv = get_option_value(m_time_processor.machine_limits.machine_max_jerk_x, id);
if (scv <= 0.0f || acceleration <= 0.0f)
return 0.0f;
return sqr(scv) * (std::sqrt(2.0f) - 1.0f) / acceleration;
}
// Marlin 2 plans with junction deviation only when M205 J > 0; classic jerk leaves it at 0.
if (m_flavor == gcfMarlinFirmware)
return get_option_value(m_time_processor.machine_limits.machine_max_junction_deviation, id);
return 0.0f;
}
float GCodeProcessor::calc_junction_acceleration(const TimeBlock& block, const Vec4f& junction_unit_vec,
PrintEstimatedStatistics::ETimeMode mode) const
{
float junction_acceleration = block.acceleration;
for (unsigned char a = X; a <= E; ++a) {
if (junction_unit_vec[a] == 0.0f)
continue;
const float axis_max_acceleration = get_axis_max_acceleration(mode, static_cast<Axis>(a), m_machine_config_idx);
if (axis_max_acceleration > 0.0f)
junction_acceleration = std::min(junction_acceleration, std::abs(axis_max_acceleration / junction_unit_vec[a]));
}
return junction_acceleration;
}
// Ported from PrusaSlicer (src/libslic3r/GCode/GCodeProcessor.cpp).
float GCodeProcessor::calc_vmax_junction_deviation(const TimeBlock& block, const TimeMachine::State& prev,
const TimeMachine::State& curr, bool has_prev_move,
PrintEstimatedStatistics::ETimeMode mode) const
{
const float junction_deviation = get_junction_deviation(mode, block.acceleration);
if (junction_deviation <= 0.0f)
return -1.0f; // classic jerk machine, the caller keeps its own computation
if (!has_prev_move)
return 0.0f; // starts from rest, the planner raises this on the reverse pass
// -1 for a straight continuation, +1 for a full reversal. Half angle identity, no acos()/sin().
float junction_cos_theta = (-prev.jd_unit_vec).dot(curr.jd_unit_vec);
if (junction_cos_theta > 0.999999f)
return 0.0f; // the path doubles back, the machine has to stop
junction_cos_theta = std::max(junction_cos_theta, -0.999999f); // guards the division below
const float sin_theta_d2 = std::sqrt(0.5f * (1.0f - junction_cos_theta)); // always positive
const Vec4f junction_vec = curr.jd_unit_vec - prev.jd_unit_vec;
const float junction_vec_norm = junction_vec.norm();
const Vec4f junction_unit_vec = (junction_vec_norm > 0.0f) ? Vec4f(junction_vec / junction_vec_norm)
: Vec4f(0.0f, 0.0f, 0.0f, 0.0f);
const float junction_acceleration = calc_junction_acceleration(block, junction_unit_vec, mode);
float vmax_junction_sqr = (junction_acceleration * junction_deviation * sin_theta_d2) / (1.0f - sin_theta_d2);
// Marlin's JD_HANDLE_SMALL_SEGMENTS: a short move through a shallow corner is treated as an arc and
// capped by the centripetal acceleration it needs. Klipper has no equivalent.
if (m_flavor != gcfKlipper && block.distance < 1.0f && junction_cos_theta < -0.7071067812f) {
// Fast acos(-t), max. error +-0.033rad. MinMax polynomial by W. Randolph Franklin:
// https://wrf.ecse.rpi.edu/Research/Short_Notes/arcsin/onlyelem.html
const float neg = junction_cos_theta < 0.0f ? -1.0f : 1.0f;
const float t = neg * junction_cos_theta;
const float asinx = 0.032843707f + t * (-1.451838349f + t * (29.66153956f + t * (-131.1123477f +
t * (262.8130562f + t * (-242.7199627f + t * (84.31466202f))))));
const float junction_theta = float(0.5 * M_PI) + neg * asinx; // acos(-t), bottoms out at 0.033
vmax_junction_sqr = std::min(vmax_junction_sqr, (block.distance * junction_acceleration) / junction_theta);
}
// Never faster than either of the two moves the junction joins.
vmax_junction_sqr = std::min(vmax_junction_sqr, std::min(sqr(block.feedrate_profile.cruise), sqr(prev.feedrate)));
return std::sqrt(vmax_junction_sqr);
}
float GCodeProcessor::get_axis_max_jerk(PrintEstimatedStatistics::ETimeMode mode, Axis axis) const
{
const size_t id = static_cast<size_t>(mode);

View File

@@ -637,6 +637,10 @@ class Print;
//For line move, there are same. For arc move, there are different.
Vec3f enter_direction;
Vec3f exit_direction;
// Orca: move direction over all four axes, scaled by 1 / block.distance. Used by
// calc_vmax_junction_deviation(), which needs E to see extrusion-rate changes
// between collinear moves the way Marlin and Klipper do.
Vec4f jd_unit_vec;
void reset();
};
@@ -1488,6 +1492,16 @@ class Print;
float get_axis_max_acceleration(PrintEstimatedStatistics::ETimeMode mode, Axis axis, int machine_idx) const;
float get_axis_max_jerk_with_jd(PrintEstimatedStatistics::ETimeMode mode, Axis axis, float acceleration) const;
float get_axis_max_jerk_with_jd(PrintEstimatedStatistics::ETimeMode mode, Axis axis) const;
// Orca: junction deviation for a block at the given acceleration, 0 for a classic jerk machine.
float get_junction_deviation(PrintEstimatedStatistics::ETimeMode mode, float acceleration) const;
// Orca: acceleration along the junction direction, clamped by the per axis limits.
float calc_junction_acceleration(const TimeBlock& block, const Vec4f& junction_unit_vec,
PrintEstimatedStatistics::ETimeMode mode) const;
// Orca: entry speed from the junction deviation model, which limits a corner by its angle alone
// and is therefore isotropic, unlike per axis jerk. Negative means classic jerk applies instead.
float calc_vmax_junction_deviation(const TimeBlock& block, const TimeMachine::State& prev,
const TimeMachine::State& curr, bool has_prev_move,
PrintEstimatedStatistics::ETimeMode mode) const;
float get_axis_max_jerk(PrintEstimatedStatistics::ETimeMode mode, Axis axis) const;
Vec3f get_xyz_max_jerk(PrintEstimatedStatistics::ETimeMode mode) const;
float get_retract_acceleration(PrintEstimatedStatistics::ETimeMode mode) const;

View File

@@ -7,9 +7,14 @@
#include "test_utils.hpp"
#include <cmath>
#include <fstream>
#include <iomanip>
#include <map>
#include <memory>
#include <sstream>
#include <string>
#include <vector>
using namespace Slic3r;
using Catch::Matchers::WithinAbs;
@@ -418,3 +423,162 @@ TEST_CASE("Per-slot machine limits follow the active nozzle", "[GCodeTiming][Mul
REQUIRE_THAT(times[2], Catch::Matchers::WithinRel(101.0 / 200.0, 0.10));
}
}
// Junction planning decides the speeds the "actual speed" / "actual flow" preview shows. Per-axis
// jerk limits a corner by the largest single-axis component of the velocity change, allowing sqrt(2)
// more speed on a diagonal than on an axis -- a four-lobed ripple around every circle. Klipper and
// Marlin 2 with M205 J plan with junction deviation instead, which sees only the corner angle.
namespace {
// One acceleration everywhere and axis limits far above it, so only the junction model under test
// can slow a corner down.
FullPrintConfig make_junction_config(GCodeFlavor flavor, double corner_velocity, double junction_deviation)
{
FullPrintConfig config;
config.gcode_flavor.value = flavor;
config.filament_diameter.values = {1.75};
config.filament_map.values = {1};
const std::vector<double> accel = {1000.0, 1000.0};
const std::vector<double> axis = {20000.0, 20000.0};
const std::vector<double> speed = {500.0, 500.0};
config.machine_max_acceleration_extruding.values = accel;
config.machine_max_acceleration_travel.values = accel;
config.machine_max_acceleration_retracting.values = accel;
config.machine_max_acceleration_x.values = axis;
config.machine_max_acceleration_y.values = axis;
config.machine_max_acceleration_z.values = axis;
config.machine_max_acceleration_e.values = axis;
config.machine_max_speed_x.values = speed;
config.machine_max_speed_y.values = speed;
config.machine_max_speed_z.values = speed;
config.machine_max_speed_e.values = speed;
// Klipper reads this as the square corner velocity, Marlin as classic jerk.
config.machine_max_jerk_x.values = {corner_velocity, corner_velocity};
config.machine_max_jerk_y.values = {corner_velocity, corner_velocity};
config.machine_max_jerk_z.values = {corner_velocity, corner_velocity};
// Kept out of the way so it never binds in the classic-jerk comparisons.
config.machine_max_jerk_e.values = {100.0, 100.0};
config.machine_max_junction_deviation.values = {junction_deviation, junction_deviation};
config.machine_min_extruding_rate.values = {0.0, 0.0};
config.machine_min_travel_rate.values = {0.0, 0.0};
return config;
}
constexpr double junction_x = 60.0;
constexpr double junction_y = 60.0;
// Two 40mm travels meeting at (junction_x, junction_y) with the given turn, rotated by `orientation`.
// 40mm is long enough to reach the commanded 150mm/s and brake back to any corner speed these tests
// produce. Travels (no E) keep the junction vector purely geometric, as the formulas below assume.
std::string corner_gcode(double turn_deg, double orientation_deg)
{
const double len = 40.0;
const double a_in = orientation_deg * M_PI / 180.0;
const double a_out = (orientation_deg + turn_deg) * M_PI / 180.0;
std::ostringstream os;
os << std::fixed << std::setprecision(4)
<< "M83\n"
<< "G1 Z0.2 F1200\n"
<< "G1 X" << junction_x - len * std::cos(a_in) << " Y" << junction_y - len * std::sin(a_in) << " F6000\n"
<< "G1 X" << junction_x << " Y" << junction_y << " F9000\n"
<< "G1 X" << junction_x + len * std::cos(a_out) << " Y" << junction_y + len * std::sin(a_out) << " F9000\n";
return os.str();
}
// Speed allowed through the corner: the vertex ending the incoming move carries that block's exit
// speed, and the vertices the actual-speed pass inserts are all strictly interior.
double corner_speed(const GCodeProcessorResult& r)
{
for (const auto& mv : r.moves)
if (mv.type == EMoveType::Travel &&
std::abs(mv.position.x() - junction_x) < 1e-3 &&
std::abs(mv.position.y() - junction_y) < 1e-3)
return mv.actual_feedrate;
return -1.0;
}
double planned_corner_speed(GCodeFlavor flavor, double corner_velocity, double junction_deviation,
double turn_deg, double orientation_deg = 0.0)
{
GCodeProcessor proc;
run_processor(proc, make_junction_config(flavor, corner_velocity, junction_deviation),
corner_gcode(turn_deg, orientation_deg).c_str());
return corner_speed(proc.get_result());
}
} // namespace
TEST_CASE("Klipper corners are planned with junction deviation derived from the square corner velocity",
"[GCodeTiming][JunctionDeviation]")
{
// jd = scv^2 * (sqrt(2) - 1) / max_accel, then v^2 = jd * accel * sin(t/2) / (1 - sin(t/2)).
// The acceleration cancels: the corner speed depends only on the scv and the angle.
const double scv = 5.0;
SECTION("a right angle is taken at exactly the square corner velocity") {
// sin(t/2) = sqrt(0.5) at 90 degrees, so v == scv -- the definition of the square corner
// velocity, and what makes the mapping above the right one.
REQUIRE_THAT(planned_corner_speed(gcfKlipper, scv, 0.0, 90.0), Catch::Matchers::WithinRel(scv, 0.02));
}
SECTION("a shallow corner is taken far faster than the per-axis jerk model allows") {
// 6 degrees: sin(t/2) = cos(3 deg), so v = 5 * sqrt((sqrt(2) - 1) * 728.68) = 86.9mm/s. Per-axis
// jerk ignores the angle and caps the velocity *change* (2v*sin(3 deg)), giving 47.8mm/s.
const double jd_speed = planned_corner_speed(gcfKlipper, scv, 0.0, 6.0);
const double jerk_speed = planned_corner_speed(gcfMarlinLegacy, scv, 0.0, 6.0);
REQUIRE_THAT(jd_speed, Catch::Matchers::WithinRel(86.87, 0.02));
REQUIRE_THAT(jerk_speed, Catch::Matchers::WithinRel(47.75, 0.02));
}
}
TEST_CASE("Junction deviation limits a corner by its angle alone, not by its orientation",
"[GCodeTiming][JunctionDeviation]")
{
// The four-lobed ripple on circular walls is per-axis jerk being anisotropic: a velocity change
// lying on an axis gets sqrt(2) less headroom than the same change on the diagonal.
const double scv = 5.0;
const double turn = 6.0;
SECTION("Klipper plans both orientations identically") {
const double on_axis = planned_corner_speed(gcfKlipper, scv, 0.0, turn, 0.0);
const double diagonal = planned_corner_speed(gcfKlipper, scv, 0.0, turn, 45.0);
REQUIRE(on_axis > 0.0);
REQUIRE_THAT(diagonal, Catch::Matchers::WithinRel(on_axis, 0.02));
}
SECTION("the classic jerk model keeps its orientation dependence") {
const double on_axis = planned_corner_speed(gcfMarlinLegacy, scv, 0.0, turn, 0.0);
const double diagonal = planned_corner_speed(gcfMarlinLegacy, scv, 0.0, turn, 45.0);
REQUIRE(on_axis > 0.0);
REQUIRE(diagonal / on_axis > 1.2);
}
}
TEST_CASE("Junction deviation is only used where the firmware actually plans with it",
"[GCodeTiming][JunctionDeviation]")
{
const double jerk = 5.0;
SECTION("Marlin 2 with M205 J disabled keeps the classic jerk planning") {
// machine_max_junction_deviation == 0 is how a Marlin 2 printer says it runs classic jerk.
const double classic = planned_corner_speed(gcfMarlinLegacy, jerk, 0.0, 90.0);
REQUIRE(classic > 0.0);
REQUIRE_THAT(planned_corner_speed(gcfMarlinFirmware, jerk, 0.0, 90.0),
Catch::Matchers::WithinRel(classic, 1e-4));
}
SECTION("Marlin 2 with M205 J enabled switches to junction deviation") {
// sqrt(1000 * 0.05 * 2.4142136) = 11.0mm/s, independent of the jerk values it no longer reads.
REQUIRE_THAT(planned_corner_speed(gcfMarlinFirmware, jerk, 0.05, 90.0),
Catch::Matchers::WithinRel(10.99, 0.02));
}
SECTION("machines without junction deviation are untouched by the jerk values it would ignore") {
// A flavor that never enters the junction deviation path must ignore the setting entirely.
const double without = planned_corner_speed(gcfMarlinLegacy, jerk, 0.0, 90.0);
REQUIRE_THAT(planned_corner_speed(gcfMarlinLegacy, jerk, 0.05, 90.0),
Catch::Matchers::WithinRel(without, 1e-4));
}
}