# Wolfram Function Repository

Instant-use add-on functions for the Wolfram Language

Function Repository Resource:

Numerically solve for isolated roots of a square system of polynomial equations using monodromy

Contributed by:
Aravind Baskar

ResourceFunction["NSolveByMonodromy"][ finds isolated roots of the polynomial equation | |

ResourceFunction["NSolveByMonodromy"][{ finds isolated roots of the system of polynomial equations |

ResourceFunction["NSolveByMonodromy"] attempts to numerically approximate many isolated roots of a square system of polynomial equations through an iterative numerical continuation algorithm using monodromy.

ResourceFunction["NSolveByMonodromy"] deals primarily with polynomial equations but can handle some analytic transcendental equations with limited success.

The system *eqns* can be of the following form:

lhs==rhs | equations |

expr | functions implicitly set to zero |

A single variable or a list of variables matching the dimension of *eqns* must be specified.

ResourceFunction["NSolveByMonodromy"] finds an initial nonsingular root using FindRoot via a random guess to kick-start the first monodromy iteration. Alternatively, an user-defined initial point may be provided through the option "InitialPoint".

ResourceFunction["NSolveByMonodromy"] takes the following options:

"MaxRoots" | Infinity | Maximum number of roots to return |

"Radius" | 100 | Radius parameter for the construction of monodromy loops |

"InitialPoint" | {} | User-defined starting guess for the first loop |

"MaxIterations" | 10 | Maximum number of iterations of monodromy loops |

"MaxSteps" | 1000 | Maximum number of steps allowed in a single path |

"MinStepSize" | 10^{-16} | Minimum size of a single step used during path tracking |

"Parallel" | False | Allow parallel computation of monodromy paths |

Numerically solve a univariate polynomial equation in *x*:

In[1]:= |

Out[2]= |

Numerically solve a multivariate system of polynomial equations:

In[3]:= |

Out[6]= |

Create a real-valued polynomial objective function *f* in three variables *x*, *y*, *z* using random coefficient parameters:

In[7]:= |

Find stationary points of the objective *f*:

In[8]:= |

Out[11]= |

Select the real valued ones:

In[12]:= |

Out[13]= |

Function values at these critical points:

In[14]:= |

Out[14]= |

Sorted critical points according to the function value:

In[15]:= |

Out[16]= |

Two of these (namely, #1 and #3) are local minima with all positive eigenvalues of the local Hessian, while the rest are saddle points:

In[17]:= |

Out[18]= |

Local minima and saddle points can be also discerned via the following ContourPlot3D visualization:

In[19]:= |

Out[20]= |

Default execution may not find all isolated roots:

In[21]:= |

Out[22]= |

All the roots may be obtained by increasing or decreasing "Radius" parameter:

In[23]:= |

Out[25]= |

"MaxRoots" can be used to compute a subset of the total roots for unwieldy polynomial systems:

In[26]:= |

Generate a random large polynomial system of dimension 5 and total degree 100,000:

In[27]:= |

Compute 50 roots:

In[28]:= |

Out[8]= |

In some cases transcendental functions can be minimized using "MaxRoots" termination criterion:

In[29]:= |

Visualize the real valued stationary points:

In[30]:= |

Out[31]= |

Direct kinematics of a Gough-Stewart parallel 6-degree of freedom platform:

In[32]:= |

Out[36]= |

Katsura-5 system from magnetism in physics:

In[37]:= |

Out[41]= |

A stationary chemical kinetics problem:

In[42]:= |

Out[46]= |

NSolveByMonodromy can find many or all of the isolated roots with or without an initial point, e.g., a bifurcation problem:

In[47]:= |

Out[51]= |

FindRoot attempts to find one root at a time based on an initial point:

In[52]:= |

Out[52]= |

NSolveByMonodromy compares favorably against NSolve in some large systems, e.g., 6R inverse position kinematics:

In[53]:= |

In[54]:= |

Out[17]= |

NSolve default execution is 4 times slower and returns duplicate instances requiring an additional filtration step to remove them:

In[55]:= |

Out[56]= |

In[57]:= |

Out[57]= |

NSolve execution in high-precision takes even longer than the default execution:

In[58]:= |

Out[58]= |

Another example, a reduced 8-dimensional system arising in economics:

In[59]:= |

Out[63]= |

NSolve default execution misses two isolated roots and returns some duplicate instances instead:

In[64]:= |

Out[65]= |

In[66]:= |

Out[66]= |

NSolve in high-precision returns the correct set of roots in a comparable time with NSolveByMonodromy:

In[67]:= |

Out[68]= |

In[69]:= |

Out[69]= |

NSolveByMonodromy does not accept a search range as input, unlike NSolve:

In[70]:= |

Select the real valued ones:

In[71]:= |

Out[71]= |

NSolve solves this system if a search range is provided:

In[72]:= |

Out[72]= |

In[73]:= |

Out[73]= |

NSolveByMonodromy can only solve square systems:

In[74]:= |

Out[74]= |

NSolve computes some witness points for such under-determined systems:

In[75]:= |

Out[75]= |

NSolveByMonodromy cannot solve non-analytic systems:

In[76]:= |

Out[76]= |

NSolve solves this system using inverse functions:

In[77]:= |

Out[77]= |

NSolveByMonodromy may fail to proceed beyond initial point in systems with special monomial structures and/or parameters:

In[78]:= |

Out[81]= |

NSolve finds the approximation of all the 20 roots for the Wilkinson's polynomial:

In[82]:= |

Out[82]= |

Static equilibria of a weight suspended by two massless cables, each of different given lengths, in a plane:

In[83]:= |

This 10-dimensional system admits 12 nonsingular isolated roots in a generic case, as found by NSolveByMonodromy:

In[84]:= |

Out[85]= |

NSolve default execution takes a long time:

In[86]:= |

Out[86]= |

The presence of singular and positive-dimensional sets impacts the results of NSolve:

In[87]:= |

Out[87]= |

Select the real valued ones in *x*, *y*:

In[88]:= |

Out[89]= |

Visualize the real equilibria, in which Configuration 1 is the only stable one with both cables under tension:

In[90]:= |

Out[90]= |

- Monodromy of z^n + b z + 1 = 0
- The Monodromy Group of an Algebraic Function
- Jan Verschelde demos
- Inria SAGA test suite
- ExampleSystems in Macaulay2
- Al-Roomi benchmarks

Wolfram Language 13.0 (December 2021) or above

- 1.0.0 – 30 June 2023

This work is licensed under a Creative Commons Attribution 4.0 International License