Skip to content

Instantly share code, notes, and snippets.

View tkmharris's full-sized avatar
💭
508 Loop Detected

Tom Harris tkmharris

💭
508 Loop Detected
View GitHub Profile
@tkmharris
tkmharris / akiyama-tanigawa-bernoulli.hs
Last active May 26, 2021 20:38
Short Haskell implementation of the Akiyama-Tanigawa algorithm for computing Bernoulli numbers.
{-- Implementation of the Akiyama-Tanigawa algorithm
for computing Bernoulli numbers. --}
{-- Implement the algorithm by:
1. Begin with list of reciprocals of positive integers.
2. Write function to get next row and iteratie it to get list of lists.
3. Take the head of each list in the list to get Bernoulli nums.
Everything evaluates lazily, so we can do this nicely with infinite lists.
See tkmh.space/flotsam/haskell-akiyama-tanigawa for more.--}
@tkmharris
tkmharris / next_prime_bias.hs
Created October 17, 2018 18:45
Haskell function for looking at next prime equivalence class bias, as described in https://arxiv.org/abs/1603.03720
{--
We provide a function 'nextPrimeBias' that when called as
nextPrimeBias a b n r
returns the probability that, among the first r primes, a prime that is equivalent to a mod n is followed by a prime that is equivalent to b mod n.
We see that these probabilities are not evenly distributed among equivalence classes, as described in https://arxiv.org/abs/1603.03720 :
nextPrimeBias 7 1 10 100000 = 0.25736558
nextPrimeBias 7 3 10 100000 = 0.27695382
nextPrimeBias 7 7 10 100000 = 0.144993
nextPrimeBias 7 9 10 100000 = 0.3206876
@tkmharris
tkmharris / mordell.sage
Last active October 14, 2018 15:12
Get integer solutions of the Mordell equation y^2 = x^3 - d for certain special values of d.
## Mordell's equation
'''
Get integer solutions of the Mordell equation
y^2 = x^3 - d
for certain special values of d.
See mordell.txt for details.
'''
def mordell_test(d, x,y):
@tkmharris
tkmharris / hurwitz_radon.py
Last active October 13, 2018 21:13
Calculate the Hurwitz-Radon number of a positive integer
def hurwitz_radon(n):
'''Takes a positve integer n; returns the Hurwitz-Radon number p(n).
p(n) defined as follows:
Let n = (2^v)m with m odd and write v = 4a + b, 0<= b < 4.
Then p(n) = 8a + 2^b.'''
if not isinstance(n, int):
raise TypeError("hurwitz_radon only defined for positive integers")
if n <= 0:
raise ValueError("hurwitz_radon only defined for positive integers")
# helper to get largest power of 2 divinding n