diff --git a/core/math/linalg/general.odin b/core/math/linalg/general.odin index 8c4f2954a..4a0150972 100644 --- a/core/math/linalg/general.odin +++ b/core/math/linalg/general.odin @@ -417,6 +417,13 @@ adjugate :: proc{ matrix4x4_adjugate, } +cofactor :: proc{ + matrix1x1_cofactor, + matrix2x2_cofactor, + matrix3x3_cofactor, + matrix4x4_cofactor, +} + inverse_transpose :: proc{ matrix1x1_inverse_transpose, matrix2x2_inverse_transpose, @@ -479,9 +486,9 @@ matrix3x3_determinant :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (det: T) # } @(require_results) matrix4x4_determinant :: proc "contextless" (m: $M/matrix[4, 4]$T) -> (det: T) #no_bounds_check { - a := adjugate(m) + c := cofactor(m) for i in 0..<4 { - det += m[0, i] * a[0, i] + det += m[0, i] * c[0, i] } return } @@ -497,6 +504,47 @@ matrix1x1_adjugate :: proc "contextless" (x: $M/matrix[1, 1]$T) -> (y: M) #no_bo @(require_results) matrix2x2_adjugate :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) #no_bounds_check { + y[0, 0] = +x[1, 1] + y[0, 1] = -x[0, 1] + y[1, 0] = -x[1, 0] + y[1, 1] = +x[0, 0] + return +} + +@(require_results) +matrix3x3_adjugate :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { + y[0, 0] = +(m[1, 1] * m[2, 2] - m[2, 1] * m[1, 2]) + y[1, 0] = -(m[1, 0] * m[2, 2] - m[2, 0] * m[1, 2]) + y[2, 0] = +(m[1, 0] * m[2, 1] - m[2, 0] * m[1, 1]) + y[0, 1] = -(m[0, 1] * m[2, 2] - m[2, 1] * m[0, 2]) + y[1, 1] = +(m[0, 0] * m[2, 2] - m[2, 0] * m[0, 2]) + y[2, 1] = -(m[0, 0] * m[2, 1] - m[2, 0] * m[0, 1]) + y[0, 2] = +(m[0, 1] * m[1, 2] - m[1, 1] * m[0, 2]) + y[1, 2] = -(m[0, 0] * m[1, 2] - m[1, 0] * m[0, 2]) + y[2, 2] = +(m[0, 0] * m[1, 1] - m[1, 0] * m[0, 1]) + return +} + +@(require_results) +matrix4x4_adjugate :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { + for i in 0..<4 { + for j in 0..<4 { + sign: T = 1 if (i + j) % 2 == 0 else -1 + y[i, j] = sign * matrix_minor(x, j, i) + } + } + return +} + + +@(require_results) +matrix1x1_cofactor :: proc "contextless" (x: $M/matrix[1, 1]$T) -> (y: M) #no_bounds_check { + y = x + return +} + +@(require_results) +matrix2x2_cofactor :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) #no_bounds_check { y[0, 0] = +x[1, 1] y[0, 1] = -x[1, 0] y[1, 0] = -x[0, 1] @@ -505,7 +553,7 @@ matrix2x2_adjugate :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) #no_bo } @(require_results) -matrix3x3_adjugate :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { +matrix3x3_cofactor :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { y[0, 0] = +(m[1, 1] * m[2, 2] - m[2, 1] * m[1, 2]) y[0, 1] = -(m[1, 0] * m[2, 2] - m[2, 0] * m[1, 2]) y[0, 2] = +(m[1, 0] * m[2, 1] - m[2, 0] * m[1, 1]) @@ -520,7 +568,7 @@ matrix3x3_adjugate :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) #no_bo @(require_results) -matrix4x4_adjugate :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { +matrix4x4_cofactor :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { for i in 0..<4 { for j in 0..<4 { sign: T = 1 if (i + j) % 2 == 0 else -1 @@ -556,19 +604,19 @@ matrix2x2_inverse_transpose :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: @(require_results) matrix3x3_inverse_transpose :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d := determinant(x) when intrinsics.type_is_integer(T) { for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[i, j] / d + y[i, j] = c[i, j] / d } } } else { id := 1/d for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[i, j] * id + y[i, j] = c[i, j] * id } } } @@ -577,22 +625,22 @@ matrix3x3_inverse_transpose :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: @(require_results) matrix4x4_inverse_transpose :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d: T for i in 0..<4 { - d += x[0, i] * a[0, i] + d += x[0, i] * c[0, i] } when intrinsics.type_is_integer(T) { for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[i, j] / d + y[i, j] = c[i, j] / d } } } else { id := 1/d for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[i, j] * id + y[i, j] = c[i, j] * id } } } @@ -625,19 +673,19 @@ matrix2x2_inverse :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) #no_bou @(require_results) matrix3x3_inverse :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d := determinant(x) when intrinsics.type_is_integer(T) { for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[j, i] / d + y[i, j] = c[j, i] / d } } } else { id := 1/d for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[j, i] * id + y[i, j] = c[j, i] * id } } } @@ -646,22 +694,22 @@ matrix3x3_inverse :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bou @(require_results) matrix4x4_inverse :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d: T for i in 0..<4 { - d += x[0, i] * a[0, i] + d += x[0, i] * c[0, i] } when intrinsics.type_is_integer(T) { for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[j, i] / d + y[i, j] = c[j, i] / d } } } else { id := 1/d for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[j, i] * id + y[i, j] = c[j, i] * id } } } diff --git a/core/math/linalg/glsl/linalg_glsl.odin b/core/math/linalg/glsl/linalg_glsl.odin index ca61891cb..bd2cf416a 100644 --- a/core/math/linalg/glsl/linalg_glsl.odin +++ b/core/math/linalg/glsl/linalg_glsl.odin @@ -1882,6 +1882,13 @@ adjugate :: proc{ adjugate_matrix4x4, } +cofactor :: proc{ + cofactor_matrix1x1, + cofactor_matrix2x2, + cofactor_matrix3x3, + cofactor_matrix4x4, +} + inverse_transpose :: proc{ inverse_transpose_matrix1x1, inverse_transpose_matrix2x2, @@ -1944,9 +1951,9 @@ determinant_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (det: T) { } @(require_results) determinant_matrix4x4 :: proc "contextless" (m: $M/matrix[4, 4]$T) -> (det: T) { - a := adjugate(m) + c := cofactor(m) #no_bounds_check for i in 0..<4 { - det += m[0, i] * a[0, i] + det += m[0, i] * c[0, i] } return } @@ -1962,6 +1969,47 @@ adjugate_matrix1x1 :: proc "contextless" (x: $M/matrix[1, 1]$T) -> (y: M) { @(require_results) adjugate_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) { + y[0, 0] = +x[1, 1] + y[0, 1] = -x[0, 1] + y[1, 0] = -x[1, 0] + y[1, 1] = +x[0, 0] + return +} + +@(require_results) +adjugate_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) { + y[0, 0] = +(m[1, 1] * m[2, 2] - m[2, 1] * m[1, 2]) + y[1, 0] = -(m[1, 0] * m[2, 2] - m[2, 0] * m[1, 2]) + y[2, 0] = +(m[1, 0] * m[2, 1] - m[2, 0] * m[1, 1]) + y[0, 1] = -(m[0, 1] * m[2, 2] - m[2, 1] * m[0, 2]) + y[1, 1] = +(m[0, 0] * m[2, 2] - m[2, 0] * m[0, 2]) + y[2, 1] = -(m[0, 0] * m[2, 1] - m[2, 0] * m[0, 1]) + y[0, 2] = +(m[0, 1] * m[1, 2] - m[1, 1] * m[0, 2]) + y[1, 2] = -(m[0, 0] * m[1, 2] - m[1, 0] * m[0, 2]) + y[2, 2] = +(m[0, 0] * m[1, 1] - m[1, 0] * m[0, 1]) + return +} + +@(require_results) +adjugate_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) { + for i in 0..<4 { + for j in 0..<4 { + sign: T = 1 if (i + j) % 2 == 0 else -1 + y[i, j] = sign * matrix_minor(x, j, i) + } + } + return +} + + +@(require_results) +cofactor_matrix1x1 :: proc "contextless" (x: $M/matrix[1, 1]$T) -> (y: M) { + y = x + return +} + +@(require_results) +cofactor_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) { y[0, 0] = +x[1, 1] y[0, 1] = -x[1, 0] y[1, 0] = -x[0, 1] @@ -1970,7 +2018,7 @@ adjugate_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) { } @(require_results) -adjugate_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) { +cofactor_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) { y[0, 0] = +(m[1, 1] * m[2, 2] - m[2, 1] * m[1, 2]) y[0, 1] = -(m[1, 0] * m[2, 2] - m[2, 0] * m[1, 2]) y[0, 2] = +(m[1, 0] * m[2, 1] - m[2, 0] * m[1, 1]) @@ -1985,7 +2033,7 @@ adjugate_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) { @(require_results) -adjugate_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) { +cofactor_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) { for i in 0..<4 { for j in 0..<4 { sign: T = 1 if (i + j) % 2 == 0 else -1 @@ -2021,19 +2069,19 @@ inverse_transpose_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: @(require_results) inverse_transpose_matrix3x3 :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d := determinant(x) when intrinsics.type_is_integer(T) { for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[i, j] / d + y[i, j] = c[i, j] / d } } } else { id := 1/d for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[i, j] * id + y[i, j] = c[i, j] * id } } } @@ -2042,22 +2090,22 @@ inverse_transpose_matrix3x3 :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: @(require_results) inverse_transpose_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d: T for i in 0..<4 { - d += x[0, i] * a[0, i] + d += x[0, i] * c[0, i] } when intrinsics.type_is_integer(T) { for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[i, j] / d + y[i, j] = c[i, j] / d } } } else { id := 1/d for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[i, j] * id + y[i, j] = c[i, j] * id } } } @@ -2090,19 +2138,19 @@ inverse_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) { @(require_results) inverse_matrix3x3 :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d := determinant(x) when intrinsics.type_is_integer(T) { for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[j, i] / d + y[i, j] = c[j, i] / d } } } else { id := 1/d for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[j, i] * id + y[i, j] = c[j, i] * id } } } @@ -2111,22 +2159,22 @@ inverse_matrix3x3 :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bou @(require_results) inverse_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d: T for i in 0..<4 { - d += x[0, i] * a[0, i] + d += x[0, i] * c[0, i] } when intrinsics.type_is_integer(T) { for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[j, i] / d + y[i, j] = c[j, i] / d } } } else { id := 1/d for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[j, i] * id + y[i, j] = c[j, i] * id } } } diff --git a/core/math/linalg/hlsl/linalg_hlsl.odin b/core/math/linalg/hlsl/linalg_hlsl.odin index a89fdddd3..cca70f9c8 100644 --- a/core/math/linalg/hlsl/linalg_hlsl.odin +++ b/core/math/linalg/hlsl/linalg_hlsl.odin @@ -1514,6 +1514,13 @@ adjugate :: proc{ adjugate_matrix4x4, } +cofactor :: proc{ + cofactor_matrix1x1, + cofactor_matrix2x2, + cofactor_matrix3x3, + cofactor_matrix4x4, +} + inverse_transpose :: proc{ inverse_transpose_matrix1x1, inverse_transpose_matrix2x2, @@ -1568,9 +1575,9 @@ determinant_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (det: T) { } @(require_results) determinant_matrix4x4 :: proc "contextless" (m: $M/matrix[4, 4]$T) -> (det: T) { - a := adjugate(m) + c := cofactor(m) #no_bounds_check for i in 0..<4 { - det += m[0, i] * a[0, i] + det += m[0, i] * c[0, i] } return } @@ -1586,6 +1593,47 @@ adjugate_matrix1x1 :: proc "contextless" (x: $M/matrix[1, 1]$T) -> (y: M) { @(require_results) adjugate_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) { + y[0, 0] = +x[1, 1] + y[0, 1] = -x[0, 1] + y[1, 0] = -x[1, 0] + y[1, 1] = +x[0, 0] + return +} + +@(require_results) +adjugate_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) { + y[0, 0] = +(m[1, 1] * m[2, 2] - m[2, 1] * m[1, 2]) + y[1, 0] = -(m[1, 0] * m[2, 2] - m[2, 0] * m[1, 2]) + y[2, 0] = +(m[1, 0] * m[2, 1] - m[2, 0] * m[1, 1]) + y[0, 1] = -(m[0, 1] * m[2, 2] - m[2, 1] * m[0, 2]) + y[1, 1] = +(m[0, 0] * m[2, 2] - m[2, 0] * m[0, 2]) + y[2, 1] = -(m[0, 0] * m[2, 1] - m[2, 0] * m[0, 1]) + y[0, 2] = +(m[0, 1] * m[1, 2] - m[1, 1] * m[0, 2]) + y[1, 2] = -(m[0, 0] * m[1, 2] - m[1, 0] * m[0, 2]) + y[2, 2] = +(m[0, 0] * m[1, 1] - m[1, 0] * m[0, 1]) + return +} + +@(require_results) +adjugate_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) { + for i in 0..<4 { + for j in 0..<4 { + sign: T = 1 if (i + j) % 2 == 0 else -1 + y[i, j] = sign * matrix_minor(x, j, i) + } + } + return +} + + +@(require_results) +cofactor_matrix1x1 :: proc "contextless" (x: $M/matrix[1, 1]$T) -> (y: M) { + y = x + return +} + +@(require_results) +cofactor_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) { y[0, 0] = +x[1, 1] y[0, 1] = -x[1, 0] y[1, 0] = -x[0, 1] @@ -1594,7 +1642,7 @@ adjugate_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) { } @(require_results) -adjugate_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) { +cofactor_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) { y[0, 0] = +(m[1, 1] * m[2, 2] - m[2, 1] * m[1, 2]) y[0, 1] = -(m[1, 0] * m[2, 2] - m[2, 0] * m[1, 2]) y[0, 2] = +(m[1, 0] * m[2, 1] - m[2, 0] * m[1, 1]) @@ -1609,7 +1657,7 @@ adjugate_matrix3x3 :: proc "contextless" (m: $M/matrix[3, 3]$T) -> (y: M) { @(require_results) -adjugate_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) { +cofactor_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) { for i in 0..<4 { for j in 0..<4 { sign: T = 1 if (i + j) % 2 == 0 else -1 @@ -1645,19 +1693,19 @@ inverse_transpose_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: @(require_results) inverse_transpose_matrix3x3 :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d := determinant(x) when intrinsics.type_is_integer(T) { for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[i, j] / d + y[i, j] = c[i, j] / d } } } else { id := 1/d for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[i, j] * id + y[i, j] = c[i, j] * id } } } @@ -1666,22 +1714,22 @@ inverse_transpose_matrix3x3 :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: @(require_results) inverse_transpose_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d: T for i in 0..<4 { - d += x[0, i] * a[0, i] + d += x[0, i] * c[0, i] } when intrinsics.type_is_integer(T) { for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[i, j] / d + y[i, j] = c[i, j] / d } } } else { id := 1/d for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[i, j] * id + y[i, j] = c[i, j] * id } } } @@ -1714,19 +1762,19 @@ inverse_matrix2x2 :: proc "contextless" (x: $M/matrix[2, 2]$T) -> (y: M) { @(require_results) inverse_matrix3x3 :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d := determinant(x) when intrinsics.type_is_integer(T) { for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[j, i] / d + y[i, j] = c[j, i] / d } } } else { id := 1/d for i in 0..<3 { for j in 0..<3 { - y[i, j] = a[j, i] * id + y[i, j] = c[j, i] * id } } } @@ -1735,22 +1783,22 @@ inverse_matrix3x3 :: proc "contextless" (x: $M/matrix[3, 3]$T) -> (y: M) #no_bou @(require_results) inverse_matrix4x4 :: proc "contextless" (x: $M/matrix[4, 4]$T) -> (y: M) #no_bounds_check { - a := adjugate(x) + c := cofactor(x) d: T for i in 0..<4 { - d += x[0, i] * a[0, i] + d += x[0, i] * c[0, i] } when intrinsics.type_is_integer(T) { for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[j, i] / d + y[i, j] = c[j, i] / d } } } else { id := 1/d for i in 0..<4 { for j in 0..<4 { - y[i, j] = a[j, i] * id + y[i, j] = c[j, i] * id } } } diff --git a/tests/issues/run.bat b/tests/issues/run.bat index dcea3d483..7ed43205d 100644 --- a/tests/issues/run.bat +++ b/tests/issues/run.bat @@ -16,6 +16,7 @@ set COMMON=-define:ODIN_TEST_FANCY=false -file -vet -strict-style ..\..\..\odin test ..\test_issue_2637.odin %COMMON% || exit /b ..\..\..\odin test ..\test_issue_2666.odin %COMMON% || exit /b ..\..\..\odin test ..\test_issue_4210.odin %COMMON% || exit /b +..\..\..\odin test ..\test_issue_4584.odin %COMMON% || exit /b @echo off diff --git a/tests/issues/run.sh b/tests/issues/run.sh index c3bc00e24..54543980e 100755 --- a/tests/issues/run.sh +++ b/tests/issues/run.sh @@ -17,6 +17,7 @@ $ODIN test ../test_issue_2615.odin $COMMON $ODIN test ../test_issue_2637.odin $COMMON $ODIN test ../test_issue_2666.odin $COMMON $ODIN test ../test_issue_4210.odin $COMMON +$ODIN test ../test_issue_4584.odin $COMMON if [[ $($ODIN build ../test_issue_2395.odin $COMMON 2>&1 >/dev/null | grep -c "Error:") -eq 2 ]] ; then echo "SUCCESSFUL 1/1" else diff --git a/tests/issues/test_issue_4584.odin b/tests/issues/test_issue_4584.odin new file mode 100644 index 000000000..9eea23420 --- /dev/null +++ b/tests/issues/test_issue_4584.odin @@ -0,0 +1,198 @@ +// Tests issue #4584 https://github.com/odin-lang/Odin/issues/4584 +package test_issues + +import "core:testing" +import "core:log" +import "core:math/linalg" +import glm "core:math/linalg/glsl" +import hlm "core:math/linalg/hlsl" + +@test +test_adjugate_2x2 :: proc(t: ^testing.T) { + I := linalg.identity(matrix[2,2]int) + m := matrix[2,2]int { + -3, 2, + -1, 0, + } + expected := matrix[2,2]int { + 0, -2, + 1, -3, + } + testing.expect_value(t, linalg.adjugate(m), expected) + testing.expect_value(t, linalg.determinant(m), 2) + testing.expect_value(t, linalg.adjugate(m) * m, 2 * I) + testing.expect_value(t, m * linalg.adjugate(m), 2 * I) + + testing.expect_value(t, glm.adjugate(m), expected) + testing.expect_value(t, glm.determinant(m), 2) + testing.expect_value(t, glm.adjugate(m) * m, 2 * I) + testing.expect_value(t, m * glm.adjugate(m), 2 * I) + + testing.expect_value(t, hlm.adjugate(m), expected) + testing.expect_value(t, hlm.determinant(m), 2) + testing.expect_value(t, hlm.adjugate(m) * m, 2 * I) + testing.expect_value(t, m * hlm.adjugate(m), 2 * I) +} + +@test +test_adjugate_3x3 :: proc(t: ^testing.T) { + I := linalg.identity(matrix[3,3]int) + m := matrix[3,3]int { + -3, 2, -5, + -1, 0, -2, + 3, -4, 1, + } + expected := matrix[3,3]int { + -8, 18, -4, + -5, 12, -1, + 4, -6, 2, + } + testing.expect_value(t, linalg.adjugate(m), expected) + testing.expect_value(t, linalg.determinant(m), -6) + testing.expect_value(t, linalg.adjugate(m) * m, -6 * I) + testing.expect_value(t, m * linalg.adjugate(m), -6 * I) + + testing.expect_value(t, glm.adjugate(m), expected) + testing.expect_value(t, glm.determinant(m), -6) + testing.expect_value(t, glm.adjugate(m) * m, -6 * I) + testing.expect_value(t, m * glm.adjugate(m), -6 * I) + + testing.expect_value(t, hlm.adjugate(m), expected) + testing.expect_value(t, hlm.determinant(m), -6) + testing.expect_value(t, hlm.adjugate(m) * m, -6 * I) + testing.expect_value(t, m * hlm.adjugate(m), -6 * I) +} + +@test +test_adjugate_4x4 :: proc(t: ^testing.T) { + I := linalg.identity(matrix[4,4]int) + m := matrix[4,4]int { + -3, 2, -5, 1, + -1, 0, -2, 2, + 3, -4, 1, 3, + 4, 5, 6, 7, + } + expected := matrix[4,4]int { + -144, 266, -92, -16, + -57, 92, -5, -16, + 105, -142, 55, 2, + 33, -96, 9, -6, + } + testing.expect_value(t, linalg.adjugate(m), expected) + testing.expect_value(t, linalg.determinant(m), -174) + testing.expect_value(t, linalg.adjugate(m) * m, -174 * I) + testing.expect_value(t, m * linalg.adjugate(m), -174 * I) + + testing.expect_value(t, glm.adjugate(m), expected) + testing.expect_value(t, glm.determinant(m), -174) + testing.expect_value(t, glm.adjugate(m) * m, -174 * I) + testing.expect_value(t, m * glm.adjugate(m), -174 * I) + + testing.expect_value(t, hlm.adjugate(m), expected) + testing.expect_value(t, hlm.determinant(m), -174) + testing.expect_value(t, hlm.adjugate(m) * m, -174 * I) + testing.expect_value(t, m * hlm.adjugate(m), -174 * I) +} + +@test +test_inverse_regression_2x2 :: proc(t: ^testing.T) { + I := linalg.identity(matrix[2,2]f32) + m := matrix[2,2]f32 { + -3, 2, + -1, 0, + } + expected := matrix[2,2]f32 { + 0.0, -1.0, + 1.0/2.0, -3.0/2.0, + } + expect_float_matrix_value(t, linalg.inverse(m), expected) + expect_float_matrix_value(t, linalg.inverse_transpose(m), linalg.transpose(expected)) + expect_float_matrix_value(t, linalg.inverse(m) * m, I) + expect_float_matrix_value(t, m * linalg.inverse(m), I) + + expect_float_matrix_value(t, glm.inverse(m), expected) + expect_float_matrix_value(t, glm.inverse_transpose(m), glm.transpose(expected)) + expect_float_matrix_value(t, glm.inverse(m) * m, I) + expect_float_matrix_value(t, m * glm.inverse(m), I) + + expect_float_matrix_value(t, hlm.inverse(m), expected) + expect_float_matrix_value(t, hlm.inverse_transpose(m), hlm.transpose(expected)) + expect_float_matrix_value(t, hlm.inverse(m) * m, I) + expect_float_matrix_value(t, m * hlm.inverse(m), I) +} + +@test +test_inverse_regression_3x3 :: proc(t: ^testing.T) { + I := linalg.identity(matrix[3,3]f32) + m := matrix[3,3]f32 { + -3, 2, -5, + -1, 0, -2, + 3, -4, 1, + } + expected := matrix[3,3]f32 { + 4.0/3.0, -3.0, 2.0/3.0, + 5.0/6.0, -2.0, 1.0/6.0, + -2.0/3.0, 1.0, -1.0/3.0, + } + expect_float_matrix_value(t, linalg.inverse(m), expected) + expect_float_matrix_value(t, linalg.inverse_transpose(m), linalg.transpose(expected)) + expect_float_matrix_value(t, linalg.inverse(m) * m, I) + expect_float_matrix_value(t, m * linalg.inverse(m), I) + + expect_float_matrix_value(t, glm.inverse(m), expected) + expect_float_matrix_value(t, glm.inverse_transpose(m), glm.transpose(expected)) + expect_float_matrix_value(t, glm.inverse(m) * m, I) + expect_float_matrix_value(t, m * glm.inverse(m), I) + + expect_float_matrix_value(t, hlm.inverse(m), expected) + expect_float_matrix_value(t, hlm.inverse_transpose(m), hlm.transpose(expected)) + expect_float_matrix_value(t, hlm.inverse(m) * m, I) + expect_float_matrix_value(t, m * hlm.inverse(m), I) +} + +@test +test_inverse_regression_4x4 :: proc(t: ^testing.T) { + I := linalg.identity(matrix[4,4]f32) + m := matrix[4,4]f32 { + -3, 2, -5, 1, + -1, 0, -2, 2, + 3, -4, 1, 3, + 4, 5, 6, 7, + } + expected := matrix[4,4]f32 { + 24.0/29.0, -133.0/87.0, 46.0/87.0, 8.0/87.0, + 19.0/58.0, -46.0/87.0, 5.0/174.0, 8.0/87.0, + -35.0/58.0, 71.0/87.0, -55.0/174.0, -1.0/87.0, + -11.0/58.0, 16.0/29.0, -3.0/58.0, 1.0/29.0, + } + expect_float_matrix_value(t, linalg.inverse(m), expected) + expect_float_matrix_value(t, linalg.inverse_transpose(m), linalg.transpose(expected)) + expect_float_matrix_value(t, linalg.inverse(m) * m, I) + expect_float_matrix_value(t, m * linalg.inverse(m), I) + + expect_float_matrix_value(t, glm.inverse(m), expected) + expect_float_matrix_value(t, glm.inverse_transpose(m), glm.transpose(expected)) + expect_float_matrix_value(t, glm.inverse(m) * m, I) + expect_float_matrix_value(t, m * glm.inverse(m), I) + + expect_float_matrix_value(t, hlm.inverse(m), expected) + expect_float_matrix_value(t, hlm.inverse_transpose(m), hlm.transpose(expected)) + expect_float_matrix_value(t, hlm.inverse(m) * m, I) + expect_float_matrix_value(t, m * hlm.inverse(m), I) +} + +@(private="file") +expect_float_matrix_value :: proc(t: ^testing.T, value, expected: $M/matrix[$N, N]f32, loc := #caller_location, value_expr := #caller_expression(value)) -> bool { + ok := true + outer: for i in 0.. 1e-6 { + ok = false + break outer + } + } + } + if !ok do log.errorf("expected %v to be %v, got %v", value_expr, expected, value, location=loc) + return ok +} \ No newline at end of file