My planets are crazy - a problem with Java

This morning I wanted to write a useless program, and I ended up with this extremely minimal astronomical simulator in processing . The version of the MCVE code attached at the end of the message.

Then I added several planets, leaned back in my chair and prepared to watch as some planets revolve around. However, a problem arose.

The problem is that two initially static planets that are too close to each other tend to throw each other at ultra-high speed, and they will both disappear from the screen forever. This is clearly against the principle of momentum. Planets do not. I begin to doubt that something is wrong with my implementation of Newton's Law. But at some other times, planets work great, provided that they maintain a “safe” distance.

If you want to study this problem, I have carefully included two settings for you. The first one adds two planets and is perfectly stable. (You can do this by commenting on each line that has “Three.”) You can test the second, just enter the default code. This planet is going crazy.

This problem seems like a mystery that comes from some area of ​​computer science that I have never studied (although I am a noob). For the sake of my hair, I would really appreciate it if anyone could explain.

Start of code:

PVector planetOneLocation = new PVector(300, 200); PVector planetOneSpeed = new PVector(0, -.1); float planetOneMass = 1; PVector planetTwoLocation = new PVector(100, 200); PVector planetTwoSpeed = new PVector(0, .1); float planetTwoMass = 1; PVector planetThreeLocation = new PVector(200, 200); PVector planetThreeSpeed = new PVector(0, 0); float planetThreeMass = 10; float g = 5; void setup() { size(500, 500); } void draw() { updatePlanetOne(); updatePlanetTwo(); updatePlanetThree(); planetOneLocation.add(planetOneSpeed); planetTwoLocation.add(planetTwoSpeed); planetThreeLocation.add(planetThreeSpeed); background(0); ellipse(planetOneLocation.x, planetOneLocation.y, 10*planetOneMass, 10*planetOneMass); ellipse(planetTwoLocation.x, planetTwoLocation.y, 10*planetTwoMass, 10*planetTwoMass); ellipse(planetThreeLocation.x, planetThreeLocation.y, 10*planetThreeMass, 10*planetThreeMass); } void updatePlanetOne() { PVector accDir = PVector.sub(planetTwoLocation, planetOneLocation); float a = g * planetTwoMass / (accDir.mag() * accDir.mag()); accDir.normalize(); PVector acceleration = accDir.mult(a); planetOneSpeed.add(acceleration); accDir = PVector.sub(planetThreeLocation, planetOneLocation); a = g * planetThreeMass / (accDir.mag() * accDir.mag()); accDir.normalize(); acceleration = accDir.mult(a); planetOneSpeed.add(acceleration); } void updatePlanetTwo() { PVector accDir = PVector.sub(planetOneLocation, planetTwoLocation); float a = g * planetOneMass / (accDir.mag() * accDir.mag()); accDir.normalize(); PVector acceleration = accDir.mult(a); planetTwoSpeed.add(acceleration); accDir = PVector.sub(planetThreeLocation, planetTwoLocation); a = g * planetThreeMass / (accDir.mag() * accDir.mag()); accDir.normalize(); acceleration = accDir.mult(a); planetTwoSpeed.add(acceleration); } void updatePlanetThree() { PVector accDir = PVector.sub(planetOneLocation, planetThreeLocation); float a = g * planetOneMass / (accDir.mag() * accDir.mag()); accDir.normalize(); PVector acceleration = accDir.mult(a); planetThreeSpeed.add(acceleration); accDir = PVector.sub(planetTwoLocation, planetThreeLocation); a = g * planetTwoMass / (accDir.mag() * accDir.mag()); accDir.normalize(); acceleration = accDir.mult(a); planetThreeSpeed.add(acceleration); } 

Update : after some effort, I change the float variables to double. But my planets are shrinking from the screen. I think, besides the double / floating issue, there are actually some resolution issues. I did not determine any time step, and this also contributes to inaccurate behavior, especially when the speed is high.

Update 2 : setting the time step helps a lot. Some settings that went crazy without any time steps now work fine. But as long as there is a possibility that the center of the two planets will be very close, there will be chances that the system will begin to marvel again. To solve this problem, you need a better integrator.

Updata 3 . In response to @ kevin-workman, I put my beautiful code here to replace the project source code. Added the same third planet in the original post and updated the corresponding mathematics of Newton. The mv.add(p.speed.mult(p.mass)); his tests looks even if mv.add(p.speed.mult(p.mass)); commented out, the third planet is still going crazy (now that I have changed the demo code to the minimum version, there is no such line anymore). I think the error introduced in mult() contributes, but the specifically unstable integrator also plays an important role.

+7
processing physics
source share
2 answers

This problem has nothing to do with float vs double precision. float has more than enough precision for this, and actually Processing uses a float for everyone by default, so any double values ​​you are trying to use will simply be lost.

All your problems are caused by this line:

 for (Planet p: planets) { mv.add(p.speed.mult(p.mass)); } 

In particular, this bit:

 p.speed.mult(p.mass) 

This line multiplies each planet’s speed by its mass.

You might think that p.speed.mult(p.mass) will leave the original p.speed PVector unchanged, but it is not. PVector not immutable, so calling the mult() function changes the underlying instance of PVector .

Your first two planets have mass of 1 , so this line does not affect them. But your third planet has mass of 10 , which means that you multiply your speed by 10 for every single frame.

You can verify this by simply changing the mass of one of your original planets to 10 or even 2 , or by changing the mass of the third planet to 1 .

So, to fix your main problem, just get rid of this line or change it so that it does not change p.speed PVector .

Further information can be found in the processing link for PVector#mult() :

Multiplies a vector by a scalar. The version of the method that uses float acts directly on the vector on which it is called (as in the first example above). The versions that receive both the PVector and the float, since the arguments are static methods, and each returns a new PVector i.e. the result of a multiplication operation.

Basically, you can change this line as follows:

 PVector.mult(p.speed, p.mass) 

This will leave p.speed unchanged and return a copy with the multiplied value.

Taking a step back, you will have another problem, because your planets can be arbitrarily close to each other. In other words, their distance can approach (or equal to) zero. This does not happen in real life, and if so, you can bet that everything will go crazy. This way you will have crazy “gravity” if everything is too close. You can limit their distances.

If you have additional questions, it will be easier for you to help you if you are working with this MCVE and not publishing the whole project:

 PVector planetOneLocation = new PVector(300, 200); PVector planetOneSpeed = new PVector(0, -.1); float planetOneMass = 1; PVector planetTwoLocation = new PVector(100, 200); PVector planetTwoSpeed = new PVector(0, .1); float planetTwoMass = 10; float g = 5; void setup() { size(500, 500); } void draw() { updatePlanetOne(); updatePlanetTwo(); planetOneLocation.add(planetOneSpeed); planetTwoLocation.add(planetTwoSpeed); background(0); ellipse(planetOneLocation.x, planetOneLocation.y, 10*planetOneMass, 10*planetOneMass); ellipse(planetTwoLocation.x, planetTwoLocation.y, 10*planetTwoMass, 10*planetTwoMass); } void updatePlanetOne() { PVector accDir = PVector.sub(planetTwoLocation, planetOneLocation); float a = g * planetOneMass / (accDir.mag() * accDir.mag()); accDir.normalize(); PVector acceleration = accDir.mult(a); planetOneSpeed.add(acceleration); } void updatePlanetTwo() { PVector accDir = PVector.sub(planetOneLocation, planetTwoLocation); float a = g * planetTwoMass / (accDir.mag() * accDir.mag()); accDir.normalize(); PVector acceleration = accDir.mult(a); planetTwoSpeed.add(acceleration); } 
+10
source share

The problem is the integrating part. One of them is the Euler integrator, and it is not very accurate. Vertex integration is commonly used to model physics.

Quotation from Ilmari Karonen Answer , Verlet can be implemented as:

 acceleration = force(time, position) / mass; time += timestep; position += timestep * (velocity + timestep * acceleration / 2); newAcceleration = force(time, position) / mass; velocity += timestep * (acceleration + newAcceleration) / 2; 
0
source share

All Articles