diff --git a/linalg/main.c b/linalg/main.c index a06c9be..7c713cf 100644 --- a/linalg/main.c +++ b/linalg/main.c @@ -7,6 +7,8 @@ * Simple linear algebra mat/vec operations * * TODO: Normalizations + * TODO: Ensure suitable test coverade + * TODO: Mat3 adjoint should use a cofactor and a transpose function */ /** @@ -200,11 +202,36 @@ bool vec3_approx_eq(const Vec3 *a, const Vec3 *b, float epsilon); */ bool vec2_approx_eq(const Vec2 *a, const Vec2 *b, float epsilon); +/** + * @brief Computes the adjugate (adjoint) of a 2x2 matrix. + * + * The adjugate of a 2x2 matrix is obtained by swapping the diagonal elements + * and negating the off-diagonal elements. + * + * @param m Pointer to the 2x2 matrix. + * @return The adjugate of the input matrix. + */ +Mat2 mat2_adj(const Mat2 *m); + +/** + * @brief Computes the adjugate (adjoint) of a 3x3 matrix. + * + * The adjugate of a 3x3 matrix is the transpose of its cofactor matrix. + * It is used in computing the inverse of the matrix. + * + * @param m Pointer to the 3x3 matrix. + * @return The adjugate of the input matrix. + */ +Mat3 mat3_adj(const Mat3 *m); + #define MAT2_AT(m, row, col) ((m)->arr[(col) * 2 + (row)]) #define MAT3_AT(m, row, col) ((m)->arr[(col) * 3 + (row)]) /* Header end... */ +/* Row major */ +#define MAT2_DET(a, b, c, d) ((a) * (d) - (b) * (c)) + float mat3_det(const Mat3 *m) { float m00 = MAT3_AT(m, 0, 0); float m01 = MAT3_AT(m, 0, 1); @@ -361,6 +388,45 @@ bool vec3_approx_eq(const Vec3 *a, const Vec3 *b, float epsilon) { (fabsf(a->z - b->z) <= epsilon); } +Mat2 mat2_adj(const Mat2 *m) { + Mat2 res; + + MAT2_AT(&res, 0, 0) = MAT2_AT(m, 1, 1); + MAT2_AT(&res, 0, 1) = -MAT2_AT(m, 0, 1); + MAT2_AT(&res, 1, 0) = -MAT2_AT(m, 1, 0); + MAT2_AT(&res, 1, 1) = MAT2_AT(m, 0, 0); + + return res; +} + +Mat3 mat3_adj(const Mat3 *m) { + Mat3 res; + + const float a = MAT3_AT(m, 0, 0); + const float b = MAT3_AT(m, 0, 1); + const float c = MAT3_AT(m, 0, 2); + const float d = MAT3_AT(m, 1, 0); + const float e = MAT3_AT(m, 1, 1); + const float f = MAT3_AT(m, 1, 2); + const float g = MAT3_AT(m, 2, 0); + const float h = MAT3_AT(m, 2, 1); + const float i = MAT3_AT(m, 2, 2); + + MAT3_AT(&res, 0, 0) = MAT2_DET(e, f, h, i); + MAT3_AT(&res, 1, 0) = -MAT2_DET(d, f, g, i); + MAT3_AT(&res, 2, 0) = MAT2_DET(d, e, g, h); + + MAT3_AT(&res, 0, 1) = -MAT2_DET(b, c, h, i); + MAT3_AT(&res, 1, 1) = MAT2_DET(a, c, g, i); + MAT3_AT(&res, 2, 1) = -MAT2_DET(a, b, g, h); + + MAT3_AT(&res, 0, 2) = MAT2_DET(b, c, e, f); + MAT3_AT(&res, 1, 2) = -MAT2_DET(a, c, d, f); + MAT3_AT(&res, 2, 2) = MAT2_DET(a, b, d, e); + + return res; +} + /* Implem end */ int main(void) {