diff --git a/src/vec3.zig b/src/vec3.zig index 9f5dc1c..05db9d8 100644 --- a/src/vec3.zig +++ b/src/vec3.zig @@ -32,6 +32,10 @@ pub fn vec3Dot(a: Vec3F64, b: Vec3F64) f64 { return @reduce(.Add, a.data * b.data); } +pub fn vec3MulScalar(a: Vec3F64, b: f64) Vec3F64 { + return .init(a.x() * b, a.y() * b, a.z() * b); +} + pub fn vec3Cross(a: Vec3F64, b: Vec3F64) Vec3F64 { const x_component = a.y() * b.z() - a.z() * b.y(); const y_component = a.z() * b.x() - a.x() * b.z(); @@ -83,6 +87,45 @@ pub fn mat3Mul(a: Mat3F64, b: Mat3F64) Mat3F64 { return .{ .row1 = c_row1, .row2 = c_row2, .row3 = c_row3 }; } +pub fn mat3Det(a: Mat3F64) f64 { + return a.row1.x() * (a.row2.y() * a.row3.z() - a.row2.z() * a.row3.y()) - + a.row1.y() * (a.row2.x() * a.row3.z() - a.row2.z() * a.row3.x()) + + a.row1.y() * (a.row2.x() * a.row3.y() - a.row2.y() * a.row3.x()); +} + +pub fn mat3Inv(a: Mat3F64) Mat3F64 { + const inv_det = 1.0 / mat3Det(a); + const row1 = vec3MulScalar( + .init( + a.row2.y() * a.row3.z() - a.row2.z() * a.row3.y(), + a.row1.z() * a.row3.y() - a.row1.y() * a.row3.z(), + a.row1.y() * a.row2.z() - a.row1.z() * a.row2.y(), + ), + inv_det, + ); + const row2 = vec3MulScalar( + .init( + a.row2.z() * a.row3.x() + a.row2.x() * a.row3.z(), + a.row1.x() * a.row3.z() - a.row1.z() * a.row3.x(), + a.row1.z() * a.row2.x() - a.row1.x() * a.row2.z(), + ), + inv_det, + ); + const row3 = vec3MulScalar( + .init( + a.row2.x() * a.row3.y() + a.row2.y() * a.row3.x(), + a.row1.y() * a.row3.x() - a.row1.x() * a.row3.y(), + a.row1.x() * a.row2.y() - a.row1.y() * a.row2.x(), + ), + inv_det, + ); + return .{ + .row1 = row1, + .row2 = row2, + .row3 = row3, + }; +} + pub fn vec3MulMat3(a: Mat3F64, b: Vec3F64) Vec3F64 { return .init( vec3Dot(a.row1, b), @@ -267,6 +310,36 @@ test "Mat3 mul works" { try std.testing.expect(c.row3.z() == 9); } +test "mat determinant works" { + const a: Mat3F64 = .{ + .row1 = Vec3F64.init(1, 2, 3), + .row2 = Vec3F64.init(0, 1, 5), + .row3 = Vec3F64.init(0, 0, 9), + }; + try std.testing.expectApproxEqAbs(9, mat3Det(a), 1e-12); +} + +test "mat inverse works" { + const a: Mat3F64 = .{ + .row1 = Vec3F64.init(1, 2, 3), + .row2 = Vec3F64.init(0, 1, 5), + .row3 = Vec3F64.init(0, 0, 9), + }; + + const ainv = mat3Inv(a); + std.debug.print("{}", .{ainv}); + std.debug.print("{}", .{ainv.row2.z()}); + + try std.testing.expectApproxEqAbs(1, ainv.row1.x(), 1e-12); + try std.testing.expectApproxEqAbs(-2, ainv.row1.y(), 1e-12); + try std.testing.expectApproxEqAbs(7.0 / 9.0, ainv.row1.z(), 1e-12); + try std.testing.expectApproxEqAbs(0, ainv.row2.x(), 1e-12); + try std.testing.expectApproxEqAbs(1, ainv.row2.y(), 1e-12); + try std.testing.expectApproxEqAbs(-5.0 / 9.0, ainv.row2.z(), 1e-12); + try std.testing.expectApproxEqAbs(0, ainv.row3.x(), 1e-12); + try std.testing.expectApproxEqAbs(0, ainv.row3.y(), 1e-12); + try std.testing.expectApproxEqAbs(1.0 / 9.0, ainv.row3.z(), 1e-12); +} test "vec3 mul mat3 works" { const b: Mat3F64 = .{ .row1 = Vec3F64.init(0, 1, 0),