{-# LANGUAGE AllowAmbiguousTypes #-} {-# LANGUAGE BangPatterns #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE NamedFieldPuns #-} {-# LANGUAGE ScopedTypeVariables #-} {-# LANGUAGE TypeApplications #-} -- | The normalization procedure of the flip rewrite system. module Moonlight.Triangulation.Internal.DcelOperations.Normalize ( drainLegalization ) where import Control.Monad (when) import Control.Monad.ST (ST) import Data.Bits ((.&.), shiftR) import Data.STRef (readSTRef) import qualified Data.Vector.Unboxed.Mutable as MUV import Moonlight.Triangulation.Internal.DcelOperations.CandidateArena ( genericCandidate , growLegalizationArena , starCandidate ) import Moonlight.Triangulation.Internal.DcelOperations.FlipRewrite (applyFlip) import Moonlight.Triangulation.Internal.DcelOperations.FlipRule ( LegalizationLaw (..) , diagonalFires ) import Moonlight.Triangulation.Internal.DcelOperations.Twin (reverseIndex) import Moonlight.Triangulation.Internal.Mutable ( MutableDcel (..) , readConstraint , readFace , readNext , readOrigin , readPointX , readPointY , readPrevious ) import Moonlight.Triangulation.Internal.OperationState ( Counter (..) , OperationState , legalizationArena , storeLegalizationArena ) import Moonlight.Triangulation.Internal.PackedIndex (packIndex) import Moonlight.Triangulation.Internal.Probe (KnownProbe (..)) -- | The normalization procedure of a confluent terminating rewrite system, and -- one canonical legalization engine because a normalization procedure is what -- it is. -- -- The objects are the triangulations of a fixed point set. The single rule is -- the Lawson flip: a locally illegal diagonal is replaced by the other -- diagonal of its quadrilateral. A normal form is a mesh with no illegal -- diagonal left to fire on. -- -- /Termination/ is by the lifted-paraboloid potential. Send each point to -- @(x, y, x² + y²)@ and read a triangulation as a piecewise-linear surface -- over the point set; @illegalDiagonal@ is exactly the test that the flip -- lowers that surface, so every rewrite strictly decreases it, and a finite -- point set has finitely many triangulations. Exact cocircularity is the one -- case where the surface does not move — the four lifted points are coplanar -- and both diagonals give the same surface — so there the potential is the -- diagonal's own key order, and the rule fires only downward in it. Without -- that tie-break a cocircular quadrilateral flips forever. -- -- /Confluence/ is Delaunay's theorem, in its strong form: a triangulation with -- no locally illegal diagonal is globally Delaunay. Local normality is thus -- global normality, the normal form is unique, and every rewrite order reaches -- it. That is what licenses the arena below to be a LIFO stack rather than a -- priority queue, and it is why callers may seed it in whatever order is -- cheapest to produce — the fan first, or the hull turns first, or both -- interleaved — without any of them changing the mesh that comes out. -- -- All callers differ only in how they seed the arena; topology mutation and -- propagation have exactly one owner. -- The stack top, the maximum top, and the flip count are strict loop -- variables, returned once when the drain finishes — the mesh reports nothing -- per candidate, and 'applyFlip' reports nothing at all. The phantom -- 'KnownProbe' parameter counts popped candidates for the instrumented lane -- and is erased everywhere else. -- -- A popped candidate is read once. Turning it against the star vertex, judging -- it, and rewriting it are three questions about the same two half-edge -- records, and the apex the turn looks for is the apex the judgement needs, so -- one pass over the quadrilateral answers all three. drainLegalization :: forall p s vertex directed undirected face . KnownProbe p => MutableDcel s vertex directed undirected face -> OperationState s -> Int -> Int -> LegalizationLaw -> ST s (Int, Int) drainLegalization mutable operation seededTop starVertex law = do -- No constraint can appear during a drain, so a mesh holding none at entry -- never needs the per-candidate protection read. constrained <- readSTRef (mdConstraintCount mutable) initialArena <- legalizationArena operation -- A star candidate is turned so that the inserted vertex is the apex -- opposite the diagonal, so that apex is the same vertex on every star -- candidate the drain pops and its two coordinates are read once here -- instead of once per candidate. starX <- if starVertex < 0 then pure 0 else readPointX mutable starVertex starY <- if starVertex < 0 then pure 0 else readPointY mutable starVertex let !guarded = constrained /= 0 -- Everything the cavity fan did not create is pinned. The border loop's -- legality is already decided by the outside triangle it keeps, and -- testing it against a neighbourhood that is still inverted could only -- produce a spurious verdict. 'ValidMesh' pins nothing, which is the -- floor no undirected index falls below. !floorPair = case law of ValidMesh -> 0 CavityRepair floorEdge -> floorEdge loop !arena !top !maxTop !flips !candidates | top <= 0 = pure (flips, maxTop, candidates, arena) | otherwise = do let !nextTop = top - 1 packedWord <- MUV.unsafeRead arena nextTop let !packed = fromIntegral packedWord :: Int !rawEdge = packed `shiftR` 1 !isStar = packed .&. 1 == 0 eligible <- if rawEdge `shiftR` 1 < floorPair then pure False else if guarded then not <$> readConstraint mutable rawEdge else pure True if not eligible then loop arena nextTop maxTop flips (probeBump @p candidates) else do let !rawTwin = reverseIndex rawEdge rawFace <- readFace mutable rawEdge rawTwinFace <- readFace mutable rawTwin if rawFace == 0 || rawTwinFace == 0 then loop arena nextTop maxTop flips (probeBump @p candidates) else do rawBefore <- readPrevious mutable rawEdge rawTwinBefore <- readPrevious mutable rawTwin rawApex <- readOrigin mutable rawBefore rawTwinApex <- readOrigin mutable rawTwinBefore -- A star candidate names an undirected edge and has to be -- turned so the inserted vertex is the apex opposite it; a -- generic one already names the direction to test. Both -- want the same two apexes, so the one read that settles -- the turn is also the one that supplies the quadrilateral. let !turned = isStar && rawApex /= starVertex if turned && rawTwinApex /= starVertex then loop arena nextTop maxTop flips (probeBump @p candidates) else do let !edge = if turned then rawTwin else rawEdge !twin = if turned then rawEdge else rawTwin !leftFace = if turned then rawTwinFace else rawFace !rightFace = if turned then rawFace else rawTwinFace !edgePrevious = if turned then rawTwinBefore else rawBefore !twinPrevious = if turned then rawBefore else rawTwinBefore !c = if turned then rawTwinApex else rawApex !d = if turned then rawApex else rawTwinApex a <- readOrigin mutable edge b <- readOrigin mutable twin ax <- readPointX mutable a ay <- readPointY mutable a bx <- readPointX mutable b by <- readPointY mutable b cx <- if isStar then pure starX else readPointX mutable c cy <- if isStar then pure starY else readPointY mutable c dx <- readPointX mutable d dy <- readPointY mutable d if not (diagonalFires law ax ay bx by cx cy dx dy) then loop arena nextTop maxTop flips (probeBump @p candidates) else do -- The neighbourhood is complete before the rewrite -- consumes it, and the pushes carry the same edges -- in the same order they always did. edgeNext <- readNext mutable edge twinNext <- readNext mutable twin applyFlip mutable edge twin edgeNext edgePrevious twinNext twinPrevious leftFace rightFace a b c d if isStar then do let !addedTop = nextTop + 2 grown <- growLegalizationArena arena addedTop MUV.unsafeWrite grown nextTop (packIndex (starCandidate twinPrevious)) MUV.unsafeWrite grown (nextTop + 1) (packIndex (starCandidate twinNext)) loop grown addedTop (max maxTop addedTop) (flips + 1) (probeBump @p candidates) else do let !addedTop = nextTop + 4 grown <- growLegalizationArena arena addedTop MUV.unsafeWrite grown nextTop (packIndex (genericCandidate edgeNext)) MUV.unsafeWrite grown (nextTop + 1) (packIndex (genericCandidate edgePrevious)) MUV.unsafeWrite grown (nextTop + 2) (packIndex (genericCandidate twinNext)) MUV.unsafeWrite grown (nextTop + 3) (packIndex (genericCandidate twinPrevious)) loop grown addedTop (max maxTop addedTop) (flips + 1) (probeBump @p candidates) (flips, maxTop, candidates, finalArena) <- loop initialArena seededTop seededTop 0 (probeZero @p) when (MUV.length finalArena /= MUV.length initialArena) $ storeLegalizationArena operation finalArena probeCharge @p operation CounterDiagLegalizationCandidates candidates pure (flips, maxTop) {-# INLINE drainLegalization #-}