Factorial Using FFT

Question

I'm trying to implement a program in C that calculates the factorial of a very large n (up to a million), using fft and binary splitting method.

I've implemented a simple library to represent arbitrary precision integer. To calculate the fft and ifft, i use twofft.c and four1.c routines from "Numerical Recipes in C"

Up to a certain n, all goes right, but when the numbers (floating arrays) are too big, the ifft (calculate with four1),after normalization and rounding, has values that are wrong.

For example, if i have two number with 2000 digits that ends with 40 zeros, and i have to multiply them each other (using fft), when i calculate the ifft, some ending zeros become "one". this happens because when i rounded one of this "zeros", (0,50009 for examples), they became "one". Now, i don't know if is my implementation wrong or if i have to rounding this numebrs in a different way. I've tried to use both binary split method and prime factorization, but for n >= 9000, the result is wrong.

there is a way to resolve this? thanks for your attention and sorry for my bad english.

Answer 1

How do you represent arbitrary precision integers?

I mean what type are you actually using?

Can you please show us your code?

If you feel really lazy you can clone this project i've made few months ago: https://github.com/nomadster/ESP

Edit:

By further reading your post i suppose by this statement

"this happens because when i rounded one of this "zeros", (0,50009 for examples), they became "one""

that you are still unaware of the fact that fft multiplication only works when the roundoff error is smaller than 0.5. So it seems to me (if and only if i've correctly interpreted your cryptic message) that you are using a floating point type that doesn't have the required precision.