gkeee
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?
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.
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
gkeee
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.
gkeee
5
I have realized your last paragraph now. Sorry. Yes, ı want to solve by broadcasting inside.
It looks like you’re doing some discrete Fourier transforms as well. Have you seen the FFTW package?