funcs.h

Go to the documentation of this file.

// Copyright (C) 2011 CNRS - Ecole Polytechnique - INRIA.
//
// This file is part of TIFA.
//
// TIFA is free software; you can redistribute it and/or modify it under the 
// terms of the GNU Lesser General Public License as published by the Free
// Software Foundation; either version 2.1 of the License, or (at your option)
// any later version.
//
// TIFA 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 Lesser General Public License for
// more details.
//
// You should have received a copy of the GNU Lesser General Public License
// along with this library; if not, write to the Free Software Foundation,
// Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301 USA.

#if !defined(_TIFA_FUNCS_H_)

#define _TIFA_FUNCS_H_

#ifdef __cplusplus
extern "C" {
#endif

#include <limits.h>
#include <inttypes.h>
#include <math.h>
#include <gmp.h>

#include "macros.h"
#include "array.h"
#include "tifa_config.h"

#define LARGEST_MULTIPLIER 97

#define BITSIZE_LARGEST_MULTIPLIER 7

#define MAX_IPRIME_IN_MULT_CALC 31

/*
 *-----------------------------------------------------------------------------
 *                      Structures and typedefs
 *-----------------------------------------------------------------------------
 */

struct struct_mult_data_t {
    uint32_t multiplier;
    uint32_t count;
    double   sum_inv_pi;
};

typedef struct struct_mult_data_t mult_data_t;

/*
 *-----------------------------------------------------------------------------
 *                      Miscellaneous functions
 *-----------------------------------------------------------------------------
 */

uint32_t most_significant_bit(uint32_t n);

inline static uint32_t floor_log2(uint32_t n) {
    return most_significant_bit(n);
}

inline static uint32_t ceil_log2(uint32_t n) {
    uint32_t e = most_significant_bit(n);
    if (! IS_POWER_OF_2_UI(n)) {
        e++;
    }
    return e;
}

inline static uint32_t ceil_log2_mp_limb(mp_limb_t limb) {
#if GMP_LIMB_BITS == 64
    uint64_t tmp = limb & GMP_NUMB_MASK;
    if (tmp >> 32) {
        return 32 + ceil_log2((uint32_t)(tmp >> 32));
    } else {
        return ceil_log2((uint32_t)tmp);
    }
#else
    uint32_t tmp = (uint32_t)limb & GMP_NUMB_MASK;
    return ceil_log2((uint32_t)tmp);
#endif
}

void find_coprime_base(mpz_array_t* const base, const mpz_t n,
                       const mpz_array_t* const factors);

/*
 *-----------------------------------------------------------------------------
 *                      Number theoretical functions
 *-----------------------------------------------------------------------------
 */

#define NO_SQRT_MOD_P (UINT32_MAX)

#define NO_SQRT_MOD_P2 (ULONG_MAX)

int8_t kronecker_ui(uint32_t a, uint32_t b);

uint32_t powm(uint32_t base, uint32_t power, uint32_t modulus);

uint32_t sqrtm(uint32_t a, uint32_t p);

extern const unsigned short qres_mod_221[221] MAYBE_UNUSED;

extern const unsigned short qres_mod_256[256] MAYBE_UNUSED;

extern const unsigned short qres_mod_315[315] MAYBE_UNUSED;

inline static unsigned long int is_square(unsigned long int x) {
    //
    // Basic perfect square detection test.
    //
    // See for example algorithm 1.7.3 from the book "A Course in Computational
    // Algebraic Number Theory" by Henri Cohen, Springer-Verlag 1993.
    //
    // The description given in this book has been adapted to use larger
    // tables in exchange of a slight performance boost (about 30% on
    // Opteron 250).
    //
    if (qres_mod_256[x & 255] == 0) {
        //
        // Get rid of about 82.8% of non squares
        //
        return 0;
    }
    if (qres_mod_315[x % 315] == 0) {
        //
        // Get rid of about 84.8% of non squares
        //
        return 0;
    }
    if (qres_mod_221[x % 221] == 0) {
        //
        // Get rid of about 71.5% of non squares
        //
        return 0;
    }
    unsigned long int root = (unsigned long int)sqrt(x);
    if ((root * root) == x) {
        return root;
    }
    return 0;
}

bool is_prime(uint32_t n);

unsigned long int gcd_ulint(unsigned long int a, unsigned long int b);

unsigned long int modinv_ui(unsigned long int n, unsigned long int p);

unsigned long int sqrtm_p2(uint32_t a, uint32_t p);

uint32_t ks_multiplier(const mpz_t n, const uint32_t size_base);

/*
 *-----------------------------------------------------------------------------
 *                          Hash functions
 *-----------------------------------------------------------------------------
 */

uint32_t hash_rj_32(const void* const keyptr);

uint32_t hash_pjw(const void* const keyptr);

uint32_t hash_sfh_ph(const void* const keyptr);

/*
 *-----------------------------------------------------------------------------
 *                          Comparison functions
 *-----------------------------------------------------------------------------
 */

int mpz_cmp_func(const void* const mpza, const void* const mpzb);

int uint32_cmp_func(const void* const uinta, const void* const uintb);

int string_cmp_func(const void* const stra, const void* const strb);

int cmp_mult_data(const void* mda, const void* mdb);

/*
 *-----------------------------------------------------------------------------
 *                        Combinatorial functions
 *-----------------------------------------------------------------------------
 */

uint32_t n_choose_k(uint8_t n, uint8_t k);

void next_subset_lex(uint32_t n, uint32_t k, uint32_t* subset, bool* end);

void unrank_subset_lex(uint32_t n, uint32_t k, uint32_t r, uint32_t* subset);

void unrank_subset_revdoor(uint32_t n, uint32_t k, uint32_t r,
                           uint32_t* subset);

/*
 *-----------------------------------------------------------------------------
 *                      Random number generation functions
 *-----------------------------------------------------------------------------
 */

void tifa_srand(uint32_t seed);

uint32_t tifa_rand();

#ifdef __cplusplus
}
#endif

#endif

---