PolynomialDomain.java
package algebra.domains;
import java.io.Serializable;
import algebra.Utilities;
import algebra.data.Polynomial;
import lang.Functions.Function2;
public class PolynomialDomain<T extends Serializable, S extends Ring<T>> extends CompositeDomain<Polynomial<T>, T, S> {
// Helper functions for creating domains by extension.
// Fully typed versions, for extensions in up to three indeterminates.
public static <T extends Serializable, S extends Ring<T>> PolynomialDomain<T, S> extend(S subDomain, String varName) {
return new PolynomialDomain<T, S>(varName, subDomain);
}
public static <T extends Serializable, S extends Ring<T>> PolynomialDomain<Polynomial<T>, PolynomialDomain<T, S>> extend(S subDomain, String v1, String v2) {
return extend(extend(subDomain, v2), v1);
}
public static <T extends Serializable, S extends Ring<T>> PolynomialDomain<Polynomial<Polynomial<T>>, PolynomialDomain<Polynomial<T>, PolynomialDomain<T, S>>> extend(S subDomain, String v1, String v2, String v3) {
return extend(extend(subDomain, v2, v3), v1);
}
// General extension, returning a polynomial domain (with unspecified sub-domains)
// when there is at least one indeterminate.
public static PolynomialDomain extend(Ring subDomain, String var1, String... vars) {
Ring result = subDomain;
for (int i = vars.length - 1; i >= 0; i--) {
result = extend(result, vars[i]);
}
return extend(result, var1);
}
public final String varName;
public PolynomialDomain(String varName, S subDomain) {
super((Class) Polynomial.class, subDomain);
this.varName = varName;
}
protected boolean objectEquals(Polynomial<T> a, Polynomial<T> b) {
return a.deg == b.deg && subDomain.equals(a.lc, b.lc) && equals(a.red, b.red);
}
// Note, the canonical form for zero is null; the restriction that leading
// coefficients are non-zero (because the rep is sparse) dictates that zero cannot
// consistently be represented by an instance of the Polynomial class.
public boolean isZero(Polynomial<T> a) {
return a == zero();
}
public Polynomial<T> polynomial(T lc, int deg, Polynomial<T> red) {
if (isZero(lc)) {
return red;
}
return new Polynomial<T>(lc, deg, red);
}
public Polynomial<T> polynomial(T lc, int deg) {
return polynomial(lc, deg, Polynomial.ZERO);
}
public Polynomial<T> polynomial(T lc) {
return polynomial(lc, 0);
}
public Polynomial<T> polynomial(T[] a) {
Polynomial<T> result = zero();
for (int i = 0; i < a.length; i++) {
result = polynomial(a[i], i, result);
}
return result;
}
public Polynomial<T> promote(T a) {
return polynomial(a);
}
public Polynomial<T> generator() {
return polynomial(subDomain.one(), 1);
}
public Polynomial<T> neg(Polynomial<T> a) {
if (a == zero()) {
return a;
}
return polynomial(neg(a.lc), a.deg, neg(a.red));
}
public Polynomial<T> sum(Polynomial<T> a, Polynomial<T> b) {
if (a == zero()) {
return b;
}
if (b == zero()) {
return a;
}
int n1 = a.deg;
int n2 = b.deg;
if (n1 < n2) {
return polynomial(b.lc, b.deg, sum(a, b.red));
}
if (n2 < n1) {
return polynomial(a.lc, a.deg, sum(a.red, b));
}
return polynomial(sum(a.lc, b.lc), a.deg, sum(a.red, b.red));
}
// Multiply by an element of the coeficient domain.
public Polynomial<T> prd(T u, Polynomial<T> v) {
if (v == zero()) {
return v;
}
return polynomial(prd(u, v.lc), v.deg, prd(u, v.red));
}
public Polynomial<T> prd(Polynomial<T> u, Polynomial<T> v) {
if (u == zero()) {
return u;
}
if (v == zero()) {
return v;
}
return polynomial(prd(u.lc, v.lc), u.deg + v.deg, sum(prd(u.red, v), prd(polynomial(u.lc, u.deg), v.red)));
}
// Divide by an element of the coeficient domain.
public Polynomial<T> div(Polynomial<T> u, T v) {
if (u == zero()) {
return u;
}
return polynomial(div(u.lc, v), u.deg, div(u.red, v));
}
public <R> R divide(Function2<Polynomial<T>, Polynomial<T>, R> f, Polynomial<T> u, Polynomial<T> v) {
if (v == zero()) {
throw new ArithmeticException("Divide by zero");
}
if (u == zero() || u.deg < v.deg) {
return f.of(zero(), u);
}
int delta = u.deg - v.deg;
T r = div(u.lc, v.lc);
QuoAndRem<Polynomial<T>> tmp = divide(QUO_AND_REM, sub(u.red, prd(polynomial(r, delta), v.red)), v);
return f.of(polynomial(r, delta, tmp.quo), tmp.rem);
}
public Polynomial<T> quo(Polynomial<T> a, Polynomial<T> b) {
return divide(QUO, a, b);
}
public Polynomial<T> rem(Polynomial<T> a, Polynomial<T> b) {
return divide(REM, a, b);
}
// This will normally throw when the coef domain is not a field (as in the case of
// multivatrate polynomials).
public Polynomial<T> div(Polynomial<T> a, Polynomial<T> b) throws ArithmeticException {
return divide(DIVIDE_OR_FAIL, a, b);
}
// u.deg >= v.deg
private T subresultant_rec(Polynomial<T> u, Polynomial<T> v, T g, T h) {
if (v == zero()) {
return h;
}
int delta = u.deg - v.deg;
T nh = div(pow(v.lc, delta), pow(h, delta - 1));
// for now, leave out out the optimisation below
// if (v.deg == 0) { return nh; }
Polynomial<T> pseudoRem = rem(prd(pow(v.lc, delta + 1), u), v);
return subresultant_rec(v, div(pseudoRem, prd(g, pow(h, delta))), v.lc, nh);
}
public T subresultant(Polynomial<T> u, Polynomial<T> v) {
T one = subDomain.one();
return degree(u) > degree(v) ? subresultant_rec(u, v, one, one) : subresultant_rec(v, u, one, one);
}
public Polynomial<T> parse(String s) {
return Utilities.polynomial(this, s);
}
public int degree(Polynomial<T> p) {
return p == zero() ? 0 : p.deg;
}
public T[] coeficients(Polynomial<T> p) {
T[] result = new VectorDomain<T>(subDomain, degree(p) + 1).zero(); // need to fill with zeros as rep is sparse.
for (Polynomial<T> r = p; r != zero(); r = r.red) {
result[r.deg] = r.lc;
}
return result;
}
// inefficient
public String toString(Polynomial<T> a) {
if (a == zero()) {
return "0";
}
String coefs = toString(a.lc);
String coef = coefs.equals("1") && a.deg != 0 ? "" : subDomain instanceof CompositeDomain ? "(" + coefs + ")" : coefs;
String pow = a.deg == 0 ? "" : a.deg == 1 ? varName : varName + "^" + a.deg;
String mul = coef != "" && pow != "" ? "*" : "";
String rest = a.red == zero() ? "" : " + " + toString(a.red);
return coef + mul + pow + rest;
}
}