stack space in gcd operations

[David Harvey]

On Aug 16, 2009, at 4:07 PM, Niels Möller wrote:

> David Harvey <dmharvey at cims.nyu.edu> writes:
>
>> The code below creates a new thread with the requested amount of
>> stack space and attempts an n-limb mpz_invert operation. On OS X
>> this fails (bus error) e.g. when the stack is 448 KB and n = 2800,
>
>> It is reasonable for GMP to require so much stack space for xgcd?
>
> I don't think so. Direct stack allocation (automatic C variables and
> stuff) ought to be O(log(n)) with a fairly small constant factor,
> without reading the code, I'd guess around a hundred bytes or smaller
> per factor of two in n.
>
> Then there's allocation by TMP_ALLOC, which tries to allocate "small"
> (currently, up to 64 KB, it seems) chunks on the stack, and large
> chunks on the heap. For gcdext, all storage is allocated up front in
> one chunk by mpn_gcdext, based on various ITCH-macros.
>
> This allocation does not include temporary storage for the udnerlying
> multiplications, though, since the current multiplication functions
> don't use an itch/scratch interface. But I think that storage should
> be limited to, say, 2n, or 44 KB if you have 64-bit limbs.
>
> So maybe you could try to instrument the TMP_ALLOC calls to see what
> sizes are requested, and then check if they're allocated from the heap
> or the stack. Or you could try reconfiguring with --disable-alloca and
> see if that helps.

I tried intercepting a few TMP_ALLOC calls (using tal-debug.c) to see  
what actually happens.

Here is an example with n = 8000.

mpz/invert.c:46: allocate approx n limbs = 64KB (this is just under  
the 64KB limit, so goes on the stack).
mpz/invert.c:47: allocate approx n limbs = 64KB
mpz/gcdext.c:66: allocate approx n limbs = 64KB
mpz/gcdext.c:67: allocate approx n limbs = 64KB
mpz/gcdext.c:76: allocate approx n limbs = 64KB
mpz/gcdext.c:77: allocate approx n limbs = 64KB
mpn/tdiv_qr.c:369: allocate approx n limbs = 64KB (which are freed  
immediately afterwards!)

So far the peak is already 448 KB, and we haven't even got to the  
guts of the GCD yet!!!

The rest of the trace is rather long and boring, but if my code is  
correct, the peak for the whole mpz_invert call is 689544 bytes.

david