Numerically Calculating the solution of very complex equations

Question

I wanted to confirm a question of my own, and I figured out if there is a solution of the following equations such that every variable is real and $x,y\ge 0$, my question could be partially verified.

However, after two hours of running Mathematica, I thought that this would not give me the solution, and aborted it. Is there a good way to calculate a numerical solution of an equation as complex as this? I tried using other kinds of equations which the solution is known, but it didn't work either.

I really want to calculate a solution, or prove that this has a solution. Any ideas?

Edit. I add a Mathematica code for solving the equations:

t0813[x_, y_, z_, a_, b_, c_, d_, e_, f_, g_, h_, i_, j_] := 
a x^3 + b y^3 + c z^3 + d x^2 y + e x y^2 + f y ^2 z + g y z^2 + 
h z^2 x + i z x^2 + j x y z 

NSolve[{t0813[1, 1, 1, a, b, c, d, e, f, g, h, i, j] == 0, 
  t0813[1, x, y, a, b, c, d, e, f, g, h, i, j] == 0, 
  t0813[1, y, x, a, b, c, d, e, f, g, h, i, j] == 0, 
  t0813[x, 1, y, a, b, c, d, e, f, g, h, i, j] == 0, 
  t0813[y, 1, x, a, b, c, d, e, f, g, h, i, j] == 0, 
  t0813[x, y, 1, a, b, c, d, e, f, g, h, i, j] == 0, 
  t0813[y, x, 1, a, b, c, d, e, f, g, h, i, j] == 0, 
  a^2 + b^2 + c^2 == 1, 
  2 a d + 2 b e + 2 b f + 2 c g + 2 c h + 2 a i == 0.01 - 700/102, 
  2 b d + d^2 + 2 a e + e^2 + 2 c f + f^2 + 2 b g + g^2 + 2 a h + 
  h^2 + 2 c i + i^2 == 1124/102, 
  2 e f + 2 g h + 2 d i + 2 a j + 2 b j + 2 c j == 738/102, 
  2 a b + 2 a c + 2 b c + 2 d e + 2 f g + 2 h i == -564/102}, {a, b, 
  c, d, e, f, g, h, i, j, x, y}]

Answer 1

Your system has the following form $$ \begin{cases} \mathbf F_7(x,y,\mathbf p) = 0\\ \mathbf G_5(\mathbf p) = 0\\ \end{cases} $$ with $\mathbf p = (a,b,c,d,e,f,g,h,i,j), \mathbf r = (x, y, z)$. $\mathbf F_7$ are seven equations linear in $\mathbf p$ and $\mathbf G_5$ are nonlinear equations in $\mathbf p$. We have $$ T(\mathbf r,\mathbf p) = ax^3 + by^3 + cz^3 + dx^2 y + e xy^2 + f y^2 z + g y z^2 + h z^2 x + i z x^2 + j x y z, $$ thus $$ \mathbf F_7(x,y, \mathbf p) = \underbrace{\begin{pmatrix} 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 \\ 1 & x^3 & y^3 & x & x^2 & x^2 y & x y^2 & y^2 & y & x y \\ 1 & y^3 & x^3 & y & y^2 & x y^2 & x^2 y & x^2 & x & x y \\ x^3 & 1 & y^3 & x^2 & x & y & y^2 & x y^2 & x^2 y & x y \\ y^3 & 1 & x^3 & y^2 & y & x & x^2 & x^2 y & x y^2 & x y \\ x^3 & y^3 & 1 & x^2 y & x y^2 & y^2 & y & x & x^2 & x y \\ y^3 & x^3 & 1 & x y^2 & x^2 y & x^2 & x & y & y^2 & x y \end{pmatrix}}_{\mathbf A}\mathbf{p}. $$ Let's consider case when $0 < x < y, x \neq 1, y \neq 1$. Let $\mathbf B$ be the following submatrix of $\mathbf A$: $$ \mathbf B = \begin{pmatrix} 1 & 1 & 1 & 1 & 1 & 1 & 1 \\ x & x^2 & x^2 y & x y^2 & y^2 & y & x y \\ y & y^2 & x y^2 & x^2 y & x^2 & x & x y \\ x^2 & x & y & y^2 & x y^2 & x^2 y & x y \\ y^2 & y & x & x^2 & x^2 y & x y^2 & x y \\ x^2 y & x y^2 & y^2 & y & x & x^2 & x y \\ x y^2 & x^2 y & x^2 & x & y & y^2 & x y \end{pmatrix}\\ \operatorname{det} \mathbf B = (x-1)^3 (y-1)^3 (y-x)^3 (x+y+1)^2 (x y+x+y)^2 \left((x+1) y^2+(x-6) x y+x (x+1)+y\right) \neq 0. $$ Now we can express $d,e,f,g,h,i,j$ in terms of $x,y,a,b,c$ $$ \begin{pmatrix} d\\e\\f\\g\\h\\i\\j \end{pmatrix} = \frac{1}{\operatorname{det}(\mathbf B)}\mathbf{M}\begin{pmatrix} a\\b\\c \end{pmatrix} \qquad \mathbf{M} \equiv -\operatorname{det}(\mathbf B)\mathbf{B}^{-1} \begin{pmatrix} 1 & 1 & 1\\ 1 & x^3 & y^3 \\ 1 & y^3 & x^3 \\ x^3 & 1 & y^3 \\ y^3 & 1 & x^3 \\ x^3 & y^3 & 1 \\ y^3 & x^3 & 1 \end{pmatrix}. $$ or in terms of $(\xi, \eta, \zeta) = \frac{(a,b,c)}{\operatorname{det}(\mathbf B)}$ $$ \mathbf p = \mathbf {Rw} \qquad \mathbf R(x,y) \equiv \begin{pmatrix} \operatorname{det}(\mathbf B) \mathbf I_3\\ \mathbf{M} \end{pmatrix} \quad \mathbf{w} \equiv \begin{pmatrix} \xi\\ \eta\\ \zeta \end{pmatrix} $$

Now we've completely eliminated equations $\mathbf F_7$ and ended up with five equations in five unknowns $x,y,\xi,\eta,\zeta$. Since every equation in $\mathbf G_5(\mathbf p) = 0$ can be written as $$ \mathbf p^\top \mathbf Q_i \mathbf p = \alpha_i, \quad \mathbf Q_i^\top = \mathbf Q_i \quad i = 1,2,3,4,5. $$ we can rewrite them in terms of $\mathbf w$: $$ \mathbf w^\top \underbrace{\mathbf R^\top \mathbf Q_i \mathbf R}_{\mathbf S_i} \mathbf w = \alpha_i $$ Every $\mathbf S_i$ matrix has form $$ \mathbf S_i = \begin{pmatrix} a_i & b_i & b_i\\ b_i & a_i & b_i\\ b_i & b_i & a_i \end{pmatrix} $$ I think this is due to it should stay the same when $\xi, \eta, \zeta$ are permuted, thus the equations reduce to $$ a_i(x,y) (\xi^2 + \eta^2 + \zeta^2) + b_i(x,y) (2\xi \eta + 2\xi \zeta + 2\eta \zeta) = \alpha_i, \quad i = 1,2,3,4,5 $$ Introducing $\sigma = \xi^2 + \eta^2 + \zeta^2, \tau = 2\xi \eta + 2\xi \zeta + 2\eta \zeta$ we see that $$ a_i(x,y) \sigma + b_i(x,y) \tau = \alpha_i, \quad i = 1,2,3,4,5. $$ is a linear system with five equations in two unknowns $\sigma, \tau$, thus solvable only if $$ \operatorname{rank} \begin{pmatrix} a_1(x,y) & b_1(x,y)\\ \vdots&\vdots\\ a_5(x,y) & b_5(x,y) \end{pmatrix} = \operatorname{rank} \begin{pmatrix} a_1(x,y) & b_1(x,y) & \alpha_1\\ \vdots&\vdots&\vdots\\ a_5(x,y) & b_5(x,y) & \alpha_5 \end{pmatrix} $$ Since the left matrix seems to have rank $2$, we need to solve $$\left\{ \begin{gathered} \begin{vmatrix} a_1(x,y) & b_1(x,y) & \alpha_1\\ a_2(x,y) & b_2(x,y) & \alpha_2\\ a_3(x,y) & b_3(x,y) & \alpha_3 \end{vmatrix} = 0\\ \begin{vmatrix} a_2(x,y) & b_2(x,y) & \alpha_2\\ a_3(x,y) & b_3(x,y) & \alpha_3\\ a_4(x,y) & b_4(x,y) & \alpha_4 \end{vmatrix} = 0\\ \begin{vmatrix} a_3(x,y) & b_3(x,y) & \alpha_3\\ a_4(x,y) & b_4(x,y) & \alpha_4\\ a_5(x,y) & b_5(x,y) & \alpha_5 \end{vmatrix} = 0\end{gathered}\right. $$ for $x,y$.

So either you're extremely lucky with $\alpha_i$ so the over-determined polynomial system has roots, in which case there's a solution $(\sigma, \tau)$ and the problem is to find the intersection of a sphere and a hyperboloid $$ \begin{cases} \xi^2 + \eta^2 + \zeta^2 = \sigma\\ (\xi+\eta+\zeta)^2 - (\xi^2 + \eta^2 + \zeta^2) = 2\xi\eta + 2\xi\zeta + 2\eta\zeta = \tau \end{cases}\\ \begin{cases} \xi^2 + \eta^2 + \zeta^2 = \sigma\\ \xi+\eta+\zeta = \sqrt{\tau + \sigma} \end{cases} $$ which itself may have no solutions, one solution or infinitely many solutions.

But For arbitrary $\alpha_i$ the system for $\sigma, \tau$ is inconsistent whatever $x,y$ you take.