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: f2d5a6ee380c
Choose a base ref
...
head repository: mockingbirdnest/Principia
Failed to load repositories. Confirm that selected head ref is valid, then try again.
Loading
compare: 74e4b070b628
Choose a head ref
  • 13 commits
  • 6 files changed
  • 2 contributors

Commits on May 10, 2021

  1. Copy the full SHA
    53b4c55 View commit details
  2. Copy the full SHA
    c2ba900 View commit details
  3. Copy the full SHA
    d1d9efd View commit details
  4. Copy the full SHA
    3e774f2 View commit details
  5. Cache expensive begin and end calls in

    DiscreteTrajectory iteration benchmarks.
    rnlahaye committed May 10, 2021
    Copy the full SHA
    fc81d60 View commit details
  6. Avoid interpolation in DiscreteTrajecotry::Evaluate*

    when unneccessary.
    rnlahaye committed May 10, 2021
    Copy the full SHA
    dbc80dc View commit details
  7. Optimizations for DiscreteTrajectory front, back,

    and begin.
    rnlahaye committed May 10, 2021
    Copy the full SHA
    f52362a View commit details
  8. Reconfigure Forkable iterator reference conversion.

    The conversion is now done by a static method instead of a constructor.
    There is also some documentation now.
    Also fixed inadvertently broken forkable_test.
    rnlahaye committed May 10, 2021
    Copy the full SHA
    89ca90a View commit details

Commits on May 11, 2021

  1. Copy the full SHA
    d368f1b View commit details

Commits on May 13, 2021

  1. Pull request comments.

    rnlahaye committed May 13, 2021
    Copy the full SHA
    d225a6b View commit details

Commits on May 14, 2021

  1. Copy the full SHA
    1201fda View commit details
  2. Lint.

    rnlahaye committed May 14, 2021
    Copy the full SHA
    4fcbf13 View commit details
  3. Merge pull request #2981 from rnlahaye/discrete-trajectory-evaluate

    Optimizations for DiscreteTrajectory
    pleroy authored May 14, 2021
    Copy the full SHA
    74e4b07 View commit details
Showing with 186 additions and 55 deletions.
  1. +103 −22 benchmarks/discrete_trajectory.cpp
  2. +5 −4 physics/discrete_trajectory.hpp
  3. +32 −13 physics/discrete_trajectory_body.hpp
  4. +4 −0 physics/forkable.hpp
  5. +25 −6 physics/forkable_body.hpp
  6. +17 −10 physics/forkable_test.cpp
125 changes: 103 additions & 22 deletions benchmarks/discrete_trajectory.cpp
Original file line number Diff line number Diff line change
@@ -7,29 +7,57 @@
#include "benchmark/benchmark.h"
#include "geometry/frame.hpp"
#include "ksp_plugin/frames.hpp"
#include "quantities/elementary_functions.hpp"
#include "quantities/named_quantities.hpp"
#include "quantities/quantities.hpp"

namespace principia {
namespace physics {

using base::make_not_null_unique;
using base::not_null;
using geometry::Displacement;
using geometry::Frame;
using geometry::Handedness;
using geometry::Inertial;
using geometry::Instant;
using geometry::Velocity;
using ksp_plugin::World;
using quantities::AngularFrequency;
using quantities::Cos;
using quantities::Sin;
using quantities::Time;
using quantities::si::Metre;
using quantities::si::Radian;
using quantities::si::Second;

namespace {

// Creates a trajectory with the given number of steps.
not_null<std::unique_ptr<DiscreteTrajectory<World>>> CreateTrajectory(
// Creates a motionless trajectory with the given number of steps.
not_null<std::unique_ptr<DiscreteTrajectory<World>>> CreateMotionlessTrajectory(
int const steps) {
auto trajectory = make_not_null_unique<DiscreteTrajectory<World>>();
Instant t;
for (int i = 0; i < steps; i++, t += 1 * Second) {
trajectory->Append(t, {World::origin, World::unmoving});
for (Time t = 0 * Second; t < steps * Second; t += 1 * Second) {
trajectory->Append(Instant() + t, {World::origin, World::unmoving});
}
return trajectory;
}

// Creates a circular trajectory with the given number of steps in the XZ plane.
not_null<std::unique_ptr<DiscreteTrajectory<World>>> CreateCircularTrajectory(
int const steps) {
auto trajectory = make_not_null_unique<DiscreteTrajectory<World>>();
constexpr AngularFrequency ω = 1 * Radian / Second;
for (Time t = 0 * Second; t < steps * Second; t += 1 * Second) {
double const sin_ωt = Sin(ω * t);
double const cos_ωt = Cos(ω * t);
trajectory->Append(
Instant() + t,
{World::origin +
Displacement<World>({cos_ωt * Metre, 0 * Metre, sin_ωt * Metre}),
Velocity<World>({-sin_ωt * Metre / Second,
0 * Metre / Second,
cos_ωt * Metre / Second})});
}
return trajectory;
}
@@ -52,7 +80,7 @@ not_null<DiscreteTrajectory<World>*> ForkAt(DiscreteTrajectory<World>& parent,

void BM_DiscreteTrajectoryFront(benchmark::State& state) {
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateTrajectory(4);
CreateMotionlessTrajectory(4);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

@@ -63,7 +91,7 @@ void BM_DiscreteTrajectoryFront(benchmark::State& state) {

void BM_DiscreteTrajectoryBack(benchmark::State& state) {
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateTrajectory(4);
CreateMotionlessTrajectory(4);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

@@ -74,7 +102,7 @@ void BM_DiscreteTrajectoryBack(benchmark::State& state) {

void BM_DiscreteTrajectoryBegin(benchmark::State& state) {
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateTrajectory(4);
CreateMotionlessTrajectory(4);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

@@ -85,7 +113,7 @@ void BM_DiscreteTrajectoryBegin(benchmark::State& state) {

void BM_DiscreteTrajectoryEnd(benchmark::State& state) {
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateTrajectory(4);
CreateMotionlessTrajectory(4);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

@@ -94,78 +122,131 @@ void BM_DiscreteTrajectoryEnd(benchmark::State& state) {
}
}

void BM_DiscreteTrajectoryTMin(benchmark::State& state) {
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateMotionlessTrajectory(4);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

for (auto _ : state) {
fork->t_min();
}
}

void BM_DiscreteTrajectoryTMax(benchmark::State& state) {
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateMotionlessTrajectory(4);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

for (auto _ : state) {
fork->t_max();
}
}

void BM_DiscreteTrajectoryCreateDestroy(benchmark::State& state) {
int const steps = state.range(0);
for (auto _ : state) {
CreateTrajectory(steps);
CreateMotionlessTrajectory(steps);
}
}

void BM_DiscreteTrajectoryIterate(benchmark::State& state) {
int const steps = state.range(0);
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateTrajectory(steps);
CreateMotionlessTrajectory(steps);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

auto const begin = fork->begin();
auto const end = fork->end();
for (auto _ : state) {
for (auto it = fork->begin(); it != fork->end(); ++it) {
for (auto it = begin; it != end; ++it) {
}
}
}

void BM_DiscreteTrajectoryReverseIterate(benchmark::State& state) {
int const steps = state.range(0);
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateTrajectory(steps);
CreateMotionlessTrajectory(steps);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

auto const begin = fork->begin();
auto const end = fork->end();
for (auto _ : state) {
for (auto it = fork->end(); it != fork->begin(); --it) {
for (auto it = end; it != begin; --it) {
}
}
}

void BM_DiscreteTrajectoryFind(benchmark::State& state) {
int const steps = state.range(0);
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateTrajectory(steps);
CreateMotionlessTrajectory(steps);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

for (auto _ : state) {
// These times are in different segments of the fork.
fork->Find(Instant() + 333 * Second);
fork->Find(Instant() + 667 * Second);
fork->Find(Instant() + 833 * Second);
fork->Find(Instant() + 1.0 / 3.0 * steps * Second);
fork->Find(Instant() + 2.0 / 3.0 * steps * Second);
fork->Find(Instant() + 5.0 / 6.0 * steps * Second);
}
}

void BM_DiscreteTrajectoryLowerBound(benchmark::State& state) {
int const steps = state.range(0);
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateTrajectory(steps);
CreateMotionlessTrajectory(steps);
not_null<DiscreteTrajectory<World>*> const fork =
ForkAt(*ForkAt(*trajectory, 0.5), 0.75);

for (auto _ : state) {
// These times are in different segments of the fork.
fork->LowerBound(Instant() + 333 * Second);
fork->LowerBound(Instant() + 667 * Second);
fork->LowerBound(Instant() + 833 * Second);
fork->LowerBound(Instant() + 1.0 / 3.0 * steps * Second);
fork->LowerBound(Instant() + 2.0 / 3.0 * steps * Second);
fork->LowerBound(Instant() + 5.0 / 6.0 * steps * Second);
}
}

void BM_DiscreteTrajectoryEvaluateDegreesOfFreedomExact(
benchmark::State& state) {
int const steps = state.range(0);
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateCircularTrajectory(4);

Instant const t = Instant() + 1 * Second;
for (auto _ : state) {
trajectory->EvaluateDegreesOfFreedom(t);
}
}

void BM_DiscreteTrajectoryEvaluateDegreesOfFreedomInterpolated(
benchmark::State& state) {
not_null<std::unique_ptr<DiscreteTrajectory<World>>> const trajectory =
CreateCircularTrajectory(4);

Instant const t = Instant() + 1.5 * Second;
for (auto _ : state) {
trajectory->EvaluateDegreesOfFreedom(t);
}
}

BENCHMARK(BM_DiscreteTrajectoryFront);
BENCHMARK(BM_DiscreteTrajectoryBack);
BENCHMARK(BM_DiscreteTrajectoryBegin);
BENCHMARK(BM_DiscreteTrajectoryEnd);
BENCHMARK(BM_DiscreteTrajectoryTMin);
BENCHMARK(BM_DiscreteTrajectoryTMax);
BENCHMARK(BM_DiscreteTrajectoryCreateDestroy)->Range(8, 1024);
BENCHMARK(BM_DiscreteTrajectoryIterate)->Range(8, 1024);
BENCHMARK(BM_DiscreteTrajectoryReverseIterate)->Range(8, 1024);
BENCHMARK(BM_DiscreteTrajectoryFind)->Range(8, 1024);
BENCHMARK(BM_DiscreteTrajectoryLowerBound)->Range(8, 1024);
BENCHMARK(BM_DiscreteTrajectoryEvaluateDegreesOfFreedomExact);
BENCHMARK(BM_DiscreteTrajectoryEvaluateDegreesOfFreedomInterpolated);

} // namespace physics
} // namespace principia
9 changes: 5 additions & 4 deletions physics/discrete_trajectory.hpp
Original file line number Diff line number Diff line change
@@ -46,6 +46,9 @@ class DiscreteTrajectoryIterator
DiscreteTrajectoryTraits<Frame>> {
public:
struct reference {
explicit reference(
typename DiscreteTrajectoryTraits<Frame>::TimelineConstIterator it);

Instant const& time;
DegreesOfFreedom<Frame> const& degrees_of_freedom;
};
@@ -256,11 +259,9 @@ class DiscreteTrajectory : public Forkable<DiscreteTrajectory<Frame>,
std::vector<DiscreteTrajectory<Frame>**> const& forks);

// Returns the Hermite interpolation for the left-open, right-closed
// trajectory segment containing the given |time|, or, if |time| is |t_min()|,
// returns a first-degree polynomial which should be evaluated only at
// |t_min()|.
// trajectory segment bounded above by |upper|.
Hermite3<Instant, Position<Frame>> GetInterpolation(
Instant const& time) const;
Iterator const& upper) const;

Timeline timeline_;

45 changes: 32 additions & 13 deletions physics/discrete_trajectory_body.hpp
Original file line number Diff line number Diff line change
@@ -31,18 +31,23 @@ Instant const& DiscreteTrajectoryTraits<Frame>::time(
return it->first;
}

template<typename Frame>
DiscreteTrajectoryIterator<Frame>::reference::reference(
typename DiscreteTrajectoryTraits<Frame>::TimelineConstIterator it)
: time(it->first), degrees_of_freedom(it->second) {}

template<typename Frame>
typename DiscreteTrajectoryIterator<Frame>::reference
DiscreteTrajectoryIterator<Frame>::operator*() const {
auto const& it = this->current();
return {it->first, it->second};
return reference(it);
}

template<typename Frame>
std::optional<typename DiscreteTrajectoryIterator<Frame>::reference>
DiscreteTrajectoryIterator<Frame>::operator->() const {
auto const& it = this->current();
return std::make_optional<reference>({it->first, it->second});
return std::make_optional<reference>(it);
}

template<typename Frame>
@@ -311,19 +316,37 @@ Instant DiscreteTrajectory<Frame>::t_max() const {
template<typename Frame>
Position<Frame> DiscreteTrajectory<Frame>::EvaluatePosition(
Instant const& time) const {
return GetInterpolation(time).Evaluate(time);
auto const iter = this->LowerBound(time);
if (iter->time == time) {
return iter->degrees_of_freedom.position();
}
CHECK_LT(t_min(), time);
CHECK_GT(t_max(), time);
return GetInterpolation(iter).Evaluate(time);
}

template<typename Frame>
Velocity<Frame> DiscreteTrajectory<Frame>::EvaluateVelocity(
Instant const& time) const {;
return GetInterpolation(time).EvaluateDerivative(time);
Instant const& time) const {
auto const iter = this->LowerBound(time);
if (iter->time == time) {
return iter->degrees_of_freedom.velocity();
}
CHECK_LT(t_min(), time);
CHECK_GT(t_max(), time);
return GetInterpolation(iter).EvaluateDerivative(time);
}

template<typename Frame>
DegreesOfFreedom<Frame> DiscreteTrajectory<Frame>::EvaluateDegreesOfFreedom(
Instant const& time) const {
auto const interpolation = GetInterpolation(time);
auto const iter = this->LowerBound(time);
if (iter->time == time) {
return iter->degrees_of_freedom;
}
CHECK_LT(t_min(), time);
CHECK_GT(t_max(), time);
auto const interpolation = GetInterpolation(iter);
return {interpolation.Evaluate(time), interpolation.EvaluateDerivative(time)};
}

@@ -638,13 +661,9 @@ void DiscreteTrajectory<Frame>::FillSubTreeFromMessage(

template<typename Frame>
Hermite3<Instant, Position<Frame>> DiscreteTrajectory<Frame>::GetInterpolation(
Instant const& time) const {
CHECK_LE(t_min(), time);
CHECK_GE(t_max(), time);
// This is the upper bound of the interval upon which we will do the
// interpolation.
auto const upper = this->LowerBound(time);
auto const lower = upper == this->begin() ? upper : --Iterator{upper};
Iterator const& upper) const {
CHECK(upper != this->begin());
auto const lower = --Iterator{upper};
return Hermite3<Instant, Position<Frame>>{
{lower->time, upper->time},
{lower->degrees_of_freedom.position(),
4 changes: 4 additions & 0 deletions physics/forkable.hpp
Original file line number Diff line number Diff line change
@@ -39,6 +39,10 @@ using geometry::Instant;
// TimelineConstIterator must be an STL-like iterator in the timeline of
// Tr4jectory. |time()| must return the corresponding time. Iterators must be
// STL-like and *must*not* be invalidated when the trajectory changes.
//
// The It3rator class must export declarations similar to the following:
// using reference = ...;
// where reference can be constructed from a TimelineConstIterator.

template<typename Tr4jectory, typename It3rator, typename Traits>
class Forkable;
Loading