// Unit vector in the direction of the collision

double ax=dx/distance, ay=dy/distance;

// Projection of the velocities in these axes

double va1=(vx1*ax+vy1*ay), vb1=(-vx1*ay+vy1*ax);

double va2=(vx2*ax+vy2*ay), vb2=(-vx2*ay+vy2*ax);

// New velocities in these axes (after collision): ed<=1, for elastic collision ed=1

double vaP1=va1 + (1+ed)*(va2-va1)/(1+mass1/mass2);

double vaP2=va2 + (1+ed)*(va1-va2)/(1+mass2/mass1);

// Undo the projections

vx1=vaP1*ax-vb1*ay; vy1=vaP1*ay+vb1*ax;// new vx,vy for ball 1 after collision

vx2=vaP2*ax-vb2*ay; vy2=vaP2*ay+vb2*ax;// new vx,vy for ball 2 after collision

In the above code you only deal with x, and y direction, but i want to try to also include a z direction as well and

I do not know how to do the calculation with x, y, and z direction since the ball I'm working on has an x, y, and z position vector.