datenwolf/linmath.h

a lean linear math library, aimed at graphics programming. Supports vec3, vec4, mat4x4 and quaternions

linmath.h -- A small library for linear math as required for computer graphics

linmath.h provides the most used types required for programming computer graphics:

• vec3 -- 3 element vector of floats
• vec4 -- 4 element vector of floats (4th component used for homogenous computations)
• mat4x4 -- 4 by 4 elements matrix, computations are done in column major order
• quat -- quaternion

The types are deliberately named like the types in GLSL. In fact they are meant to be used for the client side computations and passing to same typed GLSL uniforms.

linmath.h's Issues

mat4x4_invert implementation

Hi @datenwolf,

While noodling around on my own C math library, I stumbled across linmath.h and was very interested in your implementation of mat4x4_invert (see below).

linmath.h/linmath.h

Lines 278 to 318 in 0538757

	LINMATH_H_FUNC void mat4x4_invert(mat4x4 T, mat4x4 const M)
	{
	float s[6];
	float c[6];
	s[0] = M[0][0]*M[1][1] - M[1][0]*M[0][1];
	s[1] = M[0][0]*M[1][2] - M[1][0]*M[0][2];
	s[2] = M[0][0]*M[1][3] - M[1][0]*M[0][3];
	s[3] = M[0][1]*M[1][2] - M[1][1]*M[0][2];
	s[4] = M[0][1]*M[1][3] - M[1][1]*M[0][3];
	s[5] = M[0][2]*M[1][3] - M[1][2]*M[0][3];
	c[0] = M[2][0]*M[3][1] - M[3][0]*M[2][1];
	c[1] = M[2][0]*M[3][2] - M[3][0]*M[2][2];
	c[2] = M[2][0]*M[3][3] - M[3][0]*M[2][3];
	c[3] = M[2][1]*M[3][2] - M[3][1]*M[2][2];
	c[4] = M[2][1]*M[3][3] - M[3][1]*M[2][3];
	c[5] = M[2][2]*M[3][3] - M[3][2]*M[2][3];
	/* Assumes it is invertible */
	float idet = 1.0f/( s[0]*c[5]-s[1]*c[4]+s[2]*c[3]+s[3]*c[2]-s[4]*c[1]+s[5]*c[0] );
	T[0][0] = ( M[1][1] * c[5] - M[1][2] * c[4] + M[1][3] * c[3]) * idet;
	T[0][1] = (-M[0][1] * c[5] + M[0][2] * c[4] - M[0][3] * c[3]) * idet;
	T[0][2] = ( M[3][1] * s[5] - M[3][2] * s[4] + M[3][3] * s[3]) * idet;
	T[0][3] = (-M[2][1] * s[5] + M[2][2] * s[4] - M[2][3] * s[3]) * idet;
	T[1][0] = (-M[1][0] * c[5] + M[1][2] * c[2] - M[1][3] * c[1]) * idet;
	T[1][1] = ( M[0][0] * c[5] - M[0][2] * c[2] + M[0][3] * c[1]) * idet;
	T[1][2] = (-M[3][0] * s[5] + M[3][2] * s[2] - M[3][3] * s[1]) * idet;
	T[1][3] = ( M[2][0] * s[5] - M[2][2] * s[2] + M[2][3] * s[1]) * idet;
	T[2][0] = ( M[1][0] * c[4] - M[1][1] * c[2] + M[1][3] * c[0]) * idet;
	T[2][1] = (-M[0][0] * c[4] + M[0][1] * c[2] - M[0][3] * c[0]) * idet;
	T[2][2] = ( M[3][0] * s[4] - M[3][1] * s[2] + M[3][3] * s[0]) * idet;
	T[2][3] = (-M[2][0] * s[4] + M[2][1] * s[2] - M[2][3] * s[0]) * idet;
	T[3][0] = (-M[1][0] * c[3] + M[1][1] * c[1] - M[1][2] * c[0]) * idet;
	T[3][1] = ( M[0][0] * c[3] - M[0][1] * c[1] + M[0][2] * c[0]) * idet;
	T[3][2] = (-M[3][0] * s[3] + M[3][1] * s[1] - M[3][2] * s[0]) * idet;
	T[3][3] = ( M[2][0] * s[3] - M[2][1] * s[1] + M[2][2] * s[0]) * idet;
	}

The implementation is quite a bit more efficient than the more traditional mat4x4 inverse functions I've seen and implemented. Is the implementation above derived from a Guauss-Jordan elimination approach? I'm just curious how you arrived at it and if there are any corner cases to watch out for.

I see it was added as part of this commit - 0f1050f which does mention Guauss-Jordan elimination and iterative Jacobian (which I hadn't heard of before).

Thanks very much for your time,

Tom

quat_mul_vec3 causes SIGSEGV (Address boundary error)

It seems that the quat_mul_vec3 function takes in a vec3 as a result parameter and then passes it on to several other functions expecting a quat, resulting in an Address boundary error.

Readme tyops

"linmath.h provides the most used types required programming computer graphice:"
maybe better:
"linmath.h provides the most used types required for programming computer graphics:"

STOP. CONVERTING. DOUBLES. TO. FLOATS.

So uhh yeah, whenever I try to use linmath, I get 12 warnings all telling me that "Severity Code Description Project File Line Column Suppression State
Warning C4244 'initializing': conversion from 'double' to 'float', possible loss of data". And guess where they are? LinMath! pls fix

indexing error

The indexing for the q variable in quat_from_mat4x4 function is incorrect. I recommend using #defines for the indexes for the whole file, such as:
if(r < 1e-6) {
q[QUAT_S] = 1.f;
q[QUAT_X] = q[QUAT_Y] = q[QUAT_Z] = 0.f;
return;
}

q[QUAT_S] = r/2.f;
q[QUAT_X] = (M[p[0]][p[1]] - M[p[1]][p[0]])/(2.f*r);
q[QUAT_Y] = (M[p[2]][p[0]] - M[p[0]][p[2]])/(2.f*r);
q[QUAT_Z] = (M[p[2]][p[1]] - M[p[1]][p[2]])/(2.f*r);

Great library, by the way.

Bug in mat4x4_perspective

Hello again.

In mat4x4_perspective, you have confused the row/column order for the first four assignments. This probably went undetected because m[0][0] is the same when inverted, and the rest of the assignments are just "=0.0f;"

The issue arises when you supply an undefined matrix as target, i which case random undefined memory will occupy the values instead, leading to strange viewport bugs.

mat4x4_scale_aniso leaves target W vector uninitialized

The mat4x4_scale_aniso function should copy the W vector from the reference matrix to the target matrix. As it stands:

mat4x4 id;
mat4x4_identity(id);

mat4x4 scaled;
mat4x4_scale_aniso(scaled, id, 0.1, 0.4, 0.9); // scales x, y, z, doesn't copy w

Found this while running some really basic tests on the library.

quat_from_mat4 not working

Unfortunately, the function doesn't give correct results.

Here's a proposal to change it (taken from https://www.euclideanspace.com/maths/geometry/rotations/conversions/matrixToQuaternion/index.htm):

float tr = M[0][0] + M[1][1] + M[2][2];

  if (tr > 0) { 
    float S = sqrt(tr+1.0) * 2; // S=4*qw 
    q[3] = 0.25 * S;
    q[0] = (M[1][2] - M[2][1]) / S;
    q[1] = (M[2][0] - M[0][2]) / S; 
    q[2] = (M[0][1] - M[1][0]) / S; 
  } else if ((M[0][0] > M[1][1])&(M[0][0] > M[2][2])) { 
    float S = sqrt(1.0 + M[0][0] - M[1][1] - M[2][2]) * 2; // S=4*qx 
    q[3] = (M[1][2] - M[2][1]) / S;
    q[0] = 0.25 * S;
    q[1] = (M[1][0] + M[0][1]) / S; 
    q[2] = (M[2][0] + M[0][2]) / S; 
  } else if (M[1][1] > M[2][2]) { 
    float S = sqrt(1.0 + M[1][1] - M[0][0] - M[2][2]) * 2; // S=4*qy
    q[3] = (M[2][0] - M[0][2]) / S;
    q[0] = (M[1][0] + M[0][1]) / S; 
    q[1] = 0.25 * S;
    q[2] = (M[2][1] + M[1][2]) / S; 
  } else { 
    float S = sqrt(1.0 + M[2][2] - M[0][0] - M[1][1]) * 2; // S=4*qz
    q[3] = (M[0][1] - M[1][0]) / S;
    q[0] = (M[2][0] + M[0][2]) / S;
    q[1] = (M[2][1] + M[1][2]) / S;
    q[2] = 0.25 * S;
  }

acos() instead of acosf() was used

linmath.h/linmath.h

Line 599 in 0538757

float const angle = acos(vec3_mul_inner(a_, b_)) * s;

Usage of acos() function instead of acosf() makes compiler do a float-to-double conversion and back without needs.

It's kind of weird that it's written in C99

Hi,

nice library! But it's kind of weird that this is written in C99. Makes it quite non idiomatic in *nix terms. It also bloats the tools size, because of all the macro expansions.

warning: pointers to arrays with different qualifiers are incompatible in ISO C

after commit 0538757

linmath.h: In function ‘mat4x4_mul’:
linmath.h:174:16: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
174 | mat4x4_dup(M, temp);
| ^~~~
linmath.h: In function ‘mat4x4_translate_in_place’:
linmath.h:198:17: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
198 | mat4x4_row(r, M, i);
| ^
linmath.h: In function ‘mat4x4_rotate’:
linmath.h:225:19: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
225 | mat4x4_scale(S, S, s);
| ^
linmath.h:229:17: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
229 | mat4x4_sub(C, C, T);
| ^
linmath.h:229:20: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
229 | mat4x4_sub(C, C, T);
| ^
linmath.h:231:19: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
231 | mat4x4_scale(C, C, c);
| ^
linmath.h:233:17: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
233 | mat4x4_add(T, T, C);
| ^
linmath.h:233:20: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
233 | mat4x4_add(T, T, C);
| ^
linmath.h:234:17: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
234 | mat4x4_add(T, T, S);
| ^
linmath.h:234:20: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
234 | mat4x4_add(T, T, S);
| ^
linmath.h:237:20: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
237 | mat4x4_mul(R, M, T);
| ^
linmath.h: In function ‘mat4x4_rotate_X’:
linmath.h:252:19: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
252 | mat4x4_mul(Q, M, R);
| ^
linmath.h: In function ‘mat4x4_rotate_Y’:
linmath.h:264:19: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
264 | mat4x4_mul(Q, M, R);
| ^
linmath.h: In function ‘mat4x4_rotate_Z’:
linmath.h:276:19: warning: pointers to arrays with different qualifiers are incompatible in ISO C [-Wpedantic]
276 | mat4x4_mul(Q, M, R);

warning: return type defaults to int

warning: return type defaults to ‘int’ [-Wimplicit-int]
LINMATH_H_FUNC mat4x4_arcball(mat4x4 R, mat4x4 M, vec2 _a, vec2 _b, float s)

in commit "transplanted arcball function"

static inline void mat4x4_arcball(mat4x4 R, mat4x4 M, vec2 _a, vec2 _b, float s)
was changed to:
LINMATH_H_FUNC mat4x4_arcball(mat4x4 R, mat4x4 M, vec2 _a, vec2 _b, float s)

missing void

LINMATH_H_FUNC void mat4x4_arcball(mat4x4 R, mat4x4 M, vec2 _a, vec2 _b, float s)

Give parameters/variables better names

So for example, in a duplication function, name the params source and dest, instead of M and N.
Or name a column parameter 'column' instead of 'r'.

Reusing variables in quat_mul_vec3 for vector operations produces bad results

The new implementation of quat_mul_vec3 introduced in commit 06b68f6 incorrectly calculates 't' and 'u'. The problem's caused by the variables being used as both input and output to vec3_mul_cross, which causes the variables to be corrupted as the cross product is taken. Unlike with the simple vector operations, input variables to vec3_mul_cross must be distinct from the output variable.

I've suggested a fix for it in PR #22.

static inline void vec3_mul_cross(vec3 r, vec3 const a, vec3 const b) {
    r[0] = a[1]*b[2] - a[2]*b[1];
    r[1] = a[2]*b[0] - a[0]*b[2];
    r[2] = a[0]*b[1] - a[1]*b[0];
}

Needs a test suite

Every function of linmath.h should be complemented by a matching suite of tests. While those can't assure correctness, they can act as a canary toward behavior changes.

Usage double numbers instead of float

By default C compiler treats floating point numbers as double, so "1." is always converted to float, "1.f" should be used instead.

linmath.h/linmath.h

Line 582 in 0538757

z_a = sqrtf(1. - vec2_mul_inner(a, a));

Issue with constructing vectors

Hiya, Loving this module so far, but I have an issue constructing a 3 dimension vector.
I have a large array of floats that represent 3d vertices. I wish to loop over them on a per vertex basis so i constructed the following piece of code.

for(int i=0; i<n; i++) {
   vec3 v = VEC3(array[i+0], array[i+1], array[i+2]);
}

But this gives me the following compile error (using g++)

linmath.h:148:30: error: array must be initialized with a brace-enclosed initializer
  #define VEC3(x,y,z)   (vec3){x,y,z}
                              ^~~~~~~
vertexbuf.h:76:17: note: in expansion of macro ‘VEC3’
       vec3 p1 = VEC3(grid[i+0], grid[i+1], grid[i+2]);

Am I doing something wrong here?

mat4x4_mul_vec4 Correction

I believe Line 184:
r[j] += M[i][j] * v[i];

Should be:
r[j] += M[j][i] * v[i];

Three small issues

1:
mat4x4_col() and mat4x4_row() are defined twice, 115-119, and 266-272, respectively. Interestingly enough, two uses a loop, one uses an unrolled loop, and one uses memcpy().

2:
At 183, you are missing a semicolon.

3:
At 71: vec4_mul_inner( vec4 a, vec4 b ),
"a" and "b" should be const, also you are passing const variables to the function elsewhere, eg. "mat4x4_translate_in_place()".

mat4x4_rotate_z, order of operation intended?

I tried this function with a translation Matrix as input for M. When I use the resulting matrix it first rotates around the origin and after that, it translates, but given the functions description, I would have thought that it would take the transformation and add the operation "rotate around the origin" afterwards.
Is that behaviour intended?

static inline void mat4x4_rotate_Z(mat4x4 Q, mat4x4 M, float angle)
{
	float s = sinf(angle);
	float c = cosf(angle);
	mat4x4 R = {
		{   c,   s, 0.f, 0.f},
		{  -s,   c, 0.f, 0.f},
		{ 0.f, 0.f, 1.f, 0.f},
		{ 0.f, 0.f, 0.f, 1.f}
	};
	mat4x4_mul(Q, M, R);
}

To get this behaviour I expected I had to change :
mat4x4_mul(Q, M, R); to mat4x4_mul(Q, R, M);

Edit: 21.04.21 fix code formatting

Bogus assignment in quat_from_mat4x4 suggests broken numerics

This line very likely is indicative of a deeper lying problem.

linmath.h/linmath.h

Line 558 in 3eef828

m = r;

using mat4x4_scale as a MODEL uniform

can there be a way to use mat4x4_scale to build a scale matrix with the w-component set to 1? that way we can just send it on up to the vertex shader as our scale matrix.

Examples or docs

Hi, thanks for the library! I have searched for some example code using the library (or docs). Does anyone know of anything?

Thanks,
Erik

Usage of For Loops

Have you benchmarked these functions using the existing short-range For loops versus unrolling them out? From the vantage point of the compiled product, using Floops of only 3 and 4 increments seems woefully verbose.

Compiling with Visual Studio 2012

Hi,

I was looking for a small but efficient GLSL friendly C library, and found yours very useful. I just wanted to point out couple of things regarding the compilation with the Visual Studio 2012 Express, as C code:

• 'inline' needs to be substituted with '__inline'
• a couple of matrix declarations just need to be shifted before the first function call take place ('mat4x4 S ...;' and 'mat4x4 C;' in 'mat4x4_rotate' fun)

redundant #include <string.h>

after commit 0538757

redundant #include <string.h>

#include <string.h>
#include <math.h>
#include <string.h>