Python中快速求解与两已知圆相切的圆坐标的高效方法

你的问题本质是求两个已知圆的交点（待求圆的圆心是这两个圆的交点，这两个圆的圆心分别是(c1.x,c1.y)、(c2.x,c2.y)，半径分别为r1+r、r2+r）。用SymPy的符号求解效率极低，因为符号计算需要处理通用场景，而这个问题有直接的解析解，用数值计算实现能把单次计算耗时降到微秒级。

推导过程 #

1. 先将原方程平方消去根号：

   (c1.x - a)² + (c1.y - b)² = (r1 + r)²  --- (1)
   (c2.x - a)² + (c2.y - b)² = (r2 + r)²  --- (2)
   

2. 用(1)-(2)消去二次项，得到线性方程：
   展开后整理为：

   2*(c2.x - c1.x)*a + 2*(c2.y - c1.y)*b = (c2.x² + c2.y² - (r2 + r)²) - (c1.x² + c1.y² - (r1 + r)²)
   

   令：
   • dx = c2.x - c1.x
   • dy = c2.y - c1.y
   • d_sq = dx² + dy²
   • r_sum1 = r1 + r
   • r_sum2 = r2 + r
   • C = (c2.x² + c2.y² - r_sum2²) - (c1.x² + c1.y² - r_sum1²)
     线性方程简化为：

   2*dx*a + 2*dy*b = C
   

3. 从线性方程中解出一个变量（比如b），代入方程(1)，得到关于a的一元二次方程，求解后即可得到两个可能的(a,b)解（对应待求圆在两个已知圆的不同侧相切）。

Python实现代码 #

直接用数值计算实现上述推导，无需符号库：

import math

def find_tangent_circle_center(c1_x, c1_y, r1, c2_x, c2_y, r2, r):
    dx = c2_x - c1_x
    dy = c2_y - c1_y
    d_sq = dx**2 + dy**2
    if d_sq == 0:
        raise ValueError("两个已知圆的圆心重合，无法求解")
    
    r_sum1 = r1 + r
    r_sum2 = r2 + r
    
    # 计算线性方程的系数
    A = 2 * dx
    B = 2 * dy
    C_val = (c2_x**2 + c2_y**2 - r_sum2**2) - (c1_x**2 + c1_y**2 - r_sum1**2)
    
    # 根据系数绝对值大小选择解的变量，避免除零或精度损失
    if abs(B) > abs(A):
        # 解出a = (C_val - B*b)/A，代入方程(1)得到关于b的二次方程
        E = A**2 + B**2
        F = 2 * (B*(c1_x*A - C_val) - A**2 * c1_y)
        G = (C_val - A*c1_x)**2 + A**2 * c1_y**2 - A**2 * r_sum1**2
        
        discriminant = F**2 - 4*E*G
        if discriminant < 0:
            return []  # 无实根（两圆不相交）
        
        sqrt_d = math.sqrt(discriminant)
        b1 = (-F + sqrt_d) / (2*E)
        b2 = (-F - sqrt_d) / (2*E)
        
        a1 = (C_val - B*b1) / A
        a2 = (C_val - B*b2) / A
        
        return [(a1, b1), (a2, b2)]
    else:
        # 解出b = (C_val - A*a)/B，代入方程(1)得到关于a的二次方程
        E = A**2 + B**2
        F = 2 * (A*(c1_y*B - C_val) - B**2 * c1_x)
        G = (C_val - B*c1_y)**2 + B**2 * c1_x**2 - B**2 * r_sum1**2
        
        discriminant = F**2 - 4*E*G
        if discriminant < 0:
            return []
        
        sqrt_d = math.sqrt(discriminant)
        a1 = (-F + sqrt_d) / (2*E)
        a2 = (-F - sqrt_d) / (2*E)
        
        b1 = (C_val - A*a1) / B
        b2 = (C_val - A*a2) / B
        
        return [(a1, b1), (a2, b2)]

效率说明 #

这个实现完全基于数值运算，没有符号推导的开销，单次计算耗时在微秒级别，比SymPy的solve()快几个数量级，完全能满足你每次迭代250次的需求。

注意事项 #

• 函数返回两个解，对应待求圆在两个已知圆连线的两侧相切，你可以根据圆打包的边界条件（比如在方形区域内）选择合适的解。
• 如果判别式小于0，说明两圆没有交点（无法相切），返回空列表，这种情况需要在你的打包算法中做异常处理。

内容的提问来源于stack exchange，提问作者Kai