How do I vectorize adding array elements based on indices?

% I'm building a form of spiking neural network. Long story, but the reason
%  for the inquiry is rather brief, and (I thought) simple.

% Background, if it helps. The network has neurons and synapses. Each synapse
%  has a train (capturing all of the spikes) and destination (the index of the
%  neuron to which the spike is added). Each neuron has a potential, to which
%  input spikes are added (basically, 1s or 0s). 
% So, we enter an optimization problem...

% When the program is "fully operational," the element sizes of the potential
%  and train/destination can be upwards of 100,000 and 600,000,000, 
% respectively! I would like to use vectorization, maybe even GPU Arrays, 
% but I can't seem to figure out how to efficiently vectorize without the 
% results being wrong!

% take the following arrays:
>> potential = [ 1 , 2 , 3 , 4 , 5 ];
>> train      = [ 8 , 7 , 8 , 7 , 8 , 7 , 8 , 7 ];
>> destination= [ 1 , 1 , 1 , 2 , 2 , 2 , 3 , 4 ];
% run the following for-loop
>> for i = 1 : 8
potential( destination(i) ) = potential( destination(i) ) + mod( train(i) , 2 );
end
% arrive at a solution
>> potential
     3    3    4    4    5

% try to vectorize: the data is wrong
>> potential(destination) = potential(destination) + mod(train,2)
     1    3    3    5    5

% try to use "arrayfun": not even the right size, anymore!
>> potential = arrayfun( @(i) potential(destination(i))+mod(train(i),2), destination)
     1     1     1     2     2     2     1     3

[Stephen23]

Your example output values do not match what your code actually produces:

>> potential   = [1,2,3,4,5];
>> train       = [8,7,8,7,8,7,8,7];
>> destination = [1,1,1,2,2,2,3,4];
>> for k = 1:8, potential(destination(k)) = potential(destination(k)) + mod(train(k),2); end
>> potential
potential =
     2     4     3     5     5

Answer 1

[Stephen23]

Using accumarray:

>> potential = [1,2,3,4,5];
>> train     = [8,7,8,7,8,7,8,7];
>> destination = [1,1,1,2,2,2,3,4];
>> potential(:) + accumarray(destination(:),mod(train(:),2),[numel(potential),1])
ans =
     2
     4
     3
     5
     5