From 118980e02e59ff31871df59dce257075394f3533 Mon Sep 17 00:00:00 2001 From: Aaditya Dhruv Date: Sun, 25 Jan 2026 15:10:37 -0600 Subject: wip --- include/cglm/simd/sse2/affine.h | 115 ++++++++ include/cglm/simd/sse2/mat2.h | 48 ++++ include/cglm/simd/sse2/mat3.h | 76 ++++++ include/cglm/simd/sse2/mat4.h | 573 ++++++++++++++++++++++++++++++++++++++++ include/cglm/simd/sse2/quat.h | 54 ++++ 5 files changed, 866 insertions(+) create mode 100644 include/cglm/simd/sse2/affine.h create mode 100644 include/cglm/simd/sse2/mat2.h create mode 100644 include/cglm/simd/sse2/mat3.h create mode 100644 include/cglm/simd/sse2/mat4.h create mode 100644 include/cglm/simd/sse2/quat.h (limited to 'include/cglm/simd/sse2') diff --git a/include/cglm/simd/sse2/affine.h b/include/cglm/simd/sse2/affine.h new file mode 100644 index 0000000..0619995 --- /dev/null +++ b/include/cglm/simd/sse2/affine.h @@ -0,0 +1,115 @@ +/* + * Copyright (c), Recep Aslantas. + * + * MIT License (MIT), http://opensource.org/licenses/MIT + * Full license can be found in the LICENSE file + */ + +#ifndef cglm_affine_mat_sse2_h +#define cglm_affine_mat_sse2_h +#if defined( __SSE__ ) || defined( __SSE2__ ) + +#include "../../common.h" +#include "../intrin.h" + +CGLM_INLINE +void +glm_mul_sse2(mat4 m1, mat4 m2, mat4 dest) { + /* D = R * L (Column-Major) */ + glmm_128 l, r0, r1, r2, r3, v0, v1, v2, v3; + + l = glmm_load(m1[0]); + r0 = glmm_load(m2[0]); + r1 = glmm_load(m2[1]); + r2 = glmm_load(m2[2]); + r3 = glmm_load(m2[3]); + + v0 = _mm_mul_ps(glmm_splat_x(r0), l); + v1 = _mm_mul_ps(glmm_splat_x(r1), l); + v2 = _mm_mul_ps(glmm_splat_x(r2), l); + v3 = _mm_mul_ps(glmm_splat_x(r3), l); + + l = glmm_load(m1[1]); + v0 = glmm_fmadd(glmm_splat_y(r0), l, v0); + v1 = glmm_fmadd(glmm_splat_y(r1), l, v1); + v2 = glmm_fmadd(glmm_splat_y(r2), l, v2); + v3 = glmm_fmadd(glmm_splat_y(r3), l, v3); + + l = glmm_load(m1[2]); + v0 = glmm_fmadd(glmm_splat_z(r0), l, v0); + v1 = glmm_fmadd(glmm_splat_z(r1), l, v1); + v2 = glmm_fmadd(glmm_splat_z(r2), l, v2); + v3 = glmm_fmadd(glmm_splat_z(r3), l, v3); + + l = glmm_load(m1[3]); + v3 = glmm_fmadd(glmm_splat_w(r3), l, v3); + + glmm_store(dest[0], v0); + glmm_store(dest[1], v1); + glmm_store(dest[2], v2); + glmm_store(dest[3], v3); +} + +CGLM_INLINE +void +glm_mul_rot_sse2(mat4 m1, mat4 m2, mat4 dest) { + /* D = R * L (Column-Major) */ + + glmm_128 l, r0, r1, r2, v0, v1, v2; + + l = glmm_load(m1[0]); + r0 = glmm_load(m2[0]); + r1 = glmm_load(m2[1]); + r2 = glmm_load(m2[2]); + + v0 = _mm_mul_ps(glmm_splat_x(r0), l); + v1 = _mm_mul_ps(glmm_splat_x(r1), l); + v2 = _mm_mul_ps(glmm_splat_x(r2), l); + + l = glmm_load(m1[1]); + v0 = glmm_fmadd(glmm_splat_y(r0), l, v0); + v1 = glmm_fmadd(glmm_splat_y(r1), l, v1); + v2 = glmm_fmadd(glmm_splat_y(r2), l, v2); + + l = glmm_load(m1[2]); + v0 = glmm_fmadd(glmm_splat_z(r0), l, v0); + v1 = glmm_fmadd(glmm_splat_z(r1), l, v1); + v2 = glmm_fmadd(glmm_splat_z(r2), l, v2); + + glmm_store(dest[0], v0); + glmm_store(dest[1], v1); + glmm_store(dest[2], v2); + glmm_store(dest[3], glmm_load(m1[3])); +} + +CGLM_INLINE +void +glm_inv_tr_sse2(mat4 mat) { + __m128 r0, r1, r2, r3, x0, x1, x2, x3, x4, x5; + + r0 = glmm_load(mat[0]); + r1 = glmm_load(mat[1]); + r2 = glmm_load(mat[2]); + r3 = glmm_load(mat[3]); + x1 = _mm_set_ps(1.0f, 0.0f, 0.0f, 0.0f); + + _MM_TRANSPOSE4_PS(r0, r1, r2, x1); + + x2 = glmm_shuff1(r3, 0, 0, 0, 0); + x3 = glmm_shuff1(r3, 1, 1, 1, 1); + x4 = glmm_shuff1(r3, 2, 2, 2, 2); + x5 = glmm_float32x4_SIGNMASK_NEG; + + x0 = glmm_fmadd(r0, x2, glmm_fmadd(r1, x3, _mm_mul_ps(r2, x4))); + x0 = _mm_xor_ps(x0, x5); + + x0 = _mm_add_ps(x0, x1); + + glmm_store(mat[0], r0); + glmm_store(mat[1], r1); + glmm_store(mat[2], r2); + glmm_store(mat[3], x0); +} + +#endif +#endif /* cglm_affine_mat_sse2_h */ diff --git a/include/cglm/simd/sse2/mat2.h b/include/cglm/simd/sse2/mat2.h new file mode 100644 index 0000000..31b3a29 --- /dev/null +++ b/include/cglm/simd/sse2/mat2.h @@ -0,0 +1,48 @@ +/* + * Copyright (c), Recep Aslantas. + * + * MIT License (MIT), http://opensource.org/licenses/MIT + * Full license can be found in the LICENSE file + */ + +#ifndef cglm_mat2_sse_h +#define cglm_mat2_sse_h +#if defined( __SSE__ ) || defined( __SSE2__ ) + +#include "../../common.h" +#include "../intrin.h" + +CGLM_INLINE +void +glm_mat2_mul_sse2(mat2 m1, mat2 m2, mat2 dest) { + __m128 x0, x1, x2, x3, x4; + + x1 = glmm_load(m1[0]); /* d c b a */ + x2 = glmm_load(m2[0]); /* h g f e */ + + x3 = glmm_shuff1(x2, 2, 2, 0, 0); + x4 = glmm_shuff1(x2, 3, 3, 1, 1); + x0 = _mm_movelh_ps(x1, x1); + x2 = _mm_movehl_ps(x1, x1); + + /* + dest[0][0] = a * e + c * f; + dest[0][1] = b * e + d * f; + dest[1][0] = a * g + c * h; + dest[1][1] = b * g + d * h; + */ + x0 = glmm_fmadd(x0, x3, _mm_mul_ps(x2, x4)); + + glmm_store(dest[0], x0); +} + +CGLM_INLINE +void +glm_mat2_transp_sse2(mat2 m, mat2 dest) { + /* d c b a */ + /* d b c a */ + glmm_store(dest[0], glmm_shuff1(glmm_load(m[0]), 3, 1, 2, 0)); +} + +#endif +#endif /* cglm_mat2_sse_h */ diff --git a/include/cglm/simd/sse2/mat3.h b/include/cglm/simd/sse2/mat3.h new file mode 100644 index 0000000..f07320c --- /dev/null +++ b/include/cglm/simd/sse2/mat3.h @@ -0,0 +1,76 @@ +/* + * Copyright (c), Recep Aslantas. + * + * MIT License (MIT), http://opensource.org/licenses/MIT + * Full license can be found in the LICENSE file + */ + +#ifndef cglm_mat3_sse_h +#define cglm_mat3_sse_h +#if defined( __SSE__ ) || defined( __SSE2__ ) + +#include "../../common.h" +#include "../intrin.h" + +CGLM_INLINE +void +glm_mat3_mul_sse2(mat3 m1, mat3 m2, mat3 dest) { + __m128 l0, l1, l2, r0, r1, r2, x0, x1, x2, x3, x4, x5, x6, x7, x8, x9; + + l0 = _mm_loadu_ps(m1[0]); + l1 = _mm_loadu_ps(&m1[1][1]); + + r0 = _mm_loadu_ps(m2[0]); + r1 = _mm_loadu_ps(&m2[1][1]); + + x8 = glmm_shuff1(l0, 0, 2, 1, 0); /* a00 a02 a01 a00 */ + x1 = glmm_shuff1(r0, 3, 0, 0, 0); /* b10 b00 b00 b00 */ + x2 = _mm_shuffle_ps(l0, l1, _MM_SHUFFLE(1, 0, 3, 3)); /* a12 a11 a10 a10 */ + x3 = _mm_shuffle_ps(r0, r1, _MM_SHUFFLE(2, 0, 3, 1)); /* b20 b11 b10 b01 */ + x0 = _mm_mul_ps(x8, x1); + + x6 = glmm_shuff1(l0, 1, 0, 2, 1); /* a01 a00 a02 a01 */ + x7 = glmm_shuff1(x3, 3, 3, 1, 1); /* b20 b20 b10 b10 */ + l2 = _mm_load_ss(&m1[2][2]); + r2 = _mm_load_ss(&m2[2][2]); + x1 = _mm_mul_ps(x6, x7); + l2 = glmm_shuff1(l2, 0, 0, 1, 0); /* a22 a22 0.f a22 */ + r2 = glmm_shuff1(r2, 0, 0, 1, 0); /* b22 b22 0.f b22 */ + + x4 = glmm_shuff1(x2, 0, 3, 2, 0); /* a10 a12 a11 a10 */ + x5 = glmm_shuff1(x2, 2, 0, 3, 2); /* a11 a10 a12 a11 */ + x6 = glmm_shuff1(x3, 2, 0, 0, 0); /* b11 b01 b01 b01 */ + x2 = glmm_shuff1(r1, 3, 3, 0, 0); /* b21 b21 b11 b11 */ + + x8 = _mm_unpackhi_ps(x8, x4); /* a10 a00 a12 a02 */ + x9 = _mm_unpackhi_ps(x7, x2); /* b21 b20 b21 b20 */ + + x0 = glmm_fmadd(x4, x6, x0); + x1 = glmm_fmadd(x5, x2, x1); + + x2 = _mm_movehl_ps(l2, l1); /* a22 a22 a21 a20 */ + x3 = glmm_shuff1(x2, 0, 2, 1, 0); /* a20 a22 a21 a20 */ + x2 = glmm_shuff1(x2, 1, 0, 2, 1); /* a21 a20 a22 a21 */ + x4 = _mm_shuffle_ps(r0, r1, _MM_SHUFFLE(1, 1, 2, 2)); /* b12 b12 b02 b02 */ + + x5 = glmm_shuff1(x4, 3, 0, 0, 0); /* b12 b02 b02 b02 */ + x4 = _mm_movehl_ps(r2, x4); /* b22 b22 b12 b12 */ + x0 = glmm_fmadd(x3, x5, x0); + x1 = glmm_fmadd(x2, x4, x1); + + /* + Dot Product : dest[2][2] = a02 * b20 + + a12 * b21 + + a22 * b22 + + 0 * 00 */ + x2 = _mm_movelh_ps(x8, l2); /* 0.f a22 a12 a02 */ + x3 = _mm_movelh_ps(x9, r2); /* 0.f b22 b21 b20 */ + x2 = glmm_vdots(x2, x3); + + _mm_storeu_ps(&dest[0][0], x0); + _mm_storeu_ps(&dest[1][1], x1); + _mm_store_ss (&dest[2][2], x2); +} + +#endif +#endif /* cglm_mat3_sse_h */ diff --git a/include/cglm/simd/sse2/mat4.h b/include/cglm/simd/sse2/mat4.h new file mode 100644 index 0000000..2127e72 --- /dev/null +++ b/include/cglm/simd/sse2/mat4.h @@ -0,0 +1,573 @@ +/* + * Copyright (c), Recep Aslantas. + * + * MIT License (MIT), http://opensource.org/licenses/MIT + * Full license can be found in the LICENSE file + */ + +#ifndef cglm_mat_sse_h +#define cglm_mat_sse_h +#if defined( __SSE__ ) || defined( __SSE2__ ) + +#include "../../common.h" +#include "../intrin.h" + +#define glm_mat4_inv_precise_sse2(mat, dest) glm_mat4_inv_sse2(mat, dest) + +CGLM_INLINE +void +glm_mat4_scale_sse2(mat4 m, float s) { + __m128 x0; + x0 = glmm_set1(s); + + glmm_store(m[0], _mm_mul_ps(glmm_load(m[0]), x0)); + glmm_store(m[1], _mm_mul_ps(glmm_load(m[1]), x0)); + glmm_store(m[2], _mm_mul_ps(glmm_load(m[2]), x0)); + glmm_store(m[3], _mm_mul_ps(glmm_load(m[3]), x0)); +} + +CGLM_INLINE +void +glm_mat4_transp_sse2(mat4 m, mat4 dest) { + __m128 r0, r1, r2, r3; + + r0 = glmm_load(m[0]); + r1 = glmm_load(m[1]); + r2 = glmm_load(m[2]); + r3 = glmm_load(m[3]); + + _MM_TRANSPOSE4_PS(r0, r1, r2, r3); + + glmm_store(dest[0], r0); + glmm_store(dest[1], r1); + glmm_store(dest[2], r2); + glmm_store(dest[3], r3); +} + +CGLM_INLINE +void +glm_mat4_mul_sse2(mat4 m1, mat4 m2, mat4 dest) { + /* D = R * L (Column-Major) */ + + glmm_128 l, r0, r1, r2, r3, v0, v1, v2, v3; + + l = glmm_load(m1[0]); + r0 = glmm_load(m2[0]); + r1 = glmm_load(m2[1]); + r2 = glmm_load(m2[2]); + r3 = glmm_load(m2[3]); + + v0 = _mm_mul_ps(glmm_splat_x(r0), l); + v1 = _mm_mul_ps(glmm_splat_x(r1), l); + v2 = _mm_mul_ps(glmm_splat_x(r2), l); + v3 = _mm_mul_ps(glmm_splat_x(r3), l); + + l = glmm_load(m1[1]); + v0 = glmm_fmadd(glmm_splat_y(r0), l, v0); + v1 = glmm_fmadd(glmm_splat_y(r1), l, v1); + v2 = glmm_fmadd(glmm_splat_y(r2), l, v2); + v3 = glmm_fmadd(glmm_splat_y(r3), l, v3); + + l = glmm_load(m1[2]); + v0 = glmm_fmadd(glmm_splat_z(r0), l, v0); + v1 = glmm_fmadd(glmm_splat_z(r1), l, v1); + v2 = glmm_fmadd(glmm_splat_z(r2), l, v2); + v3 = glmm_fmadd(glmm_splat_z(r3), l, v3); + + l = glmm_load(m1[3]); + v0 = glmm_fmadd(glmm_splat_w(r0), l, v0); + v1 = glmm_fmadd(glmm_splat_w(r1), l, v1); + v2 = glmm_fmadd(glmm_splat_w(r2), l, v2); + v3 = glmm_fmadd(glmm_splat_w(r3), l, v3); + + glmm_store(dest[0], v0); + glmm_store(dest[1], v1); + glmm_store(dest[2], v2); + glmm_store(dest[3], v3); +} + +CGLM_INLINE +void +glm_mat4_mulv_sse2(mat4 m, vec4 v, vec4 dest) { + __m128 x0, x1, m0, m1, m2, m3, v0, v1, v2, v3; + + m0 = glmm_load(m[0]); + m1 = glmm_load(m[1]); + m2 = glmm_load(m[2]); + m3 = glmm_load(m[3]); + + x0 = glmm_load(v); + v0 = glmm_splat_x(x0); + v1 = glmm_splat_y(x0); + v2 = glmm_splat_z(x0); + v3 = glmm_splat_w(x0); + + x1 = _mm_mul_ps(m3, v3); + x1 = glmm_fmadd(m2, v2, x1); + x1 = glmm_fmadd(m1, v1, x1); + x1 = glmm_fmadd(m0, v0, x1); + + glmm_store(dest, x1); +} + +CGLM_INLINE +float +glm_mat4_det_sse2(mat4 mat) { + __m128 r0, r1, r2, r3, x0, x1, x2; + + /* 127 <- 0, [square] det(A) = det(At) */ + r0 = glmm_load(mat[0]); /* d c b a */ + r1 = glmm_load(mat[1]); /* h g f e */ + r2 = glmm_load(mat[2]); /* l k j i */ + r3 = glmm_load(mat[3]); /* p o n m */ + + /* + t[1] = j * p - n * l; + t[2] = j * o - n * k; + t[3] = i * p - m * l; + t[4] = i * o - m * k; + */ + x0 = glmm_fnmadd(glmm_shuff1(r3, 0, 0, 1, 1), glmm_shuff1(r2, 2, 3, 2, 3), + _mm_mul_ps(glmm_shuff1(r2, 0, 0, 1, 1), + glmm_shuff1(r3, 2, 3, 2, 3))); + /* + t[0] = k * p - o * l; + t[0] = k * p - o * l; + t[5] = i * n - m * j; + t[5] = i * n - m * j; + */ + x1 = glmm_fnmadd(glmm_shuff1(r3, 0, 0, 2, 2), glmm_shuff1(r2, 1, 1, 3, 3), + _mm_mul_ps(glmm_shuff1(r2, 0, 0, 2, 2), + glmm_shuff1(r3, 1, 1, 3, 3))); + + /* + a * (f * t[0] - g * t[1] + h * t[2]) + - b * (e * t[0] - g * t[3] + h * t[4]) + + c * (e * t[1] - f * t[3] + h * t[5]) + - d * (e * t[2] - f * t[4] + g * t[5]) + */ + x2 = glmm_fnmadd(glmm_shuff1(r1, 1, 1, 2, 2), glmm_shuff1(x0, 3, 2, 2, 0), + _mm_mul_ps(glmm_shuff1(r1, 0, 0, 0, 1), + _mm_shuffle_ps(x1, x0, _MM_SHUFFLE(1, 0, 0, 0)))); + x2 = glmm_fmadd(glmm_shuff1(r1, 2, 3, 3, 3), + _mm_shuffle_ps(x0, x1, _MM_SHUFFLE(2, 2, 3, 1)), + x2); + + x2 = _mm_xor_ps(x2, glmm_float32x4_SIGNMASK_NPNP); + + return glmm_hadd(_mm_mul_ps(x2, r0)); +} + +CGLM_INLINE +void +glm_mat4_inv_fast_sse2(mat4 mat, mat4 dest) { + __m128 r0, r1, r2, r3, + v0, v1, v2, v3, + t0, t1, t2, t3, t4, t5, + x0, x1, x2, x3, x4, x5, x6, x7, x8, x9; + + /* x8 = _mm_set_ps(-0.f, 0.f, -0.f, 0.f); */ + x8 = glmm_float32x4_SIGNMASK_NPNP; + x9 = glmm_shuff1(x8, 2, 1, 2, 1); + + /* 127 <- 0 */ + r0 = glmm_load(mat[0]); /* d c b a */ + r1 = glmm_load(mat[1]); /* h g f e */ + r2 = glmm_load(mat[2]); /* l k j i */ + r3 = glmm_load(mat[3]); /* p o n m */ + + x0 = _mm_movehl_ps(r3, r2); /* p o l k */ + x3 = _mm_movelh_ps(r2, r3); /* n m j i */ + x1 = glmm_shuff1(x0, 1, 3, 3 ,3); /* l p p p */ + x2 = glmm_shuff1(x0, 0, 2, 2, 2); /* k o o o */ + x4 = glmm_shuff1(x3, 1, 3, 3, 3); /* j n n n */ + x7 = glmm_shuff1(x3, 0, 2, 2, 2); /* i m m m */ + + x6 = _mm_shuffle_ps(r2, r1, _MM_SHUFFLE(0, 0, 0, 0)); /* e e i i */ + x5 = _mm_shuffle_ps(r2, r1, _MM_SHUFFLE(1, 1, 1, 1)); /* f f j j */ + x3 = _mm_shuffle_ps(r2, r1, _MM_SHUFFLE(2, 2, 2, 2)); /* g g k k */ + x0 = _mm_shuffle_ps(r2, r1, _MM_SHUFFLE(3, 3, 3, 3)); /* h h l l */ + + t0 = _mm_mul_ps(x3, x1); + t1 = _mm_mul_ps(x5, x1); + t2 = _mm_mul_ps(x5, x2); + t3 = _mm_mul_ps(x6, x1); + t4 = _mm_mul_ps(x6, x2); + t5 = _mm_mul_ps(x6, x4); + + /* t1[0] = k * p - o * l; + t1[0] = k * p - o * l; + t2[0] = g * p - o * h; + t3[0] = g * l - k * h; */ + t0 = glmm_fnmadd(x2, x0, t0); + + /* t1[1] = j * p - n * l; + t1[1] = j * p - n * l; + t2[1] = f * p - n * h; + t3[1] = f * l - j * h; */ + t1 = glmm_fnmadd(x4, x0, t1); + + /* t1[2] = j * o - n * k + t1[2] = j * o - n * k; + t2[2] = f * o - n * g; + t3[2] = f * k - j * g; */ + t2 = glmm_fnmadd(x4, x3, t2); + + /* t1[3] = i * p - m * l; + t1[3] = i * p - m * l; + t2[3] = e * p - m * h; + t3[3] = e * l - i * h; */ + t3 = glmm_fnmadd(x7, x0, t3); + + /* t1[4] = i * o - m * k; + t1[4] = i * o - m * k; + t2[4] = e * o - m * g; + t3[4] = e * k - i * g; */ + t4 = glmm_fnmadd(x7, x3, t4); + + /* t1[5] = i * n - m * j; + t1[5] = i * n - m * j; + t2[5] = e * n - m * f; + t3[5] = e * j - i * f; */ + t5 = glmm_fnmadd(x7, x5, t5); + + x4 = _mm_movelh_ps(r0, r1); /* f e b a */ + x5 = _mm_movehl_ps(r1, r0); /* h g d c */ + + x0 = glmm_shuff1(x4, 0, 0, 0, 2); /* a a a e */ + x1 = glmm_shuff1(x4, 1, 1, 1, 3); /* b b b f */ + x2 = glmm_shuff1(x5, 0, 0, 0, 2); /* c c c g */ + x3 = glmm_shuff1(x5, 1, 1, 1, 3); /* d d d h */ + + v2 = _mm_mul_ps(x0, t1); + v1 = _mm_mul_ps(x0, t0); + v3 = _mm_mul_ps(x0, t2); + v0 = _mm_mul_ps(x1, t0); + + v2 = glmm_fnmadd(x1, t3, v2); + v3 = glmm_fnmadd(x1, t4, v3); + v0 = glmm_fnmadd(x2, t1, v0); + v1 = glmm_fnmadd(x2, t3, v1); + + v3 = glmm_fmadd(x2, t5, v3); + v0 = glmm_fmadd(x3, t2, v0); + v2 = glmm_fmadd(x3, t5, v2); + v1 = glmm_fmadd(x3, t4, v1); + + /* + dest[0][0] = f * t1[0] - g * t1[1] + h * t1[2]; + dest[0][1] =-(b * t1[0] - c * t1[1] + d * t1[2]); + dest[0][2] = b * t2[0] - c * t2[1] + d * t2[2]; + dest[0][3] =-(b * t3[0] - c * t3[1] + d * t3[2]); */ + v0 = _mm_xor_ps(v0, x8); + + /* + dest[2][0] = e * t1[1] - f * t1[3] + h * t1[5]; + dest[2][1] =-(a * t1[1] - b * t1[3] + d * t1[5]); + dest[2][2] = a * t2[1] - b * t2[3] + d * t2[5]; + dest[2][3] =-(a * t3[1] - b * t3[3] + d * t3[5]);*/ + v2 = _mm_xor_ps(v2, x8); + + /* + dest[1][0] =-(e * t1[0] - g * t1[3] + h * t1[4]); + dest[1][1] = a * t1[0] - c * t1[3] + d * t1[4]; + dest[1][2] =-(a * t2[0] - c * t2[3] + d * t2[4]); + dest[1][3] = a * t3[0] - c * t3[3] + d * t3[4]; */ + v1 = _mm_xor_ps(v1, x9); + + /* + dest[3][0] =-(e * t1[2] - f * t1[4] + g * t1[5]); + dest[3][1] = a * t1[2] - b * t1[4] + c * t1[5]; + dest[3][2] =-(a * t2[2] - b * t2[4] + c * t2[5]); + dest[3][3] = a * t3[2] - b * t3[4] + c * t3[5]; */ + v3 = _mm_xor_ps(v3, x9); + + /* determinant */ + x0 = _mm_shuffle_ps(v0, v1, _MM_SHUFFLE(0, 0, 0, 0)); + x1 = _mm_shuffle_ps(v2, v3, _MM_SHUFFLE(0, 0, 0, 0)); + x0 = _mm_shuffle_ps(x0, x1, _MM_SHUFFLE(2, 0, 2, 0)); + + x0 = _mm_rcp_ps(glmm_vhadd(_mm_mul_ps(x0, r0))); + + glmm_store(dest[0], _mm_mul_ps(v0, x0)); + glmm_store(dest[1], _mm_mul_ps(v1, x0)); + glmm_store(dest[2], _mm_mul_ps(v2, x0)); + glmm_store(dest[3], _mm_mul_ps(v3, x0)); +} + +/* old one */ +#if 0 +CGLM_INLINE +void +glm_mat4_inv_sse2(mat4 mat, mat4 dest) { + __m128 r0, r1, r2, r3, + v0, v1, v2, v3, + t0, t1, t2, t3, t4, t5, + x0, x1, x2, x3, x4, x5, x6, x7, x8, x9; + + /* x8 = _mm_set_ps(-0.f, 0.f, -0.f, 0.f); */ + x8 = glmm_float32x4_SIGNMASK_NPNP; + x9 = glmm_shuff1(x8, 2, 1, 2, 1); + + /* 127 <- 0 */ + r0 = glmm_load(mat[0]); /* d c b a */ + r1 = glmm_load(mat[1]); /* h g f e */ + r2 = glmm_load(mat[2]); /* l k j i */ + r3 = glmm_load(mat[3]); /* p o n m */ + + x0 = _mm_movehl_ps(r3, r2); /* p o l k */ + x3 = _mm_movelh_ps(r2, r3); /* n m j i */ + x1 = glmm_shuff1(x0, 1, 3, 3 ,3); /* l p p p */ + x2 = glmm_shuff1(x0, 0, 2, 2, 2); /* k o o o */ + x4 = glmm_shuff1(x3, 1, 3, 3, 3); /* j n n n */ + x7 = glmm_shuff1(x3, 0, 2, 2, 2); /* i m m m */ + + x6 = _mm_shuffle_ps(r2, r1, _MM_SHUFFLE(0, 0, 0, 0)); /* e e i i */ + x5 = _mm_shuffle_ps(r2, r1, _MM_SHUFFLE(1, 1, 1, 1)); /* f f j j */ + x3 = _mm_shuffle_ps(r2, r1, _MM_SHUFFLE(2, 2, 2, 2)); /* g g k k */ + x0 = _mm_shuffle_ps(r2, r1, _MM_SHUFFLE(3, 3, 3, 3)); /* h h l l */ + + t0 = _mm_mul_ps(x3, x1); + t1 = _mm_mul_ps(x5, x1); + t2 = _mm_mul_ps(x5, x2); + t3 = _mm_mul_ps(x6, x1); + t4 = _mm_mul_ps(x6, x2); + t5 = _mm_mul_ps(x6, x4); + + /* t1[0] = k * p - o * l; + t1[0] = k * p - o * l; + t2[0] = g * p - o * h; + t3[0] = g * l - k * h; */ + t0 = glmm_fnmadd(x2, x0, t0); + + /* t1[1] = j * p - n * l; + t1[1] = j * p - n * l; + t2[1] = f * p - n * h; + t3[1] = f * l - j * h; */ + t1 = glmm_fnmadd(x4, x0, t1); + + /* t1[2] = j * o - n * k + t1[2] = j * o - n * k; + t2[2] = f * o - n * g; + t3[2] = f * k - j * g; */ + t2 = glmm_fnmadd(x4, x3, t2); + + /* t1[3] = i * p - m * l; + t1[3] = i * p - m * l; + t2[3] = e * p - m * h; + t3[3] = e * l - i * h; */ + t3 = glmm_fnmadd(x7, x0, t3); + + /* t1[4] = i * o - m * k; + t1[4] = i * o - m * k; + t2[4] = e * o - m * g; + t3[4] = e * k - i * g; */ + t4 = glmm_fnmadd(x7, x3, t4); + + /* t1[5] = i * n - m * j; + t1[5] = i * n - m * j; + t2[5] = e * n - m * f; + t3[5] = e * j - i * f; */ + t5 = glmm_fnmadd(x7, x5, t5); + + x4 = _mm_movelh_ps(r0, r1); /* f e b a */ + x5 = _mm_movehl_ps(r1, r0); /* h g d c */ + + x0 = glmm_shuff1(x4, 0, 0, 0, 2); /* a a a e */ + x1 = glmm_shuff1(x4, 1, 1, 1, 3); /* b b b f */ + x2 = glmm_shuff1(x5, 0, 0, 0, 2); /* c c c g */ + x3 = glmm_shuff1(x5, 1, 1, 1, 3); /* d d d h */ + + v2 = _mm_mul_ps(x0, t1); + v1 = _mm_mul_ps(x0, t0); + v3 = _mm_mul_ps(x0, t2); + v0 = _mm_mul_ps(x1, t0); + + v2 = glmm_fnmadd(x1, t3, v2); + v3 = glmm_fnmadd(x1, t4, v3); + v0 = glmm_fnmadd(x2, t1, v0); + v1 = glmm_fnmadd(x2, t3, v1); + + v3 = glmm_fmadd(x2, t5, v3); + v0 = glmm_fmadd(x3, t2, v0); + v2 = glmm_fmadd(x3, t5, v2); + v1 = glmm_fmadd(x3, t4, v1); + + /* + dest[0][0] = f * t1[0] - g * t1[1] + h * t1[2]; + dest[0][1] =-(b * t1[0] - c * t1[1] + d * t1[2]); + dest[0][2] = b * t2[0] - c * t2[1] + d * t2[2]; + dest[0][3] =-(b * t3[0] - c * t3[1] + d * t3[2]); */ + v0 = _mm_xor_ps(v0, x8); + + /* + dest[2][0] = e * t1[1] - f * t1[3] + h * t1[5]; + dest[2][1] =-(a * t1[1] - b * t1[3] + d * t1[5]); + dest[2][2] = a * t2[1] - b * t2[3] + d * t2[5]; + dest[2][3] =-(a * t3[1] - b * t3[3] + d * t3[5]);*/ + v2 = _mm_xor_ps(v2, x8); + + /* + dest[1][0] =-(e * t1[0] - g * t1[3] + h * t1[4]); + dest[1][1] = a * t1[0] - c * t1[3] + d * t1[4]; + dest[1][2] =-(a * t2[0] - c * t2[3] + d * t2[4]); + dest[1][3] = a * t3[0] - c * t3[3] + d * t3[4]; */ + v1 = _mm_xor_ps(v1, x9); + + /* + dest[3][0] =-(e * t1[2] - f * t1[4] + g * t1[5]); + dest[3][1] = a * t1[2] - b * t1[4] + c * t1[5]; + dest[3][2] =-(a * t2[2] - b * t2[4] + c * t2[5]); + dest[3][3] = a * t3[2] - b * t3[4] + c * t3[5]; */ + v3 = _mm_xor_ps(v3, x9); + + /* determinant */ + x0 = _mm_shuffle_ps(v0, v1, _MM_SHUFFLE(0, 0, 0, 0)); + x1 = _mm_shuffle_ps(v2, v3, _MM_SHUFFLE(0, 0, 0, 0)); + x0 = _mm_shuffle_ps(x0, x1, _MM_SHUFFLE(2, 0, 2, 0)); + + x0 = _mm_div_ps(glmm_set1(1.0f), glmm_vhadd(_mm_mul_ps(x0, r0))); + + glmm_store(dest[0], _mm_mul_ps(v0, x0)); + glmm_store(dest[1], _mm_mul_ps(v1, x0)); + glmm_store(dest[2], _mm_mul_ps(v2, x0)); + glmm_store(dest[3], _mm_mul_ps(v3, x0)); +} +#endif + +CGLM_INLINE +void +glm_mat4_inv_sse2(mat4 mat, mat4 dest) { + __m128 r0, r1, r2, r3, s1, s2, + v0, v1, v2, v3, v4, v5, + t0, t1, t2, + x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x13; + + /* s1 = _mm_set_ps(-0.f, 0.f, -0.f, 0.f); */ + s1 = glmm_float32x4_SIGNMASK_NPNP; + s2 = glmm_shuff1(s1, 2, 1, 2, 1); + + /* 127 <- 0 */ + r1 = glmm_load(mat[1]); /* h g f e */ + r0 = glmm_load(mat[0]); /* d c b a */ + r3 = glmm_load(mat[3]); /* p o n m */ + r2 = glmm_load(mat[2]); /* l k j i */ + + x4 = _mm_unpackhi_ps(r0, r2); /* l d k c */ + x5 = _mm_unpacklo_ps(r0, r2); /* j b i a */ + x6 = _mm_unpackhi_ps(r1, r3); /* p h o g */ + x7 = _mm_unpacklo_ps(r1, r3); /* n f m e */ + + x0 = _mm_unpackhi_ps(x7, x5); /* j n b f */ + x1 = _mm_unpacklo_ps(x7, x5); /* i m a e */ + x2 = _mm_unpackhi_ps(x6, x4); /* l p d h */ + x3 = _mm_unpacklo_ps(x6, x4); /* k o c g */ + + /* c2 = c * h - d * g c12 = a * g - c * e c8 = a * f - b * e + c1 = k * p - l * o c11 = i * o - k * m c7 = i * n - j * m + c4 = a * h - d * e c6 = b * h - d * f c10 = b * g - c * f + c3 = i * p - l * m c5 = j * p - l * n c9 = j * o - k * n */ + + x8 = _mm_shuffle_ps(x0, x3, _MM_SHUFFLE(3, 1, 3, 1)); /* k c j b */ + x9 = _mm_shuffle_ps(x0, x3, _MM_SHUFFLE(2, 0, 2, 0)); /* o g n f */ + + x10 = glmm_shuff1(x2, 2, 0, 2, 0); /* p h p h */ + x11 = glmm_shuff1(x2, 3, 1, 3, 1); /* l d l d */ + +#if 0 /* TODO measure both */ + x12 = _mm_shuffle_ps(x4, x5, _MM_SHUFFLE(1, 0, 1, 0)); /* i a k c */ + x13 = _mm_shuffle_ps(x6, x7, _MM_SHUFFLE(1, 0, 1, 0)); /* m e o g */ +#else + x12 = _mm_movelh_ps(x4, x5); /* i a k c */ + x13 = _mm_movelh_ps(x6, x7); /* m e o g */ +#endif + + t0 = _mm_mul_ps(x12, x10); + t1 = _mm_mul_ps(x5, x6); + t2 = _mm_mul_ps(x5, x9); + + t0 = glmm_fnmadd(x11, x13, t0); + t1 = glmm_fnmadd(x4, x7, t1); + t2 = glmm_fnmadd(x8, x7, t2); + + /* det */ + /* v0: c3 * c10 + c4 * c9 + c1 * c8 + c2 * c7 */ + /* v1: c5 * c12 + c6 * c11 */ + + v5 = glmm_set1_rval(1.0f); + v0 = glmm_shuff1(t2, 2, 3, 0, 1); + v1 = glmm_shuff1(t1, 0, 1, 2, 3); + v0 = _mm_mul_ps(t0, v0); + v1 = _mm_mul_ps(t1, v1); + v2 = glmm_shuff1(v1, 1, 0, 0, 1); + v3 = glmm_shuff1(v0, 0, 1, 2, 3); + v1 = _mm_add_ps(v1, v2); + v0 = _mm_add_ps(v0, v3); + v2 = glmm_shuff1(v0, 1, 0, 0, 1); + v0 = _mm_add_ps(v0, v2); + + v0 = _mm_sub_ps(v0, v1); /* det */ + v0 = _mm_div_ps(v5, v0); /* idt */ + + /* multiply t0,t1,t2 by idt to reduce 1mul below: 2eor+4mul vs 3mul+4eor */ + t0 = _mm_mul_ps(t0, v0); + t1 = _mm_mul_ps(t1, v0); + t2 = _mm_mul_ps(t2, v0); + + v0 = glmm_shuff1(t0, 0, 0, 1, 1); /* c2 c2 c1 c1 */ + v1 = glmm_shuff1(t0, 2, 2, 3, 3); /* c4 c4 c3 c3 */ + v2 = glmm_shuff1(t1, 0, 0, 1, 1); /* c12 c12 c11 c11 */ + v3 = glmm_shuff1(t1, 2, 2, 3, 3); /* c6 c6 c5 c5 */ + v4 = glmm_shuff1(t2, 0, 0, 1, 1); /* c8 c8 c7 c7 */ + v5 = glmm_shuff1(t2, 2, 2, 3, 3); /* c10 c10 c9 c9 */ + + /* result */ + + /* dest[0][0] = (f * c1 - g * c5 + h * c9) * idt; + dest[0][1] = (b * c1 - c * c5 + d * c9) * ndt; + dest[0][2] = (n * c2 - o * c6 + p * c10) * idt; + dest[0][3] = (j * c2 - k * c6 + l * c10) * ndt; + + dest[1][0] = (e * c1 - g * c3 + h * c11) * ndt; + dest[1][1] = (a * c1 - c * c3 + d * c11) * idt; + dest[1][2] = (m * c2 - o * c4 + p * c12) * ndt; + dest[1][3] = (i * c2 - k * c4 + l * c12) * idt; + + dest[2][0] = (e * c5 - f * c3 + h * c7) * idt; + dest[2][1] = (a * c5 - b * c3 + d * c7) * ndt; + dest[2][2] = (m * c6 - n * c4 + p * c8) * idt; + dest[2][3] = (i * c6 - j * c4 + l * c8) * ndt; + + dest[3][0] = (e * c9 - f * c11 + g * c7) * ndt; + dest[3][1] = (a * c9 - b * c11 + c * c7) * idt; + dest[3][2] = (m * c10 - n * c12 + o * c8) * ndt; + dest[3][3] = (i * c10 - j * c12 + k * c8) * idt; */ + + r0 = _mm_mul_ps(x0, v0); + r1 = _mm_mul_ps(x1, v0); + r2 = _mm_mul_ps(x1, v3); + r3 = _mm_mul_ps(x1, v5); + + r0 = glmm_fnmadd(x3, v3, r0); + r1 = glmm_fnmadd(x3, v1, r1); + r2 = glmm_fnmadd(x0, v1, r2); + r3 = glmm_fnmadd(x0, v2, r3); + + r0 = glmm_fmadd(x2, v5, r0); + r1 = glmm_fmadd(x2, v2, r1); + r2 = glmm_fmadd(x2, v4, r2); + r3 = glmm_fmadd(x3, v4, r3); + + /* 4xor may be fastart then 4mul, see above */ + r0 = _mm_xor_ps(r0, s1); + r1 = _mm_xor_ps(r1, s2); + r2 = _mm_xor_ps(r2, s1); + r3 = _mm_xor_ps(r3, s2); + + glmm_store(dest[0], r0); + glmm_store(dest[1], r1); + glmm_store(dest[2], r2); + glmm_store(dest[3], r3); +} +#endif +#endif /* cglm_mat_sse_h */ diff --git a/include/cglm/simd/sse2/quat.h b/include/cglm/simd/sse2/quat.h new file mode 100644 index 0000000..def0fe2 --- /dev/null +++ b/include/cglm/simd/sse2/quat.h @@ -0,0 +1,54 @@ +/* + * Copyright (c), Recep Aslantas. + * + * MIT License (MIT), http://opensource.org/licenses/MIT + * Full license can be found in the LICENSE file + */ + +#ifndef cglm_quat_simd_h +#define cglm_quat_simd_h +#if defined( __SSE__ ) || defined( __SSE2__ ) + +#include "../../common.h" +#include "../intrin.h" + +CGLM_INLINE +void +glm_quat_mul_sse2(versor p, versor q, versor dest) { + /* + + (a1 b2 + b1 a2 + c1 d2 − d1 c2)i + + (a1 c2 − b1 d2 + c1 a2 + d1 b2)j + + (a1 d2 + b1 c2 − c1 b2 + d1 a2)k + a1 a2 − b1 b2 − c1 c2 − d1 d2 + */ + + __m128 xp, xq, x1, x2, x3, r, x, y, z; + + xp = glmm_load(p); /* 3 2 1 0 */ + xq = glmm_load(q); + x1 = glmm_float32x4_SIGNMASK_NPNP; /* TODO: _mm_set1_ss() + shuff ? */ + r = _mm_mul_ps(glmm_splat_w(xp), xq); + + x2 = _mm_unpackhi_ps(x1, x1); + x3 = glmm_shuff1(x1, 3, 2, 0, 1); + x = glmm_splat_x(xp); + y = glmm_splat_y(xp); + z = glmm_splat_z(xp); + + x = _mm_xor_ps(x, x1); + y = _mm_xor_ps(y, x2); + z = _mm_xor_ps(z, x3); + + x1 = glmm_shuff1(xq, 0, 1, 2, 3); + x2 = glmm_shuff1(xq, 1, 0, 3, 2); + x3 = glmm_shuff1(xq, 2, 3, 0, 1); + + r = glmm_fmadd(x, x1, r); + r = glmm_fmadd(y, x2, r); + r = glmm_fmadd(z, x3, r); + + glmm_store(dest, r); +} + +#endif +#endif /* cglm_quat_simd_h */ -- cgit