Some gravity particles

[User137]

Problem is in that each thread will access to same particle data at the same time. I cannot split them to isolated index groups, because each particle is compared to every other. So many threads could apply forces to same particle at the same time. But it's just a simple sum operation, the calculation does not depend on value of the movement vector, so it might not matter. I do not know. I could do the "// Move particles" part in main thread after the forces are updated, cause i think that's the fastest part of it. I would also put main thread waiting for those to finish before continuing. It will go to same 1-4 fps while counting it anyway, it's hardly controllable at all.

Extra treat, shaders:

#version 120
attribute vec3 in_position;
attribute vec4 in_color;
varying vec4 color;
uniform mat4 pmv;
void main() {
  color = in_color;
  gl_Position = pmv * vec4(in_position, 1.0);
}

#version 120
varying vec4 color;
void main() {
  gl_FragColor = color;
}

edit: Using 4 physics threads now on the time-critical part, and CPU use went up to 90% or so. I don't at least see any differences in the visual outcome, sure it counts it faster.

[Dan]

you should not approach this task with a brute force of the CPU, because you will surely get stuck.

I can make a couple of suggestions to improve the algorithm:

1 - use GPU to do your calculations.
2 - cluster your space (each cluster can contain hundreds of particles). calculate the total mass and the center of gravity of each cluster. attract individual particles with each other within a cluster (or you might create sub-clusters if you have too many particles). attract individual particles within a cluster to the center of mass of the other clusters. this will significantly decrease the time complexity of your calculations (at the moment you have O(N) = N * N).

[User137]

Can i do clustering, if start situation is static mass that is mix of both gravity types? I think at beginning there would then be as many clusters as there are particles. I'm also not very convinced yet, because cluster's gravity center location is a living thing that depends on each of its children. I guess it would be more of a heuristic simulation.

Another funny sidenote about this simulator, is that currently the radius of the particle mass is 10. But if i reduce it down to 0.001 or so - all particles in very tiny spot next to eachother, they literally explode in a "Big Bang". (But i've seen that before, didn't come as a surprise...)

I think the clustering might work, but with different start situation. Still it would be very rough estimation. Say 50-100 clusters with random mass (mass being amount of assigned particles). Then each particle only orbiting the 1 cluster, and each cluster orbiting all the other clusters. Problem becomes of deciding when 1 particle switches to orbit another cluster, and i'm not sure i could solve it. Firstly the borderline between 2 clusters would be calculated with Newton's equations, spot between them where gravity force is equal. Wouldn't i be doing that calculation for each particle per each cluster all the time? Or each cluster would count and store its nearest neighbours every once in a while.

Anyway that's too much work for simple test.

[Dan]

there might be one problem with using clusters, and that is the diminishing gravity over distance. I am not sure if it is possible to simulate this effect. I might give it a try myself.

[SilverWarior]

Problem is in that each thread will access to same particle data at the same time.

Unles your threads are only reading the data and not changing it there should be no problem. Problem would only ocur if one of the threads would be changing some values of one particle while other tread would be reading it (data corruption).

I could do the "// Move particles" part in main thread after the forces are updated, cause i think that's the fastest part of it. I would also put main thread waiting for those to finish before continuing.

First of all never and I mean NEVER put the min thread on hold for other threads to finish. Why? Doing so you prevent your application to process Windows messages or System messages if you are using some other platform. And if you do so then your application would be considered as being hanged application by the OS and your OS might simply go and forcefully kill it. You definitly don't want this.
If you use double-buffering approach you actually dont even need main thread to do any moving.

Let me try to explain you how would I approach to this problem:
First I would save information for each particle in class object with double-buffering so I that I have acces to the variables using properties. Why using properties?
Becouse with getter and setter methods of properties you can easily implement double-buffering.
The class code for particle wouls look something like this:

type
  TParticle = class(TObject)
  private
    PSecondBuffer: ^Boolean;
    FPosition1: T3DPoint;
    FPosition2: T3DPoint;
    FForce1: T3DVector;
    FForce2: T3DVector;
  protected
    function GetPosition: T3DPosition;
    procedure SetPosition(AValue: T3DPoint);
    function GetForce: T3DVector;
    procedure SetForce (AValue: T3DVector);
  public
    constructor Create(AOwner: TObject; var SecondBuffer: Boolean); override;
    property Position: T3DPoint read GetPosition write SetPosition;
    property Force: T3DVector read GetForce write SetForce;
  end;

//We use this global variable to determine whether we will read from primary or secondary buffer
var SecondBuffer: Boolean; 

implementation


{ TParticle }
constructor TParticle.Create(AOwner: TObject; var SecondBuffer: Boolean);
begin
    inherited;
    //We assing pointer to external boalean variable for controling whether
    //we read from primary or secondary buffer
    FSecondBuffer := @SecondBuffer;
end;

function TParticle.GetForce: T3DVector;
begin
    //If the second buffer is true then we read from second force variable
    //else we read from first force variable
    if FSecondBuffer = True then Result := FForce2 //Secondarry buffer
    else Result := FForce1; //Primary buffer
end;

function TParticle.GetPosition: T3DPosition;
begin
    //If the second buffer is true then we read from second postiion variable
    //else we read from first position variable
    if FSecondBuffer = True then Result := FPosition2 //Secondarry buffer
    else Result := FPosition1; //Primary buffer
end;

procedure TParticle.SetForce(AValue: T3DVector);
begin
    //With the difference of getter methods we now if second buffer is true
    //write data to first force variable to avoid changing the data that might
    //be read by other threads
    //We never change the same variable as it is available for reading
    if FSecondBuffer = True then FForce1 := AValue //Primary buffer
    esle FForce2 := AValue; //Secondarry buffer
end;

procedure TParticle.SetPosition(AValue: T3DPoint);
begin
    //We handle position variable changing same as above.
    if FSecondBuffer = ture then FPosition1 := AValue //Primary buffer
    else FPosition2 := AValue; //Secondarry buffer
end;

As you see each class has two internal variables for every property (value) that might be changed by other threads during curent cycle while still preserving acces to current data without it being corupted by other threads.
For contgrolling this we are actually using external boolean variable. This means that with a single call we can change, which internal variable would be used for reading and which for writing, for even a few thousands of class objects if we need it. This is probably the best solution for controlling class double-buffering.
NOTE1: If you need per class double-bufering control you need as many global variables as you have distict classes.
NOTE2: You could also save this information into the class itself but then you will have to update this information fora each class which could result in several thousands calls.

[Darkhog]

Cool demo. Could you share sources and exe file so other people can play with it?

[User137]

This is the only part where thread collisions can happen. It is possible that some force vectors are not added to movement, but is that likely?

            movement.x:=movement.x+vec.x;
            movement.y:=movement.y+vec.y;
            movement.z:=movement.z+vec.z;
            star[j].movement.x:=star[j].movement.x-vec.x;
            star[j].movement.y:=star[j].movement.y-vec.y;
            star[j].movement.z:=star[j].movement.z-vec.z;

Also the mainloop's part for threads is like:

    i2:=count div 4;
    i3:=i2*2;
    i4:=i2*3;
    TPhysicsThread.Create(self, 0, i2-1);
    TPhysicsThread.Create(self, i2, i3-1);
    TPhysicsThread.Create(self, i3, i4-1);
    TPhysicsThread.Create(self, i4, count-2);
    while threads>0 do begin
      sleep(1); // It takes usually over 140ms, so sleep here is not unnecessary
      Application.ProcessMessages;
    end;
...
constructor TPhysicsThread.Create(parent: TGame; const first, last: cardinal);
...
Thread execute:
...
    for i:=FFirst to FLast do
      with star[i] do
        for j:=i+1 to count-1 do begin

So it will not stop responding. I know it could be done better though, even so that fps wouldn't be affected. I also realize now that all the threads aren't getting equal workload at all. Especially the first thread gets biggest work, and 4th one the smallest.

Win32-binaries and source code here: https://docs.google.com/file/d/0B7FI...it?usp=sharing
Requires Lazarus and nxpascal to compile, possibly SVN version. All the required files are here:
https://code.google.com/p/nxpascal/s...%2Ftrunk%2Fsrc

[SilverWarior]

Also the mainloop's part for threads is like:

    i2:=count div 4;
    i3:=i2*2;
    i4:=i2*3;
    TPhysicsThread.Create(self, 0, i2-1);
    TPhysicsThread.Create(self, i2, i3-1);
    TPhysicsThread.Create(self, i3, i4-1);
    TPhysicsThread.Create(self, i4, count-2);
    while threads>0 do begin
      sleep(1); // It takes usually over 140ms, so sleep here is not unnecessary
      Application.ProcessMessages;
    end;
...
constructor TPhysicsThread.Create(parent: TGame; const first, last: cardinal);
...
Thread execute:
...
    for i:=FFirst to FLast do
      with star[i] do
        for j:=i+1 to count-1 do begin

So it will not stop responding. I know it could be done better though, even so that fps wouldn't be affected. I also realize now that all the threads aren't getting equal workload at all. Especially the first thread gets biggest work, and 4th one the smallest.

First of all get that Application.ProcessMessages out of that loop. It is not good to call Application.Process messages to often. And the reason why your main thread has so much load is becouse of Application.ProcessMessavges. It is quite resource consuming call. And all that processing time that is used by Application.ProcessMessages is taken away from one of your other threads making everything slower. And since you sad that you only have dual core the impact on overal speed can be quite big.
Wanna see how much Application.ProcessMessages affect your application. Start the simulation with setting sthat will cause low FPS, rotate the view and then pause the simulation. You will se that after you pause the simulation the overal view will quickly rotate to proper position while when being at low FPS it can be falling beind a litle.

[Darkhog]

I disagree that taking away ProcessMessages is solution. If you do that, application will "hang" (won't respond to paint request, button clicks, etc.). Though I agree it shouldn't be called every iteration. Maybe every tenth iteration or so.

[SilverWarior]

No the application still works OK but the performance gain is much smaler than I imagined it would be.
Yes I'm actually testing the source code to see where there may be major bottlenecks.