First, let me state two lemmas that demand tedious computations. Let $B(x, r)$ denote the circle centered at $x$ with radius $r$.
Lemma 1: Let $B(x,r)$ be a circle contained in $B(0, 1)$. Suppose we sample two points $p_1, p_2$ inside $B(0, 1)$, and $B(x', r')$ is the circle with $p_1p_2$ as diameter. Then we have
$$\mathbb{P}(x' \in x + dx, r' \in r + dr) = \frac{8}{\pi} r dxdr.$$
Lemma 2: Let $B(x,r)$ be a circle contained in $B(0, 1)$. Suppose we sample three points $p_1, p_2, p_3$ inside $B(0, 1)$, and $B(x', r')$ is the circumcircle of $p_1p_2p_3$. Then we have
$$\mathbb{P}(x' \in x + dx, r' \in r + dr) = \frac{24}{\pi} r^3 dxdr.$$
Furthermore, conditioned on this happening, the probability that $p_1p_2p_3$ is acute is exactly $1/2$.
Given these two lemmas, let's see how to compute the probability in question. Let $p_1, \cdots, p_n$ be the $n$ points we selected, and $C$ is the smallest circle containing them. For each $i < j < k$, let $C_{ijk}$ denote the circumcircle of three points $p_i, p_j, p_k$. For each $i < j$, let $D_{ij}$ denote the circle with diameter $p_i, p_j$. Let $E$ denote the event that $C$ is contained in $B(0, 1)$.
First, a geometric statement.
Claim: Suppose no four of $p_i$ are concyclic, which happens with probability $1$. Then exactly one of the following scenarios happen.
There exists unique $1 \leq i < j < k \leq n$ such that $p_i, p_j, p_k$ form an acute triangle and $C_{ijk}$ contains all the points $p_1, \cdots, p_n$. In this case, $C = C_{ijk}$.
There exists unique $1 \leq i < j \leq n$ such that $D_{ij}$ contains all the points $p_1, \cdots, p_n$. In this case, $C = D_{ij}$.
Proof: This is not hard to show, and is listed on wikipedia.
Let $E_1$ be the event that $E$ happens and we are in scenario $1$. Let $E_2$ be the event that $E$ happens and we are in scenario $2$.
We first compute the probability that $E_1$ happens. It is
$$\mathbb{P}(E_1) = \sum_{1 \leq i < j < k \leq n} \mathbb{P}(\forall \ell \neq i,j,k, p_\ell \in C_{ijk} , C_{ijk} \subset B(0, 1), p_ip_jp_k \text{acute}).$$
Conditioned on $C_{ijk} = B(x, r)$, Lemma 2 shows that this happens with probability $\frac{1}{2}r^{2(n - 3)} \mathbb{1}_{|x| + r \leq 1}$. Lemma 2 also tells us the distribution of $(x, r)$. Integrating over $(x, r)$, we conclude that
$$\mathbb{P}(\forall \ell \neq i,j,k, p_\ell \in C_{ijk} , C_{ijk} \subset B(0, 1), p_ip_jp_k \text{acute}) = \int_{B(0, 1)}\int_0^{1 - |x|} r^{2(n - 3)} \cdot \frac{12}{\pi} r^3 dr dx.$$
Thus we have
$$\mathbb{P}(E_1) = \binom{n}{3}\int_{B(0, 1)}\int_0^{1 - |x|} r^{2(n - 3)} \cdot \frac{12}{\pi} r^3 dr dx.$$
We can first integrate the $x$-variable to get
$$\int_{B(0, 1)}\int_0^{1 - |x|} r^{2(n - 3)} \cdot \frac{12}{\pi} r^3dr dx = 12 \int_0^1 r^{2n - 3}(1 - r)^2 dr.$$
Note that
$$\int_0^1 r^{2n - 3}(1 - r)^2 dr = \frac{(2n - 3)! * 2!}{(2n)!} = \frac{2}{2n * (2n - 1) * (2n - 2)}.$$
So we conclude that
$$\mathbb{P}(E_1) = \frac{n - 2}{2n - 1}.$$
We next compute the probability that $E_2$ happens. It is
$$\mathbb{P}(E_2) = \sum_{1 \leq i < j \leq n} \mathbb{P}(\forall \ell \neq i,j, p_\ell \in D_{ij} , D_{ij} \subset B(0, 1)).$$
Conditioned on $D_{ij} = B(x, r)$, this happens with probability $r^{2(n - 2)} \mathbb{1}_{|x| + r \leq 1}$. Lemma 1 tells us the distribution of $(x, r)$. So we conclude that the probability that
$$\mathbb{P}(\forall \ell \neq i,j, p_\ell \in D_{ij} , D_{ij} \subset B(0, 1)) = \int_{B(0, 1)}\int_0^{1 - |x|} r^{2(n - 2)} \cdot \frac{8}{\pi} r dr dx.$$
So
$$\mathbb{P}(E_2) = \binom{n}{2}\int_{B(0, 1)}\int_0^{1 - |x|} r^{2(n - 2)} \cdot \frac{8}{\pi} r dr dx.$$
We compute that
$$\int_{B(0, 1)}\int_0^{1 - |x|} r^{2(n - 2)} \cdot \frac{8}{\pi} r dr dx = 8\int_0^1 r^{2n - 3} (1 - r)^2 dr = 8 \frac{(2n - 3)! 2!}{(2n)!}.$$
So we conclude that
$$\mathbb{P}(E_2) = 8 \binom{n}{2} \frac{(2n - 3)! 2!}{(2n)!} = \frac{2}{2n - 1}.$$
Finally, we get
$$\mathbb{P}(E) = \mathbb{P}(E_1) + \mathbb{P}(E_2) = \boxed{\frac{n}{2n - 1}}.$$
The proofs of the two lemmas are really not very interesting. The main tricks are some coordinate changes. Let's look at Lemma 1 for example. The trick is to make the coordinate change $p_1 = x + r (\cos \theta, \sin \theta), p_2 = x + r (-\cos \theta, -\sin \theta)$. One can compute the Jacobian of the coordinate change as something like
$$J = \begin{bmatrix} 1 & 0 & 1 & 0 \\
0 & 1 & 0 & 1 \\
\cos \theta & \sin \theta & -\cos \theta & -\sin \theta \\
-r\sin \theta & r\cos \theta & r\sin \theta & -r\cos \theta \\
\end{bmatrix}.$$
And we can compute that $|\det J| = 4r$. As $p_1, p_2$ has density function $\frac{1}{\pi^2} \mathbf{1}_{p_1, p_2 \in B(0, 1)}$, the new coordinate system $(x, r, \theta)$ has density function
$$\frac{4r}{\pi^2} \mathbf{1}_{p_1, p_2 \in B(0, 1)}$$
The second term can be dropped as it is always $1$ in the neighbor of $(x, r)$. To get the density of $(x, r)$ you can integrate in the $\theta$ variable to conclude that the density of $(x, r)$ is
$$\frac{8r}{\pi}$$
as desired.
The proof of Lemma 2 is analogous, except you can use the more complicated coordinate change from $(p_1, p_2, p_3)$ to $(x, r, \theta_1, \theta_2, \theta_3)$
$$p_1 = x + r (\cos \theta_1, \sin \theta_1), p_2 = x + r (\cos \theta_2, \sin \theta_2), p_3 = x + r (\cos \theta_3, \sin \theta_3).$$
The Jacobian $J$ is now $6$ dimensional, and Mathematica tells me that its determinant is
$$|\det J| = r^3|\sin(\theta_1 - \theta_2) + \sin(\theta_2 - \theta_3) + \sin(\theta_3 - \theta_1)|.$$
So we just need to integrate this in $\theta_{1,2,3}$! Unfortunately, Mathematica failed to do this integration, but I imagine you can do this by hand and get the desired Lemma.
Best Answer
This is equivalent to the problem of packing equal circles in a large circle. Suppose the smallest enclosing circle in your problem has radius $r$, then the smallest circle in which $n$ unit circles can be packed has radius $2r+1$.
Some solutions for the latter problem are given here and here.