Source code

001/*
002 * $Id: GreatestCommonDivisorPrimitive.java 5872 2018-07-20 16:01:46Z kredel $
003 */
004
005package edu.jas.fd;
006
007
008import org.apache.logging.log4j.Logger;
009import org.apache.logging.log4j.LogManager; 
010
011import edu.jas.poly.GenPolynomial;
012import edu.jas.poly.GenSolvablePolynomial;
013import edu.jas.poly.PolyUtil;
014import edu.jas.structure.GcdRingElem;
015import edu.jas.structure.RingFactory;
016
017
018/**
019 * (Non-unique) factorization domain greatest common divisor common algorithms
020 * with primitive polynomial remainder sequence.
021 * @param <C> coefficient type
022 * @author Heinz Kredel
023 */
024
025public class GreatestCommonDivisorPrimitive<C extends GcdRingElem<C>> extends
026                GreatestCommonDivisorAbstract<C> {
027
028
029    private static final Logger logger = LogManager.getLogger(GreatestCommonDivisorPrimitive.class);
030
031
032    private static final boolean debug = logger.isDebugEnabled();
033
034
035    /**
036     * Constructor.
037     * @param cf coefficient ring.
038     */
039    public GreatestCommonDivisorPrimitive(RingFactory<C> cf) {
040        super(cf);
041    }
042
043
044    /**
045     * Univariate GenSolvablePolynomial greatest common divisor. Uses
046     * pseudoRemainder for remainder.
047     * @param P univariate GenSolvablePolynomial.
048     * @param S univariate GenSolvablePolynomial.
049     * @return gcd(P,S) with P = P'*gcd(P,S) and S = S'*gcd(P,S).
050     */
051    @Override
052    public GenSolvablePolynomial<C> leftBaseGcd(GenSolvablePolynomial<C> P, GenSolvablePolynomial<C> S) {
053        if (S == null || S.isZERO()) {
054            return P;
055        }
056        if (P == null || P.isZERO()) {
057            return S;
058        }
059        if (P.ring.nvar > 1) {
060            throw new IllegalArgumentException(this.getClass().getName() + " no univariate polynomial");
061        }
062        boolean field = P.ring.coFac.isField();
063        long e = P.degree(0);
064        long f = S.degree(0);
065        GenSolvablePolynomial<C> q;
066        GenSolvablePolynomial<C> r;
067        if (f > e) {
068            r = P;
069            q = S;
070            long g = f;
071            f = e;
072            e = g;
073        } else {
074            q = P;
075            r = S;
076        }
077        if (debug) {
078            logger.debug("degrees: e = " + e + ", f = " + f);
079        }
080        C c;
081        if (field) {
082            r = r.monic();
083            q = q.monic();
084            c = P.ring.getONECoefficient();
085        } else {
086            r = (GenSolvablePolynomial<C>) r.abs();
087            q = (GenSolvablePolynomial<C>) q.abs();
088            C a = rightBaseContent(r);
089            C b = rightBaseContent(q);
090            r = divide(r, a); // indirection
091            q = divide(q, b); // indirection
092            c = gcd(a, b); // indirection
093        }
094        //System.out.println("baseCont: gcd(cont) = " + b);
095        if (r.isONE()) {
096            return r.multiply(c);
097        }
098        if (q.isONE()) {
099            return q.multiply(c);
100        }
101        GenSolvablePolynomial<C> x;
102        //System.out.println("baseGCD: q = " + q);
103        //System.out.println("baseGCD: r = " + r);
104        while (!r.isZERO()) {
105            x = FDUtil.<C> leftBaseSparsePseudoRemainder(q, r);
106            q = r;
107            if (field) {
108                r = x.monic();
109            } else {
110                r = leftBasePrimitivePart(x);
111            }
112            //System.out.println("baseGCD: q = " + q);
113            //System.out.println("baseGCD: r = " + r);
114        }
115        return (GenSolvablePolynomial<C>) (q.multiply(c)).abs();
116    }
117
118
119    /**
120     * Univariate GenSolvablePolynomial right greatest common divisor. Uses
121     * pseudoRemainder for remainder.
122     * @param P univariate GenSolvablePolynomial.
123     * @param S univariate GenSolvablePolynomial.
124     * @return gcd(P,S) with P = gcd(P,S)*P' and S = gcd(P,S)*S'.
125     */
126    @Override
127    public GenSolvablePolynomial<C> rightBaseGcd(GenSolvablePolynomial<C> P, GenSolvablePolynomial<C> S) {
128        if (S == null || S.isZERO()) {
129            return P;
130        }
131        if (P == null || P.isZERO()) {
132            return S;
133        }
134        if (P.ring.nvar > 1) {
135            throw new IllegalArgumentException(this.getClass().getName() + " no univariate polynomial");
136        }
137        boolean field = P.ring.coFac.isField();
138        long e = P.degree(0);
139        long f = S.degree(0);
140        GenSolvablePolynomial<C> q;
141        GenSolvablePolynomial<C> r;
142        if (f > e) {
143            r = P;
144            q = S;
145            long g = f;
146            f = e;
147            e = g;
148        } else {
149            q = P;
150            r = S;
151        }
152        if (debug) {
153            logger.debug("degrees: e = " + e + ", f = " + f);
154        }
155        C c;
156        if (field) {
157            r = r.monic();
158            q = q.monic();
159            c = P.ring.getONECoefficient();
160        } else {
161            r = (GenSolvablePolynomial<C>) r.abs();
162            q = (GenSolvablePolynomial<C>) q.abs();
163            C a = leftBaseContent(r);
164            C b = leftBaseContent(q);
165            r = divide(r, a); // indirection
166            q = divide(q, b); // indirection
167            c = gcd(a, b); // indirection
168        }
169        //System.out.println("baseCont: gcd(cont) = " + b);
170        if (r.isONE()) {
171            return r.multiply(c);
172        }
173        if (q.isONE()) {
174            return q.multiply(c);
175        }
176        GenSolvablePolynomial<C> x;
177        //System.out.println("baseGCD: q = " + q);
178        //System.out.println("baseGCD: r = " + r);
179        while (!r.isZERO()) {
180            x = FDUtil.<C> rightBaseSparsePseudoRemainder(q, r);
181            q = r;
182            if (field) {
183                r = x.monic();
184            } else {
185                r = rightBasePrimitivePart(x);
186            }
187            //System.out.println("baseGCD: q = " + q);
188            //System.out.println("baseGCD: r = " + r);
189        }
190        return (GenSolvablePolynomial<C>) (q.multiplyLeft(c)).abs();
191    }
192
193
194    /**
195     * Univariate GenSolvablePolynomial left recursive greatest common divisor.
196     * Uses pseudoRemainder for remainder.
197     * @param P univariate recursive GenSolvablePolynomial.
198     * @param S univariate recursive GenSolvablePolynomial.
199     * @return gcd(P,S) with P = P'*gcd(P,S)*p and S = S'*gcd(P,S)*s, where
200     *         deg_main(p) = deg_main(s) == 0.
201     */
202    @Override
203    public GenSolvablePolynomial<GenPolynomial<C>> leftRecursiveUnivariateGcd(
204                    GenSolvablePolynomial<GenPolynomial<C>> P, GenSolvablePolynomial<GenPolynomial<C>> S) {
205        if (S == null || S.isZERO()) {
206            return P;
207        }
208        if (P == null || P.isZERO()) {
209            return S;
210        }
211        if (P.ring.nvar > 1) {
212            throw new IllegalArgumentException("no univariate polynomial");
213        }
214        boolean field = P.leadingBaseCoefficient().ring.coFac.isField();
215        long e = P.degree(0);
216        long f = S.degree(0);
217        GenSolvablePolynomial<GenPolynomial<C>> q, r, x, qs, rs;
218        if (f > e) {
219            r = P;
220            q = S;
221            long g = f;
222            f = e;
223            e = g;
224        } else if (f < e) {
225            q = P;
226            r = S;
227        } else { // f == e
228            if (P.leadingBaseCoefficient().degree() > S.leadingBaseCoefficient().degree()) {
229                q = P;
230                r = S;
231            } else {
232                r = P;
233                q = S;
234            }
235        }
236        if (debug) {
237            logger.debug("degrees: e = " + e + ", f = " + f);
238        }
239        if (field) {
240            r = PolyUtil.<C> monic(r);
241            q = PolyUtil.<C> monic(q);
242        } else {
243            r = (GenSolvablePolynomial<GenPolynomial<C>>) r.abs();
244            q = (GenSolvablePolynomial<GenPolynomial<C>>) q.abs();
245        }
246        GenSolvablePolynomial<C> a = rightRecursiveContent(r);
247        rs = FDUtil.<C> recursiveLeftDivide(r, a);
248        if (debug) {
249            logger.info("recCont a = " + a + ", r = " + r);
250            logger.info("recCont r/a = " + rs + ", r%a = " + r.subtract(rs.multiply(a)));
251            if (!r.equals(rs.multiply(a))) {
252                System.out.println("recGcd, r         = " + r);
253                System.out.println("recGcd, cont(r)   = " + a);
254                System.out.println("recGcd, pp(r)     = " + rs);
255                System.out.println("recGcd, pp(r)c(r) = " + rs.multiply(a));
256                System.out.println("recGcd, c(r)pp(r) = " + rs.multiplyLeft(a));
257                throw new RuntimeException("recGcd, pp: not divisible");
258            }
259        }
260        r = rs;
261        GenSolvablePolynomial<C> b = rightRecursiveContent(q);
262        qs = FDUtil.<C> recursiveLeftDivide(q, b);
263        if (debug) {
264            logger.info("recCont b = " + b + ", q = " + q);
265            logger.info("recCont q/b = " + qs + ", q%b = " + q.subtract(qs.multiply(b)));
266            if (!q.equals(qs.multiply(b))) {
267                System.out.println("recGcd, q         = " + q);
268                System.out.println("recGcd, cont(q)   = " + b);
269                System.out.println("recGcd, pp(q)     = " + qs);
270                System.out.println("recGcd, pp(q)c(q) = " + qs.multiply(b));
271                System.out.println("recGcd, c(q)pp(q) = " + qs.multiplyLeft(b));
272                throw new RuntimeException("recGcd, pp: not divisible");
273            }
274        }
275        q = qs;
276        //no: GenSolvablePolynomial<C> c = leftGcd(a, b); // go to recursion
277        GenSolvablePolynomial<C> c = rightGcd(a, b); // go to recursion
278        logger.info("Gcd(contents) c = " + c + ", a = " + a + ", b = " + b);
279        if (r.isONE()) {
280            return r.multiply(c);
281        }
282        if (q.isONE()) {
283            return q.multiply(c);
284        }
285        if (debug) {
286            logger.info("r.ring = " + r.ring.toScript());
287        }
288        while (!r.isZERO()) {
289            x = FDUtil.<C> recursiveSparsePseudoRemainder(q, r);
290            q = r;
291            r = rightRecursivePrimitivePart(x);
292            if (field) {
293                r = PolyUtil.<C> monic(r);
294            }
295        }
296        if (debug) {
297            logger.info("gcd(pp) = " + q + ", ring = " + P.ring.toScript());
298        }
299        q = (GenSolvablePolynomial<GenPolynomial<C>>) q.multiply(c).abs();
300        return q;
301    }
302
303
304    /**
305     * Univariate GenSolvablePolynomial right recursive greatest common divisor.
306     * Uses pseudoRemainder for remainder.
307     * @param P univariate recursive GenSolvablePolynomial.
308     * @param S univariate recursive GenSolvablePolynomial.
309     * @return gcd(P,S) with P = p*gcd(P,S)*P' and S = s*gcd(P,S)*S', where
310     *         deg_main(p) = deg_main(s) == 0.
311     */
312    @Override
313    public GenSolvablePolynomial<GenPolynomial<C>> rightRecursiveUnivariateGcd(
314                    GenSolvablePolynomial<GenPolynomial<C>> P, GenSolvablePolynomial<GenPolynomial<C>> S) {
315        if (S == null || S.isZERO()) {
316            return P;
317        }
318        if (P == null || P.isZERO()) {
319            return S;
320        }
321        if (P.ring.nvar > 1) {
322            throw new IllegalArgumentException("no univariate polynomial");
323        }
324        boolean field = P.leadingBaseCoefficient().ring.coFac.isField();
325        long e = P.degree(0);
326        long f = S.degree(0);
327        GenSolvablePolynomial<GenPolynomial<C>> q, r, x, qs, rs;
328        if (f > e) {
329            r = P;
330            q = S;
331            long g = f;
332            f = e;
333            e = g;
334        } else if (f < e) {
335            q = P;
336            r = S;
337        } else { // f == e
338            if (P.leadingBaseCoefficient().degree() > S.leadingBaseCoefficient().degree()) {
339                q = P;
340                r = S;
341            } else {
342                r = P;
343                q = S;
344            }
345        }
346        if (debug) {
347            logger.debug("RI-degrees: e = " + e + ", f = " + f);
348        }
349        if (field) {
350            r = PolyUtil.<C> monic(r);
351            q = PolyUtil.<C> monic(q);
352        } else {
353            r = (GenSolvablePolynomial<GenPolynomial<C>>) r.abs();
354            q = (GenSolvablePolynomial<GenPolynomial<C>>) q.abs();
355        }
356        GenSolvablePolynomial<C> a = leftRecursiveContent(r);
357        rs = FDUtil.<C> recursiveRightDivide(r, a);
358        if (debug) {
359            logger.info("RI-recCont a = " + a + ", r = " + r);
360            logger.info("RI-recCont r/a = " + r + ", r%a = " + r.subtract(rs.multiplyLeft(a)));
361            if (!r.equals(rs.multiplyLeft(a))) {
362                System.out.println("RI-recGcd, r         = " + r);
363                System.out.println("RI-recGcd, cont(r)   = " + a);
364                System.out.println("RI-recGcd, pp(r)     = " + rs);
365                System.out.println("RI-recGcd, c(r)pp(r) = " + rs.multiplyLeft(a));
366                throw new RuntimeException("RI-recGcd, pp: not divisible");
367            }
368        }
369        r = rs;
370        GenSolvablePolynomial<C> b = leftRecursiveContent(q);
371        qs = FDUtil.<C> recursiveRightDivide(q, b);
372        if (debug) {
373            logger.info("RI-recCont b = " + b + ", q = " + q);
374            logger.info("RI-recCont q/b = " + qs + ", q%b = " + q.subtract(qs.multiplyLeft(b)));
375            if (!q.equals(qs.multiplyLeft(b))) {
376                System.out.println("RI-recGcd, q         = " + q);
377                System.out.println("RI-recGcd, cont(q)   = " + b);
378                System.out.println("RI-recGcd, pp(q)     = " + qs);
379                System.out.println("RI-recGcd, c(q)pp(q) = " + qs.multiplyLeft(b));
380                throw new RuntimeException("RI-recGcd, pp: not divisible");
381            }
382        }
383        q = qs;
384        //no: GenSolvablePolynomial<C> c = rightGcd(a, b); // go to recursion
385        GenSolvablePolynomial<C> c = leftGcd(a, b); // go to recursion
386        logger.info("RI-Gcd(contents) c = " + c + ", a = " + a + ", b = " + b);
387        if (r.isONE()) {
388            return r.multiplyLeft(c);
389        }
390        if (q.isONE()) {
391            return q.multiplyLeft(c);
392        }
393        if (debug) {
394            logger.info("RI-r.ring = " + r.ring.toScript());
395        }
396        while (!r.isZERO()) {
397            x = FDUtil.<C> recursiveRightSparsePseudoRemainder(q, r);
398            q = r;
399            r = leftRecursivePrimitivePart(x);
400            if (field) {
401                r = PolyUtil.<C> monic(r);
402            }
403        }
404        if (debug) {
405            logger.info("RI-gcd(pp) = " + q + ", ring = " + P.ring.toScript());
406        }
407        q = (GenSolvablePolynomial<GenPolynomial<C>>) q.multiplyLeft(c).abs();
408        return q;
409    }
410
411}