jimbelk

Assuming there's no friction on the table, the chances ought to be pretty high that they'll all eventually go in.

Hi! This seems like a neat website -- kind of like early reddit.

This is a good question, and I sometimes do put too much time into answering questions on Math Stack exchange. Right now, though, I'm fairly happy with the time I put into this simulation.

First of all, there's the fact that I did this on a Friday evening, after I came home from work. (I spent about six hours on it, and maybe another hour or two on Saturday morning). There are probably more productive uses of my time, but it was fun, and the most likely alternative was probably something like "read reddit" as opposed to "do something productive". It's nice to have leisure-time activities that are at least marginally related to my career.

Speaking of which, here's a few reasons that this was helpful for me:

1. I got a little better at programming in Mathematica. I use Mathematica a lot for both my teaching and research, so it always helps to practice.

2. I learned about a whole field of physics called contact mechanics, which I otherwise wouldn't have been aware of. It's part of my job as a math professor to know about random things like this.

3. I got some experience with modeling in applied mechanics problems. This may be helpful for me in the future -- as a professor at Bard, I regularly have to supervise senior math majors who want to work on research-level problems. The next time a senior who's interested in applied mechanics wants a research problem, I'll have one available.

That's how it's helpful for me, but I suppose it's less clear that this was helpful for the rest of the world. I guess, at the very least, it certainly seems to have entertained a bunch of people, and maybe someone learned something about modeling or differential equations or physics.

Hi!

Sorry about that. I really ought to upgrade to 9.

Apparently Mathematica added some new features in version 9 that are intended to better handle differential equations with discontinuities (see http://reference.wolfram.com/mathematica/tutorial/NDSolveWhe...). Apparently this is getting tripped up at the boundary between the two behaviors.

One simple thing that might work is to change the force law to either

  force[disp1_, disp2_, dir_] := -k Abs[(disp1 - disp2).dir ]^(3/2) dir* If[(disp1 - disp2).dir > 0, 1, 0]

or

  force[disp1_, disp2_, dir_] := -k Abs[(disp1 - disp2).dir ]^(3/2) dir* Piecewise[{{1, (disp1 - disp2).dir > 0}}, 0]

If these don't fix it, then you would need to play around with the options for NDSolve.