Source code

001/*
002 * $Id: PolyGBUtil.java 4049 2012-07-25 17:10:49Z kredel $
003 */
004
005package edu.jas.gbufd;
006
007
008import java.util.ArrayList;
009import java.util.List;
010
011import org.apache.log4j.Logger;
012
013import edu.jas.gb.GroebnerBaseAbstract;
014import edu.jas.poly.ExpVector;
015import edu.jas.poly.GenPolynomial;
016import edu.jas.poly.GenPolynomialRing;
017import edu.jas.poly.PolyUtil;
018import edu.jas.structure.GcdRingElem;
019import edu.jas.structure.RingElem;
020
021
022/**
023 * Package gbufd utilities.
024 * @author Heinz Kredel
025 */
026
027public class PolyGBUtil {
028
029
030    private static final Logger logger = Logger.getLogger(PolyGBUtil.class);
031
032
033    private static boolean debug = logger.isDebugEnabled();
034
035
036    /**
037     * Test for resultant.
038     * @param A generic polynomial.
039     * @param B generic polynomial.
040     * @param r generic polynomial.
041     * @return true if res(A,B) isContained in ideal(A,B), else false.
042     */
043    public static <C extends GcdRingElem<C>> 
044      boolean isResultant(GenPolynomial<C> A, GenPolynomial<C> B, GenPolynomial<C> r) {
045        if (r == null || r.isZERO()) {
046            return true;
047        }
048        GroebnerBaseAbstract<C> bb = GBFactory.<C> getImplementation(r.ring.coFac);
049        List<GenPolynomial<C>> F = new ArrayList<GenPolynomial<C>>(2);
050        F.add(A);
051        F.add(B);
052        List<GenPolynomial<C>> G = bb.GB(F);
053        //System.out.println("G = " + G);
054        GenPolynomial<C> n = bb.red.normalform(G, r);
055        //System.out.println("n = " + n);
056        return n.isZERO();
057    }
058
059
060    /**
061     * Top pseudo reduction wrt the main variables.
062     * @param P generic polynomial.
063     * @param A list of generic polynomials sorted according to appearing main variables.
064     * @return top pseudo remainder of P wrt. A for the appearing variables.
065     */
066    public static <C extends RingElem<C>> 
067      GenPolynomial<C> topPseudoRemainder(List<GenPolynomial<C>> A, GenPolynomial<C> P) {
068        if (A == null || A.isEmpty()) {
069            return P.monic();
070        }
071        if (P.isZERO()) {
072            return P;
073        }
074        //System.out.println("remainder, P = " + P);
075        GenPolynomialRing<C> pfac = A.get(0).ring;
076        if (pfac.nvar <= 1) { // recursion base 
077            GenPolynomial<C> R = PolyUtil.<C> baseSparsePseudoRemainder(P, A.get(0));
078            return R.monic();
079        }
080        // select polynomials according to the main variable
081        GenPolynomialRing<GenPolynomial<C>> rfac = pfac.recursive(1);
082        GenPolynomial<C> Q = A.get(0); // wrong, must eventually search polynomial
083        GenPolynomial<GenPolynomial<C>> qr = PolyUtil.<C> recursive(rfac, Q);
084        GenPolynomial<GenPolynomial<C>> pr = PolyUtil.<C> recursive(rfac, P);
085        GenPolynomial<GenPolynomial<C>> rr;
086        if (qr.isONE()) {
087            return P.ring.getZERO();
088        }
089        if (qr.degree(0) > 0) {
090            rr = PolyUtil.<C> recursiveSparsePseudoRemainder(pr, qr);
091            //System.out.println("remainder, pr = " + pr);
092            //System.out.println("remainder, qr = " + qr);
093            //System.out.println("remainder, rr = " + rr);
094        } else {
095            rr = pr;
096        }
097        if (rr.degree(0) > 0) {
098            GenPolynomial<C> R = PolyUtil.<C> distribute(pfac, rr);
099            return R.monic();
100            // not further reduced wrt. other variables = top-reduction only
101        }
102        List<GenPolynomial<C>> zeroDeg = zeroDegrees(A);
103        GenPolynomial<C> R = topPseudoRemainder(zeroDeg, rr.leadingBaseCoefficient());
104        R = R.extend(pfac, 0, 0L);
105        return R.monic();
106    }
107
108
109    /**
110     * Top coefficient pseudo remainder of the leading coefficient of P wrt A in the main variables.
111     * @param P generic polynomial in n+1 variables.
112     * @param A list of generic polynomials in n variables sorted according to appearing main variables.
113     * @return pseudo remainder of the leading coefficient of P wrt A.
114     */
115    public static <C extends RingElem<C>> 
116      GenPolynomial<C> topCoefficientPseudoRemainder(List<GenPolynomial<C>> A, GenPolynomial<C> P) {
117        if (A == null || A.isEmpty()) {
118            return P.monic();
119        }
120        if (P.isZERO()) {
121            return P;
122        }
123        GenPolynomialRing<C> pfac = P.ring;
124        GenPolynomialRing<C> pfac1 = A.get(0).ring;
125        if (pfac1.nvar <= 1) { // recursion base 
126            GenPolynomial<C> a = A.get(0);
127            GenPolynomialRing<GenPolynomial<C>> rfac = pfac.recursive(pfac.nvar - 1);
128            GenPolynomial<GenPolynomial<C>> pr = PolyUtil.<C> recursive(rfac, P);
129            // ldcf(P,x_m) = q a + r 
130            GenPolynomial<GenPolynomial<C>> rr = PolyGBUtil.<C> coefficientPseudoRemainderBase(pr, a);
131            GenPolynomial<C> R = PolyUtil.<C> distribute(pfac, rr);
132            return R.monic();
133        }
134        // select polynomials according to the main variable
135        GenPolynomialRing<GenPolynomial<C>> rfac1 = pfac1.recursive(1);
136        int nv = pfac.nvar - pfac1.nvar;
137        GenPolynomialRing<GenPolynomial<C>> rfac = pfac.recursive(1 + nv);
138        GenPolynomialRing<GenPolynomial<GenPolynomial<C>>> rfac2 = rfac.recursive(nv);
139        if (debug) {
140            logger.info("rfac =" + rfac);
141        }
142        GenPolynomial<GenPolynomial<C>> pr = PolyUtil.<C> recursive(rfac, P);
143        GenPolynomial<GenPolynomial<GenPolynomial<C>>> pr2 = PolyUtil.<GenPolynomial<C>> recursive(rfac2, pr);
144        //System.out.println("recursion, pr2 = " + pr2);
145        GenPolynomial<C> Q = A.get(0);
146        GenPolynomial<GenPolynomial<C>> qr = PolyUtil.<C> recursive(rfac1, Q);
147        GenPolynomial<GenPolynomial<GenPolynomial<C>>> rr;
148        if (qr.isONE()) {
149            return P.ring.getZERO();
150        }
151        if (qr.degree(0) > 0) {
152            // pseudo remainder:  ldcf(P,x_m) = a q + r 
153            rr = PolyGBUtil.<C> coefficientPseudoRemainder(pr2, qr);
154            //System.out.println("recursion, qr  = " + qr);
155            //System.out.println("recursion, pr  = " + pr2);
156            //System.out.println("recursion, rr  = " + rr);
157        } else {
158            rr = pr2;
159        }
160        // reduction wrt. the other variables
161        List<GenPolynomial<C>> zeroDeg = zeroDegrees(A);
162        GenPolynomial<GenPolynomial<C>> Rr = PolyUtil.<GenPolynomial<C>> distribute(rfac, rr);
163        GenPolynomial<C> R = PolyUtil.<C> distribute(pfac, Rr);
164        R = topCoefficientPseudoRemainder(zeroDeg, R);
165        return R.monic();
166    }
167
168
169    /**
170     * Polynomial leading coefficient pseudo remainder.
171     * @param P generic polynomial in n+1 variables.
172     * @param A generic polynomial in n variables.
173     * @return pseudo remainder of the leading coefficient of P wrt A, with
174     *         ldcf(A)<sup>m'</sup> P = quotient * A + remainder.
175     */
176    public static <C extends RingElem<C>> 
177      GenPolynomial<GenPolynomial<GenPolynomial<C>>> coefficientPseudoRemainder(
178                    GenPolynomial<GenPolynomial<GenPolynomial<C>>> P, GenPolynomial<GenPolynomial<C>> A) {
179        if (A == null || A.isZERO()) { // findbugs
180            throw new ArithmeticException(P + " division by zero " + A);
181        }
182        if (A.isONE()) {
183            return P.ring.getZERO();
184        }
185        if (P.isZERO() || P.isONE()) {
186            return P;
187        }
188        GenPolynomialRing<GenPolynomial<GenPolynomial<C>>> pfac = P.ring;
189        GenPolynomialRing<GenPolynomial<C>> afac = A.ring; // == pfac.coFac
190        GenPolynomial<GenPolynomial<GenPolynomial<C>>> r = P;
191        GenPolynomial<GenPolynomial<C>> h;
192        GenPolynomial<GenPolynomial<GenPolynomial<C>>> hr;
193        GenPolynomial<GenPolynomial<C>> ldcf = P.leadingBaseCoefficient();
194        long m = ldcf.degree(0);
195        long n = A.degree(0);
196        GenPolynomial<C> c = A.leadingBaseCoefficient();
197        GenPolynomial<GenPolynomial<C>> cc = afac.getZERO().sum(c);
198        //System.out.println("cc = " + cc);
199        ExpVector e = A.leadingExpVector();
200        for (long i = m; i >= n; i--) {
201            if (r.isZERO()) {
202                return r;
203            }
204            GenPolynomial<GenPolynomial<C>> p = r.leadingBaseCoefficient();
205            ExpVector g = r.leadingExpVector();
206            long k = p.degree(0);
207            if (i == k) {
208                GenPolynomial<C> pl = p.leadingBaseCoefficient();
209                ExpVector f = p.leadingExpVector();
210                f = f.subtract(e);
211                r = r.multiply(cc); // coeff cc
212                h = A.multiply(pl, f); // coeff ac
213                hr = new GenPolynomial<GenPolynomial<GenPolynomial<C>>>(pfac, h, g);
214                r = r.subtract(hr);
215            } else {
216                r = r.multiply(cc);
217            }
218            //System.out.println("r = " + r);
219        }
220        if (r.degree(0) < P.degree(0)) { // recursion for degree
221            r = coefficientPseudoRemainder(r, A);
222        }
223        return r;
224    }
225
226
227    /**
228     * Polynomial leading coefficient pseudo remainder, base case.
229     * @param P generic polynomial in 1+1 variables.
230     * @param A generic polynomial in 1 variable.
231     * @return pseudo remainder of the leading coefficient of P wrt. A, with
232     *         ldcf(A)<sup>m'</sup> P = quotient * A + remainder.
233     */
234    public static <C extends RingElem<C>> 
235      GenPolynomial<GenPolynomial<C>> coefficientPseudoRemainderBase(
236                    GenPolynomial<GenPolynomial<C>> P, GenPolynomial<C> A) {
237        if (A == null || A.isZERO()) { // findbugs
238            throw new ArithmeticException(P + " division by zero " + A);
239        }
240        if (A.isONE()) {
241            return P.ring.getZERO();
242        }
243        if (P.isZERO() || P.isONE()) {
244            return P;
245        }
246        GenPolynomialRing<GenPolynomial<C>> pfac = P.ring;
247        GenPolynomialRing<C> afac = A.ring; // == pfac.coFac
248        GenPolynomial<GenPolynomial<C>> r = P;
249        GenPolynomial<C> h;
250        GenPolynomial<GenPolynomial<C>> hr;
251        GenPolynomial<C> ldcf = P.leadingBaseCoefficient();
252        long m = ldcf.degree(0);
253        long n = A.degree(0);
254        C c = A.leadingBaseCoefficient();
255        GenPolynomial<C> cc = afac.getZERO().sum(c);
256        //System.out.println("cc = " + cc);
257        ExpVector e = A.leadingExpVector();
258        for (long i = m; i >= n; i--) {
259            if (r.isZERO()) {
260                return r;
261            }
262            GenPolynomial<C> p = r.leadingBaseCoefficient();
263            ExpVector g = r.leadingExpVector();
264            long k = p.degree(0);
265            if (i == k) {
266                C pl = p.leadingBaseCoefficient();
267                ExpVector f = p.leadingExpVector();
268                f = f.subtract(e);
269                r = r.multiply(cc); // coeff cc
270                h = A.multiply(pl, f); // coeff ac
271                hr = new GenPolynomial<GenPolynomial<C>>(pfac, h, g);
272                r = r.subtract(hr);
273            } else {
274                r = r.multiply(cc);
275            }
276            //System.out.println("r = " + r);
277        }
278        if (r.degree(0) < P.degree(0)) { // recursion for degree
279            r = coefficientPseudoRemainderBase(r, A);
280        }
281        return r;
282    }
283
284
285    /**
286     * Extract polynomials with degree zero in the main variable.
287     * @param A list of generic polynomials in n variables.
288     * @return Z = [a_i] with deg(a_i,x_n) = 0 and in n-1 variables.
289     */
290    public static <C extends RingElem<C>> List<GenPolynomial<C>> zeroDegrees(List<GenPolynomial<C>> A) {
291        if (A == null || A.isEmpty()) {
292            return A;
293        }
294        GenPolynomialRing<C> pfac = A.get(0).ring;
295        GenPolynomialRing<GenPolynomial<C>> rfac = pfac.recursive(1);
296        List<GenPolynomial<C>> zeroDeg = new ArrayList<GenPolynomial<C>>(A.size());
297        for (int i = 0; i < A.size(); i++) {
298            GenPolynomial<C> q = A.get(i);
299            GenPolynomial<GenPolynomial<C>> fr = PolyUtil.<C> recursive(rfac, q);
300            if (fr.degree(0) == 0) {
301                zeroDeg.add(fr.leadingBaseCoefficient());
302            }
303        }
304        return zeroDeg;
305    }
306
307}