mp-units icon indicating copy to clipboard operation
mp-units copied to clipboard

Implement simple vector and tensor types

Open mpusz opened this issue 1 year ago • 14 comments

To not depend on external libraries and to provide correct, fast, and simple examples of LA usage with quantities we need simple vector and tensor types (up to 3 dimensions) together with the following operations:

  • a + b - addition where both arguments should be of the same quantity kind and character
  • a - b - subtraction where both arguments should be of the same quantity kind and character
  • a % b - modulo where both arguments should be of the same quantity kind and character
  • a * b - multiplication where one of the arguments has to be a scalar
  • a / b - division where the divisor has to be scalar
  • a ⋅ b - dot product of two vectors
  • a × b - cross product of two vectors
  • |a| - magnitude of a vector
  • a ⊗ b - tensor product of two vectors or tensors
  • a ⋅ b - inner product of two tensors
  • a ⋅ b - inner product of tensor and vector
  • a : b - scalar product of two tensors

mpusz avatar Jul 12 '24 13:07 mpusz

Hi. I have a few questions regarding this issue:

  1. Should tensors of arbitrary rank be supported, or will the rank be capped, since most "normal" physics does not require more than rank-4 tensors? (Even those are rather rare. The only two rank-4 tensors I could think of off the top of my head were the elasticity tensor and the Riemann curvature tensor)
  2. When you speak of inner products between tensors or a tensor and a vector, do you actually mean contractions? So, e.g., if a_ij, and b_jk denote the components of two rank-2 tensors is (a ⋅ b)_ik = Σ_j a_ij b_jk?
  3. What do you mean with "scalar product of two tensors"? I have never heard of such a thing. The way I remember it from my physics studies is that "scalar products" is a synonym for "inner products" which (for non-complex fields) are positive definite, symmetric, bilinear forms. That means, they map two vectors to a positive scalar value, not tensors of other rank. I think I heard people use the term inner product between tensors to talk about what I know as contractions, hence my second question. However, I really don't know what a scalar product of two tensors should then be.

PatrickKa avatar Jul 14 '24 13:07 PatrickKa

Ad 3. Seems like I was a bit off with my terminology. According to Wikipedia, "scalar product" is a synonym for "dot product". This does not help my confusinon though, since the dot product is just a special case of an inner product.

PatrickKa avatar Jul 14 '24 13:07 PatrickKa

Should tensors of arbitrary rank be supported, or will the rank be capped, since most "normal" physics does not require more than rank-4 tensors? (Even those are rather rare. The only two rank-4 tensors I could think of off the top of my head were the elasticity tensor and the Riemann curvature tensor)

Yes, you are right. ISO 80000 part 2 only provides operations on tensors of the second and fourth order.

mpusz avatar Jul 14 '24 15:07 mpusz

When you speak of inner products between tensors or a tensor and a vector, do you actually mean contractions? So, e.g., if a_ij, and b_jk denote the components of two rank-2 tensors is (a ⋅ b)_ik = Σ_j a_ij b_jk?

Yes, it is how ISO 80000-2 defines them. Please contact me via email, and I will send you more details if needed.

mpusz avatar Jul 14 '24 15:07 mpusz

What do you mean with "scalar product of two tensors"?

ISO 80000-2 defines a scalar product of two tensors T and S of the second order as $T : S = \sum_i\sum_jT_{ij} S_{ij}$ and the result is a scalar quantity.

mpusz avatar Jul 14 '24 15:07 mpusz

Thank you for clarifying. It seems like I should have just read ISO 80000-2 😅. I won't have time in the foreseeable future to work on this issue, but if it is still open when I have less on my plate, I'd be happy to help with it.

Also, I didn't know that I can use $\LaTeX$ here. That's cool.

PatrickKa avatar Jul 15 '24 06:07 PatrickKa

I've got a couple of questions/comments on this issue.

To not depend on external libraries and to provide correct, fast, and simple examples of LA usage with quantities we need simple vector and tensor types (up to 3 dimensions) together with the following operations:

Why artificially limit to 3 dimensions?

  • a % b - modulo where both arguments should be of the same quantity kind and character

What does it even mean to take the "modulo" of vectors or tensors?

  • a × b - cross product of two vectors

We could do this, for 3 dimensions only. But the cross product operation is a hacky kludge, mathematically. It's just a disguised wedge product, which would be valid for arbitrary dimensions.

If you're giving different priorities to these operations, then I'd put this one in a lower priority tier.

chiphogg avatar Jul 15 '24 14:07 chiphogg

Why artificially limit to 3 dimensions?

My understanding is that ISQ talks only about Cartesian (orthonormal) coordinates in ordinary space. This even means that we do not have to provide 1- and 2-dimensional vectors to support ISQ. I proposed to support them because it should be quite easy and would allow us to model specific scenarios for some quantities.

What does it even mean to take the "modulo" of vectors or tensors?

You are right; this operation probably makes sense for scalars only.

We could do this, for 3 dimensions only.

And this is why I said that we need 3-dimensional vectors and matrices/tensors.

But the cross product operation is a hacky kludge, mathematically. It's just a disguised wedge product, which would be valid for arbitrary dimensions.

Agree, but this is how ISQ defines its quantities. And we need to support it as well as all others with the same priority.

mpusz avatar Jul 15 '24 14:07 mpusz

To be clear, I do not want to invent a new full-blown LA library here. We need two simple types working on three components (possibly also less) and a bunch of operations on them. That is it. If someone wants to do more advanced stuff, there are always Eigen, Blaze, and other beautiful libraries out there. We will also provide some interop examples with them, if possible.

mpusz avatar Jul 15 '24 14:07 mpusz

Those types are only meant to be used as representation types of vectors and tensor quantities in this library and will mostly be useful for our unit testing and examples. We may also expose them in the library interface to our users, but we will not force their usage. Alternatives should also be easy to use.

Those types will not be used to model huge vectors or matrices of quantities (even if homogeneous). For that, we already have many LA libraries.

mpusz avatar Jul 15 '24 14:07 mpusz

If we're going that route, I suggest not even making the dimensionality configurable at all. Just keep it simple and use 3D vectors.

chiphogg avatar Jul 15 '24 15:07 chiphogg

I'm used to handle 6 degree of freedom (3 distances & 3 orientaton angles) in a 3D carthesian world, hence I'm often use 6D-vectors. E.G. for 'locations' ([m],[m],[m],[rad],[rad],[rad]) and for 'forces' ([N],[N],[N],[Nm],[Nm],[Nm]) etc. At present without any unit control. I realize that this can also be set up differently, e.g. as a struct of 2 3D vectors. But it worries me to have to do this in order to benefit from your great work. Don't let that slow you down, I just wanted to report...

Spammed avatar Jul 25 '24 16:07 Spammed

@Spammed, a full "linear algebra + quantities" solution would be able to handle this. We would need interfaces to provide a separate unit for each index, and facilities for correctly combining them in matrix multiplication.

And it's known how to do this --- but there just isn't an open source library that provides this yet. Daniel Withopf began presenting one solution publicly back in 2021, but it's a closed source library that's internal to Bosch. Also in 2021, I implemented a proof of concept for wrapping Eigen with this kind of interface, and showed zero runtime performance loss. I presented this idea in this segment of my CppCon talk, but I didn't get a chance to implement actual usable end-user interfaces, so this still isn't implemented in Au yet (see aurora-opensource/au#70 to track progress, but know that it will be a while before we start).

All of which is to say: what you want is possible and desirable, even with a single vector (not a struct of 2 3D vectors), but the stopgap solution in this present issue won't give it to you, and I don't know of any open source solutions that are planned with a known timeframe (this is planned for Au, but there is no timeframe, and we surely won't start on it in 2024).

chiphogg avatar Jul 27 '24 15:07 chiphogg

Thank you very much for your response and for addressing and summarizing the problem. I greatly appreciate the work of all of you and know that your joint deliberations are leading step by step to a better future.

Spammed avatar Jul 28 '24 18:07 Spammed

Hi guys! I came across this thread and I’d love to contribute to mp-units as my first issue. This will also be my first contribution to the project, and I’m still learning modern C++, but I’ll do my best and I’m really excited to learn from the process. 🙂

Jullija avatar Oct 07 '25 09:10 Jullija

Welcome @Jullija! It is great to have you on board! 🙌

You can start from the Introduction chapter of our docs. Contributing guide is particularly useful here.

If you encounter any issues or have questions, please don't hesitate to ask.

mpusz avatar Oct 07 '25 17:10 mpusz