# The existence of solution of a multivariate polynomial

Given a multivariate polynomial of degree 4 on $$v$$ (real) variables such that there are $$N$$ terms in the polynomial. Is there an efficient algorithm to verify if it has a solution?

Not aware of any algorithm for this problem. Can someone please comment.

I think the solution itself can be calculated by Newton Raphson algorithm (but I am unaware of its running time complexity).

No. The problem is NP-hard.

It is NP-hard to test whether a system of quadratic equations has any solution. See Reduction from vertex-cover to system of quadratic equations, Determine whether the system of equations is underdetermined or overdetermined (related: 3-SAT and Systems of Nonlinear Modular Equations).

Let $$f_1(x)=0$$, ..., $$f_k(x)=0$$ be a system of quadratic equations, where each $$f_i$$ is a multivariate quadratic polynomial on the variables. Define

$$g(x) = f_1(x)^2 + \dots + f_k(x)^2.$$

Then the single equation $$g(x)=0$$ has a solution if and only if the system of quadratic equations $$f_1(x)=0$$, ..., $$f_k(x)=0$$ has a solution. Moreover, $$g$$ has degree 4.

This reduction proves that your problem is NP-hard.

• thank you. Not much background in mathematics so a follow up: Does that mean the worst case run time of Newton Raphson or bisection method is exponential w.r.t the size of problem?
– xyz
Mar 14, 2023 at 17:14
• @xyz, I believe the bisection method can only be used for univariate polynomials. I'm not familiar with the complexity-theory worst-case running time of Newton Raphson. If the above argument is right, presumably it must imply that its worst-case running time can be exponential in some settings, but I am not able to independent verify that or explain why it might be true.
– D.W.
Mar 14, 2023 at 18:04
• Just as a comment, it might be more efficient in practice to compute the cylindrical algebraic decomposition of the polynomial, despite that algorithm being $O(2^{2^v})$ in the worst case. Mar 15, 2023 at 6:46

Mainly there are two methods, estimation of Eigenvalues of a matrix is the command method we use.

    #include <iostream>
#include <Eigen/Core>
#include <Eigen/Geometry>
#include <Eigen/Dense>
#include <unsupported/Eigen/Polynomials>
#include <cmath>
#include <vector>
#include <complex>

template<typename T>
class RootFinder {
public:
RootFinder() = default;
~RootFinder() = default;

inline std::vector<std::complex<T>> calculateEigenValuesOfMatrix(Eigen::MatrixXd A){
Eigen::EigenSolver<Eigen::MatrixXd> ES(A);
std::vector<std::complex<T>> eigen_values;
for(int i =0;i<ES.eigenvalues().rows();++i){
eigen_values.push_back(ES.eigenvalues()[i]);
}
return eigen_values;
};

inline std::vector<std::complex<T>> calculateRootsOfPoly(Eigen::VectorXd coeffs){
Eigen::PolynomialSolver<T, Eigen::Dynamic> solver;
solver.compute(coeffs);
const typename Eigen::PolynomialSolver<T, Eigen::Dynamic>::RootsType & r = solver.roots();
std::vector<std::complex<T>> roots;
for(int i =0;i<r.rows();++i)
{
roots.push_back(r[i]);
}
return roots;
};

private:

};


Here some of example how to use this method to find root of any higher order polynomial https://github.com/GPrathap/motion_planning/blob/main/motion_planning_alg/src/main.cpp