Filling coordinates from distance with pyomo solver

Question

I have a dataset with the distance between points and some coordinates. I want to "fill" probables coordinates thanks to the distances I know.

For this I am trying to solve the following optimisation problem :

Minimize the sum on all the edges of (distance_known(ij)²- ((xj-xi)²+(yj-yi)²)).

I wrote the following script in python : (I use an example with a 6-nodes graph) I seems to not work as I want by returning WARNING: Loading a SolverResults object with a warning status into model.name="FillingCoordinates"; - termination condition: infeasible - message from solver: Ipopt 3.11.1\x3a Converged to a locally infeasible point. Problem may be infeasible. I do not understand why ? May someone help me with this ?

L = np.array([[1,2,3],[1,3,np.sqrt(13)],[1,4,2],[1,5,np.sqrt(8)],[1,6,5],[2,3,2],[3,4,3],[4,5,2],[5,6,np.sqrt(13)]])
P1=[1,0,0]
P2=[4,0,2]
P3=[3,-3,2]
P4=[4,0,2]
P5=[5,2,2]
P6=[6,5,0]
Constr = np.array([P1,P2,P3,P4,P6])

model = pyo.ConcreteModel()

model.name = "FillingCoordinates"
model.lines = pyo.Set(initialize=(i for i in range(len(L))))


model.nodes = pyo.Set(initialize=(i+1 for i in range(6)))

model.x = pyo.Var(model.nodes, initialize = 0, domain=pyo.Integers)
model.y = pyo.Var(model.nodes, initialize = 0, domain=pyo.Integers)



model.constraint_x = pyo.ConstraintList()
model.constraint_y = pyo.ConstraintList()

for c in Constr:
    model.constraint_x.add(model.x[c[0]]==c[1])
    model.constraint_y.add(model.x[c[0]] == c[2])



for node in model.nodes:
    model.constraint_y.add(model.y[node]>=-1)
    model.constraint_y.add(model.y[node]<=3)
    model.constraint_x.add(model.x[node] >= -4)
    model.constraint_x.add(model.x[node] <= 6)



def func_objective(model):
    objective_expr = sum([L[line][2]**2
                          - ((model.x[L[line][0]]-model.x[L[line][1]])**2
                          + (model.y[L[line][0]]-model.y[L[line][1]])**2)
                          for line in model.lines])
    
    return objective_expr



model.objective = pyo.Objective(rule=func_objective,sense=pyo.minimize)

solver = pyo.SolverFactory('ipopt')



solver.solve(model, tee=True)