divisors.hpp

Go to the documentation of this file.

/* -*- coding: latin-1 -*- */
/** \file common/sigmaodd/divisors.hpp (January 7, 2018)
 * \brief
 * Functions in link with divisor notion:
 * sum of divisors, factorization, GCD, coprime, ...
 *
 * GPLv3 --- Copyright (C) 2017, 2018 Olivier Pirson
 * http://www.opimedia.be/
 */

#ifndef PROGS_SRC_COMMON_SIGMAODD_DIVISORS_HPP_
#define PROGS_SRC_COMMON_SIGMAODD_DIVISORS_HPP_

// \cond
#include <cstdint>

#include <iostream>
#include <utility>
#include <vector>
// \endcond

#include "helper.hpp"


namespace sigmaodd {

  /* ********
   * Struct *
   **********/
  /** \brief
   * Structure to represent a factor with its exponent.
   */
  struct FactorExp {
    nat_type factor;
    unsigned int exp;

    constexpr
    bool
    operator==(const FactorExp &other) const;

    constexpr
    nat_type
    value() const;
  };



  /* ************
   * Prototypes *
   **************/
  /** \brief
   * Return "factor^exponent" representation of factor_exp.
   */
  std::ostream& operator<<(std::ostream &out, const FactorExp &factor_exp);


  /** \brief
   * Heuristic to return a list of factors with their exponents
   * such that all factors are coprimes
   * and the product == a*b.
   * If failed to have a list at least two factors then return an empty list.
   *
   * @param a > 1
   * @param b > 1
   */
  std::vector<FactorExp>
  coprime_factor_exps(nat_type a, nat_type b);


  /** \brief
   * Divide n by d until the result is not divisible by d,
   * and return (result, number of divisions).
   *
   * @param n != 0
   * @param d > 1
   *
   * @return (result, number of divisions)
   */
  std::pair<nat_type, unsigned int>
  divide_until_nb(nat_type n, nat_type d);


  /** \brief
   * Return n divided by 2 until the result is odd.
   *
   * @param n != 0
   */
  inline
  nat_type
  divide_until_odd(nat_type n);


  /** \brief
   * Divide n by 2 until the result is odd,
   * and return (result, number of divisions).
   *
   * @param n != 0
   *
   * @return (result, number of divisions)
   */
  inline
  std::pair<nat_type, unsigned int>
  divide_until_odd_nb(nat_type n);


  /** \brief
   * Return the first divisor of n > 1 (or 1 if n <= 1)
   */
  nat_type
  first_divisor(nat_type n);


  /** \brief
   * Return the Greatest Common Divisor of a and b.
   *
   * If a and b == 0,
   * then return 0.
   *
   * Return the Greatest Common Divisor of a and b
   */
  constexpr
  nat_type
  gcd(nat_type a, nat_type b);


  /** \brief
   * Return true iff a and b are coprime (relatively prime).
   */
  constexpr
  nat_type
  is_coprime(nat_type a, nat_type b);


  /** \brief
   * Return true iff d divide n,
   * i.e. if n is divisible by d.
   *
   * @param d != 0
   * @param n
   */
  constexpr
  bool
  is_divide(nat_type d, nat_type n);


  /** \brief
   * Return true iff at least one of
   * 3, 7, 31 or 127
   * is an unitary divisor of n.
   *
   * These numbers are Mersenne primes.
   *
   * List of prime Mersenne numbers:
   * https://oeis.org/A000668
   */
  constexpr
  bool
  is_first_mersenne_prime_unitary_divide(nat_type n);


  /** \brief
   * Return true iff
   * is_first_mersenne_prime_unitary_divide(n) or is_square(n).
   */
  constexpr
  bool
  is_first_mersenne_prime_unitary_divide_or_square(nat_type n);


  /** \brief
   * Return true iff at least one of
   * 3, 7, 31, 127, 8191, 131071 or 524287
   * is an unitary divisor of n.
   *
   * These numbers are Mersenne primes.
   *
   * List of prime Mersenne numbers:
   * https://oeis.org/A000668
   */
  constexpr
  bool
  is_little_mersenne_prime_unitary_divide(nat_type n);


  /** \brief
   * Return true iff at least one of
   * 3, 7, 31, 127, 8191, 131071, 524287, 2147483647 or 2305843009213693951
   * is an unitary divisor of n.
   *
   * These numbers are all the Mersenne primes < 2^64.
   *
   * List of prime Mersenne numbers:
   * https://oeis.org/A000668
   */
  constexpr
  bool
  is_mersenne_prime_unitary_divide(nat_type n);


  /** \brief
   * Return true iff n is a perfect square.
   */
  bool
  is_square(nat_type n);


  /** \brief
   * Return true iff (d divide n) and (d NOT divide n/d).
   *
   * @param d != 0
   * @param n
   */
  constexpr
  bool
  is_unitary_divide(nat_type d, nat_type n);


  /** \brief
   * Return pollard_rho(n, rand() % n, floor square root of n).
   *
   * @param n != 0
   */
  nat_type
  pollard_rho(nat_type n);


  /** \brief
   * Heuristic to find a proper divisor of n with the Pollard's rho heuristic.
   *
   * random is the "random" number used by the heuristic
   * and max_iteration the maximum number of iterations to find the divisor.
   *
   * If a divisor is find
   * then return it,
   * else return 0.
   *
   * @param n != 0
   * @param random < n
   * @param max_iteration
   *
   * @return 0 or 1 < (divisor of n) < n
   */
  nat_type
  pollard_rho(nat_type n, nat_type random, nat_type max_iteration);


  /** \brief
   * Try pollard_rho(n) a maximum of nb_tries times
   * and return 0 or the first divisor find.
   *
   * @param n != 0
   * @param nb_tries
   *
   * @return 0 or 1 < (divisor of n) < n
   */
  nat_type
  pollard_rho_repeat(nat_type n, nat_type nb_tries);


  /** \brief
   * Calculates the sum of all divisors of n by the factorization method
   * and returns it.
   *
   * @param n != 0
   *
   * @return the sum of all divisors of n
   */
  nat_type
  sum_divisors__factorize(nat_type n);


  /** \brief
   * Calculates the sum of all divisors of n by the naive method
   * and returns it.
   *
   * @param n != 0
   *
   * @return the sum of all divisors of n
   */
  nat_type
  sum_divisors__naive(nat_type n);


  /** \brief
   * Calculates the sum of odd divisors of n by the factorization method
   * and returns it.
   *
   * https://oeis.org/A000593
   *
   * @param n != 0
   *
   * @return the sum of all odd divisors of n
   */
  nat_type
  sum_odd_divisors__factorize(nat_type n);


  /** \brief
   * Calculates the sum of odd divisors of n by the naive method
   * and returns it.
   *
   * @param n != 0
   *
   * @return the sum of all odd divisors of n
   */
  nat_type
  sum_odd_divisors__naive(nat_type n);


  /** \brief
   * Return an upper bound of sum odd divisors.
   *
   * @param n != 0
   */
  nat_type
  sum_odd_divisors_upper_bound(nat_type n);  // ???


  /** \brief
   * Return an upper bound of sum odd divisors of  n
   * using the DeTemple inequalities.
   *
   * @param n != 0
   */
  nat_type
  sum_odd_divisors_upper_bound__DeTemple(nat_type n);


  /** \brief
   * Return an upper bound of sum odd divisors of  n
   * using the DeTemple inequalities. ???
   *
   * @param n != 0
   */
  nat_type
  sum_odd_divisors_upper_bound__DeTemple_ceil(nat_type n);


  /** \brief
   * Return an upper bound of sum odd divisors of n
   * using the DeTemple inequalities
   * and the knowledge that n have no divisors <= k.
   *
   * @param n != 0
   * @param k odd <= sqrt(n)
   */
  nat_type
  sum_odd_divisors_upper_bound__DeTemple(nat_type n, nat_type k);  // ???

}  // namespace sigmaodd

#include "divisors__inline.hpp"

#endif  // PROGS_SRC_COMMON_SIGMAODD_DIVISORS_HPP_