For a recent project I had to generate a stimulus set consisting of random non-overlapping circles in MATLAB. My solution was to keep a list of the outside of a gradually increasing outer ring of circles, and add new ones, expanding the ring. This requirement is similar to circle-packing, except I needed to have gaps between the circles.

If you use especially large gaps, or the ratio between your maximum and minimum circle sizes is great, I detected the occasional overlap. But for the parameters I needed, I tested over 10000 tests of 400 circles, with no overlaps.

Here is how it works:

1 2 3 4 5 6 | % Parameters for the circles setup.min_circle_radius = 15; setup.max_circle_radius = 20; setup.min_gap = 10; setup.max_gap = 20; setup.n_circles = 400; |

First I set up some parameters – minimum and maximum circle radiuses and gaps, plus the number of circles I want to generate.

1 2 3 4 5 6 7 8 9 10 11 12 | function [ circles ] = CirclePacking( setup ) % Start with a circle in the middle first_circle = generate_random_radius_circle( setup ); first_circle.position = [ 0, 0 ]; % Add second circle at random angle second_circle = generate_random_radius_circle( setup ); gap = randi([setup.min_gap, setup.max_gap]); total_distance = first_circle.radius + gap + second_circle.radius; [x, y] = pol2cart(2 * pi * rand, total_distance); second_circle.position = [x, y]; |

Next I begin with two circles, one in the centre, and another at a random angle from the first.

1 2 3 | function [circle] = generate_random_radius_circle( setup ) circle.radius = randi([setup.min_circle_radius, setup.max_circle_radius]); end |

I represent circles as structs, first defined solely by radius, but then by position. All circles are first generated using this function which makes a random radius.

Back to the main function:

1 2 3 4 5 6 7 8 9 | first_circle.before = 2; first_circle.after = 2; second_circle.before = 1; second_circle.after = 1; % Setup doubly linked list first_circle.id = 1; circles(1) = first_circle; % have to set first circle manually circles = save_circle(second_circle, circles); |

I set up a doubly linked list, that is, each circle keeps track of the circle before and after it. With just two, before and after are the same.

1 2 3 4 5 6 | function [circles, circle] = save_circle(circle, circles) if (~isfield(circle, 'id')) circle.id = size(circles,2) + 1; end circles(circle.id) = circle; end |

‘save_circle’ adds or saves a circle to the main circles array. I couldn’t work out after a google search how to initialise an empty array of structs, so I didn’t waste time doing so.

Back to main function:

1 2 3 4 | current_position = 1; for i = 1:setup.n_circles [circles, current_position] = add_new_circle( setup, circles, current_position ); end |

To finish the main function, I just call ‘add_new_circle’, maintaining a variable ‘current_position’ which points to a particular circle in the circles array.

This is the end of the main function, most of what happens is going on inside ‘add_new_circle’, but before I go into it, I will show a few more important helper functions:

1 2 3 4 5 6 7 8 9 10 | function circles_overlap = circles_overlap(circle1, circle2, min_gap) dist = edist(circle1.position, circle2.position); circles_overlap = dist < circle1.radius + circle2.radius + min_gap; end function edist = edist(a, b) % euclidean distance edist = sqrt((a(1)-b(1))^2+(a(2)-b(2))^2); end |

The above functions are self explanatory, just checking whether circles overlap, and a function for euclidean distance.

1 2 3 4 5 6 7 | function next_circle = next_circle(circle, circles) next_circle = circles(circle.after); end function previous_circle = previous_circle(circle, circles) previous_circle = circles(circle.before); end |

The above functions allow me to get the next and previous circles – this is not necessary really, but I find it’s easier on the eyes using these functions.

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 | function new_circle = poisiton_new_circle(new_circle, circle1, circle2, gap1, gap2) distance1 = new_circle.radius + gap1 + circle1.radius; distance2 = new_circle.radius + gap2 + circle2.radius; distance_between = edist(circle1.position, circle2.position); % Calculate angle to new circle angle = get_law_cosines_angle_a(distance1, distance2, distance_between); % Calculate angle of first two circles vector_first_to_second = circle2.position - circle1.position; angle_first_to_second = atan2(vector_first_to_second(2), vector_first_to_second(1)); angle_of_new_from_first = mod(angle_first_to_second + angle, 2 * pi); [x, y] = pol2cart(angle_of_new_from_first, distance1); new_circle.position = [x, y] + circle1.position; end function angle = get_law_cosines_angle_a(ab,bc,ca) num = (ab * ab) + (ca * ca) - (bc * bc); den = 2 * ab * ca; angle = acos(num / den); end |

The ‘position_new_circle’ function puts a circle next to circle1 and circle2, taking into account the gaps passed in too. It always results in the anticlockwise order circle1, new_circle, circle2. It is using the law of cosines to do this – I won’t explain the maths here.

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 | function [circles, first_circle_index] = add_new_circle( setup, circles, first_circle_index ) % circles must contain 2 or more circles gap_with_first = randi([setup.min_gap, setup.max_gap]); gap_with_second = randi([setup.min_gap, setup.max_gap]); first_circle = circles(first_circle_index); second_circle = next_circle(first_circle, circles); new_circle = generate_random_radius_circle(setup); % Check backwards for overlaps for i = 1:8; new_circle = poisiton_new_circle(new_circle, first_circle, second_circle, gap_with_first, gap_with_second); before_circle = previous_circle(first_circle, circles); number_to_check = min(size(circles,2) ,8); overlaps = false; for j=1:number_to_check % Check if circle overlaps if (before_circle.id == second_circle.id) break; end if (circles_overlap(new_circle, before_circle, setup.min_gap)) first_circle = before_circle; first_circle_index = first_circle.id; overlaps = true; break; end before_circle = previous_circle(before_circle, circles); end if (~overlaps) break; end end % Update linked list circles = update_list(circles, first_circle, second_circle, new_circle); end |

This is the main part of the algorithm that works out where to put a new circle. It starts by setting the gaps, and taking the key circle, and the circle after that. ‘first_circle_index’ indexes the key circles, all circles are generated around this circle, until overlaps. Often, a new circle will overlap one of the previous circles in the chain. When this occurs I make the overlapping circle the new key circle, and generate circles around that.

I found that checking 8 circles back was more than enough (on the idea that you can fit roughly 6 circles around another circle – if they are similar size) – but this can be modified to check all circles, at a performance hit.

So after positioning a circle between the key circle and the next one, I check circles before the key (making sure I don’t check the same circle). If the circles overlap, I set the key circle to the one that overlaps.

If there are no overlaps, then I update the linked list with that circle.

1 2 3 4 5 6 7 8 9 10 11 | function circles = update_list(circles, circle1, circle2, new_circle) new_circle.before = circle1.id; new_circle.after = circle2.id; [circles, new_circle] = save_circle(new_circle, circles); circle1.after = new_circle.id; circles = save_circle(circle1, circles); circle2.before = new_circle.id; circles = save_circle(circle2, circles); end |

Note this can very obviously be massively sped up – but remembering that ‘premature optimization is the root of all evil’ – I have not yet needed to so.