Changes for version 0.74 - 2026-03-20
- API Changes
- is_pseudoprime, is_euler_pseudoprime, and is_strong_pseudoprime will use an implicit base of 2 rather than a missing base error. This also means is_pseudoprime($n, @bases) works properly.
- is_strong_pseudoprime passes for bases equal to 0 mod n. This matches what the GMP code always did, and means primes return 1 for any base.
- divisor_sum(0) returns 0. divisors(0) returns an empty list.
- divisors takes an optional second argument to prune the returned list.
- consecutive_integer_lcm(0) was returning 0. Now returns 1, matching the original GMP API as well as the OEIS series.
- some error messages with invalid arguments have changed.
- Callers must 'use Math::Prime::Util::MemFree;' if they want to use it.
- bernfrac was documented to only accept non-negative inputs. This is now enforced.
- bernfrac and harmfrac return values are only bigints if needed.
- The real functions no longer modify their behavior based on whether the "bignum" package was loaded, whether in scope or not.
- ADDED
- powersum(n,k) sum of k-th powers of positive ints <= n
- sumpowerful(n[,k]) sum of k-powerful numbers <= n
- sumliouville(n) L(n), sum of liouville(1..n)
- sumtotient(n) sum of euler_phi(1..n)
- prime_bigomega(n) # of factors (with multiplicity)
- prime_omega(n) # of factors (distinct)
- powint(a, b) signed integer a^b
- mulint(a, b) signed integer a*b
- addint(a, b) signed integer a+b
- subint(a, b) signed integer a-b
- add1int(n) signed integer n+1
- sub1int(n) signed integer n-1
- divint(a, b) signed integer a/b (floor)
- modint(a, b) signed integer a%b (floor)
- cdivint(a, b) signed integer a/b (ceiling)
- divrem(a, b) Euclidean quotient and remainder
- fdivrem(a, b) Floored quotient and remainder
- cdivrem(a, b) Ceiling quotient and remainder
- tdivrem(a, b) Truncated quotient and remainder
- lshiftint(n, k) left shift n by k bits
- rshiftint(n, k) right shift n by k bits (truncating)
- rashiftint(n, k) right shift n by k bits (flooring)
- absint(n) integer absolute value
- negint(n) integer negation
- cmpint(a, b) integer comparison (like <=>)
- signint(a, b) integer sign (-1,0,1)
- random_safe_prime for n-bit safe primes
- inverse_li_nv(x) float inverse logarithmic integral
- is_gaussian_prime(a,b) is a+bi a Gaussian prime
- is_lucky(n) predicate for lucky number sieve
- is_happy(n) if n happy number, returns height of n
- is_happy(n,base,exponent) as above but with given base and power
- is_smooth(n,k) is n a k-smooth number
- is_rough(n,k) is n a k-rough number
- is_practical(n) is n a practical number
- is_delicate_prime(n[,b]) is n a digitally delicate prime (base b)
- is_sum_of_squares(n[,k]) can n be a sum of k squares (dflt k=2)
- is_congruent_number(n) is n a "congruent number"
- is_perfect_number(n) is n equal to sum of its proper divisors
- is_omega_prime(k,n) is n divisible by exactly k primes
- is_almost_prime(k,n) does n have exactly k prime factors
- is_chen_prime(n) is n prime and n+2 prime or semiprime
- is_powerfree(n[,k]) is n a k-powerfree number (default k=2)
- is_powerful(n[,k]) is n a k-powerful number (default k=2)
- is_perfect_power(n) is n a perfect power (1,4,8,9,16,25,..)
- is_odd(n) is n an odd integer (not divisible by 2)
- is_even(n) is n an even integer (divisible by 2)
- is_divisible(n,d) is n exactly divisible by d
- is_congruent(n,c,d) is n congruent to c mod d
- is_qr(a,n) is a a quadratic non-residue mod n
- is_cyclic(n) does n have exactly one group of order n
- forsquarefreeint {} lo,hi loop over square-free integers
- almost_primes(k,[start,]end) array ref of k-almost-primes
- almost_prime_count(k,n) count of integers with exactly k factors
- almost_prime_count_approx(k,n) fast approximate k-almost-prime count
- almost_prime_count_lower(k,n) fast k-almost-prime count lower bound
- almost_prime_count_upper(k,n) fast k-almost-prime count upper bound
- nth_almost_prime(k,n) The nth number with exactly k factors
- nth_almost_prime_approx(k,n) fast approximate nth k-almost-prime
- nth_almost_prime_lower(k,n) fast nth k-almost-prime lower bound
- nth_almost_prime_upper(k,n) fast nth k-almost-prime upper bound
- foralmostprimes {} k,a,b loop over k-almost-primes
- omega_primes(k,[start,]end) array ref of k-omega-primes
- omega_prime_count(k,n) count nums divisible by exactly k primes
- nth_omega_prime(k,n) The nth number dvsbl by exactly k primes
- powerful_numbers([lo,]hi[,k]) array ref of k-powerful from lo to hi
- powerful_count(n[,k]) count of k-powerful numbers <= n
- nth_powerful(n[,k]) the nth k-powerful number (default k=2)
- powerfree_count(n[,k]) count of k-powerfree numbers <= n
- nth_powerfree(n[,k]) The nth k-powerfree number (default k=2)
- powerfree_sum(n[,k]) sum of k-powerfree numbers <= n
- powerfree_part(n[,k]) remove excess powers from n
- powerfree_part_sum(n[,k]) sum of k-powerfree parts for 1 to n
- squarefree_kernel(n) integer radical of |n|
- perfect_power_count([lo,] hi) count of perfect powers
- smooth_count(n,k) count of k-smooth numbers <= n
- rough_count(n,k) count of k-rough numbers <= n
- allsqrtmod(a,n) all square roots of a (mod n)
- qnr(n) least quadratic non-residue
- subfactorial(n) count of derangements of n objects
- fubini(n) Fubini (Ordered Bell) number
- falling_factorial(x,n) falling factorial
- rising_factorial(x,n) rising factorial
- binomialmod(n,k,m) fast binomial(n,k) mod m
- negmod(a, n) -a mod n
- submod(a, b, n) a - b mod n
- muladdmod(a, b, c, n) a * b + c mod n
- mulsubmod(a, b, c, n) a * b - c mod n
- rootmod(a, k, n) modular k-th root
- allrootmod(a,k,n) all k-th roots of a (mod n)
- cornacchia(d,n) find x,y for x^2 + d * y^2 = n
- vecequal(\@a, \@b) compare two array refs for equality
- tozeckendorf(n) Zeckendorf binary string from n
- fromzeckendorf(str) n from Zeckendorf binary string
- lucasu(p,q,k) U(p,q)_k
- lucasv(p,q,k) V(p,q)_k
- lucasuv(p,q,k) U(p,q)_k, V(p,q)_k
- lucasumod(p,q,k,n) U(p,q)_k mod n
- lucasvmod(p,q,k,n) V(p,q)_k mod n
- lucasuvmod(p,q,k,n) (U(p,q)_k mod n, V(p,q)_k mod n)
- pisano_period(n) period of modular Fibonacci sequence
- prime_powers([start,] end) array ref of prime powers
- next_prime_power(n) next prime power: p > n
- prev_prime_power(n) previous prime power: p < n
- prime_power_count([start,] end) count of prime powers
- prime_power_count_approx(n) fast approximate prime power count
- prime_power_count_lower(n) fast prime power count lower bound
- prime_power_count_upper(n) fast prime power count upper bound
- nth_prime_power(n) the nth prime power
- nth_prime_power_approx(n) fast approximate nth prime power
- nth_prime_power_lower(n) fast nth prime power lower bound
- nth_prime_power_upper(n) fast nth prime power upper bound
- next_perfect_power(n) next perfect power: p > n
- prev_perfect_power(n) previous perfect power: p < n
- perfect_power_count([beg,] end) count of perfect powers
- perfect_power_count_approx(n) fast approximate perfect power count
- perfect_power_count_lower(n) fast perfect power count lower bound
- perfect_power_count_upper(n) fast perfect power count upper bound
- nth_perfect_power(n) the nth perfect power
- nth_perfect_power_approx(n) fast approximate nth perfect power
- nth_perfect_power_lower(n) fast nth perfect power lower bound
- nth_perfect_power_upper(n) fast nth perfect power upper bound
- next_chen_prime(n) next Chen prime: p > n
- lucky_numbers([start],] end) array ref of lucky numbers
- lucky_count([start,] end) count of lucky numbers
- lucky_count_approx(n) fast approximate lucky count
- lucky_count_lower(n) fast lucky count lower bound
- lucky_count_upper(n) fast lucky count upper bound
- nth_lucky(n) the nth lucky number
- nth_lucky_approx(n) fast approximate nth lucky number
- nth_lucky_lower(n) fast nth lucky number lower bound
- nth_lucky_upper(n) fast nth lucky number upper bound
- chinese2([a1,m1],[a2,m2],...) CRT returning (solution,modulus)
- frobenius_number(...) Frobenius number of a set
- vecmex(...) least non-negative value not in list
- vecpmex(...) least positive value not in list
- vecuniq(...) remove duplicates from list of integers
- vecfreq(...) return hash of item => count from list
- vecsingleton(...) remove all items that aren't unique
- vecslide { ... } @list calls block for each pair in list
- vecsort(@L) numerically sort integer list
- vecsorti(\@L) in-place numeric sort a list ref
- vecsample(k,@list) return k random elements of list
- toset(...) convert to int set (unique sorted aref)
- setinsert(\@A,$v) insert integer v into integer set A
- setinsert(\@A,\@B) insert list B into integer set A
- setremove(\@A,$v) remove integer v from integer set A
- setremove(\@A,\@B) remove list B values from integer set A
- setinvert(\@A,$v) if v is in set A, remove, otherwise add
- setinvert(\@A,\@B) invert for all values in integer list B
- setcontains(\@A,...) are list values all in int set A
- setcontains(\@A,\@B) is int set B a subset of int set A
- setcontainsany(\@A,...) are any list values in int set A
- setcontainsany(\@A,\@B) is any value in B in int set A
- setbinop {...} \@A[,\@B] apply op to cross product of set(s)
- sumset(\@A[,\@B]) new set from a+b {a:A,b:B}
- setunion(\@A,\@B) union of two integer lists
- setintersect(\@A,\@B) intersection of two integer lists
- setminus(\@A,\@B) difference of two integer lists
- setdelta(\@A,\@B) symmetric difference of two int lists
- is_sidon_set(\@L) is integer list L a Sidon set
- is_sumfree_set(\@L) is integer list L a sum-free set
- set_is_disjoint(\@A,\@B) is set B disjoint from set A
- set_is_equal(\@A,\@B) is set B equal to set A
- set_is_subset(\@A,\@B) is set B a subset of set A
- set_is_proper_subset(\@A,\@B) is set B a proper subset of set A
- set_is_superset(\@A,\@B) is set B a superet of set A
- set_is_proper_superset(\@A,\@B) is set B a proper superet of set A
- set_is_proper_intersection(\@A,\@B) "" a proper intersection ""
- contfrac(n,d) list of continued fraction for n/d
- from_contfrac(@A) return (p,q) rational from cfrac list
- next_calkin_wilf(n,d) next breadth-first CW rational
- next_stern_brocot(n,d) next breadth-first SB rational
- calkin_wilf_n(n,d) index of breadth-first CW rational
- stern_brocot_n(n,d) index of breadth-first SB rational
- nth_calkin_wilf(n) CW rational at breadth-first index n
- nth_stern_brocot(n) SB rational at breadth-first index n
- nth_stern_diatomic(n) Stern's Diatomic series; fusc(n)
- farey(n) list of Farey sequence order n
- farey(n,k) k'th entry of Farey sequence order n
- next_farey(n,[p,q]) next order-n rational after p/q
- farey_rank(n,[p,q]) number of F_n less than p/q
- minimal_goldbach_pair(n) least prime p where n-p is also prime
- goldbach_pair_count(n) count of how many prime pairs sum to n
- goldbach_pairs(n) array of all p where p and n-p are prime
- FIXES
- nth_ramanujan_* returns undef with input of 0.
- PP code for is_polygonal with n=0 was wrong. Fixed and refactored.
- Some large results that used GMP are properly objectified.
- When using GMP, native int results are no longer objectified.
- AKS was modding a with r, since 2012.
- is_totient(2^k) with k >= 33, found by Trizen in 2019.
- forsquarefree and forfactored were missing a destroy. From Trizen.
- todigitstring with non-standard bases and long strings. From Trizen.
- Fix v0.62 breaking extended precision long doubles on some machines. Big difference in epsilon-level precision for LambertW, Li, Zeta.
- sieve_range documentation for depth didn't exactly match the XS code.
- is_extra_strong_lucas_pseudoprime(5) was returning 0.
- lucas_sequence wasn't right for some outlier cases. This did not impact any primality tests or other internal code.
- gcd(-n) and lcd(-n) would return -n instead of n.
- All mod functions are more consistent. Like Pari, we use mod abs(n). If n=0, we return undef. If n=1, we return 0.
- invmod(0,1) = 0 instead of undef. Both cases make sense, but Pari, Mathematica, SAGE, and Math::BigInt all return 0.
- sqrtmod was not solving some cases, e.g. sqrtmod(4,8).
- chinese uses abs modulos, and values are all pre-modded. This fixes negative inputs, e.g. chinese([-4,17],[17,-19])
- The 70-rt-bignum test was horrendously slow with bignum 0.60+. It was bypassing the validation that sanitized the input.
- Fix logint for base > n. Github #75.
- lcm with empty list returns 1 instead of 0. Github #73.
- fix an XS stack issue with calling other routines inside forprimes
- is_pseudoprime and is_euler_pseudoprime have consistent XS, PP, and GMP behaviour. Specifically, single digit inputs and bases a multiple of n.
- kronecker in XS with 63-bit signed a and 64-bit unsigned n fixed.
- forfactored and forsquarefree were incorrect near 64-bit boundary.
- trial_factor etc. take an SV rather than silently coerce to a UV.
- PP trial_factor corrects for MPU::GMP's different API.
- PP is_perfect_square was wrong for 18446744073709551521 on 64-bit.
- FUNCTIONALITY AND PERFORMANCE
- LMO prime count is 1.2-1.5x faster.
- refactor objectification. Doesn't force. Enable on Perl < 5.9.
- save a nanosecond or two in 32-bit primality testing.
- Rewrite Ramanujan prime bounds using Axler (2017). For large values this is much improved. Big speedup for large exact nth/count.
- Approx Ramanujan prime nth and count are much more accurate.
- Internally, approx prime counts are faster, as we don't need to compute RiemannR to extreme precision.
- Faster and tighter nth_prime_upper bounds.
- Rewritten C version of Mertens based on pseudocode from Trizen. Over 100x faster at 2^32 and grows O(n^(2/3)) vs O(n).
- fromdigits for bigints uses subquadratic algorithm, thanks from Trizen. It also will call GMP if possible. The large example runs 100x faster.
- ExponentialIntegral with quadmath is much more precise across the range.
- rewrote C factorialmod. 1.1-2x faster for small inputs, 3-4x for large. The PP non-GMP version also has optimizations added.
- sum_primes faster for some inputs.
- Rewrote legendre_phi and removed code duplication. Much faster. Also makes _legendre_pi(n) about 4x faster.
- Rewrote prime count caching and removed code duplication. Uses less memory at large sizes and is faster.
- Error message "must be a positive integer" changed to "must be a non-negative integer" when zero is allowed.
- More accurate semiprime_count_approx and nth_semiprime_approx.
- Switched from Kahan to Neumaier sum in some real computations.
- is_primitive_root faster. Much faster for p^k and 2p^k with k >= 2.
- znprimroot faster. 2-20x faster for p^k and 2p^k with k >= 2.
- native integer input parsing is a little faster. A big deal for fast functions and long lists.
- rootint, is_power, {next,prev}_perfect_power up to 2x faster.
- inverse_totient memory leak fixed, memory use reduced, slightly faster.
- forsetproduct handles tied input arrays.
- stirling first kind without GMP is faster for some inputs.
- OTHER
- Legendre, Meissel, and Lehmer prime counting got rewritten to use our common code. 3x less code. No longer buildable standalone. Included by default, though still not called internally.
- lucas_sequence(n,p,q,k) deprecated. Use lucasuvmod(p,q,k,n) instead.
- Documentation more consistent with 'positive' vs 'non-negative'.
- valuation(n,k) now will error if k < 2. This follows Pari and SAGE. undef is returned for n = 0. Arguable "Inf" would be preferable.
- Added ability to select the bigint library, and rewrote large parts of the Perl backend to support other bigint classes (e.g. Math::GMPz and Math::GMP) rather than being specific to Math::BigInt.
- Older versions of Test::More with Math::GMPz and Math::GMP have issues comparing big integers. Lots of the tests will stringify the output so is(x,y,txt) will work even if x or y is an object.
Modules
Utilities related to prime numbers, including fast sieves and factoring
Pure Perl ChaCha20 CSPRNG
Elliptic curve operations for affine points
Elliptic curve operations for projective points
Get a good random seed
An auto-free object for Math::Prime::Util
Pure Perl version of Math::Prime::Util
PP front end for Math::Prime::Util
Primality proofs and certificates
A tied array for primes
An object iterator for primes
Generate random primes
Perl Big Float versions of Riemann Zeta and R functions
Number theory utilities
Examples
- examples/README
- examples/abundant.pl
- examples/csrand-gmp.pl
- examples/csrand.pl
- examples/fibprime-mce.pl
- examples/fibprime-serial.pl
- examples/fibprime-threads.pl
- examples/find_mr_bases.pl
- examples/findomegaseq.c
- examples/inverse_totient.pl
- examples/ktuplet-threads.pl
- examples/ktuplet.pl
- examples/numseqs.pl
- examples/porter.pl
- examples/project_euler_010.pl
- examples/project_euler_021.pl
- examples/project_euler_037.pl
- examples/project_euler_047.pl
- examples/project_euler_049.pl
- examples/project_euler_069.pl
- examples/project_euler_070.pl
- examples/project_euler_072.pl
- examples/project_euler_095.pl
- examples/project_euler_131.pl
- examples/project_euler_142.pl
- examples/project_euler_193.pl
- examples/project_euler_211.pl
- examples/project_euler_214.pl
- examples/project_euler_342.pl
- examples/project_euler_357.pl
- examples/sophie_germain.pl
- examples/twin_primes.pl
- examples/verify-cert.pl
- examples/verify-gmp-ecpp-cert.pl
- examples/verify-primegaps.pl
- examples/verify-sage-ecpp-cert.pl