But this is irrelevant. I was pointing out the flaw in the high-level strategy:

while (condition)
   add_point(P, Q)

which is a serial one-by-one element addition. If you really aim for higher speed you need to rethink what you are doing, not pretend that some black-box "addPoints(a, b)" can do wonders for you. The black-box does whatever the formula to add 2 points does, and has no knowledge of your context (adding many points).

Here's a clue:

z = p[0].x - jp[0].x
for (i = 1 .. n) {
  d = p[i ].x - jp[i ].x
  z *= d
  q[i ] = z
}

t = inv(z)
for (i = n - 1 ... 0) {
   xd = i > 0 ? t * q[i-1] : t
   t *= q[i ]
  // finish addition..
}

This trades 3(N-1) multiplications with a single inversion
Depending on size of N you can get a speedup up to the factor of around cost(1 inversion) / cost(3 multiplications).