Sorry for the confusion. I was trying to get a simple example to work so I
wouldn't get distracted by details. The @everywhere did the trick.
This is the fastest parallel version of the code that I was able to get
working. However, I easily run into memory limitations (8GB RAM) as I
increase the variable 'limit'.
function expensive_hat(S::Array{Complex{Float64},2},limit::Int64)
samples = size(S,2);
Ly = int64(limit/4);
Lx = int64(limit/2);
@everywhere add_three_tuple(x,y) = (x[1]+y[1], x[2]+y[2], x[3]+y[3])
tuphat = @parallel (add_three_tuple) for k = 1:samples
A = B = zeros(Float64,Ly,Lx);
C = complex(zeros(Float64,Ly,Lx));
for j = 1:Lx, i = 1:Ly
Sx = S[i,k];
Sy = S[j,k];
Sxy = S[i+j,k];
Sxyc = conj(Sxy);
A[i,j] += abs2(Sy*Sx);
B[i,j] += abs2(sqrt(Sxyc*Sxy));
C[i,j] += Sxyc*Sy*Sx;
end
(A, B, C)
end
D = tuphat[1].*tuphat[2]/samples/samples;
E = tuphat[3]/samples;
return (D, E);
end
srand(77);
data = rand(24000,64);
limit = 4000;
S = rfft(data,1)./24000;
@time (D, E) = expensive_hat(S,limit);
On Tuesday, February 4, 2014 3:01:17 PM UTC-5, David Salamon wrote:
>
> huh.
>
> maybe @everywhere in front of the function definition? I'm not sure
>
>
> On Tue, Feb 4, 2014 at 10:53 AM, Alex C <[email protected]
> <javascript:>>wrote:
>
>> Thanks for the hint. Getting rid of 'mx' and 'my' definitely helps.
>>
>> I couldn't figure out how to implement the parallel version of tuple
>> adding. This is what I've got. It crashes my Julia Studio console when I
>> try to run it. What am I missing?
>>
>> add_two_tuple(x,y) = (x[1]+y[1], x[2]+y[2], x[3]+y[3])
>>
>>
>> @parallel (add_two_tuple) for i = 1:100
>>
>> x = (1,2,3)
>>
>> end
>>
>>
>> I also don't have the SharedArray function (running Julia 0.2.0) so I
>> couldn't implement the other alternative.
>>
>> Alex
>>
>>
>> On Monday, February 3, 2014 11:45:57 AM UTC-5, David Salamon wrote:
>>
>>> You're not out of the no-slicing woods yet. Looks like you can get rid
>>> of `mx` and `my`
>>>
>>> for i=1:limit, j=1:int64(limit/2)
>>> end
>>>
>>>
>>>
>>> As far as parallelizing, you could define:
>>> three_tup_add(a, b, c) = (a[1] + b[1] + c[1], a[2] + b[2] + c[2], a[3] +
>>> b[3] + c[3])
>>>
>>> and then do a @parallel (three_tup_add) over your sample index?
>>>
>>> for that matter, why not compute the two parts of the answer directly
>>> rather than going via A, B, and C?
>>>
>>>
>>>
>>>
>>>
>>>
>>> On Mon, Feb 3, 2014 at 8:11 AM, Alex C <[email protected]> wrote:
>>>
>>>> Thanks. I've re-written the function to minimize the amount of copying
>>>> (i.e. slicing) that is required. But now, I'm befuddled as to how to
>>>> parallelize this function using Julia. Any suggestions?
>>>>
>>>> Alex
>>>>
>>>> function expensive_hat(S::Array{Complex{Float64},2},
>>>> mx::Array{Int64,2}, my::Array{Int64,2})
>>>>
>>>> samples = 64
>>>> A = zeros(size(mx));
>>>> B = zeros(size(mx));
>>>> C = zeros(size(mx));
>>>>
>>>> for i = 1:samples
>>>> Si = S[:,i];
>>>> Sx = Si[mx];
>>>> Sy = Si[my];
>>>> Sxy = Si[mx+my];
>>>> Sxyc = conj(Sxy);
>>>>
>>>> A += abs2(Sy .* Sx);
>>>> B += abs2(sqrt(Sxyc .* Sxy));
>>>> C += Sxyc .* Sy .* Sx;
>>>> end
>>>>
>>>> ans = (A .* B ./ samples ./ samples, C./samples)
>>>> return ans
>>>>
>>>> end
>>>>
>>>> data = rand(24000,64);
>>>> limit = 2000;
>>>>
>>>> ix = int64([1:limit/2]);
>>>> iy = ix[1:end/2];
>>>> mg = zeros(Int64,length(iy),length(ix));
>>>> mx = broadcast(+,ix',mg);
>>>> my = broadcast(+,iy,mg);
>>>> S = rfft(data,1)./24000;
>>>>
>>>> @time (AB, C) = expensive_hat(S,mx,my);
>>>>
>>>>
>>>
>>>
>
