Putting N hard spheres randomly in given volume

I need to put $N$ spheres with given radius $R$ randomly in a Volume $[-0.5,0.5]^3$, without any overlap of spheres.

If I choose values so that all the spheres will occupy ~57% of the total volume, I find it difficult to get results.
The basic scheme to place the $i$-th particle after $i-1$ particles have been placed:

However, the more particles have been placed, the less likely my random generator will find a position with enough space AND the more time the overlap check needs. For $N=500$, $R=0.5/7.7$, the algorithm is barely able to place 340 particles (overlap is called billions of times at this point).

C++ code:

using namespace std;

//this function places N hard spheres of radius R on random coordinates in a [-0.5,0.5]^3 box

vector <particle> random_positions_hard_spheres (int N, double R)
vector <particle> positions;

static seed_seq seed_sequence 100, 78639, 193, 55555, 3, 2348089, 6474367, 264441578 ;
static mt19937 gen (seed_sequence);
uniform_real_distribution<double> random_double(-0.5, 0.5);
double x,y,z;
double d=2*R;

for(size_t i=0; i<N; ++i)

x=random_double(gen);
y=random_double(gen);
z=random_double(gen);

positions.push_back(x,y,z);

//find position with no overlap
int j=0;
while(overlap_index(positions, i, d))

positions[i]=random_double(gen), random_double(gen), random_double(gen);
++j;

cout <<"overlap was called "+to_string(j)+" times"<<endl;
cout <<"particle "+to_string(i)+" was placed"<<endl;



return positions;



the overlap function checks whether particle $i$ overlaps with any previously placed particles.

1 Answer
1

Welcome to Computational Science StackExchange!

What you are trying to do is fairly difficult, and will need a more sophisticated method. That volume fraction of 0.57 lies just within the crystalline branch of the hard sphere phase diagram. The generation of dense random packings is an active area of current research.

In principle, a Monte Carlo simulation, or a collision-by-collision molecular dynamics program, can "equilibrate" an initial configuration prepared by randomly placing spheres which are significantly smaller than your desired size. Then, the idea is to steadily increase the sphere diameter, during the course of the simulation, until the desired volume fraction is reached.

A Monte Carlo code would not be too difficult to write, based on what you have already:

The alternative, of steadily increasing the diameter in the course of a molecular dynamics program, is one of the standard computational approaches to "random packing". Due to Lubachevsky and Stillinger, an example of its use can be found on arXiv (also published in Phys Rev Lett, 84 2064 (2000)).

The problem, of course, is that in principle, at that volume fraction, the system should be stable in the crystalline state, not in a randomly packed state. In practice, this may not be a problem: the system is likely to get "jammed" in a randomly packed state, at a volume fraction that depends on the compression rate, as discussed in the Phys Rev Lett.

MD programs for hard spheres can be found online, for example DynamO, or you can find guidance on how to write one in various books and papers. There are some standard tricks to improve the program speed, basically by making lists of near neighbours of the spheres: these would be worthwhile, but would not give huge speedups, for a system size of order 500.

Thanks for contributing an answer to Computational Science Stack Exchange!

But avoid …

Use MathJax to format equations. MathJax reference.

To learn more, see our tips on writing great answers.

Some of your past answers have not been well-received, and you're in danger of being blocked from answering.

Please pay close attention to the following guidance:

But avoid …

To learn more, see our tips on writing great answers.

Required, but never shown

Required, but never shown

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.