Solving position trilateration with scipy? 4 equations, 3 variables

Question

I'm trying to find the coordinates of a point in cartesian space knowing only how far it is from other known points. I originally tried scipy.optimize.fsolve with 3 distance formulas and quickly learned that 4 distances are necessary because 3 spheres may overlap in up to 2 different locations.

scipy.optimize.fsolve does not seem to work when I have more equations than I have variables. With fsolve, I tried:

from scipy.optimize import fsolve


def func1(C):
    func =     [(atom14[0]  - C[0])**2 + (atom14[1]  - C[1])**2 + (atom14[2]  - C[2])**2 - Cdist1**2]
    func.append((atom106[0] - C[0])**2 + (atom106[1] - C[1])**2 + (atom106[2] - C[2])**2 - Cdist2**2)
    func.append((atom125[0] - C[0])**2 + (atom125[1] - C[1])**2 + (atom125[2] - C[2])**2 - Cdist3**2)
    func.append((atom76[0]  - C[0])**2 + (atom76[1]  - C[1])**2 + (atom76[2]  - C[2])**2 - Cdist4**2)
    return func

solve = fsolve(func1, [0, 0, 0])

where atomxxx[0:3] are known x, y, z cartesian coordinates of each of the 4 other points in space and Cdist# are the distance of the desired point from each.

I am overwhelmed by the other optimization options, many of which require different syntax I am unfamiliar with (new to this library). What can I use to solve this 3 dimensional problem with at least 4 equations?

Thanks