Skip to content

Instantly share code, notes, and snippets.

@nathanfarlow
Created April 28, 2024 01:56
Show Gist options
  • Save nathanfarlow/fe1c399362160a76f0cc8ac9cc27e714 to your computer and use it in GitHub Desktop.
Save nathanfarlow/fe1c399362160a76f0cc8ac9cc27e714 to your computer and use it in GitHub Desktop.
def solve(f, inverse, period, target, bits):
"""Find some integer n such that f(n) ≈ target where f is periodic and
invertible"""
inverse = inverse(target.n(bits))
period = period.n(bits)
basis1 = [1 * 2**bits, 0, 0, inverse * 2**bits]
basis2 = [0, 1, 0, period * 2**bits]
basis3 = [0, 0, 1, -1 * 2**bits]
# Any integer combination of the basis vectors with coefficients a, b, and
# c will look like [a * 2^bits, b, c, 2^bits * (a * atan(target) + b * pi -
# c)]
# LLL gives us short basis vectors, so we'll weigh the entries depending on
# how important it is that the entries are small. For example, we weigh the
# last term (a * atan(target) + b * pi - c) extremely high since it is
# critical that this equation = 0 for the solution to make sense in the
# first place. Similarly, we weigh the `a` term equally high since in our
# original formulation, `a` needs to be 1 since `a` is the coefficient on
# the atan term. Then we can take the best c term of the LLL-reduced basis
# vectors to be our solution. Somehow our solution magically appears.
L = matrix(QQ, [basis1, basis2, basis3])
L = L.LLL()
dist = lambda n: abs(f(n) - target).n(bits)
c = min(L[:, 2].list(), key=dist)
n = round(c)
print(f"n = {n}")
print(f"f(n) = {f(n).n(bits)}")
print(f"error = {dist(n).n(bits)}")
print(f"number of correct digits: {int(-log(dist(n), 10).n())}")
if __name__ == "__main__":
print(f"Solving for n such that tan(n) ≈ e")
solve(f=tan, inverse=atan, period=pi, target=e, bits=256)
print()
print(f"Solving for n such that sin(n) ≈ e - 2")
solve(f=sin, inverse=asin, period=2 * pi, target=e - 2, bits=256)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment