Path generation for non-intersecting disc movement on a plane

Question

What I'm looking for

I have 300 or fewer discs of equal radius on a plane. At time 0 each disc is at a position. At time 1 each disc is at a potentially different position. I'm looking to generate a 2D path for each disc for times between 0 and 1 such that the discs do not intersect and the paths are relatively efficient (short) and of low curvature if possible. (for example, straight lines are preferable to squiggly lines)

• Lower computation time is generally more important than exactness of solution. (for example, a little intersection is okay, and I don't necessarily need an optimal result)
• However, discs shouldn't teleport through each other, stop or slow abruptly, or change direction abruptly -- the "smoother" the better. Only exception is time 0 and 1.
• Paths can be expressed in a sampled form or piecewise linear nature (or better) -- I'm not worried about having truly smooth paths via splines. (I can approximate that if I so need.)

What I've tried

You can see a demo of my best attempt (via Javascript + WebGL). Be warned, it will load slowly on older computers due to the computations involved. It appears to work in Firefox/Chrome/IE11 under Windows.

In this demo I've represented each disc as an "elastic band" in 3D (that is, each disc has a position at each time) and ran a simple game-style physics engine that resolves constraints and treats each point in time like a mass with springs to the previous/next time. ('Time' in this case is just the third dimension.)

This actually works pretty well for small N (<20), but in common test cases (for example, start with discs arranged in circle, move each disc to the opposite point on the circle) this fails to generate convincing paths since the constraints and elasticity propagate slowly throughout the springs. (for example, if I slice time into 100 discrete levels, tension in the elastic bands only propagates one level per each simulation cycle) This makes good solutions require many (>10000) iterations, and that is tediously slow for my application. It also fails to reasonably resolve many N>40 cases, but this may be simply because I can't feasibly run enough iterations.

What else I've tried

My initial attempt was a hill-climber that started with straight-line paths which were gradually mutated. Solutions which measured better than the currently best solution replaced the currently best solution. Better measurements resulted from the amount of intersection (that is, completely overlapping measured worse than just grazing) and the length of the paths (shorter paths were better).

This produced some surprisingly good results, but unreliably, likely getting stuck in local minima very often. It was extremely slow for N>20. I tried applying a few techniques (simulated annealing, a genetic algorithms approach, etc) in an attempt to get around the local minima issue, but I never had much success.

What I'm trying

I'm optimizing the "elastic band" model so that tension and constraints propagate much more quickly in the time dimension. This would save a good deal of needed iterations in many cases, however in highly-constrained scenarios (for example, many discs trying to cross the same location) an untenable amount of iterations would still be required. I'm no expert on how to solve constraints or propagate springs more quickly (I've tried reading a few papers on non-stretchable cloth simulation, but I haven't been able to figure out if they apply), so I'd be interested in if there's a good way to go about this.

Ideas on the table

• Spektre has implemented a very fast RTS-style unit movement algorithm that works admirably well. It's fast and elegant, however it suffers from RTS-movement style problems: sudden direction changes, units can stop abruptly to resolve collisions. Additionally, units do not all arrive at their destination at the same time, which is essentially an abrupt stop. This may be a good heuristic to make viable non-smooth paths after which the paths could be resampled in time and a "smoothing" algorithm could be run (much like the one used in my demo.)
• Ashkan Kzme has suggested that the problem may be related to network flows. It would appear that the minimum cost flow problem could work, as long as space and time could be discritized in a reasonable manner, and the running times could be kept down. The advantage here is that it's a well studied set of problems, but sudden velocity changes would still be an issue and some sort of "smoothing" post-steps may be desirable. The stumbling block I'm currently having is deciding on a network representation of space-time that wouldn't result in discs teleporting through each other.
• Jay Kominek posted an answer that uses a nonlinear optimizer to optimize quadratic Bezier curves with some promising results.

Hypersensitive

Your question is not specific enough. Obviously you could pick a point and have all the circle rotate around that point. Or you could have all circles move outward in a straight line at a rate proportional to their distance from the center. Many other similar things could be invented. You need to more carefully define what kind of motion you are trying to achieve. – Fictitious

@TylerDurden: Moving in a semicircle is a good analytical solution for that specific problem, but I think the OP meant the crcle setup as an illustrative test case for a more general algorithm where the discs and their destinations may be arbitrary. I don't think that the question as is is ill-defined. – Neanderthal

"However, discs shouldn't ... move sporadically" - Is the velocity along the path constant or can it change? If so, you could probably get away with chosing simpler paths with a non-constant velocity profile. For example if you have to discs that cross each othzer when going straight, you could move one disc slowly first in the first half of the time period and faster thereafter and the other way round for the second disc, so that the paths cross, but the discs don't meet. – Neanderthal

I've changed all of the "e.g." to "for example" to help make it clearer for Tyler. – Hypersensitive

@M Oehm, a constant velocity would be preferable but isn't strictly necessary. Lower accelerations and jerks (in the physics sense) would be preferable. Above all else, it's more important that the paths are visually pleasing over time. Graceful, even. – Hypersensitive

Hmm, what kind of run-time are you hoping for? – Ogdoad

@Jay Kominek, preferably "real time" speeds -- computed paths could be visualized as beginning/ending configurations are modified by a user. To be a bit more lenient, if I could generate a solution in less than half a second, I could live with that. That said, I'll take what I can get -- it may not be possible to reasonable achieve any of this in half a second. – Hypersensitive

Added a web demo of my best attempt for reference. (See the "What I've tried" section.) – Hypersensitive

@Hypersensitive perhaps i'm not getting your question completely, but is there a way you can use the maximum separate paths algorithm, which is solved using Network Flow? if you find enough separate paths along a modeled graph, then you are good to go, ofcourse assuming that you have definite paths in your model. – Logan

@Ashkan Kzme, I'm having a little trouble finding an algorithm or topic by the name "maximum separate paths". Could you clarify what you're thinking of? – Hypersensitive

@Hypersensitive there is a polynomial time algorithm for finding the maximum number of separate paths from s to t in a graph, and the algorithm uses Network Flow. I think your case is something like the "Airline Scheduling Problem", read this : en.wikipedia.org/wiki/… , and this : en.wikipedia.org/wiki/Maximum_flow_problem#Airline_scheduling – Logan

@Ashkan Kzme, Network flows seem like an interesting approach for this problem. I can see the potential efficacy of the airline scheduling problem if we discretize time and space. I'm going to need to think about this some more when I have more time, but this seems like a solid tack. – Hypersensitive

@Hypersensitive please do and if you find it correct, tell me to post it as an answer. – Logan

Let us continue this discussion in chat. – Logan

RTS games have to solve a more complicated version of this problem, with obstacles and different sized units and formations. Also they need to do it in real time. You should look into how they do it. – Pendulous

Just thinking out loud, I'd try to view the problem as a static problem in 3D space, in which each disk becomes a tube connecting the start to the end. Then I'd apply collision detection to the tubes and attempt to resolve the collisions. – Aliunde

@Andy Newman, That's pretty much the spirit of my existing solution. This works pretty well for a low number of discs, albeit too slowly. For a large number of disks this can fail pretty spectacularly. – Hypersensitive

I don't have the time to look into it deeply, but the Spring RTS (opensource) has some pretty good handling of collision-avoiding, formation-preserving pathing; see videos at youtube.com/watch?v=6578hTKM1jE . Source is available here: github.com/spring/spring/tree/develop/rts/Sim/Path . Both a quad-tree approach and a heatmap/flowmap approach are implemented and available to Spring RTS -based games – Service

Answer 1

Have played with this for fun a bit and here the result:

Algorithm:

1. process each disc
2. set speed as constant*destination_vector
   • multiplicative constant a
   • and limit the speed to constant v afterwards
3. test if new iterated position does not conflict any other disc
4. if it does rotate the speed in one direction by some angle step ang
5. loop until free direction found or full circle covered

6. if no free direction found mark disc as stuck

   This is how it looks like for circle to inverse circle path:

   This is how it looks like for random to random path:

   stuck disc are yellow (none in these cases) and not moving discs are at destination already. This can also get stuck if there is no path like if disc already in destination circles another ones destination. To avoid that you need also change the colliding disc also ... You can play with the ang,a,v constants to make different appearance and also you could try random direction of angle rotation to avoid that swirling/twister movement

Here the source code I used (C++):

//---------------------------------------------------------------------------
const int    discs =23;     // number of discs
const double disc_r=5;      // disc radius
const double disc_dd=4.0*disc_r*disc_r;
struct _disc
    {
    double x,y,vx,vy;       // actual position
    double x1,y1;           // destination
    bool _stuck;            // is currently stuck?
    };
_disc disc[discs];          // discs array
//---------------------------------------------------------------------------
void disc_generate0(double x,double y,double r)     // circle position to inverse circle destination
    {
    int i;
    _disc *p;
    double a,da;
    for (p=disc,a=0,da=2.0*M_PI/double(discs),i=0;i<discs;a+=da,i++,p++)
        {
        p->x =x+(r*cos(a));
        p->y =y+(r*sin(a));
        p->x1=x-(r*cos(a));
        p->y1=y-(r*sin(a));
        p->vx=0.0;
        p->vy=0.0;
        p->_stuck=false;
        }
    }
//---------------------------------------------------------------------------
void disc_generate1(double x,double y,double r)     // random position to random destination
    {
    int i,j;
    _disc *p,*q;
    double a,da;
    Randomize();
    for (p=disc,a=0,da=2.0*M_PI/double(discs),i=0;i<discs;a+=da,i++,p++)
        {
        for (j=-1;j<0;)
            {
            p->x=x+(2.0*Random(r))-r;
            p->y=y+(2.0*Random(r))-r;
            for (q=disc,j=0;j<discs;j++,q++)
             if (i!=j)
              if (((q->x-p->x)*(q->x-p->x))+((q->y-p->y)*(q->y-p->y))<disc_dd)
               { j=-1; break; }
            }
        for (j=-1;j<0;)
            {
            p->x1=x+(2.0*Random(r))-r;
            p->y1=y+(2.0*Random(r))-r;
            for (q=disc,j=0;j<discs;j++,q++)
             if (i!=j)
              if (((q->x1-p->x1)*(q->x1-p->x1))+((q->y1-p->y1)*(q->y1-p->y1))<disc_dd)
               { j=-1; break; }
            }
        p->vx=0.0;
        p->vy=0.0;
        p->_stuck=false;
        }
    }
//---------------------------------------------------------------------------
void disc_iterate(double dt)                    // iterate positions
    {
    int i,j,k;
    _disc *p,*q;
    double v=25.0,a=10.0,x,y;
    const double ang=10.0*M_PI/180.0,ca=cos(ang),sa=sin(ang);
    const int n=double(2.0*M_PI/ang);
    for (p=disc,i=0;i<discs;i++,p++)
        {
        p->vx=a*(p->x1-p->x); if (p->vx>+v) p->vx=+v; if (p->vx<-v) p->vx=-v;
        p->vy=a*(p->y1-p->y); if (p->vy>+v) p->vy=+v; if (p->vy<-v) p->vy=-v;
        x=p->x; p->x+=(p->vx*dt);
        y=p->y; p->y+=(p->vy*dt);
        p->_stuck=false;
        for (k=0,q=disc,j=0;j<discs;j++,q++)
         if (i!=j)
          if (((q->x-p->x)*(q->x-p->x))+((q->y-p->y)*(q->y-p->y))<disc_dd)
            {
            k++; if (k>=n) { p->x=x; p->y=y; p->_stuck=true; break; }
            p->x=+(p->vx*ca)+(p->vy*sa); p->vx=p->x;
            p->y=-(p->vx*sa)+(p->vy*ca); p->vy=p->y;
            p->x=x+(p->vx*dt);
            p->y=y+(p->vy*dt);
            j=-1; q=disc-1;
            }
        }
    }
//---------------------------------------------------------------------------

Usage is simple:

1. call generate0/1 with center and radius of your plane where discs will be placed
2. call iterate (dt is time elapsed in seconds)
3. draw the scene

if you want to change this to use t=<0,1>

1. loop iterate until all disc at destination or timeout
2. remember any change in speed for each disc in a list need the position or speed vector and time it occur
3. after loop rescale the discs list all to the range of <0,1>
4. render/animate the rescaled lists

[Notes]

My test is running in real time but I did not apply the <0,1> range and have not too many discs. So you need to test if this is fast enough for your setup.

To speed up you can:

• enlarge the angle step
• test the collision after rotation against last collided disc and only when free test the rest...
• segmentate the disc into (overlapping by radius) regions handle each region separately
• also I think some field approach here could speed up things like create field map once in a while for better determine the obstacle avoidance direction

[edit1] some tweaks to avoid infinite oscillations around obstacle

For more discs some of them get stuck bouncing around already stopped disc. To avoid that just change the ang step direction once in a while this is the result:

you can see the oscillating bouncing before finish

this is the changed source:

void disc_iterate(double dt)                    // iterate positions
    {
    int i,j,k;
    static int cnt=0;
    _disc *p,*q;
    double v=25.0,a=10.0,x,y;
    const double ang=10.0*M_PI/180.0,ca=cos(ang),sa=sin(ang);
    const int n=double(2.0*M_PI/ang);
    // process discs
    for (p=disc,i=0;i<discs;i++,p++)
        {
        // compute and limit speed
        p->vx=a*(p->x1-p->x); if (p->vx>+v) p->vx=+v; if (p->vx<-v) p->vx=-v;
        p->vy=a*(p->y1-p->y); if (p->vy>+v) p->vy=+v; if (p->vy<-v) p->vy=-v;
        // stroe old and compute new position
        x=p->x; p->x+=(p->vx*dt);
        y=p->y; p->y+=(p->vy*dt);
        p->_stuck=false;
        // test if coliding
        for (k=0,q=disc,j=0;j<discs;j++,q++)
         if (i!=j)
          if (((q->x-p->x)*(q->x-p->x))+((q->y-p->y)*(q->y-p->y))<disc_dd)
            {
            k++; if (k>=n) { p->x=x; p->y=y; p->_stuck=true; break; }   // if full circle covered? stop
            if (int(cnt&128))   // change the rotation direction every 128 iterations
                {
                // rotate +ang
                p->x=+(p->vx*ca)+(p->vy*sa); p->vx=p->x;
                p->y=-(p->vx*sa)+(p->vy*ca); p->vy=p->y;
                }
            else{
                //rotate -ang
                p->x=+(p->vx*ca)-(p->vy*sa); p->vx=p->x;
                p->y=+(p->vx*sa)+(p->vy*ca); p->vy=p->y;
                }
            // update new position and test from the start again
            p->x=x+(p->vx*dt);
            p->y=y+(p->vy*dt);
            j=-1; q=disc-1;
            }
        }
    cnt++;
    }

Cyclostome

@Hypersensitive finished editing, PS the constants are set for view 256x256 pixels and disc radius = 5 pixels – Cyclostome

This has a lot going for it: actual code, fast, simple to implement. This follows Strilanc's suggestion -- look at how RTS games treat unit movement. For RTSes it's not necessarily undesirable for units to stop entirely, and arrival times vary greatly. These attributes don't map to this problem very well, however, since velocity and direction tends to change sporadically. That said, this algorithm could be run, and then the optimization steps of the algorithm in my demo could be applied to help smooth out the paths. I've edited the question to reflect this answer. – Hypersensitive

Answer 2

It isn't perfect, but my best idea has been to move the discs along quadratic Bezier curves. That means you've got just 2 free variables per disc that you're trying to find values for.

At that point, you can "plug" an error function into a nonlinear optimizer. Longer you're willing to wait, the better your solution will be, in terms of discs avoiding each other.

Only one actual hit:

Doesn't bother displaying hits, the discs actually start overlapped:

I've produced a full example, but the key is the error function to be minimized, which I reproduce here:

double errorf(unsigned n, const double *pts, double *grad,
              void *data)
{
  problem_t *setup = (problem_t *)data;
  double error = 0.0;

  for(int step=0; step<setup->steps; step++) {
    double t = (1.0+step) / (1.0+setup->steps);
    for(int i=0; i<setup->N; i++)
      quadbezier(&setup->starts[2*i],
                 &pts[2*i],
                 &setup->stops[2*i],
                 t,
                 &setup->scratch[2*i]);

    for(int i=0; i<setup->N; i++)
      for(int j=i+1; j<setup->N; j++) {
        double d = distance(&setup->scratch[2*i],
                            &setup->scratch[2*j]);
        d /= RADIUS;
        error += (1.0/d) * (1.0/d);
      }
  }

  return error / setup->steps;
}

Ignore n, grad and data. setup describes the specific problem being optimized, number of discs, and where they start and stop. quadbezier does the Bezier curve interpolation, placing its answer into ->scratch. We check ->steps points part way along the path, and measure how close the discs are to one another at each step. To make the optimization problem smoother, it doesn't have a hard switch when the discs start touching, it just tries to keep them all as far apart from one another as possible.

Completely compilable code, Makefile and some Python for turning a bunch of quadratic bezier curves into a series of images is available at https://github.com/jkominek/discs

Performance is a bit sluggish on huge numbers of points, but there are a number of options for improvement.

1. If the user is making minor tweaks to the starting and finishing positions, then after every tweak, rerun the optimization in the background, using the previous solution as the new starting point. Fixing up a close solution should be faster than recreating it from scratch every time.
2. Parallelize the n^2 loop over all points.
3. Check to see if other optimization algorithms will do better on this data. Right now it starts with a global optimization pass, and then does a local optimization pass. There are algorithms which already "know" how to do that sort of thing, and are probably smarter about it.
4. If you can figure out how to compute the gradient function for free or close to, I'm sure it would be worth it to do so, and switch to algorithms that can make use of the gradient information. It might be worth it even if the gradient isn't cheap.
5. Replace the whole steps thing with a suboptimization that finds the t at which the two discs are closest, and then uses that distance for the error. Figuring out the gradient for that suboptimization should be much easier.
6. Better data structures for the intermediate points, so you don't perform a bunch of unnecessary distance calculations for discs that are very far apart.
7. Probably more?

Ogdoad

looks interesting but your example is colliding/intersecting the discs – Cyclostome

giving it a fraction of a second, instead of a few milliseconds, and it finds a solution to the first one with only one hit. the second example has the discs starting overlapped, it's a demonstration that this handles completely degenerate input without weirdness or special cases. This guarantees smoothness without any postprocessing step (which would presumably also have to try and avoid introducing intersections). – Ogdoad

So far, this may be the most interesting approach for my problem so far presented since it addresses all of my preferences to some degree. I'm not sure I'd pick Bezier curves, but the general approach still applies. – Hypersensitive

Answer 3

The usual solution for this kind of problem is to use what is called a "heat map" (or "influence map"). For every point in the field, you compute a "heat" value. The disks move towards high values and away from cold values. Heat maps are good for your type of problem because they are very simple to program, yet can generate sophisticated, AI-like behavior.

For example, imagine just two disks. If your heat map rule is equi-radial, then the disks will just move towards each other, then back away, oscillating back and forth. If your rule randomizes intensity on different radials, then the behavior will be chaotic. You can also make the rule depend on velocity in which case disks will accelerate and decelerate as they move around.

Generally, speaking the heat map rule should make areas "hotter" at they approach some optimal distance from a disk. Places that are too near a disk, or too far away get "colder". By changing this optimal distance you can determine how close the disks congregate together.

Here are a couple of articles with example code showing how to use heat maps:

http://haufler.org/2012/05/26/beating-the-scribd-ai-challenge-implementing-traits-through-heuristics-part-1/

http://www.gamedev.net/page/resources/_/technical/artificial-intelligence/the-core-mechanics-of-influence-mapping-r2799

Game AI Pro, Volume 2, chapter on Heat Maps

Fictitious

Answer 4

I don't have enough rep to comment yet, so sorry for the non-answer. But to the RTS angle, RTS's generally use the A* algorithm for path finding. Is there a reason you're insisting on using a physics-based model?

Secondly, your attempt you linked that operates rather smoothly, but with the acceleration in the middle, behaves how I initially thought. Since your model treats it as a rubber band, it basically is looking for which way to rotate for the shortest path to the desired location.

If you arent worried about a physical approach, I would attempt as follows: Try to move directly toward the target. if it collides, it should attempt to roll clockwise around its most recent collision until it is in a position on the vector at 90 degrees to the vector from current location to the target location.

If we assume a test case of 5 in a row at the top of a box and five in a row at the bottom, they will move directly toward each other until they collide. The entire top row will slide to the right until they fall over the edge of the bottom row as it moves to the left and floats over the edge of the top row. (Think of what the whiskey and water shot glass trick looks like when it starts)

Since the motion is not determined by a potential energy stored in the spring which will accelerate the object during a rotation, you have complete control over how the speed changes during the simulation.

In a circular test like you have above, if all disks are initialized with the same speed, the entire clump will go to the middle, collide and twist as a unit for approximately a quarter turn at which point they will break away and head for their goal.

If the timing is lightly randomized, I think you'll get the behavior you're looking for.

I hope this helps.

Selfconsequence