How can I find the position of the center point of a circle tangent to two ellipses?

Hello everyone,

Given the information about two ellipses, I want to know how to find a circle with a fixed radius tangent to the ellipse.

When we know information about ellipses E1 and E2, we want to know the center point of circle C.

I tried several methods, but with no success.

If you look at the link below, the ellipse is expressed as a matrix and the straight line tangent to the two ellipses is calculated using solve. I wonder if the circle can be obtained in the same way.

https://kr.mathworks.com/matlabcentral/answers/635464-how-to-find-common-tangents-between-2-ellipses

[Catalytic]

Even if it could, how would that be useful?

Answer 1

[Matt J]

I would probably solve for the circle parameters v1,v2 plus the two points of tangency x1,y1,x2,y2. That's 6 unknowns, and you already have 4 equations. You can get 2 more equations from the tangency conditions. The conditions are,

=0

and similarly for (x2,y2).

[Engineer_Park]

Thanks a lot. Can you explain a little bit more about the formula?

[Matt J]

The gradients of the ellipse equations are perpendicular to the tangents of the ellipse. So, if the gradient of both the circle and the ellipse are parallel to one another at a point, they must share a tangent line at that point.

[Engineer_Park]

Thank you.

but, I don't know how to apply the above equation.

You mean finding the Jacobian matrix and its determinant?

If so, can I ask if we need to find the determinants of ellipse E and circle C respectively, or how to combine them?

[Matt J]

Yes, you can think of it as setting to zero the Jacobian determinants of [E1;C] and [E2;C]. Each 2x2 determinant will give you an additional equation. The Symbolic Math Toolbox supports both jacobian() and det(). It should be a pretty direct translation into code.

https://www.mathworks.com/help/symbolic/sym.jacobian.html

https://www.mathworks.com/help/symbolic/sym.det.html?searchHighlight=det&s_tid=srchtitle_det_2

[Engineer_Park]

I'm truly grateful for your kindness.

I wrote the code below.

but it didn't work when Matlab finds roots (last line)

What did I do wrong..?

%% Ellipse information

% E1 coeff
% ax^2 + bxy + cy^2 + dx + ey + f=0
a=0.06137;
b=-0.61943;
c=6.23450;
d=0;
e=0;
f=-0.28671;
% 
% % E2 coeff
% % a1x^2 + b1xy + c1y^2 + d1x + e1y + f1=0
a1=0.16785;
b1=1.00923;
c1=6.20813;
d1=0.74001;
e1=3.73117;
f1=0.1492;

%% Ellipse plot
syms x y
% ellipse3=(x*cs1+y*sn1)^2/(a_1^2)+(-x*sn1+y*cs1)^2/(b_1^2)-1==0
ellipse1=a*x^2 + b*x*y + c*y^2 + d*x + e*y + f == 0
ellipse2=a1*x^2 + b1*x*y + c1*y^2 + d1*x + e1*y + f1 == 0

figure(1)
fimplicit(ellipse1)
% fimplicit(ellipse3)
hold on
fimplicit(ellipse2)

%% Find circles
r=50; % circle radius
% x0, y0 : Points of tangency (E1, C)
% x1, y1 : Points of tangency (E2 ,C)
syms x0 y0 x1 y1 v1 v2 % 6 unknowns

equ1=a*x0^2 + b*x0*y0 + c*y0^2 + d*x0 + e*y0 + f == 0 % E1(x0,y0)
equ2=a1*x1^2 + b1*x1*y1 + c1*y1^2 + d1*x1 + e1*y1 + f1 == 0 % E2(x1,y1)
equ3=(x0-v1)^2+(y0-v2)^2-r^2==0 % C(x0,y0)
equ4=(x1-v1)^2+(y1-v2)^2-r^2==0 % C(x1,y1)

% Jacobian Matrix
Jacob_E1=jacobian(a*x0^2 + b*x0*y0 + c*y0^2 + d*x0 + e*y0 + f,[x0,y0]); % Jacobian Matrix of E1(x0,y0)
Jacob_E2=jacobian(a1*x1^2 + b1*x1*y1 + c1*y1^2 + d1*x1 + e1*y1 + f1,[x1,y1]); % Jacobian Matrix of E2(x1,y1)
Jacob_C_xy=jacobian((x0-v1)^2+(y0-v2)^2-r^2,[x0,y0]); % Jacobian Matrix of C(x0,y0)
Jacob_C_x1y1=jacobian((x1-v1)^2+(y1-v2)^2-r^2,[x1,y1]); % Jacobian Matrix of C(x1,y1)

Jacob_combine_E1C=[Jacob_E1 ; Jacob_C_xy]; % jacobian of E1(x0,y0),C(x0,y0)
Jacob_combine_E2C=[Jacob_E2 ; Jacob_C_x1y1]; % jacobian of E2(x1,y1),C(x1,y1)
Det_Jacob_E1C=det(Jacob_combine_E1C); % Determinant of [E1(x0,y0) ; C(x0,y0)]
Det_Jacob_E2C=det(Jacob_combine_E2C); % Determinant of [E2(x0,y0) ; C(x0,y0)]

equ5=Det_Jacob_E1C==0; % Determinant of (jacbian(E1) ; jacobian(C))
equ6=Det_Jacob_E2C==0; % Determinant of (jacobian(E2) ; jacobian(C))

sol=solve([equ1,equ2,equ3,equ4,equ5,equ6]); % find roots