Eigen: 3D vector kinematics¶
A body has a position, a velocity, and an acceleration, each a 3D vector. Where is it one timestep later, and how far from the origin is it?
// Position must be meters, velocity m/s, acceleration m/s^2, and the timestep seconds.
Eigen::Vector3d advanced_position(const Eigen::Vector3d &x_m,
const Eigen::Vector3d &v_mps,
const Eigen::Vector3d &a_mps2,
double dt_s) {
return x_m + v_mps * dt_s + 0.5 * a_mps2 * dt_s * dt_s;
}
int main() {
const Eigen::Vector3d x_m {0.0, 0.0, 120.0};
// The velocity arrives as 72 km/h downrange, so keep a second copy of it, scaled into m/s.
const Eigen::Vector3d v_kph {72.0, 0.0, 0.0};
const Eigen::Vector3d v_mps = v_kph * (1000.0 / 3600.0);
const Eigen::Vector3d a_mps2 {0.0, 0.0, -9.80665};
// The timestep is 250 ms; convert that by hand too.
const Eigen::Vector3d x_new_m = advanced_position(x_m, v_mps, a_mps2, 250.0 / 1000.0);
std::cout << x_new_m.transpose() << " m" << '\n'; // Unit label typed by hand.
std::cout << x_new_m.norm() << " m" << '\n'; // ...and again, nothing checks it.
}
Includes and usings
#include "au/au.hh"
#include <iostream>
#include "Eigen/Core"
#include "au/compatibility/eigen.hh"
#include "au/io.hh"
#include "au/units/hours.hh"
#include "au/constants/standard_gravity.hh"
#include "au/units/meters.hh"
#include "au/units/seconds.hh"
// This is a `.cc` file, so we import the names we use, one at a time. See the "Namespaces and
// includes" discussion page for why we do this rather than `using namespace au;`.
using au::kilo;
using au::Meters;
using au::milli;
using au::STANDARD_GRAVITY;
using au::norm;
using au::Quantity;
using au::QuantityD;
using au::Seconds;
using au::transpose;
using au::UnitPower;
using au::UnitQuotient;
using au::symbols::h;
using au::symbols::m;
using au::symbols::s;
// Symbols for the prefixed units we use. A prefix applier turns an existing symbol into one for
// the prefixed unit, which is the most readable of the three ways to spell this.
constexpr auto km = kilo(m);
constexpr auto ms = milli(s);
// Aliases for the vector quantity types, so the signature below reads well. A type alias
// introduces one name we chose, so it is fine at namespace scope even in a header.
using Position = Quantity<Meters, Eigen::Vector3d>;
using Velocity = Quantity<UnitQuotient<Meters, Seconds>, Eigen::Vector3d>;
using Acceleration = Quantity<UnitQuotient<Meters, UnitPower<Seconds, 2>>, Eigen::Vector3d>;
// The types state the units. Nothing to remember; nothing to convert.
Position advanced_position(const Position &x,
const Velocity &v,
const Acceleration &a,
QuantityD<Seconds> dt) {
return x + v * dt + 0.5 * a * dt * dt;
}
int main() {
const auto x = Eigen::Vector3d{0.0, 0.0, 120.0} * m;
// Any units of the right dimension will do: the conversion is generated at compile time.
const auto v = Eigen::Vector3d{72.0, 0.0, 0.0} * km / h;
const auto a = Eigen::Vector3d{0.0, 0.0, -1.0} * STANDARD_GRAVITY;
const Position x_new = advanced_position(x, v, a, 250.0 * ms);
std::cout << transpose(x_new) << '\n';
std::cout << norm(x_new) << '\n';
}
Note
The two tabs are aligned for comparison: blank lines where one version needs fewer statements, and extra spaces so that corresponding expressions sit in the same column. Neither is a spelling we’d recommend writing — they’re here so that flipping between the tabs shows only the real differences.
Both programs print the same two lines1:
What’s happening¶
The physics is the same on both sides: x + v \, \Delta t + \frac{1}{2} a \, \Delta t^2. What changes is who is responsible for the units. The raw code has three different unit conversions:
km/htom/sfor the velocity,mstosfor the timestep, andg_0tom/s^2for the acceleration.
The first two are at least visible in the source. The third one isn’t, because the raw version
doesn’t convert anything. Instead, it uses -9.80665 directly as a magic number, leaving users to
guess the intent.
All of those conversions vanish from the source code in the Au version, because the library automatically generates the correct conversion factors — at compile time.
The names simplify, too: unit-suffixed names such as x_m, which force the human to keep track of
the units, get replaced by the simpler x. In fact, for velocity, we get even more
simplification: both v_mps and v_kph get replaced by a single v. Its units happen to be km
/ h, but we don’t need to worry about that; we know the library will produce any necessary
conversions. These simpler names really pay off in the advanced_position() function body: when
the suffixes vanish, the underlying physics shows through more clearly.
The output lines simplify for the same reason. The raw version types the unit label by hand —
<< " m", twice, with nothing checking that it still matches what the number means. The Au version
streams the quantities themselves, and the label comes from the type: norm(x_new) is a length, so
it prints m, and transpose(x_new) is a whole vector, so it prints Eigen’s formatting of the
elements followed by the one unit they all share.
Unit symbols and constants¶
This example leans on Unit symbols, such as m and km, and Constants, such as
STANDARD_GRAVITY. They both have the same effect here: when you multiply or divide by them,
they change the units, but not the underlying stored value. If the input is already
a Quantity, you get another Quantity; and if it’s not, then it becomes one.
Unit symbols are a handy, concise way to annotate your variables with their units. Writing 12.34
* m / s has exactly the same effect as (meters / second)(12.34); it’s just a little shorter.
Again, keep in mind that unit symbols and constants do not change the underlying value. So, if
you’re following our Eigen safety guide, these do not count as “operations” that create risk for
dangling references. That’s why the variable assignments here are perfectly safe, even without
eval().
More nuance on lifetime risk
To be clear: we mean that multiplying by symbols or constants doesn’t add lifetime risk. We don’t mean unit symbols and constants preclude lifetime risk. If there is pre-existing lifetime risk, these won’t magically remove it.
Consider this example. Suppose we have two utility functions that return Eigen::Vector3d
instances:
The following example is guaranteed to dangle:
The sum holds references to its operands, which in this case are the temporary objects v1()
and v2(), neither of which survives past the semicolon at the end of the line. m doesn’t
make this safe, but it’s also not the root of the problem: the following simpler example is also
guaranteed to dangle!
We hope this discussion clarifies how unit symbols and constants relate to lifetime risk. For a fuller treatment of this topic, we recommend that all users read and understand the Eigen safety guide before using Au with Eigen.
Alias names¶
In this example, we went with Position, Velocity, and Acceleration. This is a fine approach,
but not the only possible one. If you want to use mixed units in your interfaces, you could also
define a custom “rep-named alias” for Eigen::Vector3d:
Then you could write QuantityV3<Meters> for a position, QuantityV3<KilometersPerHour> (after
defining a suitable KilometersPerHour alias), and so on.
This question is mostly a matter of taste; Au is safe either way.
Eigen safety¶
Eigen’s famously fast performance comes in part from lazy evaluation. The equally famous cost of this speed is an elevated risk of object lifetime bugs. We have a whole Eigen safety guide devoted to this topic in general. We’ll hit the highlights relevant to this example here.
First, it’s important to appreciate that Au has “risk parity” with Eigen. This means that when you
add Au to Eigen, you still have all the same risks, but you don’t get new ones. The safety
guide explains the details, but the upshot is that you should still be looking for the same warning
signs as raw Eigen, and you’ll still use the same strategies to mitigate the issues (even if some of
the particulars might change, such as eval() being a free function instead of a member function).
In this example, there are three places where expression templates occur, and thus three places that may carry object lifetime risk.
-
The arguments
vandathat we pass toadvanced_position()both need unit conversions.- These are safe because they’re assigned to a concrete type:
VelocityandAcceleration, respectively.
- These are safe because they’re assigned to a concrete type:
-
The return value of
advanced_position(),x + v * dt + 0.5 * a * dt * dt, is also an expression template.- This is safe because the return value is a concrete type:
Position. Evaluation happens when we convert to the concrete type.
- This is safe because the return value is a concrete type:
-
The
transpose(x_new)that we stream on the last line is an expression template too: Eigen’s “view” functions are operations, just like arithmetic.- This one is not assigned to a concrete type. It’s safe for the other reason: we consume it
inside the same full expression, so
x_newcannot have died or changed in the meantime.
- This one is not assigned to a concrete type. It’s safe for the other reason: we consume it
inside the same full expression, so
So, using concrete types for the input parameters and the return value automatically guarantees lifetime safety at those boundaries — and an expression you compute and consume on the spot is safe, because nothing gets deferred past its inputs.
Summary¶
Au’s Eigen support makes it easier to get your units right robustly, and often makes your code easier to read. Lifetime safety is the one thing it leaves exactly as it found it, for better and for worse, so make sure you’re familiar with the Eigen safety guide before you start using Au with Eigen.
Related reading¶
- Eigen how-to guide, for creating and using Eigen-backed quantities.
- Eigen safety, on expression templates and object lifetime.
- Eigen compatibility reference, for the full list of free functions.
- Unit symbols, including the prefix-applier form used for
km. - Constants, such as the
STANDARD_GRAVITYused here. - Element access, for reading and writing one component of a vector quantity.
-
All of our examples get compiled and run in CI, and we check that they produce the same output. ↩