Vector and matrix algebra cookbook

Vector and matrix algebra cookbook#

How do I write cross/dot/matmul without CUDA math calls?

You call the ordinary-looking member function directly on the vector or matrix expression — a.cross(b), a.dot(b), M * v, a.norm() — and the library builds the expression tree for you; there is no CUDA math intrinsic to reach for anywhere in this surface, on host or device.

One rank-generic operator set#

aether has no separate “vector algebra” and “matrix algebra” module: every member below is declared once, on aether::Expression itself (aether/expr/Expression.h), and works on ANY conforming expression — a Vec3d, a Mat33d, a batched sample read off a View — because the whole surface is generic over element_extents’s rank rather than hard-coded to “3-vector” or “3x3 matrix”. Writing a + s*b on a Mat33d uses the exact same operator+/operator* as Expression templates showed for a Vec3d.

Vectors#

  • a.cross(b) — cross product (rank-1, 3-component operands only; aether/expr/nodes/Geometric.h).

  • a.dot(b) — inner product, returning a plain scalar (element_type, not an expression) rather than something you assign further.

  • a.norm() / a.squaredNorm() / a.cubedNorm() and their reciprocal counterparts a.rNorm() / a.rSquaredNorm() / a.rCubedNorm() — the reciprocal forms exist because a division is usually more expensive than a multiply, so a formula that needs 1/norm asks for it directly instead of computing norm() and dividing.

  • a.maxNorm() — the largest-magnitude component; a.sum() — the component sum.

  • a.unitVector() — a unit-length copy of a (no in-place normalize(); this is always a new expression).

  • a.segment<Off,Len>() / a.head<N>() / a.tail<N>() — a read-only, lazy view of Len consecutive components starting at Off (aether/expr/nodes/Structural.h). These are READ views only in the current surface — there is no a.segment<Off,Len>() = expr write-back yet.

Matrices#

aether/expr/nodes/Product.h adds the matrix-shaped operations:

  • M * v — matrix-vector product (MatVec); A * B — matrix-matrix product (MatMat); outer(u, v) — outer product (Outer), all reached through natural call-site syntax, not a named method.

  • a.cwiseMul(b) — component-wise (Hadamard) product, rank-generic like every arithmetic operator above.

  • a.matDot(b) — the Frobenius inner product of two same-shaped matrices (sum of the elementwise product), the matrix analogue of dot().

  • a.transpose(), a.row<R>(), a.col<C>(), a.block<R0,C0,BR,BC>() — the same read-only structural slicing segment/head/tail offer for vectors, extended to matrices (aether/expr/nodes/Structural.h).

  • a.trace(), a.det(), a.inverse() — small closed-form matrix algebra (aether/expr/nodes/Inverse.h), deliberately limited to 2x2 and 3x3: there is no general LU decomposition or linear solve in this surface.

An assignment like dest = M * v also gets an internal shortcut: rather than recomputing the whole dot-product chain once per output component, the reused operand (v for M * v, both matrices for A * B) is materialized into a register-resident Item once and every output component reads from that cached copy — a structural detail you never have to spell out yourself, since it happens automatically at the assignment step Expression templates describes.

A runnable example#

The identity y × z = x, z × x = y, x × y = z for the standard basis vectors, checked directly against .cross():

Vec3d yCrossZ = y.cross(z);
Vec3d zCrossX = z.cross(x);
Vec3d xCrossY = x.cross(y);
EXPECT_DOUBLE_EQ(x(0), yCrossZ(0));
EXPECT_DOUBLE_EQ(y(1), zCrossX(1));
EXPECT_DOUBLE_EQ(z(2), xCrossY(2));

Here x, y and z are Vec3d values standing for the unit basis vectors (1,0,0), (0,1,0) and (0,0,1). Quaternions picks up where this page leaves off, for rotating a vector rather than combining two of them.