diff --git a/brush-strokes/src/lib/Math/Bezier/Stroke.hs b/brush-strokes/src/lib/Math/Bezier/Stroke.hs index e06591d..dc9b90a 100644 --- a/brush-strokes/src/lib/Math/Bezier/Stroke.hs +++ b/brush-strokes/src/lib/Math/Bezier/Stroke.hs @@ -1234,10 +1234,10 @@ findCuspsIn opts boxStrokeData initBoxes = , _D22_dy = T ( 𝕀 ( ℝ1 ee_s_lo ) ( ℝ1 ee_s_hi ) ) } } = ( boxStrokeData t `Seq.index` i ) s -- Ξ» = βˆ‚E/βˆ‚t / βˆ‚E/βˆ‚s - Ξ»1 = 𝕀 ee_t_lo ee_t_hi `extendedDivide` 𝕀 ee_s_lo ee_s_hi + Ξ»1 = 𝕀 ee_t_lo ee_t_hi ⊘ 𝕀 ee_s_lo ee_s_hi -- Ξ» = u / v - Ξ»2 = 𝕀 ux_lo ux_hi `extendedDivide` 𝕀 vx_lo vx_hi - Ξ»3 = 𝕀 uy_lo uy_hi `extendedDivide` 𝕀 vy_lo vy_hi + Ξ»2 = 𝕀 ux_lo ux_hi ⊘ 𝕀 vx_lo vx_hi + Ξ»3 = 𝕀 uy_lo uy_hi ⊘ 𝕀 vy_lo vy_hi Ξ» = [ 𝕀 ( recip -0 ) ( recip 0 ) ] `intersectMany` Ξ»1 `intersectMany` Ξ»2 diff --git a/brush-strokes/src/lib/Math/Bezier/Stroke/EnvelopeEquation.hs b/brush-strokes/src/lib/Math/Bezier/Stroke/EnvelopeEquation.hs index 7c5482a..933f828 100644 --- a/brush-strokes/src/lib/Math/Bezier/Stroke/EnvelopeEquation.hs +++ b/brush-strokes/src/lib/Math/Bezier/Stroke/EnvelopeEquation.hs @@ -451,9 +451,9 @@ evaluateCubic bez t = let inf_bez = fmap inf bez sup_bez = fmap sup bez mins = fmap (Cubic.bezier @( T Double ) inf_bez) - $ inf t :| ( sup t : filter ( `inside` t ) ( Cubic.extrema inf_bez ) ) + $ inf t :| ( sup t : filter ( ∈ t ) ( Cubic.extrema inf_bez ) ) maxs = fmap (Cubic.bezier @( T Double ) sup_bez) - $ inf t :| ( sup t : filter ( `inside` t ) ( Cubic.extrema sup_bez ) ) + $ inf t :| ( sup t : filter ( ∈ t ) ( Cubic.extrema sup_bez ) ) in 𝕀 ( minimum mins ) ( maximum maxs ) -- | Evaluate a quadratic BΓ©zier curve, when both the coefficients and the @@ -464,7 +464,7 @@ evaluateQuadratic bez t = let inf_bez = fmap inf bez sup_bez = fmap sup bez mins = fmap (Quadratic.bezier @( T Double ) inf_bez) - $ inf t :| ( sup t : filter ( `inside` t ) ( Quadratic.extrema inf_bez ) ) + $ inf t :| ( sup t : filter ( ∈ t ) ( Quadratic.extrema inf_bez ) ) maxs = fmap (Quadratic.bezier @( T Double ) sup_bez) - $ inf t :| ( sup t : filter ( `inside` t ) ( Quadratic.extrema sup_bez ) ) + $ inf t :| ( sup t : filter ( ∈ t ) ( Quadratic.extrema sup_bez ) ) in 𝕀 ( minimum mins ) ( maximum maxs ) diff --git a/brush-strokes/src/lib/Math/Interval.hs b/brush-strokes/src/lib/Math/Interval.hs index cb80187..32575ec 100644 --- a/brush-strokes/src/lib/Math/Interval.hs +++ b/brush-strokes/src/lib/Math/Interval.hs @@ -11,10 +11,9 @@ module Math.Interval , 𝕀ℝ , isCanonical , singleton, nonDecreasing - , inside - , aabb - , extendedDivide, extendedRecip - , intersect, bisect + , (∈), (βŠ–), (∩), (⊘) + , extendedRecip + , aabb, bisect ) where @@ -72,9 +71,10 @@ nonDecreasing :: ( a -> b ) -> 𝕀 a -> 𝕀 b nonDecreasing f ( 𝕀 lo hi ) = 𝕀 ( f lo ) ( f hi ) -- | Does the given value lie inside the specified interval? -inside :: Ord a => a -> 𝕀 a -> Bool -inside x ( 𝕀 lo hi ) = x >= lo && x <= hi -{-# INLINEABLE inside #-} +(∈) :: Ord a => a -> 𝕀 a -> Bool +x ∈ 𝕀 lo hi = x >= lo && x <= hi +{-# INLINEABLE (∈) #-} +infix 4 ∈ -- | Is this interval canonical, i.e. it consists of either 1 or 2 floating point -- values only? @@ -91,8 +91,8 @@ midpoint ( 𝕀 x_inf x_sup ) = 0.5 * ( x_inf + x_sup ) -- -- Returns whether the first interval is a strict subset of the second interval -- (or the intersection is a single point). -intersect :: 𝕀 Double -> 𝕀 Double -> [ ( 𝕀 Double, Bool ) ] -intersect ( 𝕀 lo1 hi1 ) ( 𝕀 lo2 hi2 ) +(∩) :: 𝕀 Double -> 𝕀 Double -> [ ( 𝕀 Double, Bool ) ] +𝕀 lo1 hi1 ∩ 𝕀 lo2 hi2 | lo > hi = [ ] | otherwise @@ -100,6 +100,7 @@ intersect ( 𝕀 lo1 hi1 ) ( 𝕀 lo2 hi2 ) where lo = max lo1 lo2 hi = min hi1 hi2 +infix 3 ∩ -- | Bisect an interval. -- @@ -113,6 +114,15 @@ bisect x@( 𝕀 x_inf x_sup ) = 𝕀 x_inf x_mid NE.:| [ 𝕀 x_mid x_sup ] where x_mid = midpoint x +infixl 6 βŠ– +(βŠ–) :: ( Ring a, Ord a ) => 𝕀 a -> 𝕀 a -> 𝕀 a +(βŠ–) a@( 𝕀 lo1 hi1 ) b@( 𝕀 lo2 hi2 ) + | width a >= width b + = 𝕀 ( lo1 - lo2 ) ( hi1 - hi2 ) + | otherwise + = 𝕀 ( hi1 - hi2 ) ( lo1 - lo2 ) +{-# INLINEABLE (βŠ–) #-} + deriving via ViaAbelianGroup ( T ( 𝕀 Double ) ) instance Semigroup ( T ( 𝕀 Double ) ) deriving via ViaAbelianGroup ( T ( 𝕀 Double ) ) @@ -128,8 +138,10 @@ instance Torsor ( T ( 𝕀 Double ) ) ( 𝕀 Double ) where ------------------------------------------------------------------------------- -- Extended division -extendedDivide :: 𝕀 Double -> 𝕀 Double -> [ 𝕀 Double ] -extendedDivide x y = nub $ map ( x * ) ( extendedRecip y ) +-- | Extended division. +(⊘) :: 𝕀 Double -> 𝕀 Double -> [ 𝕀 Double ] +x ⊘ y = nub $ map ( x * ) ( extendedRecip y ) +infixl 7 ⊘ extendedRecip :: 𝕀 Double -> [ 𝕀 Double ] extendedRecip x@( 𝕀 lo hi ) diff --git a/brush-strokes/src/lib/Math/Root/Isolation.hs b/brush-strokes/src/lib/Math/Root/Isolation.hs index b21d15b..2decabc 100644 --- a/brush-strokes/src/lib/Math/Root/Isolation.hs +++ b/brush-strokes/src/lib/Math/Root/Isolation.hs @@ -192,7 +192,7 @@ isolateRootsIn ( RootIsolationOptions { rootIsolationAlgorithms } ) [ RootIsolationTree ( Box n ) ] go history cand | -- Check the range of the equations contains zero. - not $ ( unT ( origin @Double ) `inside` iRange ) + not $ ( unT ( origin @Double ) ∈ iRange ) -- Box doesn't contain a solution: discard it. = return [] | otherwise diff --git a/brush-strokes/src/lib/Math/Root/Isolation/Bisection.hs b/brush-strokes/src/lib/Math/Root/Isolation/Bisection.hs index f51927d..52f95d0 100644 --- a/brush-strokes/src/lib/Math/Root/Isolation/Bisection.hs +++ b/brush-strokes/src/lib/Math/Root/Isolation/Bisection.hs @@ -106,7 +106,7 @@ defaultBisectionOptions minWidth _Ξ΅_eq box = -- box(0)-consistency let iRange' :: Box d iRange' = eqs box' `monIndex` zeroMonomial - in unT ( origin @Double ) `inside` iRange' + in unT ( origin @Double ) ∈ iRange' -- box(1)-consistency --let box1Options = Box1Options _Ξ΅_eq ( toList $ universe @n ) ( toList $ universe @d ) @@ -117,7 +117,7 @@ defaultBisectionOptions minWidth _Ξ΅_eq box = -- box'' = makeBox2Consistent _minWidth box2Options eqs box' -- iRange'' :: Box d -- iRange'' = eqs box'' `monIndex` zeroMonomial - --in unT ( origin @Double ) `inside` iRange'' + --in unT ( origin @Double ) ∈ iRange'' , fallbackBisectionCoord = \ _thisRoundHist _prevRoundsHist eqs possibleCoordChoices -> let datPerCoord = diff --git a/brush-strokes/src/lib/Math/Root/Isolation/GaussSeidel.hs b/brush-strokes/src/lib/Math/Root/Isolation/GaussSeidel.hs index aca5411..040e76c 100644 --- a/brush-strokes/src/lib/Math/Root/Isolation/GaussSeidel.hs +++ b/brush-strokes/src/lib/Math/Root/Isolation/GaussSeidel.hs @@ -208,17 +208,17 @@ gaussSeidelStep as b ( T x0 ) = coerce $ x_i = x `index` i a_ii = ( as ! i ) `index` i -- Take a shortcut before performing the division if possible. - if | not $ 0 `inside` ( s - a_ii * x_i ) + if | not $ 0 ∈ ( s - a_ii * x_i ) -- No solutions: don't bother performing a division. -> [ ] - | 0 `inside` s && 0 `inside` a_ii + | 0 ∈ s && 0 ∈ a_ii -- The division would produce [-oo,+oo]: don't do anything. -> [ ( x, False ) ] -- Otherwise, perform the division. | otherwise -> do - x_i'0 <- s `extendedDivide` a_ii - ( x_i', sub_i ) <- x_i'0 `intersect` x_i + x_i'0 <- s ⊘ a_ii + ( x_i', sub_i ) <- x_i'0 ∩ x_i return $ ( set i x_i' x, sub_i && contraction ) {-# INLINEABLE gaussSeidelStep #-} @@ -244,16 +244,16 @@ gaussSeidelStep_Complete as b ( T x0 ) = coerce $ do ( x', subs ) <- fromComponents \ j -> do let x_j = x `index` j a_ij = ( as ! j ) `index` i - s_j = s `ominus` ( a_ij * x_j ) + s_j = s βŠ– ( a_ij * x_j ) -- Shortcut division if possible (see gaussSeidelStep for commentary). - if | not $ 0 `inside` ( s_j - a_ii * x_i ) + if | not $ 0 ∈ ( s_j - a_ii * x_i ) -> [ ] - | 0 `inside` s_j && 0 `inside` a_ij + | 0 ∈ s_j && 0 ∈ a_ij -> [ ( x_j, False ) ] | otherwise -> do - x_j'0 <- s_j `extendedDivide` a_ij - ( x_j', sub_j ) <- x_j'0 `intersect` x_j + x_j'0 <- s_j ⊘ a_ij + ( x_j', sub_j ) <- x_j'0 ∩ x_j return $ ( x_j', sub_j ) return ( x', (||) <$> contractions <*> subs ) return ( x', and subs ) @@ -269,14 +269,6 @@ fromComponents f = do -- TODO: this could be more efficient. {-# INLINEABLE fromComponents #-} -infixl 6 `ominus` -ominus :: 𝕀 Double -> 𝕀 Double -> 𝕀 Double -ominus a@( 𝕀 lo1 hi1 ) b@( 𝕀 lo2 hi2 ) - | width a >= width b - = 𝕀 ( lo1 - lo2 ) ( hi1 - hi2 ) - | otherwise - = 𝕀 ( hi1 - hi2 ) ( lo1 - lo2 ) - -- | The midpoint of a box. boxMidpoint :: Representable Double ( ℝ n ) => 𝕀ℝ n -> ℝ n boxMidpoint box = diff --git a/brush-strokes/src/lib/Math/Root/Isolation/Narrowing.hs b/brush-strokes/src/lib/Math/Root/Isolation/Narrowing.hs index 9ff0dd5..a41ace9 100644 --- a/brush-strokes/src/lib/Math/Root/Isolation/Narrowing.hs +++ b/brush-strokes/src/lib/Math/Root/Isolation/Narrowing.hs @@ -295,9 +295,9 @@ leftNarrow Ξ΅_eq Ξ΅_bis ff' = left_narrow else let x = 𝕀 ( sup x_left ) x_sup ( _f_x, f'_x ) = ff' x - x's = do Ξ΄ <- f_x_left `extendedDivide` f'_x + x's = do Ξ΄ <- f_x_left ⊘ f'_x let x_new = x_left - Ξ΄ - map fst $ x_new `intersect` x + map fst $ x_new ∩ x in if | null x's -> [] @@ -328,9 +328,9 @@ rightNarrow Ξ΅_eq Ξ΅_bis ff' = right_narrow else let x = 𝕀 x_inf ( inf x_right ) ( _f_x, f'_x ) = ff' x - x's = do Ξ΄ <- f_x_right `extendedDivide` f'_x + x's = do Ξ΄ <- f_x_right ⊘ f'_x let x_new = x_right - Ξ΄ - map fst $ x_new `intersect` x + map fst $ x_new ∩ x in if | null x's -> [] @@ -387,7 +387,7 @@ leftShave Ξ΅_eq Ξ΅_bis go Ξ³ x_sup x_left = let ( f_x_left, _f'_x_left ) = ff' x_left in - if 0 `inside` f_x_left + if 0 ∈ f_x_left -- Box-consistency achieved; finish. then [ 𝕀 ( inf x_left ) x_sup ] else @@ -399,7 +399,7 @@ leftShave Ξ΅_eq Ξ΅_bis -- NB: this always uses the initial width (to "avoid asymptotic behaviour" according to the paper) ( f_guess, f'_guess ) = ff' guess x_minus_guess = 𝕀 ( min x_sup ( succFP $ sup guess ) ) x_sup - in if not ( 0 `inside` f_guess ) + in if not ( 0 ∈ f_guess ) then -- We successfully shaved "guess" off; go round again after removing it. -- TODO: here we could go back to the top with a new "w" maybe? @@ -410,9 +410,9 @@ leftShave Ξ΅_eq Ξ΅_bis -- "guesses'" refining where the function can be zero. let guesses' :: [ ( 𝕀 Double ) ] guesses' = do - Ξ΄ <- f_x_left `extendedDivide` f'_guess + Ξ΄ <- f_x_left ⊘ f'_guess let guess' = singleton ( inf guess ) - Ξ΄ - map fst $ guess' `intersect` guess + map fst $ guess' ∩ guess w_guess = width guess w_guesses' | null guesses' @@ -453,30 +453,30 @@ sbc Ξ΅ ff' = go | otherwise = let x_mid = 0.5 * ( x_l + x_r ) ( left_done, left_todo ) - | 0 `inside` ( fst $ ff' ( 𝕀 x_l x_lp ) ) + | 0 ∈ ( fst $ ff' ( 𝕀 x_l x_lp ) ) = ( True, [ ] ) - | not $ 0 `inside` ( fst $ ff' i_l ) + | not $ 0 ∈ ( fst $ ff' i_l ) = ( False, [ ] ) | otherwise = let xls = do let l = 𝕀 x_lp x_lp - Ξ΄ <- fst ( ff' l ) `extendedDivide` snd ( ff' i_l ) - map fst $ ( l - Ξ΄ ) `intersect` i_l + Ξ΄ <- fst ( ff' l ) ⊘ snd ( ff' i_l ) + map fst $ ( l - Ξ΄ ) ∩ i_l in ( False, xls ) where x_lp = min ( succFP x_l ) x_mid i_l = 𝕀 x_lp x_mid ( right_done, right_todo ) - | 0 `inside` ( fst $ ff' ( 𝕀 x_rm x_r ) ) + | 0 ∈ ( fst $ ff' ( 𝕀 x_rm x_r ) ) = ( True, [ ] ) - | not $ 0 `inside` ( fst $ ff' i_r ) + | not $ 0 ∈ ( fst $ ff' i_r ) = ( False, [ ] ) | otherwise = let xrs = do let r = 𝕀 x_rm x_rm - Ξ΄ <- fst ( ff' r ) `extendedDivide` snd ( ff' i_r ) - map fst $ ( r - Ξ΄ ) `intersect` i_r + Ξ΄ <- fst ( ff' r ) ⊘ snd ( ff' i_r ) + map fst $ ( r - Ξ΄ ) ∩ i_r in ( False, xrs ) where x_rm = max ( prevFP x_r ) x_mid i_r = 𝕀 x_mid x_rm