# Cubic polynomial smoothly connecting two circles

Given two circles with radii $$R_L$$ and $$R_R$$ and centers at $$(-(R_L+a),\,0)$$ and $$(R_R+a,\,0)$$, respectively, find a cubic polynomial $$p(x)=b+cx^2+dx^3$$ that smoothly connects the two circles.

$$b$$ is a parameter so $$p(0)=b$$ and the linear term of the polynomial is omitted because we want $$\frac{\mathrm{d}p}{\mathrm{d}x}\Big|_{x=0}=0$$.

My attempt at a solution. Let $$\mathrm{C}_{L,R}$$ the equations for the upper half of the $$L,R$$ circles. I formulate two equations relating $$\mathrm{C}_{L,R}$$ and $$p$$, and two equations relating $$\frac{\mathrm{d}}{\mathrm{d}x}\mathrm{C}_{L,R}$$ and $$\frac{\mathrm{d}p}{\mathrm{d}x}$$. Let $$x_{L,R}$$ be the points where $$p(x)$$ and $$\mathrm{C}_{L,R}(x)$$ intersect, then:

$$\mathrm{C}_L(x_L)-p(x_L) = 0$$ $$\mathrm{C}_R(x_R)-p(x_R) = 0$$ $$\frac{\mathrm{d}}{\mathrm{d}x}\mathrm{C}_L(x_L) - \frac{\mathrm{d}p}{\mathrm{d}x}(x_L) = 0$$ $$\frac{\mathrm{d}}{\mathrm{d}x}\mathrm{C}_R(x_R) - \frac{\mathrm{d}p}{\mathrm{d}x}(x_R) = 0$$ so I have a system of four nonlinear equations with four unknowns $$(x_L,\,x_R,\,c,\,d)$$. I coded a simple Newton's method for the system and it works well for some combinations of parameters $$(a,\,b,\,R_L,\,R_R)$$ when the initial guess is close enough, especially when $$|R_L-R_R|$$ is not too large and I use a constant damping for the Newton iterations. I can find initial guesses that I think are good via a graphical interface I coded. However, as $$|R_L-R_R|$$ gets larger the solver fails spectacularly to converge even with very close initial guesses and very small damping. (I should add that I'm actually taking the square of the equations to avoid square roots of negative numbers during the Newton iterations).

My question is threefold:

a) what other method or modification can I use to make the solver more stable?

b) this problem seems to me like it should be solved somewhere, do you know a reference?

c) more generally, is there a reason this should fail as horribly as it does when $$|R_L-R_R|>>1$$?

• I would try parametrizing this with only two unknowns $x_{L,R}$, fitting a cubic polynomial $a_0+a_1x+a_2x^2+a_3x^3$ to the points and tangents on the circles, and solving for $a_0=b,a_1=0$ instead. – user856 Oct 14 '19 at 16:48
• Thank you. I'm not sure I follow. If I did that, wouldn't $a_2$ and $a_3$ be unknowns as well? If I had two unknowns and four equations I'm not even sure what I could do. Also, setting $a_0=b$ and $a_1=0$ is exactly what I'm doing. By fitting do you mean I should try some optimisation technique? – mvaldez Oct 14 '19 at 17:03
• Do you mean that parameter $b$ si fixed : all curves must pass through point $(0,b)$ ? – Jean Marie Oct 14 '19 at 17:05
• yes, the parameters $(a,\,b,\,R_L,\,R_R)$ are fixed so the curve I want to find must pass through $(0,\,b)$ and have vanishing derivative there – mvaldez Oct 14 '19 at 17:07
• I mean treat only the points $x_{L,R}$ as unknowns. Compute a cubic polynomial $a_0+a_1x+a_2x^2+a_3x^3$ passing through them, which can be done in closed form. Then you get $c=a_2,d=a_3$. But you may not have $a_0=b,a_1=0$, so you need to choose $x_{L,R}$ to satisfy them. Thus you have two nonlinear equations in two unknowns. – user856 Oct 14 '19 at 17:09

Here is a method giving all solutions in a deterministic way.

Take a look at the following figure : Fig. 1 : The exhaustive set of 8 solution curves in the case $$a=1, \ b=1, \ R_L=3, \ R_R=1$$. Please note that some tangencies are internal to the circles.

How has this result been obtained ? By using a Computer Algebra System with a rather simple system of 4 polynomial equations in the 4 unknowns $$c,d,x_L,x_R$$ where (your notations) the two last parameters are the abscissas of the tangency points of the cubic curve with equation

$$y=f(x)=b+cx^2+dx^3$$

with the Left and Right circle resp.

Here they are:

$$\begin{cases}(x_L+R_L+a)^2+f(x_L)^2&=&R_L^2& \ \ (i)\\ (x_R-R_R-a)^2+f(x_R)^2&=&R_R^2& \ \ (ii)\\ \dfrac{f(x_L)}{x_L+R_L+a}&=&- \dfrac{1}{f'(x_L)}& \ \ (iii)\\ \dfrac{f(x_R)}{x_R-R_R-a}&=&- \dfrac{1}{f'(x_R)}& \ \ (iv)\\ \end{cases}\tag{1}$$

The first two equations express that $$M_L=(x_L,f(x_L))$$ and $$M_R=(x_R,f(x_R))$$ belong to their resp. circles.

Let us now explain the third equation. Let $$C_L(-(R_L+a),0)$$ denotes the center of the first circle. The slope of the tangent in $$M_L$$ is $$f'(x_L)$$. $$\vec{C_LM_L}=\binom{x_L+R_L+a}{f(x_L)-0}$$, being orthogonal to this tangent, has a slope $$-\dfrac{1}{f'(x_L)}$$ (see there). Similar reasoning for the fourth equation.

Remarks :

1) I haven't found the limitations you mention in the case of a large gap between $$R_L$$ and $$R_R$$. For example when $$R_L=100$$ and $$R_R=1$$, (with $$a=1$$ and $$b=0$$), one finds 14 solutions...

2) The third equation can be written : $$f(x_L)f'(x_L)+x_L+R_L+a=0$$ : this avoids possible division by zero. Same thing for the fourth equation.

Another case (Fig. 2), this time with $$R_L=R_R$$ displaying some spurious solutions under the form of ... second degree curves, i.e., parabolas (one can indeed have $$d=0$$) : Fig. 2.

Still another case, with four solutions (one of them looks to be a parabola but isn't): Fig. 3.

Here is the Matlab program that has given figure 1
(running time : around 30 seconds on my rather slow computer). Please note the "isreal" conditions (we want the different unknowns to be real) : in fact in the first case, there are $$40$$ solutions, $$32$$ of them being spurious, i.e., with complex coefficients...

Final remark : In fact, there exists a different way to solve system (1) by doing it in two steps ; first, by expressing the fact that equations (i) and (iii) have a common root $$x_L$$ giving (using a "resultant" ) a first (non-linear) constraint between $$c$$ and $$d$$. Doing the same for equations (ii) and (iv) for common root $$x_R$$, gives a second constraint. Then, in a second step, solve the system of 2 (non-linear...) equations in 2 unknowns $$c$$ and $$d$$.

syms  xL xR c d;%symbolic variables declaration
a=1;b=1;m=max([RL,RR]);p=2*m+a;
axis([-p+3,p+3,-p-3,p+3]);axis equal;
f=@(x)(b+c*x*x+d*x*x*x);%cubic function
fp=@(x)(2*c*x+3*d*x*x);%its derivative
%The system of 4 equations in 4 unknowns ; sol. in  [XL,XR,C,D]
[XL,XR,C,D]=solve(...
(xL+RL+a)^2+f(xL)^2==RL^2,
(xR-RR-a)^2+f(xR)^2==RR^2,...
xL+f(xL)*fp(xL)==-RL-a,...
xR+f(xR)*fp(xR)==RR+a,...
xL,xR,c,d);% the 4 variables
plot([-p,p],[0,0]);\$x axis
plot([0,0],[-p,p]);%y axis
UC=exp(i*(0:0.01:2*pi));%unit circle
plot(-(RL+a)+RL*UC,'k');%circle C_L
plot((RR+a)+RR*UC,'k');%circle C_R
x=-p:0.01:p;%range of values for x
for k=1:length(XL)
c=C(k);d=D(k);xL=XL(k);xR=XR(k);
if isreal(c) && isreal(d) && isreal(xL) && isreal(xR)
plot(x,b+c*x.*x+d*x.*x.*x,'color',rand(1,3));hold on;
end;
end;

• Thank you Jean Marie, your answer was very helpful. This new set of equations seem to be much more stable – mvaldez Oct 22 '19 at 20:45
• Thanks. I have added a small remark (remark 2). – Jean Marie Oct 23 '19 at 7:06
• See as well the final remark I just added. – Jean Marie Oct 23 '19 at 7:42