You pointed out one error in your own answer, here's another. I haven't figured out the adjustments, but presumably they cancel.
You've done
∫₀¹ ∫₀¹⁻ᶜ (2π r)³ dr dc
However,
- the `dr` integral is assuming that the radii are uniformly likely in [0, r]
- the `dc` integral is assuming that the centers are uniformly likely in [0, c]
You need to wait these integrals by the conditional probability distributions.
Sorry, this was incorrect. In this integral we're just calculating the geometrical question, what is is the 5d volume of the collection of 3 points whose circumcenter is on {(x, 0) : x in [0, 1]}. So the "overall" probability of the problem has nothing to do with it.