Solving intersection point of two circles using sympy in python
I'm trying to find the intersection point of two functions, which are in this case circle1 and circle2 with a radius of r. When I try to find the intersection point, it takes Sympy forever to calculate.
Here is an image with the intersection:
I think this is because there are an infinite number of intersection points, when theta_1 and theta_2 are raised. Thus I guess sympy keeps solving it until infinity.
Unfortunately, solve or solveset dont allow to solve an equation for two variables within a given a range.
The code I used is:
def SYMPY():
theta_1, theta_2, r = symbols('theta_1, theta_2, r', real=True)
r = 10
x1 = r*cos(theta_1) + 5
x2 = r*cos(theta_2)
y1 = r*sin(theta_1)
y2 = r*sin(theta_2)
test = sp.Eq((x1-x2)+(y1-y2),0)
ans = sp.solve(test, theta_1, theta_2)
print(ans)
SYMPY()
Does anybody know how to tackle this problem?
I know Sympy also has the circle() method, but I would like to calculate it using the equations above. This is because I have another (large) equation that I would like to solve, that also involves two circles and an intersection point.
Thanks, Mudi :)
Solution 1:
To test whether (x1,y1) is the same point as (x2,y2), Eq((x1-x2)+(y1-y2),0)doesn't do the job. You need both Eq(x1,x2) and Eq(y1,y2). (Note that you shouldn't make r a symbolic variable, as you assign a constant to it. Also note that (theta_1, theta_2) should be written as a tuple (or a list)).
For some reason, unclear to me, sympy doesn't find a solution:
from sympy import symbols, sin, cos, Eq, solve
theta_1, theta_2 = symbols('theta_1, theta_2', real=True)
r = 10
x1 = r * cos(theta_1) + 5
x2 = r * cos(theta_2)
y1 = r * sin(theta_1)
y2 = r * sin(theta_2)
test = [Eq(x1, x2), Eq(y1, y2)]
ans = solve(test, (theta_1, theta_2))
print(ans) # [] ???
An alternative is to use the direct circle equations. The angles can be found using the arc tangent.
from sympy import symbols, sin, cos, atan2, Eq, solve
x1, y1 = symbols('x1 y1', real=True)
r = 10
xm1, ym1 = 5, 0 # center of the first circle
xm2, ym2 = 0, 0 # center of the second circle
test = [Eq((x1 - xm1) ** 2 + (y1 - ym1) ** 2, r ** 2),
Eq((x1 - xm2) ** 2 + (y1 - ym2) ** 2, r ** 2), ]
ans = solve(test, (x1, y1))
for xa, ya in ans:
theta1 = atan2(xa - xm1, ya - ym1)
theta2 = atan2(xa - xm2, ya - ym2)
print(f'theta1: {theta1}, (in radians: {theta1.evalf():.6f})')
print(f'theta2: {theta2}, (in radians: {theta2.evalf():.6f})')
This finds two solutions for (x1, y1):
[(5/2, -5*sqrt(15)/2),
(5/2, 5*sqrt(15)/2)]
And corresponding angles:
theta1: -pi + atan(sqrt(15)/15), (in radians: -2.888912)
theta2: pi - atan(sqrt(15)/15), (in radians: 2.888912)
theta1: -atan(sqrt(15)/15), (in radians: -0.252680)
theta2: atan(sqrt(15)/15), (in radians: 0.252680)
Solution 2:
The answer from JohanC is good but I thought I'd add that you can actually get the answer in fully general terms for arbitarry symbols if you use Cartesian coordinates rather than polar:
In [88]: x1, x2, y1, y2, xi, yi, r = symbols('x1, x2, y1, y2, xi, yi, r', real=True)
In [89]: eq1 = Eq((x1 - xi)**2 + (y1 - yi)**2, r**2)
In [90]: eq2 = Eq((x2 - xi)**2 + (y2 - yi)**2, r**2)
In [91]: solve([eq1, eq2], [xi, yi])
Out[91]:
⎡⎛ ⎛ ____________________________________________
⎢⎜ ⎜ ╱ ⎛ 2 2 2 2
⎢⎜ 2 2 2 2 ⎜y₁ y₂ ╲╱ -⎝x₁ - 2⋅x₁⋅x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂
⎢⎜x₁ - x₂ + y₁ - y₂ + (-2⋅y₁ + 2⋅y₂)⋅⎜── + ── - ───────────────────────────────────────────────
⎢⎜ ⎜2 2 ⎛ 2
⎢⎜ ⎝ 2⋅⎝x₁ - 2⋅x₁⋅
⎢⎜─────────────────────────────────────────────────────────────────────────────────────────────────
⎢⎜ 2⋅(x₁ - x₂)
⎣⎝
_______________________________________________________ ⎞
⎞ ⎛ 2 2 2 2 2⎞ ⎟
⎠⋅⎝- 4⋅r + x₁ - 2⋅x₁⋅x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠ ⋅(x₁ - x₂)⎟
─────────────────────────────────────────────────────────────────⎟ __________________
2 2 2⎞ ⎟ ╱ ⎛ 2
x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠ ⎠ y₁ y₂ ╲╱ -⎝x₁ - 2⋅x₁⋅x₂ +
──────────────────────────────────────────────────────────────────, ── + ── - ─────────────────────
2 2
⎞ ⎛
⎟ ⎜
⎟ ⎜ 2
_________________________________________________________________________________ ⎟ ⎜x₁
2 2 2⎞ ⎛ 2 2 2 2 2⎞ ⎟ ⎜
x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠⋅⎝- 4⋅r + x₁ - 2⋅x₁⋅x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠ ⋅(x₁ - x₂)⎟ ⎜
───────────────────────────────────────────────────────────────────────────────────────────⎟, ⎜────
⎛ 2 2 2 2⎞ ⎟ ⎜
2⋅⎝x₁ - 2⋅x₁⋅x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠ ⎠ ⎝
⎛ __________________________________________________
⎜ ╱ ⎛ 2 2 2 2⎞ ⎛
2 2 2 ⎜y₁ y₂ ╲╱ -⎝x₁ - 2⋅x₁⋅x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠⋅⎝- 4
- x₂ + y₁ - y₂ + (-2⋅y₁ + 2⋅y₂)⋅⎜── + ── + ─────────────────────────────────────────────────────
⎜2 2 ⎛ 2
⎝ 2⋅⎝x₁ - 2⋅x₁⋅x₂ + x
───────────────────────────────────────────────────────────────────────────────────────────────────
2⋅(x₁ - x₂)
_________________________________________________ ⎞
2 2 2 2 2⎞ ⎟
⋅r + x₁ - 2⋅x₁⋅x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠ ⋅(x₁ - x₂)⎟
───────────────────────────────────────────────────────────⎟ ________________________
2 2 2⎞ ⎟ ╱ ⎛ 2 2
₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠ ⎠ y₁ y₂ ╲╱ -⎝x₁ - 2⋅x₁⋅x₂ + x₂ +
────────────────────────────────────────────────────────────, ── + ── + ───────────────────────────
2 2
⎞⎤
⎟⎥
⎟⎥
___________________________________________________________________________ ⎟⎥
2 2⎞ ⎛ 2 2 2 2 2⎞ ⎟⎥
y₁ - 2⋅y₁⋅y₂ + y₂ ⎠⋅⎝- 4⋅r + x₁ - 2⋅x₁⋅x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠ ⋅(x₁ - x₂)⎟⎥
─────────────────────────────────────────────────────────────────────────────────────⎟⎥
⎛ 2 2 2 2⎞ ⎟⎥
2⋅⎝x₁ - 2⋅x₁⋅x₂ + x₂ + y₁ - 2⋅y₁⋅y₂ + y₂ ⎠ ⎠⎦