unit test ellipse implicitization

stla · May 2, 2024 · 747a9b0 · 747a9b0
1 parent 9f3e7d1
commit 747a9b0
Show file tree

Hide file tree

Showing 3 changed files with 36 additions and 36 deletions.
diff --git a/scripts/script30.hs b/scripts/script30.hs
@@ -1,7 +1,4 @@
 import Math.Algebra.Hspray
-import qualified Data.Foldable       as DF
-import qualified Data.HashMap.Strict           as HM
-import qualified Data.Sequence           as S
 
 cost = qlone 1
 sint = qlone 2
@@ -15,17 +12,12 @@ n_equations = length equations
 coordinates = [qlone (m + i) | i <- [1 .. n_equations]]
 generators = relations ++ zipWith (^-^) equations coordinates 
 gb = groebnerBasis generators True
-gb' = drop 1 gb
-f :: Exponents -> Bool
-f expnts = S.null expnts || DF.all (0 ==) (S.take n_variables expnts)
+-- gb' = drop 1 gb
 isfree :: QSpray -> Bool
-isfree spray = DF.all f (allExponents spray)
-results = filter isfree gb'
-dropXis = HM.mapKeys 
-            (\(Powers exps n) -> 
-                Powers (S.drop n_variables exps) (n - n_variables))
-results' = map dropXis results
-showResults = map (prettyQSprayXYZ ["a", "b"]) results'
+isfree spray = not $ or (map (involvesVariable spray) [1 .. n_variables]) 
+results = filter isfree gb
+results' = map (dropVariables n_variables) results 
+showResults = map (prettyQSprayXYZ ["a", "b", "x", "y"]) results'
 
 {- sprays = generators
 j0       = length sprays

diff --git a/src/Math/Algebra/Hspray.hs b/src/Math/Algebra/Hspray.hs
@@ -2269,18 +2269,15 @@ psPolynomial n k
 -- (use the function with the same name in the __jackpolynomials__ package 
 -- if you need efficiency)
 isSymmetricSpray :: forall a. (AlgField.C a, Eq a) => Spray a -> Bool
-isSymmetricSpray spray = check1 && check2 
+isSymmetricSpray spray = check
   where
     n = numberOfVariables spray
     indices = [1 .. n]
     gPolys  = [esPolynomial n i ^-^ lone (n + i) | i <- indices]
     gbasis  = groebner0 gPolys
     spray'  = removeConstantTerm spray
     g       = sprayDivisionRemainder spray' gbasis
-    gpowers = HM.keys g
-    check1  = minimum (map nvariables gpowers) > n
-    expnts  = map exponents gpowers
-    check2  = DF.all (DF.all (0 ==)) (map (S.take n) expnts) 
+    check   = not $ any (involvesVariable g) [1 .. n]
 
 -- | Whether a spray can be written as a polynomial of a given list of sprays;
 -- this polynomial is returned if this is true
@@ -2303,7 +2300,7 @@ isPolynomialOf spray sprays
     n   = maximum $ map numberOfVariables sprays
     result
       | nov > n   = (False, Nothing)
-      | otherwise = (checks, poly)
+      | otherwise = (check, poly)
         where
           m            = length sprays
           yPolys       = [lone (n + i) | i <- [1 .. m]]  
@@ -2312,20 +2309,10 @@ isPolynomialOf spray sprays
           constantTerm = getConstantTerm spray
           spray'       = removeConstantTerm spray
           g            = sprayDivisionRemainder spray' gbasis0
-          gpowers      = HM.keys g
-          check1       = minimum (map nvariables gpowers) > n
-          check2       = 
-            DF.all (DF.all (0 ==)) (map (S.take n . exponents) gpowers)
-          checks       = check1 && check2
-          poly         = if checks
-            then Just g''
+          check        = not $ any (involvesVariable g) [1 .. n]
+          poly         = if check
+            then Just (constantTerm +> dropVariables n g)
             else Nothing
-          g' = dropXis g
-          g'' = if constantTerm == AlgAdd.zero 
-            then g' 
-            else addTerm g' (nullPowers, constantTerm)
-          dropXis = HM.mapKeys f
-          f (Powers expnnts nv) = Powers (S.drop n expnnts) (nv - n)
 
 
 -- resultant ------------------------------------------------------------------

diff --git a/tests/Main.hs b/tests/Main.hs
@@ -95,7 +95,8 @@ import           Math.Algebra.Hspray            ( Spray,
                                                   qlone',
                                                   qmonomial,
                                                   isHomogeneousSpray,
-                                                  psPolynomial
+                                                  psPolynomial,
+                                                  prettyQSprayXYZ
                                                 )
 import           MathObj.Matrix                 ( fromRows )
 import qualified MathObj.Matrix                 as MathMatrix
@@ -584,9 +585,29 @@ main = defaultMain $ testGroup
         z = lone 3 :: Spray Rational
       assertEqual "" 
         (isPolynomialOf x [x ^+^ y^*^z, y, z]) 
-        (True, Just $ x ^-^ y^*^z),
-
-    testCase "power sum polynomials" $ do
+        (True, Just $ x ^-^ y^*^z)
+
+    , testCase "Groebner implicitization ellipse" $ do
+      let
+        cost = qlone 1
+        sint = qlone 2
+        n_variables = 2
+        a = qlone 3
+        b = qlone 4
+        equations = [a ^*^ cost, b ^*^ sint]
+        relations = [cost^**^2 ^+^ sint^**^2 ^-^ unitSpray]
+        m = maximum (map numberOfVariables equations)
+        coordinates = [qlone (m + i) | i <- [1 .. length equations]]
+        generators = relations ++ zipWith (^-^) equations coordinates 
+        gb = groebnerBasis generators True
+        isfree :: QSpray -> Bool
+        isfree spray = not $ any (involvesVariable spray) [1 .. n_variables]
+        results = filter isfree gb
+        results' = map (dropVariables n_variables) results 
+        showResults = map (prettyQSprayXYZ ["a", "b", "x", "y"]) results'
+      assertEqual "" showResults ["a^2.b^2 - a^2.y^2 - b^2.x^2"]
+
+    , testCase "power sum polynomials" $ do
       let
         x = lone 1 :: Spray Rational
         y = lone 2 :: Spray Rational