moonlight-triangulation-1.4.0.4: test/algebra/Moonlight/Triangulation/PowerDiagramSpec.hs
-- | Exact bounded power-cell laws through the public owner.
module Moonlight.Triangulation.PowerDiagramSpec (tests) where
import Control.Monad (unless)
import Data.Foldable (traverse_)
import qualified Data.List as List
import Data.List.NonEmpty (NonEmpty (..))
import qualified Data.List.NonEmpty as NonEmpty
import qualified Data.Map.Strict as Map
import qualified Data.Set as Set
import Moonlight.Triangulation
( ExactPoint
, ExactRational
, ExactVector (..)
, BoundedPowerDiagram
, ConvexPolygon
, Point (..)
, PowerCellDisposition (..)
, PowerDualEdge (..)
, PowerDiagramError (..)
, PowerSite
, PowerWeight
, PowerWeightError (..)
, RegularSiteDisposition (..)
, RegularEdge
, boundedPowerDiagram
, convexPolygon
, convexPolygonPoints
, exactPointCoordinates
, exactPointFromPoint
, exactRayDirection
, exactRayOrigin
, exactSegmentEndpoints
, planarLayerRegions
, powerCellDisposition
, powerCellDispositions
, powerDiagramActiveBoundaries
, powerDiagramCoincidentEquivalentCells
, powerDiagramCoincidentDominatedCells
, powerDiagramEmptyCells
, powerDiagramInputSites
, powerDiagramLowerDimensionalCells
, powerDiagramPlanarLayer
, powerDiagramPublishedCells
, powerDiagramOracleCells
, powerDiagramMaximumCellConstraints
, powerDiagramRegularEdges
, powerDiagramSubmittedSiteConstraints
, powerSite
, powerSiteLabel
, powerSitePosition
, powerSiteWeight
, powerWeight
, powerWeightExact
, regularEdgeDual
, regularEdgeLabels
, regularEdges
, regularFaceDualPoint
, regularFaces
, regularNeighbours
, regularSiteDisposition
, regularTriangulation
, regularTriangulationEdges
, regularTriangulationFaces
, regularTriangulationVisibleSites
, translateExactPoint
)
import Moonlight.Triangulation.Exact
( ExactClipDisposition (..)
, ExactClosedHalfPlane
, exactAffineLine
, exactClipRetainedPolygon
, exactClosedHalfPlane
, exactRetainedPolygon
, exactRetainedPolygonPoints
)
import Support (assertEqual, integerPoint, requireRight)
tests :: IO ()
tests = do
testEqualWeightsAndCommonShift
testThreeSiteExactPartitionCoverage
testRegularTriangleDualRays
testRegularBoundedDualSegment
testRegularCoplanarUpperFacet
testRegularCollinearClassification
testRegularLowerDimensionalAndHiddenSites
testSparsePowerMatchesCompleteOracle
testEqualCoincidentResolution
testDominantCoincidentResolution
testLowerDimensionalCell
testDistinctEmptyCell
testDyadicNearParallelBisectors
testBinary64PrecisionSites
testDuplicateLabels
testNonFiniteWeight
putStrLn "power diagram: ok"
testRegularTriangleDualRays :: IO ()
testRegularTriangleDualRays = do
zero <- admittedWeight 0
firstSite <- admittedSite "first" (Point 0 0) zero
secondSite <- admittedSite "second" (Point 2 0) zero
thirdSite <- admittedSite "third" (Point 0 2) zero
let sites = firstSite :| [secondSite, thirdSite]
(regular, receipt) <-
requireRight "three-site regular topology" (regularTriangulation sites)
assertEqual "three visible regular sites" 3 (regularTriangulationVisibleSites receipt)
assertEqual "one regular face" 1 (length (regularFaces regular))
assertEqual "three regular boundary edges" 3 (length (regularEdges regular))
dualPoint <-
case regularFaces regular of
[face] -> pure (regularFaceDualPoint face)
faces -> fail ("three-site topology expected one face, got " <> show (length faces))
traverse_ (assertBoundaryDualRay sites dualPoint) (regularEdges regular)
traverse_
(\site ->
let label = powerSiteLabel site
in do
assertEqual
("visible regular disposition for " <> label)
(Just RegularSiteVisible)
(regularSiteDisposition label regular)
assertEqual
("two regular neighbours for " <> label)
2
(Set.size (regularNeighbours label regular)))
sites
testRegularBoundedDualSegment :: IO ()
testRegularBoundedDualSegment = do
zero <- admittedWeight 0
sites <-
traverse
(\(label, point) -> admittedSite label point zero)
( ("south-west", Point 0 0)
:| [ ("south-east", Point 4 0)
, ("north-east", Point 3 3)
, ("north-west", Point 0 4)
]
)
(regular, receipt) <-
requireRight "four-site regular topology" (regularTriangulation sites)
assertEqual
"regular planar Euler equation"
1
( regularTriangulationVisibleSites receipt
- regularTriangulationEdges receipt
+ regularTriangulationFaces receipt
)
let faceDuals = Set.fromList (fmap regularFaceDualPoint (regularFaces regular))
bounded =
[ segment
| edge <- regularEdges regular
, BoundedPowerDual segment <- [regularEdgeDual edge]
]
case bounded of
[segment] ->
let (firstEndpoint, secondEndpoint) = exactSegmentEndpoints segment
in unless
(Set.member firstEndpoint faceDuals && Set.member secondEndpoint faceDuals)
(fail "bounded regular dual does not join its two incident face duals")
segments ->
fail ("four-site topology expected one bounded dual, got " <> show (length segments))
testRegularCoplanarUpperFacet :: IO ()
testRegularCoplanarUpperFacet = do
sites <-
traverse
prepareOracleSite
( ("a-bottom", Point 0.25 0.25, -0.875)
:| [ ("b-south-west", Point 0 0, 0)
, ("c-south-east", Point 1 0, 1)
, ("d-center", Point 0.5 0.5, 0.5)
, ("e-north-east", Point 1 1, 2)
, ("f-north-west", Point 0 1, 1)
]
)
(regular, receipt) <-
requireRight "coplanar upper regular facet" (regularTriangulation sites)
traverse_
(\label ->
assertEqual
("upper-facet vertex remains visible: " <> label)
(Just RegularSiteVisible)
(regularSiteDisposition label regular))
["b-south-west", "c-south-east", "e-north-east", "f-north-west"]
assertEqual
"upper-facet interior generator remains lower-dimensional"
(Just RegularSiteLowerDimensional)
(regularSiteDisposition "d-center" regular)
assertEqual
"strictly lower lifted generator remains hidden"
(Just RegularSiteHidden)
(regularSiteDisposition "a-bottom" regular)
assertEqual "coplanar upper facet has four visible vertices" 4 (regularTriangulationVisibleSites receipt)
assertEqual "coplanar upper facet receives one deterministic diagonal" 2 (regularTriangulationFaces receipt)
testRegularLowerDimensionalAndHiddenSites :: IO ()
testRegularLowerDimensionalAndHiddenSites = do
domain <- squareDomain
zero <- admittedWeight 0
lowerWeight <- admittedWeight (-1.5)
hiddenWeight <- admittedWeight (-2)
firstSite <- admittedSite "first" (Point 0 0) zero
secondSite <- admittedSite "second" (Point 2 0) zero
thirdSite <- admittedSite "third" (Point 0 2) zero
lowerSite <- admittedSite "center" (Point 0.5 0.5) lowerWeight
hiddenSite <- admittedSite "center" (Point 0.5 0.5) hiddenWeight
let lowerSites = firstSite :| [secondSite, thirdSite, lowerSite]
hiddenSites = firstSite :| [secondSite, thirdSite, hiddenSite]
(lowerRegular, _) <-
requireRight "coplanar regular topology" (regularTriangulation lowerSites)
assertEqual
"coplanar interior generator remains lower-dimensional"
(Just RegularSiteLowerDimensional)
(regularSiteDisposition "center" lowerRegular)
(lowerDiagram, lowerReceipt) <-
requireRight "bounded lower-dimensional power cell" (boundedPowerDiagram domain lowerSites)
case powerCellDisposition "center" lowerDiagram of
Just (LowerDimensionalPowerCell _) -> pure ()
other -> fail ("expected lower-dimensional center cell, got " <> dispositionTag other)
assertEqual "one lower-dimensional oracle cell" 1 (powerDiagramOracleCells lowerReceipt)
(hiddenRegular, _) <-
requireRight "hidden regular topology" (regularTriangulation hiddenSites)
assertEqual
"strictly interior lifted generator is hidden"
(Just RegularSiteHidden)
(regularSiteDisposition "center" hiddenRegular)
(hiddenDiagram, hiddenReceipt) <-
requireRight "bounded hidden power cell" (boundedPowerDiagram domain hiddenSites)
assertEqual "hidden generator has empty cell" (Just EmptyPowerCell) (powerCellDisposition "center" hiddenDiagram)
assertEqual "hidden generator needs no HPI oracle" 0 (powerDiagramOracleCells hiddenReceipt)
unless
( powerDiagramMaximumCellConstraints hiddenReceipt <= powerDiagramRegularEdges hiddenReceipt
)
(fail "power construction retained more axes than the regular graph")
testRegularCollinearClassification :: IO ()
testRegularCollinearClassification = do
sites <- traverse prepareCollinearSite (0 :| [1 .. 8])
(regular, _) <-
requireRight "collinear regular topology" (regularTriangulation sites)
traverse_
(\index ->
assertEqual
("collinear disposition for " <> show index)
(Just (if even index then RegularSiteVisible else RegularSiteHidden))
(regularSiteDisposition (show index) regular))
([0 .. 8] :: [Int])
where
prepareCollinearSite :: Int -> IO (PowerSite String)
prepareCollinearSite index = do
weight <- admittedWeight (if odd index then -(2 / 256) else 0)
admittedSite (show index) (Point (fromIntegral index / 16) 0) weight
testSparsePowerMatchesCompleteOracle :: IO ()
testSparsePowerMatchesCompleteOracle = do
domain <- squareDomain
sites <-
traverse
prepareOracleSite
( ("a", Point 1.0 1.0, 0.25)
:| [ ("b", Point 3.0 0.8, -0.125)
, ("c", Point 5.2 1.4, 0.375)
, ("d", Point 8.5 1.0, 0.0)
, ("e", Point 1.7 3.8, -0.25)
, ("f", Point 4.1 4.4, 0.125)
, ("g", Point 7.6 3.5, -0.375)
, ("h", Point 9.1 5.4, 0.25)
, ("i", Point 1.0 7.8, 0.0)
, ("j", Point 3.6 8.9, 0.5)
, ("k", Point 6.4 7.5, -0.125)
, ("l", Point 8.8 9.0, 0.375)
]
)
(diagram, receipt) <-
requireRight "sparse regular-neighbour power diagram" (boundedPowerDiagram domain sites)
traverse_
(\site -> do
oracle <- completeOraclePowerCell domain sites site
assertEqual
("regular-neighbour cell equals complete HPI oracle for " <> powerSiteLabel site)
(Just oracle)
(powerCellDisposition (powerSiteLabel site) diagram))
sites
let siteCount = NonEmpty.length sites
completeConstraintCount = siteCount * (siteCount - 1)
unless (powerDiagramSubmittedSiteConstraints receipt < completeConstraintCount) $
fail "regular-neighbour construction did not eliminate the complete pairwise cell schedule"
unless (powerDiagramMaximumCellConstraints receipt <= powerDiagramRegularEdges receipt) $
fail "regular-neighbour construction retained a non-topological radical axis"
prepareOracleSite
:: (String, Point, Double)
-> IO (PowerSite String)
prepareOracleSite (label, point, weightValue) = do
weight <- admittedWeight weightValue
admittedSite label point weight
data OracleScore = OracleScore
{ oracleXCoefficient :: !ExactRational
, oracleYCoefficient :: !ExactRational
, oracleConstant :: !ExactRational
}
oracleScore :: PowerSite label -> IO OracleScore
oracleScore site = do
point <- requireRight "oracle source point" (exactPointFromPoint (powerSitePosition site))
let (coordinateX, coordinateY) = exactPointCoordinates point
pure
OracleScore
{ oracleXCoefficient = 2 * coordinateX
, oracleYCoefficient = 2 * coordinateY
, oracleConstant =
powerWeightExact (powerSiteWeight site)
- coordinateX * coordinateX
- coordinateY * coordinateY
}
completeOraclePowerCell
:: ConvexPolygon
-> NonEmpty (PowerSite String)
-> PowerSite String
-> IO (PowerCellDisposition String)
completeOraclePowerCell domain sites owner = do
retainedDomain <-
requireRight "oracle retained domain" (exactRetainedPolygon (convexPolygonPoints domain))
ownerScore <- oracleScore owner
halfPlanes <-
traverse
(oracleWinningHalfPlane ownerScore)
( filter
((/= powerSiteLabel owner) . powerSiteLabel)
(NonEmpty.toList sites)
)
(disposition, _) <-
requireRight
("complete HPI oracle for " <> powerSiteLabel owner)
(exactClipRetainedPolygon retainedDomain halfPlanes)
case disposition of
ExactClipFullDimensional retained ->
PublishedPowerCell
<$> requireRight
"oracle convex publication"
(convexPolygon (exactRetainedPolygonPoints retained))
ExactClipLowerDimensional points ->
pure (LowerDimensionalPowerCell points)
ExactClipEmpty -> pure EmptyPowerCell
oracleWinningHalfPlane
:: OracleScore
-> PowerSite label
-> IO ExactClosedHalfPlane
oracleWinningHalfPlane owner competitor = do
competitorScore <- oracleScore competitor
line <-
requireRight
"oracle radical axis"
( exactAffineLine
(oracleXCoefficient owner - oracleXCoefficient competitorScore)
(oracleYCoefficient owner - oracleYCoefficient competitorScore)
(oracleConstant owner - oracleConstant competitorScore)
)
pure (exactClosedHalfPlane line)
assertBoundaryDualRay
:: NonEmpty (PowerSite String)
-> ExactPoint
-> RegularEdge String
-> IO ()
assertBoundaryDualRay sites expectedOrigin edge =
case regularEdgeDual edge of
UnboundedPowerDual ray -> do
assertEqual "regular ray starts at incident face dual" expectedOrigin (exactRayOrigin ray)
let ExactVector directionX directionY = exactRayDirection ray
unless (directionX /= 0 || directionY /= 0) $
fail ("regular edge has a zero dual-ray direction: " <> show (regularEdgeLabels edge))
let sample = translateExactPoint expectedOrigin (exactRayDirection ray)
(firstLabel, secondLabel) = regularEdgeLabels edge
firstSite <- requireSite firstLabel sites
secondSite <- requireSite secondLabel sites
firstDistance <- powerDistance sample firstSite
secondDistance <- powerDistance sample secondSite
assertEqual "regular ray remains on its radical axis" firstDistance secondDistance
traverse_
(\competitor -> do
competitorDistance <- powerDistance sample competitor
unless (firstDistance <= competitorDistance) $
fail ("regular ray points outside the common winning cone: " <> show (firstLabel, secondLabel)))
sites
other -> fail ("regular triangle boundary expected a ray, got " <> show other)
requireSite :: Eq label => label -> NonEmpty (PowerSite label) -> IO (PowerSite label)
requireSite label sites =
case List.find ((== label) . powerSiteLabel) (NonEmpty.toList sites) of
Nothing -> fail "regular topology references a missing source site"
Just site -> pure site
testEqualWeightsAndCommonShift :: IO ()
testEqualWeightsAndCommonShift = do
domain <- squareDomain
zero <- admittedWeight 0
shifted <- admittedWeight 7
zeroSites <-
( :| )
<$> admittedSite "left" (Point 2 5) zero
<*> traverse (uncurry3 admittedSite) [("right", Point 8 5, zero)]
shiftedSites <-
( :| )
<$> admittedSite "left" (Point 2 5) shifted
<*> traverse (uncurry3 admittedSite) [("right", Point 8 5, shifted)]
(diagram, receipt) <- requireRight "equal-weight power diagram" (boundedPowerDiagram domain zeroSites)
(shiftedDiagram, _) <- requireRight "common-shift power diagram" (boundedPowerDiagram domain shiftedSites)
assertEqual "common additive power shift" (powerCellDispositions diagram) (powerCellDispositions shiftedDiagram)
assertPermutationInvariant "equal-weight power diagram" domain zeroSites diagram
assertPermutationInvariant "common-shift power diagram" domain shiftedSites shiftedDiagram
expectedLeft <-
requireRight
"expected left power cell"
(convexPolygon (integerPoint 0 0 :| [integerPoint 5 0, integerPoint 5 10, integerPoint 0 10]))
expectedRight <-
requireRight
"expected right power cell"
(convexPolygon (integerPoint 5 0 :| [integerPoint 10 0, integerPoint 10 10, integerPoint 5 10]))
assertEqual "equal-weight left bisector cell" (Just (PublishedPowerCell expectedLeft)) (powerCellDisposition "left" diagram)
assertEqual "equal-weight right bisector cell" (Just (PublishedPowerCell expectedRight)) (powerCellDisposition "right" diagram)
assertEqual "receipt input sites" 2 (powerDiagramInputSites receipt)
assertEqual
"one submitted site constraint per cell"
2
(powerDiagramSubmittedSiteConstraints receipt)
assertEqual
"parallel domain and site boundaries coalesce"
8
(powerDiagramActiveBoundaries receipt)
assertEqual "receipt published cells" 2 (powerDiagramPublishedCells receipt)
assertEqual "receipt empty cells" 0 (powerDiagramEmptyCells receipt)
traverse_ (assertPublishedVerticesWin zeroSites) (powerCellDispositions diagram)
layer <- requireRight "power diagram planar layer" (powerDiagramPlanarLayer "outside" diagram)
assertEqual "derived layer contains both cells" ["left", "right"] (Map.keys (planarLayerRegions layer))
assertOutsideCollision diagram
testThreeSiteExactPartitionCoverage :: IO ()
testThreeSiteExactPartitionCoverage = do
domain <- squareDomain
weight <- admittedWeight 0
left <- admittedSite "left" (Point 2 5) weight
middle <- admittedSite "middle" (Point 5 5) weight
right <- admittedSite "right" (Point 8 5) weight
let sites = left :| [middle, right]
(diagram, receipt) <-
requireRight "three-site exact partition" (boundedPowerDiagram domain sites)
assertEqual "three-site partition publishes every cell" 3 (powerDiagramPublishedCells receipt)
assertEqual
"three-site exact cell areas cover the square"
200
(publishedTwiceAreaSum diagram)
traverse_
(\permutation -> do
(permuted, _) <-
requireRight
"permuted three-site exact partition"
(boundedPowerDiagram domain permutation)
assertEqual
"three-site partition is permutation invariant"
(powerCellDispositions diagram)
(powerCellDispositions permuted)
assertEqual
"permuted exact cell areas cover the square"
200
(publishedTwiceAreaSum permuted))
(nonEmptyPermutations sites)
nonEmptyPermutations :: NonEmpty value -> [NonEmpty value]
nonEmptyPermutations =
foldMap (maybe [] pure . NonEmpty.nonEmpty)
. List.permutations
. NonEmpty.toList
assertPermutationInvariant
:: String
-> ConvexPolygon
-> NonEmpty (PowerSite String)
-> BoundedPowerDiagram String
-> IO ()
assertPermutationInvariant label domain sites expected =
traverse_
(\permutation -> do
(permuted, _) <-
requireRight
(label <> " permutation")
(boundedPowerDiagram domain permutation)
assertEqual
(label <> " is permutation invariant")
(powerCellDispositions expected)
(powerCellDispositions permuted))
(nonEmptyPermutations sites)
publishedTwiceAreaSum :: BoundedPowerDiagram label -> ExactRational
publishedTwiceAreaSum =
List.foldl'
(\total (_, disposition) -> case disposition of
PublishedPowerCell polygon -> total + exactPolygonTwiceArea (convexPolygonPoints polygon)
LowerDimensionalPowerCell _ -> total
EmptyPowerCell -> total
CoincidentEquivalentTo _ -> total
CoincidentDominatedBy _ -> total)
0
. powerCellDispositions
exactPolygonTwiceArea :: NonEmpty ExactPoint -> ExactRational
exactPolygonTwiceArea points =
abs
( List.foldl'
(\twiceArea (firstPoint, secondPoint) ->
let (firstX, firstY) = exactPointCoordinates firstPoint
(secondX, secondY) = exactPointCoordinates secondPoint
in twiceArea + firstX * secondY - firstY * secondX)
0
(cyclicPairs points)
)
cyclicPairs :: NonEmpty value -> [(value, value)]
cyclicPairs (firstValue :| remainingValues) =
zip
(firstValue : remainingValues)
(remainingValues <> [firstValue])
testEqualCoincidentResolution :: IO ()
testEqualCoincidentResolution = do
domain <- squareDomain
weight <- admittedWeight 3
a <- admittedSite "a" (Point 4 4) weight
b <- admittedSite "b" (Point 4 4) weight
(diagram, receipt) <- requireRight "equivalent coincident power sites" (boundedPowerDiagram domain (b :| [a]))
assertPublished "coincident canonical representative" (powerCellDisposition "a" diagram)
assertEqual "equal-function coincident alias" (Just (CoincidentEquivalentTo "a")) (powerCellDisposition "b" diagram)
assertEqual "equivalent coincidence receipt" 1 (powerDiagramCoincidentEquivalentCells receipt)
assertEqual "equivalent coincidence is not dominance" 0 (powerDiagramCoincidentDominatedCells receipt)
(permuted, _) <- requireRight "permuted equivalent sites" (boundedPowerDiagram domain (a :| [b]))
assertEqual "equivalent resolution permutation invariance" (powerCellDispositions diagram) (powerCellDispositions permuted)
testDominantCoincidentResolution :: IO ()
testDominantCoincidentResolution = do
domain <- squareDomain
high <- admittedWeight 3
low <- admittedWeight 2
winner <- admittedSite "winner" (Point 4 4) high
dominated <- admittedSite "dominated" (Point 4 4) low
(diagram, receipt) <- requireRight "dominant coincident power sites" (boundedPowerDiagram domain (dominated :| [winner]))
assertPublished "dominant coincident representative" (powerCellDisposition "winner" diagram)
assertEqual "lower-weight coincident site" (Just (CoincidentDominatedBy "winner")) (powerCellDisposition "dominated" diagram)
assertEqual "dominant coincidence receipt" 1 (powerDiagramCoincidentDominatedCells receipt)
assertEqual "dominance is not equivalence" 0 (powerDiagramCoincidentEquivalentCells receipt)
assertPermutationInvariant "dominant coincident power sites" domain (dominated :| [winner]) diagram
testLowerDimensionalCell :: IO ()
testLowerDimensionalCell = do
domain <- squareDomain
zero <- admittedWeight 0
suppressed <- admittedWeight (-1)
left <- admittedSite "left" (Point 0 5) zero
right <- admittedSite "right" (Point 2 5) zero
middle <- admittedSite "middle" (Point 1 5) suppressed
let sites = left :| [right, middle]
(diagram, receipt) <- requireRight "one-dimensional bounded power cell" (boundedPowerDiagram domain sites)
case powerCellDisposition "middle" diagram of
Just (LowerDimensionalPowerCell points) ->
assertEqual
"middle cell is the exact x=1 segment"
(List.sort [integerPoint 1 0, integerPoint 1 10])
(List.sort (NonEmpty.toList points))
other -> fail ("one-dimensional power cell: expected retained segment, got " <> dispositionTag other)
assertEqual "lower-dimensional receipt" 1 (powerDiagramLowerDimensionalCells receipt)
assertEqual "lower-dimensional cell is not empty" 0 (powerDiagramEmptyCells receipt)
assertEqual "two full-dimensional neighbours" 2 (powerDiagramPublishedCells receipt)
traverse_ (assertPublishedVerticesWin sites) (powerCellDispositions diagram)
layer <- requireRight "lower-dimensional derived layer" (powerDiagramPlanarLayer "outside" diagram)
assertEqual "derived layer omits the one-dimensional cell" ["left", "right"] (Map.keys (planarLayerRegions layer))
assertPermutationInvariant "one-dimensional bounded power cell" domain sites diagram
testDistinctEmptyCell :: IO ()
testDistinctEmptyCell = do
domain <- squareDomain
ordinary <- admittedWeight 0
suppressed <- admittedWeight (-1000)
winner <- admittedSite "winner" (Point 0 0) ordinary
hidden <- admittedSite "hidden" (Point 5 5) suppressed
(diagram, receipt) <- requireRight "distinct empty power cell" (boundedPowerDiagram domain (winner :| [hidden]))
assertPublished "dominant distinct site" (powerCellDisposition "winner" diagram)
assertEqual "distinct site may have empty bounded cell" (Just EmptyPowerCell) (powerCellDisposition "hidden" diagram)
assertEqual "empty disposition receipt" 1 (powerDiagramEmptyCells receipt)
assertPermutationInvariant "distinct empty power cell" domain (winner :| [hidden]) diagram
testDyadicNearParallelBisectors :: IO ()
testDyadicNearParallelBisectors = do
domain <- squareDomain
weight <- admittedWeight 0
let epsilon = 2 ** (-20) :: Double
origin <- admittedSite "origin" (Point 0 0) weight
horizontal <- admittedSite "horizontal" (Point 2 0) weight
tilted <- admittedSite "tilted" (Point 2 epsilon) weight
let sites = origin :| [horizontal, tilted]
(diagram, receipt) <-
requireRight
"dyadic near-parallel power bisectors"
(boundedPowerDiagram domain sites)
assertEqual "near-parallel cells retain full dimension" 3 (powerDiagramPublishedCells receipt)
traverse_ (assertPublishedVerticesWin sites) (powerCellDispositions diagram)
assertPermutationInvariant "dyadic near-parallel power bisectors" domain sites diagram
testBinary64PrecisionSites :: IO ()
testBinary64PrecisionSites = do
domain <- squareDomain
firstWeight <- admittedWeight 0.2
secondWeight <- admittedWeight (-0.3)
thirdWeight <- admittedWeight 0.7
firstSite <- admittedSite "first" (Point 0.1 0.3) firstWeight
secondSite <- admittedSite "second" (Point 9.7 0.2) secondWeight
thirdSite <- admittedSite "third" (Point 4.9 9.6) thirdWeight
let sites = firstSite :| [secondSite, thirdSite]
(diagram, receipt) <-
requireRight
"exact binary64 power sites"
(boundedPowerDiagram domain sites)
assertEqual "binary64 fixture publishes every cell" 3 (powerDiagramPublishedCells receipt)
traverse_ (assertPublishedVerticesWin sites) (powerCellDispositions diagram)
assertPermutationInvariant "exact binary64 power sites" domain sites diagram
testDuplicateLabels :: IO ()
testDuplicateLabels = do
domain <- squareDomain
weight <- admittedWeight 0
firstSite <- admittedSite "duplicate" (Point 1 1) weight
secondSite <- admittedSite "duplicate" (Point 9 9) weight
case boundedPowerDiagram domain (firstSite :| [secondSite]) of
Left (DuplicatePowerSiteLabel "duplicate") -> pure ()
other -> fail ("duplicate power-site label: expected typed refusal, got " <> show other)
testNonFiniteWeight :: IO ()
testNonFiniteWeight =
case powerWeight (0 / 0) of
Left (PowerWeightNonFinite _) -> pure ()
other -> fail ("non-finite power weight: expected typed refusal, got " <> show other)
squareDomain :: IO ConvexPolygon
squareDomain =
requireRight
"square clipping domain"
(convexPolygon (integerPoint 0 0 :| [integerPoint 10 0, integerPoint 10 10, integerPoint 0 10]))
admittedWeight :: Double -> IO PowerWeight
admittedWeight value =
requireRight
"finite power weight"
(powerWeight value)
admittedSite :: String -> Point -> PowerWeight -> IO (PowerSite String)
admittedSite label point weight = requireRight "admitted power site" (powerSite label point weight)
uncurry3 :: (a -> b -> c -> result) -> (a, b, c) -> result
uncurry3 function (firstValue, secondValue, thirdValue) = function firstValue secondValue thirdValue
assertPublished :: String -> Maybe (PowerCellDisposition label) -> IO ()
assertPublished _ (Just (PublishedPowerCell _)) = pure ()
assertPublished label other = fail (label <> ": expected published cell, got " <> dispositionTag other)
dispositionTag :: Maybe (PowerCellDisposition label) -> String
dispositionTag Nothing = "missing"
dispositionTag (Just (PublishedPowerCell _)) = "published"
dispositionTag (Just (LowerDimensionalPowerCell _)) = "lower-dimensional"
dispositionTag (Just EmptyPowerCell) = "empty"
dispositionTag (Just (CoincidentEquivalentTo _)) = "coincident-equivalent"
dispositionTag (Just (CoincidentDominatedBy _)) = "coincident-dominated"
assertOutsideCollision :: BoundedPowerDiagram String -> IO ()
assertOutsideCollision diagram =
case powerDiagramPlanarLayer "left" diagram of
Left (PowerDiagramOutsideLabelCollides "left") -> pure ()
other -> fail ("power layer outside-label collision: expected refusal, got " <> show other)
assertPublishedVerticesWin
:: NonEmpty (PowerSite String)
-> (String, PowerCellDisposition String)
-> IO ()
assertPublishedVerticesWin sites (ownerLabel, disposition) =
case disposition of
PublishedPowerCell polygon ->
assertCellPointsWin ownerLabel sites (convexPolygonPoints polygon)
LowerDimensionalPowerCell points ->
assertCellPointsWin ownerLabel sites points
EmptyPowerCell -> pure ()
CoincidentEquivalentTo _ -> pure ()
CoincidentDominatedBy _ -> pure ()
assertCellPointsWin
:: String
-> NonEmpty (PowerSite String)
-> NonEmpty ExactPoint
-> IO ()
assertCellPointsWin ownerLabel sites points =
case lookupOwner ownerLabel (NonEmpty.toList sites) of
Nothing -> fail ("power-cell owner missing: " <> ownerLabel)
Just owner ->
traverse_
(\point -> traverse_ (assertOwnerWinsAt point owner) sites)
points
lookupOwner :: Eq label => label -> [PowerSite label] -> Maybe (PowerSite label)
lookupOwner label = List.find ((== label) . powerSiteLabel)
assertOwnerWinsAt
:: ExactPoint
-> PowerSite String
-> PowerSite String
-> IO ()
assertOwnerWinsAt point owner competitor = do
ownerValue <- powerDistance point owner
competitorValue <- powerDistance point competitor
unless (ownerValue <= competitorValue) $
fail
( "power-cell vertex violates source inequality: "
<> show (powerSiteLabel owner, powerSiteLabel competitor, ownerValue, competitorValue)
)
powerDistance :: ExactPoint -> PowerSite label -> IO ExactRational
powerDistance point site = do
exactSite <- requireRight "site exact position" (exactPointFromPoint (powerSitePosition site))
let (x, y) = exactPointCoordinates point
(siteX, siteY) = exactPointCoordinates exactSite
deltaX = x - siteX
deltaY = y - siteY
pure (deltaX * deltaX + deltaY * deltaY - powerWeightExact (powerSiteWeight site))