Re: Fast calculation of pair correlation function
- To: mathgroup at smc.vnet.net
- Subject: [mg98645] Re: Fast calculation of pair correlation function
- From: Ray Koopman <koopman at sfu.ca>
- Date: Wed, 15 Apr 2009 04:59:33 -0400 (EDT)
- References: <gs1ntd$7ij$1@smc.vnet.net>
On Apr 14, 3:18 am, markus <markusg.p... at googlemail.com> wrote: > Hi, > > I have a quite large list (length ~ 10,000) of coordinates (positions > of particles from a simulation), e.g. data={{x1,y1,z1}, > {x2,y2,z2},...}, and I am looking for fast way to calculate the "pair > correlation" or "radial distribution function", which tells me, for > each value of r, how many particles are separated by a distance r (+- > epsilon). > Usually, this can be calculated by something like: > > Do[Do[histogr[[ Ceiling[Norm[ data[[i]]-data[[j]] ] / delta] ]], {j, i > +1, Length[data]}], {i, 1, Length[data]-1}] > > where histogr is just a list whose index corresponds to the distance > r. > > Unfortunately I have found that for lists of length >1,000, this way > of calculating is very slow (it can take ~minutes), compared to the > equivalent C code, which proceeds in a few seconds. I have also tried > to "Compile" the function, but the speed does not increase... > > Does anybody know some fast way of calculating the pair distribution > function in Mathematica? This will do 12 000 particles in 55 sec on an old 1.6 GHz Mac G5. The time scales roughly as n^2. Note that it wants the x, y, and z vectors passed separately. It returns a list of frequency counts for bins whose upper bounds are w, 2*w, ..., m@w, where m*w is at least twice the distance from the centroid to the point farthest from the centroid. func = Compile[{{x,_Real,1},{y,_Real,1},{z,_Real,1},{w,_Real}}, Module[{n,m,f}, n = Length@x; m = Ceiling[2./w * Sqrt@Max[(x-Tr@x/n)^2 + (y-Tr@y/n)^2 + (z-Tr@z/n)^2]]; f = Table[0,{m}]; Do[Scan[f[[#]]++&, Ceiling[Sqrt[ (Drop[x,i]-x[[i]])^2 + (Drop[y,i]-y[[i]])^2 + (Drop[z,i]-z[[i]])^2]/w]], {i,n-1}]; f]]