Notes
The explicit form of the lower bound in this answer may be wrong. In the updated numerical results from Mark, it can produce a configuration with a $r$ smaller than the lower bound here. It looks like the constraint on one of the points (the most likely candidate is point $F$) is redundant. We are essentially back to square one.... C'est la vie....
Please consider this as a supplement to Mark H's answer.
In Mark's answer, the centers are roughly symmetrical in the horizontal direction.
If one compute the Voronoi diagram associated with these centers, one obtain a figure
looks roughly what is shown below.
$\hspace0.75in$ 
The points $A, B, C, D$ are those $4$ centers in the first quadrant.
The orange lines are the boundary of the Voronoi cells. The points
$E, F, G, H, I$ are the intersections of the boundaries of these cells and the semi-circle. Let $O = (0,0)$ and $X = (1,0)$.
The key observation is the distances (those illustrated in pink) among these points
are all roughly equal. If the algorithm in Mark's answer converge to some configuration of centers, that configuration should be a local minimum of the maximum-minimum distance functional. This means these distances should be equal to each other exactly. i.e.
$$|AX| = |AE| = |AF| = |AG| = |BE| = |BF| = |BH|\\
= |CF| = |CG| = |CH| = |CI| = |DH| = |DI|$$
To find such a configuration, we first relax the constraint that $D$ lies on $y$-axis. We assume $A$ is located near what's in Mark's answer.
Let $r$ be the common values of above $13$ distances. Let $A = (1-u,v)$ and $B = (0,w)$. The condition $|AX| = r$ leads to $r^2 = u^2 + v^2$.
We then proceed to express the positions of $E, F, G, C, H, I$ (in that order) and finally $D$ in terms of these 3 variables $u,v,w$.
Aside from the formula of $I$ and $D$, they are not that horrible.
- $|AE| = r$ implies $E = (1-2u,0)$.
- $|BE| = r$ leads to
$$w^2 + (1-2u)^2 - u^2 - v^2 = 0\tag{*1}$$
- $|AE| = |BE| = |AF| = |BF|$ implies $F = A + B - E = (u,v+w)$.
- $|AX| = |GX|$ and $|OX| = |OG|$ implies
$$G = \verb/Refl/(A,X) =
\left(\frac{2(1-u)^2}{v^2+(1-u)^2}-1,\frac{2(1-u)v}{v^2+(1-u)^2}\right)$$
where $\displaystyle\;\verb/Refl/(\vec{U},\vec{V}) = 2\frac{\vec{V}\cdot\vec{U}}{|\vec{U}|^2}\vec{U} - \vec{V}$ maps point $V$ to its mirror image with respect to $OU$.
- $|AF| = |AG| = |CF| = |CG|$ implies
$$C = F + G - A = \left(\frac{2(1-u)^2}{v^2+(1-u)^2}+2u-2,w+\frac{2(1-u)v}{v^2+(1-u)^2}\right)$$
$|BF| = |BH| = |CF| = |CH|$ implies
$$H = B + C - F =
\left(\frac{2(1-u)^2}{v^2+(1-u)^2}+u-2,w+ \frac{2(1-u)v}{v^2+(1-u)^2}-v\right)$$
$|CG| = |CI|$ and $|OG| = |OI|$ implies $ I = \verb/Refl/(C,G) = $ a horrible mess!
- $|CH| = |CI| = |DH| = |DI|$ implies $D = H + I - C = $ another horrible mess!
If we put back the constraint that $D$ lies on the $y$-axis, we obtain
$$((u-1)v^2+u^3-u^2-u+1)w^2 + (4u^2-8u+4)vw\\
+ (4u^3-4u^2-4u+4)v^2 + 4u^5-12u^4+12u^3-4u^2\\
= 0\tag{*2}$$
We can eliminate $w$ by computing the resultant between the two polynomials in $(*1)$ and $(*2)$.
The resultant has the form $(u-1)^2f(u,v)$ where $f(u,v)$ is a polynomial of degree $8$ in $u, v$.
We can simplify this expression a little bit by a change of variable.
Let $(u,v) = (rs, r\sqrt{1-s^2})$,
the condition becomes
$$\begin{align}
g(r,s) \stackrel{def}{=} &\;\; f(rs,r\sqrt{1-s^2})\\
= &\;\;16r^6s^4+(-16r^7-32r^5-16r^3)s^3+(24r^6+48r^4+24r^2)s^2\\
&\;\; + (8r^7-24r^5-40r^3-8r)s+r^8-12r^6+22r^4+4r^2+1\\
= &\;\; 0
\end{align}
$$
To complete our task, we need to find the point along the curve $g(r,s) = 0$
with minimal $r$. At that point, the tangent vector of the curve will be pointing in the $s$-direction. This means $\frac{\partial g}{\partial s}(r,s) = 0$ at that particular point.
We can eliminate $s$ by computing the resultant between $g(r,s)$ and $\frac{\partial g}{\partial s}(r,s)$. The resultant is reasonably simple
$$-1048576 r^{18} (r^2-1)^{12} (27r^8 - 324 r^6 - 270 r^4-36 r^2 + 11)$$
Eliminating the impossible values $0$ and $1$ for $r$, the local extremium
of maximum-minimum distance functional should be a root of the quartic equation:
$$27r^8 - 324 r^6 - 270 r^4-36 r^2 + 11 = 0$$
Solving this equation and compare the numerical values we find from Mark's answer, the $r$ we seek equals to
$$r = \sqrt{ 3 + 2\sqrt{3} - 2\sqrt{5 + \sqrt{27} - \frac{1}{\sqrt{27}}}} \approx 0.3719887399641683$$
and at the least, this is a local minimum for the maximum-minimum distance functional. With this, we can numerically back out the parameters $u,v,w$
$$\begin{cases}
u &\approx 0.3259601005065833\\
v &\approx 0.1792362562035586\\
w &\approx 0.1312100460994327
\end{cases}$$
and the location of the centers
$$\left\{\begin{array}{lcll}
x_5 \leftrightarrow A & \approx & (0.6740398994934167,& 0.1792362562035586),\\
x_1 \leftrightarrow B & \approx & (0,& 0.1312100460994327),\\
x_3 \leftrightarrow C & \approx & (0.5198397097356646, & 0.6279149145866445),\\
x_2 \leftrightarrow D & \approx & (0,& 0.7661472706667444)
\end{array}\right.$$
