fix(math): correct select logic in quaternion and svd

Fixed conditional selection logic in quaternion and SVD math functions by swapping select argument order for correctness. Fixed LookRotationSafe and normalizesafe to return valid quaternions. Corrected SVD helper functions for proper value swapping and safe reciprocal. Added unit tests for matrix, reflection, projection, refraction, quaternion normalization, LookRotationSafe, and SVD operations. Incremented project version to 1.3.3. Minor formatting and using directive updates.
This commit is contained in:
2026-04-07 22:18:55 +09:00
parent 7ea7e8aa9e
commit b08662b77d
7 changed files with 191 additions and 15 deletions

View File

@@ -7,7 +7,7 @@
<AllowUnsafeBlocks>True</AllowUnsafeBlocks> <AllowUnsafeBlocks>True</AllowUnsafeBlocks>
<GeneratePackageOnBuild>true</GeneratePackageOnBuild> <GeneratePackageOnBuild>true</GeneratePackageOnBuild>
<Authors>Misaki</Authors> <Authors>Misaki</Authors>
<AssemblyVersion>1.3.2</AssemblyVersion> <AssemblyVersion>1.3.3</AssemblyVersion>
<Version>$(AssemblyVersion)</Version> <Version>$(AssemblyVersion)</Version>
<PackageProjectUrl>https://git.personalnas.com/Misaki/Misaki.HighPerformance.git</PackageProjectUrl> <PackageProjectUrl>https://git.personalnas.com/Misaki/Misaki.HighPerformance.git</PackageProjectUrl>
<RepositoryUrl>https://git.personalnas.com/Misaki/Misaki.HighPerformance.git</RepositoryUrl> <RepositoryUrl>https://git.personalnas.com/Misaki/Misaki.HighPerformance.git</RepositoryUrl>

View File

@@ -60,7 +60,7 @@ public partial struct quaternion : IEquatable<quaternion>
var u_mask = new uint4((int)u_sign >> 31); var u_mask = new uint4((int)u_sign >> 31);
var t_mask = new uint4(asint(t) >> 31); var t_mask = new uint4(asint(t) >> 31);
float tr = 1.0f + abs(u.x); var tr = 1.0f + abs(u.x);
var sign_flips = new uint4(0x00000000, 0x80000000, 0x80000000, 0x80000000) ^ (u_mask & new uint4(0x00000000, 0x80000000, 0x00000000, 0x80000000)) ^ (t_mask & new uint4(0x80000000, 0x80000000, 0x80000000, 0x00000000)); var sign_flips = new uint4(0x00000000, 0x80000000, 0x80000000, 0x80000000) ^ (u_mask & new uint4(0x00000000, 0x80000000, 0x00000000, 0x80000000)) ^ (t_mask & new uint4(0x80000000, 0x80000000, 0x80000000, 0x00000000));
@@ -428,7 +428,7 @@ public partial struct quaternion : IEquatable<quaternion>
var mx = max(max(forwardLengthSq, upLengthSq), tLengthSq); var mx = max(max(forwardLengthSq, upLengthSq), tLengthSq);
var accept = mn > 1e-35f && mx < 1e35f && isfinite(forwardLengthSq) && isfinite(upLengthSq) && isfinite(tLengthSq); var accept = mn > 1e-35f && mx < 1e35f && isfinite(forwardLengthSq) && isfinite(upLengthSq) && isfinite(tLengthSq);
return new quaternion(select(float4.unitW, new quaternion(new float3x3(t, cross(forward, t), forward)).value, accept)); return new quaternion(select(new quaternion(new float3x3(t, cross(forward, t), forward)).value, float4.unitW, accept));
} }
/// <summary>Returns true if the quaternion is equal to a given quaternion, false otherwise.</summary> /// <summary>Returns true if the quaternion is equal to a given quaternion, false otherwise.</summary>
@@ -559,7 +559,7 @@ public static partial class math
{ {
var x = q.value; var x = q.value;
var len = dot(x, x); var len = dot(x, x);
return quaternion(select(Mathematics.quaternion.identity.value, x * rsqrt(len), len > FLT_MIN_NORMAL)); return quaternion(select(x * rsqrt(len), Mathematics.quaternion.identity.value, len > FLT_MIN_NORMAL));
} }
/// <summary> /// <summary>
@@ -574,7 +574,7 @@ public static partial class math
{ {
var x = q.value; var x = q.value;
var len = dot(x, x); var len = dot(x, x);
return quaternion(select(defaultvalue.value, x * rsqrt(len), len > FLT_MIN_NORMAL)); return quaternion(select(x * rsqrt(len), defaultvalue.value, len > FLT_MIN_NORMAL));
} }
/// <summary>Returns the natural exponent of a quaternion. Assumes w is zero.</summary> /// <summary>Returns the natural exponent of a quaternion. Assumes w is zero.</summary>
@@ -757,7 +757,7 @@ public static partial class math
{ {
i = adj(m, out var det); i = adj(m, out var det);
var c = abs(det) > epsilon; var c = abs(det) > epsilon;
var detInv = select(new float3(1f), rcp(det), c); var detInv = select(rcp(det), new float3(1f), c);
i = scaleMul(detInv, i); i = scaleMul(detInv, i);
return c; return c;
} }

View File

@@ -1,4 +1,4 @@
using System.Runtime.CompilerServices; using System.Runtime.CompilerServices;
namespace Misaki.HighPerformance.Mathematics; namespace Misaki.HighPerformance.Mathematics;
@@ -18,23 +18,23 @@ public static class svd
private static void condSwap(bool c, ref float x, ref float y) private static void condSwap(bool c, ref float x, ref float y)
{ {
var tmp = x; var tmp = x;
x = math.select(x, y, c); x = math.select(y, x, c);
y = math.select(y, tmp, c); y = math.select(tmp, y, c);
} }
[MethodImpl(MethodImplOptions.AggressiveInlining)] [MethodImpl(MethodImplOptions.AggressiveInlining)]
private static void condNegSwap(bool c, ref float3 x, ref float3 y) private static void condNegSwap(bool c, ref float3 x, ref float3 y)
{ {
var tmp = -x; var tmp = -x;
x = math.select(x, y, c); x = math.select(y, x, c);
y = math.select(y, tmp, c); y = math.select(tmp, y, c);
} }
[MethodImpl(MethodImplOptions.AggressiveInlining)] [MethodImpl(MethodImplOptions.AggressiveInlining)]
private static quaternion condNegSwapQuat(bool c, quaternion q, float4 mask) private static quaternion condNegSwapQuat(bool c, quaternion q, float4 mask)
{ {
const float halfSqrt2 = 0.707106781186548f; const float halfSqrt2 = 0.707106781186548f;
return math.mul(q, math.select(quaternion.identity.value, mask * halfSqrt2, c)); return math.mul(q, math.select(mask * halfSqrt2, quaternion.identity.value, c));
} }
[MethodImpl(MethodImplOptions.AggressiveInlining)] [MethodImpl(MethodImplOptions.AggressiveInlining)]
@@ -68,7 +68,7 @@ public static class svd
var ch = 2f * (pq.x - pq.y); // approx cos(a/2) var ch = 2f * (pq.x - pq.y); // approx cos(a/2)
var sh = pq.z; // approx sin(a/2) var sh = pq.z; // approx sin(a/2)
var r = math.select(new float4(s8, s8, s8, c8), new float4(sh, sh, sh, ch), g * sh * sh < ch * ch) * mask; var r = math.select(new float4(sh, sh, sh, ch), new float4(s8, s8, s8, c8), g * sh * sh < ch * ch) * mask;
return math.normalize(r); return math.normalize(r);
} }
@@ -76,7 +76,7 @@ public static class svd
private static quaternion qrGivensQuat(float2 pq, float4 mask) private static quaternion qrGivensQuat(float2 pq, float4 mask)
{ {
var l = math.sqrt(pq.x * pq.x + pq.y * pq.y); var l = math.sqrt(pq.x * pq.x + pq.y * pq.y);
var sh = math.select(0f, pq.y, l > EPSILON_NORMAL_SQRT); var sh = math.select(pq.y, 0f, l > EPSILON_NORMAL_SQRT);
var ch = math.abs(pq.x) + math.max(l, EPSILON_NORMAL_SQRT); var ch = math.abs(pq.x) + math.max(l, EPSILON_NORMAL_SQRT);
condSwap(pq.x < 0f, ref sh, ref ch); condSwap(pq.x < 0f, ref sh, ref ch);
@@ -148,7 +148,7 @@ public static class svd
[MethodImpl(MethodImplOptions.AggressiveInlining)] [MethodImpl(MethodImplOptions.AggressiveInlining)]
private static float3 rcpsafe(float3 x, float epsilon = EPSILON_RCP) => private static float3 rcpsafe(float3 x, float epsilon = EPSILON_RCP) =>
math.select(math.rcp(x), float3.zero, math.abs(x) < epsilon); math.select(float3.zero, math.rcp(x), math.abs(x) < epsilon);
[MethodImpl(MethodImplOptions.AggressiveInlining)] [MethodImpl(MethodImplOptions.AggressiveInlining)]
public static float3x3 svdInverse(float3x3 a) public static float3x3 svdInverse(float3x3 a)

View File

@@ -6,7 +6,9 @@ using Misaki.HighPerformance.Test.Benchmark;
var faceDirection = math.normalize(float3.zero - new float3(0.0f, 0.0f, 5.0f)); var faceDirection = math.normalize(float3.zero - new float3(0.0f, 0.0f, 5.0f));
var test = quaternion.LookRotation(faceDirection, math.up()); var test = quaternion.LookRotation(faceDirection, math.up());
var test2 = quaternion.LookRotationSafe(faceDirection, math.up());
var rotation = quaternion.EulerXYZ(new float3(0, math.radians(180.0f), 0)); var rotation = quaternion.EulerXYZ(new float3(0, math.radians(180.0f), 0));
Console.WriteLine(test); Console.WriteLine(test);
Console.WriteLine(test2);
Console.WriteLine(rotation); Console.WriteLine(rotation);

View File

@@ -1,4 +1,5 @@
using Misaki.HighPerformance.Mathematics; using Misaki.HighPerformance.Mathematics;
using System.Numerics;
namespace Misaki.HighPerformance.Test.UnitTest.Mathematics; namespace Misaki.HighPerformance.Test.UnitTest.Mathematics;
@@ -288,4 +289,64 @@ public class TestMathFunctions
Assert.AreEqual(5f, result.y, 1e-6f); // condition is false, so select b.y Assert.AreEqual(5f, result.y, 1e-6f); // condition is false, so select b.y
Assert.AreEqual(3f, result.z, 1e-6f); // condition is true, so select a.z Assert.AreEqual(3f, result.z, 1e-6f); // condition is true, so select a.z
} }
[TestMethod]
public void TestMatrixFunctions()
{
// Test determinant and inverse of a 2x2 matrix
var m = new float2x2(new float2(1f, 2f), new float2(3f, 4f));
var det = math.determinant(m);
Assert.AreEqual(-2f, det, 1e-6f);
var inv = math.inverse(m);
var expectedInv = new float2x2(new float2(-2f, 1f), new float2(1.5f, -0.5f));
Assert.AreEqual(expectedInv.c0.x, inv.c0.x, 1e-6f);
Assert.AreEqual(expectedInv.c0.y, inv.c0.y, 1e-6f);
Assert.AreEqual(expectedInv.c1.x, inv.c1.x, 1e-6f);
Assert.AreEqual(expectedInv.c1.y, inv.c1.y, 1e-6f);
}
[TestMethod]
public void TestRflections()
{
// Test reflect function
var incident = new float3(1f, -1f, 0f);
var normal = math.normalize(new float3(0f, 1f, 0f));
var reflected = math.reflect(incident, normal);
// The reflected vector should be (1, 1, 0)
Assert.AreEqual(1f, reflected.x, 1e-6f);
Assert.AreEqual(1f, reflected.y, 1e-6f);
Assert.AreEqual(0f, reflected.z, 1e-6f);
}
[TestMethod]
public void TestProjections()
{
// Test project function
var vector = new float3(1f, 2f, 3f);
var onNormal = math.normalize(new float3(0f, 1f, 0f));
var projected = math.project(vector, onNormal);
// The projection of (1,2,3) onto the Y-axis should be (0,2,0)
Assert.AreEqual(0f, projected.x, 1e-6f);
Assert.AreEqual(2f, projected.y, 1e-6f);
Assert.AreEqual(0f, projected.z, 1e-6f);
}
[TestMethod]
public void TestRefraction()
{
// Test refract function
var incident = math.normalize(new float3(1f, -1f, 0f));
var normal = math.normalize(new float3(0f, 1f, 0f));
var eta = 0.5f; // Refractive index ratio
var refracted = math.refract(incident, normal, eta);
// The refracted vector should be approximately (0.707, -0.707, 0)
Assert.AreEqual(0.3535534f, refracted.x, 1e-6f);
Assert.AreEqual(-0.9354144f, refracted.y, 1e-6f);
Assert.AreEqual(0f, refracted.z, 1e-6f);
}
} }

View File

@@ -151,6 +151,18 @@ public class TestQuaternion
Assert.AreEqual(1f, math.length(normalized.value), 1e-6f); Assert.AreEqual(1f, math.length(normalized.value), 1e-6f);
} }
[TestMethod]
public void TestQuaternionNormalizeSafe()
{
// Create non-unit quaternion
var q = new quaternion(1f, 2f, 3f, 4f);
var normalized = math.normalizesafe(q);
// Should have length 1
Assert.AreEqual(1f, math.length(normalized.value), 1e-6f);
}
[TestMethod] [TestMethod]
public void TestMatrixFromQuaternion() public void TestMatrixFromQuaternion()
{ {
@@ -240,4 +252,27 @@ public class TestQuaternion
Assert.AreEqual(1f, rotatedUp.y, 1e-6f); Assert.AreEqual(1f, rotatedUp.y, 1e-6f);
Assert.AreEqual(0f, rotatedUp.z, 1e-6f); Assert.AreEqual(0f, rotatedUp.z, 1e-6f);
} }
[TestMethod]
public void TestLookRotationSafe()
{
var forward = new float3(0f, 0f, -1f);
var up = new float3(0f, 1f, 0f);
var q = quaternion.LookRotationSafe(forward, up);
// Should be unit quaternion
Assert.AreEqual(1f, math.length(q.value), 1e-6f);
// Should rotate (0,0,1) to approximately (0,0,1)
var rotatedForward = math.mul(q, forward);
Assert.AreEqual(0f, rotatedForward.x, 1e-6f);
Assert.AreEqual(0f, rotatedForward.y, 1e-6f);
Assert.AreEqual(1f, rotatedForward.z, 1e-6f);
// Should rotate (0,1,0) to approximately (0,1,0)
var rotatedUp = math.mul(q, up);
Assert.AreEqual(0f, rotatedUp.x, 1e-6f);
Assert.AreEqual(1f, rotatedUp.y, 1e-6f);
Assert.AreEqual(0f, rotatedUp.z, 1e-6f);
}
} }

View File

@@ -0,0 +1,78 @@
using System;
using System.Collections.Generic;
using System.Text;
using Microsoft.VisualStudio.TestTools.UnitTesting;
using Misaki.HighPerformance.Mathematics;
namespace Misaki.HighPerformance.Test.UnitTest.Mathematics;
[TestClass]
public class TestSVD
{
private const float TOLERANCE = 1e-5f;
private static void AssertFloat3x3Equal(float3x3 expected, float3x3 actual, float tolerance = TOLERANCE)
{
Assert.AreEqual(expected.c0.x, actual.c0.x, tolerance, "c0.x mismatch");
Assert.AreEqual(expected.c0.y, actual.c0.y, tolerance, "c0.y mismatch");
Assert.AreEqual(expected.c0.z, actual.c0.z, tolerance, "c0.z mismatch");
Assert.AreEqual(expected.c1.x, actual.c1.x, tolerance, "c1.x mismatch");
Assert.AreEqual(expected.c1.y, actual.c1.y, tolerance, "c1.y mismatch");
Assert.AreEqual(expected.c1.z, actual.c1.z, tolerance, "c1.z mismatch");
Assert.AreEqual(expected.c2.x, actual.c2.x, tolerance, "c2.x mismatch");
Assert.AreEqual(expected.c2.y, actual.c2.y, tolerance, "c2.y mismatch");
Assert.AreEqual(expected.c2.z, actual.c2.z, tolerance, "c2.z mismatch");
}
private static float3x3 CreateDiagonal(float3 s)
{
return new float3x3(new float3(s.x, 0, 0), new float3(0, s.y, 0), new float3(0, 0, s.z));
}
[TestMethod]
public void TestSVDInverse_Identity()
{
var identity = float3x3.identity;
var inverse = svd.svdInverse(identity);
AssertFloat3x3Equal(identity, inverse);
}
[TestMethod]
public void TestSVDInverse_Diagonal()
{
var scale = new float3(2f, 0.5f, 4f);
var a = CreateDiagonal(scale);
var expected = CreateDiagonal(1f / scale);
var actual = svd.svdInverse(a);
AssertFloat3x3Equal(expected, actual);
}
[TestMethod]
public void TestSVDRotation_PureRotation()
{
var q = quaternion.AxisAngle(math.normalize(new float3(1f, 2f, 3f)), 1.2f);
var a = new float3x3(q);
var result = svd.svdRotation(a);
// SVD rotation should extract the rotation part.
// For a pure rotation matrix, result should be the same as input q (or -q)
var dot = math.dot(q.value, result.value);
Assert.IsTrue(math.abs(dot) > 0.999f, $"Rotation mismatch: dot product {dot}");
}
[TestMethod]
public void TestSVDInverse_RotationAndScale()
{
var q = quaternion.AxisAngle(math.normalize(new float3(0.5f, -0.2f, 0.8f)), 0.5f);
var r = new float3x3(q);
var scale = new float3(1.5f, 0.8f, 2.0f);
// A = R * S
var a = math.mulScale(r, scale);
var inverse = svd.svdInverse(a);
// Check inverse * A = Identity
var identityResult = math.mul(inverse, a);
AssertFloat3x3Equal(float3x3.identity, identityResult, 1e-4f);
}
}