Geek Answers Handbook

Infinite Summation in Python

I have a function, since I have to do infinite summation over (over all integers) numerically. Summation does not always have to converge, since I can change the internal parameters. The function looks like this:

m(g, x, q0) = sum(abs(g(x - n*q0))^2 for n in Integers) m(g, q0) = minimize(m(g, x, q0) for x in [0, q0])

using pseudo code Pythonic

Using the Scipy integration methods, I just flooring n and integrated, as with fixed x,

 m(g, z, q0) = integrate.quad(lambda n: abs(g(x - int(n)*q0))**2, -inf, +inf)[0]

This works very well, but then I need to do the optimization on x as a function of x, and then do another summation on what gives the integral of the optimization of the integral. To a large extent this takes a very long time.

Do you know that it is better to do the summation faster? Manual coding seemed slow to him.

I am currently working with

 g(x) = (2/sqrt(3))*pi**(-0.25)*(1 - x**2)*exp(-x**2/2)

but the solution must be general

Document from Wavelet Conversion, Time Localization, and Signal Analysis by Daubechies (IEEE 1990)

thanks

+5
source share
3 answers

Thanks to all the helpful commentary, I wrote my own adder, which seems to work pretty fast. Someone has recommendations to make it better, I will gladly take them.

I will test this for the problem I'm working on, and as soon as it shows success, I will claim that it is functional.

 def integers(blk_size=100): x = arange(0, blk_size) while True: yield x yield -x -1 x += blk_size # # For convergent summation # on not necessarily finite sequences # processes in blocks which can be any size # shape that the function can handle # def converge_sum(f, x_strm, eps=1e-5, axis=0): total = sum(f(x_strm.next()), axis=axis) for x_blk in x_strm: diff = sum(f(x_blk), axis=axis) if abs(linalg.norm(diff)) <= eps: # Converged return total + diff else: total += diff
+4
source

g(x) almost certainly your bottleneck. A very quick and dirty solution would be to vectorize it to work with an array of integers, and then use np.trapz to evaluate the integral using the trapezoid rule:

 import numpy as np # appropriate range and step size depends on how accurate you need to be and how # quickly the sum converges xmin = -1000000 xmax = 1000000 dx = 1 x = np.arange(xmin, xmax + dx, dx) gx = (2 / np.sqrt(3)) * np.pi**(-0.25)*(1 - x**2) * np.exp(-x**2 / 2) sum_gx = np.trapz(gx, x, dx)

Alternatively, you can overwrite g(x) with Cython or numba to speed it up.

+2
source

There Numba chance significantly improves speed - http://numba.pydata.org

It is a little painful to install, but very easy to use. Take a look at: https://jakevdp.imtqy.com/blog/2015/02/24/optimizing-python-with-numpy-and-numba/

+1
source

More articles:


All Articles