ericlengyel / terathon-math-library

C++ math library for 2D/3D/4D vector, matrix, quaternion, and geometric algebra.

Terathon Math Library

This is a C++ math library containing classes for vectors, matrices, quaternions, and elements of projective geometric algebra. The specific classes are the following:

• Vector2D – A 2D vector (x, y) that extends to four dimensions as (x, y, 0, 0).
• Vector3D – A 3D vector (x, y, z) that extends to four dimensions as (x, y, z, 0).
• Vector4D – A 4D vector (x, y, z, w).
• Point2D – A 2D point (x, y) that extends to four dimensions as (x, y, 0, 1).
• Point3D – A 3D point (x, y, z) that extends to four dimensions as (x, y, z, 1).
• Bivector3D – A 3D bivector x e₂₃ + y e₃₁ + z e₁₂.
• Matrix2D – A 2×2 matrix.
• Matrix3D – A 3×3 matrix.
• Matrix4D – A 4×4 matrix.
• Transform2D – A 3×3 matrix with fourth row always (0, 0, 1).
• Transform3D – A 4×4 matrix with fourth row always (0, 0, 0, 1).
• Quaternion – A conventional quaternion xi + yj + zk + w.
• DualNum – A dual number s + tε.

2D rigid geometric algebra

• FlatPoint2D – A 2D flat point x e₁ + y e₂ + z e₃.
• Line2D – A 2D line x e₂₃ + y e₃₁ + z e₁₂.
• Motor2D – A 2D motion operator Q_x e₁ + Q_y e₂ + Q_z e₃ + Q_w 𝟙.
• Flector2D - A 2D reflection operator F_x e₂₃ + F_y e₃₁ + F_z e₁₂ + F_w 1.

3D rigid geometric algebra

• FlatPoint3D – A 3D flat point x e₁ + y e₂ + z e₃ + w e₄.
• Line3D – A 3D line l_vx e₄₁ + l_vy e₄₂ + l_vz e₄₃ + l_mx e₂₃ + l_my e₃₁ + l_mz e₁₂.
• Plane3D – A 3D plane x e₂₃₄ + y e₃₁₄ + z e₁₂₄ + w e₃₂₁.
• Motor3D – A 3D motion operator Q_vx e₄₁ + Q_vy e₄₂ + Q_vz e₄₃ + Q_vw 𝟙 + Q_mx e₂₃ + Q_my e₃₁ + Q_mz e₁₂ + Q_mw 1.
• Flector3D - A 3D reflection operator F_px e₁ + F_py e₂ + F_pz e₃ + F_pw e₄ + F_gx e₄₂₃ + F_gy e₄₃₁ + F_gz e₄₁₂ + F_gw e₃₂₁.

2D conformal geometric algebra

• RoundPoint2D – A 2D round point x e₁ + y e₂ + z e₃ + w e₄.
• Dipole2D – A 2D dipole d_gx e₂₃ + d_gy e₃₁ + d_gz e₁₂ + d_px e₄₁ + d_py e₄₂ + d_pz e₄₃.
• Circle2D – A 2D circle w e₃₂₁ + x e₄₂₃ + y e₄₃₁ + z e₄₁₂.

3D conformal geometric algebra

• RoundPoint3D – A 3D round point x e₁ + y e₂ + z e₃ + w e₄ + u e₅.
• Dipole3D – A 3D dipole d_vx e₄₁ + d_vy e₄₂ + d_vz e₄₃ + d_mx e₂₃ + d_my e₃₁ + d_mz e₁₂ + d_px e₁₅ + d_py e₂₅ + d_pz e₃₅ + d_pw e₄₅.
• Circle3D – A 3D circle c_gx e₄₂₃ + c_gy e₄₃₁ + c_gz e₄₁₂ + c_gw e₃₂₁ + c_vx e₄₁₅ + c_vy e₄₂₅ + c_vz e₄₃₅ + c_mx e₂₃₅ + c_my e₃₁₅ + c_mz e₁₂₅.
• Sphere3D – A 3D sphere u e₁₂₃₄ + x e₄₂₃₅ + y e₄₃₁₅ + z e₄₁₂₅ + w e₃₂₁₅.

Component Swizzling

Vector components can be swizzled using shading-language syntax. As an example, the following expressions are all valid for a Vector3D object v:

• v.x – The x component of v.
• v.xy – A 2D vector having the x and y components of v.
• v.yzx – A 3D vector having the components of v in the order (y, z, x).

Support for repeated components in a swizzle can be enabled by defining TERATHON_SWIZZLE_REPEAT. This is disabled by default because the large number of additional swizzling possibilities increases compile times substantially. Swizzles with repeated components are always const so that it's not possible to assign to them.

Rows, columns, and submatrices can be extracted from matrix objects using a similar syntax. As an example, the following expressions are all valid for a Matrix3D object m:

• m.m12 – The (1,2) entry of m.
• m.row0 – The first row of m.
• m.col1 – The second column of m.
• m.matrix2D – The upper-left 2×2 submatrix of m.
• m.transpose – The transpose of m.

All of the above are generally free operations, with no copying, when their results are consumed by an expression. For more information, see Eric Lengyel's 2018 GDC talk Linear Algebra Upgraded.

Geometric Algebra

The ^ operator is overloaded for cases in which the wedge or antiwedge product can be applied between vectors, bivectors, flat points, lines, planes, round points, dipoles, circles, and spheres. (Note that ^ has lower precedence than just about everything else, so parentheses will be necessary.)

The library does not provide operators that directly calculate the geometric product and antiproduct because they would tend to generate inefficient code and produce intermediate results having unnecessary types when something like the sandwich product Q ⟇ p ⟇ ~Q appears in an expression. Instead, there are Transform() functions that take some object p for the first parameter and the motor Q with which to transform it for the second parameter.

See Eric Lengyel's Projective Geometric Algebra website for more information about operations among these types.

API Documentation

There is API documentation embedded in the header files. The formatted equivalent can be found in the C4 Engine documentation.

Licensing

Separate proprietary licenses are available from Terathon Software. Please send an email with details about your particular use case if you are interested.

terathon-math-library's Issues

wrong formula in doc for Vector2D Reject

I think the formula in the documentation is wrong.

It should not be (a - a.b)b, but rather a - (a.b)b

Terathon-Math-Library/TSVector2D.h

Line 124 in ed1cbdf

//# Returns (<b>a</b>&#x202F;&minus;&#x202F;<b>a</b>&#x202F;&sdot;&#x202F;<b>b</b>)<b>b</b>, which is the rejection of $a$ from $b$ under the assumption that the magnitude of $b$ is one.

The actual code seems fine, though:

Terathon-Math-Library/TSVector2D.h

Line 399 in ed1cbdf

return (a - b * Dot(a, b));

TSMath.h is missing

the TSMath.h header seems to be missing, it's being called from TSAlgebra.h , Thank you

Unit tests

Any plans on adding unit tests to the library? I think it could prevent issues like these:

#5

TSBasic.h is missing.

Looks like TSBasic.h is missing (included in TSMatrix3D.cpp)

add topic tags

I suggest adding topics such as geometric-algebra, quaternion in the About section

Terathon::Arctan() too imprecise

Hi,

We (Masterbrain Bytes/Cattle and Crops) have camera code that uses polar coordinates to handle input, and we also want the camera to move (and turn) with the vehicle the player is currently driving.
So we're using code like

float zrot = Arctan( vehicleTransform.m10, vehicleTransform.m00 );
azimuth += zrot;

This generally works, but we noticed that when just driving in a circle, at certain angles the camera stuttered.
When replacing Arctan() with atan2f() those problems went away.

This could've been it, but this issue somehow nerdsniped me so I started investigating and found out that Arctan(float x) has errors of up to 0.004459 (compared to atan(x) in double precision).

When looking up how a fast atan() can be implemented, I stumbled upon the idea of just doing a linear interpolation between values of a lookup table.
So I implemented it (with the existing lookup-table), just interpolating between table[i] and table[i+1], using b - i as the lerp factor, and to my big surprise this works really well - the maximal error seems to be below 0.000020, i.e. over 200 times lower.
In addition to that, it even seems to be about twice as fast, and it works with our camera code.

I tested this extensively both with values between +/- 10'000 (in steps of 0.00001) and values between +/- 1.5 (in steps of 0.00000001) and it indeed seems to work quite well.

Here's my test program (tested on both Linux and Windows), including both your and my Arctan(float) implementation: atantest-upd.cpp.txt (Github doesn't allow attaching .cpp files so I had to append .txt).
UPDATE: Replaced the file with an updated one with one more element in table; also added function with 130 element lookup table.
Feel free to use it (also in your commercial code).

Here's my test results from both Linux (on a AMD Ryzen 2700X) and Win10 (on an Intel Haswell i7-4771): atanresults.txt

Cheers,
Daniel