I realized today that while 3D is great, is not necessary to simulate the gravity particle effect, or should you call it "anti-gravity". This can be done in 2D, with big gain in performance, and better results to look at.

Now after testing it a little, it didn't seem to work somehow, or maybe i'm missing something. I also think that i'm able to rid of Sqrt() completely.
Code:
    for i:=FFirst to FLast do
      with star[i] do
        for j:=i+1 to count-1 do begin
          xx:=v.x-star[j].v.x;
          yy:=v.y-star[j].v.y;
          d:=xx*xx+yy*yy;
          if d>0 then begin
            d:=(G_FORCE*G_FORCE)/(d*d*d); // Compiler should optimize the constant multiplication to 1 value
            // If not, just remove multiplication and change G_FORCE value
            if star[j].positive<>positive then d:=-d;
            vec.x:=(star[j].v.x-v.x)*d;
            vec.y:=(star[j].v.y-v.y)*d;
            movement.x:=movement.x+vec.x;
            movement.y:=movement.y+vec.y;
            star[j].movement.x:=star[j].movement.x-vec.x;
            star[j].movement.y:=star[j].movement.y-vec.y;
          end;
        end;
The particles scatter endlessly, with no clear groupings done as was in 3D. I have tried with many different parameters. I wouldn't want to make boundaries to the "universe", where particles would bounce backwards, but i might have to. I also tried with every particle being same polarity, and they still did not group. It's not because of removal of sqrt(), that was almost last step in my tests, improvement to calculation time i'd say.