Anonymous user
Cyclotomic polynomial: Difference between revisions
no edit summary
(Created page with "{{task|Cyclotomic Polynomial}} The nth cyclotomic polynomial, for any positive integer n, is the unique irreducible polynomial of largest degree w...") |
No edit summary |
||
Line 9:
* The sequence [[oeis:A013594|A013594]] with the smallest order of cyclotomic polynomial containing n or -n as a coefficient.
<br><br>
=={{header|Java}}==
<lang java>
import java.util.ArrayList;
import java.util.Collections;
import java.util.Comparator;
import java.util.HashMap;
import java.util.List;
import java.util.Map;
import java.util.TreeMap;
public class CyclotomicPolynomial {
@SuppressWarnings("unused")
private static int divisions = 0;
private static int algorithm = 2;
public static void main(String[] args) throws Exception {
System.out.println("Task 1: cyclotomic polynomials for n <= 30:");
for ( int i = 1 ; i <= 30 ; i++ ) {
Polynomial p = cyclotomicPolynomial(i);
System.out.printf("CP[%d] = %s%n", i, p);
}
System.out.println("Task 2: Smallest cyclotomic polynomial with n or -n as a coefficient:");
int n = 0;
for ( int i = 1 ; i <= 10 ; i++ ) {
while ( true ) {
n++;
Polynomial cyclo = cyclotomicPolynomial(n);
if ( cyclo.hasCoefficientAbs(i) ) {
System.out.printf("CP[%d] has coefficient with magnitude = %d%n", n, i);
n--;
break;
}
}
}
}
private static final Map<Integer, Polynomial> COMPUTED = new HashMap<>();
private static Polynomial cyclotomicPolynomial(int n) {
if ( COMPUTED.containsKey(n) ) {
return COMPUTED.get(n);
}
//System.out.println("COMPUTE: n = " + n);
if ( n == 1 ) {
// Polynomial: x - 1
Polynomial p = new Polynomial(1, 1, -1, 0);
COMPUTED.put(1, p);
return p;
}
Map<Integer,Integer> factors = getFactors(n);
if ( factors.containsKey(n) ) {
// n prime
List<Term> termList = new ArrayList<>();
for ( int index = 0 ; index < n ; index++ ) {
termList.add(new Term(1, index));
}
Polynomial cyclo = new Polynomial(termList);
COMPUTED.put(n, cyclo);
return cyclo;
}
else if ( factors.size() == 2 && factors.containsKey(2) && factors.get(2) == 1 && factors.containsKey(n/2) && factors.get(n/2) == 1 ) {
// n = 2p
int prime = n/2;
List<Term> termList = new ArrayList<>();
int coeff = -1;
for ( int index = 0 ; index < prime ; index++ ) {
coeff *= -1;
termList.add(new Term(coeff, index));
}
Polynomial cyclo = new Polynomial(termList);
COMPUTED.put(n, cyclo);
return cyclo;
}
else if ( factors.size() == 1 && factors.containsKey(2) ) {
// n = 2^h
int h = factors.get(2);
List<Term> termList = new ArrayList<>();
termList.add(new Term(1, (int) Math.pow(2, h-1)));
termList.add(new Term(1, 0));
Polynomial cyclo = new Polynomial(termList);
COMPUTED.put(n, cyclo);
return cyclo;
}
else if ( factors.size() == 1 && ! factors.containsKey(n) ) {
// n = p^k
int p = 0;
for ( int prime : factors.keySet() ) {
p = prime;
}
int k = factors.get(p);
List<Term> termList = new ArrayList<>();
for ( int index = 0 ; index < p ; index++ ) {
termList.add(new Term(1, index * (int) Math.pow(p, k-1)));
}
Polynomial cyclo = new Polynomial(termList);
COMPUTED.put(n, cyclo);
return cyclo;
}
else if ( factors.size() == 2 && factors.containsKey(2) ) {
// n = 2^h * p^k
int p = 0;
for ( int prime : factors.keySet() ) {
if ( prime != 2 ) {
p = prime;
}
}
List<Term> termList = new ArrayList<>();
int coeff = -1;
int twoExp = (int) Math.pow(2, factors.get(2)-1);
int k = factors.get(p);
for ( int index = 0 ; index < p ; index++ ) {
coeff *= -1;
termList.add(new Term(coeff, index * twoExp * (int) Math.pow(p, k-1)));
}
Polynomial cyclo = new Polynomial(termList);
COMPUTED.put(n, cyclo);
return cyclo;
}
else if ( factors.containsKey(2) && ((n/2) % 2 == 1) && (n/2) > 1 ) {
// CP(2m)[x] = CP(-m)[x], n odd integer > 1
Polynomial cycloDiv2 = cyclotomicPolynomial(n/2);
List<Term> termList = new ArrayList<>();
for ( Term term : cycloDiv2.polynomialTerms ) {
termList.add(term.exponent % 2 == 0 ? term : term.negate());
}
Polynomial cyclo = new Polynomial(termList);
COMPUTED.put(n, cyclo);
return cyclo;
}
// General Case
if ( algorithm == 0 ) {
// Slow - uses basic definition.
List<Integer> divisors = getDivisors(n);
// Polynomial: ( x^n - 1 )
Polynomial cyclo = new Polynomial(1, n, -1, 0);
for ( int i : divisors ) {
Polynomial p = cyclotomicPolynomial(i);
cyclo = cyclo.divide(p);
}
COMPUTED.put(n, cyclo);
return cyclo;
}
else if ( algorithm == 1 ) {
// Faster. Remove Max divisor (and all divisors of max divisor) - only one divide for all divisors of Max Divisor
List<Integer> divisors = getDivisors(n);
int maxDivisor = Integer.MIN_VALUE;
for ( int div : divisors ) {
maxDivisor = Math.max(maxDivisor, div);
}
List<Integer> divisorsExceptMax = new ArrayList<Integer>();
for ( int div : divisors ) {
if ( maxDivisor % div != 0 ) {
divisorsExceptMax.add(div);
}
}
// Polynomial: ( x^n - 1 ) / ( x^m - 1 ), where m is the max divisor
Polynomial cyclo = new Polynomial(1, n, -1, 0).divide(new Polynomial(1, maxDivisor, -1, 0));
for ( int i : divisorsExceptMax ) {
Polynomial p = cyclotomicPolynomial(i);
cyclo = cyclo.divide(p);
}
COMPUTED.put(n, cyclo);
return cyclo;
}
else if ( algorithm == 2 ) {
// Fastest
// Let p ; q be primes such that p does not divide n, and q q divides n.
// Then CP(np)[x] = CP(n)[x^p] / CP(n)[x]
int m = 1;
Polynomial cyclo = cyclotomicPolynomial(m);
List<Integer> primes = new ArrayList<>(factors.keySet());
Collections.sort(primes);
for ( int prime : primes ) {
// CP(m)[x]
Polynomial cycloM = cyclo;
// Compute CP(m)[x^p].
List<Term> termList = new ArrayList<>();
for ( Term t : cycloM.polynomialTerms ) {
termList.add(new Term(t.coefficient, t.exponent * prime));
}
cyclo = new Polynomial(termList).divide(cycloM);
m = m * prime;
}
// Now, m is the largest square free divisor of n
int s = n / m;
// Compute CP(n)[x] = CP(m)[x^s]
List<Term> termList = new ArrayList<>();
for ( Term t : cyclo.polynomialTerms ) {
termList.add(new Term(t.coefficient, t.exponent * s));
}
cyclo = new Polynomial(termList);
COMPUTED.put(n, cyclo);
return cyclo;
}
else {
throw new RuntimeException("ERROR 103: Invalid algorithm.");
}
}
private static final List<Integer> getDivisors(int number) {
List<Integer> divisors = new ArrayList<Integer>();
long sqrt = (long) Math.sqrt(number);
for ( int i = 1 ; i <= sqrt ; i++ ) {
if ( number % i == 0 ) {
divisors.add(i);
int div = number / i;
if ( div != i && div != number ) {
divisors.add(div);
}
}
}
return divisors;
}
private static final Map<Integer,Map<Integer,Integer>> allFactors = new TreeMap<Integer,Map<Integer,Integer>>();
static {
Map<Integer,Integer> factors = new TreeMap<Integer,Integer>();
factors.put(2, 1);
allFactors.put(2, factors);
}
public static Integer MAX_ALL_FACTORS = 100000;
public static final Map<Integer,Integer> getFactors(Integer number) {
if ( allFactors.containsKey(number) ) {
return allFactors.get(number);
}
Map<Integer,Integer> factors = new TreeMap<Integer,Integer>();
if ( number % 2 == 0 ) {
Map<Integer,Integer> factorsdDivTwo = getFactors(number/2);
factors.putAll(factorsdDivTwo);
factors.merge(2, 1, (v1, v2) -> v1 + v2);
if ( number < MAX_ALL_FACTORS )
allFactors.put(number, factors);
return factors;
}
boolean prime = true;
long sqrt = (long) Math.sqrt(number);
for ( int i = 3 ; i <= sqrt ; i += 2 ) {
if ( number % i == 0 ) {
prime = false;
factors.putAll(getFactors(number/i));
factors.merge(i, 1, (v1, v2) -> v1 + v2);
if ( number < MAX_ALL_FACTORS )
allFactors.put(number, factors);
return factors;
}
}
if ( prime ) {
factors.put(number, 1);
if ( number < MAX_ALL_FACTORS )
allFactors.put(number, factors);
}
return factors;
}
private static final class Polynomial {
private List<Term> polynomialTerms;
// Format - coeff, exp, coeff, exp, (repeating in pairs) . . .
public Polynomial(int ... values) {
if ( values.length % 2 != 0 ) {
throw new IllegalArgumentException("ERROR 102: Polynomial constructor. Length must be even. Length = " + values.length);
}
polynomialTerms = new ArrayList<>();
for ( int i = 0 ; i < values.length ; i += 2 ) {
Term t = new Term(values[i], values[i+1]);
polynomialTerms.add(t);
}
Collections.sort(polynomialTerms, new TermSorter());
}
public Polynomial() {
// zero
polynomialTerms = new ArrayList<>();
polynomialTerms.add(new Term(0,0));
}
private boolean hasCoefficientAbs(int coeff) {
for ( Term term : polynomialTerms ) {
if ( Math.abs(term.coefficient) == coeff ) {
return true;
}
}
return false;
}
private Polynomial(List<Term> termList) {
if ( termList.size() == 0 ) {
// zero
termList.add(new Term(0,0));
}
else {
// Remove zero terms if needed
for ( int i = 0 ; i < termList.size() ; i++ ) {
if ( termList.get(i).coefficient == 0 ) {
termList.remove(i);
}
}
}
if ( termList.size() == 0 ) {
// zero
termList.add(new Term(0,0));
}
polynomialTerms = termList;
Collections.sort(polynomialTerms, new TermSorter());
}
public Polynomial divide(Polynomial v) {
//long start = System.currentTimeMillis();
divisions++;
Polynomial q = new Polynomial();
Polynomial r = this;
long lcv = v.leadingCoefficient();
long dv = v.degree();
while ( r.degree() >= v.degree() ) {
long lcr = r.leadingCoefficient();
long s = lcr / lcv; // Integer division
Term term = new Term(s, r.degree() - dv);
q = q.add(term);
r = r.add(v.multiply(term.negate()));
}
//long end = System.currentTimeMillis();
//System.out.printf("Divide: Elapsed = %d, Degree 1 = %d, Degree 2 = %d%n", (end-start), this.degree(), v.degree());
return q;
}
public Polynomial add(Polynomial polynomial) {
List<Term> termList = new ArrayList<>();
int thisCount = polynomialTerms.size();
int polyCount = polynomial.polynomialTerms.size();
while ( thisCount > 0 || polyCount > 0 ) {
Term thisTerm = thisCount == 0 ? null : polynomialTerms.get(thisCount-1);
Term polyTerm = polyCount == 0 ? null : polynomial.polynomialTerms.get(polyCount-1);
if ( thisTerm == null ) {
termList.add(polyTerm.clone());
polyCount--;
}
else if (polyTerm == null ) {
termList.add(thisTerm.clone());
thisCount--;
}
else if ( thisTerm.degree() == polyTerm.degree() ) {
Term t = thisTerm.add(polyTerm);
if ( t.coefficient != 0 ) {
termList.add(t);
}
thisCount--;
polyCount--;
}
else if ( thisTerm.degree() < polyTerm.degree() ) {
termList.add(thisTerm.clone());
thisCount--;
}
else {
termList.add(polyTerm.clone());
polyCount--;
}
}
return new Polynomial(termList);
}
public Polynomial add(Term term) {
List<Term> termList = new ArrayList<>();
boolean added = false;
for ( int index = 0 ; index < polynomialTerms.size() ; index++ ) {
Term currentTerm = polynomialTerms.get(index);
if ( currentTerm.exponent == term.exponent ) {
added = true;
if ( currentTerm.coefficient + term.coefficient != 0 ) {
termList.add(currentTerm.add(term));
}
}
else {
termList.add(currentTerm.clone());
}
}
if ( ! added ) {
termList.add(term.clone());
}
return new Polynomial(termList);
}
public Polynomial multiply(Term term) {
List<Term> termList = new ArrayList<>();
for ( int index = 0 ; index < polynomialTerms.size() ; index++ ) {
Term currentTerm = polynomialTerms.get(index);
termList.add(currentTerm.clone().multiply(term));
}
return new Polynomial(termList);
}
public Polynomial clone() {
List<Term> clone = new ArrayList<>();
for ( Term t : polynomialTerms ) {
clone.add(new Term(t.coefficient, t.exponent));
}
return new Polynomial(clone);
}
public long leadingCoefficient() {
// long coefficient = 0;
// long degree = Integer.MIN_VALUE;
// for ( Term t : polynomialTerms ) {
// if ( t.degree() > degree ) {
// coefficient = t.coefficient;
// degree = t.degree();
// }
// }
return polynomialTerms.get(0).coefficient;
}
public long degree() {
// long degree = Integer.MIN_VALUE;
// for ( Term t : polynomialTerms ) {
// if ( t.degree() > degree ) {
// degree = t.degree();
// }
// }
return polynomialTerms.get(0).exponent;
}
@Override
public String toString() {
StringBuilder sb = new StringBuilder();
boolean first = true;
for ( Term term : polynomialTerms ) {
if ( first ) {
sb.append(term);
first = false;
}
else {
sb.append(" ");
if ( term.coefficient > 0 ) {
sb.append("+ ");
sb.append(term);
}
else {
sb.append("- ");
sb.append(term.negate());
}
}
}
return sb.toString();
}
}
private static final class TermSorter implements Comparator<Term> {
@Override
public int compare(Term o1, Term o2) {
return (int) (o2.exponent - o1.exponent);
}
}
// Note: Cyclotomic Polynomials have small coefficients. Not appropriate for general polynomial usage.
private static final class Term {
long coefficient;
long exponent;
public Term(long c, long e) {
coefficient = c;
exponent = e;
}
public Term clone() {
return new Term(coefficient, exponent);
}
public Term multiply(Term term) {
return new Term(coefficient * term.coefficient, exponent + term.exponent);
}
public Term add(Term term) {
if ( exponent != term.exponent ) {
throw new RuntimeException("ERROR 102: Exponents not equal.");
}
return new Term(coefficient + term.coefficient, exponent);
}
public Term negate() {
return new Term(-coefficient, exponent);
}
public long degree() {
return exponent;
}
@Override
public String toString() {
if ( coefficient == 0 ) {
return "0";
}
if ( exponent == 0 ) {
return "" + coefficient;
}
if ( coefficient == 1 ) {
if ( exponent == 1 ) {
return "x";
}
else {
return "x^" + exponent;
}
}
if ( exponent == 1 ) {
return coefficient + "x";
}
return coefficient + "x^" + exponent;
}
}
}
</lang>
{{out}}
<pre>
Task 1: cyclotomic polynomials for n <= 30:
CP[1] = x - 1
CP[2] = x + 1
CP[3] = x^2 + x + 1
CP[4] = x^2 + 1
CP[5] = x^4 + x^3 + x^2 + x + 1
CP[6] = x^2 - x + 1
CP[7] = x^6 + x^5 + x^4 + x^3 + x^2 + x + 1
CP[8] = x^4 + 1
CP[9] = x^6 + x^3 + 1
CP[10] = x^4 - x^3 + x^2 - x + 1
CP[11] = x^10 + x^9 + x^8 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + x + 1
CP[12] = x^4 - x^2 + 1
CP[13] = x^12 + x^11 + x^10 + x^9 + x^8 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + x + 1
CP[14] = x^6 - x^5 + x^4 - x^3 + x^2 - x + 1
CP[15] = x^8 - x^7 + x^5 - x^4 + x^3 - x + 1
CP[16] = x^8 + 1
CP[17] = x^16 + x^15 + x^14 + x^13 + x^12 + x^11 + x^10 + x^9 + x^8 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + x + 1
CP[18] = x^6 - x^3 + 1
CP[19] = x^18 + x^17 + x^16 + x^15 + x^14 + x^13 + x^12 + x^11 + x^10 + x^9 + x^8 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + x + 1
CP[20] = x^8 - x^6 + x^4 - x^2 + 1
CP[21] = x^12 - x^11 + x^9 - x^8 + x^6 - x^4 + x^3 - x + 1
CP[22] = x^10 - x^9 + x^8 - x^7 + x^6 - x^5 + x^4 - x^3 + x^2 - x + 1
CP[23] = x^22 + x^21 + x^20 + x^19 + x^18 + x^17 + x^16 + x^15 + x^14 + x^13 + x^12 + x^11 + x^10 + x^9 + x^8 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + x + 1
CP[24] = x^8 - x^4 + 1
CP[25] = x^20 + x^15 + x^10 + x^5 + 1
CP[26] = x^12 - x^11 + x^10 - x^9 + x^8 - x^7 + x^6 - x^5 + x^4 - x^3 + x^2 - x + 1
CP[27] = x^18 + x^9 + 1
CP[28] = x^12 - x^10 + x^8 - x^6 + x^4 - x^2 + 1
CP[29] = x^28 + x^27 + x^26 + x^25 + x^24 + x^23 + x^22 + x^21 + x^20 + x^19 + x^18 + x^17 + x^16 + x^15 + x^14 + x^13 + x^12 + x^11 + x^10 + x^9 + x^8 + x^7 + x^6 + x^5 + x^4 + x^3 + x^2 + x + 1
CP[30] = x^8 + x^7 - x^5 - x^4 - x^3 + x + 1
Task 2: Smallest cyclotomic polynomial with n or -n as a coefficient:
CP[1] has coefficient with magnitude = 1
CP[105] has coefficient with magnitude = 2
CP[385] has coefficient with magnitude = 3
CP[1365] has coefficient with magnitude = 4
CP[1785] has coefficient with magnitude = 5
CP[2805] has coefficient with magnitude = 6
CP[3135] has coefficient with magnitude = 7
CP[6545] has coefficient with magnitude = 8
CP[6545] has coefficient with magnitude = 9
CP[10465] has coefficient with magnitude = 10
</pre>
| |||