tinker: NonLinearArithmatic

src/Theory/LinearArithmatic/Constraint.hs

@@ -47,6 +47,16 @@ data AffineExpr = AffineExpr
   { getConstant :: Rational  
   , getCoeffMap :: LinearExpr }
 
+instance Num AffineExpr where
+  AffineExpr constA coeffsA + AffineExpr constB coeffsB = 
+    AffineExpr (constA+constB) (M.unionWith (+) coeffsA coeffsB)
+  (*) = error "AffineExpr not closed under multiplication"
+  abs = error "TODO"
+  signum = error "TODO"
+  fromInteger n = AffineExpr (fromInteger n) M.empty
+  negate (AffineExpr constant coeffs) = 
+    AffineExpr (negate constant) (M.map negate coeffs)
+
 instance Show AffineExpr where
   show (AffineExpr constant coeff_map) =
     let 

src/Theory/NonLinearRealArithmatic/CAD.hs

@@ -116,8 +116,8 @@ split = \case
       , Closed end   end 
       ] 
 
-isolateRoots :: forall a. (Fractional a, Ord a, Show a) => UnivariatePolynomial a -> [Interval a]
-isolateRoots polynomial = concatMap go $ maybeToList $ cauchyBound polynomial
+isolateRootsIn :: forall a. (Fractional a, Ord a, Show a) => Interval a -> UnivariatePolynomial a -> [Interval a]
+isolateRootsIn initial_interval polynomial = go initial_interval
   where
     sturm_seq = sturmSequence polynomial
 
@@ -129,3 +129,9 @@ isolateRoots polynomial = concatMap go $ maybeToList $ cauchyBound polynomial
       | otherwise = undefined
       where
         root_count = countRootsIn sturm_seq interval
+
+isolateRoots :: forall a. (Fractional a, Ord a, Show a) => UnivariatePolynomial a -> [Interval a]
+isolateRoots polynomial =
+  case cauchyBound polynomial of
+    Nothing       -> []
+    Just interval -> isolateRootsIn interval polynomial

src/Theory/NonLinearRealArithmatic/Constraint.hs

@@ -1,4 +1,3 @@
-
 module Theory.NonLinearRealArithmatic.Constraint 
   ( Constraint
   , ConstraintRelation(..)

src/Theory/NonLinearRealArithmatic/Polynomial.hs

@@ -33,6 +33,7 @@ import Data.Maybe ( maybeToList, fromMaybe )
 import Prelude hiding (map)
 import Control.Arrow ((&&&))
 import Control.Exception (assert)
+import qualified Theory.LinearArithmatic.Constraint as LIA
 
 -- | Map of variables to integer exponents.
 type Monomial = IntMap Int
@@ -256,3 +257,33 @@ class Num a => Assignable a where
   eval :: Assignment a -> Polynomial a -> a
   eval assignment polynomial = sum (evalTerm assignment <$> getTerms polynomial)
 
+type UnivariatePolynomial a = [ (Int, a) ]
+
+{-| 
+  Extracts list of exponent/coefficient pairs, sorted by exponent. 
+  Returns Nothing, if the polynomial is multivariate.
+-}
+toUnivariate :: Polynomial a -> Maybe (UnivariatePolynomial a)
+toUnivariate polynomial =
+  case varsIn polynomial of
+    [ var ] -> Just 
+      $ List.sortOn fst 
+      $ fmap (exponentOf var . getMonomial &&& getCoeff)
+      $ getTerms polynomial
+    zero_or_two_or_more_vars -> Nothing
+
+-- | If possible, converts polynomial to affine expression. Otherwise returns `Nothing`.
+toAffineExpr :: Polynomial Rational -> Maybe LIA.AffineExpr
+toAffineExpr (Polynomial terms) = 
+  let
+    from_term :: Term Rational -> Maybe LIA.AffineExpr
+    from_term (Term coeff monomial) = 
+      case M.toList monomial of
+        -- constant term:
+        [] -> Just $ LIA.AffineExpr coeff M.empty
+        -- single variable with exponent one:
+        [ (1, var) ] -> Just $ LIA.AffineExpr 0 (M.singleton var coeff)
+        -- more than one variable and/or exponents greater one:
+        _non_affine_term -> Nothing
+  in
+    sum <$> traverse from_term terms

src/Theory/NonLinearRealArithmatic/Subtropical.hs

@@ -50,7 +50,7 @@ frame polynomial = undefined -- TODO
       $ L.filter ((/= 0) . getCoeff) (getTerms polynomial)
 
 findDominatingDirection :: (Num a, Ord a) => Polynomial a -> Maybe (Assignment Int)
-findDominatingDirection terms = undefined
+findDominatingDirection terms = error "TODO: compute using Simplex"
   where
     pos_terms = filter ((> 0) . getCoeff) (getTerms terms)
 
@@ -132,7 +132,7 @@ intermediateRoot polynomial neg_sol pos_sol =
         Just [ (0, c), (1, b), (2, a) ] ->
           [ (-b + sqrt (b^2 - 4*a*c)) / (2*a)
           , (-b - sqrt (b^2 - 4*a*c)) / (2*a) ]
-        -- TODO: higher degree polynomial ==> use bisection?
+        -- TODO: higher degree polynomial ==> use CAD
         Just higher_degree_polynomial -> error "TODO: subtropical does not support higher degree polynomials yet"
 
     t_solution :: Assignment a