Iterations

At the start of each iteration, we have an ellipse point $P(\varphi_n)$ and the corresponding inclination $\varphi_n$, the signed distance $k_n\ne0$ between $P(\varphi_n)$ and the test point $A$ and the value of $t_n=\tan{\zeta_n/2}$ where $\zeta_n$ is the slope of the line joining $P(\varphi_n)$ to the test point.

If the estimate $P(\varphi_n)$ were perfect, we would have $\varphi_n=\varphi$ and $k_n=h$. If $h\ne0$ the circle centered on the test point $A$ with radius $\vert k_n\vert$ would be tangent to the ellipse. In our iterative process, the estimate is not perfect and the circle intersects the ellipse at least at one other point $P(\tilde{\varphi}_n)$, $\tilde{\varphi}_n\ne\varphi_n$. Figures 3 and 4 show these circle/ellipse intersections in two different situations. As shown by these figures, there can be up to three intersection points in addition to the already known point $P(\varphi_n)$.

Figure 3: interior point

Figure 4: exterior point

The additional intersection points can be computed by introducing $\tau=\tan{\zeta/2}$ in equation (5) an multiplying by $(1+\tau^2)^2$. We get:

$$\begin{align*}\begin{aligned}[t]k^2 &[(1+\tau^2)^2-f(2-f)(1-\tau^2)^2]\\ -2k &[(... ...\tau^2)]\\ + &[(1-f)^2(r^2-a_e^2)+z^2](1+\tau^2)^2 = 0 \end{aligned}\end{align*}$$

$$\Leftrightarrow \phantom{+}\alpha\tau^4+\beta\tau^3+\gamma\tau^2+\delta\tau+\varepsilon = 0$$

$$\quad \left\{\begin{aligned}\alpha &= (1-f)^2((r+k)^2-a_e^2)+z^2\\ \bet... ... \delta &= -4kz\\ \varepsilon &=(1-f)^2((r-k)^2-a_e^2)+z^2 \end{aligned}\right.$$

This is a quartic equation for which we already know one root: $t_n=\tan{\zeta_n/2}$. We can therefore reduce this equation to a normalized cubic polynomial by a simple coefficients identification (the value of the constant term $\varepsilon$ of the quartic is of course not needed for this degree reduction):

$$\tau^3+a_1\tau^2+a_2\tau+a_3=0 \quad \left\{\begin{aligned}a_1 &= \frac{\beta}{\alpha} + t_n\\... ...t_n\\ a_3 &= \frac{\delta}{\alpha} + a_2 t_n \end{aligned}\right.$$

The other circle/ellipse intersections are given by the real roots of equation (8). In order to compute them, we compute first the discriminant $D$:

$$Q=\frac{3a_2-a_1^2}{9} \qquad R=\frac{9a_1a_2-27a_3-2a_1^3}{54} \qquad D=Q^3+R^2$$

If $D$ is positive, the following expressions compute the two real numbers $S$ et $T$ and allow to deduce the unique real root $\tilde{t}_{n,1}$. We do not need the two complex roots $\tilde{t}_{n,2}$ and $\tilde{t}_{n,3}$. This corresponds to the case depicted in figure 4. It also occurs for interior points when we are close to the solution.

$$S=\sqrt[3]{R+\sqrt{Q^3+R^2}} \qquad T=\sqrt[3]{R-\sqrt{Q^3+R^2}}$$

$$\left\{\begin{aligned}\tilde{t}_{n,1} &= S+T-\frac{a_1}{3}\\ \til... ...&= -\frac{S+T}{2}-\frac{a_1}{3}-i\frac{\sqrt{3}(S-T)}{2}\\ \end{aligned}\right.$$ (9)

If $D$ is negative, the three roots are all real ones. This corresponds to the case depicted in figure 3. In this case, the trigonometric expressions of the root are simpler to use:

$$\left\{\begin{aligned}\tilde{t}_{n,1} &= 2\sqrt{-Q}\cos{\frac{\th... ...} &= 2\sqrt{-Q}\cos{\frac{\theta+4\pi}{3}}-\frac{a_1}{3}\\ \end{aligned}\right. \quad \cos\theta = \frac{R}{\sqrt{-Q^3}}$$ (10)

For each $\tilde{t}_{n,i}$ representing a circle/ellipse intersection, we compute $\tilde{\zeta}_{n,i}$ and $P(\tilde{\varphi}_{n,i})$ using equation (4). Since this point belongs to the ellipse we deduce $\tilde{\varphi}_{n,i}$ using equation (2).

When there are three roots, only one of them² leads to a point on the same half plane ($z<0$ or $z>0$) as the test point $A$ and the already known intersection point $P(\varphi_n)$. We select this root and disgard the two other ones. This selection is done by comparing the signs of $\varphi_n$ and $\tilde{\varphi}_{n,i}$.