I have a system of $n$ rational functions of degree $d$, $f_1, \ldots, f_n$ in variables $x_1, \ldots, x_n$, i.e. $\forall i, f_i = p_i/q_i$ where $p_i, q_i$ are multivariate polynomials of degree at most $d$. I would like to find a good approximations to the roots of this system. When I searched for this, the most popular method seems to be the multidimensional Newton iteration. However, it is known that Newton iteration has bad global convergence. I have two questions:
Is there a method that has super-linear convergence rate, and has good global convergence? I'd be interested even in methods even in the univariate setting.
Does the fact that my system is made up of rational functions help/hinder Newton iteration? I understand that $f_i$'s have the same roots as $p_i$, but since the multidimensional version needs the evaluation of a Jacobian, does it cause problems?
I'm interested only in the theory and don't care about actual performance/implementation.
Roots of rational polynomials
$$f(x)=\frac{p(x)}{q(x)} $$ Any root of $p(x)$ is also a root of of $f(x)$ independent of $q(x)$. A root of $q$ is denoted a pole, however not relevant for this question.
Univariant
$$f(x)=\frac{p(x)}{q(x)}=\frac{a_nx^n....a_0x^0}{b_mx^m....+1} $$
Regarding an univariant polynomial, do not use gradient or stochastic methods. Build the companion matrix (this is a matrix such that the characteristic polynom is $p(x)$) and find the eigenvalues. The pure real eigenvalues are the roots of $p(x)$. A polynomial of degree d has d roots in $\mathbb{C}$ and at most d roots in $\mathbb{R}$, Fundamental theorem of algebra. If your are looking at domain on $\mathbb{R}$, [a,b], you can simply keep all the eigenvalues which fullfill: $Im(\lambda)= 0 ∧ \lambda>a ∧ \lambda<b $. Those are then the roots of your polynomials in said domain. Make sure you don't return more than d roots.
Multivariant
$$f(x_1,...,x_i)=\frac{p(x_1,...,x_i)}{q(x_1,...,x_i)}$$
As you can have a multiplicity of roots, global convergence to a single root is likley by Newton or Levenberg-Marquardt, depending on the initial guess. Make sure you don't compute the gradients nummerically, as they can be easily computed analytically.
CMA-ES is very likely to globaly converge to one root, however not recommended for polynomials. I just mentions this method because it is my favorite.
If you want to find all roots of $p(x_1,...,x_i)$ you can add a radial basis function to $p(x_1,...,x_i)$ for every root you find. Then your globally convergent method converges to the next root, as the radial basis function add nonzero value to the root, thereby it is no longer a root. Or you start a few 100 newton iterations with random initial guess and then you are likely to get most of your roots.
However, you can also find a root in a one dimensional subspace, $x_1$ for example, (treating all other variables as constants) similar as the univariant case. If the root is complex, indicating a local minima you should continue in an orthogonal subspace, $x_2$ for example.