0

I am trying to find a solution to the following system where f and g are R^2 -> R^2 functions:

f(x1,x2) = (y1,y2)
g(y1,y2) = (x1,x2)

I tried solving it using scipy.optimize.fsolve as follows:

def eqm(vars):
    x1,x2,y1,y2 = vars
    eq1 = f([x1, x2])[0] - y1
    eq2 = f([x1, x2])[1] - y2
    eq3 = g([y1, y2])[0] - x1
    eq4 = g([y1, y2])[1] - x2
    return [eq1, eq2, eq3, eq4]

fsolve(eqm, x0 = [1,0.5,1,0.5])

Although it is returning an output, it does not seem to be a correct one as it does not seem to satisfy the two conditions, and seems to vary a lot with the x0 specified. Also getting a warning: 'The iteration is not making good progress, as measured by the improvement from the last ten iterations.' I do know for a fact that a unique solution exists, which I have obtained algebraically.

Not sure what is going on and if there is a simpler way of solving it, especially using just two equations instead of splitting up into 4. Something like:



def equations(vars):
    X,Y = vars
    eq1 = f(X)-Y
    eq2 = g(Y)-X
    return [eq1, eq2]

fsolve(equations, x0 =[[1,0.5],[1,0.5]])

Suggestions on other modules e.g. sympy are also welcome!

McLovin
  • 1
  • 1
  • 2
    Welcome to SO. Can you please provide `f` and `g` such that your question contains a [minimal **reproducible** example](https://stackoverflow.com/help/minimal-reproducible-example)? – joni Apr 17 '22 at 17:44
  • Thanks, I have added. But I am afraid it is not too 'minimal' – McLovin May 04 '22 at 20:50

2 Answers2

1

First, I recommend working with numpy arrays since manipulating these is simpler than lists.

I've slighlty rewritten your code:

import scipy.optimize as opt

def f(x):
  return x
def g(x):
  return x

def func(vars):
  input = np.array(vars)
  eq1 = f(input[:2]) - input[2:]
  eq2 = g(input[2:]) - input[:2]
  return np.concatenate([eq1, eq2])

root = opt.fsolve(func, [1, 1, 0., 1.2])

print(root)
print(func(root)) # should be close to zeros

What you have should work correctly, so I believe there is something wrong with the equations you're using. If you provide those, I can try to see what may be wrong.

wizzzz1
  • 21
  • 2
  • Thanks for your comment. I have edited to include the equations. The equations themselves are being derived from two optimization problems. Also, trying to think if a different method for the optimization might yield better results. – McLovin May 04 '22 at 20:44
0

This seems to be more of a problem of numerical mathematics than Python coding. Your functions may have "ugly" behavior around the solution, may be strongly non-linear or contain singularities. We cannot help further without seeing the functions. One thing you might try is to instead solve a system

g(f(x)) - x = 0

and simplify g(f(x)) as much as possible analytically. Then calculate y = f(x) after solving the equation.

Martin Sipka
  • 56
  • 3
  • Thanks for your comment. I have included the specific functions in the last edit. The functions are basically the outcome of two optimization problems. Do you think the optimization method being used could be leading to potentially 'ugly' behavior? – McLovin May 04 '22 at 20:45