I haven't implemented the algorithm myself, but Gino Van den Bergen does write about the answer to your question in his book "Collision Detection in Interactive 3D Environments", iirc. Also, iirc, EPA will always only return a single point of intersection -- not multiple points.

My hunch would be that if you can grab the penetration depth through your implementation, then the best point of contact to represent the collision would be the deepest point. You should be able to just pull that right out of your EPA data once the algorithm terminates. If not then you can possibly just query a support point along the intersection normal.

IMHO EPA has too many numerical instability for closing contacts. Also, you'll need a contact caching mechanism later (which doesn't have the same quality as an single-shot contact generation method such SAT + Sutherland-Hodgman according to my research).

Here you can find an old naive GJK-SAT implementation I wrote. No collision skins added.

To get the closest points, I've followed this EPA approach. Assumming you have a tetahedron that encloses the origin of the MD at the end, you expand the polytope removing and adding vertex as described in the link. The data structure he uses is not the best one (he's not using structures such the half-edge), and the insert/removing operations on the polytope may be slow most of the time (cylinders, spheres, curved objects).

The data structure he uses is not the best one (such the half-edge)

What is bad about the half-edge data structure and what data structure do you recommend instead?

The data structure he uses is not the best one (such the half-edge)

What is bad about the half-edge data structure and what data structure do you recommend instead?

Oh! I was suppose to mean "he's not using structures such the half-edge" (tried to connect the "best-one-such")).

Edited.