Support matrices, vectors, and linear algebra

[chiphogg]

This was mentioned in our CppCon 2021 talk, but it hasn't yet actually made its way to the library. This is likely to be frustrating for many who found the library via the talk, so I wanted to give some backstory and share the future plans.

I had worked out the core details from the "units" side as far back as July 2021 (which is why it made it into the talk). However, before implementing a production version, I found Daniel Withopf's work (seen most recently in Daniel's CppCon 2022 talk). Besides the fact that it was essentially an existence proof for a production-quality solution, the most exciting aspect was that it also handled frames of reference. This was consistent with our internal experience that when it comes to matrices and vectors, supporting units is nice to have, but supporting reference frames adds much more value. At this point, reference frame support became an essential feature, which raised the bar for the production API design, to the point where it fell by the wayside.

In the long term, this may be for the best. In that intervening time, I wrote the new foundations and more fluent interfaces for the units library --- essentially, the parts that became Au. The matrix APIs we build on top of this stronger foundation will be much better and more usable than what we would have built in 2021. However, in the short term, it means we don't have any solution to offer, which I know is frustrating.

For clarity: this feature is still considered "in scope" for the library. I hope we will be able to provide an open source matrix solution that handles both units and reference frames, so the world can benefit. However, 2023 is a critical year for Aurora's product roadmap, so it's unlikely to happen in this year.

Meanwhile, if this is a feature you're interested in, please feel free to "vote" by reacting with :+1: to help us gauge community interest!

[catskul]

Is there any advice you can give for people using this in the short term? E.g. is there a sub-standard way to do this that we might implement on top of the base libarary in the mean time?

[chiphogg]

I wish I had a good answer for this. Here are some of the best ideas I have right now.

1. Use a strong convention. This only works for some use cases, but it's what we do internally at Aurora, for example. We use the convention that vectors are always expressed in base SI units. Since it's so universal, it basically never leads to confusion. While it's less "slick", and you do need to do manual conversions between scalar and vector types, it actually does work better than I would have guessed! Of course, for users whose applications aren't well suited to base SI units, this is a non-starter.

2. Use wg21/linear_algebra. I got this from the corresponding mp-units webpage https://mpusz.github.io/units/use_cases/linear_algebra.html. I haven't actually used the library myself. However, since they have some examples on the mp-units webpage --- and, since it shouldn't matter which units library we're using --- then I would expect this to work out of the box. Of course, it would be a "linear-algebra-on-units" solution, so it would never support mixed-unit types.

3. Use Eigen and submit patches. We used a "patched Eigen" approach for years at Uber ATG. It wasn't too bad. The approach is to use a fork of Eigen instead of using "base" Eigen directly. Then, try writing what you want to write, and whenever you get a compiler error with Eigen, figure out what assumption it's making and fix it in your local fork. Since it's in an open source context, it would be good to go the extra mile and submit those patches as upstream PRs to Eigen. In this case, the fork wouldn't be intended to diverge permanently; you'd hope the patches get accepted so that over time, more and more of Eigen could work out of the box. Again, this is a "linear-algebra-on-units" approach, which has shortfalls, but can be fine for people who don't need mixed-unit containers.

I think this thread at mp-units https://github.com/mpusz/units/issues/301#issuecomment-943664661 may also be relevant.

Sorry I don't have a better answer for you! I would really love to dig into this task, but we're in the home stretch this year for building Aurora's first product, and that's taking up my mental energy right now.

On Tue, Apr 18, 2023 at 4:57 PM catskul @.***> wrote:

Is there any advice you can give for people using this in the short term? E.g. is there a sub-standard way to do this that we might implement on top of the base libarary in the mean time?

— Reply to this email directly, view it on GitHub https://github.com/aurora-opensource/au/issues/70#issuecomment-1513785714, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAN4IYAIFVI7ZPXKFAZGRR3XB3525ANCNFSM6AAAAAAT6BTXQ4 . You are receiving this because you authored the thread.Message ID: @.***>

[chiphogg]

(Also: I'm on vacation and I don't have my work laptop, so if the formatting isn't the best in these responses, then that's why!)

On Mon, Apr 24, 2023 at 4:15 PM Charles Hogg @.***> wrote:

I wish I had a good answer for this. Here are some of the best ideas I have right now.

1. Use a strong convention. This only works for some use cases, but it's what we do internally at Aurora, for example. We use the convention that vectors are always expressed in base SI units. Since it's so universal, it basically never leads to confusion. While it's less "slick", and you do need to do manual conversions between scalar and vector types, it actually does work better than I would have guessed! Of course, for users whose applications aren't well suited to base SI units, this is a non-starter.

2. Use wg21/linear_algebra. I got this from the corresponding mp-units webpage https://mpusz.github.io/units/use_cases/linear_algebra.html. I haven't actually used the library myself. However, since they have some examples on the mp-units webpage --- and, since it shouldn't matter which units library we're using --- then I would expect this to work out of the box. Of course, it would be a "linear-algebra-on-units" solution, so it would never support mixed-unit types.

3. Use Eigen and submit patches. We used a "patched Eigen" approach for years at Uber ATG. It wasn't too bad. The approach is to use a fork of Eigen instead of using "base" Eigen directly. Then, try writing what you want to write, and whenever you get a compiler error with Eigen, figure out what assumption it's making and fix it in your local fork. Since it's in an open source context, it would be good to go the extra mile and submit those patches as upstream PRs to Eigen. In this case, the fork wouldn't be intended to diverge permanently; you'd hope the patches get accepted so that over time, more and more of Eigen could work out of the box. Again, this is a "linear-algebra-on-units" approach, which has shortfalls, but can be fine for people who don't need mixed-unit containers.

I think this thread at mp-units https://github.com/mpusz/units/issues/301#issuecomment-943664661 may also be relevant.

Sorry I don't have a better answer for you! I would really love to dig into this task, but we're in the home stretch this year for building Aurora's first product, and that's taking up my mental energy right now.

On Tue, Apr 18, 2023 at 4:57 PM catskul @.***> wrote:

Is there any advice you can give for people using this in the short term? E.g. is there a sub-standard way to do this that we might implement on top of the base libarary in the mean time?

— Reply to this email directly, view it on GitHub https://github.com/aurora-opensource/au/issues/70#issuecomment-1513785714, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAN4IYAIFVI7ZPXKFAZGRR3XB3525ANCNFSM6AAAAAAT6BTXQ4 . You are receiving this because you authored the thread.Message ID: @.***>

[connorjak]

For my application (aircraft / spacecraft dynamics & systems modeling), this feature would be most useful as a utility to do conversions of homogenous-unit vector and matrix quantities, retaining the dimensional analysis of the vector/matrix as a whole instead of getting weighted down by individually assigning units to indices. If I need heterogenous-unit data structures, I'm probably best served by making a struct.

Typical use cases:

1. Data parsing program intaking scalars, vectors, and matrices from a file as double precision values with an arbitrary mix of units, then passing that data to various math models that use SI.
2. Realtime or post-sim analysis software converting units from SI to a user-selected unit.

The above use-cases seem like they're most easily accomplished by 3 (insert units in Eigen template types).

There are a few notable things I have run into in implementing arbitrary transformations between Reference Frames. Representing the state of a body (location, rotation, velocity, ang velocity, acceleration, ang accel, etc.) requires two defined frames: "With Respect To" and "Expressed In". "With Respect To" defines the observer, and "Expressed In" defines the basis vectors of the coordinate system that is used to resolve the physical quantity into 3 axes. Transforming an acceleration from an initial assumed frame to a given "with respect to" target frame requires performing 2nd derivative Transport Theorem through all intermediary frames in a frame hierarchy to retain correctness through rotating/accelerating frames. Then, transforming to an "expressed in" frame requires a rotation of the vector through all those intermediary frames. There is a large amount of dynamics math that we did in getting this right; automated reference frame transformations may be out of scope for Au.

[chiphogg]

Thanks for this additional input on use cases!

For my application (aircraft / spacecraft dynamics & systems modeling), this feature would be most useful as a utility to do conversions of homogenous-unit vector and matrix quantities, retaining the dimensional analysis of the vector/matrix as a whole instead of getting weighted down by individually assigning units to indices. If I need heterogenous-unit data structures, I'm probably best served by making a struct.

I hope the APIs for heterogeneous units wouldn't be very burdensome. Ideally, it'd be possible to define it once in a readable way with a concise name, and then just use that name in most callsites.

That said, of course the homogeneous case should be much simpler. My hope is that we'll end up being able to call auto x = meters(Eigen::Vector3d{1., 2., 3.}), and also retrieve the underlying vector with x.in(meters). Before I could contemplate doing that, I'd need to make the QuantityMaker APIs a little more generic, and I'd have to make sure that extra indirection wouldn't hurt compile times too much.

Typical use cases:

1. Data parsing program intaking scalars, vectors, and matrices from a file as double precision values with an arbitrary mix of units, then passing that data to various math models that use SI.
2. Realtime or post-sim analysis software converting units from SI to a user-selected unit.

The above use-cases seem like they're most easily accomplished by 3 (insert units in Eigen template types).

The good news about that is that you're not blocked on me! :grin: I'd suggest simply doing this, then seeing what breaks, figuring out how to patch Eigen, and putting together a PR. It shouldn't take too many iterations to get something that's reasonably usable for simple use cases like this.

There are a few notable things I have run into in implementing arbitrary transformations between Reference Frames. Representing the state of a body (location, rotation, velocity, ang velocity, acceleration, ang accel, etc.) requires two defined frames: "With Respect To" and "Expressed In". "With Respect To" defines the observer, and "Expressed In" defines the basis vectors of the coordinate system that is used to resolve the physical quantity into 3 axes. Transforming an acceleration from an initial assumed frame to a given "with respect to" target frame requires performing 2nd derivative Transport Theorem through all intermediary frames in a frame hierarchy to retain correctness through rotating/accelerating frames. Then, transforming to an "expressed in" frame requires a rotation of the vector through all those intermediary frames. There is a large amount of dynamics math that we did in getting this right; automated reference frame transformations may be out of scope for Au.

I agree this'd be out of scope for Au, but the frame annotations would still make it a lot easier to make a robust solution. You'd need to encode your transformation in an object or function that takes two frame parameters. The details of that transformation would be whatever they need to be, but by declaring the source and destination frames, you could get compile-time checking to make sure you're applying it to objects in the correct frames.

[catskul]

Just wanted to check in to see if anything might be likely to change in the short or medium term w.r.t. availability of matrix/linear algebra support.

[chiphogg]

It's painful for me to say this, but: no, we're not likely to be able to work on this in the short or medium term. Reasons include:

• high amount of design effort needed for interfaces
• high amount of effort on implementation and testing
• inability to get the practical production experience we'd need to validate this

I'll elaborate on that last part a bit. We want to make sure the things we add are really useful and usable in production code. Adopting a feature like this would change the way we use our APIs far more than adopting other Au features. Right now, where we're at as a company is strongly focusing on delivering our Aurora Driver Ready milestone, and preparing for the subsequent commercial launch in 2024. We won't be able to afford the distraction and opportunity cost of changing our interfaces to try out a promising but untested new approach.

It's definitely a bummer --- I would really enjoy being able to dedicate to this feature the time and attention it needs in order to be good. It's an itch I've been longing to scratch for a long time.

We do take ":+1:" votes into account in choosing what to work on, but we also take into account the size and complexity of the task. That's why the second-place issue, #43 (and parts of #90 it depends on), are likely to get done far sooner. It was simple enough that I was able to do basically all of the hard parts on the plane during my recent business travel. This is true for the rest of milestone 0.4.0 as well: every individual issue tagged there is small and self contained, and thus easy to fit in.

The best stopgap strategy IMO is still to use a fork of Eigen (say), so that you can patch issues as they come up. Based on our experience at ATG, the first issues that come up tend to be the most impactful, so fixing them can quickly get many use cases working. Plus, only one project really has to do this if they're submitting their changes upstream, because if Eigen accepts them, then everyone will be able to benefit.

[catskul]

Thank you for the in-depth reply and the advice. Your effort and responsiveness is impressive and appreciated.

[mpusz]

I hope I do not overstep here, as @chiphogg is working with me on a proposal to standardize the physical quantities and units library in C++29 based on mp-units.

In the upcoming releases, we will work on improving the linear algebra support for vector and tensor quantities. A short summary can be found here: https://github.com/mpusz/mp-units/issues/301#issuecomment-2225025754 and more details here:

• https://github.com/mpusz/mp-units/issues/463
• https://github.com/mpusz/mp-units/issues/439

Please provide feedback and describe your intended use cases so we can better adjust the design to your needs.

[chiphogg]

No overstepping, @mpusz --- you're always welcome to share news and insights from your library here. :slightly_smiling_face: I'm looking forward to seeing how these ideas turn out in practice.