< Summary

Information
Class: FixedMathSharp.Geometry.BoxCylinderAnalyticFeatures
Assembly: FixedMathSharp
File(s): /home/runner/work/FixedMathSharp/FixedMathSharp/src/FixedMathSharp/Geometry/Wide/OrientedBox/BoxCylinderAnalyticFeatures.cs
Line coverage
100%
Covered lines: 197
Uncovered lines: 0
Coverable lines: 197
Total lines: 328
Line coverage: 100%
Branch coverage
100%
Covered branches: 86
Total branches: 86
Branch coverage: 100%
Method coverage

Feature is only available for sponsors

Upgrade to PRO version

Metrics

MethodBranch coverage Crap Score Cyclomatic complexity Line coverage
TryGetBest(...)100%3434100%
IsEdgeConeDominated(...)100%1010100%
KeepAxis(...)100%1414100%
BuildRankingAxis(...)100%1616100%
RestoreRawGapDenominator(...)100%11100%
HasVertexRimSeparation(...)100%1212100%
WritePlaneDirection(...)100%11100%
Unit(...)100%11100%
.ctor(...)100%11100%

File(s)

/home/runner/work/FixedMathSharp/FixedMathSharp/src/FixedMathSharp/Geometry/Wide/OrientedBox/BoxCylinderAnalyticFeatures.cs

#LineLine coverage
 1//=======================================================================
 2// BoxCylinderAnalyticFeatures.cs
 3//=======================================================================
 4// MIT License, Copyright (c) 2024-present David Oravsky (mrdav30)
 5// See LICENSE file in the project root for full license information.
 6//=======================================================================
 7
 8using System;
 9using static FixedMathSharp.Geometry.CylinderContactAlgebra;
 10
 11namespace FixedMathSharp.Geometry;
 12
 13/// <summary>
 14/// Exact rational support directions and vertex/rim separation certificates
 15/// for a box and a finite cylinder. Coordinates stay in the authored box frame.
 16/// </summary>
 17/// <remarks>
 18/// The support sphere splits at box-edge planes, the cylinder's side plane
 19/// and its two poles. This owner handles their intersections and side arcs;
 20/// BoxCylinderEdgeContacts handles the remaining edge/rim arcs. Inside an
 21/// open vertex/rim patch, zero meridional curvature rules out a positive
 22/// minimum. Its outward residual can still separate, so that certificate
 23/// must not be omitted merely because it cannot supply penetration depth.
 24/// </remarks>
 25internal static class BoxCylinderAnalyticFeatures
 26{
 27    private const int Words = ConvexContactCandidate.Words;
 28
 29    internal static bool TryGetBest(in BoxCylinderGeometry geometry,
 30        Span<ulong> bestValues, Span<int> bestSigns, out int bestGapSign,
 31        out Vector3d localSupportSigns)
 32    {
 118333        var selection = new Selection(bestValues, bestSigns);
 118334        Span<ulong> direction = stackalloc ulong[3 * Words];
 118335        Span<int> directionSigns = stackalloc int[3];
 118336        WriteDirection(geometry.Axis, direction, directionSigns);
 118337        KeepAxis(geometry, direction, directionSigns, ref selection);
 944238        for (int axis = 0; axis < 3 && selection.GapSign >= 0; axis++)
 39        {
 353840            WideAxis3 unit = Unit(axis);
 353841            WriteDirection(unit, direction, directionSigns);
 353842            KeepAxis(geometry, direction, directionSigns, ref selection);
 353843            WriteDirection(WideAxis3.Cross(unit, geometry.Axis), direction, directionSigns);
 353844            KeepAxis(geometry, direction, directionSigns, ref selection);
 45        }
 46
 118347        if (!geometry.Radius.IsZero)
 48        {
 49            // On a cylinder's side-normal plane, each box vertex contributes
 50            // a linear term plus a circle. Its radial direction and the plane's
 51            // box-edge boundaries exhaust the possible minima.
 2106252            for (int corner = 0; corner < 8 && selection.GapSign >= 0; corner++)
 53            {
 936154                WideAxis3 offset = geometry.GetCapVertexOffset(corner, 0);
 936155                ProjectPerpendicular(geometry.Axis, offset, direction, directionSigns);
 936156                KeepAxis(geometry, direction, directionSigns, ref selection);
 5615257                for (int cap = -1; cap <= 1; cap += 2)
 58                {
 1872059                    if (HasVertexRimSeparation(geometry, corner, cap))
 60                    {
 561                        RestoreRawGapDenominator(geometry, bestValues, direction);
 562                        bestGapSign = -1;
 563                        localSupportSigns = default;
 564                        return false;
 65                    }
 66                }
 67            }
 68
 69            // The edge quartic eliminates sqrt(Q) using K=q-p*t. K=0 is a
 70            // rational direction, not permission to divide by zero. Include
 71            // it here using the full support function, independent of a cone.
 72            // When p=q=0 the projected disk's two principal axes suffice.
 932673            for (int edge = 0; edge < 3 && selection.GapSign >= 0; edge++)
 74            {
 349375                int j = (edge + 1) % 3;
 349376                int k = (edge + 2) % 3;
 349377                Signed320 aj = BoxCylinderGeometry.GetComponent(geometry.Axis, j);
 349378                Signed320 ak = BoxCylinderGeometry.GetComponent(geometry.Axis, k);
 349379                WritePlaneDirection(j, k, aj, ak, direction, directionSigns);
 349380                KeepAxis(geometry, direction, directionSigns, ref selection);
 349381                WritePlaneDirection(j, k, WideArithmetic.Negate(ak), aj, direction, directionSigns);
 349382                KeepAxis(geometry, direction, directionSigns, ref selection);
 3487083                for (int cone = 0; cone < 4 && selection.GapSign >= 0; cone++)
 84                {
 1394285                    int corner = ((cone & 1) != 0 ? 1 << j : 0)
 1394286                        | ((cone & 2) != 0 ? 1 << k : 0);
 8365287                    for (int cap = -1; cap <= 1; cap += 2)
 88                    {
 2788489                        WideAxis3 c = geometry.GetCapVertexOffset(corner, cap);
 2788490                        Signed320 cj = BoxCylinderGeometry.GetComponent(c, j);
 2788491                        Signed320 ck = BoxCylinderGeometry.GetComponent(c, k);
 2788492                        int sj = (cone & 1) != 0 ? 1 : -1;
 2788493                        int sk = (cone & 2) != 0 ? 1 : -1;
 2788494                        if (cj.Sign * sj * ck.Sign * sk > 0)
 95                        {
 2283696                            WritePlaneDirection(j, k, cj, ck, direction, directionSigns);
 2283697                            KeepAxis(geometry, direction, directionSigns, ref selection);
 98                        }
 99                    }
 100                }
 101            }
 102        }
 1178103        RestoreRawGapDenominator(geometry, bestValues, direction);
 1178104        bestGapSign = selection.GapSign;
 1178105        localSupportSigns = selection.LocalSigns;
 1178106        return bestGapSign >= 0;
 107    }
 108
 109    internal static bool IsEdgeConeDominated(in BoxCylinderGeometry geometry,
 110        int edge, int first, int second, int firstSign, int secondSign, int cap,
 111        in ConvexContactCandidate best)
 112    {
 2464113        int corner = (firstSign > 0 ? 1 << first : 0) | (secondSign > 0 ? 1 << second : 0);
 2464114        WideAxis3 offset = geometry.GetCapVertexOffset(corner, cap);
 2464115        Signed320 cu = BoxCylinderGeometry.GetComponent(offset, first);
 2464116        Signed320 cv = BoxCylinderGeometry.GetComponent(offset, second);
 2464117        if (cu.Sign * firstSign < 0 || cv.Sign * secondSign < 0)
 572118            return false;
 119
 120        // In this signed unit quadrant, cu*nu+cv*nv >= min(cu,cv).
 121        // The projected disk support is at least R*|a_edge|/sqrt(U).
 122        // Their sum bounds the entire cone, including both parameter charts.
 123        // It is positive here: radius and the edge's axis component are nonzero.
 1892124        Span<ulong> values = stackalloc ulong[ConvexContactCandidate.Slots * Words];
 1892125        Span<int> signs = stackalloc int[ConvexContactCandidate.Slots];
 3784126        values.Clear(); signs.Clear();
 3784127        values[0] = 1; signs[0] = 1; // Comparison-only candidate: arbitrary unit normal.
 1892128        Span<ulong> work = stackalloc ulong[4 * Words];
 5676129        Span<ulong> m = Slot(work, 0), other = Slot(work, 1), re = Slot(work, 2);
 1892130        Span<ulong> product = Slot(work, 3);
 3784131        Import(cu, m); Import(cv, other);
 1892132        if (WideArithmetic.CompareMagnitudeEqualLength(m, other) > 0)
 935133            other.CopyTo(m);
 1892134        Import(geometry.Radius, other);
 1892135        Import(BoxCylinderGeometry.GetComponent(geometry.Axis, edge), product);
 1892136        WideArithmetic.MultiplyMagnitudes(other, product, re);
 1892137        Import(geometry.Axis.SquaredLength, Slot(values, 9));
 138        // b²=(A+B*sqrt(U))/D, A=m²U+(Re)², B=2mRe, D=S²U.
 139        // Conservative geometry bounds put A below 817 bits and D below 682.
 1892140        WideArithmetic.MultiplyMagnitudes(m, m, product);
 1892141        WideArithmetic.MultiplyMagnitudes(product, Slot(values, 9), Slot(values, 7));
 1892142        WideArithmetic.MultiplyMagnitudes(re, re, product);
 1892143        WideArithmetic.AddMagnitudeInto(product, Slot(values, 7));
 1892144        WideArithmetic.MultiplyMagnitudes(m, re, product);
 1892145        WideArithmetic.AddEqualMagnitudes(product, product, Slot(values, 8));
 1892146        Import(Signed320.ExtendValue(geometry.RawScale), other);
 1892147        WideArithmetic.MultiplyMagnitudes(other, other, product);
 1892148        WideArithmetic.MultiplyMagnitudes(product, Slot(values, 9), Slot(values, 10));
 1892149        signs[7] = 1;
 1892150        signs[8] = IsZero(m) ? 0 : 1;
 1892151        var bound = new ConvexContactCandidate(values, signs, 1);
 1892152        return WideConvexPrismRelations.CompareConvexContactCandidates(bound, best) >= 0;
 153    }
 154
 155    private static void KeepAxis(in BoxCylinderGeometry geometry,
 156        scoped ReadOnlySpan<ulong> direction, scoped ReadOnlySpan<int> directionSigns, ref Selection selection)
 157    {
 47442158        if (selection.HasValue && selection.GapSign < 0)
 9159            return;
 47433160        if (directionSigns[0] == 0 && directionSigns[1] == 0 && directionSigns[2] == 0)
 3129161            return;
 44304162        Span<ulong> values = stackalloc ulong[ConvexContactCandidate.Slots * Words];
 44304163        Span<int> signs = stackalloc int[ConvexContactCandidate.Slots];
 44304164        BuildRankingAxis(geometry, direction, directionSigns, values, signs, out int gapSign, out int orientation);
 44304165        var candidate = new ConvexContactCandidate(values, signs, gapSign);
 44304166        if (selection.HasValue && WideConvexPrismRelations.CompareConvexContactCandidates(candidate,
 44304167                new ConvexContactCandidate(selection.Values, selection.Signs, selection.GapSign)) >= 0)
 42769168            return;
 1535169        values.CopyTo(selection.Values);
 1535170        signs.CopyTo(selection.Signs);
 1535171        selection.GapSign = gapSign;
 1535172        selection.HasValue = true;
 1535173        selection.LocalSigns = new Vector3d((Fixed64)(orientation * directionSigns[0]),
 1535174            (Fixed64)(orientation * directionSigns[1]), (Fixed64)(orientation * directionSigns[2]));
 1535175    }
 176
 177    private static void BuildRankingAxis(in BoxCylinderGeometry geometry,
 178        ReadOnlySpan<ulong> direction, ReadOnlySpan<int> directionSigns,
 179        Span<ulong> values, Span<int> signs, out int gapSign, out int orientation)
 180    {
 44304181        values.Clear();
 44304182        signs.Clear();
 44304183        Span<ulong> work = stackalloc ulong[9 * Words];
 132912184        Span<ulong> u = Slot(work, 0), n2 = Slot(work, 1), axial = Slot(work, 2);
 132912185        Span<ulong> rational = Slot(work, 3), temporary = Slot(work, 4), product = Slot(work, 5);
 88608186        Span<ulong> p = Slot(work, 6), q = Slot(work, 7);
 44304187        Span<ulong> radius = Slot(work, 8);
 44304188        Import(geometry.Axis.SquaredLength, u);
 44304189        SumSquares(direction, n2);
 44304190        Dot(geometry.CenterDifference, direction, directionSigns, rational, out int centerSign);
 44304191        orientation = centerSign < 0 ? -1 : 1;
 44304192        int rationalSign = IsZero(rational) ? 0 : -1;
 44304193        Dot(geometry.Half, direction, directionSigns, temporary, out _);
 44304194        Add(temporary, 1, rational, ref rationalSign);
 354432195        for (int component = 0; component < 3; component++)
 196        {
 132912197            Import(BoxCylinderGeometry.GetComponent(geometry.HalfExtents, component), temporary);
 132912198            WideArithmetic.MultiplyMagnitudes(temporary, Slot(direction, component), product);
 132912199            Add(product, 1, rational, ref rationalSign);
 200        }
 44304201        if (IsZero(rational))
 51202            rationalSign = 0;
 203
 204        // T=h_box(n)+|H.n|-|d.n|, U=a², V=U*n²-(a.n)² >= 0.
 205        // g²=(U*T²+R²*V+2*T*R*sqrt(U*V))/(U*S²*n²). Cancel U
 206        // before construction instead of carrying it through every ranking
 207        // product; R>=0 permits moving it outside the radical exactly.
 208        // Directions <2^577 and coordinates <2^237 give |T|<2^817,
 209        // U<2^342, R<2^235, S<2^170, V<2^1498 give A<2^1977,
 210        // B<2^1053, C<2^1840 and D<2^1838. All fit forty words,
 211        // including coefficient-sum carries.
 212        // All directions share positive U*S². Rank with denominator n²,
 213        // preserving signed order and exact ties; restore raw units once
 214        // before either TryGetBest exit. These private scaled values never
 215        // enter cone pruning, edge-root ranking or public depth rounding.
 44304216        Dot(geometry.Axis, direction, directionSigns, axial, out _);
 44304217        WideArithmetic.MultiplyMagnitudes(axial, axial, temporary);
 44304218        WideArithmetic.MultiplyMagnitudes(u, n2, product);
 44304219        WideArithmetic.SubtractEqualMagnitudes(product, temporary, q);
 44304220        Import(geometry.Radius, radius);
 44304221        WideArithmetic.MultiplyMagnitudes(radius, radius, temporary);
 44304222        WideArithmetic.MultiplyMagnitudes(temporary, q, product);
 44304223        WideArithmetic.MultiplyMagnitudes(rational, rational, p);
 44304224        WideArithmetic.MultiplyMagnitudes(p, u, temporary);
 44304225        WideArithmetic.AddEqualMagnitudes(temporary, product, Slot(values, 7));
 44304226        gapSign = rationalSign < 0 ? WideArithmetic.CompareMagnitudeEqualLength(product, temporary)
 44304227            : rationalSign > 0 || !IsZero(product) ? 1 : 0;
 44304228        WideArithmetic.MultiplyMagnitudes(rational, radius, p);
 44304229        WideArithmetic.AddEqualMagnitudes(p, p, Slot(values, 8));
 44304230        WideArithmetic.MultiplyMagnitudes(u, q, Slot(values, 9));
 44304231        n2.CopyTo(Slot(values, 10));
 44304232        signs[7] = IsZero(Slot(values, 7)) ? 0 : 1;
 44304233        signs[8] = IsZero(p) ? 0 : rationalSign;
 44304234        CylinderContactAlgebra.WriteWorldDirection(geometry.WorldBasis, direction,
 44304235            directionSigns, values, signs, orientation);
 44304236    }
 237
 238    private static void RestoreRawGapDenominator(in BoxCylinderGeometry geometry,
 239        Span<ulong> values, Span<ulong> scratch)
 240    {
 241        // The direction workspace is dead at either exit. Reuse its three
 242        // slots to reconstruct the original D=U*S²*n² (<2^1838).
 3549243        Span<ulong> u = Slot(scratch, 0), scale = Slot(scratch, 1), product = Slot(scratch, 2);
 1183244        Import(geometry.Axis.SquaredLength, u);
 1183245        Import(Signed320.ExtendValue(geometry.RawScale), scale);
 1183246        WideArithmetic.MultiplyMagnitudes(scale, scale, product);
 1183247        WideArithmetic.MultiplyMagnitudes(u, product, scale);
 1183248        WideArithmetic.MultiplyMagnitudes(Slot(values, 10), scale, product);
 1183249        product.CopyTo(Slot(values, 10));
 1183250    }
 251
 252    private static bool HasVertexRimSeparation(in BoxCylinderGeometry geometry, int corner, int cap)
 253    {
 18720254        WideAxis3 c = geometry.GetCapVertexOffset(corner, cap);
 18720255        Signed576 axial = WideAxis3.Dot(geometry.Axis, c);
 18720256        if (axial.Sign * cap > 0)
 13838257            return false;
 4882258        Span<ulong> work = stackalloc ulong[10 * Words];
 14646259        Span<ulong> u = Slot(work, 0), q = Slot(work, 1), radius = Slot(work, 2);
 14646260        Span<ulong> first = Slot(work, 3), second = Slot(work, 4), root = Slot(work, 5);
 9764261        Span<ulong> radial = work.Slice(6 * Words, 3 * Words), coordinate = Slot(work, 9);
 4882262        Span<int> radialSigns = stackalloc int[3];
 4882263        Import(geometry.Axis.SquaredLength, u);
 4882264        Import(c.SquaredLength, first);
 4882265        WideArithmetic.MultiplyMagnitudes(u, first, q);
 4882266        Import(axial, first);
 4882267        WideArithmetic.MultiplyMagnitudes(first, first, second);
 4882268        WideArithmetic.SubtractEqualMagnitudes(q, second, q);
 4882269        if (IsZero(q))
 4270            return false; // A cap pole already supplies this normal.
 4878271        Import(geometry.Radius, radius);
 4878272        WideArithmetic.MultiplyMagnitudes(radius, radius, first);
 4878273        WideArithmetic.MultiplyMagnitudes(first, u, second);
 4878274        int radialComparison = WideArithmetic.CompareMagnitudeEqualLength(q, second);
 275        // At equality the residual is purely axial (or zero). The cylinder
 276        // pole was tested first, including separation; this is not a new rim.
 4878277        if (radialComparison <= 0)
 224278            return false;
 4654279        WideArithmetic.MultiplyMagnitudes(u, q, root);
 4654280        ProjectPerpendicular(geometry.Axis, c, radial, radialSigns);
 281        // A smooth vertex/rim patch has zero meridional curvature. A positive
 282        // support gap therefore has a negative second variation there and
 283        // cannot be the minimum penetration. Only its outward residual can
 284        // separate: n=R*P_a(c)-c*sqrt(U*q). Test its exact support cone.
 9522285        for (int component = 0; component < 3; component++)
 286        {
 4756287            WideArithmetic.MultiplyMagnitudes(radius, Slot(radial, component), first);
 4756288            Signed320 value = BoxCylinderGeometry.GetComponent(c, component);
 4756289            Import(value, coordinate);
 4756290            int sign = WideConvexPrismRelations.GetConvexContactCandidateQuadraticSign(
 4756291                first, radialSigns[component], coordinate, -value.Sign, root);
 4756292            int supportSign = (corner & (1 << component)) != 0 ? 1 : -1;
 4756293            if (sign * supportSign < 0)
 4649294                return false;
 295        }
 5296        return true;
 297    }
 298
 299    private static void WritePlaneDirection(int j, int k, Signed320 u, Signed320 v,
 300        Span<ulong> direction, Span<int> signs)
 301    {
 29822302        direction.Clear();
 29822303        signs.Clear();
 59644304        Import(u, Slot(direction, j)); signs[j] = u.Sign;
 59644305        Import(v, Slot(direction, k)); signs[k] = v.Sign;
 29822306    }
 307
 3538308    private static WideAxis3 Unit(int axis) => new(
 3538309        Signed320.ExtendValue(Signed192.Signed(axis == 0 ? 1 : 0)),
 3538310        Signed320.ExtendValue(Signed192.Signed(axis == 1 ? 1 : 0)),
 3538311        Signed320.ExtendValue(Signed192.Signed(axis == 2 ? 1 : 0)));
 312    private ref struct Selection
 313    {
 314        internal Span<ulong> Values;
 315        internal Span<int> Signs;
 316        internal int GapSign;
 317        internal bool HasValue;
 318        internal Vector3d LocalSigns;
 319        internal Selection(Span<ulong> values, Span<int> signs)
 320        {
 1183321            Values = values;
 1183322            Signs = signs;
 1183323            GapSign = 0;
 1183324            HasValue = false;
 1183325            LocalSigns = default;
 1183326        }
 327    }
 328}