Okay, try that program in post #59 first, then try this. It adds another loop, so a two body system might take a little longer than that last one, but a multi-body system should run more quickly by shortening the main loop, exponentially with a greater number of bodies, without affecting the rest of your program. By the way, I just noticed you've got h=0 for that second loop which I've changed to h=1 and colored in blue. It wouldn't make a difference to your program values, but it causes the loop to run one extra time unnecessarily, so that'll save you more time right there.
Code:
``` for (int h=1;h<=NumObjects;h++) {
objvx[h] = objvx[h] * ObjMass[h];
objvy[h] = objvy[h] * ObjMass[h];
objvz[h] = objvz[h] * ObjMass[h];
}

for (int k=1;k<=NumObjects;k++) {
for (int j=k+1;j<=NumObjects;j++) {
dx = objx[j] - objx[k];
dy = objy[j] - objy[k];
dz = objz[j] - objz[k];
D = sqrt(dx*dx+dy*dy+dz*dz);
if  (((objsize[k] + objsize[j])/2) > D) {
Collisions = Collisions + 1;
CollisionObjectA[Collisions] = j;
CollisionObjectB[Collisions] = k;
}
f = ObjMass[j] * ObjMass[k] * 398600440000000/D/D/D;
if (f > 0) {
objvx[j] = objvx[j] - dx*f;
objvy[j] = objvy[j] - dy*f;
objvz[j] = objvz[j] - dz*f;
objvx[k] = objvx[k] + dx*f;
objvy[k] = objvy[k] + dy*f;
objvz[k] = objvz[k] + dz*f;
}

}
}
for (int h=1;h<=NumObjects;h++) {
objvx[h] = objvx[h]/ObjMass[h];
objvy[h] = objvy[h]/ObjMass[h];
objvz[h] = objvz[h]/ObjMass[h];
objx[h] = objx[h] + objvx[h];
objy[h] = objy[h] + objvy[h];
objz[h] = objz[h] + objvz[h];
}```