6
$\begingroup$

Mathematica can use Reduce to arrive at an answer in the following equation:

Clear[log1, log2, log3, n, k, s, x]
s = 1;
k = 0;
log1 = Sum[0/(1*n + 1)^s, {n, 0, k}]
log2 = Sum[1/(2*n + 1)^s - 1/(2*n + 2)^s, {n, 0, k}]
log3 = Sum[1/(3*n + 1)^s + 1/(3*n + 2)^s - 2/(3*n + 3)^s, {n, 0, k}]
Reduce[1/(E^log1)^x + 1/(E^log2)^x + 1/(E^log3)^x == 0, x]

Setting k=1 it still manages to solve it, and setting k=2 it gives no error message but it calculates forever.

Setting k=20 gives a first clue to what the mathematics behind it is when I get the error message:

PolynomialGCD::lrgexp: Exponent is out of bounds for function PolynomialGCD.

So how does Reduce work in solving:

Reduce[1/(E^log1)^x + 1/(E^log2)^x + 1/(E^log3)^x == 0, x]

?

Notice that log1, log2, log3 are rational numbers used as approximations of logarithms and not exact logarithms since then the equation would be unsolvable.

My ultimate goal is to understand roots and polynomials of this form completely, but I don't know where to start really. The word Galois theory comes to mind but I don't know what that is either.

Edit 19.11.2014: For integer values of log3 and letting log1 and log2 be constant, the following program;

Clear[log1, log2, log3, n, k, s, x]
 s = 1;
k = 0;
log1 = 0;
log2 = 1;
Do[
 log3 = q;
 Print[{"value of q", 
   q, Reduce[1/(E^log1)^x + 1/(E^log2)^x + 1/(E^log3)^x == 0, x]}], {q,
   3, 32}]

seems give the largest exponent of the roots equal to this sequence in the OEIS:

https://oeis.org/A061800

Some theory might also be given by this MathWorld page:

http://mathworld.wolfram.com/PolynomialRoots.html

Edit 20.11.2014:

Clear[nn, a, x, lev, sol]
nn = 71
a = Range[nn]*0;
Monitor[Do[
  sol = Flatten[{ToRules[
      Simplify[Reduce[1/(E^0)^x + 1/(E^1)^x + 1/(E^n)^x == 0, x], 
        C[1] \[Element] Integers] /. {C[1] -> 0}]}];
  lev = Level[Exp[x /. sol[[Length[sol]]]], -1];
  a[[n]] = lev[[Length[lev] - 1]], {n, 3, nn}], n]
a

{0, 0, 3, 4, 3, 6, 7, 6, 9, 10, 9, 12, 13, 12, 15, 16, 15, 18, 19, \

18, 21, 22, 21, 24, 25, 24, 27, 28, 27, 30, 31, 30, 33, 34, 33, 36, \

37, 36, 39, 40, 39, 42, 43, 42, 45, 46, 45, 48, 49, 48, 51, 52, 51, \

54, 55, 54, 57, 58, 57, 60, 61, 60, 63, 64, 63, 66, 67, 66, 69, 70, \

69}

Edit 21.11.2014:

Clear[nn, a, x, lev, sol]
nn = 70
a = Range[nn]*0;
Monitor[Do[
  sol = Flatten[{ToRules[
      Simplify[Reduce[1/(E^0)^x + 1/(E^2)^x + 1/(E^n)^x == 0, x], 
        C[1] \[Element] Integers] /. {C[1] -> 0}]}];
  lev = Level[
    Exp[If[Mod[n, 2] == 0, 2 x, x] /. sol[[Length[sol] - 1]]], {-1}];
  a[[n]] = lev[[Length[lev] - 1]];, {n, 3, nn}], n]
a

{0, 0, 2, ( 2 I)/3, 4, 2, 4, 3, 8, 2, 10, 5, 10, 6, 14, 5, 16, 8, 16, 9, 20, 8, \

22, 11, 22, 12, 26, 11, 28, 14, 28, 15, 32, 14, 34, 17, 34, 18, 38, \

17, 40, 20, 40, 21, 44, 20, 46, 23, 46, 24, 50, 23, 52, 26, 52, 27, \

56, 26, 58, 29, 58, 30, 62, 29, 64, 32, 64, 33, 68, 32}

Table[-1 + n - If[Mod[n, 2] == 0, n/2, 0] - 
  If[Mod[n, 3] == 1, 2, 0], {n, 1, nn}]

{-2, 0, 2, -1, 4, 2, 4, 3, 8, 2, 10, 5, 10, 6, 14, 5, 16, 8, 16, 9, \

20, 8, 22, 11, 22, 12, 26, 11, 28, 14, 28, 15, 32, 14, 34, 17, 34, \

18, 38, 17, 40, 20, 40, 21, 44, 20, 46, 23, 46, 24, 50, 23, 52, 26, \

52, 27, 56, 26, 58, 29, 58, 30, 62, 29, 64, 32, 64, 33, 68, 32}

Clear[nn, a, x, lev, sol]
nn = 80
a = Range[nn]*0;
Monitor[Do[
  sol = Flatten[{ToRules[
      Simplify[Reduce[1/(E^0)^x + 1/(E^3)^x + 1/(E^n)^x == 0, x], 
        C[1] \[Element] Integers] /. {C[1] -> 0}]}];
  lev = Level[
    Exp[If[Mod[n, 3] == 0, 3 x, x] /. sol[[Length[sol] - 1]]], {-1}];
  a[[n]] = lev[[Length[lev] - 1]];, {n, 3, nn}], n]
a

{0, 0, x, 3, 4, -((2 I)/ 3), 6, 7, 2, 9, 10, 3, 12, 13, 2, 15, 16, 5, 18, 19, 6, 21, 22, 5, \

24, 25, 8, 27, 28, 9, 30, 31, 8, 33, 34, 11, 36, 37, 12, 39, 40, 11, \

42, 43, 14, 45, 46, 15, 48, 49, 14, 51, 52, 17, 54, 55, 18, 57, 58, \

17, 60, 61, 20, 63, 64, 21, 66, 67, 20, 69, 70, 23, 72, 73, 24, 75, \

76, 23, 78, 79}

Table[-1 + n - If[Mod[n, 3] == 0, n - n/3, 0] - 
  If[Mod[n, 9] == 6, 2, 0], {n, 1, nn}]

{0, 1, 0, 3, 4, -1, 6, 7, 2, 9, 10, 3, 12, 13, 2, 15, 16, 5, 18, 19, \

6, 21, 22, 5, 24, 25, 8, 27, 28, 9, 30, 31, 8, 33, 34, 11, 36, 37, \

12, 39, 40, 11, 42, 43, 14, 45, 46, 15, 48, 49, 14, 51, 52, 17, 54, \

55, 18, 57, 58, 17, 60, 61, 20, 63, 64, 21, 66, 67, 20, 69, 70, 23, \

72, 73, 24, 75, 76, 23, 78, 79}

Edit 19.7.2021:

There is simple case that appears to follow a clear rule described in terms of the von Mangoldt function:

(*start*)
k = 5;
$MaxRootDegree = 1000
TableForm[
 A = Table[
   s /. Last[Solve[(Sum[(-1)^n*Exp[m/n]^s, {n, 1, k}]) == 0, s]], {m, 
    1, k}]]
(*end*)

has the same output as:

(*start*)
m = 5;
n = Times @@ Exp[MangoldtLambda[Range[m]]];
l = Length[Divisors[n]];
d = Take[Divisors[n], {l - m + 1, l}];
p = d - First[d];
$MaxRootDegree = Max[p] + 1;
TableForm[
 Table[Reverse[d][[v]]*
   Log[Root[Sum[(-1)^(j + 1)*#1^p[[j]], {j, 1, m}] &, Max[p]]], {v, 1,
    m}]]
(*end*)

Output:

{
 {60 Log[Root[1 - #1^3 + #1^8 - #1^18 + #1^48 &, 48]]},
 {30 Log[Root[1 - #1^3 + #1^8 - #1^18 + #1^48 &, 48]]},
 {20 Log[Root[1 - #1^3 + #1^8 - #1^18 + #1^48 &, 48]]},
 {15 Log[Root[1 - #1^3 + #1^8 - #1^18 + #1^48 &, 48]]},
 {12 Log[Root[1 - #1^3 + #1^8 - #1^18 + #1^48 &, 48]]}
}

Also found this link: https://reference.wolfram.com/language/tutorial/SomeNotesOnInternalImplementation.html
at:
https://scicomp.stackexchange.com/a/26473/1899

Improve this question
$\endgroup$
5
  • 1
    $\begingroup$ Related: math.stackexchange.com/questions/979420/… $\endgroup$
    – Mats Granvik
    Oct 19, 2014 at 7:22
  • $\begingroup$ NSolve instead of Reduce will work, but it looks like you want an exact answer? $\endgroup$
    – user1722
    Oct 19, 2014 at 17:11
  • 7
    $\begingroup$ (1) Reduce is looking to create a "polynomial" is some exponential of x. It will be represented internally in a dense manner at least for some preprocessing (involving polynomial gcds). When the exponents of those numerators gets large this polynomial will be too big and PolynomialGCD code will give up with that message. $\endgroup$
    – Daniel Lichtblau
    Oct 19, 2014 at 19:00
  • 4
    $\begingroup$ (2) "Galois theory" is two words. $\endgroup$
    – Daniel Lichtblau
    Oct 19, 2014 at 19:01
  • 5
    $\begingroup$ (3) I think this is a reasonable question even if I don't have a particularly good answer. Granted, as worded it might require knowledge of Solve internals. A modest reinterpretation is "Can anyone explain why this is apparently out of range for solving?" and that seems like a fair thing to request (and a reasonable variant to address). $\endgroup$
    – Daniel Lichtblau
    Oct 19, 2014 at 19:04

1 Answer 1

Reset to default
2
$\begingroup$

Edit 1.11.2021: Jeffrey Stopple answered this at Mathoverflow on July 26 2021: https://mathoverflow.net/a/398389/25104

If one simplifies the rational approximations of the logarithms one seems to be able to spot the relation and pattern as given below "(2) alternating polynomial" in this program:

(*start*)
Clear[number, log1, log2, log3, log4, s, x, X, h, n, k];
number = 10;
log1 = 0;
log2 = Round[Log[2]*number]/number;
log3 = Round[Log[3]*number]/number;
log4 = Round[Log[4]*number]/number;
Clear[s];
x = s /. FindRoot[
   1 - 1/(E^(log2))^s + 1/(E^(log3))^s - 1/(E^(log4))^s == 0, {s, 
    1 - I*20}, WorkingPrecision -> 100]
1 - 1/(E^(Log[2]))^s + 1/(E^(Log[3]))^s - 1/(E^(Log[4]))^s
s = x;
"(1) alternating series"
1 - 1/(E^(Log[2]))^s + 1/(E^(Log[3]))^s - 1/(E^(Log[4]))^s
1 - 1/(E^(log2))^s + 1/(E^(log3))^s - 1/(E^(log4))^s
Clear[x]
1 - x^Round[Log[2]*number] + x^Round[Log[3]*number] - 
 x^Round[Log[4]*number]
"(2) alternating polynomial"
X = Exp[-s/number]
1 - X^Round[Log[2]*number] + X^Round[Log[3]*number] - 
 X^Round[Log[4]*number]
(*end*)
h = {log1, log2, log3, log4}*number
"(3) Characteristic polynomial"
Expand[CharacteristicPolynomial[
    Table[Table[(-I)^(n + k - 0)*x^((h[[n]] + 1)/2)*
       x^((h[[k]] + 1)/2), {k, 1, 4}], {n, 1, 4}], x]/x^4 - 1]
(*end*)

The more interesting case is to use the identity:

1/a^(b + I*c) = 1/a^b*(Cos[c*Log[1/a]] + I*Sin[c*Log[1/a]])

as in this program:

(*start*)
Clear[number, log1, log2, log3, log4, s, x, X, h, n, k, a, b, c];
number = 7;
log1 = 0
log2 = Round[Log[2]*number]/number
log3 = Round[Log[3]*number]/number
log4 = Round[Log[4]*number]/number
Reduce[1 - 1/(E^(log2))^b*(Cos[c*(-log2)] + I*Sin[c*(-log2)]) + 
   1/(E^(log3))^b*(Cos[c*(-log3)] + I*Sin[c*(-log3)]) - 
   1/(E^(log4))^b*(Cos[c*(-log4)] + I*Sin[c*(-log4)]) == 0, {b, c}]
(*end*)

The output is long. The variable to vary is number = 7.

Clear[m, n, k, s]
TableForm[
 Table[Expand[
   FullSimplify[
    Det[A = Table[
        Table[(-1)^
           m*((-I)^(k + n)/(k^(s/2)*n^(s/2)) - If[n == k, 1, 0]), {k, 
          1, m}], {n, 1, m}]] - 1]], {m, 1, 8}]]
TableForm[A]
s = I + 1/2;
TableForm[N[A]]

The general form for the zeros of the truncated Dirichlet Eta that Mathematica gives in unfactored form are:

(*start*)
Clear[number, log1, log2, log3, log4, s, x, X, h, n, k];
$MaxRootDegree = 2000
number = 100;
integer = 0;
m = 200;
s = Table[
   number*(2*Pi*I*integer - 
      Log[Root[Sum[(-1)^(n + 1)*#1^Round[Log[n]*number], {n, 1, m}] &,
         k]]), {k, 1, Round[Log[m]*number]}];
Block[{$MaxExtraPrecision = 100}, 
 N[Sum[(-1)^(n + 1)*1/(E^(Round[Log[n]*number]/number))^s, {n, 1, m}]]]
N[Round[%, 10^-2]]
N[s]
Sort[Re[%]]
ListLinePlot[%]
(*end*)

Approximate Dirichlet eta function zeros

(*start*)
Clear[number, log1, log2, log3, log4, s, x, X, h, n, k];
$MaxRootDegree = 2000
number = 100;
integer = 0;
m = 2000;
s = Table[
   number*(2*Pi*I*integer - 
      Log[Root[Sum[(-1)^(n + 1)*#1^Round[Log[n]*number], {n, 1, m}] &,
         k]]), {k, 50, 50 + 5 + Round[Log[m]*number]*0}];
N[%, 12]
(*end*)

number*(2*Pi*I*integer - 
   Log[Root[Sum[(-1)^(n + 1)*#1^Round[Log[n]*number], {n, 1, m}] &, 
     k]])

Checking the expression for the zeros of the truncated Dirichlet eta function by feeding it into the truncated Dirichlet eta function:

(*later 9 8 2021*)(*start*)
Clear[number, log1, log2, log3, log4, s, x, X, h, n, k];
number = 120;
integer = 0;
m = 100;
$MaxRootDegree = Round[Log[m]*number] + 10
$MaxExtraPrecision = Round[Log[m]*number] + 50
Sum[(-1)^(n + 1)*1/(E^(Log[n]))^s, {n, 1, m}]
s = ParallelTable[
   number*(If[k == 1, Pi*I, 0] + 2*Pi*I*integer - 
      Log[If[k == 1, -1, 1]*
        Root[Sum[(-1)^(n + 1)*#1^Round[Log[n]*number], {n, 1, m}] &, 
         k]]), {k, 1, Round[Log[m]*number]}];
Block[{$MaxExtraPrecision = 300}, 
 N[Sum[(-1)^(n + 1)*1/(E^(Round[Log[n]*number]/number))^s, {n, 1, m}],
   20]]
N[Round[%]]
N[s]
Sort[Re[%]]
ListLinePlot[%]
(*end*)

Edit: 1.11.2021

With the program:

(*start*)
(*Mathematica 8*)
(*Dirichlet eta function*)
Clear[m, s, n, k];
n = 100;(*Increase n=100 for better precision*)
$MaxRootDegree = Round[Log[n]*n] + 10;
Round[Log[n]*n]
Monitor[s = 
   Table[n*(-Log[
        Root[Sum[(-1)^(k + 1)*#1^Round[Log[k]*n], {k, 1, n}] &, 
         m]]), {m, Round[Log[n]*n/2], Round[Log[n]*n]}];, n]
"Plot of the zeros in the complex plane:"
ListPlot[Table[{Re[s[[n]]], Im[s[[n]]]}, {n, 1, Length[s]}]]
(*end*)

One gets the picture of the -n*Log[polynomial roots] of the Dirichlet eta function:

n equal 100 Dirichlet eta zeros plot in the complex plane

where there is a tendency for real parts of the Dirichlet eta zeros to cluster around $1/2$ and $1$.

Setting n=300 it becomes even clearer:

n equal 300 Dirichlet eta zeros plot in the complex plane

For what values of c does the following sum converge?:

Clear[n, k, x]
nn = 10^3;
c = 10;
x = N[Exp[-ZetaZero[1]/c], 100]
Plot[Re[Sum[(-1)^(k + 1)*x^(Log[k]*c), {k, 1, n}]], {n, 1, nn}]
Plot[Im[Sum[(-1)^(k + 1)*x^(Log[k]*c), {k, 1, n}]], {n, 1, nn}]
Sum[(-1)^(k + 1)*x^(Log[k]*c), {k, 1, Infinity}]

convergence plot

0.*10^-98 + 0.*10^-98 I

Improve this answer
$\endgroup$
3

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.