CODEWITHSUNDEEP
algorithmsum
user19565054
user19565054

Reputation:

Algorithm to sum reciprocals

Are there any fast algorithms to evaluate a sum of the form (a * n + b) / (c * n + d) for a, b, c, d fixed, and n ranging from 1 to around 10^14 or so?

Obviously, summing each term individually won't work due to the size of the sum.

Edit: An algorithm to sum 1 / (c * n + d) would suffice, since you can split the fraction up and sum each numerator in O(1) time.

Upvotes: 0

Views: 113

Answers (2)

user20962667
user20962667

Reputation:

Here is a approximation function that is especially accurate for larger values of x as it was explained in here.

f = lambda x: 7/(6*x**15) - 691/(2730*x**13) + 5/(66*x**11) - 1/(30*x**9) + 1/(42*x**7) - 1/(30*x**5) + 1/(6*x**3) + 1/(2*x**2) + 1/x
# or 
f = lambda x: 0.97354060901812177/x**2 + 1/(x + 0.4849142940227510) # relative low accuracy and optimized for absolute error
eulers_mascheroni = 0.5772156649015328606065120900824024310421
print(f((eulers_mascheroni+(12*5)+7))) # obviously use an another language if you had speed bottlenecks

Altough you can just use a library in most languages that wil do that for you in the best implementation possible

In case of python you can use something like scipy or mpmath(which is also for floating point arithmetics)

from mpmath import *
mp.dps = 50; mp.pretty = True
print(psi(1,eulers_mascheroni+(12*5)+7))

If you were to use it check out there docs they have a great explanations of each function.

from scipy import special
from numpy import arange
x = [*map(lambda c,n,d: c * n + d + eulers_mascheroni,
    range(172),arange(24,42,.04),range(96))]
print(special.polygamma(1, x))

Upvotes: 1

user1196549
user1196549

Reputation:

You can reduce the summation to that of the inverses of n + α, which yields a shifted Harmonic number, corresponding to the Digamma function. The value is asymptotic to ln(n) + Γ (Euler's constant), with a correction for the missing/additional initial terms from 1 to floor(α).

For an approximation of the Digamma function, check https://en.wikipedia.org/wiki/Digamma_function#Computation_and_approximation

Upvotes: 1

Related Questions