using System;
using System.Collections.Generic;
using OpenNest.Engine.Jobs;
using OpenNest.Geometry;
using OpenNest.Math;
namespace OpenNest.Engine.Qwen38FlashNext.Engine;
using Math = System.Math;
/// A committed placement: model, orientation, and origin position on the sheet.
internal readonly struct PlacedPart
{
public PlacedPart(PartModel model, OrientationModel orientation, double x, double y)
{
Model = model;
Orientation = orientation;
X = x;
Y = y;
}
public PartModel Model { get; }
public OrientationModel Orientation { get; }
public double X { get; }
public double Y { get; }
}
internal readonly struct PlacementResult
{
public PlacementResult(PlacedPart part)
{
Part = part;
}
public PlacedPart Part { get; }
}
///
/// One sheet packed by this engine's own algorithm: bottom-left greedy insertion of
/// convex NFP corner candidates.
///
/// For each candidate orientation the packer builds a convex No-Fit-Polygon per
/// already-placed part as placedHull (+) disk(spacing) (+) reflect(candidateHull) -
/// a superset of the true NFP because hulls ignore concavities and cutouts, so an
/// anchor outside every NFP plus inside the anchor work-box is always legal ("strict"
/// certification). An anchor inside an NFP is still accepted when the exact material
/// gate says the parts clear: that gate inflates the placed part's material by the
/// spacing (holes shrunk, closed holes treated solid) and tests it against the
/// candidate's raw material with holes subtracted - the same inflation rule the
/// benchmark validator uses, so interlocking concave parts are recovered without ever
/// accepting an overlap or a spacing violation.
///
///
/// Candidate anchors are the corner points of the feasible region: the four anchor
/// work-box corners, every NFP vertex, and every NFP-edge/box-line crossing (slides).
/// Candidates are tried in ascending bottom-left order and the first legal one wins.
///
///
internal sealed class SheetPacker
{
private const int MaxHullVertices = 40;
private readonly double _workLeft;
private readonly double _workBottom;
private readonly double _workRight;
private readonly double _workTop;
// Per-placed-part caches (indexed by placement order).
private readonly List _inflatedBounds = new();
private readonly List<(Polygon Perimeter, List Holes)> _placedGate = new();
///
/// Fast-path outlines of the placed parts' gate material in world coordinates,
/// parallel to (O(1) translation of the shared per-
/// orientation template). A hit certifies the two
/// outer shells - hence both materials - are clear and skips the exact
/// gate, which at fine flattening triangulates thousands
/// of edges per call. Null when the outline has no usable ring.
///
private readonly List _placedGateFast = new();
// NFP caches: (placedIndex, orientationId) -> forbidden-anchor contour.
private readonly Dictionary<(int, int), ConvexContour?> _nfpCache = new();
private readonly Dictionary _orientationIds = new();
///
/// Per-orientation cache of NFP/NFP valley anchors (anchors touching two placed
/// parts at once). A committed part's NFP never changes and
/// only grows, so each (i, j) pair is intersected exactly once per orientation
/// instead of once per candidate enumeration - the re-sweep was the dominant cost
/// on crowded sheets (O(placed^2 * edges^2) per insert attempt).
///
private sealed class ValleyCache
{
public int BuiltThrough;
public readonly List<(double X, double Y)> Valleys = new();
}
private readonly Dictionary _valleyCaches = new();
private readonly Dictionary _reflectedHulls = new();
private readonly Dictionary _placedHullDisk = new();
private ConvexContour? _disk;
// Uniform spatial grid over placed parts' inflated bounds: IsLegal only tests the
// parts whose cells touch the candidate's cells, so legality stays near-constant
// as a sheet fills instead of scanning every placed part.
private readonly double _cellSize;
private readonly int _gridCols;
private readonly int _gridRows;
private readonly List[] _grid;
private SheetPacker(NestPlateStock stock, PartPreparation prep, int stockIndex)
{
Stock = stock;
StockIndex = stockIndex;
Preparation = prep;
Spacing = stock.PartSpacing;
var left = stock.Quadrant is 1 or 4 ? 0.0 : -stock.Size.Length;
var bottom = stock.Quadrant is 1 or 2 ? 0.0 : -stock.Size.Width;
_workLeft = left + stock.EdgeSpacing.Left;
_workBottom = bottom + stock.EdgeSpacing.Bottom;
_workRight = left + stock.Size.Length - stock.EdgeSpacing.Right;
_workTop = bottom + stock.Size.Width - stock.EdgeSpacing.Top;
WorkWidth = _workRight - _workLeft;
WorkHeight = _workTop - _workBottom;
// Cells roughly the size of a mid-range part: a candidate usually touches 2-6.
_cellSize = System.Math.Max(1.0, System.Math.Min(WorkWidth, WorkHeight) / 12.0);
_gridCols = System.Math.Max(1, (int)System.Math.Ceiling(WorkWidth / _cellSize));
_gridRows = System.Math.Max(1, (int)System.Math.Ceiling(WorkHeight / _cellSize));
_grid = new List[_gridCols * _gridRows];
for (var i = 0; i < _grid.Length; i++)
_grid[i] = new List();
}
public static SheetPacker Create(NestPlateStock stock, PartPreparation prep, int stockIndex) =>
new(stock, prep, stockIndex);
public NestPlateStock Stock { get; }
public int StockIndex { get; }
public PartPreparation Preparation { get; }
public double Spacing { get; }
public double WorkWidth { get; }
public double WorkHeight { get; }
public List Placed { get; } = new();
public bool IsFull => Placed.Count >= MaxPartsPerSheet;
///
/// Safety cap on parts per sheet: real sheets never exceed this, and it bounds
/// per-insert NFP work and the validator's area budget on pathological jobs.
///
public const int MaxPartsPerSheet = 500;
/// True when the part's bounds can never fit this sheet in any orientation.
public bool CanEverFit(PartModel model)
{
foreach (var angle in PartPreparation.CandidateAngles(model))
{
var orientation = Preparation.Oriented(model, angle, 0);
if (orientation.Width <= WorkWidth + 1e-9 && orientation.Height <= WorkHeight + 1e-9)
return true;
}
return false;
}
internal long DiagInsertAttempts;
internal long DiagCandidateChecks;
internal long DiagGateCalls;
internal long DiagConvexRejections;
internal long DiagFastClears;
internal long DiagFastNull;
internal long DiagFastOverlapWithHoles;
internal long DiagFastOverlaps;
internal long DiagFastUnknowns;
public string DiagStats() =>
$"inserts={DiagInsertAttempts} checks={DiagCandidateChecks} gates={DiagGateCalls} " +
$"convexRej={DiagConvexRejections} fastClear={DiagFastClears} fastOver={DiagFastOverlaps} fastUnk={DiagFastUnknowns} fastNull={DiagFastNull} fastOverHoles={DiagFastOverlapWithHoles} triNull={DiagTriNull} triNullOut={DiagTriNullOut} triFallback={DiagTriFallback}";
///
/// Greedily insert an instance: best (bottom-left) legal corner over all candidate
/// orientations. Returns false (and changes nothing) when no legal position exists.
///
public bool TryInsert(PartModel model, out PlacementResult result)
{
result = default;
var bestScore = double.MaxValue;
PlacedPart? best = null;
DiagInsertAttempts++;
foreach (var angle in PartPreparation.CandidateAngles(model))
{
var orientation = Preparation.Oriented(model, angle, Spacing);
if (orientation.Width > WorkWidth + 1e-9 || orientation.Height > WorkHeight + 1e-9)
continue;
foreach (var (x, y) in OrderedCandidates(orientation))
{
var score = Score(orientation, x, y);
if (score >= bestScore)
continue; // no later candidate (same sort) can beat it
if (!IsLegal(orientation, x, y))
continue;
bestScore = score;
best = new PlacedPart(model, orientation, x, y);
break; // first legal in ascending-score order is this orientation's best
}
}
if (best == null)
return false;
Commit(best.Value);
result = new PlacementResult(best.Value);
return true;
}
private double Score(OrientationModel orientation, double x, double y) =>
x + orientation.MinX + (y + orientation.MinY) * 1.0001;
///
/// Corner candidates in deterministic ascending bottom-left order: anchor work-box
/// corners, NFP vertices, and NFP-edge/box-line crossings.
///
private List<(double x, double y)> OrderedCandidates(OrientationModel orientation)
{
var boxLeft = _workLeft - orientation.MinX;
var boxRight = _workRight - orientation.MaxX;
var boxBottom = _workBottom - orientation.MinY;
var boxTop = _workTop - orientation.MaxY;
var seen = new HashSet<(long, long)>();
var candidates = new List<(double, double)>(128);
// Math.Clamp throws when min > max, and a part that fits the work area to
// within floating-point noise can invert the anchor box by ~1e-14. Order the
// bounds so a degenerate box collapses to its single legal point.
var anchorMinX = Math.Min(boxLeft, boxRight);
var anchorMaxX = Math.Max(boxLeft, boxRight);
var anchorMinY = Math.Min(boxBottom, boxTop);
var anchorMaxY = Math.Max(boxBottom, boxTop);
void Add(double x, double y)
{
if (x < anchorMinX - 1e-9 || x > anchorMaxX + 1e-9 || y < anchorMinY - 1e-9 || y > anchorMaxY + 1e-9)
return;
x = Math.Clamp(x, anchorMinX, anchorMaxX);
y = Math.Clamp(y, anchorMinY, anchorMaxY);
if (!seen.Add(((long)Math.Round(x * 1e6), (long)Math.Round(y * 1e6))))
return;
candidates.Add((x, y));
}
Add(boxLeft, boxBottom);
Add(boxRight, boxBottom);
Add(boxLeft, boxTop);
Add(boxRight, boxTop);
for (var i = 0; i < Placed.Count; i++)
{
var nfp = NfpFor(i, orientation);
if (nfp == null)
continue;
var n = nfp.Count;
for (var v = 0; v < n; v++)
Add(nfp.X(v), nfp.Y(v));
// Slides: NFP edges crossing the anchor box border lines.
for (var v = 0; v < n; v++)
{
var ax = nfp.X(v);
var ay = nfp.Y(v);
var bx = nfp.X((v + 1) % n);
var by = nfp.Y((v + 1) % n);
CrossLine(ax, ay, bx, by, boxLeft, true, Add);
CrossLine(ax, ay, bx, by, boxRight, true, Add);
CrossLine(ax, ay, bx, by, boxBottom, false, Add);
CrossLine(ax, ay, bx, by, boxTop, false, Add);
}
}
// Valleys between two neighbors: NFP/NFP edge intersections are the anchors
// where the candidate touches two placed parts at once - the classic
// bottom-left stable corners the single-NFP candidates cannot produce.
foreach (var (vx, vy) in ValleysFor(orientation))
Add(vx, vy);
candidates.Sort(
(p, q) =>
{
var byY = p.Item2.CompareTo(q.Item2);
return byY != 0 ? byY : p.Item1.CompareTo(q.Item1);
}
);
return candidates;
}
///
/// Cached NFP/NFP valley anchors for one orientation, extended in place with the
/// pairs involving placements committed since the last call. Each (i, j) pair is
/// intersected once per orientation for the packer's lifetime.
///
private List<(double X, double Y)> ValleysFor(OrientationModel orientation)
{
if (!_valleyCaches.TryGetValue(orientation, out var cache))
{
cache = new ValleyCache();
_valleyCaches[orientation] = cache;
}
var count = Placed.Count;
for (var j = cache.BuiltThrough; j < count; j++)
{
var nfpB = NfpFor(j, orientation);
if (nfpB == null)
continue;
for (var i = 0; i < j; i++)
{
var nfpA = NfpFor(i, orientation);
if (nfpA == null || !nfpA.Bounds.Intersects(nfpB.Bounds))
continue;
var na = nfpA.Count;
var nb = nfpB.Count;
for (var va = 0; va < na; va++)
{
var a0x = nfpA.X(va);
var a0y = nfpA.Y(va);
var a1x = nfpA.X((va + 1) % na);
var a1y = nfpA.Y((va + 1) % na);
for (var vb = 0; vb < nb; vb++)
{
var b0x = nfpB.X(vb);
var b0y = nfpB.Y(vb);
var b1x = nfpB.X((vb + 1) % nb);
var b1y = nfpB.Y((vb + 1) % nb);
if (
Math.Max(a0x, a1x) < Math.Min(b0x, b1x)
|| Math.Max(b0x, b1x) < Math.Min(a0x, a1x)
|| Math.Max(a0y, a1y) < Math.Min(b0y, b1y)
|| Math.Max(b0y, b1y) < Math.Min(a0y, a1y)
)
continue;
var r = SegmentIntersect(
a0x, a0y, a1x, a1y,
b0x, b0y, b1x, b1y
);
if (r.HasValue)
cache.Valleys.Add((r.Value.X, r.Value.Y));
}
}
}
}
cache.BuiltThrough = count;
return cache.Valleys;
}
/// Proper or endpoint intersection of two segments, if any.
private static Vector? SegmentIntersect(
double ax,
double ay,
double bx,
double by,
double cx,
double cy,
double dx,
double dy
)
{
var rx = bx - ax;
var ry = by - ay;
var sx = dx - cx;
var sy = dy - cy;
var denom = rx * sy - ry * sx;
if (Math.Abs(denom) < 1e-12)
return null; // parallel
var t = ((cx - ax) * sy - (cy - ay) * sx) / denom;
var u = ((cx - ax) * ry - (cy - ay) * rx) / denom;
if (t < -1e-9 || t > 1 + 1e-9 || u < -1e-9 || u > 1 + 1e-9)
return null;
return new Vector(ax + t * rx, ay + t * ry);
}
private static void CrossLine(
double ax,
double ay,
double bx,
double by,
double at,
bool vertical,
Action add
)
{
var (ua, ub) = vertical ? (ax, bx) : (ay, by);
if (ua == ub)
return;
var t = (at - ua) / (ub - ua);
if (t < 0 || t > 1)
return;
var along = vertical ? ay + (by - ay) * t : ax + (bx - ax) * t;
if (vertical)
add(at, along);
else
add(along, at);
}
private void Commit(PlacedPart part)
{
Placed.Add(part);
var pad = Spacing;
_inflatedBounds.Add(
new Bounds(
part.X + part.Orientation.MinX - pad,
part.Y + part.Orientation.MinY - pad,
part.X + part.Orientation.MaxX + pad,
part.Y + part.Orientation.MaxY + pad
)
);
// Gate geometry: material inflated by spacing (holes shrunk) when positive,
// raw material at zero spacing; already in world coordinates.
var gatePerimeter = part.Orientation.InflatedPerimeter ?? part.Orientation.Perimeter;
var gateHoles = part.Orientation.InflatedPerimeter != null
? part.Orientation.InflatedHoles
: part.Orientation.Holes;
var worldPerimeter = (Polygon)gatePerimeter.Clone();
worldPerimeter.Offset(part.X, part.Y);
worldPerimeter.UpdateBounds();
var worldHoles = new List(gateHoles.Count);
foreach (var hole in gateHoles)
{
var h = (Polygon)hole.Clone();
h.Offset(part.X, part.Y);
h.UpdateBounds();
worldHoles.Add(h);
}
_placedGate.Add((worldPerimeter, worldHoles));
_placedGateFast.Add(part.Orientation.GateFast?.Translated(part.X, part.Y));
GridAdd(Placed.Count - 1, _inflatedBounds[^1]);
}
// ---- uniform spatial grid (cell -> placed indices) --------------------------
private void GridAdd(int placedIndex, in Bounds bounds)
{
var c0 = System.Math.Clamp(
(int)System.Math.Floor((bounds.MinX - _workLeft) / _cellSize),
0,
_gridCols - 1
);
var c1 = System.Math.Clamp(
(int)System.Math.Floor((bounds.MaxX - _workLeft) / _cellSize),
0,
_gridCols - 1
);
var r0 = System.Math.Clamp(
(int)System.Math.Floor((bounds.MinY - _workBottom) / _cellSize),
0,
_gridRows - 1
);
var r1 = System.Math.Clamp(
(int)System.Math.Floor((bounds.MaxY - _workBottom) / _cellSize),
0,
_gridRows - 1
);
for (var r = r0; r <= r1; r++)
for (var c = c0; c <= c1; c++)
_grid[r * _gridCols + c].Add(placedIndex);
}
private readonly HashSet _nearScratch = new();
private HashSet Near(in Bounds bounds)
{
_nearScratch.Clear();
var c0 = System.Math.Clamp(
(int)System.Math.Floor((bounds.MinX - _workLeft) / _cellSize),
0,
_gridCols - 1
);
var c1 = System.Math.Clamp(
(int)System.Math.Floor((bounds.MaxX - _workLeft) / _cellSize),
0,
_gridCols - 1
);
var r0 = System.Math.Clamp(
(int)System.Math.Floor((bounds.MinY - _workBottom) / _cellSize),
0,
_gridRows - 1
);
var r1 = System.Math.Clamp(
(int)System.Math.Floor((bounds.MaxY - _workBottom) / _cellSize),
0,
_gridRows - 1
);
for (var r = r0; r <= r1; r++)
for (var c = c0; c <= c1; c++)
foreach (var index in _grid[r * _gridCols + c])
_nearScratch.Add(index);
return _nearScratch;
}
///
/// Legality of one anchor. Outside every overlapping NFP is strict certification;
/// inside one still passes when the exact material gate clears (interlocking
/// concaves and cutouts that the convex NFP cannot represent).
///
private bool IsLegal(OrientationModel orientation, double x, double y)
{
if (
x + orientation.MinX < _workLeft - 1e-9
|| x + orientation.MaxX > _workRight + 1e-9
|| y + orientation.MinY < _workBottom - 1e-9
|| y + orientation.MaxY > _workTop + 1e-9
)
return false;
var pad = Spacing;
var candidate = new Bounds(
x + orientation.MinX - pad,
y + orientation.MinY - pad,
x + orientation.MaxX + pad,
y + orientation.MaxY + pad
);
// World-space candidate material, built at most once per anchor and only when
// a convex-NFP hit actually needs the exact gate; freed with the anchor.
(Polygon Perimeter, List Holes)? gate = null;
foreach (var i in Near(candidate))
{
var bounds = _inflatedBounds[i];
if (
candidate.MinX >= bounds.MaxX
|| candidate.MaxX <= bounds.MinX
|| candidate.MinY >= bounds.MaxY
|| candidate.MaxY <= bounds.MinY
)
continue;
// Fast rejection: outside the hull-based NFP the placed and candidate
// HULLS are at least spacing apart, and hulls contain materials, so the
// materials clear - a valid certification for any shape, holed or
// concave. Inside the NFP decides nothing by itself (the hull sum
// over-approximates for concaves and holes), but when both materials are
// convex solids with uncapped hulls the sum is exact (modulo the
// circumscribed disk's chord error, which only ever rejects a hair too
// much), so interior means overlap. Everything else pays the exact
// material gate.
DiagCandidateChecks++;
var nfp = NfpFor(i, orientation);
if (nfp == null)
{
if (!TryPairVerdict(orientation, x, y, i, out var nullNfpOverlap))
{
DiagGateCalls++;
gate ??= BuildCandidateGate(orientation, x, y);
nullNfpOverlap = MaterialOverlap(gate.Value, orientation, x, y, i);
}
if (nullNfpOverlap)
return false;
continue;
}
if (!nfp.ContainsPoint(x, y))
continue; // outside the conservative forbidden sum: certified clear
if (
orientation.IsConvexSolid
&& Placed[i].Orientation.IsConvexSolid
&& orientation.Hull.Count <= MaxHullVertices
&& Placed[i].Orientation.Hull.Count <= MaxHullVertices
)
{
DiagConvexRejections++;
return false; // exact convex-convex NFP interior: overlap
}
// Cheap world-bbox test against the placed gate material before paying
// for candidate gate construction or the clipper.
if (
!_placedGate[i]
.Perimeter.BoundingBox
.Intersects(orientation.Perimeter.BoundingBox.Translate(x, y))
)
continue;
// Fast shell relation against the placed gate outline: a certified clear
// skips the exact gate entirely (no Polygon clones, no triangulation), a
// certified overlap rejects without it. Hole-bearing pairs and touches fall
// through to the exact gate.
if (TryPairVerdict(orientation, x, y, i, out var fastOverlap))
{
if (fastOverlap)
return false;
continue;
}
DiagGateCalls++;
gate ??= BuildCandidateGate(orientation, x, y);
if (MaterialOverlap(gate.Value, orientation, x, y, i))
return false;
}
return true;
}
///
/// Fast outer-shell relation for one (candidate, placed) pair against the placed
/// part's spacing-inflated gate outline, deciding whether the exact material gate
/// must run. A CERTIFIED verdict skips it: disjoint shells mean no material overlap
/// (holes only remove material), and a shell crossing/containment between two
/// hole-free polygons IS a positive-area material overlap. Hole-bearing pairs whose
/// shells overlap and near-degenerate touches fall through to the exact gate.
///
private bool TryPairVerdict(
OrientationModel orientation,
double x,
double y,
int placedIndex,
out bool overlap
)
{
overlap = false;
var placedFast = _placedGateFast[placedIndex];
var candidateFast = orientation.PerimeterFast;
if (placedFast == null || candidateFast == null)
{
DiagFastNull++;
return false;
}
var relation = FastPoly.Relate(candidateFast.Translated(x, y), placedFast);
if (relation == FastPoly.FastRelation.Overlap)
DiagFastOverlapWithHoles++;
switch (relation)
{
case FastPoly.FastRelation.Clear:
DiagFastClears++;
return true; // certified clear (spacing included in the placed gate)
case FastPoly.FastRelation.Overlap
when orientation.Holes.Count == 0 && _placedGate[placedIndex].Holes.Count == 0:
DiagFastOverlaps++;
overlap = true; // certified overlap: shells share area, nothing to subtract
return true;
default:
DiagFastUnknowns++;
return false; // exact gate must decide
}
}
private static readonly bool VerifyFastClear =
Environment.GetEnvironmentVariable("QWEN_VERIFY_FASTCLEAR") == "1";
internal long DiagFastClearMismatch;
internal long DiagTriMismatch;
internal long DiagTriNull;
internal long DiagTriNullOut;
internal long DiagTriFallback;
///
/// Exact clearance gate against one placed part: placed gate material (inflated by
/// spacing when positive) versus the candidate's raw material with holes subtracted.
///
private bool MaterialOverlap(
(Polygon Perimeter, List Holes) gate,
OrientationModel orientation,
double x,
double y,
int placedIndex
)
{
var placed = _placedGate[placedIndex];
if (!placed.Perimeter.BoundingBox.Intersects(orientation.Perimeter.BoundingBox.Translate(x, y)))
return false;
if (!gate.Perimeter.BoundingBox.Intersects(placed.Perimeter.BoundingBox))
return false;
var placedPart = Placed[placedIndex];
// Allocation-free path: both sides carry cached triangulations in their local
// frames, so the test clips triangles with the anchors as plain translations.
var candTris = orientation.MaterialTris;
var placedTris = placedPart.Orientation.GateTris;
if (candTris == null || placedTris == null)
DiagTriNull++;
if (candTris != null && placedTris != null)
{
var cached = candTris.HasOverlap(
placedTris, x, y, placedPart.X, placedPart.Y
);
if (cached.HasValue)
{
if (VerifyFastClear)
{
var truth = Collision.HasOverlap(
gate.Perimeter, placed.Perimeter, gate.Holes, placed.Holes
);
if (truth != cached.Value)
{
DiagTriMismatch++;
System.IO.File.AppendAllText(
"/tmp/triset_mismatch.log",
$"cached={cached.Value} truth={truth} candAng={orientation.Angle:F4} at ({x:F8},{y:F8}) " +
$"placedAng={placedPart.Orientation.Angle:F4} at ({placedPart.X:F8},{placedPart.Y:F8}) " +
$"candTris={candTris} candVerts={orientation.Perimeter.Vertices.Count} holes={orientation.Holes.Count} " +
$"placedVerts={placedPart.Orientation.Perimeter.Vertices.Count} placedHoles={placedPart.Orientation.Holes.Count}\n"
);
}
}
return cached.Value;
}
DiagTriNullOut++;
// Scratch overflow: fall through to the Polygon gate.
}
DiagTriFallback++;
return Preparation.MaterialOverlapMemo(
placedPart.Orientation,
placedPart.X,
placedPart.Y,
orientation,
x,
y,
() => Collision.HasOverlap(
gate.Perimeter,
placed.Perimeter,
gate.Holes,
placed.Holes
)
);
}
private (Polygon, List) BuildCandidateGate(OrientationModel orientation, double x, double y)
{
var perimeter = (Polygon)orientation.Perimeter.Clone();
perimeter.Offset(x, y);
perimeter.UpdateBounds();
var holes = new List(orientation.Holes.Count);
foreach (var hole in orientation.Holes)
{
var h = (Polygon)hole.Clone();
h.Offset(x, y);
h.UpdateBounds();
holes.Add(h);
}
return (perimeter, holes);
}
private int OrientationId(OrientationModel orientation)
{
if (!_orientationIds.TryGetValue(orientation, out var id))
{
id = _orientationIds.Count;
_orientationIds[orientation] = id;
}
return id;
}
///
/// Convex NFP of forbidden anchors: placedHull (+) disk(spacing) (+) reflect(candidateHull).
///
private ConvexContour? NfpFor(int placedIndex, OrientationModel orientation)
{
var key = (placedIndex, OrientationId(orientation));
if (_nfpCache.TryGetValue(key, out var cached))
return cached;
ConvexContour? result;
try
{
if (!_placedHullDisk.TryGetValue(placedIndex, out var placedDisk))
{
var placed = Placed[placedIndex];
var hull = placed.Orientation.Hull;
var capped = CapHull(hull, placed.X, placed.Y);
var placedHull = ConvexContour.FromVertices(capped);
placedDisk = Spacing > Tolerance.Epsilon
? NfpGeometry.Minkowski(placedHull, Disk())
: placedHull;
_placedHullDisk[placedIndex] = placedDisk;
}
if (!_reflectedHulls.TryGetValue(orientation, out var reflected))
{
var capped = CapHull(orientation.Hull, 0, 0);
reflected = NfpGeometry.Reflect(ConvexContour.FromVertices(capped));
_reflectedHulls[orientation] = reflected;
}
result = NfpGeometry.Minkowski(placedDisk, reflected);
}
catch (Exception)
{
// A degenerate Minkowski sum removes the fast rejection for this pair;
// the material gate still enforces correctness.
result = null;
}
_nfpCache[key] = result;
return result;
}
///
/// Bounded-size convex SUPERSET of (translated by dx/dy).
/// When the hull is dense, keep every k-th vertex, then shift each chord's
/// supporting line outward by the chord's maximum sagitta (the largest distance of
/// any dropped vertex to its chord). Every dropped vertex lies within the sagitta
/// of its chord, so the offset half-plane intersection contains the original hull
/// and the NFP built from it stays a conservative superset of the forbidden anchors.
///
private static List CapHull(List hull, double dx, double dy)
{
var n = hull.Count;
var shifted = new List(n);
for (var i = 0; i < n; i++)
shifted.Add(new Vector(hull[i].X + dx, hull[i].Y + dy));
if (n <= MaxHullVertices)
return shifted;
// Chord (v_i, v_{i+k}) for i in steps of k, with each chord's outward shift:
// the max perpendicular distance from any vertex it spans to the chord line.
var k = (int)Math.Ceiling(n / (double)MaxHullVertices);
var lines = new List<(double ax, double ay, double bx, double by, double shift)>();
for (var i = 0; i < n; i += k)
{
var a = shifted[i];
var b = shifted[(i + k) % n];
var span = Math.Min(k, n - i);
var sagitta = 0.0;
var length = Math.Sqrt((b.X - a.X) * (b.X - a.X) + (b.Y - a.Y) * (b.Y - a.Y));
if (length > 1e-12)
for (var j = 1; j < span; j++)
{
var p = shifted[i + j];
var distance = Math.Abs(Cross(a.X, a.Y, b.X, b.Y, p)) / length;
if (distance > sagitta)
sagitta = distance;
}
lines.Add((a.X, a.Y, b.X, b.Y, sagitta));
}
// Sutherland-Hodgman from a generous bounding box; the interior of each chord
// is the CCW left side, shifted outward (left) by the sagitta.
var minX = double.MaxValue;
var minY = double.MaxValue;
var maxX = double.MinValue;
var maxY = double.MinValue;
foreach (var v in shifted)
{
if (v.X < minX)
minX = v.X;
if (v.X > maxX)
maxX = v.X;
if (v.Y < minY)
minY = v.Y;
if (v.Y > maxY)
maxY = v.Y;
}
var margin = Math.Max(1.0, Math.Max(maxX - minX, maxY - minY));
var polygon = new List
{
new(minX - margin, minY - margin),
new(maxX + margin, minY - margin),
new(maxX + margin, maxY + margin),
new(minX - margin, maxY + margin),
};
foreach (var (ax, ay, bx, by, shift) in lines)
{
if (polygon.Count == 0)
return shifted; // degenerate; fall back to full hull
// Shift the line perpendicular away from the interior (CCW: interior is left).
var edgeX = bx - ax;
var edgeY = by - ay;
var length = Math.Sqrt(edgeX * edgeX + edgeY * edgeY);
if (length <= 1e-12)
continue;
var nx = edgeY / length;
var ny = -edgeX / length;
var ox = ax + nx * shift;
var oy = ay + ny * shift;
var input = polygon;
polygon = new List();
for (var i = 0; i < input.Count; i++)
{
var current = input[i];
var next = input[(i + 1) % input.Count];
var currentInside = Cross(ox, oy, ox + edgeX, oy + edgeY, current) >= 0;
var nextInside = Cross(ox, oy, ox + edgeX, oy + edgeY, next) >= 0;
if (currentInside)
{
polygon.Add(current);
if (!nextInside)
polygon.Add(Intersect(ox, oy, ox + edgeX, oy + edgeY, current, next));
}
else if (nextInside)
{
polygon.Add(Intersect(ox, oy, ox + edgeX, oy + edgeY, current, next));
}
}
}
return polygon.Count >= 3 ? polygon : shifted;
}
private static double Cross(double ax, double ay, double bx, double by, Vector p) =>
(bx - ax) * (p.Y - ay) - (by - ay) * (p.X - ax);
private static Vector Intersect(
double ax,
double ay,
double bx,
double by,
Vector p,
Vector q
)
{
var dx1 = bx - ax;
var dy1 = by - ay;
var dx2 = q.X - p.X;
var dy2 = q.Y - p.Y;
var cross = dx1 * dy2 - dy1 * dx2;
if (Math.Abs(cross) < 1e-300)
return p;
var t = ((p.X - ax) * dy2 - (p.Y - ay) * dx2) / cross;
return new Vector(ax + t * dx1, ay + t * dy1);
}
private ConvexContour Disk() =>
// Circumscribed so the polygon contains the true spacing disk: the NFP stays a
// conservative superset of the forbidden-anchor region.
_disk ??= ConvexContour.Disk(Spacing / Math.Cos(Math.PI / 24), 24);
}