Skip to content
Permalink

Comparing changes

Choose two branches to see what’s changed or to start a new pull request. If you need to, you can also or learn more about diff comparisons.

Open a pull request

Create a new pull request by comparing changes across two branches. If you need to, you can also . Learn more about diff comparisons here.
base repository: mockingbirdnest/Principia
Failed to load repositories. Confirm that selected base ref is valid, then try again.
Loading
base: afa7516505f6
Choose a base ref
...
head repository: mockingbirdnest/Principia
Failed to load repositories. Confirm that selected head ref is valid, then try again.
Loading
compare: 592bcce46ae8
Choose a head ref
  • 3 commits
  • 4 files changed
  • 1 contributor

Commits on Jul 11, 2021

  1. Copy the full SHA
    cb630ef View commit details
  2. Lint is stupid.

    pleroy committed Jul 11, 2021
    Copy the full SHA
    1cafa5c View commit details
  3. Merge pull request #3058 from pleroy/DownsamplingForks

    Restructure a bit the code of the discrete trajectory downsampling
    pleroy authored Jul 11, 2021
    Copy the full SHA
    592bcce View commit details
Showing with 112 additions and 123 deletions.
  1. +1 −1 base/macros.hpp
  2. +2 −0 physics/discrete_trajectory.hpp
  3. +45 −39 physics/discrete_trajectory_body.hpp
  4. +64 −83 physics/discrete_trajectory_test.cpp
2 changes: 1 addition & 1 deletion base/macros.hpp
Original file line number Diff line number Diff line change
@@ -165,7 +165,7 @@ inline void noreturn() { std::exit(0); }
// Lexicographic comparison (v1, v2, v3) ≥ (w1, w2, w3).
#define VERSION_GE(v1, v2, v3, w1, w2, w3) \
((v1) > (w1) || ((v1) == (w1) && (v2) > (w2)) || \
((v1) == (w1) && (v2) == (w2) && (v3) >= (w3)))
((v1) == (w1) && (v2) == (w2) && (v3) >= (w3)))

#define CLANG_VERSION_GE(major, minor, patchlevel) \
VERSION_GE(__clang_major__, \
2 changes: 2 additions & 0 deletions physics/discrete_trajectory.hpp
Original file line number Diff line number Diff line change
@@ -264,6 +264,8 @@ class DiscreteTrajectory : public Forkable<DiscreteTrajectory<Frame>,
Hermite3<Instant, Position<Frame>> GetInterpolation(
Iterator const& upper) const;

absl::Status UpdateDownsampling();

Timeline timeline_;

std::optional<Downsampling> downsampling_;
84 changes: 45 additions & 39 deletions physics/discrete_trajectory_body.hpp
Original file line number Diff line number Diff line change
@@ -202,46 +202,10 @@ absl::Status DiscreteTrajectory<Frame>::Append(
<< "Append out of order at " << time << ", last time is "
<< (--timeline_.end())->first;
if (downsampling_.has_value()) {
if (timeline_.size() == 1) {
downsampling_->SetStartOfDenseTimeline(timeline_.begin(), timeline_);
} else {
this->CheckNoForksBefore(this->back().time);
downsampling_->increment_dense_intervals(timeline_);
if (downsampling_->reached_max_dense_intervals()) {
std::vector<TimelineConstIterator> dense_iterators;
// This contains points, hence one more than intervals.
dense_iterators.reserve(downsampling_->max_dense_intervals() + 1);
for (TimelineConstIterator it =
downsampling_->start_of_dense_timeline();
it != timeline_.end();
++it) {
dense_iterators.push_back(it);
}
auto right_endpoints = FitHermiteSpline<Instant, Position<Frame>>(
dense_iterators,
[](auto&& it) -> auto&& { return it->first; },
[](auto&& it) -> auto&& { return it->second.position(); },
[](auto&& it) -> auto&& { return it->second.velocity(); },
downsampling_->tolerance());
if (!right_endpoints.ok()) {
// Note that the actual appending took place; the propagated status
// only reflects a lack of downsampling.
return right_endpoints.status();
}
if (right_endpoints->empty()) {
right_endpoints->push_back(dense_iterators.end() - 1);
}
TimelineConstIterator left = downsampling_->start_of_dense_timeline();
for (const auto& it_in_dense_iterators : right_endpoints.value()) {
TimelineConstIterator const right = *it_in_dense_iterators;
timeline_.erase(++left, right);
left = right;
}
downsampling_->SetStartOfDenseTimeline(left, timeline_);
}
}
return UpdateDownsampling();
} else {
return absl::OkStatus();
}
return absl::OkStatus();
}

template<typename Frame>
@@ -678,6 +642,48 @@ Hermite3<Instant, Position<Frame>> DiscreteTrajectory<Frame>::GetInterpolation(
upper->degrees_of_freedom.velocity()}};
}

template<typename Frame>
absl::Status DiscreteTrajectory<Frame>::UpdateDownsampling() {
if (timeline_.size() == 1) {
downsampling_->SetStartOfDenseTimeline(timeline_.begin(), timeline_);
} else {
this->CheckNoForksBefore(this->back().time);
downsampling_->increment_dense_intervals(timeline_);
if (downsampling_->reached_max_dense_intervals()) {
std::vector<TimelineConstIterator> dense_iterators;
// This contains points, hence one more than intervals.
dense_iterators.reserve(downsampling_->max_dense_intervals() + 1);
for (TimelineConstIterator it = downsampling_->start_of_dense_timeline();
it != timeline_.end();
++it) {
dense_iterators.push_back(it);
}
auto right_endpoints = FitHermiteSpline<Instant, Position<Frame>>(
dense_iterators,
[](auto&& it) -> auto&& { return it->first; },
[](auto&& it) -> auto&& { return it->second.position(); },
[](auto&& it) -> auto&& { return it->second.velocity(); },
downsampling_->tolerance());
if (!right_endpoints.ok()) {
// Note that the actual appending took place; the propagated status only
// reflects a lack of downsampling.
return right_endpoints.status();
}
if (right_endpoints->empty()) {
right_endpoints->push_back(dense_iterators.end() - 1);
}
TimelineConstIterator left = downsampling_->start_of_dense_timeline();
for (const auto& it_in_dense_iterators : right_endpoints.value()) {
TimelineConstIterator const right = *it_in_dense_iterators;
timeline_.erase(++left, right);
left = right;
}
downsampling_->SetStartOfDenseTimeline(left, timeline_);
}
}
return absl::OkStatus();
}

} // namespace internal_discrete_trajectory
} // namespace physics
} // namespace principia
147 changes: 64 additions & 83 deletions physics/discrete_trajectory_test.cpp
Original file line number Diff line number Diff line change
@@ -133,6 +133,39 @@ class DiscreteTrajectoryTest : public testing::Test {
return result;
}

DoublePrecision<Instant> AppendCircularTrajectory(
AngularFrequency const& ω,
Length const& r,
Time const& Δt,
Instant const& t1,
Instant const& t2,
DiscreteTrajectory<World>& trajectory) {
return AppendCircularTrajectory(
ω, r, Δt, DoublePrecision<Instant>(t1), t2, trajectory);
}

DoublePrecision<Instant> AppendCircularTrajectory(
AngularFrequency const& ω,
Length const& r,
Time const& Δt,
DoublePrecision<Instant> const& t1,
Instant const& t2,
DiscreteTrajectory<World>& trajectory) {
Speed const v = ω * r / Radian;
auto t = t1;
for (; t.value <= t2; t.Increment(Δt)) {
DegreesOfFreedom<World> const dof = {
World::origin + Displacement<World>{{r * Cos(ω * (t.value - t0_)),
r * Sin(ω * (t.value - t0_)),
0 * Metre}},
Velocity<World>{{-v * Sin(ω * (t.value - t0_)),
v * Cos(ω * (t.value - t0_)),
0 * Metre / Second}}};
trajectory.Append(t.value, dof);
}
return t;
}

Position<World> q1_, q2_, q3_, q4_;
Velocity<World> p1_, p2_, p3_, p4_;
DegreesOfFreedom<World> d1_, d2_, d3_, d4_;
@@ -800,20 +833,14 @@ TEST_F(DiscreteTrajectoryTest, IteratorSuccess) {

TEST_F(DiscreteTrajectoryTest, QuadrilateralCircle) {
DiscreteTrajectory<World> circle;

AngularFrequency const ω = 3 * Radian / Second;
Length const r = 2 * Metre;
Speed const v = ω * r / Radian;
Time const period = 2 * π * Radian / ω;
for (Time t; t <= period; t += period / 4) {
circle.Append(
t0_ + t,
{World::origin + Displacement<World>{{r * Cos(ω * t),
r * Sin(ω * t),
0 * Metre}},
Velocity<World>{{-v * Sin(ω * t),
v * Cos(ω * t),
0 * Metre / Second}}});
}
AppendCircularTrajectory(
ω, r, /*Δt=*/period / 4, /*t1=*/t0_, /*t2=*/t0_ + period, circle);

double max_r_error = 0;
double max_v_error = 0;
for (Time t; t < period; t += period / 32) {
@@ -839,20 +866,12 @@ TEST_F(DiscreteTrajectoryTest, Downsampling) {
/*tolerance=*/1 * Milli(Metre));
AngularFrequency const ω = 3 * Radian / Second;
Length const r = 2 * Metre;
Speed const v = ω * r / Radian;
for (auto t = DoublePrecision<Instant>(t0_);
t.value <= t0_ + 10 * Second;
t.Increment(10 * Milli(Second))) {
DegreesOfFreedom<World> const dof =
{World::origin + Displacement<World>{{r * Cos(ω * (t.value - t0_)),
r * Sin(ω * (t.value - t0_)),
0 * Metre}},
Velocity<World>{{-v * Sin(ω * (t.value - t0_)),
v * Cos(ω * (t.value - t0_)),
0 * Metre / Second}}};
circle.Append(t.value, dof);
downsampled_circle.Append(t.value, dof);
}
Time const Δt = 10 * Milli(Second);
Instant const t1 = t0_;
Instant const t2 = t0_ + 10 * Second;
AppendCircularTrajectory(ω, r, Δt, t1, t2, circle);
AppendCircularTrajectory(ω, r, Δt, t1, t2, downsampled_circle);

EXPECT_THAT(circle.Size(), Eq(1001));
EXPECT_THAT(downsampled_circle.Size(), Eq(77));
std::vector<Length> errors;
@@ -867,29 +886,18 @@ TEST_F(DiscreteTrajectoryTest, Downsampling) {

TEST_F(DiscreteTrajectoryTest, DownsamplingSerialization) {
DiscreteTrajectory<World> circle;
auto deserialized_circle = make_not_null_unique<DiscreteTrajectory<World>>();
circle.SetDownsampling(/*max_dense_intervals=*/50,
/*tolerance=*/1 * Milli(Metre));
deserialized_circle->SetDownsampling(/*max_dense_intervals=*/50,
/*tolerance=*/1 * Milli(Metre));
AngularFrequency const ω = 3 * Radian / Second;
Length const r = 2 * Metre;
Speed const v = ω * r / Radian;
auto t = DoublePrecision<Instant>(t0_);
for (; t.value <= t0_ + 5 * Second; t.Increment(10 * Milli(Second))) {
DegreesOfFreedom<World> const dof =
{World::origin + Displacement<World>{{r * Cos(ω * (t.value - t0_)),
r * Sin(ω * (t.value - t0_)),
0 * Metre}},
Velocity<World>{{-v * Sin(ω * (t.value - t0_)),
v * Cos(ω * (t.value - t0_)),
0 * Metre / Second}}};
circle.Append(t.value, dof);
deserialized_circle->Append(t.value, dof);
}
Time const Δt = 10 * Milli(Second);
Instant const t1 = t0_;
Instant const t2 = t0_ + 5 * Second;
auto const circle_tmax = AppendCircularTrajectory(ω, r, Δt, t1, t2, circle);

serialization::DiscreteTrajectory message;
deserialized_circle->WriteToMessage(&message, /*forks=*/{});
deserialized_circle =
circle.WriteToMessage(&message, /*forks=*/{});
auto deserialized_circle =
DiscreteTrajectory<World>::ReadFromMessage(message, /*forks=*/{});

// Serialization/deserialization preserves the size, the times, and nudges the
@@ -910,25 +918,15 @@ TEST_F(DiscreteTrajectoryTest, DownsamplingSerialization) {

// Appending may result in different downsampling because the positions differ
// a bit.
for (; t.value <= t0_ + 10 * Second; t.Increment(10 * Milli(Second))) {
DegreesOfFreedom<World> const dof =
{World::origin + Displacement<World>{{r * Cos(ω * (t.value - t0_)),
r * Sin(ω * (t.value - t0_)),
0 * Metre}},
Velocity<World>{{-v * Sin(ω * (t.value - t0_)),
v * Cos(ω * (t.value - t0_)),
0 * Metre / Second}}};
circle.Append(t.value, dof);
deserialized_circle->Append(t.value, dof);
}
Instant const t3 = t0_ + 10 * Second;
AppendCircularTrajectory(ω, r, Δt, circle_tmax, t3, circle);
AppendCircularTrajectory(ω, r, Δt, circle_tmax, t3, *deserialized_circle);

// Despite the difference in downsampling (and therefore in size) the two
// trajectories are still within the tolerance.
EXPECT_THAT(circle.Size(), Eq(77));
EXPECT_THAT(deserialized_circle->Size(), Eq(78));
for (auto t = DoublePrecision<Instant>(t0_);
t.value <= t0_ + 10 * Second;
t.Increment(10 * Milli(Second))) {
for (auto t = DoublePrecision<Instant>(t0_); t.value <= t3; t.Increment(Δt)) {
EXPECT_THAT(deserialized_circle->EvaluatePosition(t.value),
AbsoluteErrorFrom(circle.EvaluatePosition(t.value),
Le(0.23 * Milli(Metre))));
@@ -947,33 +945,16 @@ TEST_F(DiscreteTrajectoryTest, DownsamplingForgetAfter) {
/*tolerance=*/1 * Milli(Metre));
AngularFrequency const ω = 3 * Radian / Second;
Length const r = 2 * Metre;
Speed const v = ω * r / Radian;
auto t = t0_;
for (; t <= t0_ + 10 * Second; t += 10 * Milli(Second)) {
DegreesOfFreedom<World> const dof =
{World::origin + Displacement<World>{{r * Cos(ω * (t - t0_)),
r * Sin(ω * (t - t0_)),
0 * Metre}},
Velocity<World>{{-v * Sin(ω * (t - t0_)),
v * Cos(ω * (t - t0_)),
0 * Metre / Second}}};
circle.Append(t, dof);
forgotten_circle.Append(t, dof);
}
forgotten_circle.ForgetAfter(t0_ + 5 * Second);
t = forgotten_circle.back().time;
for (t += 10 * Milli(Second); t <= t0_ + 10 * Second;
t += 10 * Milli(Second)) {
DegreesOfFreedom<World> const dof =
{World::origin + Displacement<World>{{r * Cos(ω * (t - t0_)),
r * Sin(ω * (t - t0_)),
0 * Metre}},
Velocity<World>{{-v * Sin(ω * (t - t0_)),
v * Cos(ω * (t - t0_)),
0 * Metre / Second}}};
forgotten_circle.Append(t, dof);
}
EXPECT_THAT(circle.Size(), Eq(77));
Time const Δt = 1.0 / 128.0 * Second; // Yields exact times.
Instant const t1 = t0_;
Instant const t2 = t0_ + 10 * Second;
AppendCircularTrajectory(ω, r, Δt, t1, t2, circle);
AppendCircularTrajectory(ω, r, Δt, t1, t2, forgotten_circle);

forgotten_circle.ForgetAfter(t2);
AppendCircularTrajectory(
ω, r, Δt, forgotten_circle.back().time + Δt, t2, forgotten_circle);
EXPECT_THAT(circle.Size(), Eq(93));
EXPECT_THAT(forgotten_circle.Size(), Eq(circle.Size()));
std::vector<Length> errors;
for (auto const [time, degrees_of_freedom] : forgotten_circle) {