# Is it a mistake in the formula or inaccuracy in the calculations of numerical integration?

Formula for roots of trinomial $$\displaystyle z^m-az^n-1$$ with definite integration from paper Лахтинъ, “Выраженiе корней трехчленнаго алгебраическаго уравненiя посредствомъ опредѣленныхъ интеграловъ” (1890):

$$\displaystyle z_j=e^{2j\pi i/m}+\frac{1}{2\pi i}\left(e^{(2j+1)\pi i/m}\int_0^\infty log\left(1+a\frac{t^n}{1+t^m}e^{(2j+1)\pi in/m}\right)dt \\- e^{(2j-1)\pi i/m}\int_0^\infty log\left(1+a\frac{t^n}{1+t^m}e^{(2j-1)\pi in/m}\right)dt\right)$$

where natural $$m>n>0$$, $$j=0,1,...m-1$$ and $$a$$ is natural.

Formula work perfect for any $$a,m$$ and even $$n$$, but for odd $$n$$ there is an error in two roots.

Example code:

a = 7; m = 5; n = 1;
Print["\nEquation: z^", m, " - ", a, "*z^", n, " - 1 = 0\n"];
Print["Ordinary solution:"];
Print[z /. (z^m - a z^n - 1 // NSolve), "\n"];
Print["Solution with definite integration:"];
S = Table[
Exp[2 j Pi I/m] +
1/(2 Pi I) (Exp[(2 j + 1) Pi I/m]*
NIntegrate[
Log[1 + a t^n/(1 + t^m) Exp[(2 j + 1) Pi I n/m]], {t, 0,
Infinity}] -
Exp[(2 j - 1) Pi I/m]*
NIntegrate[
Log[1 + a t^n/(1 + t^m) Exp[(2 j - 1) Pi I n/m]], {t, 0,
Infinity}]),
{j, 0, m - 1}
];
Print[S];


Output:

Equation: z^5 - 7*z^1 - 1 = 0

Ordinary solution:

{-1.58871,-0.142866,0.0355442 -1.62852 I,0.0355442 +1.62852 I,1.66049}

Solution with definite integration:

NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in t near {t} = {1.58596}. NIntegrate obtained -0.71867+4.5404 I and 0.009643094017453568 for the integral and error estimates.
NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in t near {t} = {1.58596}. NIntegrate obtained -0.71867+4.5404 I and 0.009643094017453568 for the integral and error estimates.

{1.66049 +0. I,0.0355442 +1.62852 I,-1.58841+0.0000544415 I,-0.143161-0.0000544415 I,0.0355442 -1.62852 I}


If it is inaccuracy of numerical integration, then which parameters of NIntegrate need use for more correct integration?

## 1 Answer

From the code / computations below it seems that you just need higher values for MaxRecursion. Then only slow convergence messages are given ("NIntegrate::slwcon".)

a = 7; m = 5; n = 1;
Print["\nEquation: z^", m, " - ", a, "*z^", n, " - 1 = 0\n"];
Print["Ordinary solution:"];
NSolve[(z^m - a z^n - 1)]
sol = z /. NSolve[(z^m - a z^n - 1)]

(* During evaluation of In[88]:=
Equation: z^5 - 7*z^1 - 1 = 0

During evaluation of In[88]:= Ordinary solution: *)

(* {{z -> -1.58871}, {z -> -0.142866}, {z ->
0.0355442 - 1.62852 I}, {z -> 0.0355442 + 1.62852 I}, {z ->
1.66049}}

{-1.58871, -0.142866, 0.0355442 - 1.62852 I,
0.0355442 + 1.62852 I, 1.66049} *)

Print["Solution with definite integration:"]; S =
Table[Exp[2 j Pi I/m] +
1/(2 Pi I) (Exp[(2 j + 1) Pi I/m]*
NIntegrate[
Log[1 + a t^n/(1 + t^m) Exp[(2 j + 1) Pi I n/m]], {t, 0,
Infinity}, MaxRecursion -> 200] -
Exp[(2 j - 1) Pi I/m]*
NIntegrate[
Log[1 + a t^n/(1 + t^m) Exp[(2 j - 1) Pi I n/m]], {t, 0,
Infinity}, MaxRecursion -> 200]), {j, 0, m - 1}];

(*
During evaluation of In[93]:= Solution with definite integration:

During evaluation of In[93]:= NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small.

During evaluation of In[93]:= NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small. *)

S

(* {1.66049 + 0. I,
0.0355442 + 1.62852 I, -1.58871 + 7.8*10^-8 I, -0.142866 -
7.8*10^-8 I, 0.0355442 - 1.62852 I} *)


Here we see that you have all the values in your "ordinary solution" (using a certain tolerance):

Complement[S, sol,
SameTest -> (Abs[#1 - #2]/Norm[{#1, #2}, Infinity] < 10^-6 &)]

(* {} *)