image image image image image image image

Rings: efficient Java/Scala library for polynomial rings

Rings is an efficient implementation of univariate and multivariate polynomial algebra over arbitrary coefficient rings. It makes use of asymptotically fast algorithms for basic algebraic operations as well as for advanced methods like GCDs and polynomial factorization. Performance achieved in Rings is comparable to well known solutions like Singular/NTL/FLINT/Maple/Mathematica.

The key features of Rings include:

The full documentation is available at http://rings.readthedocs.io.

Set up

Interactive Rings shell

To taste what Rings can do, one can try interactive session with Ammonite REPL. You can install Rings.repl with Homebrew:

$ brew install PoslavskySV/rings/rings.repl

or just by typing the following commands at the prompt:

$ sudo curl -L -o /usr/local/bin/amm https://git.io/v5Tct && sudo chmod +x /usr/local/bin/amm
$ sudo curl -L -o /usr/local/bin/rings.repl https://git.io/vd7EY && chmod +x /usr/local/bin/rings.repl

Now run Rings.repl:

$ rings.repl
Loading...
Rings 2.2: efficient Java/Scala library for polynomial rings

@ implicit val ring = MultivariateRing(Z, Array("x", "y", "z"))
ring: MultivariateRing[IntZ] = MultivariateRing(Z, Array("x", "y", "z"), LEX)

@ val poly1 = ring("x + y - z").pow(8) 
poly1: MultivariatePolynomial[IntZ] = z^8-8*y*z^7+28*y^2*z^6-56*y^3*z^5+70*...

@ val poly2 = ring("x - y + z").pow(8) 
poly1: MultivariatePolynomial[IntZ] = z^8-8*y*z^7+28*y^2*z^6-56*y^3*z^5+70*...

@ Factor(poly1 - poly2)
res13: FactorDecomposition[MultivariatePolynomial[IntZ]] = 
       16*(x)*((-1)*z+y)
       *(z^4-4*y*z^3+6*y^2*z^2-4*y^3*z+y^4+6*x^2*z^2-12*x^2*y*z+6*x^2*y^2+x^4)
       *(z^2-2*y*z+y^2+x^2)

Java/Scala library

Rings are currently available for Java and Scala. To get started with Scala SBT, simply add the following dependence to your build.sbt file:

libraryDependencies += "cc.redberry" % "rings.scaladsl" % "2.2"

For using Rings solely in Java there is Maven artifact:

<dependency>
    <groupId>cc.redberry</groupId>
    <artifactId>rings</artifactId>
    <version>2.2</version>
</dependency>

Examples: algebra, GCDs, factorization, programming

Below examples can be evaluated directly in the Rings.repl. If using Rings in Scala, the following preambula will import all required things from Rings library:

import cc.redberry.rings

import rings.poly.PolynomialMethods._
import rings.scaladsl._
import syntax._

Java examples can be found in the complete documentation pages.


Do some algebra in Galois field GF(17^9):

// GF(17^9) (irreducible poly in Z/17[x] will be generated automaticaly)
implicit val ring = GF(17, 9, "x")

// some random element from ring
val a = ring.randomElement()
val b = a.pow(1000)
val c = 1 / b

assert ( b * c === 1)

// explicitly parse ring element from string
val d = ring("1 + x + x^2 + x^3 + 15*x^999")
// do some math ops
val some = a / (b + c) + a.pow(6) - a * b * c * d

Some math with multivariate polynomials from Z[x, y, z]:

// Z[x, y, z]
implicit val ring = MultivariateRing(Z, Array("x", "y", "z")) 

val (x, y, z) = ring("x", "y", "z") 

// do some math
val a = (x + y + z).pow(2) - 1 
val b = (x - y - z - 1).pow(2) + x + y + z - 1 
val c = (a + b + 1).pow(9) - a - b - 1

// reduce c modulo a and b (multivariate division with remainder)
val (div1, div2, rem) = c /%/% (a, b)

Univariate extended GCD in Z_{17}[x]:

// ring Z/17[x]
implicit val ring = UnivariateRingZp64(17, "x")

val x = ring("x")

val poly1 = 1 + x + x.pow(2) + x.pow(3)
val poly2 = 1 + 2 * x + 9 * x.pow(2)
val (gcd, s, t) = PolynomialExtendedGCD(poly1, poly2) match { case Array(gcd,s,t) => (gcd,s,t) }

println((gcd, s, t))

Multivariate GCD in Z[a, b, c]:

// ring Z[a, b, c]
implicit val ring = MultivariateRing(Z, Array("a", "b", "c"))

val poly1 = ring("-b-b*c-b^2+a+a*c+a^2")
val poly2 = ring("b^2+b^2*c+b^3+a*b^2+a^2+a^2*c+a^2*b+a^3")

val gcd   = PolynomialGCD(poly1, poly2)

println(s"gcd: ${ring show gcd}")

Factor polynomial in Z_{17}[x]:

// ring Z/17[x]
implicit val ring = UnivariateRingZp64(17, "x")x

val poly = ring("4 + 8*x + 12*x^2 + 5*x^5 - x^6 + 10*x^7 + x^8")

// factorize poly
val factors = Factor(poly)

println(factors)

Coefficient rings with arbitrary large characteristic are available:

// coefficient ring Z/1237940039285380274899124357 (the next prime to 2^100)
val modulus = Z("1267650600228229401496703205653")
val cfRing  = Zp(modulus)

// ring Z/1237940039285380274899124357[x]
implicit val ring = UnivariateRing(cfRing, "x")

val poly = ring("4 + 8*x + 12*x^2 + 5*x^5 + 16*x^6 + 27*x^7 + 18*x^8")

println(Factor(poly))

(large primes can be generated with BigPrimes.nextPrime method, see Prime numbers).


Factor polynomial in Z_2[x, y, z]:

// ring Z/2[x, y, z]
implicit val ring = MultivariateRingZp64(2, Array("x", "y", "z"))

val (x, y, z) = ring("x", "y", "z")

val factors = Factor(1 + (1 + x + y + z).pow(2) + (x + y + z).pow(4))

println(factors)

Factor polynomial in Z[a, b, c]:

// ring Z[a, b, c]
implicit val ring = MultivariateRing(Z, Array("a", "b", "c"))

val (a, b, c) = ring("a", "b", "c")

val factors = Factor(1 - (1 + a + b + c).pow(2) - (2 + a + b + c).pow(3))

println(ring show factors)

Factor polynomial in Q[x, y, z]:

// ring Q[x, y, z]
implicit val ring = MultivariateRing(Q, Array("x", "y", "z"))

val poly = ring(
  """
    |(1/6)*y*z + (1/6)*y^3*z^2 - (1/2)*y^6*z^5 - (1/2)*y^8*z^6
    |-(1/3)*x*z - (1/3)*x*y^2*z^2 + x*y^5*z^5 + x*y^7*z^6
    |+(1/9)*x^2*y^2*z - (1/3)*x^2*y^7*z^5 - (2/9)*x^3*y*z
    |+(2/3)*x^3*y^6*z^5 - (1/2)*x^6*y - (1/2)*x^6*y^3*z
    |+x^7 + x^7*y^2*z - (1/3)*x^8*y^2 + (2/3)*x^9*y
  """.stripMargin)

val factors = Factor(poly)

println(factors)

Polynomial rings over Z and Q:

// Ring Z[x]
UnivariateRing(Z, "x")
// Ring Z[x, y, z]
MultivariateRing(Z, Array("x", "y", "z"))
// Ring Q[a, b, c]
MultivariateRing(Q, Array("a", "b", "c"))

Polynomial rings over Z_p:

// Ring Z/3[x]
UnivariateRingZp64(3, "x")
// Ring Z/3[x, y, z]
MultivariateRingZp64(3, Array("x", "y", "z"))
// Ring Z/p[x, y, z] with p = 2^107 - 1 (Mersenne prime)
MultivariateRing(Zp(Z(2).pow(107) - 1), Array("x", "y", "z"))

Galois fields:

// Galois field with cardinality 7^10 
// (irreducible polynomial will be generated automatically)
GF(7, 10, "x")
// GF(7^3) generated by irreducible polynomial "1 + 3*z + z^2 + z^3"
GF(UnivariateRingZp64(7, "z")("1 + 3*z + z^2 + z^3"), "z")

Fractional fields:

// Field of fractions of univariate polynomials Z[x]
Frac(UnivariateRing(Z, "x"))
// Field of fractions of multivariate polynomials Z/19[x, y, z]
Frac(MultivariateRingZp64(19, Array("x", "y", "z")))

Ring of univariate polynomials over elements of Galois field GF(7^3)[x]:

// Elements of GF(7^3) are represented as polynomials
// over "z" modulo irreducible polynomial "1 + 3*z + z^2 + z^3"
val cfRing = GF(UnivariateRingZp64(7, "z")("1 + 3*z + z^2 + z^3"), "z")

assert(cfRing.characteristic().intValue() == 7)
assert(cfRing.cardinality().intValue() == 343)

// Ring GF(7^3)[x]
implicit val ring = UnivariateRing(cfRing, "x")

// Coefficients of polynomials in GF(7^3)[x] are elements of GF(7^3)
val poly = ring("1 - (1 - z^3) * x^6 + (1 - 2*z) * x^33 + x^66")

// factorize poly (in this examples there will be 9 factors)
val factors = Factor(poly)
println(s"${ring show factors}")

Ring of multivariate polynomials over elements of Galois field GF(7^3)[x, y, z]:

// Elements of GF(7^3) are represented as polynomials
// over "z" modulo irreducible polynomial "1 + 3*z + z^2 + z^3"
val cfRing = GF(UnivariateRingZp64(7, "z")("1 + 3*z + z^2 + z^3"), "z")
// Ring GF(7^3)[x]
implicit val ring = MultivariateRing(cfRing, Array("a", "b", "c"))

// Coefficients of polynomials in GF(7^3)[x] are elements of GF(7^3)
val poly = ring("1 - (1 - z^3) * a^6*b + (1 - 2*z) * c^33 + a^66")

Implement generic function for solving linear Diophantine equations:

/**
  * Solves equation \sum f_i s_i  = gcd(f_1, \dots, f_N) for given f_i and unknown s_i
  * @return a tuple (gcd, solution)
  */
def solveDiophantine[E](fi: Seq[E])(implicit ring: Ring[E]) =
  fi.foldLeft((ring(0), Seq.empty[E])) { case ((gcd, seq), f) =>
    val xgcd = ring.extendedGCD(gcd, f)
    (xgcd(0), seq.map(_ * xgcd(1)) :+ xgcd(2))
  }

Implement generic function for computing partial fraction decomposition:

/** Computes partial fraction decomposition of given rational */
def apart[E](frac: Rational[E]) = {
  implicit val ring: Ring[E] = frac.ring
  val factors = ring.factor(frac.denominator).map {case (f, exp) => f.pow(exp)}
  val (gcd,  nums) = solveDiophantine(factors.map(frac.denominator / _))
  val (ints, rats) = (nums zip factors)
    .map { case (num, den) => Rational(frac.numerator * num, den * gcd) }
    .flatMap(_.normal)       // extract integral parts from fractions
    .partition(_.isIntegral) // separate integrals and fractions
  rats :+ ints.foldLeft(Rational(ring(0)))(_ + _)
}

Apply that function to elements of different rings:

// partial fraction decomposition for rationals
// gives List(184/479, (-10)/13, 1/8, (-10)/47, 1)
val qFracs = apart( Q("1234213 / 2341352"))

// partial fraction decomposition for rational functions
val ufRing = Frac(UnivariateRingZp64(17, "x"))
// gives List(4/(16+x), 1/(10+x), 15/(1+x), (14*x)/(15+7*x+x^2))
val pFracs = apart( ufRing("1 / (3 - 3*x^2 - x^3 + x^5)") )

Implement Lagrange method for univariate interpolation:

    p(x) = \sum_i p(x_i) \Pi_{j \ne i} \frac{x_{\phantom{i}} - x_j}{x_i -x_j}
/** Lagrange polynomial interpolation formula */
def interpolate[Poly <: IUnivariatePolynomial[Poly], Coef]
    (points: Seq[(Coef, Coef)])
    (implicit ring: IUnivariateRing[Poly, Coef]) = {
      // implicit coefficient ring (setups algebraic operators on type Coef)
      implicit val cfRing: Ring[Coef] = ring.cfRing
      if (!cfRing.isField) throw new IllegalArgumentException
      points.indices
        .foldLeft(ring(0)) { case (sum, i) =>
          sum + points.indices
            .filter(_ != i)
            .foldLeft(ring(points(i)._2)) { case (product, j) =>
              product * (ring.`x` - points(j)._1) / (points(i)._1 - points(j)._1)
            }
        }
    }

Interpolate polynomial from Frac(Z_{13}[a,b,c])[x]:

// coefficient ring Frac(Z/13[a,b,c])
val cfRing = Frac(MultivariateRingZp64(2, Array("a", "b", "c")))
val (a, b, c) = cfRing("a", "b", "c")

implicit val ring = UnivariateRing(cfRing, "x")
// interpolate with Lagrange formula
val data = Seq(a -> b, b -> c, c -> a)
val poly = interpolate(data)
assert(data.forall { case (p, v) => poly.eval(p) == v })

Highlighted benchmarks

In the following plots performance of Rings is compared to Wolfram Mathematica 11. All tests were performed on MacBook Pro (15-inch, 2017), 3,1 GHz Intel Core i7, 16 GB 2133 MHz LPDDR3. The code of benchmarks can be found at GitHub. In all benchamrks random polynomials were used.

Rings vs Singular performance of gcd(a g, b g) for random polynomials (a, b, g) \in Z[x_1,x_2,x_3,x_4,x_5] each with 40 terms and degree 20 in each variable

Rings vs Mathematica performance of gcd(a g, b g) for random polynomials (a, b, g) \in Z[x_1,x_2,x_3,x_4,x_5] each with 40 terms and degree 20 in each variable

Rings vs Singular performance of factor(a b c) for random polynomials (a, b, c) \in Z_2[x_1,x_2,x_3,x_4,x_5,x_6,x_7] each with 20 terms and degree 10 in each variable

Rings vs Mathematica performance of factor(a b c) for random polynomials (a, b, c) \in Z[x,y,z] each with 20 terms and degree 10 in each variable

Univariate factorization performance on polynomials of the form (1 + \sum_{i = 1}^{i \leq deg} i \times x^i) in Z_{17}[x]. At small degrees the performance is identical, while at large degrees NTL and FLINT have much better asymptotic (probably due to more advanced algorithms for polynomial multiplication).

Index of algorithms implemented in Rings

Univariate polynomials

  1. Karatsuba multiplication (Sec. 8.1 in [GaGe03])

    used with some adaptations for multiplication of univariate polynomials:

  2. Half-GCD and Extended Half-GCD (Sec. 11 in [GaGe03])

    used (with adaptations inspired by [ShoNTL]) for univariate GCD:

  3. Subresultant polynomial remainder sequences (Sec. 7.3 in [GeCL92]):

  4. Modular GCD in Z[x] and Q[x] (Sec. 6.7 in [GaGe03], small primes version):

  5. Fast univariate division with Newton iteration (Sec. 9.1 in [GaGe03])

    used everywhere where multiple divisions (remainders) by the same divider are performed:

  6. Univariate square-free factorization in zero characteristic (Yun's algorithm) (Sec. 14.6 in [GaGe03]):

  7. Univariate square-free factorization in non-zero characteristic (Musser's algorithm) (Sec. 8.3 in [GeCL92], [Muss71]):

  8. Distinct-degree factorization (Sec. 14.2 in [GaGe03])

    plain version and adapted version with precomputed x-powers (used by default):

  9. Shoup's baby-step giant-step algorithm for distinct-degree factorization ([Shou95])

    used for factorization over fields with large cardinality:

  10. Univariate modular composition

    plain algorithm with Horner schema:

  11. Brent-Kung univariate modular composition ([BreK98], [Shou95]):

  12. Cantor-Zassenhaus algorithm (equal-degree splitting) (Sec. 14.3 in [GaGe03])

    both for odd and even characteristic:

  13. Univaraite linear p-adic Hensel lifting (Sec. 6.5 in [GeCL92]):

  14. Univaraite quadratic p-adic Hensel lifting (Sec. 15.4-15.5 in [GaGe03]):

  15. Univariate polynomial factorization over finite fields

    uses Musser's square free factorization followed by distinct-degree factorization (either x-powers or Shoup's algorithm) followed by Cantor-Zassenhaus equal-degree factorization:

  16. Univariate polynomial factorization over Z and Q

    uses factorization modulo small prime followed by Hensel lifting (adaptive linear/quadratic) and naive recombination:

  17. Univariate irreducibility test (Sec. 14.9 in [GaGe03]):

  18. Ben-Or’s generation of irreducible polynomials (Sec. 14.9 in [GaGe03]):

  19. Univariate polynomial interpolation

    Lagrange and Newton methods:

Multivariate polynomials

  1. Brown GCD over finite fields ([Brow71], Sec. 7.4 in [GeCL92], [Yang09]):

  2. Zippel's sparse GCD over finite fields ([Zipp79], [Zipp93], [dKMW05], [Yang09])

    both for monic (with fast Vandermonde systems) and non-monic (LINZIP) cases:

  3. Extended Zassenhaus GCD (EZ-GCD) over finite fields (Sec. 7.6 in [GeCL92], [MosY73]):

  4. Enhanced Extended Zassenhaus GCD (EEZ-GCD) over finite fields ([Wang80]):

  5. Modular GCD over Z with sparse interpolation ([Zipp79], [Zipp93], [dKMW05])

    (the same interpolation techniques as in ZippelGCD is used):

  6. Modular GCD over Z (small primes version):

  7. Kaltofen's & Monagan's generic modular GCD with sparse interpolation ([KalM99])

    used for computing multivariate GCD over finite fields of very small cardinality:

  8. Kaltofen's & Monagan's generic modular GCD with EEZ-GCD for modular images ([KalM99])

    used for computing multivariate GCD over finite fields of very small cardinality:

  9. Multivariate square-free factorization in zero characteristic (Yun's algorithm) ([LeeM13]):

  10. Multivariate square-free factorization in non-zero characteristic (Musser's algorithm) ([Muss71], Sec. 8.3 in [GeCL92]):

  11. Bernardin's fast dense multivariate Hensel lifting ([Bern99], [LeeM13])

    both for bivariate case (original Bernardin's paper) and multivariate case (Lee thesis) and both with and without precomputed leading coefficients:

  12. Sparse Hensel lifting ([Kalt85], [LeeM13])

  13. Fast dense bivariate factorization with recombination ([Bern99], [LeeM13]):

  14. Kaltofen's multivariate factorization in finite fields ([Kalt85], [LeeM13])

    modified version of original Kaltofen's algorithm for leading coefficient precomputation with square-free decomposition (instead of distinct variables decomposition) due to Lee is used; further adaptations are made to work in finite fields of very small cardinality; the resulting algorithm is close to [LeeM13], but at the same time has many differences in details:

  15. Kaltofen's multivariate factorization Z ([Kalt85], [LeeM13])

    (with the same modifications as for algorithm for finite fields):

  16. Multivariate polynomial interpolation with Newton method:

References

  • [GaGe03] J von zur Gathen and J Gerhard. Modern computer algebra (2 ed.). Cambridge University Press, 2003.
  • [ShoNTL] V Shoup. NTL: A library for doing number theory. www.shoup.net/ntl
  • [GeCL92] K O Geddes, S R Czapor, G Labahn. Algorithms for Computer Algebra. 1992.
  • [Muss71] D R Musser. Algorithms for polynomial factorization, Ph.D. Thesis, University of Wisconsin, 1971.
  • [Shou95] V Shoup. A new polynomial factorization algorithm and its implementation. J. Symb. Comput., 20(4):363–397, 1995.
  • [BreK98] R P Brent and H T Kung. Fast algorithms for manipulating formal power series. J. Assoc. Comput. Math. 25:581-595, 1978
  • [Brow71] W S Brown. On Euclid’s algorithm and the computation of polynomial greatest common divisors. J. ACM, 18(4):478–504, 1971.
  • [Zipp79] R E Zippel. Probabilistic algorithms for sparse polynomials. In Proceedings of the International Symposiumon on Symbolic and Algebraic Computation, EUROSAM ’79, pages 216–226, London, UK, UK, 1979. ringer-Verlag.
  • [Zipp93] R E Zippel. Effective Polynomial Computation. Kluwer International Series in Engineering and Computer Science. Kluwer Academic Publishers, 1993.
  • [dKMW05] J de Kleine, M B Monagan, A D Wittkopf. Algorithms for the Non-monic Case of the Sparse Modular GCD Algorithm. Proceeding of ISSAC ’05, ACM Press, pp. 124-131 , 2005.
  • [Yang09] S Yang. Computing the greatest common divisor of multivariate polynomials over finite fields. Master’s thesis, Simon Fraser University, 2009.
  • [MosY73] J Moses and D Y Y Yun, "The EZGCD Algorithm," pp. 159-166 in Proc. ACM Annual Conference, (1973).
  • [Wang80] P S Wang, "The EEZ-GCD Algorithm," ACM SIGSAMBull., 14 pp. 50-60 (1980).
  • [KalM99] E Kaltofen, M. B. Monagan. On the Genericity of the Modular Polynomial GCD Algorithm. Proceeding of ISSAC ’99, ACM Press, 59-66, 1999.
  • [Bern99] L Bernardin. Factorization of Multivariate Polynomials over Finite Fields. PhD thesis, ETH Zurich, 1999.
  • [LeeM13] M M-D Lee, Factorization of multivariate polynomials, Ph.D. thesis, University of Kaiserslautern, 2013
  • [Kalt85] E Kaltofen. Sparse Hensel lifting. In EUROCAL 85 European Conf. Comput. Algebra Proc. Vol. 2, pages 4–17, 1985.

License

Apache License, Version 2.0 http://www.apache.org/licenses/LICENSE-2.