Vectorizing a summation

New to Julia
1 
term = 0
for k in 1:N
    term_inner = 0
    for p in 0:m
        term_inner += binomial(m, p) * c[k]^p *f[k]^(m - p) * integral_Ik_recursive(p, γ)
    end
    term += exp(im * 2π * γ * (k - 1)) * term_inner
end

I want to compute this summation if m and \gamma are vectors. Namely, ı want to fill the array with its dimensional length(m) x length(\gamma) . How can i write effective code?

2

One option is to put your code into a function and then broadcast it

First, put your code into a function:

function mysum(m, γ; c, f, N)
    term = 0
    for k in 1:N
        term_inner = 0
        for p in 0:m
            term_inner += binomial(m, p) * c[k]^p *f[k]^(m - p) * integral_Ik_recursive(p, γ)
        end
        term += exp(im * 2π * γ * (k - 1)) * term_inner
    end
    return term
end

The arguments to this function are scalars m and γ, but you can apply it to vectors with broadcasting:

mysum.(m, γ'; c, f, N)

Another option is to use broadcasting inside your mysum function so that it works on vectors. This takes a little more work, e.g. initializing term based on the sizes of the vectors, using “dot” operations like @. term += exp(...), and so on.

3

On a separate note, this is a very inefficient way to sum a series like this. You are recomputing the binomial factor and powers like c[k]^p separately for each loop iteration, when in practice you can compute them much more efficiently using a recurrence relation from the previous term.

See e.g. Sinc function - and sinc neural networks - and function approximation for the former - #4 by stevengj

1 Like
4

Thank you so much. But i know this method, actually I want to matrix product instead of for loop.

I separated each term like;

binomial.(m', p)

term_c = c .^ p'

term_f = reverse(f .^ (m)', dims=2)

Ip_recursive.(p, γ')

But I couldn’t arrange.

5

I have realized your last paragraph now. Sorry. Yes, ı want to solve by broadcasting inside.

6

It looks like you’re doing some discrete Fourier transforms as well. Have you seen the FFTW package?