Source code

001/*
002 * $Id: CycloUtil.java 5871 2018-07-20 15:58:45Z kredel $
003 */
004
005package edu.jas.ufd;
006
007
008import java.util.ArrayList;
009import java.util.List;
010import java.util.Map;
011import java.util.SortedMap;
012
013import org.apache.logging.log4j.Logger;
014import org.apache.logging.log4j.LogManager; 
015
016import edu.jas.arith.BigInteger;
017import edu.jas.arith.PrimeInteger;
018import edu.jas.poly.ExpVector;
019import edu.jas.poly.GenPolynomial;
020import edu.jas.poly.GenPolynomialRing;
021import edu.jas.poly.Monomial;
022import edu.jas.structure.Power;
023
024
025/**
026 * Cyclotomic polynomial utilities. Adapted from Sympy Python codes.
027 * @author Heinz Kredel
028 */
029
030public class CycloUtil {
031
032
033    private static final Logger logger = LogManager.getLogger(CycloUtil.class);
034
035
036    //private static final boolean debug = logger.isDebugEnabled();
037
038
039    /**
040     * Compute n-th cyclotomic polynomial.
041     * @param n number of polynomial.
042     * @param ring univariate polynomial ring of cyclotomic polynomial.
043     * @return n-th cyclotomic polynomial.
044     */
045    public static GenPolynomial<BigInteger> cyclotomicPolynomial(GenPolynomialRing<BigInteger> ring, long n) {
046        GenPolynomialRing<BigInteger> pfac = ring;
047        if (pfac == null) {
048            throw new IllegalArgumentException("ring must be non null");
049        }
050        GenPolynomial<BigInteger> h = pfac.univariate(0).subtract(pfac.getONE());
051        //System.out.println("h = " + h);
052        SortedMap<Long, Integer> fac = PrimeInteger.factors(n);
053        logger.info("factors = " + fac);
054
055        for (Map.Entry<Long, Integer> m : fac.entrySet()) {
056            // h = dup_quo(dup_inflate(h, p, K), h, K)
057            // h = dup_inflate(h, p**(k - 1), K)
058            long p = m.getKey();
059            int k = m.getValue();
060            GenPolynomial<BigInteger> ih = h.inflate(p);
061            //System.out.println("ih = " + ih);
062            h = ih.divide(h);
063            //System.out.println("h = " + h);
064            h = h.inflate(Power.power(p, k - 1));
065            //System.out.println("h = " + h + ", power = " + Power.power(p, k-1));
066        }
067        return h;
068    }
069
070
071    /**
072     * Compute the factors of the n-th cyclotomic polynomial.
073     * @param n number of polynomial.
074     * @param ring univariate polynomial ring of cyclotomic polynomial.
075     * @return list of factors of n-th cyclotomic polynomial.
076     */
077    public static List<GenPolynomial<BigInteger>> cyclotomicDecompose(GenPolynomialRing<BigInteger> ring,
078                    long n) {
079        GenPolynomialRing<BigInteger> pfac = ring;
080        if (pfac == null) {
081            throw new IllegalArgumentException("ring must be non null");
082        }
083        GenPolynomial<BigInteger> q = pfac.univariate(0).subtract(pfac.getONE());
084        //System.out.println("q = " + q);
085        List<GenPolynomial<BigInteger>> H = new ArrayList<GenPolynomial<BigInteger>>();
086        H.add(q);
087
088        SortedMap<Long, Integer> fac = PrimeInteger.factors(n);
089        logger.info("factors = " + fac);
090
091        for (Map.Entry<Long, Integer> m : fac.entrySet()) {
092            //Q = [ dup_quo(dup_inflate(h, p, K), h, K) for h in H ]
093            //H.extend(Q)
094            long p = m.getKey();
095            int k = m.getValue();
096            List<GenPolynomial<BigInteger>> Q = new ArrayList<GenPolynomial<BigInteger>>();
097            for (GenPolynomial<BigInteger> h : H) {
098                GenPolynomial<BigInteger> g = h.inflate(p).divide(h);
099                Q.add(g);
100            }
101            //System.out.println("Q = " + Q);
102            H.addAll(Q);
103            //for i in xrange(1, k):
104            //    Q = [ dup_inflate(q, p, K) for q in Q ]
105            //    H.extend(Q)
106            for (int i = 1; i < k; i++) {
107                List<GenPolynomial<BigInteger>> P = new ArrayList<GenPolynomial<BigInteger>>();
108                for (GenPolynomial<BigInteger> h : Q) {
109                    GenPolynomial<BigInteger> g = h.inflate(p);
110                    P.add(g);
111                }
112                //System.out.println("P = " + P);
113                Q = P;
114                H.addAll(P);
115            }
116        }
117        return H;
118    }
119
120
121    /**
122     * Compute the factors of the cyclotomic polynomial.
123     * @param p cyclotomic polynomial.
124     * @return list of factors of cyclotomic polynomial p or emtpy list.
125     */
126    public static List<GenPolynomial<BigInteger>> cyclotomicFactors(GenPolynomial<BigInteger> p) {
127        List<GenPolynomial<BigInteger>> H = new ArrayList<GenPolynomial<BigInteger>>();
128        long n = p.degree();
129        if (n <= 0) {
130            return H;
131        }
132        BigInteger lc, tc;
133        lc = p.leadingBaseCoefficient();
134        tc = p.trailingBaseCoefficient();
135        if (!lc.isONE() || (!tc.isONE() && !tc.negate().isONE())) {
136            return H;
137        }
138        if (p.length() != 2) { // other coefficients must be zero
139            return H;
140        }
141        // only case: x**n +/- 1
142        //F = _dup_cyclotomic_decompose(n, K)
143        //if not K.is_one(tc_f):
144        //   return F
145        GenPolynomialRing<BigInteger> pfac = p.ring;
146        List<GenPolynomial<BigInteger>> F = cyclotomicDecompose(pfac, n);
147        if (!tc.isONE()) {
148            return F;
149        }
150        //else:
151        //   H = []
152        //   for h in _dup_cyclotomic_decompose(2*n, K):
153        //       if h not in F:
154        //          H.append(h)
155        //return H                         
156        H = new ArrayList<GenPolynomial<BigInteger>>();
157        List<GenPolynomial<BigInteger>> F2 = cyclotomicDecompose(pfac, 2L * n);
158        for (GenPolynomial<BigInteger> h : F2) {
159            if (!F.contains(h)) {
160                H.add(h);
161            }
162        }
163        return H;
164    }
165
166
167    /**
168     * Test for cyclotomic polynomial.
169     * @param p polynomial.
170     * @return true if p is a cyclotomic polynomial, else false.
171     */
172    public static boolean isCyclotomicPolynomial(GenPolynomial<BigInteger> p) {
173        long n = p.degree();
174        if (n <= 0) {
175            return false;
176        }
177        BigInteger lc, tc;
178        lc = p.leadingBaseCoefficient();
179        tc = p.trailingBaseCoefficient();
180        if (!lc.isONE() || (!tc.isONE() && !tc.negate().isONE())) {
181            //System.out.println("!lc.isONE() || (!tc.isONE() && !tc.negate().isONE())");
182            return false;
183        }
184        if (p.length() == 2) { // other coefficients must be zero
185            return true;
186        }
187        // ignore: if not irreducible:
188        GenPolynomialRing<BigInteger> ring = p.ring;
189        if (ring.nvar != 1) {
190            throw new IllegalArgumentException("not univariate polynomial");
191        }
192        GenPolynomial<BigInteger> g, h, f;
193        g = ring.getZERO().copy();
194        h = ring.getZERO().copy();
195        long x = n % 2;
196        for (Monomial<BigInteger> m : p) {
197            ExpVector e = m.e;
198            long d = e.getVal(0);
199            ExpVector e2 = e.subst(0, d / 2L);
200            if (d % 2 == x) {
201                g.doPutToMap(e2, m.c);
202            } else {
203                h.doPutToMap(e2, m.c);
204            }
205        }
206        //g = dup_sqr(dup_strip(g), K)
207        //h = dup_sqr(dup_strip(h), K)
208        g = g.multiply(g);
209        h = h.multiply(h);
210        //System.out.println("g = " + g);
211        //System.out.println("h = " + h);
212        //F = dup_sub(g, dup_lshift(h, 1, K), K)
213        ExpVector on = ExpVector.create(1, 0, 1L);
214        f = g.subtract(h.multiply(on)).abs();
215        //System.out.println("f = " + f + ", f==p: " + f.equals(p));
216        //if F == f: return True 
217        if (f.equals(p)) {
218            return true;
219        }
220        //g = dup_mirror(f, K)
221        //if F == g and dup_cyclotomic_p(g, K):
222        //   return True
223        g = ring.getZERO().copy();
224        for (Monomial<BigInteger> m : p) {
225            ExpVector e = m.e;
226            long d = e.getVal(0);
227            if (d % 2 == 1) {
228                g.doPutToMap(e, m.c.negate());
229            } else {
230                g.doPutToMap(e, m.c);
231            }
232        }
233        g = g.abs();
234        //System.out.println("g = " + g + ", f==g: " + f.equals(g));
235        if (f.equals(g) && isCyclotomicPolynomial(g)) {
236            return true;
237        }
238        //G = dup_sqf_part(F, K)
239        Squarefree<BigInteger> engine;
240        engine = SquarefreeFactory.getImplementation(lc);
241        GenPolynomial<BigInteger> G;
242        G = engine.squarefreePart(f);
243        //System.out.println("G = " + G + ", G^2==f: " + G.multiply(G).equals(f));
244        //if dup_sqr(G, K) == F and dup_cyclotomic_p(G, K):
245        //   return True
246        if (G.multiply(G).equals(f) && isCyclotomicPolynomial(G)) {
247            return true;
248        }
249        return false;
250    }
251
252}