int.c
/*
* Copyright (c) 2002, 2017 Jens Keiner, Stefan Kunis, Daniel Potts
*
* This program is free software; you can redistribute it and/or modify it under
* the terms of the GNU General Public License as published by the Free Software
* Foundation; either version 2 of the License, or (at your option) any later
* version.
*
* This program is distributed in the hope that it will be useful, but WITHOUT
* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
* FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
* details.
*
* You should have received a copy of the GNU General Public License along with
* this program; if not, write to the Free Software Foundation, Inc., 51
* Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
*/
#include "infft.h"
INT Y(exp2i)(const INT a)
{
return (1U << a);
}
INT Y(log2i)(const INT m)
{
/* Special case, zero or negative input. */
if (m <= 0)
return -1;
#if SIZEOF_PTRDIFF_T == 8
/* Hash map with return values based on De Bruijn sequence. */
static INT debruijn[64] =
{
0, 58, 1, 59, 47, 53, 2, 60, 39, 48, 27, 54, 33, 42, 3, 61, 51, 37, 40, 49,
18, 28, 20, 55, 30, 34, 11, 43, 14, 22, 4, 62, 57, 46, 52, 38, 26, 32, 41,
50, 36, 17, 19, 29, 10, 13, 21, 56, 45, 25, 31, 35, 16, 9, 12, 44, 24, 15,
8, 23, 7, 6, 5, 63
};
register uint64_t v = (uint64_t)(m);
/* Round down to one less than a power of 2. */
v |= v >> 1;
v |= v >> 2;
v |= v >> 4;
v |= v >> 8;
v |= v >> 16;
v |= v >> 32;
/* 0x03f6eaf2cd271461 is a hexadecimal representation of a De Bruijn
* sequence for binary words of length 6. The binary representation
* starts with 000000111111. This is required to make it work with one less
* than a power of 2 instead of an actual power of 2.
*/
return debruijn[(uint64_t)(v * 0x03f6eaf2cd271461LU) >> 58];
#elif SIZEOF_PTRDIFF_T == 4
/* Hash map with return values based on De Bruijn sequence. */
static INT debruijn[32] =
{
0, 9, 1, 10, 13, 21, 2, 29, 11, 14, 16, 18, 22, 25, 3, 30, 8, 12, 20, 28,
15, 17, 24, 7, 19, 27, 23, 6, 26, 5, 4, 31
};
register uint32_t v = (uint32_t)(m);
/* Round down to one less than a power of 2. */
v |= v >> 1;
v |= v >> 2;
v |= v >> 4;
v |= v >> 8;
v |= v >> 16;
/* 0x07C4ACDD is a hexadecimal representation of a De Bruijn sequence for
* binary words of length 5. The binary representation starts with
* 0000011111. This is required to make it work with one less than a power of
* 2 instead of an actual power of 2.
*/
return debruijn[(uint32_t)(v * 0x07C4ACDDU) >> 27];
#else
#error Incompatible size of ptrdiff_t
#endif
}
INT Y(next_power_of_2)(const INT x)
{
if (x < 0)
return -1;
else if (x < 2)
return x + 1;
else
{
uint64_t v = (uint64_t)x;
/* Subtract one, so we return the input if it is a power of two. */
v--;
/* Round down to one less than a power of two. */
v |= v >> 1;
v |= v >> 2;
v |= v >> 4;
#if SIZEOF_PTRDIFF_T >= 2
v |= v >> 8;
#endif
#if SIZEOF_PTRDIFF_T >= 4
v |= v >> 16;
#endif
#if SIZEOF_PTRDIFF_T >= 8
v |= v >> 32;
#endif
/* Add one to get power of two. */
v++;
return v;
}
}
void Y(next_power_of_2_exp)(const INT N, INT *N2, INT *t)
{
INT n,i,logn;
INT N_is_not_power_of_2=0;
if (N == 0)
{
*N2 = 1;
*t = 0;
}
else
{
n = N;
logn = 0;
while (n != 1)
{
if (n%2 == 1)
{
N_is_not_power_of_2=1;
}
n = n/2;
logn++;
}
if (!N_is_not_power_of_2)
{
logn--;
}
for (i = 0; i <= logn; i++)
{
n = n*2;
}
*N2 = n;
*t = logn+1;
}
}
void Y(next_power_of_2_exp_int)(const int N, int *N2, int *t)
{
int n,i,logn;
int N_is_not_power_of_2=0;
if (N == 0)
{
*N2 = 1;
*t = 0;
}
else
{
n = N;
logn = 0;
while (n != 1)
{
if (n%2 == 1)
{
N_is_not_power_of_2=1;
}
n = n/2;
logn++;
}
if (!N_is_not_power_of_2)
{
logn--;
}
for (i = 0; i <= logn; i++)
{
n = n*2;
}
*N2 = n;
*t = logn+1;
}
}
infft.h
Internal header file for auxiliary definitions and functions.