# Solve non-linear set of three equations using scipy

I need to solve a non-linear set of three equations using scipy.

However, I do not have any clue on which algorithm is suitable for my problem from a mathematical point of view (stability, convergence behaviour), since scipy provides a huge variety of different algorithms in the scipy.optimize module such as:

Here is my set of equations (updated version due to comments): Which algorithm is probably the best one for my problem?

• What is $t_{k-1}$? Is there some sort of initial condition for the $Q(0)$? – Bill Barth Jun 21 '15 at 14:10
• $t_{k-1}$ just symbolizes the time dependence of the temperatures $\vartheta_{foo}$ whereas $k-1$ means that this references to the value of the last calculation / discretization step resp. initial condition if we talk about the first step. I am not sure, but $Q(0)$ should be zero, since there is no heat transfered at init. – albert Jun 21 '15 at 14:29
• Yeah, I get that, but your last equation shows $t_k=1$, so you should either remove it or provide something about the size of the time steps. – Bill Barth Jun 21 '15 at 14:31
• @BillBarth: Since this is part of a self-written simulation software the step size depends on the number of discretisation elements specified by the user. In the simulation algorithm I converted the rotation of a wheel into time steps depending on the total number of discretisation points. Due to simplicity I assumed that the step width is equal to one second which is given if the wheel turns with 1 rpm and is represented by 360 discretisation segments / points. If you wish having a more physical approach I could add the appropiate units to all variables. – albert Jun 21 '15 at 14:41
• OK, replacing $t_k$ by $\Delta t$ makes this a bit clearer. – Bill Barth Jun 21 '15 at 15:00

Since your problem is small, you're probably best off trying fsolve or root. Both of these are interfaces to MINPACK and call HYBRD or HYBRJ. Since calculating a Jacobian matrix for your system shouldn't be hard (either do it by hand, or use your favorite computer algebra system, like SymPy, Sage, Maple, or Mathematica), you should supply a Jacobian matrix. The HYBR functions in MINPACK use "Powell's hybrid method", which uses Newton's method and checks if Newton steps will be descent steps by comparing against least-squares minimization (specifically, does the Newton step also decrease the sum-of-squared-residuals). If Newton steps are not descent steps, then it falls back to the gradient of the sum-of-squared-residuals.
scipy.optimize.root is a unified interface for the nonlinear solvers that you list in your third point, which are implemented in Python. On the other hand, scipy.optimize.fsolve is a wrapper around some Fortran functions (see docs). The interface is the same for both and I haven't seen any benchmarks (quick Google search), so just implement the function of the system and give them a go. Focus on correct implementation of this function and on choosing a sensible initial condition and then just try both root with its different methods and fsolve.