mersenneforum.org  

Go Back   mersenneforum.org > Great Internet Mersenne Prime Search > Math

Reply
 
Thread Tools
Old 2012-01-04, 06:53   #1
rapso
 

3×232 Posts
Default some FFT Mul questions

as I like programming, and I like math, the mersenne prime search is really interesting for me.
I've implemented my little math lib and I try to improve it, now I'm at the point where I want to optimze my toom-cook mul by the schoenhage–strassen algorithm (I mean, I'll use it for numbers beyond some value), and as a first step I try to grasp the integer multiply (I think the schoenhage–strassen algorithm is just a specialized case of the FFT mul(?)).

so, I've got most of it, but I've some open questions I can't find any answer for (feel free to save your time by just posting related links).

I've implemented my first version based on the summary of this svg:
http://en.wikipedia.org/wiki/File:In...ion_by_FFT.svg

first question is rather to satisfy my theoretical understanding: PowerMod[Length[x], -1, b]; is actually a call to find the multiplicative inverse in the ring of b(337), using the euclid algorithm is the proper way I guess? and if the result is <0, I guess it's correct to bring it within the bounds of the ring, but modulo of something negative is (in this case -42 ->295) will result in -42, is there a proper name for this operation?

2. how to choose b? I found that not all primes >324 work, actually, most don't work, at least the result seem to be wrong. I think my implementation is correct, is it expected that most primes don't work? is there a smarter way than try&error to find working primes?

3. r is 85, as it's the 8th root of unit in 337, but how do I calculate this? atm, I've just a lame loop trying all possible numbers to find the 8th root, is there an algorithm? I can find some info about nth-root calculation and about the root of unity calculation, but the n-th root? (not just for power of 2, or is the NTT transform just working on PO2 sizes, like FFT?)

( 4. the NTT transform in the summary is O(n*n), I shall be able to convert it to a fast O(n*log(n)) version, just like FT->FFT? any concerns?)

5. not really a question, but, if you have any links with good informations on this topic, I'd appreciate if you'd share them with me. I do that all because of my joy of math and I want to fully understand this topic.
  Reply With Quote
Old 2012-01-04, 14:46   #2
fivemack
(loop (#_fork))
 
fivemack's Avatar
 
Feb 2006
Cambridge, England

13·491 Posts
Default

I don't know a name for the operation 'add p until not negative then reduce mod p', but it's the one you need.

To do an NTT on X terms, there has to be an Xth root of unity mod p, which means that (by Lagrange's theorem) p must be ==1 mod X.

To find the Xth root of unity, start with a primitive root mod p, and take the (p-1)/X th power of it ... to check whether some q is a primitive root mod p, check that q^R is not =1 mod p whenever R is a factor of (p-1).
fivemack is offline   Reply With Quote
Old 2012-01-04, 15:40   #3
davieddy
 
davieddy's Avatar
 
"Lucan"
Dec 2006
England

2×3×13×83 Posts
Default

Quote:
Originally Posted by fivemack View Post
I don't know a name for the operation 'add p until not negative then reduce mod p', but it's the one you need.
Neither do I, but as many others must do, I find the process very handy in numerous situations.
Anyone care to suggest a name for it?

How about "Division"?

David

Last fiddled with by davieddy on 2012-01-04 at 15:44
davieddy is offline   Reply With Quote
Old 2012-01-04, 16:20   #4
science_man_88
 
science_man_88's Avatar
 
"Forget I exist"
Jul 2009
Dumbassville

100000110000002 Posts
Default

Quote:
Originally Posted by fivemack View Post
I don't know a name for the operation 'add p until not negative then reduce mod p', but it's the one you need.
+p*x is what I get from this wording also would p*x reduce to 0*x=0 mod p and hence add nothing to it ?
science_man_88 is offline   Reply With Quote
Old 2012-01-05, 01:50   #5
jasonp
Tribal Bullet
 
jasonp's Avatar
 
Oct 2004

24×13×17 Posts
Default

I've usually seen it referred to as 'normalization', or 'modular reduction' sometimes. Allowing negative numbers when they aren't strictly necessary, i.e. -p < x < p, is a redundant representation, and removing the redundancy so that -p/2 <= x < p/2 is called balanced representation by Crandall.
jasonp is offline   Reply With Quote
Old 2012-01-26, 14:50   #6
Maximus
 
Nov 2011

148 Posts
Default

Quote:
Originally Posted by rapso View Post
3. r is 85, as it's the 8th root of unit in 337, but how do I calculate this? atm, I've just a lame loop trying all possible numbers to find the 8th root, is there an algorithm? I can find some info about nth-root calculation and about the root of unity calculation, but the n-th root?
One way to calculate root of unity is exhaustive search.
Example:
8/2=4.
337=24*21.
Take a number, let it be 2.
Calculate 24*21 mod 337 = 1. But must be -1. So 2 is not eligible.
Let take 3.
34*21 mod 337 = 336 = -1.
3 is eligible.
Calculate r=321 mod 337 = 252.
So r4 mod 337 = 336 = -1. (And r8 mod 337 = 1). So r is 8th root of unity in 337.
Other roots of unity are:
r3 mod 337 = 2523 mod 337 = 226;
r5 mod 337 = 2525 mod 337 = 85;
r7 mod 337 = 2527 mod 337 = 111.

Quote:
Originally Posted by rapso View Post
(not just for power of 2, or is the NTT transform just working on PO2 sizes, like FFT?)
Fast NTT like FFT can work not only PO2 size, but on any size. So if size=p1*p2*p3..., then computing complexity of individual iterations is O(p1*p1), O(p2*p2), O(p3*p3), ... So PO2 is preferable.

Quote:
Originally Posted by rapso View Post
( 4. the NTT transform in the summary is O(n*n), I shall be able to convert it to a fast O(n*log(n)) version, just like FT->FFT? any concerns?)
Fast NTT is similar to FFT, but uses modular arithmetics instead of complex number arithmetics.
Maximus is offline   Reply With Quote
Old 2012-01-26, 16:36   #7
jasonp
Tribal Bullet
 
jasonp's Avatar
 
Oct 2004

1101110100002 Posts
Default

The typical method for calculating a root of unity is to start with numbers modulo p, a prime of the form k*2^n+1 where n is large. Pick a number g < p and then for each unique prime factor q of k*2^n verify that g^((p-1)/q) != 1 mod p. Any g for which you can show this is then a primitive root of order k*2^n, i.e. a root of unity that allows a Fourier-like transform of size k*2^n.

For a transform size smaller than that, say size N, just make g^((p - 1)/N) mod p the root of unity.

For power-of-two sizes a number-theoretic FFT looks a lot like an FFT over the complex numbers, but it doesn't use complex arithmetic. For non-power-of-two sizes I think that it's more complicated, you cannot just take a floating-point FFT and replace the arithmetic operations.
jasonp is offline   Reply With Quote
Old 2012-01-26, 18:59   #8
Maximus
 
Nov 2011

22·3 Posts
Default some corrections

Quote:
Originally Posted by Maximus View Post
Calculate 24*21 mod 337 = 1. But must be -1. So 2 is not eligible.
Saying more correct, number must not be equal to 1, but needs some additional calculations and longer to explain.

Quote:
Originally Posted by Maximus View Post
Fast NTT like FFT can work not only PO2 size, but on any size. So if size=p1*p2*p3..., then computing complexity of individual iterations is O(p1*p1), O(p2*p2), O(p3*p3), ... So PO2 is preferable.
I mean one step; one iteration takes O(SIZE*pi) of time.
Maximus is offline   Reply With Quote
Reply

Thread Tools


Similar Threads
Thread Thread Starter Forum Replies Last Post
Two questions: Dubslow GPU Computing 1 2011-08-05 18:22
Questions about the QS Carmichael Factoring 8 2007-04-10 11:30
Questions OmbooHankvald Prime Sierpinski Project 2 2005-08-01 20:18
LLR questions OmbooHankvald Math 6 2005-06-23 11:42
A few questions :) xtreme2k Lounge 59 2002-10-31 06:20

All times are UTC. The time now is 10:40.

Tue Apr 13 10:40:46 UTC 2021 up 5 days, 5:21, 1 user, load averages: 2.04, 1.67, 1.54

Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2021, Jelsoft Enterprises Ltd.

This forum has received and complied with 0 (zero) government requests for information.

Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.2 or any later version published by the Free Software Foundation.
A copy of the license is included in the FAQ.