How to resolve the algorithm Cyclotomic polynomial step by step in the Go programming language
Published on 12 May 2024 09:40 PM
How to resolve the algorithm Cyclotomic polynomial step by step in the Go programming language
Table of Contents
Problem Statement
The nth Cyclotomic polynomial, for any positive integer n, is the unique irreducible polynomial of largest degree with integer coefficients that is a divisor of x^n − 1, and is not a divisor of x^k − 1 for any k < n.
Let's start with the solution:
Step by Step solution about How to resolve the algorithm Cyclotomic polynomial step by step in the Go programming language
This code computes cyclotomic polynomials for a given integer n. Cyclotomic polynomials are polynomials that have the nth roots of unity as their roots.
The code uses a recursive algorithm to compute the cyclotomic polynomial for a given n.
The algorithm starts by checking if the value has been previously computed, in the event that it has, it returns the computed value.
If the value hasn't been computed, the code checks if n is equal to 1, in which case the cyclotomic polynomial is x - 1.
If n is not equal to 1, the code calls the getFactors function to get the prime factors of n.
Depending on the value of the prime factors, the code uses different algorithms to compute the cyclotomic polynomial.
The code uses three different algorithms, each one with a different level of speed and accuracy.
The first algorithm is the slowest but most accurate, the second algorithm is faster but less accurate, and the third algorithm is the fastest but least accurate.
The code uses the fastest algorithm by default, but the user can specify which algorithm to use by setting the algo constant to one of the following values:
0: Use the slowest but most accurate algorithm
1: Use the faster but less accurate algorithm
2: Use the fastest but least accurate algorithm.
Once the cyclotomic polynomial has been computed, the code returns the polynomial and stores it in a map so that it can be reused later.
The main function of the code calls the cycloPoly function to compute the cyclotomic polynomials for n <= 30 and prints the results.
The main function also calls the cycloPoly function to compute the smallest cyclotomic polynomial with n or -n as a coefficient and prints the results.
Source code in the go programming language
package main
import (
"fmt"
"log"
"math"
"sort"
"strings"
)
const (
algo = 2
maxAllFactors = 100000
)
func iabs(i int) int {
if i < 0 {
return -i
}
return i
}
type term struct{ coef, exp int }
func (t term) mul(t2 term) term {
return term{t.coef * t2.coef, t.exp + t2.exp}
}
func (t term) add(t2 term) term {
if t.exp != t2.exp {
log.Fatal("exponents unequal in term.add method")
}
return term{t.coef + t2.coef, t.exp}
}
func (t term) negate() term { return term{-t.coef, t.exp} }
func (t term) String() string {
switch {
case t.coef == 0:
return "0"
case t.exp == 0:
return fmt.Sprintf("%d", t.coef)
case t.coef == 1:
if t.exp == 1 {
return "x"
} else {
return fmt.Sprintf("x^%d", t.exp)
}
case t.exp == 1:
return fmt.Sprintf("%dx", t.coef)
}
return fmt.Sprintf("%dx^%d", t.coef, t.exp)
}
type poly struct{ terms []term }
// pass coef, exp in pairs as parameters
func newPoly(values ...int) poly {
le := len(values)
if le == 0 {
return poly{[]term{term{0, 0}}}
}
if le%2 != 0 {
log.Fatalf("odd number of parameters (%d) passed to newPoly function", le)
}
var terms []term
for i := 0; i < le; i += 2 {
terms = append(terms, term{values[i], values[i+1]})
}
p := poly{terms}.tidy()
return p
}
func (p poly) hasCoefAbs(coef int) bool {
for _, t := range p.terms {
if iabs(t.coef) == coef {
return true
}
}
return false
}
func (p poly) add(p2 poly) poly {
p3 := newPoly()
le, le2 := len(p.terms), len(p2.terms)
for le > 0 || le2 > 0 {
if le == 0 {
p3.terms = append(p3.terms, p2.terms[le2-1])
le2--
} else if le2 == 0 {
p3.terms = append(p3.terms, p.terms[le-1])
le--
} else {
t := p.terms[le-1]
t2 := p2.terms[le2-1]
if t.exp == t2.exp {
t3 := t.add(t2)
if t3.coef != 0 {
p3.terms = append(p3.terms, t3)
}
le--
le2--
} else if t.exp < t2.exp {
p3.terms = append(p3.terms, t)
le--
} else {
p3.terms = append(p3.terms, t2)
le2--
}
}
}
return p3.tidy()
}
func (p poly) addTerm(t term) poly {
q := newPoly()
added := false
for i := 0; i < len(p.terms); i++ {
ct := p.terms[i]
if ct.exp == t.exp {
added = true
if ct.coef+t.coef != 0 {
q.terms = append(q.terms, ct.add(t))
}
} else {
q.terms = append(q.terms, ct)
}
}
if !added {
q.terms = append(q.terms, t)
}
return q.tidy()
}
func (p poly) mulTerm(t term) poly {
q := newPoly()
for i := 0; i < len(p.terms); i++ {
ct := p.terms[i]
q.terms = append(q.terms, ct.mul(t))
}
return q.tidy()
}
func (p poly) div(v poly) poly {
q := newPoly()
lcv := v.leadingCoef()
dv := v.degree()
for p.degree() >= v.degree() {
lcp := p.leadingCoef()
s := lcp / lcv
t := term{s, p.degree() - dv}
q = q.addTerm(t)
p = p.add(v.mulTerm(t.negate()))
}
return q.tidy()
}
func (p poly) leadingCoef() int {
return p.terms[0].coef
}
func (p poly) degree() int {
return p.terms[0].exp
}
func (p poly) String() string {
var sb strings.Builder
first := true
for _, t := range p.terms {
if first {
sb.WriteString(t.String())
first = false
} else {
sb.WriteString(" ")
if t.coef > 0 {
sb.WriteString("+ ")
sb.WriteString(t.String())
} else {
sb.WriteString("- ")
sb.WriteString(t.negate().String())
}
}
}
return sb.String()
}
// in place descending sort by term.exp
func (p poly) sortTerms() {
sort.Slice(p.terms, func(i, j int) bool {
return p.terms[i].exp > p.terms[j].exp
})
}
// sort terms and remove any unnecesary zero terms
func (p poly) tidy() poly {
p.sortTerms()
if p.degree() == 0 {
return p
}
for i := len(p.terms) - 1; i >= 0; i-- {
if p.terms[i].coef == 0 {
copy(p.terms[i:], p.terms[i+1:])
p.terms[len(p.terms)-1] = term{0, 0}
p.terms = p.terms[:len(p.terms)-1]
}
}
if len(p.terms) == 0 {
p.terms = append(p.terms, term{0, 0})
}
return p
}
func getDivisors(n int) []int {
var divs []int
sqrt := int(math.Sqrt(float64(n)))
for i := 1; i <= sqrt; i++ {
if n%i == 0 {
divs = append(divs, i)
d := n / i
if d != i && d != n {
divs = append(divs, d)
}
}
}
return divs
}
var (
computed = make(map[int]poly)
allFactors = make(map[int]map[int]int)
)
func init() {
f := map[int]int{2: 1}
allFactors[2] = f
}
func getFactors(n int) map[int]int {
if f, ok := allFactors[n]; ok {
return f
}
factors := make(map[int]int)
if n%2 == 0 {
factorsDivTwo := getFactors(n / 2)
for k, v := range factorsDivTwo {
factors[k] = v
}
factors[2]++
if n < maxAllFactors {
allFactors[n] = factors
}
return factors
}
prime := true
sqrt := int(math.Sqrt(float64(n)))
for i := 3; i <= sqrt; i += 2 {
if n%i == 0 {
prime = false
for k, v := range getFactors(n / i) {
factors[k] = v
}
factors[i]++
if n < maxAllFactors {
allFactors[n] = factors
}
return factors
}
}
if prime {
factors[n] = 1
if n < maxAllFactors {
allFactors[n] = factors
}
}
return factors
}
func cycloPoly(n int) poly {
if p, ok := computed[n]; ok {
return p
}
if n == 1 {
// polynomial: x - 1
p := newPoly(1, 1, -1, 0)
computed[1] = p
return p
}
factors := getFactors(n)
cyclo := newPoly()
if _, ok := factors[n]; ok {
// n is prime
for i := 0; i < n; i++ {
cyclo.terms = append(cyclo.terms, term{1, i})
}
} else if len(factors) == 2 && factors[2] == 1 && factors[n/2] == 1 {
// n == 2p
prime := n / 2
coef := -1
for i := 0; i < prime; i++ {
coef *= -1
cyclo.terms = append(cyclo.terms, term{coef, i})
}
} else if len(factors) == 1 {
if h, ok := factors[2]; ok {
// n == 2^h
cyclo.terms = append(cyclo.terms, term{1, 1 << (h - 1)}, term{1, 0})
} else if _, ok := factors[n]; !ok {
// n == p ^ k
p := 0
for prime := range factors {
p = prime
}
k := factors[p]
for i := 0; i < p; i++ {
pk := int(math.Pow(float64(p), float64(k-1)))
cyclo.terms = append(cyclo.terms, term{1, i * pk})
}
}
} else if len(factors) == 2 && factors[2] != 0 {
// n = 2^h * p^k
p := 0
for prime := range factors {
if prime != 2 {
p = prime
}
}
coef := -1
twoExp := 1 << (factors[2] - 1)
k := factors[p]
for i := 0; i < p; i++ {
coef *= -1
pk := int(math.Pow(float64(p), float64(k-1)))
cyclo.terms = append(cyclo.terms, term{coef, i * twoExp * pk})
}
} else if factors[2] != 0 && ((n/2)%2 == 1) && (n/2) > 1 {
// CP(2m)[x] == CP(-m)[x], n odd integer > 1
cycloDiv2 := cycloPoly(n / 2)
for _, t := range cycloDiv2.terms {
t2 := t
if t.exp%2 != 0 {
t2 = t.negate()
}
cyclo.terms = append(cyclo.terms, t2)
}
} else if algo == 0 {
// slow - uses basic definition
divs := getDivisors(n)
// polynomial: x^n - 1
cyclo = newPoly(1, n, -1, 0)
for _, i := range divs {
p := cycloPoly(i)
cyclo = cyclo.div(p)
}
} else if algo == 1 {
// faster - remove max divisor (and all divisors of max divisor)
// only one divide for all divisors of max divisor
divs := getDivisors(n)
maxDiv := math.MinInt32
for _, d := range divs {
if d > maxDiv {
maxDiv = d
}
}
var divsExceptMax []int
for _, d := range divs {
if maxDiv%d != 0 {
divsExceptMax = append(divsExceptMax, d)
}
}
// polynomial: ( x^n - 1 ) / ( x^m - 1 ), where m is the max divisor
cyclo = newPoly(1, n, -1, 0)
cyclo = cyclo.div(newPoly(1, maxDiv, -1, 0))
for _, i := range divsExceptMax {
p := cycloPoly(i)
cyclo = cyclo.div(p)
}
} else if algo == 2 {
// fastest
// let p, q be primes such that p does not divide n, and q divides n
// then CP(np)[x] = CP(n)[x^p] / CP(n)[x]
m := 1
cyclo = cycloPoly(m)
var primes []int
for prime := range factors {
primes = append(primes, prime)
}
sort.Ints(primes)
for _, prime := range primes {
// CP(m)[x]
cycloM := cyclo
// compute CP(m)[x^p]
var terms []term
for _, t := range cycloM.terms {
terms = append(terms, term{t.coef, t.exp * prime})
}
cyclo = newPoly()
cyclo.terms = append(cyclo.terms, terms...)
cyclo = cyclo.tidy()
cyclo = cyclo.div(cycloM)
m *= prime
}
// now, m is the largest square free divisor of n
s := n / m
// Compute CP(n)[x] = CP(m)[x^s]
var terms []term
for _, t := range cyclo.terms {
terms = append(terms, term{t.coef, t.exp * s})
}
cyclo = newPoly()
cyclo.terms = append(cyclo.terms, terms...)
} else {
log.Fatal("invalid algorithm")
}
cyclo = cyclo.tidy()
computed[n] = cyclo
return cyclo
}
func main() {
fmt.Println("Task 1: cyclotomic polynomials for n <= 30:")
for i := 1; i <= 30; i++ {
p := cycloPoly(i)
fmt.Printf("CP[%2d] = %s\n", i, p)
}
fmt.Println("\nTask 2: Smallest cyclotomic polynomial with n or -n as a coefficient:")
n := 0
for i := 1; i <= 10; i++ {
for {
n++
cyclo := cycloPoly(n)
if cyclo.hasCoefAbs(i) {
fmt.Printf("CP[%d] has coefficient with magnitude = %d\n", n, i)
n--
break
}
}
}
}
You may also check:How to resolve the algorithm Sum multiples of 3 and 5 step by step in the BASIC programming language
You may also check:How to resolve the algorithm Reduced row echelon form step by step in the Phix programming language
You may also check:How to resolve the algorithm SOAP step by step in the PureBasic programming language
You may also check:How to resolve the algorithm Parsing/Shunting-yard algorithm step by step in the EchoLisp programming language
You may also check:How to resolve the algorithm Haversine formula step by step in the FOCAL programming language