RE: [SciPy-dev] fft code for testing

Hi,
-----Original Message----- From: Pearu Peterson [mailto:pearu@scipy.org] Sent: Sat 5/31/2003 11:49 AM To: scipy-dev@scipy.net Cc: Subject: RE: [SciPy-dev] fft code for testing
On Sat, 31 May 2003, Chuck Harris wrote:
Btw, does anybody know an efficient way how to test if an integer is a power-of-two in C? Currently djbfft wrapper uses switch (n) { case 2:;case 4:;... case 8192: .. } but there must be a better way.
Fooling around a bit, I think the construction
int foo(int n) { int ok = 1;
if ( n == 1 || n == 2 || n == 4 || n == 8 ); else ok = 0; return ok; }
produces the best looking assembly here.
There are 31 C int's that are power of two (assuming 32 bit machines). Though, may be only the first 24 or so are used in real applications; for instance, the size of double complex array of length 2**24 is 256MB.
So, when n is not a power-of-two, then at least 31 C int comparisons are required. I wonder if this is the lowest bound of #operations or can we do better?
Probably, but by the time we are doing 1024 point transforms the time spent in checking doesn't matter, whereas for the small numbers you can hardly beat the direct comparisons. Just like sequential searches can be faster than binary searches for small lists. I agree it's brute force, and if more information is needed, like _what_ power of two, then it no longer looks so appealing. Also, the assumption is that the number _will_ be a power of two, otherwise it is an error and we don't much care about the time. In the general case of random integers your point is well taken. The loop approach is certainly more adaptable to different sized integers and will produce a slightly smaller function. I have no good ideas of how to do this without some sort of search. More idle curiousity: int foo(int n) { int i, ok = 0; if (n <= (1<<30) ) { for(i = 1; i < n; i <<= 1); if (i == n) ok = 1; } return ok; } produces foo: pushl %ebp movl %esp, %ebp subl $8, %esp movl $0, -4(%ebp) cmpl $1073741824, 8(%ebp) jg .L2 movl $1, -8(%ebp) .L3: movl -8(%ebp), %eax cmpl 8(%ebp), %eax jl .L5 jmp .L4 .L5: leal -8(%ebp), %eax sall $1, (%eax) jmp .L3 .L4: movl -8(%ebp), %eax cmpl 8(%ebp), %eax jne .L2 movl $1, -4(%ebp) .L2: movl -4(%ebp), %eax leave ret gcc screws up the loop a bit, which is why looking at assembly output is bad for my blood pressure :) Chuck _______________________________________________ Scipy-dev mailing list Scipy-dev@scipy.net http://www.scipy.net/mailman/listinfo/scipy-dev

On Sat, 31 May 2003, Chuck Harris wrote:
So, when n is not a power-of-two, then at least 31 C int comparisons are required. I wonder if this is the lowest bound of #operations or can we do better?
Probably, but by the time we are doing 1024 point transforms the time spent in checking doesn't matter, whereas for the small numbers you can hardly beat the direct comparisons.
... which gives me an idea to factor the comparison sequence by introducing some auxiliary tests: int foo(int n) { int ok = 1; if (n<0x100) if (n==0x1 || n==0x2 || n==0x4 || n==0x8 || n==0x10 || n==0x20 || n==0x40 || n==0x80); else ok = 0; else if (n<0x10000) if (n==0x100 || n==0x200 || n==0x400 || n==0x800 || n==0x1000 || n==0x2000 || n==0x4000 || n==0x8000); else ok = 0; else if (n==0x10000 || n==0x20000 || n==0x40000 || n==0x80000 || n==0x100000 || n==0x200000 || n==0x400000 || n==0x800000 || n==0x1000000 || n==0x2000000 || n==0x4000000 || n==0x8000000 || n==0x10000000 || n==0x20000000 || n==0x40000000); else ok = 0; return ok; } and as a result, much less operations are required for random small integers. It would be a nice homework to find optimal steps for auxiliary less-tests by taking also into account the distribution of fft sizes used in real applications... ;-) Pearu
participants (2)
-
Chuck Harris
-
Pearu Peterson