On Mon, Mar 31, 2008 at 12:43 PM, John Cremona <[EMAIL PROTECTED]>
wrote:
>
> There was some discussion a week or so ago about a more efficient
> implentation of cyclotomic polynomials, which led to trac#2654. I
> have tried various alternatives of the prototype code posted there,
> finding that the speed varied enormously as a result of
> trivial-seeming changes.
>
> The key operation needed is to replace f(x) by f(x^m) for various
> positive integers m, where f is an integer polynomial. Doing
> {{{
> f = f(x**m)
> }}}
> was very slow, especially for large m. Faster is to apply this function:
> {{{
> def powsub(f,m):
> assert m>0
> if m==1: return f
> g={}
> d=f.dict()
> for i in d:
> g[m*i]=d[i]
> return f.parent()(g)
> }}}
>
> Is there a better way? And is this operation needed elsewhere, in
> which case the function should be available generally?
>
> I thought of using sparse polynomials, but the algorithm also needs
> exact polynomial division whcih (as far as I can see) is not available
> for sparse integer polys.
>
> Suggestions welcome, so we can put a decent patch in for this
> important function!
>
Hi John,
Sure a special function to compute f(x^n) quickly seems like a good idea.
By the way, I just did some cyclotomic polynomial benchmarking
in on Intel Core 2 Duo 2.6Ghz laptop. Here are the results:
n | Sage 2.11 | newcyc at #2654 | PARI | Magma
10^5 | 1.29 | 0.37 | 1.24 | 0.02
10^6 | ? | 32.89 | 123.8 | 0.20
10^7 | ? | ? | ? | 2.37
Magma crushes everything else yet again...
----
Relevant code:
{{{id=15|
///
}}}
{{{id=8|
def newcyc(n):
assert n>0
plist = n.prime_divisors()
X = PolynomialRing(ZZ,'X').gen()
f = X-1
m = n
for p in plist:
f = f(X**p)//f(X)
m //= p
return f.subs(X**m)
///
}}}
{{{id=9|
time f= newcyc(10^5)
///
CPU time: 0.37 s, Wall time: 0.37 s
}}}
{{{id=21|
time f= newcyc(10^6)
///
CPU time: 32.89 s, Wall time: 33.06 s
}}}
{{{id=22|
time f= newcyc(10^7)
///
}}}
{{{id=24|
time f = cyclotomic_polynomial(10^5)
///
CPU time: 1.29 s, Wall time: 1.29 s
}}}
{{{id=25|
time f = cyclotomic_polynomial(10^6)
///
CPU time: 129.18 s, Wall time: 129.54 s
}}}
{{{id=26|
time f = cyclotomic_polynomial(10^7)
///
}}}
{{{id=27|
///
}}}
{{{id=12|
def h(n):
print magma.eval('time f := CyclotomicPolynomial(%s);'%n)
///
}}}
{{{id=16|
h(10^5)
///
Time: 0.020
}}}
{{{id=13|
h(10^6)
///
Time: 0.200
}}}
{{{id=17|
h(10^7)
///
Time: 2.370
}}}
{{{id=7|
def g(n):
print gp.eval('gettime; f = polcyclo(%s); gettime/1000.0'%n)
///
}}}
{{{id=11|
g(10^5)
///
1.243000000000000000000000000
}}}
{{{id=18|
g(10^6)
///
123.7870000000000000000000000
}}}
{{{id=19|
g(10^7)
///
}}}
{{{id=20|
///
}}}
--~--~---------~--~----~------------~-------~--~----~
To post to this group, send email to sage-devel@googlegroups.com
To unsubscribe from this group, send email to [EMAIL PROTECTED]
For more options, visit this group at http://groups.google.com/group/sage-devel
URLs: http://www.sagemath.org
-~----------~----~----~----~------~----~------~--~---
