Skip to content

Commit 9716faf

Browse files
committed
fixed, improved, fast math methods + added benchmarks
1 parent 183a76c commit 9716faf

8 files changed

Lines changed: 885 additions & 123 deletions

Runtime/FastMath/FastCosine.cs

Lines changed: 0 additions & 62 deletions
This file was deleted.

Runtime/FastMath/FastCosine.cs.meta

Lines changed: 0 additions & 11 deletions
This file was deleted.

Runtime/FastMath/FastFunctions.cs

Lines changed: 112 additions & 36 deletions
Original file line numberDiff line numberDiff line change
@@ -4,9 +4,11 @@
44
// ** Repository : https://github.com/LTMX/Unity.mathx
55
#endregion
66

7+
using System;
78
using System.Runtime.CompilerServices;
89
using System.Runtime.InteropServices;
9-
using UnityEngine;
10+
using Unity.Burst;
11+
using Unity.Burst.Intrinsics;
1012

1113
namespace Unity.Mathematics
1214
{
@@ -44,6 +46,108 @@ private struct FloatIntUnion
4446
[MethodImpl(IL)] public static float2 fsqrt(this float2 f) => new(f.x.fsqrt(), f.y.fsqrt()); // to never simplify to new f2(f.xy.fastsqrt())
4547

4648

49+
[StructLayout(LayoutKind.Explicit)]
50+
private struct FloatV128
51+
{
52+
[FieldOffset(0)] public float f;
53+
[FieldOffset(0)] public v128 v;
54+
}
55+
56+
/// Fast reciprocal (1 / f).
57+
///
58+
/// Primary path : RCPSS hardware intrinsic — ~0.037% max relative error
59+
/// Fallback path : IEEE 754 exponent bit-hack — ~3.4% max relative error
60+
///
61+
/// Both are significantly faster than scalar division.
62+
/// Sign, zero, infinity and NaN are preserved on the intrinsic path.
63+
/// The bit-hack path handles normal floats only (see caveats below).
64+
[MethodImpl(MethodImplOptions.AggressiveInlining)]
65+
public static unsafe float fastrcp(this float f)
66+
{
67+
// if (X86.Sse.IsSseSupported)
68+
// {
69+
v128 result = X86.Sse.rcp_ss(*(v128*)&f);
70+
return *(float*)&result;
71+
// }
72+
//
73+
// return math.asfloat(0x7EF127EA - math.asint(f));
74+
}
75+
76+
[StructLayout(LayoutKind.Explicit)]
77+
private struct Float2V128
78+
{
79+
[FieldOffset(0)] public float2 f2;
80+
[FieldOffset(0)] public v128 v;
81+
}
82+
83+
[StructLayout(LayoutKind.Explicit)]
84+
private struct Float3V128
85+
{
86+
[FieldOffset(0)] public float3 f3;
87+
88+
[FieldOffset(0)] public v128 v;
89+
// float3 = 12 bytes, v128 = 16 bytes
90+
// lane 3 is padding garbage — rcp_ps will compute it anyway, we just ignore it
91+
}
92+
93+
[StructLayout(LayoutKind.Explicit)]
94+
private struct Float4V128
95+
{
96+
[FieldOffset(0)] public float4 f4;
97+
98+
[FieldOffset(0)] public v128 v;
99+
// float4 = 16 bytes = v128 exactly, perfect fit
100+
}
101+
102+
const int Magic = 0x7EF127EA;
103+
private static readonly v128 magic = new(0x7EF127EA, 0x7EF127EA, 0x7EF127EA, 0x7EF127EA);
104+
105+
[MethodImpl(MethodImplOptions.AggressiveInlining)]
106+
public static unsafe float2 fastrcp(this float2 f)
107+
{
108+
// broadcast the magic constant into all 4 lanes
109+
110+
// integer subtract on all lanes at once — no lookup table, pure ALU
111+
v128 result = X86.Sse2.sub_epi32(magic, *(v128*)&f);
112+
return *(float2*)&result;
113+
}
114+
115+
[MethodImpl(MethodImplOptions.AggressiveInlining)]
116+
public static unsafe float3 fastrcp(this float3 f)
117+
{
118+
// if (X86.Sse.IsSseSupported)
119+
// {
120+
v128 result = X86.Sse2.sub_epi32(magic, *(v128*)&f);
121+
return *(float3*)&result;
122+
// }
123+
//
124+
// return new float3(math.asfloat(Magic - math.asint(f.x)), math.asfloat(Magic - math.asint(f.y)), math.asfloat(Magic - math.asint(f.z)));
125+
}
126+
127+
[MethodImpl(MethodImplOptions.AggressiveInlining)]
128+
public static unsafe float4 fastrcp(this float4 f)
129+
{
130+
// if (X86.Sse.IsSseSupported)
131+
// {
132+
v128 result = X86.Sse2.sub_epi32(magic, *(v128*)&f);
133+
return *(float4*)&result;
134+
// }
135+
//
136+
// return new float4(math.asfloat(Magic - math.asint(f.x)), math.asfloat(Magic - math.asint(f.y)), math.asfloat(Magic - math.asint(f.z)),
137+
// math.asfloat(Magic - math.asint(f.w)));
138+
}
139+
140+
141+
// /// <inheritdoc cref="fastrcp(float)"/>
142+
// [MethodImpl(IL)] public static float fastrcp(this int x) => ((float)x).fastrcp();
143+
// /// <inheritdoc cref="fastrcp(float)"/>
144+
// [MethodImpl(IL)] public static float4 fastrcp(this float4 f) => new(f.x.fastrcp(), f.y.fastrcp(), f.z.fastrcp(), f.w.fastrcp());
145+
// /// <inheritdoc cref="fastrcp(float)"/>
146+
// [MethodImpl(IL)] public static float3 fastrcp(this float3 f) => new(f.x.fastrcp(), f.y.fastrcp(), f.z.fastrcp());
147+
// /// <inheritdoc cref="fastrcp(float)"/>
148+
// [MethodImpl(IL)] public static float2 fastrcp(this float2 f) => new(f.x.fastrcp(), f.y.fastrcp());
149+
150+
47151
/// Returns the distance between a and b (fast but low accuracy)
48152
[MethodImpl(IL)] public static float fdistance(float4 a, float4 b) => (a - b).flengthsq().fsqrt();
49153
/// <inheritdoc cref="fdistance(float4, float4)"/>
@@ -66,6 +170,13 @@ private struct FloatIntUnion
66170
/// <inheritdoc cref="math.lengthsq(float2)"/>
67171
[MethodImpl(IL)] public static float flengthsq(this float2 f) => f.fdot(f);
68172

173+
/// <inheritdoc cref="math.distancesq(float4, float4)"/>
174+
[MethodImpl(IL)] public static float fdistancesq(this float4 f, float4 f2) => flengthsq(f2 - f);
175+
/// <inheritdoc cref="math.distancesq(float4, float4)"/>
176+
[MethodImpl(IL)] public static float fdistancesq(this float3 f, float3 f2) => flengthsq(f2 - f);
177+
/// <inheritdoc cref="math.distancesq(float4, float4)"/>
178+
[MethodImpl(IL)] public static float fdistancesq(this float2 f, float2 f2) => flengthsq(f2 - f);
179+
69180
/// faster dot method removing to double casts
70181
[MethodImpl(IL)] public static float fdot(this float4 f, float4 f2) => f.x * f2.x + f.y * f2.y + f.z * f2.z + f.w * f2.w;
71182
/// <inheritdoc cref="fdot(float4,float4)"/>
@@ -103,40 +214,5 @@ public static int log2int(this int value)
103214
[MethodImpl(IL)] public static float3 fexp(float3 f) => new(fexp(f.x), fexp(f.y), fexp(f.z));
104215
/// <inheritdoc cref="fexp(float)"/>
105216
[MethodImpl(IL)] public static float4 fexp(float4 f) => new(fexp(f.x), fexp(f.y), fexp(f.z), fexp(f.w));
106-
107-
108-
#region Deprecated
109-
110-
// [BurstCompile]
111-
// [StructLayout(LayoutKind.Explicit)]
112-
// private struct FloatUInt32Union
113-
// {
114-
// [FieldOffset(0)] public float f;
115-
// [FieldOffset(0)] public uint u;
116-
// }
117-
//
118-
//
119-
// /// returns 1/x using fast math
120-
// [MethodImpl(IL)]
121-
// public static float frcp(this float x)
122-
// {
123-
// FloatUInt32Union fiu = new();
124-
// fiu.f = x;
125-
// fiu.u = (0xbe6eb3beU - fiu.u) >> 1; // pow( x, -0.5 )
126-
// return fiu.f * fiu.f; // pow( pow(x,-0.5), 2 ) = pow( x, -1 ) = 1.0 / x
127-
// }
128-
//
129-
//
130-
// /// returns 1/x using fast math
131-
// [MethodImpl(IL)]
132-
// public static float frcp(this int x)
133-
// {
134-
// FloatUInt32Union fiu = new();
135-
// fiu.f = x;
136-
// fiu.u = (0xbe6eb3beU - fiu.u) >> 1; // pow( x, -0.5 )
137-
// return fiu.f * fiu.f; // pow( pow(x,-0.5), 2 ) = pow( x, -1 ) = 1.0 / x
138-
// }
139-
140-
#endregion
141217
}
142218
}

Runtime/FastMath/FastTrigonometry.cs

Lines changed: 18 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -4,6 +4,8 @@
44
// ** Repository : https://github.com/LTMX/Unity.mathx
55
#endregion
66

7+
using System.Runtime.CompilerServices;
8+
79
namespace Unity.Mathematics
810
{
911
public static partial class mathx
@@ -12,7 +14,7 @@ public static partial class mathx
1214
private const float t2 = 0.405284735f;
1315

1416
/// Low precision sine (~14x faster) - always wrap input angle to -PI..PI
15-
public static float sfastsine(this float x){
17+
[MethodImpl(IL)] public static float sfastsine(this float x){
1618

1719
if (x < -PI) x += TAU;
1820
else if (x > PI) x -= TAU;
@@ -21,16 +23,15 @@ public static float sfastsine(this float x){
2123

2224
//sin(x + PI/2) = cos(x)
2325
/// Low precision cosine (~14x faster)
24-
public static float sfastcosine(this float x)
26+
[MethodImpl(IL)] public static float sfastcosine(this float x)
2527
{
2628
x += HPI;
2729
if (x > PI) x -= TAU;
2830
return x < 0 ? x * (t1+ t2* x) : x * (t1- t2* x);
2931
}
3032

31-
//
3233
/// High precision sine (~8x faster) - always wrap input angle to -PI..PI
33-
public static float fastsine(this float x)
34+
[MethodImpl(IL)] public static float fastsine(this float x)
3435
{
3536
if (x < -PI) x += TAU;
3637
else if (x > PI) x -= TAU;
@@ -44,17 +45,20 @@ public static float fastsine(this float x)
4445

4546
//sin(x + PI/2) = cos(x)
4647
/// High precision cosine (~8x faster) - always wrap input angle to -PI..PI
47-
public static float fastcosine(this float x) => fastsine(x + HPI);
48-
49-
// Overloads
50-
51-
//sin(x + PI/2) = cos(x)
52-
/// Low precision cosine (~14x faster)
53-
public static float sfastcosine(this int f)
48+
[MethodImpl(IL)] public static float fastcosine(this float x)
5449
{
55-
var x = f + HPI;
56-
if (x > PI) x -= TAU;
57-
return x < 0 ? x * (t1+ t2* x) : x * (t1- t2* x);
50+
x += HPI;
51+
if (x < -PI) x += TAU;
52+
else if (x > PI) x -= TAU;
53+
if (x < 0) {
54+
var s = x * (t1 + t2 * x);
55+
return s < 0 ? .225f * (s * -s - s) + s : .225f * (s * s - s) + s;
56+
}
57+
var s2 = x * (t1 - t2 * x);
58+
return s2 < 0 ? .225f * (s2 * -s2 - s2) + s2 : .225f * (s2 * s2 - s2) + s2;
5859
}
60+
/// Low precision cosine (~14x faster)
61+
[MethodImpl(IL)] public static float sfastcosine(this int f) => ((float)f).sfastcosine();
62+
5963
}
6064
}

0 commit comments

Comments
 (0)