我认为它只是缺少初始化的值

from scipy.optimize import fsolve

def equation(var, *data):
    a,b,c = data 
    x,y = var
    eq1 = x**2 + y**2 - a**2
    eq2 = (x - b)**2 + y**2 - c**2
   return [eq1, eq2]

x,y = fsolve(equation,[1,1], args=(1,1,1))
print(x,y)

从Thierry Lathuille的几何解释中获得灵感，可能真实的上不需要使用非线性求解器。
首先，方程（1）要求y^2=a^2-x^2，因此（x-b）^2+a^2-x^2=c^2，根据方程式（2）.后一个等式简化为b^2-2bx+a^2=c^2，其中xSol=（a^2+b^2-c^2）/（2b）作为x的可能解（对称性意味着x只有一个解）。利用这个解，我们可以检查a^2-xSol^2的符号。如果为负，则没有解。如果为非负，则解为（xSol，+ySol）和（xSol，-ySol），其中ySol=np·sqrt（a^2-xSol ^2）。
基于对下面代码的一些测试，上面的解决方案似乎可以工作（如果我错了，请纠正我）。

from scipy.optimize import fsolve
import numpy as np

abc = (1,1,4)

# Nonlinear optimizer

def equation(var, *data):
    a,b,c = data
    x,y = var
    eq1 = x**2 + y**2 - a**2
    eq2 = (x - b)**2 + y**2 - c**2
    return [eq1, eq2]

x,y = fsolve(equation,[1,1], args=abc)
print(x,y)

# Geometric solution

a = abc[0]; b = abc[1]; c = abc[2]
xSol = (a**2+b**2-c**2)/(2*b)
if a**2-xSol**2<0:
    print("No solution")
else:
    ySol = np.sqrt(a**2-xSol**2)
    print(xSol,ySol)
    print(xSol,-ySol)