improved prediction

.github/workflows/cmake-multi-platform.yml

@@ -707,7 +707,11 @@ jobs:
 
       - name: Build
         run: |
+          if [[ "${{ matrix.msystem }}" == "MINGW32" ]]; then
+            cmake --build cmake-build/build/${{ matrix.preset }} --target unit_tests omath --parallel 1
+          else
             cmake --build cmake-build/build/${{ matrix.preset }} --target unit_tests omath
+          fi
 
       - name: Run unit_tests.exe
         run: |

benchmark/benchmark_mat.cpp

@@ -2,11 +2,9 @@
 // Created by Vlad on 9/17/2025.
 //
 #include <benchmark/benchmark.h>
-
 #include <omath/omath.hpp>
 using namespace omath;
 
-
 void mat_float_multiplication_col_major(benchmark::State& state)
 {
     using MatType = Mat<128, 128, float, MatStoreType::COLUMN_MAJOR>;
@@ -15,9 +13,12 @@ void mat_float_multiplication_col_major(benchmark::State& state)
     a.set(3.f);
     b.set(7.f);
 
-
     for ([[maybe_unused]] const auto _ : state)
-        std::ignore = a * b;
+    {
+        benchmark::DoNotOptimize(a);
+        benchmark::DoNotOptimize(b);
+        benchmark::DoNotOptimize(a * b);
+    }
 }
 void mat_float_multiplication_row_major(benchmark::State& state)
 {
@@ -27,9 +28,12 @@ void mat_float_multiplication_row_major(benchmark::State& state)
     a.set(3.f);
     b.set(7.f);
 
-
     for ([[maybe_unused]] const auto _ : state)
-        std::ignore = a * b;
+    {
+        benchmark::DoNotOptimize(a);
+        benchmark::DoNotOptimize(b);
+        benchmark::DoNotOptimize(a * b);
+    }
 }
 
 void mat_double_multiplication_row_major(benchmark::State& state)
@@ -40,9 +44,12 @@ void mat_double_multiplication_row_major(benchmark::State& state)
     a.set(3.f);
     b.set(7.f);
 
-
     for ([[maybe_unused]] const auto _ : state)
-        std::ignore = a * b;
+    {
+        benchmark::DoNotOptimize(a);
+        benchmark::DoNotOptimize(b);
+        benchmark::DoNotOptimize(a * b);
+    }
 }
 
 void mat_double_multiplication_col_major(benchmark::State& state)
@@ -53,9 +60,12 @@ void mat_double_multiplication_col_major(benchmark::State& state)
     a.set(3.f);
     b.set(7.f);
 
-
     for ([[maybe_unused]] const auto _ : state)
-        std::ignore = a * b;
+    {
+        benchmark::DoNotOptimize(a);
+        benchmark::DoNotOptimize(b);
+        benchmark::DoNotOptimize(a * b);
+    }
 }
 
 BENCHMARK(mat_float_multiplication_col_major)->Iterations(5000);

benchmark/benchmark_projectile_pred.cpp

@@ -2,22 +2,105 @@
 // Created by Vlad on 9/18/2025.
 //
 #include <benchmark/benchmark.h>
-#include <omath/omath.hpp>
+#include <omath/projectile_prediction/proj_pred_engine_legacy.hpp>
-using namespace omath;
 
-using namespace omath::projectile_prediction;
+namespace
-
-constexpr float simulation_time_step = 1.f / 1000.f;
-constexpr float hit_distance_tolerance = 5.f;
-
-void source_engine_projectile_prediction(benchmark::State& state)
 {
-    constexpr Target<float> target{.m_origin = {100, 0, 90}, .m_velocity = {0, 0, 0}, .m_is_airborne = false};
+    using Engine = omath::projectile_prediction::ProjPredEngineLegacy<>;
-    constexpr Projectile<float> projectile = {.m_origin = {3, 2, 1}, .m_launch_speed = 5000.f, .m_gravity_scale = 0.4f};
+    using Projectile = omath::projectile_prediction::Projectile<float>;
+    using Target = omath::projectile_prediction::Target<float>;
 
-    for ([[maybe_unused]] const auto _: state)
+    struct PredictionScenario
-        std::ignore = ProjPredEngineLegacy<>(400.f, simulation_time_step, 50.f, hit_distance_tolerance)
+    {
-                              .maybe_calculate_aim_point(projectile, target);
+        Projectile projectile;
-}
+        Target target;
+        float gravity;
+        float simulation_time_step;
+        float maximum_simulation_time;
+        float distance_tolerance;
+        bool expects_solution;
+    };
 
-BENCHMARK(source_engine_projectile_prediction)->Iterations(10'000);
+    void run_prediction_benchmark(benchmark::State& state, const PredictionScenario& scenario)
+    {
+        const Engine engine(scenario.gravity, scenario.simulation_time_step, scenario.maximum_simulation_time,
+                            scenario.distance_tolerance);
+        auto projectile = scenario.projectile;
+        auto target = scenario.target;
+
+        if (engine.maybe_calculate_aim_point(projectile, target).has_value() != scenario.expects_solution)
+        {
+            state.SkipWithError("Projectile benchmark scenario returned an unexpected result");
+            return;
+        }
+
+        for ([[maybe_unused]] const auto _ : state)
+        {
+            benchmark::DoNotOptimize(projectile);
+            benchmark::DoNotOptimize(target);
+            auto result = engine.maybe_calculate_aim_point(projectile, target);
+            benchmark::DoNotOptimize(result);
+        }
+    }
+
+    void projectile_prediction_near_static_hit(benchmark::State& state)
+    {
+        constexpr PredictionScenario scenario{
+                .projectile = {.m_origin = {3.f, 2.f, 1.f}, .m_launch_speed = 5000.f, .m_gravity_scale = 0.4f},
+                .target = {.m_origin = {100.f, 0.f, 90.f}, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false},
+                .gravity = 400.f,
+                .simulation_time_step = 1.f / 1000.f,
+                .maximum_simulation_time = 50.f,
+                .distance_tolerance = 5.f,
+                .expects_solution = true,
+        };
+        run_prediction_benchmark(state, scenario);
+    }
+
+    void projectile_prediction_moving_hit(benchmark::State& state)
+    {
+        constexpr PredictionScenario scenario{
+                .projectile = {.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 3000.f, .m_gravity_scale = 1.f},
+                .target = {.m_origin = {500.f, 100.f, 0.f}, .m_velocity = {-50.f, 20.f, 0.f}, .m_is_airborne = false},
+                .gravity = 800.f,
+                .simulation_time_step = 1.f / 500.f,
+                .maximum_simulation_time = 30.f,
+                .distance_tolerance = 10.f,
+                .expects_solution = true,
+        };
+        run_prediction_benchmark(state, scenario);
+    }
+
+    void projectile_prediction_unreachable_full_scan(benchmark::State& state)
+    {
+        constexpr PredictionScenario scenario{
+                .projectile = {.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 1.f, .m_gravity_scale = 1.f},
+                .target = {.m_origin = {100'000.f, 0.f, 0.f}, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false},
+                .gravity = 9.81f,
+                .simulation_time_step = 1.f / 1000.f,
+                .maximum_simulation_time = 2.f,
+                .distance_tolerance = 5.f,
+                .expects_solution = false,
+        };
+        run_prediction_benchmark(state, scenario);
+    }
+
+    void projectile_prediction_receding_full_scan(benchmark::State& state)
+    {
+        constexpr PredictionScenario scenario{
+                .projectile = {.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 100.f, .m_gravity_scale = 0.f},
+                .target = {.m_origin = {100.f, 0.f, 0.f}, .m_velocity = {200.f, 0.f, 0.f}, .m_is_airborne = false},
+                .gravity = 9.81f,
+                .simulation_time_step = 1.f / 1000.f,
+                .maximum_simulation_time = 2.f,
+                .distance_tolerance = 0.01f,
+                .expects_solution = false,
+        };
+        run_prediction_benchmark(state, scenario);
+    }
+} // namespace
+
+BENCHMARK(projectile_prediction_near_static_hit);
+BENCHMARK(projectile_prediction_moving_hit);
+BENCHMARK(projectile_prediction_unreachable_full_scan);
+BENCHMARK(projectile_prediction_receding_full_scan);

docs/trigonometry/angle.md

@@ -3,8 +3,8 @@
 > Header: `omath/trigonometry/angle.hpp`
 > Namespace: `omath`
 > Template: `Angle<Type = float, min = 0, max = 360, flags = AngleFlags::Normalized>`
-> Requires: `std::is_arithmetic_v<Type>`
+> Requires: `std::is_floating_point_v<Type>`
-> Formatters: `std::formatter` for `char`, `wchar_t`, `char8_t` → `"{}deg"`
+> Formatters: `std::formatter` for `char` and `wchar_t` → `"{}deg"`
 
 ---
 
@@ -14,7 +14,7 @@
 
 Two behaviors via `AngleFlags`:
 
-* `AngleFlags::Normalized` (default): values are wrapped into `[min, max]` using `angles::wrap_angle`.
+* `AngleFlags::Normalized` (default): values are wrapped into `[min, max)` using `angles::wrap_angle`.
 * `AngleFlags::Clamped`: values are clamped to `[min, max]` using `std::clamp`.
 
 ---
@@ -28,12 +28,16 @@ enum class AngleFlags { Normalized = 0, Clamped = 1 };
 
 template<class Type = float, Type min = Type(0), Type max = Type(360),
          AngleFlags flags = AngleFlags::Normalized>
-requires std::is_arithmetic_v<Type>
+requires std::is_floating_point_v<Type>
 class Angle {
 public:
   // Construction
   static constexpr Angle from_degrees(const Type& deg) noexcept;
   static constexpr Angle from_radians(const Type& rad) noexcept;
+  static constexpr Angle from_asin(const Type& value) noexcept;
+  static constexpr Angle from_acos(const Type& value) noexcept;
+  static constexpr Angle from_atan(const Type& value) noexcept;
+  static constexpr Angle from_atan2(const Type& y, const Type& x) noexcept;
   constexpr Angle() noexcept;                  // 0 deg, adjusted by flags/range
 
   // Accessors / conversions (degrees stored internally)
@@ -45,10 +49,9 @@ public:
   Type sin() const noexcept;
   Type cos() const noexcept;
   Type tan() const noexcept;
-  Type atan() const noexcept;  // atan(as_radians())  (rarely used)
   Type cot() const noexcept;   // cos()/sin()  (watch sin≈0)
 
-  // Arithmetic (wraps or clamps per flags and [min,max])
+  // Arithmetic (wraps or clamps per flags and configured range)
   constexpr Angle& operator+=(const Angle&) noexcept;
   constexpr Angle& operator-=(const Angle&) noexcept;
   constexpr Angle  operator+(const Angle&) noexcept;
@@ -68,7 +71,7 @@ public:
 std::format("{}", Angle<float>::from_degrees(45)); // "45deg"
 ```
 
-Formatters exist for `char`, `wchar_t`, and `char8_t`.
+Formatters exist for `char` and `wchar_t`.
 
 ---
 
@@ -116,10 +119,9 @@ float deg = *yaw;                                // same as yaw.as_degrees()
 ## Semantics & notes
 
 * **Storage & units:** Internally stores **degrees** (`Type m_angle`). `as_radians()`/`from_radians()` use the project helpers in `omath::angles`.
-* **Arithmetic honors policy:** `operator+=`/`-=` and the binary `+`/`-` apply **wrap** or **clamp** in `[min,max]`, mirroring construction behavior.
+* **Arithmetic honors policy:** `operator+=`/`-=` and the binary `+`/`-` apply **wrap** or **clamp**, mirroring construction behavior.
-* **`atan()`**: returns `std::atan(as_radians())` (the arctangent of the *radian value*). This is mathematically unusual for an angle type and is rarely useful; prefer `tan()`/`atan2` in client code when solving geometry problems.
 * **`cot()` / `tan()` singularities:** Near multiples where `sin() ≈ 0` or `cos() ≈ 0`, results blow up. Guard in your usage if inputs can approach these points.
-* **Comparison:** `operator<=>` is defaulted. With normalization, distinct representatives can compare as expected (e.g., `-180` vs `180` in signed ranges are distinct endpoints).
+* **Comparison:** `operator<=>` is defaulted. Normalization canonicalizes the maximum endpoint to the minimum endpoint.
 * **No implicit numeric conversion:** There’s **no `operator Type()`**. Use `as_degrees()`/`as_radians()` (or `*angle`) explicitly—this intentional friction avoids unit mistakes.
 
 ---

docs/trigonometry/angles.md

@@ -4,7 +4,7 @@
 > Namespace: `omath::angles`
 > All functions are `[[nodiscard]]` and `noexcept` where applicable.
 
-A small set of constexpr-friendly utilities for converting between degrees/radians, converting horizontal/vertical field of view, and wrapping angles into a closed interval.
+A small set of constexpr-friendly utilities for converting between degrees/radians, converting horizontal/vertical field of view, and wrapping angles into a half-open interval.
 
 ---
 
@@ -29,9 +29,9 @@ template<class Type>
 requires std::is_floating_point_v<Type>
 Type vertical_fov_to_horizontal(const Type& vertical_fov, const Type& aspect) noexcept;
 
-// Wrap angle into [min, max] (any arithmetic type)
+// Wrap angle into [min, max) (floating-point types)
 template<class Type>
-requires std::is_arithmetic_v<Type>
+requires std::is_floating_point_v<Type>
 Type wrap_angle(const Type& angle, const Type& min, const Type& max) noexcept;
 ```
 
@@ -66,10 +66,10 @@ Formulas (in radians):
 
 ### Wrapping angles (or any periodic value)
 
-Wrap any numeric `angle` into `[min, max]`:
+Wrap any floating-point `angle` into `[min, max)`:
 
 ```cpp
-// Wrap degrees into [0, 360]
+// Wrap degrees into [0, 360)
 float a = omath::angles::wrap_angle(   370.0f, 0.0f, 360.0f); // 10
 float b = omath::angles::wrap_angle(   -15.0f, 0.0f, 360.0f); // 345
 // Signed range [-180,180]
@@ -83,10 +83,10 @@ float c = omath::angles::wrap_angle(  200.0f, -180.0f, 180.0f); // -160
 * **Type requirements**
 
     * Converters & FOV helpers require **floating-point** `Type`.
-    * `wrap_angle` accepts any arithmetic `Type` (floats or integers).
+    * `wrap_angle` accepts floating-point types.
 * **Aspect ratio** must be **positive** and finite. For `aspect == 0` the FOV helpers are undefined.
 * **Units**: FOV functions accept/return **degrees** but compute internally in radians.
-* **Wrapping interval**: Behavior assumes `max > min`. The result lies in the **closed interval** `[min, max]` with modulo arithmetic; if you need half-open behavior (e.g., `[min,max)`), adjust your range or post-process endpoint cases.
+* **Wrapping interval**: Behavior assumes `max > min`. The result lies in the half-open interval `[min, max)`.
 * **constexpr**: Converters are `constexpr`; FOV helpers are runtime constexpr-compatible except for `std::atan/std::tan` constraints on some standard libraries.
 
 ---
@@ -103,5 +103,5 @@ float v = horizontal_fov_to_vertical(90.0f, 16.0f/9.0f);
 float h = vertical_fov_to_horizontal(v, 16.0f/9.0f);
 assert(std::abs(h - 90.0f) < 1e-5f);
 
-assert(wrap_angle(360.0f, 0.0f, 360.0f) == 0.0f || wrap_angle(360.0f, 0.0f, 360.0f) == 360.0f);
+assert(wrap_angle(360.0f, 0.0f, 360.0f) == 0.0f);
 ```

include/omath/engines/cry_engine/traits/camera_trait.hpp

@@ -4,7 +4,6 @@
 
 #pragma once
 #include "omath/engines/cry_engine/formulas.hpp"
-#include "omath/internal/constexpr_math.hpp"
 #include "omath/projection/camera.hpp"
 namespace omath::cry_engine
 {
@@ -16,8 +15,8 @@ namespace omath::cry_engine
                                                        const Vector3<float>& look_at) noexcept
         {
             const auto direction = (look_at - cam_origin).normalized();
-            return {PitchAngle::from_radians(internal::asin(direction.z)),
+            return {PitchAngle::from_asin(direction.z), -YawAngle::from_atan2(direction.x, direction.y),
-                    YawAngle::from_radians(-internal::atan2(direction.x, direction.y)), RollAngle::from_radians(0.f)};
+                    RollAngle::from_radians(0.f)};
         }
 
         [[nodiscard("view matrix result should not be discarded")]]

include/omath/engines/frostbite_engine/traits/camera_trait.hpp

@@ -4,7 +4,6 @@
 
 #pragma once
 #include "omath/engines/frostbite_engine/formulas.hpp"
-#include "omath/internal/constexpr_math.hpp"
 #include "omath/projection/camera.hpp"
 
 namespace omath::frostbite_engine
@@ -18,8 +17,8 @@ namespace omath::frostbite_engine
         {
             const auto direction = (look_at - cam_origin).normalized();
 
-            return {PitchAngle::from_radians(-internal::asin(direction.y)),
+            return {-PitchAngle::from_asin(direction.y), YawAngle::from_atan2(direction.x, direction.z),
-                    YawAngle::from_radians(internal::atan2(direction.x, direction.z)), RollAngle::from_radians(0.f)};
+                    RollAngle::from_radians(0.f)};
         }
 
         [[nodiscard("view matrix result should not be discarded")]]

include/omath/engines/iw_engine/traits/camera_trait.hpp

@@ -4,7 +4,6 @@
 
 #pragma once
 #include "omath/engines/iw_engine/formulas.hpp"
-#include "omath/internal/constexpr_math.hpp"
 #include "omath/projection/camera.hpp"
 
 namespace omath::iw_engine
@@ -18,8 +17,8 @@ namespace omath::iw_engine
         {
             const auto direction = (look_at - cam_origin).normalized();
 
-            return {PitchAngle::from_radians(-internal::asin(direction.z)),
+            return {-PitchAngle::from_asin(direction.z), YawAngle::from_atan2(direction.y, direction.x),
-                    YawAngle::from_radians(internal::atan2(direction.y, direction.x)), RollAngle::from_radians(0.f)};
+                    RollAngle::from_radians(0.f)};
         }
 
         [[nodiscard("view matrix result should not be discarded")]]

include/omath/engines/opengl_engine/traits/camera_trait.hpp

@@ -4,7 +4,6 @@
 
 #pragma once
 #include "omath/engines/opengl_engine/formulas.hpp"
-#include "omath/internal/constexpr_math.hpp"
 #include "omath/projection/camera.hpp"
 
 namespace omath::opengl_engine
@@ -18,8 +17,8 @@ namespace omath::opengl_engine
         {
             const auto direction = (look_at - cam_origin).normalized();
 
-            return {PitchAngle::from_radians(internal::asin(direction.y)),
+            return {PitchAngle::from_asin(direction.y), -YawAngle::from_atan2(direction.x, -direction.z),
-                    YawAngle::from_radians(-internal::atan2(direction.x, -direction.z)), RollAngle::from_radians(0.f)};
+                    RollAngle::from_radians(0.f)};
         }
 
         [[nodiscard("view matrix result should not be discarded")]]

include/omath/engines/rage_engine/traits/camera_trait.hpp

@@ -4,7 +4,6 @@
 
 #pragma once
 #include "omath/engines/rage_engine/formulas.hpp"
-#include "omath/internal/constexpr_math.hpp"
 #include "omath/projection/camera.hpp"
 
 namespace omath::rage_engine
@@ -18,8 +17,8 @@ namespace omath::rage_engine
         {
             const auto direction = (look_at - cam_origin).normalized();
 
-            return {PitchAngle::from_radians(internal::asin(direction.z)),
+            return {PitchAngle::from_asin(direction.z), -YawAngle::from_atan2(direction.x, direction.y),
-                    YawAngle::from_radians(-internal::atan2(direction.x, direction.y)), RollAngle::from_radians(0.f)};
+                    RollAngle::from_radians(0.f)};
         }
 
         [[nodiscard("view matrix result should not be discarded")]]

include/omath/engines/source_engine/traits/camera_trait.hpp

@@ -4,7 +4,6 @@
 
 #pragma once
 #include "omath/engines/source_engine/formulas.hpp"
-#include "omath/internal/constexpr_math.hpp"
 #include "omath/projection/camera.hpp"
 
 namespace omath::source_engine
@@ -18,8 +17,8 @@ namespace omath::source_engine
         {
             const auto direction = (look_at - cam_origin).normalized();
 
-            return {PitchAngle::from_radians(-internal::asin(direction.z)),
+            return {-PitchAngle::from_asin(direction.z),
-                    YawAngle::from_radians(internal::atan2(direction.y, direction.x)), RollAngle::from_radians(0.f)};
+                    YawAngle::from_atan2(direction.y, direction.x), RollAngle::from_radians(0.f)};
         }
 
         [[nodiscard("view matrix result should not be discarded")]]

include/omath/engines/source_engine/traits/pred_engine_trait.hpp

@@ -6,6 +6,8 @@
 #include "omath/engines/source_engine/formulas.hpp"
 #include "omath/projectile_prediction/projectile.hpp"
 #include "omath/projectile_prediction/target.hpp"
+#include <cmath>
+#include <limits>
 #include <optional>
 
 namespace omath::source_engine
@@ -14,19 +16,58 @@ namespace omath::source_engine
     {
     public:
         [[nodiscard("projectile position result should not be discarded")]]
-        constexpr static Vector3<float> predict_projectile_position(const projectile_prediction::Projectile<float>& projectile,
+        constexpr static Vector3<float>
-                                                                    const float pitch, const float yaw,
+        predict_projectile_position(const projectile_prediction::Projectile<float>& projectile, const float pitch,
-                                                                    const float time, const float gravity) noexcept
+                                    const float yaw, const float time, const float gravity) noexcept
         {
             const auto launch_pos = projectile.m_origin + projectile.m_launch_offset;
-            auto current_pos = launch_pos
+            const auto pitch_angle = PitchAngle::from_degrees(-pitch);
-                               + forward_vector({PitchAngle::from_degrees(-pitch), YawAngle::from_degrees(yaw),
+            const auto yaw_angle = YawAngle::from_degrees(yaw);
-                                                 RollAngle::from_degrees(0)})
+            const auto pitch_cos = pitch_angle.cos();
-                                         * projectile.m_launch_speed * time;
+            // Roll is always zero here, so this is the exact first column of the rotation matrix.
+            const Vector3 forward{pitch_cos * yaw_angle.cos(), pitch_cos * yaw_angle.sin(), -pitch_angle.sin()};
+            auto current_pos = launch_pos + forward * projectile.m_launch_speed * time;
             current_pos.z -= (gravity * projectile.m_gravity_scale) * (time * time) * 0.5f;
 
             return current_pos;
         }
+
+        [[nodiscard("reachability result should not be discarded")]]
+        static constexpr bool
+        can_projectile_reach_target_at_time(const projectile_prediction::Projectile<float>& projectile,
+                                            const Vector3<float>& target_position, const float time,
+                                            const float gravity, const float distance_tolerance) noexcept
+        {
+            if (!(distance_tolerance >= 0.f))
+                return false;
+
+            // After undoing gravity, every possible projectile position is on a sphere with radius speed * time.
+            const auto launch_position = projectile.m_origin + projectile.m_launch_offset;
+            auto adjusted_delta = target_position - launch_position;
+            const auto gravity_displacement = (gravity * projectile.m_gravity_scale) * (time * time) * 0.5f;
+            adjusted_delta.z += gravity_displacement;
+
+            const auto target_distance_sqr = adjusted_delta.length_sqr();
+            const auto projectile_distance = std::abs(projectile.m_launch_speed * time);
+            const auto floating_point_margin =
+                    std::numeric_limits<float>::epsilon() * 8.f
+                    * (std::abs(target_position.x) + std::abs(target_position.y) + std::abs(target_position.z)
+                       + std::abs(launch_position.x) + std::abs(launch_position.y) + std::abs(launch_position.z)
+                       + std::abs(gravity_displacement) + projectile_distance + distance_tolerance + 1.f);
+
+            if (!std::isfinite(target_distance_sqr) || !std::isfinite(floating_point_margin)) [[unlikely]]
+                return true;
+
+            const auto conservative_tolerance = distance_tolerance + floating_point_margin;
+            const auto maximum_distance = projectile_distance + conservative_tolerance;
+            if (target_distance_sqr > maximum_distance * maximum_distance)
+                return false;
+
+            const auto minimum_distance =
+                    projectile_distance > conservative_tolerance ? projectile_distance - conservative_tolerance : 0.f;
+            return target_distance_sqr >= minimum_distance * minimum_distance;
+        }
+
         [[nodiscard("target position result should not be discarded")]]
         static constexpr Vector3<float> predict_target_position(const projectile_prediction::Target<float>& target,
                                                                 const float time, const float gravity) noexcept

include/omath/engines/unity_engine/traits/camera_trait.hpp

@@ -4,7 +4,6 @@
 
 #pragma once
 #include "omath/engines/unity_engine/formulas.hpp"
-#include "omath/internal/constexpr_math.hpp"
 #include "omath/projection/camera.hpp"
 
 namespace omath::unity_engine
@@ -18,8 +17,8 @@ namespace omath::unity_engine
         {
             const auto direction = (look_at - cam_origin).normalized();
 
-            return {PitchAngle::from_radians(-internal::asin(direction.y)),
+            return {-PitchAngle::from_asin(direction.y), YawAngle::from_atan2(direction.x, direction.z),
-                    YawAngle::from_radians(internal::atan2(direction.x, direction.z)), RollAngle::from_radians(0.f)};
+                    RollAngle::from_radians(0.f)};
         }
 
         [[nodiscard("view matrix result should not be discarded")]]

include/omath/engines/unreal_engine/traits/camera_trait.hpp

@@ -4,7 +4,6 @@
 
 #pragma once
 #include "omath/engines/unreal_engine/formulas.hpp"
-#include "omath/internal/constexpr_math.hpp"
 #include "omath/projection/camera.hpp"
 
 namespace omath::unreal_engine
@@ -18,8 +17,8 @@ namespace omath::unreal_engine
         {
             const auto direction = (look_at - cam_origin).normalized();
 
-            return {PitchAngle::from_radians(internal::asin(direction.z)),
+            return {PitchAngle::from_asin(direction.z), YawAngle::from_atan2(direction.y, direction.x),
-                    YawAngle::from_radians(internal::atan2(direction.y, direction.x)), RollAngle::from_radians(0.f)};
+                    RollAngle::from_radians(0.f)};
         }
 
         [[nodiscard("view matrix result should not be discarded")]]

include/omath/hud/gradient.hpp

@@ -9,7 +9,7 @@ namespace omath::hud
         LeftToRight,
     };
 
-    struct Gradient
+    struct Gradient final
     {
         Color top_left;
         Color top_right;

include/omath/linear_algebra/mat.hpp

@@ -186,7 +186,14 @@ namespace omath
                 else if constexpr (StoreType == MatStoreType::COLUMN_MAJOR)
                     return cache_friendly_multiply_col_major(other);
             }
+            if constexpr (!std::is_same_v<Type, float> && !std::is_same_v<Type, double>)
+            {
                 if constexpr (StoreType == MatStoreType::ROW_MAJOR)
+                    return cache_friendly_multiply_row_major(other);
+                else if constexpr (StoreType == MatStoreType::COLUMN_MAJOR)
+                    return cache_friendly_multiply_col_major(other);
+            }
+            else if constexpr (StoreType == MatStoreType::ROW_MAJOR)
                 return avx_multiply_row_major(other);
             else if constexpr (StoreType == MatStoreType::COLUMN_MAJOR)
                 return avx_multiply_col_major(other);
@@ -429,12 +436,21 @@ namespace omath
         cache_friendly_multiply_row_major(const Mat<Columns, OtherColumns, Type, MatStoreType::ROW_MAJOR>& other) const
         {
             Mat<Rows, OtherColumns, Type, MatStoreType::ROW_MAJOR> result;
+            const Type* left_data = m_data.data();
+            const Type* right_data = other.raw_array().data();
+            Type* result_data = result.raw_array().data();
+
             for (std::size_t row_index = 0; row_index < Rows; ++row_index)
+            {
+                const Type* left_row = left_data + row_index * Columns;
+                Type* result_row = result_data + row_index * OtherColumns;
                 for (std::size_t column_index = 0; column_index < Columns; ++column_index)
                 {
-                    const Type& current_number = at(row_index, column_index);
+                    const Type current_number = left_row[column_index];
+                    const Type* right_row = right_data + column_index * OtherColumns;
                     for (std::size_t other_column = 0; other_column < OtherColumns; ++other_column)
-                        result.at(row_index, other_column) += current_number * other.at(column_index, other_column);
+                        result_row[other_column] += current_number * right_row[other_column];
+                }
             }
             return result;
         }
@@ -444,12 +460,21 @@ namespace omath
                 const Mat<Columns, OtherColumns, Type, MatStoreType::COLUMN_MAJOR>& other) const
         {
             Mat<Rows, OtherColumns, Type, MatStoreType::COLUMN_MAJOR> result;
+            const Type* left_data = m_data.data();
+            const Type* right_data = other.raw_array().data();
+            Type* result_data = result.raw_array().data();
+
             for (std::size_t other_column = 0; other_column < OtherColumns; ++other_column)
+            {
+                const Type* right_column = right_data + other_column * Columns;
+                Type* result_column = result_data + other_column * Rows;
                 for (std::size_t column_index = 0; column_index < Columns; ++column_index)
                 {
-                    const Type& current_number = other.at(column_index, other_column);
+                    const Type current_number = right_column[column_index];
+                    const Type* left_column = left_data + column_index * Rows;
                     for (std::size_t row_index = 0; row_index < Rows; ++row_index)
-                        result.at(row_index, other_column) += at(row_index, column_index) * current_number;
+                        result_column[row_index] += left_column[row_index] * current_number;
+                }
             }
             return result;
         }
@@ -466,56 +491,92 @@ namespace omath
 
             if constexpr (std::is_same_v<Type, float>)
             {
-                // ReSharper disable once CppTooWideScopeInitStatement
                 constexpr std::size_t vector_size = 8;
+                constexpr std::size_t block_size = vector_size * 4;
                 for (std::size_t j = 0; j < OtherColumns; ++j)
                 {
                     auto* c_col = reinterpret_cast<float*>(result_mat_data + j * Rows);
+                    std::size_t i = 0;
+                    for (; i + block_size <= Rows; i += block_size)
+                    {
+                        __m256 cvec0 = _mm256_setzero_ps();
+                        __m256 cvec1 = _mm256_setzero_ps();
+                        __m256 cvec2 = _mm256_setzero_ps();
+                        __m256 cvec3 = _mm256_setzero_ps();
                         for (std::size_t k = 0; k < Columns; ++k)
                         {
-                        const float bkj = reinterpret_cast<const float*>(other_mat_data)[k + j * Columns];
+                            const __m256 bkj_vec = _mm256_set1_ps(other_mat_data[k + j * Columns]);
-                        const __m256 bkj_vec = _mm256_set1_ps(bkj);
+                            const auto* a_col_k = this_mat_data + k * Rows + i;
-
+                            cvec0 = _mm256_fmadd_ps(_mm256_loadu_ps(a_col_k), bkj_vec, cvec0);
-                        const auto* a_col_k = reinterpret_cast<const float*>(this_mat_data + k * Rows);
+                            cvec1 = _mm256_fmadd_ps(_mm256_loadu_ps(a_col_k + vector_size), bkj_vec, cvec1);
-
+                            cvec2 = _mm256_fmadd_ps(_mm256_loadu_ps(a_col_k + vector_size * 2), bkj_vec, cvec2);
-                        std::size_t i = 0;
+                            cvec3 = _mm256_fmadd_ps(_mm256_loadu_ps(a_col_k + vector_size * 3), bkj_vec, cvec3);
+                        }
+                        _mm256_storeu_ps(c_col + i, cvec0);
+                        _mm256_storeu_ps(c_col + i + vector_size, cvec1);
+                        _mm256_storeu_ps(c_col + i + vector_size * 2, cvec2);
+                        _mm256_storeu_ps(c_col + i + vector_size * 3, cvec3);
+                    }
                     for (; i + vector_size <= Rows; i += vector_size)
                     {
-                            __m256 cvec = _mm256_loadu_ps(c_col + i);
+                        __m256 cvec = _mm256_setzero_ps();
+                        for (std::size_t k = 0; k < Columns; ++k)
+                        {
+                            const __m256 bkj_vec = _mm256_set1_ps(other_mat_data[k + j * Columns]);
+                            const auto* a_col_k = this_mat_data + k * Rows;
                             const __m256 a_vec = _mm256_loadu_ps(a_col_k + i);
                             cvec = _mm256_fmadd_ps(a_vec, bkj_vec, cvec);
+                        }
                         _mm256_storeu_ps(c_col + i, cvec);
                     }
                     for (; i < Rows; ++i)
-                            c_col[i] += a_col_k[i] * bkj;
+                        for (std::size_t k = 0; k < Columns; ++k)
-                    }
+                            c_col[i] += this_mat_data[i + k * Rows] * other_mat_data[k + j * Columns];
                 }
             }
             else if (std::is_same_v<Type, double>)
-            { // double
+            {
-                // ReSharper disable once CppTooWideScopeInitStatement
                 constexpr std::size_t vector_size = 4;
+                constexpr std::size_t block_size = vector_size * 4;
                 for (std::size_t j = 0; j < OtherColumns; ++j)
                 {
                     auto* c_col = reinterpret_cast<double*>(result_mat_data + j * Rows);
+                    std::size_t i = 0;
+                    for (; i + block_size <= Rows; i += block_size)
+                    {
+                        __m256d cvec0 = _mm256_setzero_pd();
+                        __m256d cvec1 = _mm256_setzero_pd();
+                        __m256d cvec2 = _mm256_setzero_pd();
+                        __m256d cvec3 = _mm256_setzero_pd();
                         for (std::size_t k = 0; k < Columns; ++k)
                         {
-                        const double bkj = reinterpret_cast<const double*>(other_mat_data)[k + j * Columns];
+                            const __m256d bkj_vec = _mm256_set1_pd(other_mat_data[k + j * Columns]);
-                        const __m256d bkj_vec = _mm256_set1_pd(bkj);
+                            const auto* a_col_k = this_mat_data + k * Rows + i;
-
+                            cvec0 = _mm256_fmadd_pd(_mm256_loadu_pd(a_col_k), bkj_vec, cvec0);
-                        const auto* a_col_k = reinterpret_cast<const double*>(this_mat_data + k * Rows);
+                            cvec1 = _mm256_fmadd_pd(_mm256_loadu_pd(a_col_k + vector_size), bkj_vec, cvec1);
-
+                            cvec2 = _mm256_fmadd_pd(_mm256_loadu_pd(a_col_k + vector_size * 2), bkj_vec, cvec2);
-                        std::size_t i = 0;
+                            cvec3 = _mm256_fmadd_pd(_mm256_loadu_pd(a_col_k + vector_size * 3), bkj_vec, cvec3);
+                        }
+                        _mm256_storeu_pd(c_col + i, cvec0);
+                        _mm256_storeu_pd(c_col + i + vector_size, cvec1);
+                        _mm256_storeu_pd(c_col + i + vector_size * 2, cvec2);
+                        _mm256_storeu_pd(c_col + i + vector_size * 3, cvec3);
+                    }
                     for (; i + vector_size <= Rows; i += vector_size)
                     {
-                            __m256d cvec = _mm256_loadu_pd(c_col + i);
+                        __m256d cvec = _mm256_setzero_pd();
+                        for (std::size_t k = 0; k < Columns; ++k)
+                        {
+                            const __m256d bkj_vec = _mm256_set1_pd(other_mat_data[k + j * Columns]);
+                            const auto* a_col_k = this_mat_data + k * Rows;
                             const __m256d a_vec = _mm256_loadu_pd(a_col_k + i);
                             cvec = _mm256_fmadd_pd(a_vec, bkj_vec, cvec);
+                        }
                         _mm256_storeu_pd(c_col + i, cvec);
                     }
                     for (; i < Rows; ++i)
-                            c_col[i] += a_col_k[i] * bkj;
+                        for (std::size_t k = 0; k < Columns; ++k)
-                    }
+                            c_col[i] += this_mat_data[i + k * Rows] * other_mat_data[k + j * Columns];
                 }
             }
             else
@@ -536,56 +597,92 @@ namespace omath
 
             if constexpr (std::is_same_v<Type, float>)
             {
-                // ReSharper disable once CppTooWideScopeInitStatement
                 constexpr std::size_t vector_size = 8;
+                constexpr std::size_t block_size = vector_size * 4;
                 for (std::size_t i = 0; i < Rows; ++i)
                 {
-                    Type* c_row = result_mat_data + i * OtherColumns;
+                    auto* c_row = reinterpret_cast<float*>(result_mat_data + i * OtherColumns);
+                    std::size_t j = 0;
+                    for (; j + block_size <= OtherColumns; j += block_size)
+                    {
+                        __m256 cvec0 = _mm256_setzero_ps();
+                        __m256 cvec1 = _mm256_setzero_ps();
+                        __m256 cvec2 = _mm256_setzero_ps();
+                        __m256 cvec3 = _mm256_setzero_ps();
                         for (std::size_t k = 0; k < Columns; ++k)
                         {
-                        const auto aik = static_cast<float>(this_mat_data[i * Columns + k]);
+                            const __m256 aik_vec = _mm256_set1_ps(this_mat_data[i * Columns + k]);
-                        const __m256 aik_vec = _mm256_set1_ps(aik);
+                            const auto* b_row = other_mat_data + k * OtherColumns + j;
-                        const auto* b_row = reinterpret_cast<const float*>(other_mat_data + k * OtherColumns);
+                            cvec0 = _mm256_fmadd_ps(_mm256_loadu_ps(b_row), aik_vec, cvec0);
-
+                            cvec1 = _mm256_fmadd_ps(_mm256_loadu_ps(b_row + vector_size), aik_vec, cvec1);
-                        std::size_t j = 0;
+                            cvec2 = _mm256_fmadd_ps(_mm256_loadu_ps(b_row + vector_size * 2), aik_vec, cvec2);
+                            cvec3 = _mm256_fmadd_ps(_mm256_loadu_ps(b_row + vector_size * 3), aik_vec, cvec3);
+                        }
+                        _mm256_storeu_ps(c_row + j, cvec0);
+                        _mm256_storeu_ps(c_row + j + vector_size, cvec1);
+                        _mm256_storeu_ps(c_row + j + vector_size * 2, cvec2);
+                        _mm256_storeu_ps(c_row + j + vector_size * 3, cvec3);
+                    }
                     for (; j + vector_size <= OtherColumns; j += vector_size)
                     {
-                            __m256 cvec = _mm256_loadu_ps(c_row + j);
+                        __m256 cvec = _mm256_setzero_ps();
+                        for (std::size_t k = 0; k < Columns; ++k)
+                        {
+                            const __m256 aik_vec = _mm256_set1_ps(this_mat_data[i * Columns + k]);
+                            const auto* b_row = other_mat_data + k * OtherColumns;
                             const __m256 b_vec = _mm256_loadu_ps(b_row + j);
                             cvec = _mm256_fmadd_ps(b_vec, aik_vec, cvec);
-
+                        }
                         _mm256_storeu_ps(c_row + j, cvec);
                     }
                     for (; j < OtherColumns; ++j)
-                            c_row[j] += aik * b_row[j];
+                        for (std::size_t k = 0; k < Columns; ++k)
-                    }
+                            c_row[j] += this_mat_data[i * Columns + k] * other_mat_data[k * OtherColumns + j];
                 }
             }
             else if (std::is_same_v<Type, double>)
-            { // double
+            {
-                // ReSharper disable once CppTooWideScopeInitStatement
                 constexpr std::size_t vector_size = 4;
+                constexpr std::size_t block_size = vector_size * 4;
                 for (std::size_t i = 0; i < Rows; ++i)
                 {
-                    Type* c_row = result_mat_data + i * OtherColumns;
+                    auto* c_row = reinterpret_cast<double*>(result_mat_data + i * OtherColumns);
+                    std::size_t j = 0;
+                    for (; j + block_size <= OtherColumns; j += block_size)
+                    {
+                        __m256d cvec0 = _mm256_setzero_pd();
+                        __m256d cvec1 = _mm256_setzero_pd();
+                        __m256d cvec2 = _mm256_setzero_pd();
+                        __m256d cvec3 = _mm256_setzero_pd();
                         for (std::size_t k = 0; k < Columns; ++k)
                         {
-                        const auto aik = static_cast<double>(this_mat_data[i * Columns + k]);
+                            const __m256d aik_vec = _mm256_set1_pd(this_mat_data[i * Columns + k]);
-                        const __m256d aik_vec = _mm256_set1_pd(aik);
+                            const auto* b_row = other_mat_data + k * OtherColumns + j;
-                        const auto* b_row = reinterpret_cast<const double*>(other_mat_data + k * OtherColumns);
+                            cvec0 = _mm256_fmadd_pd(_mm256_loadu_pd(b_row), aik_vec, cvec0);
-
+                            cvec1 = _mm256_fmadd_pd(_mm256_loadu_pd(b_row + vector_size), aik_vec, cvec1);
-                        std::size_t j = 0;
+                            cvec2 = _mm256_fmadd_pd(_mm256_loadu_pd(b_row + vector_size * 2), aik_vec, cvec2);
+                            cvec3 = _mm256_fmadd_pd(_mm256_loadu_pd(b_row + vector_size * 3), aik_vec, cvec3);
+                        }
+                        _mm256_storeu_pd(c_row + j, cvec0);
+                        _mm256_storeu_pd(c_row + j + vector_size, cvec1);
+                        _mm256_storeu_pd(c_row + j + vector_size * 2, cvec2);
+                        _mm256_storeu_pd(c_row + j + vector_size * 3, cvec3);
+                    }
                     for (; j + vector_size <= OtherColumns; j += vector_size)
                     {
-                            __m256d cvec = _mm256_loadu_pd(c_row + j);
+                        __m256d cvec = _mm256_setzero_pd();
+                        for (std::size_t k = 0; k < Columns; ++k)
+                        {
+                            const __m256d aik_vec = _mm256_set1_pd(this_mat_data[i * Columns + k]);
+                            const auto* b_row = other_mat_data + k * OtherColumns;
                             const __m256d b_vec = _mm256_loadu_pd(b_row + j);
                             cvec = _mm256_fmadd_pd(b_vec, aik_vec, cvec);
-
+                        }
                         _mm256_storeu_pd(c_row + j, cvec);
                     }
                     for (; j < OtherColumns; ++j)
-                            c_row[j] += aik * b_row[j];
+                        for (std::size_t k = 0; k < Columns; ++k)
-                    }
+                            c_row[j] += this_mat_data[i * Columns + k] * other_mat_data[k * OtherColumns + j];
                 }
             }
             else

include/omath/projectile_prediction/proj_pred_engine_legacy.hpp

@@ -9,6 +9,7 @@
 #include "omath/projectile_prediction/proj_pred_engine.hpp"
 #include "omath/projectile_prediction/projectile.hpp"
 #include "omath/projectile_prediction/target.hpp"
+#include <cmath>
 #include <optional>
 
 namespace omath::projectile_prediction
@@ -17,9 +18,8 @@ namespace omath::projectile_prediction
     concept PredEngineConcept =
             requires(const Projectile<ArithmeticType>& projectile, const Target<ArithmeticType>& target,
                      const Vector3<ArithmeticType>& vec_a, const Vector3<ArithmeticType>& vec_b,
-                     Vector3<ArithmeticType> v3,
+                     Vector3<ArithmeticType> v3, ArithmeticType pitch, ArithmeticType yaw, ArithmeticType time,
-                     ArithmeticType pitch, ArithmeticType yaw, ArithmeticType time, ArithmeticType gravity,
+                     ArithmeticType gravity, std::optional<ArithmeticType> maybe_pitch) {
-                     std::optional<ArithmeticType> maybe_pitch) {
                 {
                     T::predict_projectile_position(projectile, pitch, yaw, time, gravity)
                 } -> std::same_as<Vector3<ArithmeticType>>;
@@ -44,8 +44,7 @@ namespace omath::projectile_prediction
     class ProjPredEngineLegacy final : public ProjPredEngineInterface<ArithmeticType>
     {
     public:
-        explicit ProjPredEngineLegacy(const ArithmeticType gravity_constant,
+        explicit ProjPredEngineLegacy(const ArithmeticType gravity_constant, const ArithmeticType simulation_time_step,
-                                      const ArithmeticType simulation_time_step,
                                       const ArithmeticType maximum_simulation_time,
                                       const ArithmeticType distance_tolerance)
             : m_gravity_constant(gravity_constant), m_simulation_time_step(simulation_time_step),
@@ -54,8 +53,9 @@ namespace omath::projectile_prediction
         }
 
         [[nodiscard]]
-        std::optional<Vector3<ArithmeticType>> maybe_calculate_aim_point(
+        std::optional<Vector3<ArithmeticType>>
-            const Projectile<ArithmeticType>& projectile, const Target<ArithmeticType>& target) const override
+        maybe_calculate_aim_point(const Projectile<ArithmeticType>& projectile,
+                                  const Target<ArithmeticType>& target) const override
         {
             const auto solution = find_solution(projectile, target);
             if (!solution)
@@ -66,15 +66,16 @@ namespace omath::projectile_prediction
         }
 
         [[nodiscard]]
-        std::optional<AimAngles<ArithmeticType>> maybe_calculate_aim_angles(
+        std::optional<AimAngles<ArithmeticType>>
-            const Projectile<ArithmeticType>& projectile, const Target<ArithmeticType>& target) const override
+        maybe_calculate_aim_angles(const Projectile<ArithmeticType>& projectile,
+                                   const Target<ArithmeticType>& target) const override
         {
             const auto solution = find_solution(projectile, target);
             if (!solution)
                 return std::nullopt;
 
-            const auto yaw = EngineTrait::calc_direct_yaw_angle(
+            const auto yaw = EngineTrait::calc_direct_yaw_angle(projectile.m_origin + projectile.m_launch_offset,
-                projectile.m_origin + projectile.m_launch_offset, solution->predicted_target_position);
+                                                                solution->predicted_target_position);
             return AimAngles<ArithmeticType>{solution->pitch, yaw};
         }
 
@@ -89,24 +90,40 @@ namespace omath::projectile_prediction
         std::optional<Solution> find_solution(const Projectile<ArithmeticType>& projectile,
                                               const Target<ArithmeticType>& target) const
         {
-            for (ArithmeticType time = ArithmeticType{0}; time < m_maximum_simulation_time;
+            if (!std::isfinite(m_simulation_time_step) || m_simulation_time_step <= ArithmeticType{0}
-                 time += m_simulation_time_step)
+                || !std::isfinite(m_maximum_simulation_time) || m_maximum_simulation_time < ArithmeticType{0}
+                || !std::isfinite(projectile.m_launch_speed) || projectile.m_launch_speed <= ArithmeticType{0}
+                || !(m_distance_tolerance >= ArithmeticType{0})) [[unlikely]]
+                return std::nullopt;
+
+            for (ArithmeticType time = ArithmeticType{0}; time <= m_maximum_simulation_time;)
             {
                 const auto predicted_target_position =
                         EngineTrait::predict_target_position(target, time, m_gravity_constant);
 
+                if (is_target_potentially_reachable(projectile, predicted_target_position, time))
+                {
                     const auto projectile_pitch =
                             maybe_calculate_projectile_launch_pitch_angle(projectile, predicted_target_position);
 
-                if (!projectile_pitch.has_value()) [[unlikely]]
+                    if (projectile_pitch.has_value()) [[likely]]
-                    continue;
+                    {
-
+                        const auto yaw = EngineTrait::calc_direct_yaw_angle(
-                if (!is_projectile_reached_target(predicted_target_position, projectile, projectile_pitch.value(),
+                                projectile.m_origin + projectile.m_launch_offset, predicted_target_position);
-                                                  time))
-                    continue;
 
+                        if (is_projectile_reached_target(predicted_target_position, projectile,
+                                                         projectile_pitch.value(), yaw, time))
                             return Solution{predicted_target_position, projectile_pitch.value()};
                     }
+                }
+
+                if (time == m_maximum_simulation_time)
+                    break;
+                const auto next_time = time + m_simulation_time_step;
+                if (!(next_time > time)) [[unlikely]]
+                    break;
+                time = next_time < m_maximum_simulation_time ? next_time : m_maximum_simulation_time;
+            }
             return std::nullopt;
         }
 
@@ -144,29 +161,46 @@ namespace omath::projectile_prediction
             const auto distance2d = EngineTrait::calc_vector_2d_distance(delta);
             const auto distance2d_sqr = distance2d * distance2d;
             const auto launch_speed_sqr = projectile.m_launch_speed * projectile.m_launch_speed;
+            const auto ballistic_term =
+                    bullet_gravity * distance2d_sqr
+                    + ArithmeticType{2} * EngineTrait::get_vector_height_coordinate(delta) * launch_speed_sqr;
 
-            ArithmeticType root = launch_speed_sqr * launch_speed_sqr
+            ArithmeticType root = launch_speed_sqr * launch_speed_sqr - bullet_gravity * ballistic_term;
-                                  - bullet_gravity
-                                            * (bullet_gravity * distance2d_sqr
-                                               + ArithmeticType{2} * EngineTrait::get_vector_height_coordinate(delta)
-                                                         * launch_speed_sqr);
 
             if (root < ArithmeticType{0}) [[unlikely]]
                 return std::nullopt;
 
             root = std::sqrt(root);
-            const ArithmeticType angle = std::atan((launch_speed_sqr - root) / (bullet_gravity * distance2d));
+            // This rationalized form avoids cancellation in launch_speed_sqr - root for low-angle shots.
+            const ArithmeticType angle = std::atan2(ballistic_term, distance2d * (launch_speed_sqr + root));
 
             return angles::radians_to_degrees(angle);
         }
 
         [[nodiscard]]
-        bool is_projectile_reached_target(const Vector3<ArithmeticType>& target_position,
+        bool is_target_potentially_reachable(const Projectile<ArithmeticType>& projectile,
-                                          const Projectile<ArithmeticType>& projectile,
+                                             const Vector3<ArithmeticType>& target_position,
-                                          const ArithmeticType pitch, const ArithmeticType time) const noexcept
+                                             const ArithmeticType time) const noexcept
+        {
+            if constexpr (requires {
+                              {
+                                  EngineTrait::can_projectile_reach_target_at_time(
+                                          projectile, target_position, time, m_gravity_constant, m_distance_tolerance)
+                              } -> std::same_as<bool>;
+                              requires noexcept(EngineTrait::can_projectile_reach_target_at_time(
+                                      projectile, target_position, time, m_gravity_constant, m_distance_tolerance));
+                          })
+                return EngineTrait::can_projectile_reach_target_at_time(projectile, target_position, time,
+                                                                        m_gravity_constant, m_distance_tolerance);
+
+            return true;
+        }
+
+        [[nodiscard]]
+        bool is_projectile_reached_target(const Vector3<ArithmeticType>& target_position,
+                                          const Projectile<ArithmeticType>& projectile, const ArithmeticType pitch,
+                                          const ArithmeticType yaw, const ArithmeticType time) const noexcept
         {
-            const auto yaw = EngineTrait::calc_direct_yaw_angle(
-                projectile.m_origin + projectile.m_launch_offset, target_position);
             const auto projectile_position =
                     EngineTrait::predict_projectile_position(projectile, pitch, yaw, time, m_gravity_constant);
 

include/omath/trigonometry/angle.hpp

@@ -6,7 +6,9 @@
 #include "omath/internal/constexpr_math.hpp"
 #include "omath/trigonometry/angles.hpp"
 #include <algorithm>
+#include <compare>
 #include <format>
+#include <type_traits>
 #include <utility>
 
 namespace omath
@@ -18,7 +20,7 @@ namespace omath
     };
 
     template<class Type = float, Type min = Type(0), Type max = Type(360), AngleFlags flags = AngleFlags::Normalized>
-    requires std::is_arithmetic_v<Type>
+    requires std::is_floating_point_v<Type>
     class Angle
     {
         Type m_angle;
@@ -43,7 +45,7 @@ namespace omath
         {
             return Angle{degrees};
         }
-        constexpr Angle() noexcept: m_angle(0)
+        constexpr Angle() noexcept: Angle(Type{0})
         {
         }
         [[nodiscard]]
@@ -52,6 +54,30 @@ namespace omath
             return Angle{angles::radians_to_degrees<Type>(degrees)};
         }
 
+        [[nodiscard]]
+        constexpr static Angle from_asin(const Type& value) noexcept
+        {
+            return from_radians(internal::asin(value));
+        }
+
+        [[nodiscard]]
+        constexpr static Angle from_acos(const Type& value) noexcept
+        {
+            return from_radians(internal::acos(value));
+        }
+
+        [[nodiscard]]
+        constexpr static Angle from_atan(const Type& value) noexcept
+        {
+            return from_radians(internal::atan(value));
+        }
+
+        [[nodiscard]]
+        constexpr static Angle from_atan2(const Type& y, const Type& x) noexcept
+        {
+            return from_radians(internal::atan2(y, x));
+        }
+
         [[nodiscard]]
         constexpr const Type& operator*() const noexcept
         {
@@ -88,12 +114,6 @@ namespace omath
             return internal::tan(as_radians());
         }
 
-        [[nodiscard]]
-        constexpr Type atan() const noexcept
-        {
-            return internal::atan(as_radians());
-        }
-
         [[nodiscard]]
         constexpr Type cot() const noexcept
         {
@@ -121,7 +141,8 @@ namespace omath
 
         constexpr Angle& operator-=(const Angle& other) noexcept
         {
-            return operator+=(-other);
+            *this = Angle{m_angle - other.m_angle};
+            return *this;
         }
 
         [[nodiscard]]
@@ -142,7 +163,7 @@ namespace omath
         [[nodiscard]]
         constexpr Angle operator-(const Angle& other) const noexcept
         {
-            return operator+(-other);
+            return Angle{m_angle - other.m_angle};
         }
 
         [[nodiscard]]
@@ -172,7 +193,6 @@ struct std::formatter<omath::Angle<T, MinV, MaxV, F>, char> final // NOLINT(*-dc
         return std::format_to(ctx.out(), "{}deg", a.as_degrees());
     }
 };
-
 // wchar_t formatter
 template<class T, T MinV, T MaxV, omath::AngleFlags F>
 struct std::formatter<omath::Angle<T, MinV, MaxV, F>, wchar_t> final // NOLINT(*-dcl58-cpp)
@@ -193,24 +213,3 @@ struct std::formatter<omath::Angle<T, MinV, MaxV, F>, wchar_t> final // NOLINT(*
         return std::format_to(ctx.out(), L"{}deg", a.as_degrees());
     }
 };
-
-// wchar_t formatter
-template<class T, T MinV, T MaxV, omath::AngleFlags F>
-struct std::formatter<omath::Angle<T, MinV, MaxV, F>, char8_t> final // NOLINT(*-dcl58-cpp)
-{
-    using AngleT = omath::Angle<T, MinV, MaxV, F>;
-
-    [[nodiscard]]
-    static constexpr auto parse(std::wformat_parse_context& ctx)
-    {
-        return ctx.begin();
-    }
-
-    template<class FormatContext>
-    [[nodiscard]]
-    auto format(const AngleT& a, FormatContext& ctx) const
-    {
-        static_assert(std::is_same_v<typename FormatContext::char_type, char8_t>);
-        return std::format_to(ctx.out(), u8"{}deg", a.as_degrees());
-    }
-};

include/omath/trigonometry/angles.hpp

@@ -48,10 +48,10 @@ namespace omath::angles
     }
 
     template<class Type>
-    requires std::is_arithmetic_v<Type>
+    requires std::is_floating_point_v<Type>
     [[nodiscard]] constexpr Type wrap_angle(const Type& angle, const Type& min, const Type& max) noexcept
     {
-        if (angle <= max && angle >= min)
+        if (angle < max && angle >= min)
             return angle;
 
         const Type range = max - min;

scripts/coverage-llvm.sh

@@ -177,9 +177,7 @@ if command -v genhtml >/dev/null 2>&1; then
         --title "Omath Coverage Report" \
         --show-details \
         --legend \
-        --demangle-cpp \
         --num-spaces 4 \
-        --sort \
         --function-coverage \
         --branch-coverage
     

tests/general/unit_test_angle.cpp

@@ -14,9 +14,14 @@ namespace
 
     // Handy aliases (defaults: Type=float, [0,360], Normalized)
     using Deg = Angle<float, static_cast<float>(0), static_cast<float>(360), AngleFlags::Normalized>;
+    using Fov = Angle<float, static_cast<float>(0), static_cast<float>(180), AngleFlags::Clamped>;
+    using Offset = Angle<float, static_cast<float>(10), static_cast<float>(20), AngleFlags::Clamped>;
     using Pitch = Angle<float, static_cast<float>(-90), static_cast<float>(90), AngleFlags::Clamped>;
     using Turn = Angle<float, static_cast<float>(-180), static_cast<float>(180), AngleFlags::Normalized>;
 
+    template<class Type>
+    concept SupportedAngleType = requires { typename Angle<Type>; };
+
     constexpr float k_eps = 1e-5f;
 
     constexpr bool close_to(const float actual, const float expected, const float epsilon)
@@ -36,6 +41,12 @@ TEST(UnitTestAngle, DefaultConstructor_IsZeroDegrees)
     EXPECT_FLOAT_EQ(a.as_degrees(), 0.0f);
 }
 
+TEST(UnitTestAngle, DefaultConstructor_AppliesRangePolicy)
+{
+    constexpr Offset a;
+    EXPECT_FLOAT_EQ(a.as_degrees(), 10.0f);
+}
+
 TEST(UnitTestAngle, FromDegrees_Normalized_WrapsAboveMax)
 {
     const Deg a = Deg::from_degrees(370.0f);
@@ -66,6 +77,14 @@ TEST(UnitTestAngle, FromRadians_And_AsRadians)
     EXPECT_NEAR(b.as_radians(), std::numbers::pi_v<float>, 1e-6f);
 }
 
+TEST(UnitTestAngle, FromInverseTrigonometricFunctions)
+{
+    EXPECT_NEAR(Pitch::from_asin(0.5f).as_degrees(), 30.0f, k_eps);
+    EXPECT_NEAR(Pitch::from_acos(0.5f).as_degrees(), 60.0f, k_eps);
+    EXPECT_NEAR(Pitch::from_atan(1.0f).as_degrees(), 45.0f, k_eps);
+    EXPECT_NEAR(Turn::from_atan2(-1.0f, -1.0f).as_degrees(), -135.0f, k_eps);
+}
+
 // ---------- Unary minus & deref ----------
 
 TEST(UnitTestAngle, UnaryMinus_Normalized)
@@ -101,17 +120,6 @@ TEST(UnitTestAngle, SinCosTanCot_BasicCases)
     EXPECT_NEAR(a90.cos(), 0.0f, 1e-4f);
 }
 
-TEST(UnitTestAngle, Atan_IsAtanOfRadians)
-{
-    // atan(as_radians). For 0° -> atan(0)=0.
-    const Deg a0 = Deg::from_degrees(0.0f);
-    EXPECT_NEAR(a0.atan(), 0.0f, k_eps);
-
-    const Deg a45 = Deg::from_degrees(45.0f);
-    // atan(pi/4) ≈ 0.665773...
-    EXPECT_NEAR(a45.atan(), 0.66577375f, 1e-6f);
-}
-
 // ---------- Compound arithmetic ----------
 
 TEST(UnitTestAngle, PlusEquals_Normalized_Wraps)
@@ -142,6 +150,16 @@ TEST(UnitTestAngle, MinusEquals_Clamped_Clamps)
     EXPECT_FLOAT_EQ(p.as_degrees(), -90.0f);
 }
 
+TEST(UnitTestAngle, Subtraction_ClampedNonSymmetricRange)
+{
+    Fov compound = Fov::from_degrees(90.0f);
+    compound -= Fov::from_degrees(10.0f);
+    EXPECT_FLOAT_EQ(compound.as_degrees(), 80.0f);
+
+    const Fov binary = Fov::from_degrees(90.0f) - Fov::from_degrees(10.0f);
+    EXPECT_FLOAT_EQ(binary.as_degrees(), 80.0f);
+}
+
 // ---------- Alternative ranges ----------
 
 TEST(UnitTestAngle, NormalizedRange_Neg180To180)
@@ -205,5 +223,12 @@ static_assert(close_to(Pitch::from_degrees(45.0f).tan(), 1.0f, 1e-4f),
               "Tan should be constexpr with embedded constexpr math");
 static_assert(close_to(Pitch::from_degrees(45.0f).cot(), 1.0f, 1e-4f),
               "Cot should be constexpr with embedded constexpr math");
-static_assert(close_to(Pitch::from_degrees(45.0f).atan(), 0.66577375f, 1e-6f),
+static_assert(close_to(Pitch::from_asin(0.5f).as_degrees(), 30.0f, k_eps),
-              "Atan should be constexpr with embedded constexpr math");
+              "From asin should be constexpr with embedded constexpr math");
+static_assert(close_to(Pitch::from_acos(0.5f).as_degrees(), 60.0f, k_eps),
+              "From acos should be constexpr with embedded constexpr math");
+static_assert(close_to(Pitch::from_atan(1.0f).as_degrees(), 45.0f, k_eps),
+              "From atan should be constexpr with embedded constexpr math");
+static_assert(close_to(Turn::from_atan2(-1.0f, -1.0f).as_degrees(), -135.0f, k_eps),
+              "From atan2 should be constexpr with embedded constexpr math");
+static_assert(!SupportedAngleType<int>, "Angle should only accept floating-point types");

tests/general/unit_test_angles.cpp

@@ -47,3 +47,10 @@ TEST(unit_test_angles, wrap_angle_negative_range)
 
     EXPECT_NEAR(wrapped, 270.f, 0.01f);
 }
+
+TEST(unit_test_angles, wrap_angle_maximum_maps_to_minimum)
+{
+    const float wrapped = omath::angles::wrap_angle(360.f, 0.f, 360.f);
+
+    EXPECT_FLOAT_EQ(wrapped, 0.f);
+}

tests/general/unit_test_mat.cpp

@@ -16,6 +16,21 @@ namespace
         const float diff = actual - expected;
         return (diff < 0.0f ? -diff : diff) <= epsilon;
     }
+
+    template<size_t Rows, size_t Columns, size_t OtherColumns, class Type, MatStoreType StoreType>
+    void expect_multiplication_matches_scalar_reference(const Mat<Rows, Columns, Type, StoreType>& left,
+                                                        const Mat<Columns, OtherColumns, Type, StoreType>& right)
+    {
+        const auto result = left * right;
+        for (size_t row = 0; row < Rows; ++row)
+            for (size_t column = 0; column < OtherColumns; ++column)
+            {
+                Type expected{};
+                for (size_t shared_index = 0; shared_index < Columns; ++shared_index)
+                    expected += left.at(row, shared_index) * right.at(shared_index, column);
+                EXPECT_EQ(result.at(row, column), expected);
+            }
+    }
 } // namespace
 
 class UnitTestMat : public ::testing::Test
@@ -92,6 +107,45 @@ TEST_F(UnitTestMat, Operator_Multiplication_Matrix)
     EXPECT_FLOAT_EQ(m3.at(1, 1), 22.0f);
 }
 
+TEST(UnitTestMatStandalone, Operator_Multiplication_RowMajorSimdAndTail)
+{
+    Mat<3, 5, float, MatStoreType::ROW_MAJOR> left;
+    Mat<5, 33, float, MatStoreType::ROW_MAJOR> right;
+    for (size_t row = 0; row < left.row_count(); ++row)
+        for (size_t column = 0; column < left.columns_count(); ++column)
+            left.at(row, column) = static_cast<float>(row * 3 + column + 1);
+    for (size_t row = 0; row < right.row_count(); ++row)
+        for (size_t column = 0; column < right.columns_count(); ++column)
+            right.at(row, column) = static_cast<float>((row + 1) * (column % 5 + 1));
+
+    expect_multiplication_matches_scalar_reference(left, right);
+}
+
+TEST(UnitTestMatStandalone, Operator_Multiplication_ColumnMajorSimdAndTail)
+{
+    Mat<17, 5, double, MatStoreType::COLUMN_MAJOR> left;
+    Mat<5, 3, double, MatStoreType::COLUMN_MAJOR> right;
+    for (size_t row = 0; row < left.row_count(); ++row)
+        for (size_t column = 0; column < left.columns_count(); ++column)
+            left.at(row, column) = static_cast<double>(row * 3 + column + 1);
+    for (size_t row = 0; row < right.row_count(); ++row)
+        for (size_t column = 0; column < right.columns_count(); ++column)
+            right.at(row, column) = static_cast<double>((row + 1) * (column + 1));
+
+    expect_multiplication_matches_scalar_reference(left, right);
+}
+
+TEST(UnitTestMatStandalone, Operator_Multiplication_IntegerFallsBackFromAvx)
+{
+    constexpr Mat<2, 3, int> left{{1, 2, 3}, {4, 5, 6}};
+    constexpr Mat<3, 2, int> right{{7, 8}, {9, 10}, {11, 12}};
+    constexpr auto result = left * right;
+    static_assert(result.at(0, 0) == 58);
+    static_assert(result.at(1, 1) == 154);
+
+    expect_multiplication_matches_scalar_reference(left, right);
+}
+
 TEST_F(UnitTestMat, Operator_Multiplication_Scalar)
 {
     Mat<2, 2> m3 = m2 * 2.0f;

tests/general/unit_test_pred_engine_trait.cpp

@@ -1,4 +1,5 @@
 // Tests for PredEngineTrait
+#include <array>
 #include <gtest/gtest.h>
 #include <omath/engines/source_engine/traits/pred_engine_trait.hpp>
 #include <omath/projectile_prediction/projectile.hpp>
@@ -77,6 +78,77 @@ TEST(PredEngineTrait, PredictProjectilePositionWithLaunchOffset)
     EXPECT_NEAR(pos_t1.z, -2.f - 9.81f * 0.5f, 1e-3f);
 }
 
+TEST(PredEngineTrait, PredictProjectilePositionMatchesRotationMatrix)
+{
+    constexpr Projectile projectile{
+            .m_origin = {10.f, -20.f, 30.f},
+            .m_launch_offset = {2.f, 3.f, -4.f},
+            .m_launch_speed = 750.f,
+            .m_gravity_scale = 0.6f,
+    };
+    struct TestCase
+    {
+        float pitch;
+        float yaw;
+        float time;
+    };
+    constexpr std::array test_cases{
+            TestCase{0.f, 0.f, 0.f},
+            TestCase{25.f, 45.f, 0.25f},
+            TestCase{-60.f, -135.f, 1.5f},
+            TestCase{120.f, 540.f, 2.f},
+    };
+    constexpr float gravity = 9.81f;
+
+    for (const auto& test_case : test_cases)
+    {
+        const auto launch_position = projectile.m_origin + projectile.m_launch_offset;
+        auto expected_position = launch_position
+                                 + forward_vector({PitchAngle::from_degrees(-test_case.pitch),
+                                                   YawAngle::from_degrees(test_case.yaw), RollAngle::from_degrees(0.f)})
+                                           * projectile.m_launch_speed * test_case.time;
+        expected_position.z -= gravity * projectile.m_gravity_scale * test_case.time * test_case.time * 0.5f;
+
+        const auto actual_position = PredEngineTrait::predict_projectile_position(
+                projectile, test_case.pitch, test_case.yaw, test_case.time, gravity);
+
+        EXPECT_NEAR(actual_position.x, expected_position.x, 1e-5f);
+        EXPECT_NEAR(actual_position.y, expected_position.y, 1e-5f);
+        EXPECT_NEAR(actual_position.z, expected_position.z, 1e-5f);
+    }
+}
+
+TEST(PredEngineTrait, ReachabilityCheckUsesDistanceTolerance)
+{
+    constexpr Projectile projectile{
+            .m_origin = {0.f, 0.f, 0.f},
+            .m_launch_speed = 100.f,
+            .m_gravity_scale = 0.f,
+    };
+
+    EXPECT_FALSE(PredEngineTrait::can_projectile_reach_target_at_time(projectile, {100.f, 0.f, 0.f}, 0.5f, 9.81f, 0.f));
+    EXPECT_TRUE(PredEngineTrait::can_projectile_reach_target_at_time(projectile, {100.f, 0.f, 0.f}, 1.f, 9.81f, 0.f));
+    EXPECT_TRUE(PredEngineTrait::can_projectile_reach_target_at_time(projectile, {101.f, 0.f, 0.f}, 1.f, 9.81f, 1.f));
+    EXPECT_TRUE(PredEngineTrait::can_projectile_reach_target_at_time(projectile, {99.f, 0.f, 0.f}, 1.f, 9.81f, 1.f));
+    EXPECT_FALSE(PredEngineTrait::can_projectile_reach_target_at_time(projectile, {101.f, 0.f, 0.f}, 1.f, 9.81f, 0.f));
+    EXPECT_FALSE(PredEngineTrait::can_projectile_reach_target_at_time(projectile, {98.f, 0.f, 0.f}, 1.f, 9.81f, 1.f));
+}
+
+TEST(PredEngineTrait, ReachabilityCheckIncludesFloatingPointError)
+{
+    constexpr Projectile projectile{
+            .m_origin = {1.f, 2.f, 3.f},
+            .m_launch_offset = {0.1f, -0.2f, 0.3f},
+            .m_launch_speed = 100.f,
+            .m_gravity_scale = 1.f,
+    };
+    constexpr float gravity = 9.81f;
+    constexpr float time = 0.3f;
+    const auto target_position = PredEngineTrait::predict_projectile_position(projectile, 25.f, 45.f, time, gravity);
+
+    EXPECT_TRUE(PredEngineTrait::can_projectile_reach_target_at_time(projectile, target_position, time, gravity, 0.f));
+}
+
 TEST(PredEngineTrait, ZeroLaunchOffsetMatchesOriginalBehavior)
 {
     Projectile p;

tests/general/unit_test_proj_pred_engine_legacy_more.cpp

@@ -1,4 +1,5 @@
 #include <gtest/gtest.h>
+#include <limits>
 #include <omath/projectile_prediction/proj_pred_engine_legacy.hpp>
 #include <omath/projectile_prediction/projectile.hpp>
 #include <omath/projectile_prediction/target.hpp>
@@ -15,25 +16,45 @@ struct FakeEngineZeroGravity
     {
         return t.m_origin;
     }
-    static Vector3<float> predict_projectile_position(const Projectile& /*p*/, float /*pitch*/, float /*yaw*/, float /*time*/, float /*gravity*/) noexcept
+    static Vector3<float> predict_projectile_position(const Projectile& /*p*/, float /*pitch*/, float /*yaw*/,
+                                                      float /*time*/, float /*gravity*/) noexcept
     {
         // Return a fixed point matching typical target used in the test
         return Vector3<float>{100.f, 0.f, 0.f};
     }
-    static float calc_vector_2d_distance(const Vector3<float>& v) noexcept { return std::hypot(v.x, v.y); }
+    static float calc_vector_2d_distance(const Vector3<float>& v) noexcept
-    static float get_vector_height_coordinate(const Vector3<float>& v) noexcept { return v.z; }
+    {
-    static Vector3<float> calc_viewpoint_from_angles(const Projectile& /*p*/, Vector3<float> /*v*/, std::optional<float> /*maybe_pitch*/) noexcept
+        return std::hypot(v.x, v.y);
+    }
+    static float get_vector_height_coordinate(const Vector3<float>& v) noexcept
+    {
+        return v.z;
+    }
+    static Vector3<float> calc_viewpoint_from_angles(const Projectile& /*p*/, Vector3<float> /*v*/,
+                                                     std::optional<float> /*maybe_pitch*/) noexcept
     {
         return Vector3<float>{1.f, 2.f, 3.f};
     }
-    static float calc_direct_pitch_angle(const Vector3<float>& /*a*/, const Vector3<float>& /*b*/) noexcept { return 12.5f; }
+    static float calc_direct_pitch_angle(const Vector3<float>& /*a*/, const Vector3<float>& /*b*/) noexcept
-    static float calc_direct_yaw_angle(const Vector3<float>& /*a*/, const Vector3<float>& /*b*/) noexcept { return 0.f; }
+    {
+        return 12.5f;
+    }
+    static float calc_direct_yaw_angle(const Vector3<float>& /*a*/, const Vector3<float>& /*b*/) noexcept
+    {
+        return 0.f;
+    }
+    static bool can_projectile_reach_target_at_time(const Projectile& /*projectile*/,
+                                                    const Vector3<float>& /*target_position*/, float /*time*/,
+                                                    float /*gravity*/, float /*distance_tolerance*/)
+    {
+        return false;
+    }
 };
 
 TEST(ProjPredLegacyMore, ZeroGravityUsesDirectPitchAndReturnsViewpoint)
 {
-    constexpr Projectile proj{ .m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 10.f, .m_gravity_scale = 0.f };
+    constexpr Projectile proj{.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 10.f, .m_gravity_scale = 0.f};
-    constexpr Target target{ .m_origin = {100.f, 0.f, 0.f}, .m_velocity = {0.f,0.f,0.f}, .m_is_airborne = false };
+    constexpr Target target{.m_origin = {100.f, 0.f, 0.f}, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false};
 
     using Engine = omath::projectile_prediction::ProjPredEngineLegacy<FakeEngineZeroGravity>;
     const Engine engine(9.8f, 0.1f, 5.f, 1e-3f);
@@ -117,3 +138,77 @@ TEST(ProjPredLegacyMore, AngleComputedButMissReturnsNullopt)
     const auto res = engine.maybe_calculate_aim_point(proj, target);
     EXPECT_FALSE(res.has_value());
 }
+
+TEST(ProjPredLegacyMore, IncludesMaximumSimulationTime)
+{
+    constexpr Projectile projectile{.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 10.f, .m_gravity_scale = 0.f};
+    constexpr Target target{.m_origin = {10.f, 0.f, 0.f}, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false};
+    const omath::projectile_prediction::ProjPredEngineLegacy<> engine(0.f, 0.1f, 1.f, 0.f);
+
+    const auto result = engine.maybe_calculate_aim_point(projectile, target);
+
+    ASSERT_TRUE(result.has_value());
+    EXPECT_NEAR(result->x, target.m_origin.x, 1e-6f);
+    EXPECT_NEAR(result->y, target.m_origin.y, 1e-6f);
+    EXPECT_NEAR(result->z, target.m_origin.z, 1e-6f);
+}
+
+TEST(ProjPredLegacyMore, RejectsInvalidSimulationSteps)
+{
+    constexpr Projectile projectile{.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 10.f, .m_gravity_scale = 0.f};
+    constexpr Target target{.m_origin = {10.f, 0.f, 0.f}, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false};
+
+    EXPECT_FALSE(omath::projectile_prediction::ProjPredEngineLegacy<>(0.f, 0.f, 1.f, 1.f)
+                         .maybe_calculate_aim_point(projectile, target));
+    EXPECT_FALSE(omath::projectile_prediction::ProjPredEngineLegacy<>(0.f, -0.1f, 1.f, 1.f)
+                         .maybe_calculate_aim_point(projectile, target));
+    EXPECT_FALSE(
+            omath::projectile_prediction::ProjPredEngineLegacy<>(0.f, std::numeric_limits<float>::infinity(), 1.f, 1.f)
+                    .maybe_calculate_aim_point(projectile, target));
+}
+
+TEST(ProjPredLegacyMore, RejectsInvalidProjectileSpeedAndTolerance)
+{
+    constexpr Target target{.m_origin = {1.f, 0.f, 0.f}, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false};
+    constexpr Projectile stopped_projectile{.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 0.f, .m_gravity_scale = 0.f};
+    constexpr Projectile moving_projectile{.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 10.f, .m_gravity_scale = 0.f};
+
+    EXPECT_FALSE(omath::projectile_prediction::ProjPredEngineLegacy<>(0.f, 0.1f, 1.f, 2.f)
+                         .maybe_calculate_aim_point(stopped_projectile, target));
+    EXPECT_FALSE(omath::projectile_prediction::ProjPredEngineLegacy<>(0.f, 0.1f, 1.f, -1.f)
+                         .maybe_calculate_aim_point(moving_projectile, target));
+}
+
+TEST(ProjPredLegacyMore, StablePitchFindsHighSpeedLowArc)
+{
+    constexpr Projectile projectile{.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 10'000.f, .m_gravity_scale = 1.f};
+    constexpr Target target{.m_origin = {100.f, 0.f, 0.f}, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false};
+    const omath::projectile_prediction::ProjPredEngineLegacy<> engine(9.81f, 0.01f, 0.02f, 0.0001f);
+
+    const auto result = engine.maybe_calculate_aim_angles(projectile, target);
+
+    ASSERT_TRUE(result.has_value());
+    EXPECT_GT(result->pitch, 0.f);
+    EXPECT_NEAR(result->yaw, 0.f, 1e-6f);
+}
+
+TEST(ProjPredLegacyMore, CoincidentTargetReturnsLaunchOrigin)
+{
+    constexpr Projectile projectile{.m_origin = {5.f, 4.f, 3.f}, .m_launch_speed = 100.f, .m_gravity_scale = 1.f};
+    constexpr Target target{.m_origin = projectile.m_origin, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false};
+    const omath::projectile_prediction::ProjPredEngineLegacy<> engine(9.81f, 0.01f, 0.01f, 0.f);
+
+    const auto result = engine.maybe_calculate_aim_point(projectile, target);
+
+    ASSERT_TRUE(result.has_value());
+    EXPECT_EQ(result.value(), projectile.m_origin);
+}
+
+TEST(ProjPredLegacyMore, TinyDistanceDoesNotUnderflowToAHit)
+{
+    constexpr Projectile projectile{.m_origin = {0.f, 0.f, 0.f}, .m_launch_speed = 1.f, .m_gravity_scale = 0.f};
+    constexpr Target target{.m_origin = {1e-30f, 0.f, 0.f}, .m_velocity = {0.f, 0.f, 0.f}, .m_is_airborne = false};
+    const omath::projectile_prediction::ProjPredEngineLegacy<> engine(0.f, 0.1f, 0.f, 0.f);
+
+    EXPECT_FALSE(engine.maybe_calculate_aim_point(projectile, target));
+}

vcpkg-configuration.json

@@ -1,7 +1,7 @@
 {
   "default-registry": {
     "kind": "git",
-    "baseline": "b1b19307e2d2ec1eefbdb7ea069de7d4bcd31f01",
+    "baseline": "0878b5224d4a4968940ee296a2e7fae2d3b62983",
     "repository": "https://github.com/microsoft/vcpkg"
   },
   "registries": [

vcpkg.json

@@ -21,7 +21,11 @@
       "dependencies": [
         {
           "name": "omath",
-          "features": ["imgui", "lua", "hooking"]
+          "features": [
+            "imgui",
+            "lua",
+            "hooking"
+          ]
         }
       ]
     },
@@ -50,16 +54,26 @@
         "opengl",
         {
           "name": "omath",
-          "features": ["hooking"],
+          "features": [
+            "hooking"
+          ],
           "platform": "windows & !arm & !uwp"
         },
         {
           "name": "imgui",
-          "features": ["glfw-binding", "opengl3-binding"]
+          "features": [
+            "glfw-binding",
+            "opengl3-binding"
+          ]
         },
         {
           "name": "imgui",
-          "features": ["dx9-binding", "dx11-binding", "dx12-binding", "win32-binding"],
+          "features": [
+            "dx9-binding",
+            "dx11-binding",
+            "dx12-binding",
+            "win32-binding"
+          ],
           "platform": "windows & !arm & !uwp"
         }
       ]