Last active
September 8, 2019 14:15
-
-
Save badgateway666/36efba86744926de35b284c29cf50a58 to your computer and use it in GitHub Desktop.
Sum fun. (BADUMTSS) v1
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
import random | |
from math import sqrt | |
import numpy as np | |
def harmonic_series(stop_after=-1): | |
""" | |
:param stop_after: | |
:return: | |
""" | |
sum = np.longdouble(0.0) | |
k = 1 | |
while True: | |
sum += 1/k | |
yield sum | |
if k == stop_after: | |
break | |
k += 1 | |
def alternating_harmonic_series(stop_after=-1): | |
""" | |
:param stop_after: | |
:return: | |
""" | |
sum = np.longdouble(0.0) | |
k = 1 | |
negative = False | |
while True: | |
if negative: | |
sum += (-1)/k | |
else: | |
sum += 1/k | |
yield sum | |
if k == stop_after: | |
break | |
k += 1 | |
negative = not negative | |
def general_harmonic_series(a, b, alternating=False, stop_after=-1): | |
""" | |
:param a: | |
:param b: | |
:param alternating: | |
:param stop_after: | |
:return: | |
""" | |
assert a != 0 | |
assert b/a > 0.0 | |
sum = np.longdouble(0.0) | |
n = 0 | |
negative = False | |
while True: | |
if negative: | |
sum += (-1) / ((a*n) + b) | |
else: | |
sum += 1 / ((a*n) + b) | |
yield sum | |
if n == stop_after: | |
break | |
if alternating: | |
negative = not negative | |
n += 1 | |
def leibniz_series(stop_after=-1): | |
""" | |
Approaches pi/4 | |
:param stop_after: if specified, stop after n-th element of series | |
:return: leibniz series generator object | |
""" | |
sum = np.longdouble(0.0) | |
k = 0 | |
negative = False | |
while True: | |
if negative: | |
sum += (-1) / ((2*k)+1) | |
else: | |
sum += 1 / ((2*k)+1) | |
yield sum | |
if k == stop_after: | |
break | |
k += 1 | |
negative = not negative | |
def p_series(p, stop_after=-1): | |
""" | |
:param stop_after: | |
:return: | |
""" | |
sum = np.longdouble(0.0) | |
k = 1 | |
while True: | |
sum += 1 / (k ** p) | |
yield sum | |
if k == stop_after: | |
break | |
k += 1 | |
def random_harmonic_series(stop_after=-1): | |
""" | |
:param stop_after: | |
:return: | |
""" | |
sum = np.longdouble(0.0) | |
k = 1 | |
while True: | |
sum += random.choice([-1, 1]) / k | |
yield sum | |
if k == stop_after: | |
break | |
k += 1 | |
def euler_accelerator(series): | |
""" | |
Lets a series converge faster | |
:param series: | |
:return: | |
""" | |
s0 = series.__next__() # Sn-1 | |
s1 = series.__next__() # Sn | |
s2 = series.__next__() # Sn+1 | |
while True: | |
yield s2 - (sqrt(s2 - s1)) / (s0 - 2*s1 + s2) | |
s0, s1, s2 = s1, s2, series.__next__() | |
def aitken_accelerator(series): | |
""" | |
Lets a series converge even faster | |
:param series: | |
:return: | |
""" | |
s0 = series.__next__() # Sn-1 | |
s1 = series.__next__() # Sn | |
s2 = series.__next__() # Sn+1 | |
while True: | |
yield s2 - ( ((s2 - s1) ** 2) / ( s2 - (2*s1) + s0) ) | |
s0, s1, s2 = s1, s2, series.__next__() | |
def multiplier(series,factor): | |
""" | |
Multiplies every element of a generator by a factor | |
:param series: Generator object from where to obtain element | |
:param factor: factor to multiply each element with | |
:return: multiplied element generator object | |
""" | |
for elem in series: | |
yield factor * elem | |
def converger(series, precision_exponent=25): | |
""" | |
:param series: | |
:param precision_exponent: | |
:return: | |
""" | |
precision = 10**(-precision_exponent) | |
print(f"Converging {series} to a delta precision of {precision}") | |
prev = np.inf | |
for idx, val in enumerate(series): | |
delta = abs(val - prev) | |
if delta < precision: | |
print(f"Reached target convergence precision at Step #{idx}") | |
return idx, val | |
prev = val | |
if __name__ == '__main__': | |
to_converge = [ aitken_accelerator(leibniz_series())] | |
for result in map(converger, to_converge): | |
print(f"After #{result[0]} Steps, converged to a value of {result[1]*4}/4") | |
print(f"Difference to pi: {(result[1]*4) - np.pi}") |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment