Optimizing Modulus (was Re: Would someone mind elaborating the explanation in the manual?)

[David McKen]

Hmm, that might work.

I am also trying a speed-space analysis of my program and I would
like to know how much memory is used by a mpz_t? I work in loads of
about 1000000 to 10000000 numbers a load so breaking this up like
that would probably be quite efficient. 

This chunk of code is probably the main chunk that called that mod
function so much. I only posted the first chunk because alot of the
numbers won't even get past it (that and the fact that I forgot
completly about this loop).  The size of the list of prime numbers is
variable. I know it is making some difference as the bigger I make
the list the faster the entire process gets (the upper limit
unfortunately is the amount of memory in the system, once it starts
swapping to disk things slow down). Right now I am testing using all
the prime numbers between 100 and 1,000,000 (this was used when
calculating all my past times).

--Begin
//Final optimization. Use a list of prime numbers to further weed out
//Unnesesary divisions.
break_loop = false;
curr_prime = prime_num_list.begin();
while(curr_prime != prime_num_list.end()) {
  remainder = mpz_divisible_ui_p(up_counter, *curr_prime);
  if(remainder) {
    curr_prime++;
  } else {
    break_loop = true;
    break;
  }
}
if(break_loop) {
  update_func(up_counter, up_counter, 2);
  continue;
}
--End

If the numbers get past both my hard coded checks (divisibility by
all primes less than 100, repeat of the if statement I posted
earlier) and this loop (all other primes greater then 100 less than
whatever my limit is) then I get lazy and feed it into the
divisibility function that comes with GMP (which is really slow).

Another question:
How fast is the C++ interface? OOP is usually slower but as far as I
can tell most wrappers if they are inline show just about no real
performance difference. Is GMP different?

--- Sisyphus <kalinabears at iinet.net.au> wrote:
> David McKen wrote:
> 
> > 
> > --Begin
> >       remainder = mpz_fdiv_ui (up_counter, prime_product_3_29);
> >       if(!((remainder %  3) &&
> >            (remainder %  5) &&
> >           (remainder %  7) &&
> >            (remainder % 11) &&
> >           (remainder % 13) &&
> >            (remainder % 17) &&
> >            (remainder % 19) &&
> >            (remainder % 23) &&
> >            (remainder % 29)) )
> >       {
> >         update_func (up_counter, up_counter, 2);
> >         continue;
> >       }
> > --End
> > 
> 
> I'm assuming that up_counter is odd.
> Looks to me that what you really should be doing is sieving a range
> of 
> numbers (let's say the odd numbers in the range up_counter to 
> up_counter+10000).
> 
> You then find the first number in the range that is divisible by 3,
> and 
> cross it off, and cross off every third number after it. Then find
> the 
> first number in the range that is divisible by 5 and cross it off,
> and 
> cross off every fifth number after it. And so on for 7, 11, ... 29.
> 
> By then, the only numbers that haven't been crossed off are the
> numbers 
> that are not divisible by [3..29].
> 
> You keep track of which numbers have been "crossed off" using a bit
> 
> vector (5000 bits long) initialised with all bits set to (say) 1.
> "Crossing off" simply involves changing the appropriate bit to 0.
> 
> You can then determine which numbers have been crossed off by
> examining 
> the bits in the vector.
> 
> eg bits 5, 7 and 8 set to 1 would tell you that up_counter+10, 
> up_counter+14 and up_counter+16 are not divisible by any of
> [3..29].
> 
> (In case I'm not being clear, 'up_counter' is now just the initial
> value 
> given - it doesn't get incremented.)
> 
> One inefficiency of this method is that many values are getting
> "crossed 
> off" more than once - which is wasteful. This can be addressed (by
> some 
> sort of 'wheel', I think), but I can't find a link at the moment,
> and 
> it's not something I've personally attempted. I'm not sure of the 
> savings to be made.
> 
> I find I can take (eg) a ten thousand decimal digit odd number, x,
> and 
> find all of the odd numbers in the range [x .. x+100000] not
> divisible 
> by any prime less than 2^20 - in about 2.1 seconds on my 1MHz pc.
> That includes the time it takes to generate the primes in the first
> 
> place with a Sieve of Eratosthenes. You can vary the size of the
> range 
> without greatly affecting the timing. A range of [x .. x+200000]
> adds 
> about 0.01 seconds to the time taken. The big efect is, of course,
> how 
> deeply you sieve. Sieving (in this example) for all primes less
> than 
> 2^25 takes 56 seconds. (I'm not making use of the time-saving
> technique 
> we started off discussing, btw. Something for me to consider.)
> 
> Anyway ... that was just by way of example ... I'm just a wannabe 
> programmer, so if someone asks "how come it's taking me so long ?",
> 
> don't be surprised.
> 
> And it's probably OT for this group. (I've found sci.crypt to be a
> very 
> useful newsgroup for matters pertaining to primes and factoring.)
> Sorry, 
> I can't really help with the on-topic aspects you raised.
> 
> Any thoughts, Kevin, on including a
> 'mpz_sieve_eratosthenes(max_prime)' 
> function with GMP ? I envisage a function that returns a mpz such
> that 
> if (and only if) bit x is set, then 2 * x is prime. Obviously such
> a mpz 
> does not encode the prime number '2' - but if you want to encode
> '2', 
> you double the size of the output mpz.
> 
> My sieve of eratosthenes is just a straight C impelementation,
> where the 
> output vector is, in fact, an array of ints. I don't know if there
> would 
> be significant gains in using the low level GMP functions to create
> an 
> mpz that encodes the same information.
> 
> Another issue might be that it is slower to iterate through an mpz
> bit 
> by bit than it is to iterate through an array of ints bit by bit.
> 
> Cheers,
> Rob

=====
Signed 
David Mcken

Life Sucks
Live with it

__________________________________
Do you Yahoo!?
Exclusive Video Premiere - Britney Spears
http://launch.yahoo.com/promos/britneyspears/