Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Polynomial Multiplication with NTT

  1. Oct 28, 2006 #1
    Hello everyone, I'm new to this forum, and I've registered to seek advice on this topic. What I'm attempting to do is multiply two very large numbers together, represented as polynomials, by the use of the number theory transform (ntt). I've done much research on this topic and I've had some success, but I have run into problems as well.

    So far, I've coded a (hopefully) working ntt algorithm in c++ which is based off of some of the code found in the apfloat package (http://www.apfloat.org/apfloat/). With this algorithm, I can successfully ntt a polynomial (represented as an array of numbers), and then inverse ntt the polynomial back to its original state.

    Seeing as I could do that, I proceded to the next step in my project. I ntt'd a polynomial, squared each value, and then inverse ntt'd the new polynomial. This should yield a polynomial that represents the square of the original number, as this method works with a fast fourier transform. Here's where my trouble lies. For some reason, sometimes I do get the correct polynomial as a result, and sometimes I don't. There seems to be no real pattern in the results that I get. If I use a different base for the polynomial, sometimes certain numbers will work when they didn't in another base, but no one base always yields correct results. I would like to know if the error lies in my ntt algorithm, or elsewhere in my method. I was attempting this with a modulus of 5, a length of 4, and a primitive root of 2 for the ntt, but i get similar results for a modulus of 17, a length of 16, and a primitive root of 3 as well.

    As I said before, I know that with the fast fourier transform, you can fft a polynomial, square each value, and inverse fft the new polynomial to get the square of the original number. My whole project relies on the premise that this also works for the number theory transform, which I believe it should, shouldn't it? Below is the code for my number theory transform. Any help with figuring out what is wrong would be much appreciated.
    Thank you,
    Mike

    (tested and compiles with dev c++)
    Code (Text):

    #include <iostream>
    using namespace std;



    // evaluates exponents, with a modulus
    long power(long base,long power,long mod)
    {
           long final=1;
           for (int i=0; i<power; i++)
           {
               final = ((final%mod) * (base%mod)) % mod;
           }
           return final;
    }



    // prints an array
    void print(long x[],int N)
    {
         for (int i=0; i<N; i++)
         {
             cout<<x[i]<<endl;
         }
         cout<<endl;
    }



    // number theory transform algorithm
    void ntt(long data[],long N,long mod,long pri)
    {
        long k;
        long j;
        long w;
        long wj=1;
        long wk=1;
        long temp[N];

        w = power(pri,((mod - 1)/N) % mod,mod) % mod;

        for (j=0; j<N; j++)
        {
            temp[j]=0;
            wk=1;
            for (k=0; k<N; k++)
            {
                temp[j] = ((temp[j]%mod) + ((wk%mod)*(data[k]%mod))%mod) % mod;
                wk = ((wk%mod)*(wj%mod)) % mod;
            }
            wj = ((wj%mod)*(w%mod)) % mod;
        }

        for (k=0; k<N; k++)
        {
            data[k] = temp[k] % mod;
        }

    }



    // inverse number theory transform algorithm
    void intt(long data[],long N,long mod,long pri)
    {
         long k;
         long j;
         long w;
         long wk=1;
         long wj=1;
         long temp[N];
         long rN;
         
         w = power(pri,(mod-1-((mod-1)/N)) % mod,mod) % mod;
         rN = power(N,(mod-2)%mod,mod) % mod;
         
         for (j=0; j<N; j++)
         {
             temp[j]=0;
             wk=1;
             for (k=0; k<N; k++)
             {
                temp[j] = ((temp[j]%mod) + ((wk%mod)*(data[k]%mod))%mod) % mod;
                wk = ((wk%mod)*(wj%mod)) % mod;
             }
             wj = ((wj%mod)*(w%mod)) % mod;
         }
         
         for (k=0; k<N; k++)
         {
             data[k]=((temp[k]%mod)*(rN%mod)) % mod;
         }
         
    }



    // main
    int main()
    {
        long N=16;
        long mod=17;
        long pri=3;
        long data[N];

        data[0]=3;
        data[1]=1;
        data[2]=0;
        data[3]=0;
        for (int i=4; i<N; i++)
        {
            data[i]=0;
        }
       
        print(data,N);
        ntt(data,N,mod,pri);
        print(data,N);
        for (int i=0; i<N; i++)
        {
            data[i]=data[i]*data[i];
        }
        print(data,N);
        intt(data,N,mod,pri);
        print(data,N);
       
        cout<<endl<<endl<<endl;
        system("PAUSE");
        return 0;
    }
     
     
  2. jcsd
  3. Oct 28, 2006 #2

    Hurkyl

    User Avatar
    Staff Emeritus
    Science Advisor
    Gold Member

    Can you give examples of times you get the wrong answer?
     
  4. Oct 28, 2006 #3
    Yes, here is a list of results that i get in base 16 and 17, for an ntt with a modulus of 17, length of 16, and primitive root of 3.

    Code (Text):

                   base16  base17
    num   num^2    result  result
    16     256      256    
    17     289      289     289
    18     324      324     324
    19     361      361     361
    20     400      400     400
    21     441      424     441
    22     484      450     467
    23     529      495     495
    24     576      525     542
    25     625      285     574
    26     676      319     319
    27     729      338     355
    28     784      376     376
    29     841      416     416
    30     900      441     458
    31     961      468     485
    32     1024     1024    514
    33     1089     1089    545
    34     1156     1156    1156
     
    so...
    in base 16, the results deviate from num^2 when num=21,22,23,24,25,26,27,28,29,30,31
    in base 17, the results deviate from num^2 when num=22,23,24,25,26,27,28,29,30,31,32,33

    I guess there is a pattern with the incorrect results, but i don't know what to make of it or how to fix the problem.
     
    Last edited: Oct 28, 2006
  5. Oct 29, 2006 #4

    Hurkyl

    User Avatar
    Staff Emeritus
    Science Advisor
    Gold Member

    What is "num" and "result" in that table?


    I was looking for a specific example of a polynomial (and the corresponding N, mod, and pri), and the polynomial your code thought was its square.


    By the way, did you notice that all of the differences are divisible by 17? e.g. for num = 26: 676 - 319 = 357 = 17 * 21.


    One of the things I was suspecting is that your code gives the right answer modulo 17... which is the best it can hope to do because all the arithmetic is done modulo 17.

    I had wanted to see an actual polynomial so I could test that hypothesis, but that results column suggests it too.
     
  6. Oct 29, 2006 #5

    shmoe

    User Avatar
    Science Advisor
    Homework Helper

    When you multiply two polynomials using a number theory transform mod 17, you get the product of these two polynomials mod 17 (please correct me if I'm wrong, I can't say I've ever looked at a ntt before today-keep this in mind when reading this post i could be off). For example, when squaring 21 using base 16, you used 21=1*16+5. This corresponds to squaring the polynomial x+5. The answer you want is x^2+10x+25, but mod 17 the square of this polynomial is x^2+10x+8, which explains your off by 17.

    Mod 17 is a different story, 21=1*17+4. The square of x+4 is x^2+8x+16, which is also it's square mod 17. Essentially you get problems when the square of your coefficients is 17 or greater. (other problems will occur when your transform isn't long enough to get the highest powers of the polynomial). By looking at the coefficients in your base 16 and 17 representations of "num" in your table, it should be clear why the rows your mod 17 ntt fail or succeed.

    There's a few ways around this. First you can represent your numbers by a small enough base with respect to your transform to ensure the square of the coefficients is never too large. Try squaring those numbers again using the same ntt but working base 5 (e.g. 21=4*5+1, you can use your mod 17 ntt to square 4*x+1 without problems).

    Second, you could use multiple transforms to different bases and the chinese remainder theorem to find the coeffcients of your squared polynomials exactly. for example, if you are representing your numbers in base 16 you wanted to know (x+5)^2 to find 21^2. You know it's x^2+10*x+8 mod 17 from your mod 17 ntt. Do another ntt using the primitive root 2 mod 5, you'll get (x+5)^2=x^2 mod 5 for example (okie, that is obvious without the transform). By the Chinese remainder theorem you'd then know the coefficients of (x+5)^2 mod 17*5=85, which is sufficient here.


    edit-had "mod" where I should have had "base" and a 16 that should have been a 17
     
    Last edited: Oct 29, 2006
  7. Oct 29, 2006 #6
    Thank you shmoe, I understand exactly what's going wrong now. In order to be absolutely sure that a polynomial can be squared using this method without any problems, I think that this inequality must evaluate true:

    [tex](d+1)(m-1)^2 < p[/tex]

    where d is the degree of the polynomial, m is the base used, and p is the prime modulus used in the ntt. This is because m-1 is the largest possible diget value for the given base, and d+1 times the square of the largest possible diget is the size of the largest possible coefficient when the polynomial is squared. The largest possible coefficient must be smaller than p for it to work.

    According to this, even a base of 5 would not be good enough to safely square all first degree polynomials when p=17 because (1+1)(5-1)^2 = 32 > 17. The largest base that you can use to safely square polynomials when p=17 is 3! (not factorial :tongue:) And that's only talking about first degree polynomials, i.e. only two digets in base 3.

    Thanks for the help guys, I'll let you know if I run into any more snags along my journey.
     
    Last edited: Oct 30, 2006
  8. Oct 29, 2006 #7

    Hurkyl

    User Avatar
    Staff Emeritus
    Science Advisor
    Gold Member

    Oh, whoops! I missed that your goal is to multiply integers, not polynomials. :blushing:


    Hrm, this is interesting: I seem to remember that when squaring an integer n, you usually want to work with a modulus bigger than n^2... but off hand, I don't see anything wrong with your suggestion on a lower bound for the modulus!
     
    Last edited: Oct 29, 2006
  9. Nov 15, 2006 #8
    Multiplication works perfectly now, but i'm wondering about division. If i wish to divide two sequences of numbers with an ntt, would I just ntt those two sequences, and multiply the first sequence's terms by the modular inverse of the second sequence's terms? Is division even possible to do with an ntt (i know it's possible with an ftt)? Any help would be appreciated, thanks.
     
  10. Nov 15, 2006 #9

    Hurkyl

    User Avatar
    Staff Emeritus
    Science Advisor
    Gold Member

    I can tell you what I know about polynomials. When using one of these evaluate/interpolate algorithms, you can only do exact division: remainders must be zero.

    Since you're really working with integers, but converting to polynomials, you have a very difficult problem: not only does the integer division have to be exact, but so does the polynomial division! :frown:

    Of course, I would try coding it up to see if it works; sometimes the heavens align and good things happen. :smile:
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?



Similar Discussions: Polynomial Multiplication with NTT
  1. Feral Polynomials (Replies: 3)

Loading...