diff --git a/ChangeLog.md b/ChangeLog.md new file mode 100644 --- /dev/null +++ b/ChangeLog.md @@ -0,0 +1,3 @@ +# Changelog for Ritt-Wu + +## Unreleased changes diff --git a/LICENSE b/LICENSE new file mode 100644 --- /dev/null +++ b/LICENSE @@ -0,0 +1,30 @@ +Copyright Author name here (c) 2019 + +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions are met: + + * Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + + * Redistributions in binary form must reproduce the above + copyright notice, this list of conditions and the following + disclaimer in the documentation and/or other materials provided + with the distribution. + + * Neither the name of Author name here nor the names of other + contributors may be used to endorse or promote products derived + from this software without specific prior written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR +A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT +OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, +DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY +THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT +(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. diff --git a/README.md b/README.md new file mode 100644 --- /dev/null +++ b/README.md @@ -0,0 +1,1 @@ +# Ritt-Wu diff --git a/Ritt-Wu.cabal b/Ritt-Wu.cabal new file mode 100644 --- /dev/null +++ b/Ritt-Wu.cabal @@ -0,0 +1,114 @@ +cabal-version: 1.12 + +-- This file has been generated from package.yaml by hpack version 0.31.1. +-- +-- see: https://github.com/sol/hpack +-- +-- hash: 4cc0561c5348e7e793df07bd56a8e30373c5664853c2f505a7cf6a1f9e64f2d0 + +name: Ritt-Wu +version: 0.1.0.0 +synopsis: Parallel implementation of Ritt-Wu's algorithm. +description: Please see the README on GitHub at +homepage: https://github.com/githubuser/Ritt-Wu#readme +category: Algorithm +bug-reports: https://github.com/githubuser/Ritt-Wu/issues +author: Jose Luis Seraquive Cuenca +maintainer: jose.seraquive@gmail.com +copyright: 2019 Jose Seraquive +license: BSD3 +license-file: LICENSE +build-type: Simple +extra-source-files: + README.md + ChangeLog.md + +source-repository head + type: git + location: https://github.com/githubuser/Ritt-Wu + +library + exposed-modules: + Polynomial.Monomial + Polynomial.Terms + Polynomial.Polynomial + Polynomial.Wu + other-modules: + Paths_Ritt_Wu + hs-source-dirs: + src + build-depends: + base >=4.7 && <5 + , algebra -any + --, containers >= 0.6.0.1 + -- , containers >= 0.6.0.1 + , criterion -any + , deepseq >= 1.4.4.0 + , massiv == 0.4.2.0 + , scheduler >= 1.4.2 + , sscript -any + default-language: Haskell2010 + +executable Ritt-Wu-exe + main-is: Main.hs + other-modules: + Paths_Ritt_Wu + hs-source-dirs: + app + ghc-options: -threaded -rtsopts -eventlog "-with-rtsopts=-N -s -h -i0.1" + build-depends: + Ritt-Wu + , algebra + , base >=4.7 && <5 + --, containers >= 0.6.0.1 + --, containers >= 0.6.0.1 + , criterion + , deepseq >= 1.4.4.0 + , massiv == 0.4.2.0 + , scheduler >= 1.4.2 + , sscript + default-language: Haskell2010 + +test-suite Ritt-Wu-test + type: exitcode-stdio-1.0 + main-is: Spec.hs + other-modules: + Paths_Ritt_Wu + hs-source-dirs: + test + ghc-options: -threaded -rtsopts + build-depends: + Ritt-Wu + , base >=4.7 && <5 + --, containers >= 0.6.0.1 + -- ,containers >= 0.6.0.1 + , deepseq >= 1.4.4.0 + , scheduler >= 1.4.2 + , tasty + , tasty-hunit + , tasty-quickcheck + , tasty-smallcheck + default-language: Haskell2010 + +benchmark Ritt-Wu-benchmark + type: exitcode-stdio-1.0 + hs-source-dirs: bench + main-is: Benchmarking.hs + other-modules: + Paths_Ritt_Wu + ghc-options: -threaded "-with-rtsopts=-N -s " + build-depends: + Ritt-Wu + , algebra + , base >=4.7 && <5 + --, containers >= 0.6.0.1 + --, containers >= 0.6.0.1 + , criterion + , deepseq >= 1.4.4.0 + , massiv == 0.4.2.0 + , scheduler >= 1.4.2 + , sscript + + default-language: Haskell2010 + + diff --git a/Setup.hs b/Setup.hs new file mode 100644 --- /dev/null +++ b/Setup.hs @@ -0,0 +1,2 @@ +import Distribution.Simple +main = defaultMain diff --git a/app/Main.hs b/app/Main.hs new file mode 100644 --- /dev/null +++ b/app/Main.hs @@ -0,0 +1,153 @@ +module Main where + +import Polynomial.Polynomial +import Polynomial.Monomial +import Polynomial.Terms +import Polynomial.Wu +-- -- http://symbolicdata.org/XMLResources/IntPS/Behnke.xml +-- a1 = Poly [Term(1,m[7]), Term(-1,mp[2][7])] :: Poly Rational Revlex +-- a2 = Poly [Term(1,mp[2][7]), Term(-1,mp[3][7])] :: Poly Rational Revlex +-- a3 = Poly [Term(1,mp[3][7]), Term(-1,mp[4][7])] :: Poly Rational Revlex +-- a4 = Poly [Term(1,mp[4][7]), Term(-1, mp[5][7])] :: Poly Rational Revlex +-- a5 = Poly [Term(1,m[6,1]), Term(1, mp[2,3][6,1]),Term(1, mp[3,4][6,1]), Term(1, mp[4,5][6,1]), Term(1, mp[1,5][1,6]) ] :: Poly Rational Revlex +-- ps1= a1 : a2: a3: a4: a5: [] +-- -- http://symbolicdata.org/XMLResources/IntPS/ZeroDim.example_1.xml +-- a6 = Poly [Term(1,m[2]), Term(1, mp[2][2]), Term(-10,m[0])] :: Poly Rational Revlex +-- a7 = Poly [Term(1,m[2]), Term(1, m[1,1]), Term(2,mp[2][2]), Term(-16,m[0])] :: Poly Rational Revlex +-- ps2 = a6 : a7 :[] + +--a6 = Poly [ Term(1,mp[2][1]), Term(1,mp[3][2]), Term(1,mp[4][2]), Term(-1,m[2]) ] :: Poly Rational Revlex +--a7 = Poly [Term(1,mp[2,3][1,1]), Term(1,mp[4][2]), Term(-1,m[])] :: Poly Rational Revlex +--a8 = Poly [Term(1,mp[2,3,4][1,1,1]), Term(-1,mp[2][2]), Term(-1,mp[3][2]), Term(-1,mp[4][1]), Term(1,m[])] :: Poly Rational Revlex +--ps2 = a6 : a7: a8 : [] +-- ******************************************************************************************************************************************************************************************************** +f1 = p [tp 1 [2] [2], t (-1) [1, 1], t (-1) []] :: Poly Rational Revlex +f2 = p [t 1 [2], tp (-2) [3] [1]] :: Poly Rational Revlex +f3 = p [tp 1 [3] [2], t (-1) [1, 1], t 1 []] :: Poly Rational Revlex +ps1 = parps [f2,f1, f3] +---------------------------------------------------------------------------------------------------- +a6 = + p [tp 1 [2] [2], tp 1 [3] [2], tp 1 [4] [2], t (-1) [2]] :: Poly Rational Revlex +a7 = + p [tp 1 [2, 3] [1, 1], tp 1 [4] [2], t (-1) []] :: Poly Rational Revlex +a8 = + p + [ tp 1 [2, 3, 4] [1, 1, 1] + , tp (-1) [2] [2] + , tp (-1) [3] [2] + , tp (-1) [4] [1] + , t 1 [] + ] :: Poly Rational Revlex +ps2' = a6 : a7 : a8 : [] +ps2 = parps ps2' +---------------------------------------------------------------------------------------------------- +a9 = + p + [ tp 1 [1, 10] [2, 1] + , tp 2 [1, 2, 11] [1, 1, 1] + , tp 3 [2, 12] [2, 1] + , tp 1 [1, 4] [1, 1] + , tp 2 [2, 5] [1, 1] + , tp 1 [7] [1] + ] :: Poly Rational Revlex +a10 = + p + [ tp 3 [1, 9] [2, 1] + , tp 2 [1, 2, 10] [1, 1, 1] + , tp 1 [2, 11] [2, 1] + , tp 2 [1, 3] [1, 1] + , tp 1 [2, 4] [1, 1] + , tp 1 [6] [1] + ] :: Poly Rational Revlex +a11 = + p + [ tp 1 [1, 9] [3, 1] + , tp 1 [1, 2, 10] [2, 1, 1] + , tp 1 [1, 2, 11] [1, 2, 1] + , tp 1 [2, 12] [3, 1] + , tp 1 [1, 3] [2, 1] + , tp 1 [1, 2, 4] [1, 1, 1] + , tp 1 [2, 5] [2, 1] + , tp 1 [1, 6] [1, 1] + , tp 1 [2, 7] [1, 1] + , tp 1 [8] [1] + ] :: Poly Rational Revlex +ps3' = a9 : a10 : a11 : [] +ps3 = parps ps3' +---------------------------------------------------------------------------------------------------- +a12 = p [ tp 1 [1,2,3][2,4,1], tp 1 [1,2][2,4], tp 2 [1,2,3][2,2,1], tp 6 [1,2][2,2], tp (-4) [1,2][1,3], tp 1 [2,3][4,1], tp (-1) [2][4], tp 1 [1,3][2,1], tp 1 [1][2], tp 4 [1,2][1,1], tp 2 [2,3][2,1], tp (-6) [2][2], tp 1 [3][1], t (-1) [] + ] :: Poly Rational Revlex +a13 = p [ tp 1 [1,2,3][8,3,1], tp 1 [1,2,3][8,1,1], tp 4 [1,2,3][6,3,1], tp (-8) [1,2][6,3], tp 8 [1,2][5,4], tp 4 [1,2,3][6,1,1], tp 8 [1,2][6,1], tp (-48) [1,2][5,2], tp 6 [1,2,3][4,3,1], tp 48 [1,2][4,3], tp (-8) [1,2][3,4], tp 8 [1][5], tp 6 [1,2,3][4,1,1], tp (-48) [1,2][4,1], tp (48) [1,2][3,2], tp 4 [1,2,3][2,3,1], tp (-8) [1,2][2,3], tp (-8) [1][3], tp 4 [1,2,3][2,1,1], tp (8) [1,2][2,1], tp (1) [2,3][3,1], tp 1 [2,3][1,1] + ] :: Poly Rational Revlex +ps4 = parps [a12,a13] +---------------------------------------------------------------------------------------------------- +a14 = p[ tp 35 [2][1], tp 40 [3][1], tp 25 [4][1], tp (-27) [5][1]] :: Poly Rational Revlex +a15 = p[ tp 45 [2][1], tp 35 [5][1], tp (-165) [6][1], t (-36) []] :: Poly Rational Revlex +a16 = p[ tp (-11)[5,6][1,1], tp (3)[6][2],tp (99) [1][1] ] :: Poly Rational Revlex +a17 = p[ tp 25 [2,5][1,1], tp (-165)[6][2], tp (15) [1][1], tp 30 [3][1], tp (-18) [4][1] ] :: Poly Rational Revlex +a18 = p[ tp 15 [2,4][1,1], tp 20 [3,5][1,1], tp (-9) [1][1]] :: Poly Rational Revlex +a19 = p[ tp (-11)[6][3], tp 1 [1,2][1,1], tp 2 [3,4][1,1] ] :: Poly Rational Revlex +ps5 = parps [a14,a15,a16,a17,a18,a19] +---------------------------------------------------------------------------------------------------- +a20 = p[t (-4) [1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1]] :: Poly Rational Revlex +a21 = p[tp (-4) [1][2], tp 1 [2][2], tp 1 [3][2], tp 1 [4][2], tp 4 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], t (-3) []] :: Poly Rational Revlex +a22 = p[tp (5) [1][3], tp 4 [3,4][2,1], tp 3 [2][2], tp 2 [1,4][1,1], tp 4 [1][1], tp 1 [2][1], tp 1 [3][1], tp 2 [4][1], t (-1) []] :: Poly Rational Revlex +a23 = p[tp 5 [3][4], tp 1 [4][3], tp 16 [1][2], tp 3 [2][2], tp (-4) [4][1], t (-1) []] :: Poly Rational Revlex +ps6 = parps [a20,a21,a22,a23] +---------------------------------------------------------------------------------------------------- +a24 = p[tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1]] :: Poly Rational Revlex +a25 = p[tp 1 [1,2][1,1], tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [1,5][1,1], tp 1 [4,5][1,1]] :: Poly Rational Revlex +a26 = p[tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1 [1,2,5][1,1,1], tp 1 [1,4,5][1,1,1], tp 1 [3,4,5][1,1,1]] :: Poly Rational Revlex +a27 = p[tp 1 [1,2,3,4][1,1,1,1], tp 1 [1,2,3,5][1,1,1,1], tp 1 [1,2,4,5][1,1,1,1], tp 1 [1,3,4,5][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1]] :: Poly Rational Revlex +a28 = p[t 1 [1,1,1,1,1], t (-1) []] :: Poly Rational Revlex +ps7 = parps [a24,a25,a26, a27, a28] +---------------------------------------------------------------------------------------------------- +a29 = p[tp (-2) [1,3][1,1], tp (-2) [3][2], tp (-2)[1][1], tp (8) [3][1], t (-2) []]:: Poly Rational Revlex +a30 = p[tp (-3)[1,2,3][1,1,1], tp (2) [1,3,4][1,1,1], tp (4)[3,4][2,1], tp (3)[2,3][1,1], tp 1 [1,4][1,1], tp (-7) [3,4][1,1]]:: Poly Rational Revlex +a31 = p[tp 1 [1,2][2,2], tp (-2)[1,2,4][2,1,1], t (-2) [1,1,1,1], tp 1 [1,4][2,2], tp 2 [1,3,4][1,1,2], tp 1 [3,4][2,2], tp (-2)[1,2][1,2], t (4) [1,1,1], tp 2 [1,2,4][1,1,1], tp 2 [2,3,4][1,1,1], tp (-4)[1,4][1,2], tp (-4) [3,4][1,2], tp 1 [2][2], tp 2 [1,3][1,1], tp 10 [3][2], tp (-4)[2,4][1,1], tp 4 [4][2], tp (-2)[1][1], tp (-8)[3][1], t 2 [] ]:: Poly Rational Revlex +a32 = p[t 3 [2,2], t 12 [1,1,1,1], tp (-3) [1,4][2,2], tp 6 [1,3,4][1,1,2], tp (-3) [3,4][2,2], tp (-6) [1,2][1,2], tp 12 [1,2,4][1,1,1], tp 12 [2,3,4][1,1,1], tp (-4) [1][2], tp 3 [2][2], tp 5 [3][2], tp (-12) [2,4][1,1], tp 12 [4][2], tp (-6) [3][1], t 5 [] ] :: Poly Rational Revlex +ps8 = parps [a29,a30,a31,a32] +---------------------------------------------------------------------------------------------------- +-- ojito + +---------------------------------------------------------------------------------------------------- +a39 = p[tp 1 [2,3][1,5], t 1 [1], t (-2)[] ]:: Poly Rational Revlex +a40 = p[tp 1 [1,3][5,1], tp 1 [2][1], t (-2)[]]:: Poly Rational Revlex +a41 = p[t 1 [1,5], tp 1 [3][1], t (-2) []]:: Poly Rational Revlex +ps11 = parps [a39,a40,a41] +---------------------------------------------------------------------------------------------------- +a46 = p[tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1], tp 1 [6][1]] :: Poly Rational Revlex +a47 = p[tp 1 [1,2][1,1], tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [4,5][1,1], tp 1 [1,6][1,1], tp 1 [5,6][1,1]] :: Poly Rational Revlex +a48 = p[tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1[3,4,5][1,1,1], tp 1 [1,2,6][1,1,1], tp 1 [1,5,6][1,1,1], tp 1 [4,5,6][1,1,1]] :: Poly Rational Revlex +a49 = p[tp 1 [1,2,3,4][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1], tp 1 [1,2,3,6][1,1,1,1], tp 1 [1,2,5,6][1,1,1,1], tp 1 [1,4,5,6][1,1,1,1], tp 1 [3,4,5,6][1,1,1,1]] :: Poly Rational Revlex +a50 = p[tp 1 [1,2,3,4,5][1,1,1,1,1], tp 1 [1,2,3,4,6][1,1,1,1,1], tp 1 [1,2,3,5,6][1,1,1,1,1], tp 1 [1,2,4,5,6][1,1,1,1,1], tp 1 [1,3,4,5,6][1,1,1,1,1], tp 1 [2,3,4,5,6][1,1,1,1,1]] :: Poly Rational Revlex +a51 = p[t 1 [1,1,1,1,1,1], t (-1) [] ] :: Poly Rational Revlex +ps13 = parps [a46,a47,a48,a49,a50,a51] +---------------------------------------------------------------------------------------------------- +a52 = p [t 1 [1,1],tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [1,5][1,1], tp 1 [4,5][1,1]] :: Poly Rational Revlex +a53 = p [tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1]] :: Poly Rational Revlex +a54 = p [tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1 [1,2,5][1,1,1], tp 1 [1,4,5][1,1,1], tp 1 [3,4,5][1,1,1]] :: Poly Rational Revlex +a55 = p [t 1 [1,1,1,1], tp 1 [1,2,4,5][1,1,1,1], tp 1 [1,3,4,5][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1], tp 1 [2,3,4][1,1,1] ] :: Poly Rational Revlex +a56 = p [ t 1 [1,1,1,1,1], t (-1) []] :: Poly Rational Revlex +ps14 = parps [a52, a53,a54,a55,a56] +---------------------------------------------------------------------------------------------------- +a57 = p [tp 1 [1,2,3][2,1,1], tp 1 [1,2,3][1,2,1], tp 1 [1,2,3][1,1,2], tp 1 [1,1,1][1,1,1], tp 1 [1,2][1,1], tp 1 [1,3][1,1], tp 1 [2,3][1,1]] :: Poly Rational Revlex +a58 = p [tp 1 [1,2,3][2,2,1], tp 1 [1,2,3][1,2,2], tp 1 [1,2,3][2,1,1], t 1 [1,1,1], tp 1 [2,3][1,1], t 1 [1], tp 1 [3][1] ] :: Poly Rational Revlex +a59 = p [t 1 [2,2,2], t 1 [2,2,1], t 1 [1,2,1], t 1 [1,1,1], tp 1 [1,3][1,1], tp 1 [3][1], t 1 []] :: Poly Rational Revlex +ps15 = parps [a57, a58, a59] +---------------------------------------------------------------------------------------------------- + + +main :: IO() +main = do --putStrLn "Para usar los benchmar usar cabal new-bench" + print $ charSetMPS ps1 + --charSetPfS pall + --charSetPfPr pallR + --charSetNormalPS ps2 + + + + +-- documentation: cabal haddock +-- bench: cabal new-bench -02, with the compiled file use the bgroup optins +-- execute: cabal new-run diff --git a/bench/Benchmarking.hs b/bench/Benchmarking.hs new file mode 100644 --- /dev/null +++ b/bench/Benchmarking.hs @@ -0,0 +1,326 @@ +module Main(main) where +import Criterion.Main +import Prelude as P +import Polynomial.Monomial +import Polynomial.Wu +import Polynomial.Polynomial +import Data.Massiv.Array as A + +a33 = p[tp 1 [2][4], tp 1 [1,2,3][1,2,1], t 1 [2], t (-2)[1,1], tp 1 [2][1], tp 1 [3][2] ] :: Poly Rational Revlex +a34 = p[tp 1 [1,2][1,4], tp 1 [2,3][1,4], t (-2) [2,1], t (-3) []]:: Poly Rational Revlex +a35 = p[t (-1) [3,2], t 1 [1,1,3], tp 1 [2][4], t 1 [1,2,1], t (-2)[1,1]]:: Poly Rational Revlex +ps9 = parps [a33,a34,a35] +--import Polynomial.Terms +main :: IO () +main + --Example 1 book calculator + = do + let f1 = p [tp 1 [2] [2], t (-1) [1, 1], t (-1) []] :: Poly Rational Revlex + f2 = p [t 1 [2], tp (-2) [3] [1]] :: Poly Rational Revlex + f3 = p [tp 1 [3] [2], t (-1) [1, 1], t 1 []] :: Poly Rational Revlex + ps1 = parps [f2,f1, f3] + -- -------------------------------------------------- + -- -- example 6.1 pag 15 (Wang 2004, Epsilon-A14) + a6 = + p [tp 1 [2] [2], tp 1 [3] [2], tp 1 [4] [2], t (-1) [2]] :: Poly Rational Revlex + a7 = + p [tp 1 [2, 3] [1, 1], tp 1 [4] [2], t (-1) []] :: Poly Rational Revlex + a8 = + p + [ tp 1 [2, 3, 4] [1, 1, 1] + , tp (-1) [2] [2] + , tp (-1) [3] [2] + , tp (-1) [4] [1] + , t 1 [] + ] :: Poly Rational Revlex + ps2' = a6 : a7 : a8 : [] + ps2 = parps ps2' + -- http://symbolicdata.org/XMLResources/IntPS/DiscrC2.xml + a9 = + p + [ tp 1 [1, 10] [2, 1] + , tp 2 [1, 2, 11] [1, 1, 1] + , tp 3 [2, 12] [2, 1] + , tp 1 [1, 4] [1, 1] + , tp 2 [2, 5] [1, 1] + , tp 1 [7] [1] + ] :: Poly Rational Revlex + a10 = + p + [ tp 3 [1, 9] [2, 1] + , tp 2 [1, 2, 10] [1, 1, 1] + , tp 1 [2, 11] [2, 1] + , tp 2 [1, 3] [1, 1] + , tp 1 [2, 4] [1, 1] + , tp 1 [6] [1] + ] :: Poly Rational Revlex + a11 = + p + [ tp 1 [1, 9] [3, 1] + , tp 1 [1, 2, 10] [2, 1, 1] + , tp 1 [1, 2, 11] [1, 2, 1] + , tp 1 [2, 12] [3, 1] + , tp 1 [1, 3] [2, 1] + , tp 1 [1, 2, 4] [1, 1, 1] + , tp 1 [2, 5] [2, 1] + , tp 1 [1, 6] [1, 1] + , tp 1 [2, 7] [1, 1] + , tp 1 [8] [1] + ] :: Poly Rational Revlex + ps3' = a9 : a10 : a11 : [] + ps3 = parps ps3' + --------------------------------------------------- + -- http://symbolicdata.org/XMLResources/IntPS/Schiele_1.xml + + a12 = p [ tp 1 [1,2,3][2,4,1], tp 1 [1,2][2,4], tp 2 [1,2,3][2,2,1], tp 6 [1,2][2,2], tp (-4) [1,2][1,3], tp 1 [2,3][4,1], tp (-1) [2][4], tp 1 [1,3][2,1], tp 1 [1][2], tp 4 [1,2][1,1], tp 2 [2,3][2,1], tp (-6) [2][2], tp 1 [3][1], t (-1) [] + ] :: Poly Rational Revlex + a13 = p [ tp 1 [1,2,3][8,3,1], tp 1 [1,2,3][8,1,1], tp 4 [1,2,3][6,3,1], tp (-8) [1,2][6,3], tp 8 [1,2][5,4], tp 4 [1,2,3][6,1,1], tp 8 [1,2][6,1], tp (-48) [1,2][5,2], tp 6 [1,2,3][4,3,1], tp 48 [1,2][4,3], tp (-8) [1,2][3,4], tp 8 [1][5], tp 6 [1,2,3][4,1,1], tp (-48) [1,2][4,1], tp (48) [1,2][3,2], tp 4 [1,2,3][2,3,1], tp (-8) [1,2][2,3], tp (-8) [1][3], tp 4 [1,2,3][2,1,1], tp (8) [1,2][2,1], tp (1) [2,3][3,1], tp 1 [2,3][1,1] + ] :: Poly Rational Revlex + ps4 = parps [a12,a13] + -------------------------------------------------- + -- http://symbolicdata.org/XMLResources/IntPS/Trinks.xml + a14 = p[ tp 35 [2][1], tp 40 [3][1], tp 25 [4][1], tp (-27) [5][1]] :: Poly Rational Revlex + a15 = p[ tp 45 [2][1], tp 35 [5][1], tp (-165) [6][1], t (-36) []] :: Poly Rational Revlex + a16 = p[ tp (-11)[5,6][1,1], tp (3)[6][2],tp (99) [1][1] ] :: Poly Rational Revlex + a17 = p[ tp 25 [2,5][1,1], tp (-165)[6][2], tp (15) [1][1], tp 30 [3][1], tp (-18) [4][1] ] :: Poly Rational Revlex + a18 = p[ tp 15 [2,4][1,1], tp 20 [3,5][1,1], tp (-9) [1][1]] :: Poly Rational Revlex + a19 = p[ tp (-11)[6][3], tp 1 [1,2][1,1], tp 2 [3,4][1,1] ] :: Poly Rational Revlex + ps5 = parps [a14,a15,a16,a17,a18,a19] + -- -------------------------------------------------- + -- http://symbolicdata.org/XMLResources/IntPS/ZeroDim.example_14.xml + a20 = p[t (-4) [1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1]] :: Poly Rational Revlex + a21 = p[tp (-4) [1][2], tp 1 [2][2], tp 1 [3][2], tp 1 [4][2], tp 4 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], t (-3) []] :: Poly Rational Revlex + a22 = p[tp (5) [1][3], tp 4 [3,4][2,1], tp 3 [2][2], tp 2 [1,4][1,1], tp 4 [1][1], tp 1 [2][1], tp 1 [3][1], tp 2 [4][1], t (-1) []] :: Poly Rational Revlex + a23 = p[tp 5 [3][4], tp 1 [4][3], tp 16 [1][2], tp 3 [2][2], tp (-4) [4][1], t (-1) []] :: Poly Rational Revlex + ps6 = parps [a20,a21,a22,a23] + -------------------------------------------------- + -- http://symbolicdata.org/XMLResources/IntPS/Cyclic_5.xml + -- v,w,x,y,z + -- 1,2,3,4,5 + a24 = p[tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1]] :: Poly Rational Revlex + a25 = p[tp 1 [1,2][1,1], tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [1,5][1,1], tp 1 [4,5][1,1]] :: Poly Rational Revlex + a26 = p[tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1 [1,2,5][1,1,1], tp 1 [1,4,5][1,1,1], tp 1 [3,4,5][1,1,1]] :: Poly Rational Revlex + a27 = p[tp 1 [1,2,3,4][1,1,1,1], tp 1 [1,2,3,5][1,1,1,1], tp 1 [1,2,4,5][1,1,1,1], tp 1 [1,3,4,5][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1]] :: Poly Rational Revlex + a28 = p[t 1 [1,1,1,1,1], t (-1) []] :: Poly Rational Revlex + ps7 = parps [a24,a25,a26, a27, a28] + -------------------------------------------------- + -- http://symbolicdata.org/XMLResources/IntPS/Fee_1.xml q,c,p,d a29 = + -- q,c,p,d + -- 1,2,3,4 + a29 = p[tp (-2) [1,3][1,1], tp (-2) [3][2], tp (-2)[1][1], tp (8) [3][1], t (-2) []]:: Poly Rational Revlex + a30 = p[tp (-3)[1,2,3][1,1,1], tp (2) [1,3,4][1,1,1], tp (4)[3,4][2,1], tp (3)[2,3][1,1], tp 1 [1,4][1,1], tp (-7) [3,4][1,1]]:: Poly Rational Revlex + a31 = p[tp 1 [1,2][2,2], tp (-2)[1,2,4][2,1,1], t (-2) [1,1,1,1], tp 1 [1,4][2,2], tp 2 [1,3,4][1,1,2], tp 1 [3,4][2,2], tp (-2)[1,2][1,2], t (4) [1,1,1], tp 2 [1,2,4][1,1,1], tp 2 [2,3,4][1,1,1], tp (-4)[1,4][1,2], tp (-4) [3,4][1,2], tp 1 [2][2], tp 2 [1,3][1,1], tp 10 [3][2], tp (-4)[2,4][1,1], tp 4 [4][2], tp (-2)[1][1], tp (-8)[3][1], t 2 [] ]:: Poly Rational Revlex + a32 = p[t 3 [2,2], t 12 [1,1,1,1], tp (-3) [1,4][2,2], tp 6 [1,3,4][1,1,2], tp (-3) [3,4][2,2], tp (-6) [1,2][1,2], tp 12 [1,2,4][1,1,1], tp 12 [2,3,4][1,1,1], tp (-4) [1][2], tp 3 [2][2], tp 5 [3][2], tp (-12) [2,4][1,1], tp 12 [4][2], tp (-6) [3][1], t 5 [] ] :: Poly Rational Revlex + ps8 = parps [a29,a30,a31,a32] + -------------------------------------------------- + -- http://symbolicdata.org/XMLResources/IntPS/Weispfenning-94.xml + a33 = p[tp 1 [2][4], tp 1 [1,2,3][1,2,1], t 1 [2], t (-2)[1,1], tp 1 [2][1], tp 1 [3][2] ] :: Poly Rational Revlex + a34 = p[tp 1 [1,2][1,4], tp 1 [2,3][1,4], t (-2) [2,1], t (-3) []]:: Poly Rational Revlex + a35 = p[t (-1) [3,2], t 1 [1,1,3], tp 1 [2][4], t 1 [1,2,1], t (-2)[1,1]]:: Poly Rational Revlex + ps9 = parps [a33,a34,a35] + -------------------------------------------------- + -- http://symbolicdata.org/XMLResources/IntPS/Fateman.xml + a36 = p[t 2 [3], tp 2 [2][3], tp 2 [3][3], tp 1 [4][3]]:: Poly Rational Revlex + -- 6*p^4*q+4*p^3*q^2+2*p^2*q^3+4*p*q^4+6*p^4*r+8*p^3*q*r+4*p*q^3*r+6*q^4*r+4*p^3*r^2+4*q^3*r^2+2*p^2*r^3+4*p*q*r^3+2*q^2*r^3+4*p*r^4+4*q*r^4+3*p^4*s+4*p^3*q*s+2*p*q^3*s+3*q^4*s+4*p^3*r*s+4*q^3*r*s+2*p*r^3*s+2*q*r^3*s+3*r^4*s+p^3*s^2+q^3*s^2+r^3*s^2+p^2*s^3+2*p*q*s^3+q^2*s^3+2*p*r*s^3+2*q*r*s^3+r^2*s^3+p*s^4+q*s^4+r*s^4-s^5 + a37 = p[]:: Poly Rational Revlex + a38 = p[t 2 [5], tp 2 [2][5], tp 2 [3][5], tp 1 [4][5]]:: Poly Rational Revlex + ps10 = parps [a36,a37,a38] + -- http://symbolicdata.org/XMLResources/IntPS/Sym3_5.xml + a39 = p[tp 1 [2,3][1,5], t 1 [1], t (-2)[] ]:: Poly Rational Revlex + a40 = p[tp 1 [1,3][5,1], tp 1 [2][1], t (-2)[]]:: Poly Rational Revlex + a41 = p[t 1 [1,5], tp 1 [3][1], t (-2) []]:: Poly Rational Revlex + ps11 = parps [a39,a40,a41] + -- http://symbolicdata.org/XMLResources/IntPS/Wu-90.xml -- no seguro 4 + -- http://symbolicdata.org/XMLResources/IntPS/Cyclic_6.xml + a46 = p[tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1], tp 1 [6][1]] :: Poly Rational Revlex + a47 = p[tp 1 [1,2][1,1], tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [4,5][1,1], tp 1 [1,6][1,1], tp 1 [5,6][1,1]] :: Poly Rational Revlex + a48 = p[tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1[3,4,5][1,1,1], tp 1 [1,2,6][1,1,1], tp 1 [1,5,6][1,1,1], tp 1 [4,5,6][1,1,1]] :: Poly Rational Revlex + a49 = p[tp 1 [1,2,3,4][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1], tp 1 [1,2,3,6][1,1,1,1], tp 1 [1,2,5,6][1,1,1,1], tp 1 [1,4,5,6][1,1,1,1], tp 1 [3,4,5,6][1,1,1,1]] :: Poly Rational Revlex + a50 = p[tp 1 [1,2,3,4,5][1,1,1,1,1], tp 1 [1,2,3,4,6][1,1,1,1,1], tp 1 [1,2,3,5,6][1,1,1,1,1], tp 1 [1,2,4,5,6][1,1,1,1,1], tp 1 [1,3,4,5,6][1,1,1,1,1], tp 1 [2,3,4,5,6][1,1,1,1,1]] :: Poly Rational Revlex + a51 = p[t 1 [1,1,1,1,1,1], t (-1) [] ] :: Poly Rational Revlex + ps13 = parps [a46,a47,a48,a49,a50,a51] +--http://symbolicdata.org/XMLResources/IntPS/Cyclic_5_1.xml + a52 = p [t 1 [1,1],tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [1,5][1,1], tp 1 [4,5][1,1]] :: Poly Rational Revlex + a53 = p [tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1]] :: Poly Rational Revlex + a54 = p [tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1 [1,2,5][1,1,1], tp 1 [1,4,5][1,1,1], tp 1 [3,4,5][1,1,1]] :: Poly Rational Revlex + a55 = p [t 1 [1,1,1,1], tp 1 [1,2,4,5][1,1,1,1], tp 1 [1,3,4,5][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1], tp 1 [2,3,4][1,1,1] ] :: Poly Rational Revlex + a56 = p [ t 1 [1,1,1,1,1], t (-1) []] :: Poly Rational Revlex + ps14 = parps [a52, a53,a54,a55,a56] +--http://symbolicdata.org/XMLResources/IntPS/Cyclic_7_1.xml + a57 = p [tp 1 [1,2,3][2,1,1], tp 1 [1,2,3][1,2,1], tp 1 [1,2,3][1,1,2], tp 1 [1,1,1][1,1,1], tp 1 [1,2][1,1], tp 1 [1,3][1,1], tp 1 [2,3][1,1]] :: Poly Rational Revlex + a58 = p [tp 1 [1,2,3][2,2,1], tp 1 [1,2,3][1,2,2], tp 1 [1,2,3][2,1,1], t 1 [1,1,1], tp 1 [2,3][1,1], t 1 [1], tp 1 [3][1] ] :: Poly Rational Revlex + a59 = p [t 1 [2,2,2], t 1 [2,2,1], t 1 [1,2,1], t 1 [1,1,1], tp 1 [1,3][1,1], tp 1 [3][1], t 1 []] :: Poly Rational Revlex + ps15 = parps [a57, a58, a59] +-- http://symbolicdata.org/XMLResources/IntPS/Czapor-86a.xml + a60 = p [] :: Poly Rational Revlex + a61 = p [] :: Poly Rational Revlex + a62 = p [] :: Poly Rational Revlex + ps16 = parps [a60, a61, a62] + defaultMain + -- ******************************************** + [ bgroup + "Computing CharSet" + [ + -- bgroup + -- "SP1" + -- [ + -- bench "nf" $ nf (charSetMSP) ps1 + -- --bench "whnf" $ whnf (charSetMSP) ps1 + -- ] + -- , + -- bgroup + -- "PS1" + -- [bench "nf" $ nf (charSetMPS) ps1 + -- --bench "nf" $ nf (A.sum) (fromLists' Par [[0,0,0,0,0],[0,1,2,3,4],[0,2,4,6,8]] :: Array U Ix2 Double) + -- --bench "whnf" $ whnf (charSetMPS) ps1 +-- ] + -- , + -- bgroup + -- "SP2" + -- [bench "nf" $ nf (charSetMSP) ps2 + -- --bench "whnf" $ whnf (charSetMSP) ps2 + -- ] + -- , + -- bgroup + -- "PS2" + -- [bench "nf" $ nf (charSetMPS) ps2 + -- -- bench "whnf" $ whnf (charSetMPS) ps2 + -- ] + -- , + -- bgroup + -- "SP3" + -- [bench "nf" $ nf (charSetMSP) ps3 + -- --bench "whnf" $ whnf (charSetMSP) ps3 + -- ] + -- , + -- bgroup + -- "PS3" + -- [bench "nf" $ nf (charSetMPS) ps3 + -- --bench "whnf" $ whnf (charSetMPS) ps3 + -- ] + -- , + -- bgroup + -- "SP4" + -- [bench "nf" $ nf (charSetMSP) ps4 + -- --bench "whnf" $ whnf (charSetMSP) ps4 + -- ] + -- , + -- bgroup + -- "PS4" + -- [bench "nf" $ nf (charSetMPS) ps4 + -- --bench "whnf" $ nf (charSetMPS) ps4 + -- ] + -- , + -- bgroup + -- "SP5" + -- [bench "nf" $ nf (charSetMSP) ps5 + -- --bench "whnf" $ nf (charSetMSP) ps5 + -- ] + -- , + -- bgroup + -- "PS5" + -- [bench "nf" $ nf (charSetMPS) ps5 + -- --bench "whnf" $ nf (charSetMPS) ps5 + -- ] + -- , + -- bgroup + -- "SP6" + -- [bench "nf" $ nf (charSetMSP) ps6 + -- --bench "whnf" $ nf (charSetMSP) ps6 + -- ] + -- , + -- bgroup + -- "PS6" + -- [bench "nf" $ nf (charSetMPS) ps6 + -- --bench "whnf" $ nf (charSetMPS) ps6 + -- ] + -- , + -- bgroup + -- "SP7" + -- [bench "nf" $ nf (charSetMSP) ps7 + -- --bench "whnf" $ nf (charSetMSP) ps7 + -- ] + -- , + -- bgroup + -- "PS7" + -- [bench "nf" $ nf (charSetMPS) ps7 + -- --bench "whnf" $ nf (charSetMPS) ps7 + -- ] + -- , + -- bgroup + -- "SP8" + -- [bench "nf" $ nf (charSetMSP) ps8 + -- --bench "whnf" $ nf (charSetMSP) ps8 + -- ] + -- , + bgroup + "PS8" + [bench "nf" $ nf (charSetMPS) ps8 + --bench "whnf" $ nf (charSetMPS) ps8 + ] + -- , + -- bgroup + -- "SP9" + -- [bench "nf" $ nf (charSetMSP) ps11 + -- --bench "whnf" $ nf (charSetMSP) ps11 + -- ] + -- , + -- bgroup + -- "PS9" + -- [bench "nf" $ nf (charSetMPS) ps11 + -- --bench "whnf" $ nf (charSetMPS) ps11 + -- ] + -- , + -- bgroup + -- "SP10" + -- [bench "nf" $ nf (charSetMSP) ps13 + -- --bench "whnf" $ nf (charSetMPS) ps13 + -- ] + -- , + -- bgroup + -- "PS10" + -- [bench "nf" $ nf (charSetMPS) ps13 + -- --bench "nf" $ nf (charSetMPS) ps13 + -- ] + -- , + -- bgroup + -- "SP11" + -- [bench "nf" $ nf (charSetMSP) ps14 + -- --bench "whnf" $ nf (charSetMPS) ps14 + -- ] + -- , + -- bgroup + -- "PS11" + -- [bench "nf" $ nf (charSetMPS) ps14 + -- --bench "nf" $ nf (charSetMPS) ps14 + -- ] + -- , + -- bgroup + -- "SP12" + -- [bench "nf" $ nf (charSetMSP) ps15 + -- --bench "whnf" $ nf (charSetMPS) ps15 + -- ] + -- , + -- bgroup + -- "PS12" + -- [bench "nf" $ nf (charSetMPS) ps15 + -- --bench "nf" $ nf (charSetMPS) ps15 + -- ] + + ] + ] +-- ./Ritt-Wu-benchmark +RTS -N -l +-- threadscope Ritt-Wu-benchmark.eventlog + +-- el comando anterior genera estos achivos por lo cual es innecesario la siguiente linea +-- ./Ritt-Wu-benchmark +RTS -p -RTS +-- cat Ritt-Wu-benchmark.prof > profiling.txt + +-- hp2ps -c Ritt-Wu-benchmark.hp +-- gnome-open Ritt-Wu-benchmark.ps + --return $ charSetMPS ps6 -- excelent + --return $ charSetMSP ps5 -- excelent to parallel + --return $ charSetMPS ps7 + --return $ charSetMSP ps8 --- parallel + --return $ charSetMPS ps9 -- parallel + --return $ charSetMPS ps1 diff --git a/src/Polynomial/Monomial.hs b/src/Polynomial/Monomial.hs new file mode 100644 --- /dev/null +++ b/src/Polynomial/Monomial.hs @@ -0,0 +1,178 @@ +{-# LANGUAGE DeriveGeneric #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE GeneralizedNewtypeDeriving #-} +{-# LANGUAGE MultiParamTypeClasses #-} +{-| +-- Module : Data.Massiv.Array +-- Copyright : +-- License : BSD3 +-- Maintainer : Jose Seraquive +-- Stability : experimental +-- Portability : non-portable +-} +module Polynomial.Monomial + ( + -- Types + Mon(..), + Lex(..), + Revlex(..), + -- Classes + + -- Functions + m, + mp, + + ) +where +import Data.Massiv.Array as A +import Prelude as P +import Data.Char.SScript +import Numeric.Algebra as N +import Data.Function +import Control.DeepSeq +import GHC.Generics (Generic, Generic1) + + +-- | A wrapper for monomials with a certain (monomial) order +newtype Mon ord = Mon (Array P Ix1 Int) deriving (Generic, NFData, Eq) +--instance NFData (Mon ord) where +-- rnf x = seq x () + +-- * Names for orderings. +-- We didn't choose to define one single type for ordering names for the extensibility. + +-- | Lexicographical order +data Lex = Lex +-- | Reverse lexicographical order +data Revlex = Revlex + +---multiDeg :: Mon m -> Mon m -> [Int] +-- ----------------------<< FUNCTIONS >>-------------------- +-- | Monomial with terms +m :: [Int] -> Mon ord +m [] = Mon $ A.fromList (ParN 1) [] +m [0] = Mon $ A.fromList (ParN 1) [] --error "A empty term should be represented as m[]" +m xs = Mon $ A.fromList (ParN 1) xs + +-- Function that recive the x_i position with the corresponding exp +-- | Monomial with the term position +mp :: [Int] -> [Int] -> Mon ord +mp xs xz + | lxs /= lxz = error "The size of the position and exponent doesn't correspond" + | otherwise = m $ mpRev xs xz 1 + where + lxs = length xs + lxz = length xz + +mpRev :: [Int] -> [Int] -> Int -> [Int] +mpRev [] [] _ = [] +mpRev (p:ps) (x:xs) n + | p == 0 = error "The initial monomial position is 1 not 0" + | n /= p = 0 : mpRev (p : ps) (x : xs) next + | otherwise = x : mpRev ps xs next + where + next = n P.+ 1 +----------------------------------------------------------------------------------------- +-- lex order +-- mp' :: [Int] -> [Int] -> Mon ord +-- mp' xs xz +-- | lxs /= lxz = error "The size of the position and exponent doesn't correspond" +-- | otherwise = m $ reverse (mpRev xs xz 1) +-- where +-- lxs = length xs +-- lxz = length xz +----------------------------------------------------------------------------------------- +------- <>-------------- +instance Show (Mon ord) where + show (Mon m) = formatSS $ showMon (A.toList m) 1 + +showMon :: (Num a, Eq a, Ord a, Show a) => [a] -> Int -> String +showMon [] _ = "" +showMon (x:xs) s + | x == 0 = printMon + | null xs = format + | otherwise = format ++ printMon + where + next = succ + format = "x_{" ++ show s ++ "}^{" ++ show x ++ "}" + printMon = showMon xs (next s) +--------------------------------------------------------------------------------------- +instance Multiplicative (Mon ord) where + (*) (Mon n) (Mon n') = m (quitZero $ on aux toList n n') + +aux :: [Int] -> [Int] -> [Int] +aux [] [] = [] +aux (x:xs) [] = x : aux xs [] +aux [] (y:yp) = y : aux [] yp +aux (x:xs) (y:yp) = x P.+ y : aux xs yp + +quitZero :: [Int] -> [Int] +quitZero [] = [] +quitZero xs + | last xs == 0 = quitZero $ init xs + | otherwise = xs +--------------------------------------------------------------------------------------- +instance Division (Mon ord) where + (/) (Mon n) (Mon n') = m (quitZero $ on aux' (A.toList) n n') + +aux' :: [Int] -> [Int] -> [Int] +aux' [] [] = [] +aux' (x:xs) [] = x : aux' xs [] +aux' [] (y:yp) = y : aux' [] yp +aux' (x:xs) (y:yp) = x P.- y : aux' xs yp + +---------------------------------------pal mar------------------------------------------------ +instance Additive (Mon ord) where + (+) (Mon m)(Mon m') = Mon $ verificationMl m m' +--------------------------------------------------------------------------------------- +instance Semiring (Mon ord) +instance Abelian (Mon ord) +instance Monoidal (Mon ord) where + zero = Mon empty +instance LeftModule Integer (Mon ord) where + (.*) = undefined +instance RightModule Integer (Mon ord) where + (*.) = undefined +instance LeftModule Natural (Mon ord) where + (.*) = undefined +instance RightModule Natural (Mon ord) where + (*.) = undefined +--------------------------------------------------------------------------------------- +instance Group (Mon ord) where + (-) (Mon m)(Mon m') = Mon $ verificationMl m m' + +verificationMl :: Array P Ix1 Int -> Array P Ix1 Int -> Array P Ix1 Int +verificationMl xs xz + | xs == xz = xs + | otherwise = error "The monomial doesn't match " +--------------------------------------------------------------------------------------- + -- sujeto a revision cunado se realize la diviision + -- polynomial, delay es el mas optimo porque no se realizan + -- operaciones "estrictas sobre los monomios" + -- si seq evalua sequencial es correcto (verificar) + -- deberia usar force? + + +--------------------------------------------------------------------------------------- +instance Ord (Mon Lex) where + compare (Mon m)(Mon m') = on lex' A.toList m m' + (>) (Mon m)(Mon m')= on (P.>) A.toList m m' + (<) (Mon m)(Mon m')= on (P.<) A.toList m m' + +lex' :: (Num a, Eq a, Ord a) => [a] -> [a] -> Ordering +lex' [] [] = EQ +lex' [] _ = LT +lex' _ [] = GT +lex' (x:xs) (y:ys) + | x == y = lex' xs ys + | x P.> y = GT + | otherwise = LT +--------------------------------------------------------------------------------------- +instance Ord (Mon Revlex) where + compare (Mon m)(Mon m')= on revlex' A.toList m m' + +revlex' :: (Num a, Eq a, Ord a) => [a] -> [a] -> Ordering +revlex' = on lex' P.reverse +--------------------------------------------------------------------------------------- +instance Unital (Mon ord ) where + one = undefined --Mon empty diff --git a/src/Polynomial/Polynomial.hs b/src/Polynomial/Polynomial.hs new file mode 100644 --- /dev/null +++ b/src/Polynomial/Polynomial.hs @@ -0,0 +1,405 @@ +{-# LANGUAGE DeriveGeneric #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE GeneralizedNewtypeDeriving #-} +{-# LANGUAGE MultiParamTypeClasses #-} + + +module Polynomial.Polynomial + ( + Poly(..), + class', + ld, + prem, + varDegree, + lt, + -- lm, + lv, + initOfv, + initOfv', + max', + spoly, + basicSpoly, + lcmP, + withoutFactor, + lcmT, + factor, + simP, + listElimination, + p, + orden, + tsum, + totalDeg, + numTerms + ) where +import Polynomial.Terms +import Polynomial.Monomial +import Prelude as P +import Data.Massiv.Array as A +import Control.Scheduler +import Numeric.Algebra as N +import Data.Function -- on +import Data.Char.SScript +import Data.List as L +import Control.DeepSeq +import GHC.Generics (Generic, Generic1) + + +newtype Poly t ord = Poly (Array N Ix1 (Term t ord)) deriving (Generic, NFData, Eq) +--p[Term 1 $ m[1,2,3], Term 3 $ m[4,5,6], Term 0 $ m[]] :: Poly Rational Revlex +-- Term 4 $ m[1,2,3,4] :: Term Rational Revlex +--instance NFData (Poly t ord) where +-- rnf x = seq x () + +p :: (NFData k ) => [Term k Revlex] -> Poly k Revlex +p xs = Poly $ A.fromList (Seq) xs + + +---------------------------------------------------- << FUNCTIONS >>----------------------------------------- +--- The class is independent of the monomial order, in theory but in Revlex we take the last elemens in a monomial +class' ::(NFData t ) => Poly t ord -> Int +class' (Poly p) = A.maximum' $ A.map (f) p + where + g (Mon m) = m + f (Term k mon) = A.elemsCount (g mon) +--class' xs = max' $ P.map (elemsCount . getMon . snd . getT) (getP xs) +max' :: [Int] -> Int +max' [] = error "There is no max in a empty list" +max' [x] = x +max' (x:x':xs) = + max' + ((if x >= x' + then x + else x') : + xs) +-------------------------------------------------------------- +-- leading variable +lv :: (NFData t ) => Poly t Revlex -> Int +lv = class' +-- LEX +--lv :: Poly t Lex -> String +--lv p = "x_{" ++ show ( p) ++ "}" +-------------------------------------------------------------- +-- leading degreee +ld :: (NFData t ) => Poly t Revlex -> Int +ld p@(Poly xp) = A.maximum' . findingLast $ A.computeAs N predicate + --max' . P.map last' . toListP $ A.computeAs N predicate + where + -- toListP = A.toList . A.map h + g (Mon m) = m + f (Term k mon) = A.elemsCount (g mon) + -- h (Term k mon) = A.toList (g mon) + predicate = A.filterS (\x -> (f x) == class' p) xp + k (Term k mon) = if elemsCount (g mon) /= 0 then g mon ! ((elemsCount (g mon)) P.-1) else 0 + findingLast = A.map k + + +-- ld :: Poly t Revlex -> Int +-- ld xp = +-- max' $ +-- P.map +-- last' +-- [ A.toList x +-- | x <- P.map (getMon . snd . getT) (getP xp) +-- , class' xp == elemsCount x +-- ] +-- LEX +-- ld1 :: Poly t Lex -> Int +-- ld1 xp = +-- max' $ +-- P.map +-- last' +-- [ A.toList x +-- | x <- P.map (getMon . snd . getT) (getP xp) +-- , class' xp == elemsCount x +-- ] +-- -------------------------------------------------------------- +-- -- multidegree +-- --mdP :: Poly t ord -> [Int] +-- --mdP = undefined +-- -------------------------------------------------------------- +--leading term +lt :: (NFData t, Eq t) => Poly t Revlex -> Term t Revlex +lt poly@(Poly terms)= h . i . computeAs N $ fil terms + where + fil = A.filterS (\x -> f x == ld poly) + g (Mon m) = m + f (Term k mon) = if elemsCount (g mon) /= 0 then g mon ! (elemsCount (g mon) P.-1) else 0 + i = A.quicksort + h = head . A.toList + +varDegree :: (NFData t, Eq t) => Poly t Revlex -> Int -> Int +varDegree (Poly terms ) n = maxel . A.map takeEl . computeAs N $ fil terms + where + fil = A.filterS (\x -> f x >= n) + g (Mon m) = m + f (Term k mon) = if elemsCount (g mon) /= 0 then elemsCount (g mon) else 0 + maxel = N.sum --maximum' + takeEl (Term k mon) = g mon ! (n P.- 1) + +-- lmMax' :: (Eq t )=>[Term t Revlex] -> Term t Revlex +-- lmMax' [x] = x +-- lmMax' (x:x':xs) = lmMax' +-- ((if x >= x' +-- then x +-- else x') : +-- xs) +--MULTIVARIABLE ------------------------------------------------- +-- -- --------------------------------------------------------------- +-- lm :: (Num t, Eq t ) => Poly t Revlex -> Term t Revlex +-- lm xp = Term (1, xs) +-- where +-- (x,xs) = getT $ lt xp +-- -- -------------------------------------------------------------- +-- -- -- Initial of a palynomial with respect to a variable +initOfv :: (NFData t, Eq t) => Poly t Revlex -> Poly t Revlex +initOfv p@(Poly xp) = Poly $ A.computeAs N (A.map takeInit (A.computeAs N predicate)) + where + predicate = A.filterS (\x -> f x == class' p && h x == ld p) xp + f (Term k mon) = A.elemsCount (g mon) + g (Mon m) = m + h (Term k mon) = g mon ! (A.elemsCount (g mon) P.- 1 ) + +takeInit :: Term t ord -> Term t ord +takeInit (Term k mod) + | A.elemsCount (g mod) == 1 = Term k $ m[] + | otherwise = Term k $ Mon ( A.computeAs P ( A.takeS (Sz (A.elemsCount (g mod)) P.- 1 ) (g mod) )) + where + g (Mon m) = m + +initOfv' :: (NFData t, Eq t) => Poly t Revlex -> Poly t Revlex +initOfv' p@(Poly xp) = Poly $ (A.computeAs N predicate) + where + predicate = A.filterS (\x -> f x == class' p && h x == ld p) xp + f (Term k mon) = A.elemsCount (g mon) + g (Mon m) = m + h (Term k mon) = g mon ! (A.elemsCount (g mon) P.- 1 ) + +-- -------------------------------------------------------------- + +-- lc :: (Num t, Eq t ) => Poly t Revlex -> Term t Revlex +-- lc xp = Term (x, m[0]) +-- where +-- (x,xs) = getT $ lt xp +-- -------------------------------------------------------------- +-- -- reduction of the leading monomial +-- red :: (Eq t) => Poly t Revlex -> Poly t Revlex +-- red xs = Poly [x | x <- getP xs, x /= a] +-- where +-- a = lt xs + +-- ------------------------------------------------------------ +-- --Monomial exponentation +-- -- expM :: (Num k) => Term k ord -> Int -> Term k ord +-- -- expM m i = Term (p, Mon $ A.map (P.*i) b) +-- -- where +-- -- (k,m') = getT m +-- -- p = foldl (P.*) 1 $ L.replicate i k +-- -- b = getMon m' +-- -------------------------------------------------------------- +-- --- comparison total degree (muti degree monomial) +-- mdM :: Term k ord -> [Int] +-- mdM = toList . getMon . snd . getT +-- -------------------------------------------------------------- +-- lcm para polynomios +lcmP :: (NFData t, Num t, Ord t)=> Poly t Revlex -> Poly t Revlex -> [Poly t Revlex] +lcmP a@(Poly poly) b@(Poly poly') + | A.elemsCount poly == 1 && A.elemsCount poly' == 1 = (p $ on lcmT (head . A.toList) poly poly' : []) :[] + | A.elemsCount poly' == 1 = withoutFactor a N.* (p $ lcmT (head . A.toList $ poly') (factor a) : []) :[] + | A.elemsCount poly == 1 = withoutFactor b N.* (p $ lcmT (head . A.toList $ poly) (factor b) : []) :[] + | otherwise = [p $ [on lcmT factor a b],withoutFactor a, withoutFactor b] + +factor :: (NFData t, Num t) => Poly t Revlex -> Term t Revlex +factor (Poly poly) = Term 1 $ m commList + where + f (Term k mon) = g mon + g (Mon m) = A.toList m + commList = L.foldl1' (P.zipWith min) (A.toList $ A.map f poly) + + +withoutFactor ::(NFData t, Num t) => Poly t Revlex -> Poly t Revlex +withoutFactor (Poly poly) = Poly $ A.computeAs N ( A.map f poly) + where + f m = quit m m' + m' = factor $ Poly poly + +quit :: Term k ord -> Term k ord -> Term k ord +quit (Term k mon) (Term k' mon') + | mon' == zero = Term k mon + | otherwise = Term k $ m (on quit' f mon mon') + where + f (Mon m) = A.toList m + +quit' :: [Int]->[Int]->[Int] +quit' as [] = as +quit' (a:as) (b:bs) = a P.- b : quit' as bs + + +-- lcm'' :: (Ord t, Num t )=> [Term t Revlex] -> [Term t Revlex] -> [Term t Revlex] +-- lcm'' xs xp +-- | length xs == 1 && length xp == 1 = [ on lcmT head xs xp] +-- | otherwise = undefined + +-- -------------------------------------------------------------- +lcmT :: (Ord t, Num t) => Term t Revlex -> Term t Revlex -> Term t Revlex +lcmT (Term k mon)(Term k' mon')= lcm' mon mon' + +lcm' :: (Ord t, Num t )=> Mon Revlex -> Mon Revlex -> Term t Revlex +lcm' (Mon m)(Mon m') = Term 1 $ Mon (A.computeAs P (A.zipWith max m m' )) +-------------------------------------------------------------- +-- mon :: Term t ord -> Mon ord +-- mon m = Mon $ getMon b +-- where +-- (a,b) = getT m +-- -------------------------------------------------------------- +-- -- spolynomial +basicSpoly ::(NFData t, Fractional t, Num t, Eq t, Ord t ) => Poly t Revlex -> Poly t Revlex -> Poly t Revlex +basicSpoly f g = simP x' (initOfv' f) f N.- simP x' (initOfv' g) g + where + x' = on lcmP initOfv' f g + +--length (getP f) P.> 1 && length (getP g) P.> 1 = [ withoutFactor f, withoutFactor g, Poly[lcmT (factor f) (factor g)] ] + +simP :: (NFData t, Fractional t, Num t, Eq t ) => [Poly t Revlex]->Poly t Revlex->Poly t Revlex->Poly t Revlex +simP f1 f2 f3 + | length f1 == 3 = listElimination f1 f2 f3 + | (orden . head $ f1) == orden f2 = f3 --si + | orden f3 == orden f2 = head f1 --si + | (f . head $ f1) == 1 && (f f2) == 1 = p[h f1 N./ g f2] N.* f3 --tiene sentido que se pregunte por f1 == 1 porque este puede ser un polynomio + | f f2 == 1 = b N.* p [a N./ g f2] N.* f3 + | orden b == orden b' = p[a N./ a'] N.* f3 + | otherwise = error "Option not consider in simP" + where + (a,b) = (factor . head $ f1, (withoutFactor . head) f1) + (a',b') = (factor f2, withoutFactor f2) + f (Poly poly) = A.elemsCount poly + g (Poly poly) = head . A.toList $ poly + h [f1] = g f1 + + + +orden :: (NFData t)=> Poly t Revlex -> Poly t Revlex +orden (Poly poly) = p $ P.map h ( sortBy (\( k, mon) ( k', mon') -> compare mon mon') (A.toList $ A.map f poly )) + where + h (k, xs)= Term k $ m xs + g (Mon m) = m + f (Term k mon) = (k, A.toList (g mon)) + +listElimination :: (NFData t, Fractional t, Eq t, Num t) => [Poly t Revlex] -> Poly t Revlex -> Poly t Revlex -> Poly t Revlex +listElimination [a,b,c] f2 f3 + | orden b == orden e = f3 N.* p[ (head . A.toList $ f a) N./ d] N.* c + | orden c == orden e = f3 N.* p[ (head . A.toList $ f a) N./ d] N.* b + | orden f3 == orden f2 = a N.* b N.* c + | otherwise = error "Option not consider list Elimination" + where + (d,e) = (factor f2, withoutFactor f2) + f (Poly poly) = poly + +spoly :: (NFData t, Fractional t, Ord t) => Poly t Revlex -> Poly t Revlex -> Poly t Revlex +spoly f g + | f == zero = f + | ld f >= ld g && class' g == class' f = spoly (basicSpoly f g) g + | otherwise = f +------------------------------------------------------------ +-- pseudo remainder +prem :: (NFData t, Eq t, Num t) => Poly t Revlex -> Poly t Revlex -> Poly t Revlex +prem g f + | g /= zero && class' g == class' f && ld g >= m = prem calc f + | otherwise = g + where + m = ld f + calc = (initOfv f N.* g) N.- (initOfv g N.* f N.* p [Term 1 $ mp[lv f][ld g P.- m] ]) +------ <> ------------------------------------ +instance (NFData t, Ord t, Num t, Show t) => Show (Poly t ord) where + show (Poly p) = showPoly p + +showPoly :: (NFData t, Ord t, Num t, Show t) => Array N Ix1 (Term t ord) -> String +showPoly terms = (concat . A.toList) $ A.map f terms + where + f (Term k mon) + | k P.> 0 = " + " ++ show k ++ show mon + | k P.< 0 = " - " ++ show (abs k) ++ show mon + | otherwise = error "term with empty coefficient k" + +-- --------------------------------------------------------------- +instance (NFData t, Num t, Eq t) => Additive (Poly t Revlex) where + (+) (Poly poly)(Poly poly') = + p $ + P.filter removeZero $ + P.map tsum $ + groupBy (\( k, mon) ( k', mon') -> (==) mon mon') $ + sortBy (\( k, mon) ( k', mon') -> compare mon mon') (A.toList $ A.map f (computeAs N concatp )) + where + g (Mon m) = m + f (Term k mon) = (k, A.toList (g mon)) + concatp = A.append' 1 poly poly' +--sum term in a polynomial +tsum :: (Num k) => [(k, [Int])] -> Term k ord +tsum [(k,x)] = Term k (m x) +tsum ((k, x):xs) = (Term k $ m x) N.+ tsum xs + +-- remove term 0 form a polynomial +removeZero :: (Num k, Eq k) => Term k ord -> Bool +removeZero (Term k mon) + | k /= 0 = True + | otherwise = False + + +instance (NFData t, Num t, Eq t) => Semiring (Poly t Revlex) +instance (NFData t, Num t, Eq t) => Abelian (Poly t Revlex) +instance (NFData t, Num t, Eq t) => Multiplicative (Poly t Revlex) where + (*) (Poly poly)(Poly poly') = simp $ p [ x N.* y | x <- A.toList poly, y <- A.toList poly' ] + +simp :: (NFData t, Num t, Eq t) => Poly t Revlex -> Poly t Revlex +simp (Poly poly) = + p $ + P.filter removeZero $ + P.map tsum $ + groupBy (\( k, mon) ( k', mon') -> (==) mon mon') $ + sortBy (\( k, mon) ( k', mon') -> compare mon mon') (A.toList $ A.map f poly ) + where + g (Mon m) = m + f (Term k mon) = (k, A.toList (g mon)) +--------------------------------------------------------------------------- + +instance (NFData t, Num t, Eq t) => Group (Poly t Revlex) where + (-) (Poly poly)(Poly poly') = + p $ + P.filter removeZero $ + P.map tsum $ + groupBy (\( k, mon) ( k', mon') -> (==) mon mon') $ + sortBy (\( k, mon) ( k', mon') -> compare mon mon') (A.toList $ A.map f (computeAs N concatp )) + where + g (Mon m) = m + f (Term k mon) = (k, A.toList (g mon)) + h (Term k mon) = Term (- k) mon + concatp = A.append' 1 poly (A.map h poly') + + +instance (NFData t, Num t, Eq t) => LeftModule Integer (Poly t Revlex) where + (.*) = undefined +instance (NFData t, Num t, Eq t) => RightModule Integer (Poly t Revlex) where + (*.) = undefined +instance (NFData t, Num t, Eq t) => LeftModule Natural (Poly t Revlex) where + (.*) = undefined +instance (NFData t, Num t, Eq t) => RightModule Natural (Poly t Revlex) where + (*.) = undefined + +instance (NFData t, Num t, Eq t) => Monoidal (Poly t Revlex) where + zero = p[] --Poly (Array N Ix1 m[]) --zero --p[] + + +last' :: [Int] -> Int +last' [] = 0 +last' xs = last xs + + +-- find the term with the maximum of the sum of all powers +totalDeg :: (NFData t) => Poly t Revlex -> Int +totalDeg (Poly p) = A.maximum' $ A.map (f) p + where + g (Mon m) = m + f (Term k mon) = A.sum (g mon) +-- Contador de terminos en un polinomio +numTerms :: (NFData t) => Poly t Revlex -> Int +numTerms (Poly p) = A.elemsCount p + diff --git a/src/Polynomial/Terms.hs b/src/Polynomial/Terms.hs new file mode 100644 --- /dev/null +++ b/src/Polynomial/Terms.hs @@ -0,0 +1,97 @@ +{-# LANGUAGE DeriveGeneric #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE MultiParamTypeClasses #-} +{-# LANGUAGE GeneralizedNewtypeDeriving #-} + +module Polynomial.Terms + -- Types + ( Term(..) + ) where +import Prelude as P +import Polynomial.Monomial +import Numeric.Algebra as N +import Data.Foldable +import Data.Function +import Control.DeepSeq +import Data.Massiv.Array as A +import Data.Char.SScript +import GHC.Generics (Generic) +-- | Constraint synonym for rings that can be used as polynomial coefficient. + +-- data Term k ord = +-- Term +-- { getT :: (k , Mon ord) +-- } deriving (Eq) + +data Term k ord = Term !k !(Mon ord) deriving (Eq, Generic) -- Term 4 $ m[1,2,3,4] :: Term Rational Revlex +instance (NFData k) => NFData (Term k ord) where + rnf (Term k mon) = rnf k `deepseq` rnf mon `deepseq` () +-- ----------------------<< FUNCTIONS >>-------------------- +-- ----------------------<< INSTANCES >>-------------------- +instance (Num k, Show k, Eq k) => Show (Term k ord) where + show xs = showMon xs + +showMon :: (Num k, Show k, Eq k) => Term k ord -> String +showMon (Term k mon) + | mon == zero = show k + | otherwise = (formatSS . show) k ++ show mon + +----------------------------------------------------------------------------------------- +instance (Num k) => Additive (Term k ord) where + (+) xs xp = addPol' xs xp + +addPol' :: (Num k) => Term k mon -> Term k mon -> Term k mon +addPol' (Term k mon) (Term k' mon') = Term (k P.+ k')( mon N.+ mon') +----------------------------------------------------------------------------------------- +instance (Num k) => Multiplicative (Term k ord) where + (*) xs xz = multPol' xs xz + +multPol' :: (Num k) => Term k mon -> Term k mon -> Term k mon +multPol' (Term k mon) (Term k' mon') = Term (k P.* k')( mon N.* mon') + +----------------------------------------------------------------------------------------- +instance (Fractional k, Num k) => Division (Term k ord) where + (/) xs xp = divPol' xs xp + +divPol' :: (Fractional k, Num k) => Term k mon -> Term k mon -> Term k mon +divPol' (Term k mon) (Term k' mon') = Term (k P./ k')( mon N./ mon') + +instance (Num k) => Unital (Term k ord) where + one = undefined -- Term (0, one) + +instance (Num k) => Semiring (Term k ord) +instance (Num k) => Abelian (Term k ord) +instance (Num k) => Monoidal (Term k ord) where + zero = Term 0 zero -- Term (0, zero) +instance (Num k) => LeftModule Integer (Term k ord) where + (.*) = undefined +instance (Num k) => RightModule Integer (Term k ord) where + (*.) = undefined +instance (Num k) => LeftModule Natural (Term k ord) where + (.*) = undefined +instance (Num k) => RightModule Natural (Term k ord) where + (*.) = undefined +instance (Num k) => Group (Term k ord) where + (-) xs xz = subPol' xs xz + +subPol' ::(Num k) => Term k mon -> Term k mon -> Term k mon +subPol' (Term k mon) (Term k' mon') = Term (k P.- k')( mon N.- mon') + +------------------------------------------------------------------------------- +--instance Eq (Term k ord) where +-- (==) = +-- instance (Eq k)=> Ord (Term k Lex) where +-- compare = on compare (snd . getT) +-- (<) = on (P.<) (snd . getT) +-- (>) = on (P.>) (snd . getT) + +-- λ> a = Term (2,m[5,2,3,2]) :: Term Int Lex +-- λ> b = Term (8,m[4,2,3,2]) :: Term Int Lex +-- λ> a P.> b +-- True +-- λ> + +instance (Eq k)=> Ord (Term k Revlex) where + compare (Term k mon)(Term k' mon') = compare mon mon' + (<) (Term k mon)(Term k' mon') = (P.<) mon mon' + (>) (Term k mon)(Term k' mon') = (P.>) mon mon' diff --git a/src/Polynomial/Wu.hs b/src/Polynomial/Wu.hs new file mode 100644 --- /dev/null +++ b/src/Polynomial/Wu.hs @@ -0,0 +1,660 @@ +{-# LANGUAGE DeriveGeneric #-} +{-# LANGUAGE FlexibleInstances #-} +{-# LANGUAGE GeneralizedNewtypeDeriving #-} +{-# LANGUAGE MultiParamTypeClasses #-} +module Polynomial.Wu ( +PS(..), +-- -- *Division +-- psBsPS, +-- psBsSP, +-- --par +psBsUParP, +psBsUParS, +-- --seq +-- psBsSeqPS, +-- psBsSeqSP, + +-- -- pure IO +-- --psBsSParP, +-- --psBsUParPUpio, + +-- -- * char Set +-- charSetNormalSP, +-- charSetNormalPS, + +charSetMSP, +charSetMPS, + +-- charSetMSeqSP, +-- charSetMSeqPS, + + +t, +tp, +parps + +-- basicSet + ) where +import Control.Scheduler +import Polynomial.Monomial +import Polynomial.Terms +import Polynomial.Polynomial +import Data.List as L +import Data.Function +import Data.Massiv.Array as A +import Control.Monad +import Numeric.Algebra as N +import Prelude as P +import Control.Concurrent +import Foreign.Marshal.Unsafe +import System.IO.Unsafe +import Control.DeepSeq +import GHC.Generics (Generic, Generic1) + +data Idp t ord = Idp !(Poly t ord) !Int deriving (Eq, Generic) +newtype PS t ord = PS (Array N Ix1 (Idp t ord )) deriving (Eq,Generic, NFData) + +instance (NFData t) => NFData (Idp t ord) where + rnf (Idp t ord) = rnf t `deepseq` rnf ord `deepseq` () + +-- polynomial id +idp :: Poly t ord -> Idp t ord +idp xs = Idp xs 0 + + +-- polynomials id to parallel polynomial set +ps :: (NFData t) => [Idp t ord] -> PS t ord +ps xs = PS $ A.fromList (ParN 2) xs + +--wrapper to write terms +t :: k -> [Int] -> Term k ord +t k mon = Term k $ m mon + +--wrapper to write terms with a givne position +tp :: k -> [Int] -> [Int] -> Term k ord +tp k mon mon' = Term k $ mp mon mon' + +--parallel forma of a polynomial set +parps :: (NFData t) =>[Poly t ord] -> PS t ord +parps xs = ps $ P.map idp xs + +-- how write polynomials +p1 = p [t 3 [1,2], tp 4 [2,4][1,2], t 4 [1,2] ] :: Poly Rational Revlex +p2 = p [t 3 [2], t 4 [1] ] :: Poly Rational Revlex +p3 = p [t 3 [5,6], t 4 [1] ] :: Poly Rational Revlex +p4 = p [t 3 [2,1,2], t 4 [1] ] :: Poly Rational Revlex +p5 = p [t 3 [1], t 6 [2] ] :: Poly Rational Revlex +p6 = p [t 3 [1,2,3], t 6 [2] ] :: Poly Rational Revlex +example = parps [p1,p2,p3,p4,p5,p6] + + +---------------------------------------------------------------------------------------------------- +instance (NFData t, Ord t, Num t, Show t) => Show (PS t ord) where + show xs = concat $ showPoly xs + +showPoly :: (NFData t, Ord t, Show t, Num t) => PS t ord -> [String] +showPoly (PS xs) = conv + where + conv = L.map (\ x -> showIdp x ++ "\n") (A.toList xs) + --f (Idp poly id) = show poly ++ "<<" ++show id ++ ">>," +---------------------------------------------------------------------------------------------------- +setClassPs :: (NFData t, Eq t) =>PS t ord -> PS t ord +setClassPs = sortGrupPolyClass . searchPolyClass + +searchPolyClass :: (NFData t )=> PS t ord -> PS t ord +searchPolyClass (PS ps) = PS $ A.computeAs N (A.map f ps) + where + f (Idp poly id) = Idp poly (class' poly) + + +sortGrupPolyClass :: (NFData t, Eq t) => PS t ord -> PS t ord +sortGrupPolyClass (PS xs) = PS conversion + where + conversion = A.quicksort xs + +instance (NFData t, Eq t) => Ord (Idp t ord) where + compare (Idp poly id)(Idp poly' id') = compare id id' + + +-- sucesive pseudo division of a polynomial and a polynomial set (SPOLY) +sspoly :: (NFData t, Fractional t, Num t, Eq t, Ord t) => Idp t Revlex -> PS t Revlex -> Idp t Revlex +sspoly (Idp pol id)(PS xs) = Idp (sspoly' pol (P.map f (A.toList xs))) 0 + where + f (Idp pol' id') = pol' + +sspoly' :: (NFData t, Fractional t, Ord t, Num t, Eq t) => Poly t Revlex -> [Poly t Revlex] -> Poly t Revlex +sspoly' rm [] = rm +sspoly' rm (b:bs) + | a /= p[] = sspoly' a bs + | otherwise = p[] + where + a = spoly rm b + +-- sucesive pseudo division of a polynomial and a polynomial set +sprem :: (NFData t, Num t, Eq t) => Idp t Revlex -> PS t Revlex -> Idp t Revlex +sprem (Idp pol id)(PS xs) = Idp (sprem' pol (P.map f (A.toList xs))) 0 + where + f (Idp pol' id') = pol' + +sprem' :: (NFData t, Num t, Eq t) => Poly t Revlex -> [Poly t Revlex] -> Poly t Revlex +sprem' rm [] = rm +sprem' rm (b:bs) + | a /= p[] = sprem' a bs + | otherwise = p[] + where + a = prem rm b + +-- -- reducction of two list of polynomials +red :: (NFData t, Eq t) => PS t ord -> PS t ord -> PS t ord +red (PS poly)(PS poly') = ps $ on (P.foldl (flip delete)) f poly' poly + where + f = A.toList + + +-- -- Monad IO division +-- psBsUParPUpio :: (Num t, Eq t) =>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- psBsUParPUpio ps bs = unsafePerformIO ( traverseConcurrently Par (\p -> p `deepseq` pure p) [sprem x bs | x <- ps] ) + +----------------------------------------------------------------------------- +-- take a Basic Set BS from a PS (but only lower rank ) +-- basicSet' :: (Num t, Eq t) => [[Poly t Revlex]] -> [Poly t Revlex] +-- basicSet' xs = searchLowerRank xs [] +------------------------------------------------ +-- -- a function that search a poly with lower rank and reduced x_i +-- search Poly Lower Rank and Reduced +groupClass ::(NFData t, Num t, Eq t)=> PS t ord -> [[Idp t ord]] +groupClass (PS poly) = groupBy (on (==) f) $ A.toList poly + where + f (Idp pol id) = id + +sPLRR :: (NFData t, Num t, Eq t) => [[Idp t Revlex]] -> [Idp t Revlex] -> PS t Revlex +sPLRR (x:xs) [] = sPLRR xs [head $ lDPolys x] +sPLRR [] xp = ps xp +sPLRR xs xp + | f (head xp) == 0 = error "Contradictory AS There is (are) constant inside the polynomial set" + | h /= z = sPLRR (tail xs) (h : xp) + | otherwise = sPLRR (tail xs) xp + where + f(Idp poly id) = class' poly + g = lDPolys . head + h = reducedP (g xs) (head xp) + z = idp (p[]) + +-- check if a given polynomial is reduced with respect to a set +reducedP :: (NFData t, Num t, Eq t)=> [Idp t Revlex] -> Idp t Revlex -> Idp t Revlex +reducedP [] _ = idp (p[]) +reducedP (x:xs) z + | isReduced x z = x + | otherwise = reducedP xs z +----------------------------------------------------------------- +-- izq the entering polynomial, right the forming basic set +-- give two polyomials a b check if a is reduced to b +isReduced :: (NFData t, Eq t, Num t) => Idp t Revlex -> Idp t Revlex -> Bool +isReduced (Idp poly id) (Idp poly' id') + -- degP poly (class' poly') P.< degP poly' (class' poly') = True + | prem poly poly' == poly = True + | otherwise = False +-- ----------------------------------------------------------------- +-- compare the degree of two polynomials in a given variable +degP ::(NFData t ) => Poly t Revlex -> Int -> Int +degP ( Poly poly) cls = A.maximum' . A.map h . A.computeAs N . A.filterS (\x -> f x >= cls) $ poly + where + f(Term k mon) = A.elemsCount (g mon) + g(Mon m) = m + h (Term k mon) = (g mon) ! (cls P.-1) + +---------------------------------------------------------------------------------------------------- +--search the polynomial (S) with lower degree in a list of polynomials with the same class +lDPolys :: (NFData t) => [Idp t Revlex]-> [Idp t Revlex] +lDPolys xs = head $ groupBy (on (==) i) (sortBy (on compare i) f) + where + f = L.map g xs + g (Idp poly id) = Idp poly (ld poly) + i (Idp poly id) = id + +----------------------------------------------------------------------------- +quitEmptyPoly :: (NFData t, Eq t) => PS t Revlex -> PS t Revlex +quitEmptyPoly (PS pol)= PS $ A.computeAs N $ ( A.filterS (\x -> x /= f) pol) + where + f = idp $ p[] +---------------------------------------------------------------------------- + +instance (NFData t, Ord t, Num t, Show t) => Show (Idp t ord) where + show xs = showIdp xs + +showIdp ::(NFData t, Ord t, Num t, Show t)=> Idp t ord -> String +showIdp (Idp poly id) = show poly -- ++ "<<" ++show id ++ ">>" +---------------------------------------------------------------- +-- take a basic set using the classification class +basicSet :: (NFData t, Eq t, Num t) => PS t Revlex -> PS t Revlex +basicSet xs = sPLRR classify [] + where + classify = groupClass . setClassPs $ xs + + +-- -- //////////////////////////////////////////////////////////////////////////////////////////////////// +-- --begin: char Set Normal PS--------------------------------------------------------------------------- +-- charSetNormalPS :: [Poly Rational Revlex] -> [Poly Rational Revlex] +-- charSetNormalPS ps +-- | bs == ps = ps +-- | rs == ps' = bs +-- | rs /= [] = charSetNormalPS (rs++ bs) +-- | otherwise = bs +-- where +-- rs = bsDividePsPS ps' bs +-- bs = basicSet ps +-- ps' = red bs ps +-- --------------------------------------------------------------------------------------------------- +-- bsDividePsPS ::(Num t, Eq t) =>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- bsDividePsPS ps' bs = quitEmptyPoly (psBsPS ps' bs) +-- ---------------------------------------------------------------------------------------------------- +-- psBsPS :: (Num t, Eq t)=>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- psBsPS ps bs = [ force (sprem x bs) | x <- ps] + +-- --end: char Set Normal PS-------------------------------------------------------------------------- +-- --begin: char Set Normal SP------------------------------------------------------------------------- +-- charSetNormalSP :: [Poly Rational Revlex] -> [Poly Rational Revlex] +-- charSetNormalSP ps +-- | bs == ps = ps +-- | rs == ps' = bs +-- | rs /= [] = charSetNormalSP (rs++ bs) +-- | otherwise = bs +-- where +-- rs = bsDividePsSP ps' bs +-- bs = basicSet ps +-- ps' +-- = red bs ps +-- ---------------------------------------------------------------------------------------------------- +-- bsDividePsSP ::(Fractional t, Ord t, Num t, Eq t) =>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- bsDividePsSP ps' bs = quitEmptyPoly (psBsSP ps' bs) +-- ---------------------------------------------------------------------------------------------------- +-- psBsSP :: (Fractional t, Ord t, Num t, Eq t)=>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- psBsSP ps bs = [ force (sspoly x bs) | x <- ps] + +-- ----- *****************************----------------------- +-- begin: char Set parallel division PS +--characteristic set in parallel computation using pseudo remainder +charSetMPS :: (NFData t, Eq t, Num t) => PS t Revlex -> PS t Revlex +charSetMPS ps + | bs == ps = ps + | rs == ps' = bs + | rs /= parps[] = charSetMPS (clean (PS $ A.computeAs N $ on (A.append' 1) f rs bs)) + | otherwise = bs + where + rs = clean $ bsDividePsMPS ps' bs + bs = clean $ basicSet ps + ps' = red bs ps + f(PS poly) = poly; + clean (PS poly) = PS $ A.computeAs N (A.map g poly); + g (Idp poly id)= Idp poly 0 + +-- orden' :: PS t ord -> PS t ord +-- orden' (PS polys) = A.map (f) polys +-- where +-- f (Idp poly id) = g poly () + +-- f(PS poly) = poly; clean (PS poly) = PS $ A.computeAs N (A.map g poly); g (Idp poly id)= Idp poly 0 + +bsDividePsMPS ::(NFData t, Num t, Eq t) => PS t Revlex -> PS t Revlex -> PS t Revlex +bsDividePsMPS ps' bs = quitEmptyPoly (psBsUParP ps' bs) +-- ---------------------------------------------------------------------------------------------------- +psBsUParP :: (NFData t, Num t, Eq t) => PS t Revlex -> PS t Revlex -> PS t Revlex +psBsUParP (PS ps) (bs) = PS $ A.computeAs N $ A.map (`sprem` bs) ps +-- --psBsUParP ps bs = unsafePerformIO ( traverseConcurrently (ParN 4) (\p -> p `deepseq` pure p) [sprem x bs | x <- ps] ) +-- psBsUParP ps bs = unsafePerformIO ( traverseConcurrently Par (pure $!) [sprem x bs | x <- ps] ) +-- --psBsUParP ps bs = unsafePerformIO ( traverseConcurrently (ParN $ fromIntegral (length ps)) (\p -> myThreadId >>= print >> pure p) [sprem x bs | x <- ps] ) +-- laia ps bs = unsafePerformIO ( traverseConcurrently (Par) (\p -> myThreadId >>= print >> pure p) [div x bs | x <- ps] ) +-- -- ahora me pide el nfdata pra poly (sin siquiera haberlo puesto en el bench) +--end : char Ser parallel division PS-------------------------------------------------------------------------------------------------- + +--begin: characteristic set Parallel division Spoly +charSetMSP :: (NFData t, Fractional t , Ord t, Eq t, Num t) => PS t Revlex -> PS t Revlex +charSetMSP ps + | bs == ps = ps + | rs == ps' = bs + | rs /= parps[] = charSetMSP (clean (PS $ A.computeAs N $ on (A.append' 1) f rs bs)) + | otherwise = bs + where + rs = clean $ bsDividePsMSP ps' bs + bs = clean $ basicSet ps + ps' = red bs ps + f(PS poly) = poly; + clean (PS poly) = PS $ A.computeAs N (A.map g poly); + g (Idp poly id)= Idp poly 0 + +bsDividePsMSP ::(NFData t, Fractional t, Ord t, Num t, Eq t) => PS t Revlex -> PS t Revlex -> PS t Revlex +bsDividePsMSP ps' bs = quitEmptyPoly (psBsUParS ps' bs) + +psBsUParS :: (NFData t, Fractional t, Num t, Eq t, Ord t) =>PS t Revlex -> PS t Revlex -> PS t Revlex +psBsUParS (PS ps) (bs) = PS $ A.computeAs N $ A.map (`sspoly` bs) ps + +----- *****************************----------------------- +-- -- charSet sequetial computation +-- --begin: char set sequential PS--------------------------------------------------------------------- +-- charSetMSeqPS :: [Poly Rational Revlex] -> [Poly Rational Revlex] +-- charSetMSeqPS ps +-- | bs == ps = ps +-- | rs == ps' = bs +-- | rs /= [] = charSetMSeqPS (rs ++ bs) +-- | otherwise = bs +-- where +-- rs = bsDividePsMseqPS ps' bs +-- bs = basicSet ps +-- ps' = red bs ps +-- ------------------------------------------------------------------------------------------------- +-- bsDividePsMseqPS ::(Fractional t, Ord t, Num t, Eq t) =>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- bsDividePsMseqPS ps' bs = quitEmptyPoly (psBsSeqPS ps' bs ) +-- --------------------------------------------------------------------------------------------------- +-- psBsSeqPS :: (Num t, Eq t)=>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- psBsSeqPS ps bs = unsafePerformIO ( traverseConcurrently Seq (\p -> p `deepseq` pure p) [sprem x bs | x <- ps] ) +-- --end: char set sequential PS------------------------------------------------------------- +-- --begin: char set sequential SP-------------------------------------------------------------------------------------------------- +-- charSetMSeqSP :: [Poly Rational Revlex] -> [Poly Rational Revlex] +-- charSetMSeqSP ps +-- | bs == ps = ps +-- | rs == ps' = bs +-- | rs /= [] = charSetMSeqSP (rs++ bs) +-- | otherwise = bs +-- where +-- rs = bsDividePsSP ps' bs +-- bs = basicSet ps +-- ps' = red bs ps +-- ----------------------------------------------------------------------------- +-- bsDividePsMseqSP ::(Fractional t, Ord t, Num t, Eq t) =>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- bsDividePsMseqSP ps' bs = quitEmptyPoly (psBsSeqSP ps' bs) +-- ---------------------------------------------------------------------------------------------------- +-- psBsSeqSP :: (Fractional t, Ord t, Num t, Eq t)=>[Poly t Revlex] -> [Poly t Revlex] -> [Poly t Revlex] +-- psBsSeqSP ps bs = unsafePerformIO ( traverseConcurrently Seq (\p -> p `deepseq` pure p) [sspoly x bs | x <- ps] ) +-- --end: charset sequential SP +-- -- ********************************************************************** + +-- -- parallel division pure IO +-- psBsSParP :: [Poly Rational Revlex] -> [Poly Rational Revlex] -> IO[Poly Rational Revlex] +-- psBsSParP ps bs = traverseConcurrently Par (\p -> p `deepseq` pure p) [sprem x bs | x <- ps] :: IO [Poly Rational Revlex] + +-- psBsSParS :: [Poly Rational Revlex] -> [Poly Rational Revlex] -> IO[Poly Rational Revlex] +-- psBsSParS ps bs = traverseConcurrently Par (\p -> p `deepseq` pure p) [sspoly x bs | x <- ps] :: IO [Poly Rational Revlex] + +-- -- sequential division pure IO +-- psBsSSeqP :: [Poly Rational Revlex] -> [Poly Rational Revlex] -> IO[Poly Rational Revlex] +-- psBsSSeqP ps bs = traverseConcurrently Seq (\p -> p `deepseq` pure p) [sprem x bs | x <- ps] :: IO [Poly Rational Revlex] + +-- psBsSSeqS :: [Poly Rational Revlex] -> [Poly Rational Revlex] -> IO[Poly Rational Revlex] +-- psBsSSeqS ps bs = traverseConcurrently Seq (\p -> p `deepseq` pure p) [sspoly x bs | x <- ps] :: IO [Poly Rational Revlex] + +-- psBsSParn1P :: [Poly Rational Revlex] -> [Poly Rational Revlex] -> IO[Poly Rational Revlex] +-- psBsSParn1P ps bs = traverseConcurrently (ParN 1) (\p -> p `deepseq` pure p) [sprem x bs | x <- ps] :: IO [Poly Rational Revlex] + +-- instance NFData (IO[Poly Rational Revlex]) where +-- rnf x = seq x () +-- instance NFData (Poly t Revlex) where +-- rnf x = seq x () +-- -- psByTfPM :: (Num t, Eq t) => [Poly t Revlex] -> [Poly t Revlex]-> IO [Poly t Revlex] +-- -- psByTfPM ps bs = withScheduler Par $ \x -> do +-- -- p <- ps +-- -- scheduleWork x $ pure (sprem p bs ) + +-- --psByTfP :: (Num t, Eq t) => [Poly t Revlex] -> [Poly t Revlex]-> IO [[Poly t Revlex]] +-- --psByTfP ps bs = withScheduler Par $ \x -> scheduleWork x $ pure ([(sprem x) bs | x <- ps]) -- paralelisa? + +-- --psByTfPR :: (Num t, Eq t) => [Poly t Revlex] -> [Poly t Revlex]-> IO [[Poly t Revlex]] +-- --psByTfPR ps bs = withScheduler Par $ \x -> replicateM (length ps) $ scheduleWork x $ pure ([(sprem x) bs | x <- ps]) -- replica one single work + +-- --psByTfPS :: (Num t, Eq t) => [Poly t Revlex] -> [Poly t Revlex]-> IO [Poly t Revlex] +-- --psByTfPS ps bs = withScheduler Par $ \x -> do scheduleWork x $ pure ( sprem (ps !! 0) bs ) +-- -- scheduleWork x $ pure ( sprem (ps !! 1) bs ) -- there is no way we can know the length of a set +-- -- como crear funciones de acuerdo al input + +-- -- laia :: Int -> Int -> IO [Int] +-- -- laia a b= withScheduler Par $ \x -> replicateM_ a $ scheduleWork x $ pure (L.foldl' (P.+) 0 [0 .. b]) -- scheduleId x + +-- -- laia' :: Int -> Int -> IO[Int] +-- -- laia' a b = withScheduler Par $ \x -> do replicateM a $ scheduleWork x $ pure (L.foldl' (P.+) 0 [0 .. b]) -- scheduleId x +-- -- -- scheduleId = (`scheduleWorkId` (\ i -> threadDelay 1000000 >> print i >> pure (L.foldl' (P.+) 0 [0 .. 100000] ))) + +-- -- laia'' :: Int -> IO [Int] +-- -- laia'' a = +-- -- withScheduler Par $ \x -> do +-- -- scheduleWork x $ pure (L.foldl' (P.+) 0 [0 .. div a 2]) +-- -- scheduleWork x $ pure (L.foldl' (P.+) 0 [div a 2 .. a]) +-- ---------------- ************** ----------------------------- +-- ----------------------------------------------------------------------------- +-- --red [] xp = xp +-- --red (x:xs) xp = red xs (delete x xp) +-- ----------------------------------------------------------------------------- +-- -- graphics calculator examples +-- -- p1 = Poly [Term(1,mp[2][2]), Term(-1,m[1,1]), Term(-1,m[0])] :: Poly Rational Revlex +-- -- p2 = Poly [Term(1,m[2]), Term(-2,mp[3][1])]:: Poly Rational Revlex +-- -- p3 = Poly [Term(1,mp[3][2]), Term(-1,m[1,1]), Term(1,m[0])] :: Poly Rational Revlex +-- -- p4 = Poly [Term3,mp[2][5]] :: Poly Rational Revlex +-- -- ps = [Poly [Term(1,mp[2][2]), Term(-1,m[1,1]), Term(-1,m[0])] :: Poly Rational Revlex, Poly [Term(1,m[2]), Term(-2,mp[3][1])]:: Poly Rational Revlex, Poly [Term(1,mp[3][2]), Term(-1,m[1,1]), Term(1,m[0])] :: Poly Rational Revlex] +-- --Poly [Term(3,mp[2][5])] :: Poly Rational Revlex] +-- -- ps1 = [Poly [Term(1,mp[2][2]), Term(-1,m[1,1]), Term(-1,m[0])] :: Poly Rational Revlex, Poly [Term(1,m[2]), Term(-2,mp[3][1])]:: Poly Rational Revlex, Poly [Term(1,mp[3][2]), Term(-1,m[1,1]), Term(1,m[0])] :: Poly Rational Revlex, Poly [Term(3,mp[2][5])] :: Poly Rational Revlex] + +-- -- Buchberger +-- -- f1 = Poly [Term(1,mp[1,4][1,2]), Term(1,mp[4][2]), Term(-1,mp[1,2,4][1,1,1]), Term(-1,mp[2,4][1,1]), Term(1,m[1,1]), Term(3, mp[2][1])] :: Poly Rational Revlex; f2 = Poly [Term(1,mp[1,4][1,1]), Term(1,mp[3][1]), Term(-1,m[1,1])] :: Poly Rational Revlex; f3 = Poly [Term(1,mp[3,4][1,1]), Term(-2,mp[2][2]),Term(-1,m[1,1]),Term(-1,m[0])] :: Poly Rational Revlex; f4 = prem f1 f2 +-- -- asc +-- -- f1 = Poly [Term(1,m[5])] :: Poly Rational Revlex; f2 = Poly [Term(1,m[6]), Term(1,mp[2][1])] :: Poly Rational Revlex; f4 = Poly [Term(1,m[1]),Term(1,mp[2][3])]; f3 = Poly[Term(1,m[2])] :: Poly Rational Revlex +-- -- lower = [Poly[Term(1,m[6])], Poly[Term(4,m[3])], Poly[Term(6,m[5])], Poly[Term(10,m[3])], Poly[Term(1, mp[2][3]), Term(2, m[3])], Poly[Term(5,mp[2][5])], Poly[Term(3, mp[3][2]), Term(4, mp[1][2])], Poly[Term(3, mp[3][1])]] + +-- -- c1 = [Poly[Term(1,m[6])], Poly[Term(4,m[3])], Poly[Term(6,m[5])], Poly[Term(10,m[3])] :: Poly Rational Revlex; c2 = Poly[Term(1, mp[1,2][2,1]), Term(2, m[3])], Poly[Term(5,mp[2][5])] :: Poly Rational Revlex; c3 = Poly[Term(3, mp[3][2]), Term(4, mp[1][2])], Poly[Term(3, mp[3][1])]] :: Poly Rational Revlex +-- -- c1 = [Poly[Term(1,m[6])], Poly[Term(4,m[5])], Poly[Term(6,m[7])], Poly[Term(10,m[3])]] :: [Poly Rational Revlex]; c2 = [Poly[Term(4, mp[1,2][2,1])],Poly[Term(1, mp[1,2][2,1]), Term(2, m[3])], Poly[Term(5,mp[2][5])], Poly[Term(3,mp[1,2][3,4])]] :: [Poly Rational Revlex]; c3 = [Poly[Term(3, mp[3][11]), Term(4, mp[1][2])], Poly[Term(3, mp[3][10])], Poly [Term(5,mp[1,2,3][1,0,7])]] :: [Poly Rational Revlex] +-- -- calculator example +-- -- f1 = Poly[Term(1,mp[2][2]), Term(-1, m[1,1]), Term(-1, m[0])] :: Poly Rational Revlex; f2 = Poly[Term(1,m[2]), Term(-2,mp[3][1])] :: Poly Rational Revlex; f3 = Poly[Term(1,mp[3][2]), Term(-1,m[1,1]), Term(1,m[0])] :: Poly Rational Revlex +-- -- c1 = [Poly[Term(1,m[6])], Poly[Term(4,m[3])], Poly[Term(6,m[5])], Poly[Term(10,m[3])]] :: [Poly Rational Revlex]; c2 = [Poly[Term(1, mp[1,2][2,1]), Term(2, m[3])], Poly[Term(5,mp[2][5])], Poly[Term(3,mp[1,2][1,1])]] :: [Poly Rational Revlex]; c3 = [Poly[Term(3, mp[3][2]), Term(4, mp[1][2])], Poly[Term(3, mp[3][1])]] :: [Poly Rational Revlex]; ps = concat $ c1:c2:c3:[] + + + +-- ----------------------------------------------------------------------------------------- +-- -- s1 = +-- -- Poly [Term (1, mp' [2] [2]), Term(-1, m' [1, 1]), Term (-1, m' [0])] :: Poly Rational Revlex + +-- -- s2 = Poly [Term (1, m' [2]), Term (-2, mp' [3] [1])] :: Poly Rational Revlex + +-- -- s3 = +-- -- Poly [Term (1, mp' [3] [2]), Term (-1, m' [1, 1]), Term(1, m' [0])] :: Poly Rational Revlex +-- -- pallL = [s1,s2,s3] + +f4 = p[t 3 [1], tp (-1) [2][1], t (-7) []] :: Poly Rational Revlex +f5 = p[t 2 [1], tp 3 [2][1], t (-1) []]:: Poly Rational Revlex +pall1= [f4,f5] + +f6 =p[t 1 [1,1], t 1 [1], tp 1 [2][1]] :: Poly Rational Revlex +f7 =p[t 1 [1,2], t 1 [1], tp 1 [2][1]] :: Poly Rational Revlex + +pall2 = parps [f6,f7] + +-- -- http://symbolicdata.org/XMLResources/IntPS/Behnke.xml +-- a1 = Poly [Term(1,m[7]), Term(-1,mp[2][7])] :: Poly Rational Revlex +-- a2 = Poly [Term(1,mp[2][7]), Term(-1,mp[3][7])] :: Poly Rational Revlex +-- a3 = Poly [Term(1,mp[3][7]), Term(-1,mp[4][7])] :: Poly Rational Revlex +-- a4 = Poly [Term(1,mp[4][7]), Term(-1, mp[5][7])] :: Poly Rational Revlex +-- -- a5 = Poly [Term(1,m[6,1]), Term(1, mp[2,3][6,1]),Term(1, mp[3,4][6,1]), Term(1, mp[4,5][6,1]), Term(1, mp[1,5][1,6]) ] :: Poly Rational Revlex +-- ps2= a1 : a2: a3: a4: a5: [] + +-- ps = [Poly[Term(1,m[2])]] :: [Poly Rational Revlex] +-- ps1 = [Poly[Term(1,m[2])]] :: [Poly Rational Revlex] + +-- t1 = Poly [Term(1, m[4])]:: Poly Rational Revlex +-- --t1 = Poly [Term(1, m[4]), Term(1,m[])]:: Poly Rational Revlex +-- t2 = Poly [Term(1,m[2])]::Poly Rational Revlex +-- t3 = Poly [Term(1,mp[2][2]),Term(1, m[2]) ]::Poly Rational Revlex +-- --t3 = Poly [Term(1,mp[1,2][4,2])]::Poly Rational Revlex +-- t4 = Poly [Term(1,mp[3][2])]::Poly Rational Revlex +-- pt = [t1,t2,t3,t4] +-- ******************************************************************************************************************************************************************************************************** +f1 = p [tp 1 [2] [2], t (-1) [1, 1], t (-1) []] :: Poly Rational Revlex +f2 = p [t 1 [2], tp (-2) [3] [1]] :: Poly Rational Revlex +f3 = p [tp 1 [3] [2], t (-1) [1, 1], t 1 []] :: Poly Rational Revlex +ps1 = parps [f2,f1, f3] +---------------------------------------------------------------------------------------------------- +a6 = + p [tp 1 [2] [2], tp 1 [3] [2], tp 1 [4] [2], t (-1) [2]] :: Poly Rational Revlex +a7 = + p [tp 1 [2, 3] [1, 1], tp 1 [4] [2], t (-1) []] :: Poly Rational Revlex +a8 = + p + [ tp 1 [2, 3, 4] [1, 1, 1] + , tp (-1) [2] [2] + , tp (-1) [3] [2] + , tp (-1) [4] [1] + , t 1 [] + ] :: Poly Rational Revlex +ps2' = a6 : a7 : a8 : [] +ps2 = parps ps2' +---------------------------------------------------------------------------------------------------- +a9 = + p + [ tp 1 [1, 10] [2, 1] + , tp 2 [1, 2, 11] [1, 1, 1] + , tp 3 [2, 12] [2, 1] + , tp 1 [1, 4] [1, 1] + , tp 2 [2, 5] [1, 1] + , tp 1 [7] [1] + ] :: Poly Rational Revlex +a10 = + p + [ tp 3 [1, 9] [2, 1] + , tp 2 [1, 2, 10] [1, 1, 1] + , tp 1 [2, 11] [2, 1] + , tp 2 [1, 3] [1, 1] + , tp 1 [2, 4] [1, 1] + , tp 1 [6] [1] + ] :: Poly Rational Revlex +a11 = + p + [ tp 1 [1, 9] [3, 1] + , tp 1 [1, 2, 10] [2, 1, 1] + , tp 1 [1, 2, 11] [1, 2, 1] + , tp 1 [2, 12] [3, 1] + , tp 1 [1, 3] [2, 1] + , tp 1 [1, 2, 4] [1, 1, 1] + , tp 1 [2, 5] [2, 1] + , tp 1 [1, 6] [1, 1] + , tp 1 [2, 7] [1, 1] + , tp 1 [8] [1] + ] :: Poly Rational Revlex +ps3' = a9 : a10 : a11 : [] +ps3 = parps ps3' +---------------------------------------------------------------------------------------------------- +a12 = p [ tp 1 [1,2,3][2,4,1], tp 1 [1,2][2,4], tp 2 [1,2,3][2,2,1], tp 6 [1,2][2,2], tp (-4) [1,2][1,3], tp 1 [2,3][4,1], tp (-1) [2][4], tp 1 [1,3][2,1], tp 1 [1][2], tp 4 [1,2][1,1], tp 2 [2,3][2,1], tp (-6) [2][2], tp 1 [3][1], t (-1) [] + ] :: Poly Rational Revlex +a13 = p [ tp 1 [1,2,3][8,3,1], tp 1 [1,2,3][8,1,1], tp 4 [1,2,3][6,3,1], tp (-8) [1,2][6,3], tp 8 [1,2][5,4], tp 4 [1,2,3][6,1,1], tp 8 [1,2][6,1], tp (-48) [1,2][5,2], tp 6 [1,2,3][4,3,1], tp 48 [1,2][4,3], tp (-8) [1,2][3,4], tp 8 [1][5], tp 6 [1,2,3][4,1,1], tp (-48) [1,2][4,1], tp (48) [1,2][3,2], tp 4 [1,2,3][2,3,1], tp (-8) [1,2][2,3], tp (-8) [1][3], tp 4 [1,2,3][2,1,1], tp (8) [1,2][2,1], tp (1) [2,3][3,1], tp 1 [2,3][1,1] + ] :: Poly Rational Revlex +ps4 = parps [a12,a13] +---------------------------------------------------------------------------------------------------- +a14 = p[ tp 35 [2][1], tp 40 [3][1], tp 25 [4][1], tp (-27) [5][1]] :: Poly Rational Revlex +a15 = p[ tp 45 [2][1], tp 35 [5][1], tp (-165) [6][1], t (-36) []] :: Poly Rational Revlex +a16 = p[ tp (-11)[5,6][1,1], tp (3)[6][2],tp (99) [1][1] ] :: Poly Rational Revlex +a17 = p[ tp 25 [2,5][1,1], tp (-165)[6][2], tp (15) [1][1], tp 30 [3][1], tp (-18) [4][1] ] :: Poly Rational Revlex +a18 = p[ tp 15 [2,4][1,1], tp 20 [3,5][1,1], tp (-9) [1][1]] :: Poly Rational Revlex +a19 = p[ tp (-11)[6][3], tp 1 [1,2][1,1], tp 2 [3,4][1,1] ] :: Poly Rational Revlex +ps5 = parps [a14,a15,a16,a17,a18,a19] +---------------------------------------------------------------------------------------------------- +a20 = p[t (-4) [1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1]] :: Poly Rational Revlex +a21 = p[tp (-4) [1][2], tp 1 [2][2], tp 1 [3][2], tp 1 [4][2], tp 4 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], t (-3) []] :: Poly Rational Revlex +a22 = p[tp (5) [1][3], tp 4 [3,4][2,1], tp 3 [2][2], tp 2 [1,4][1,1], tp 4 [1][1], tp 1 [2][1], tp 1 [3][1], tp 2 [4][1], t (-1) []] :: Poly Rational Revlex +a23 = p[tp 5 [3][4], tp 1 [4][3], tp 16 [1][2], tp 3 [2][2], tp (-4) [4][1], t (-1) []] :: Poly Rational Revlex +ps6 = parps [a20,a21,a22,a23] +---------------------------------------------------------------------------------------------------- +a24 = p[tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1]] :: Poly Rational Revlex +a25 = p[tp 1 [1,2][1,1], tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [1,5][1,1], tp 1 [4,5][1,1]] :: Poly Rational Revlex +a26 = p[tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1 [1,2,5][1,1,1], tp 1 [1,4,5][1,1,1], tp 1 [3,4,5][1,1,1]] :: Poly Rational Revlex +a27 = p[tp 1 [1,2,3,4][1,1,1,1], tp 1 [1,2,3,5][1,1,1,1], tp 1 [1,2,4,5][1,1,1,1], tp 1 [1,3,4,5][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1]] :: Poly Rational Revlex +a28 = p[t 1 [1,1,1,1,1], t (-1) []] :: Poly Rational Revlex +ps7 = parps [a24,a25,a26, a27, a28] +---------------------------------------------------------------------------------------------------- +a29 = p[tp (-2) [1,3][1,1], tp (-2) [3][2], tp (-2)[1][1], tp (8) [3][1], t (-2) []]:: Poly Rational Revlex +a30 = p[tp (-3)[1,2,3][1,1,1], tp (2) [1,3,4][1,1,1], tp (4)[3,4][2,1], tp (3)[2,3][1,1], tp 1 [1,4][1,1], tp (-7) [3,4][1,1]]:: Poly Rational Revlex +a31 = p[tp 1 [1,2][2,2], tp (-2)[1,2,4][2,1,1], t (-2) [1,1,1,1], tp 1 [1,4][2,2], tp 2 [1,3,4][1,1,2], tp 1 [3,4][2,2], tp (-2)[1,2][1,2], t (4) [1,1,1], tp 2 [1,2,4][1,1,1], tp 2 [2,3,4][1,1,1], tp (-4)[1,4][1,2], tp (-4) [3,4][1,2], tp 1 [2][2], tp 2 [1,3][1,1], tp 10 [3][2], tp (-4)[2,4][1,1], tp 4 [4][2], tp (-2)[1][1], tp (-8)[3][1], t 2 [] ]:: Poly Rational Revlex +a32 = p[t 3 [2,2], t 12 [1,1,1,1], tp (-3) [1,4][2,2], tp 6 [1,3,4][1,1,2], tp (-3) [3,4][2,2], tp (-6) [1,2][1,2], tp 12 [1,2,4][1,1,1], tp 12 [2,3,4][1,1,1], tp (-4) [1][2], tp 3 [2][2], tp 5 [3][2], tp (-12) [2,4][1,1], tp 12 [4][2], tp (-6) [3][1], t 5 [] ] :: Poly Rational Revlex +ps8 = parps [a29,a30,a31,a32] +---------------------------------------------------------------------------------------------------- +-- ojito + +---------------------------------------------------------------------------------------------------- +a39 = p[tp 1 [2,3][1,5], t 1 [1], t (-2)[] ]:: Poly Rational Revlex +a40 = p[tp 1 [1,3][5,1], tp 1 [2][1], t (-2)[]]:: Poly Rational Revlex +a41 = p[t 1 [1,5], tp 1 [3][1], t (-2) []]:: Poly Rational Revlex +ps11 = parps [a39,a40,a41] +---------------------------------------------------------------------------------------------------- +a46 = p[tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1], tp 1 [6][1]] :: Poly Rational Revlex +a47 = p[tp 1 [1,2][1,1], tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [4,5][1,1], tp 1 [1,6][1,1], tp 1 [5,6][1,1]] :: Poly Rational Revlex +a48 = p[tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1[3,4,5][1,1,1], tp 1 [1,2,6][1,1,1], tp 1 [1,5,6][1,1,1], tp 1 [4,5,6][1,1,1]] :: Poly Rational Revlex +a49 = p[tp 1 [1,2,3,4][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1], tp 1 [1,2,3,6][1,1,1,1], tp 1 [1,2,5,6][1,1,1,1], tp 1 [1,4,5,6][1,1,1,1], tp 1 [3,4,5,6][1,1,1,1]] :: Poly Rational Revlex +a50 = p[tp 1 [1,2,3,4,5][1,1,1,1,1], tp 1 [1,2,3,4,6][1,1,1,1,1], tp 1 [1,2,3,5,6][1,1,1,1,1], tp 1 [1,2,4,5,6][1,1,1,1,1], tp 1 [1,3,4,5,6][1,1,1,1,1], tp 1 [2,3,4,5,6][1,1,1,1,1]] :: Poly Rational Revlex +a51 = p[t 1 [1,1,1,1,1,1], t (-1) [] ] :: Poly Rational Revlex +ps13 = parps [a46,a47,a48,a49,a50,a51] +---------------------------------------------------------------------------------------------------- +a52 = p [t 1 [1,1],tp 1 [2,3][1,1], tp 1 [3,4][1,1], tp 1 [1,5][1,1], tp 1 [4,5][1,1]] :: Poly Rational Revlex +a53 = p [tp 1 [1][1], tp 1 [2][1], tp 1 [3][1], tp 1 [4][1], tp 1 [5][1]] :: Poly Rational Revlex +a54 = p [tp 1 [1,2,3][1,1,1], tp 1 [2,3,4][1,1,1], tp 1 [1,2,5][1,1,1], tp 1 [1,4,5][1,1,1], tp 1 [3,4,5][1,1,1]] :: Poly Rational Revlex +a55 = p [t 1 [1,1,1,1], tp 1 [1,2,4,5][1,1,1,1], tp 1 [1,3,4,5][1,1,1,1], tp 1 [2,3,4,5][1,1,1,1], tp 1 [2,3,4][1,1,1] ] :: Poly Rational Revlex +a56 = p [ t 1 [1,1,1,1,1], t (-1) []] :: Poly Rational Revlex +ps14 = parps [a52, a53,a54,a55,a56] +---------------------------------------------------------------------------------------------------- +a57 = p [tp 1 [1,2,3][2,1,1], tp 1 [1,2,3][1,2,1], tp 1 [1,2,3][1,1,2], tp 1 [1,1,1][1,1,1], tp 1 [1,2][1,1], tp 1 [1,3][1,1], tp 1 [2,3][1,1]] :: Poly Rational Revlex +a58 = p [tp 1 [1,2,3][2,2,1], tp 1 [1,2,3][1,2,2], tp 1 [1,2,3][2,1,1], t 1 [1,1,1], tp 1 [2,3][1,1], t 1 [1], tp 1 [3][1] ] :: Poly Rational Revlex +a59 = p [t 1 [2,2,2], t 1 [2,2,1], t 1 [1,2,1], t 1 [1,1,1], tp 1 [1,3][1,1], tp 1 [3][1], t 1 []] :: Poly Rational Revlex +ps15 = parps [a57, a58, a59] +---------------------------------------------------------------------------------------------------- + + +-- class' XE take max +-- totalDeg XE +-- number of poly + +--ld XE + +-- numTerms + +helper :: (NFData t) => PS t Revlex -> [Poly t Revlex] +helper (PS b) = L.map f (A.toList b ) + where + f (Idp poly x) = poly + + +tdp :: (NFData t, Eq t, Num t) =>PS t Revlex -> Int -> [Int] +tdp xs s = L.map (\x -> varDegree x s) $ helper (charSetMPS xs) + +tdp' :: (NFData t, Eq t, Num t, Fractional t, Ord t) =>PS t Revlex -> Int -> [Int] +tdp' xs s = L.map (\x -> varDegree x s) $ helper (charSetMSP xs) + +-- jose :: (NFData t, Eq t, Num t, Fractional t, Ord t) => PS t Revlex -> Int -> [[Int]] +-- jose xs x +-- | x /= 1 = ((tdp xs x) : jose xs (x P.-1) +-- | otherwise = ((tdp xs x)) :[] + + +cls ::(NFData t, Eq t, Num t) => PS t Revlex -> [Int] +cls p = L.map class' $ helper ( p) + +std :: (NFData t, Eq t, Num t) =>PS t Revlex -> [Int] +std p = L.map totalDeg $ helper ( p) + +np :: (NFData t, Eq t, Num t) =>PS t Revlex -> Int +np p = length $ helper (charSetMPS p) +np' :: (NFData t, Eq t, Num t, Fractional t, Ord t) =>PS t Revlex -> Int +np' p = length $ helper (charSetMSP p) +--------------------------------------------------------- +-- lde :: (NFData t, Eq t, Num t) =>PS t Revlex -> [Int] +-- lde p = L.map (\x-> varDegree x 1) $ helper (charSetMPS p) +-- lde' :: (NFData t, Eq t, Num t, Fractional t, Ord t) =>PS t Revlex -> [Int] +--lde' p = L.map (\x-> varDegree x 1) $ helper (charSetMSP p) + +nt :: (NFData t, Eq t, Num t) => PS t Revlex -> [Int] +nt p = L.map numTerms $ helper (charSetMPS p) + +nt' :: (NFData t, Eq t, Num t, Fractional t, Ord t) => PS t Revlex -> [Int] +nt' p = L.map numTerms $ helper (charSetMSP p) + + +-- cls (mas alto), td, np; class' cambiar el valor; nt +a = p[tp 2 [2][3], tp (-1) [2][2], t 1 [2,1]] :: Poly Rational Revlex +b = p [t 1 [1,2], t 1 []] :: Poly Rational Revlex + + +c = p[tp 1 [1,4][1,2], tp 1 [4][2], tp (-1) [1,2,4][1,1,1], tp (-1) [2,4][1,1], tp 1 [1,2][1,1], tp 3 [2][1] ] :: Poly Rational Revlex +d = p[tp 1 [1,4][1,1], tp 1 [3][1], tp (-1) [1,2][1,1]] :: Poly Rational Revlex +e = p[tp 1 [3,4][1,1], tp (-2) [2][2], tp (-1) [1,2][1,1], t (-1)[]] :: Poly Rational Revlex +f = parps $ c:e:d:[] diff --git a/test/Spec.hs b/test/Spec.hs new file mode 100644 --- /dev/null +++ b/test/Spec.hs @@ -0,0 +1,31 @@ +{-# LANGUAGE ScopedTypeVariables #-} + +module Main where + +import Test.Tasty.QuickCheck as QC +import Test.Tasty +import Lib + +reverseTests :: TestTree +reverseTests = testGroup "Tests over reverse" + [QC.testProperty "reverse respect lenght" $ + \(lst::[Integer])-> length (reverse' lst) == length lst] + + + + +-- import Test.Tasty +-- import Test.Tasty.HUnit as HU +-- import Lib + +-- hunitTestInsertOnLeaf :: TestTree +-- hunitTestInsertOnLeaf = HU.testCase "Insert 'a' on empty tree" $ +-- assertEqual "Insertion is wrong" (treeInsert 'a' Leaf)(Node 'a' Leaf Leaf) + +-- allTests :: TestTree +-- allTests = testGroup "Tasty Test" [ +-- testGroup "Hunit Test" [hunitTestInsertOnLeaf] +-- ] + +main :: IO() +main = defaultMain reverseTests