Problem plotting finite sum involving HurwitzZeta function

Question

I'm having problems regarding to the following code:

Z1D[\[Beta]_, \[Xi]_, \[Theta]_] = (1/\[Xi] (2/\[Beta]^2 - 
   1) + ((\[Theta] Exp[-\[Beta]])/2 - 1/2) + 
Sum[((-\[Beta])^n \[Xi]^(n/2))/n! HurwitzZeta[-(n/2), 
   1 + 1/\[Xi]], {n, 1, 80}]);

Cv0[\[Beta]_, \[Xi]_] = \[Beta]^2 D[
    Log[Z1D[\[Beta], \[Xi], 1]], {\[Beta], 2}];

CSS[\[Beta]_, a_] = 
  Cv0[\[Beta], 18] + 
   a/2 \[Beta]^2 D[\[Beta]^2/Z1D[\[Beta], 18, 1] D[
        Z1D[\[Beta], 18, 1], {\[Beta], 2}] - 
      Log[Z1D[\[Beta], 18, 1]]^2, {\[Beta], 2}];

Quiet[Plot[{CSS[\[Beta], 0], CSS[\[Beta], 5/100], CSS[\[Beta], 1/10], 
   CSS[\[Beta], 5/10]}, {\[Beta], 0, 5}, WorkingPrecision -> 50, 
  PlotPoints -> 100, MaxRecursion -> 1, PlotLegends -> "Expressions"]]

Obviously, there is some issue regarding to the HurwitzZeta function, but I cannot understand what it is. Can someone give some advice?

As a guide, my mathematical formalism ensures that the functions which I'm plotting are real-valued and positive for beta>0.

Thanks!.

Answer 1

I'm not sure if the definitions you posted are different from the ones used in your plot, or you have lingering definitions that are interfering, but the plot I get from your code is slightly different than the one you show.

Z1D[\[Beta]_, \[Xi]_, \[Theta]_] = (1/\[Xi] (2/\[Beta]^2 - 
       1) + ((\[Theta] Exp[-\[Beta]])/2 - 1/2) + 
    Sum[((-\[Beta])^n \[Xi]^(n/2))/n! HurwitzZeta[-(n/2), 
       1 + 1/\[Xi]], {n, 1, 80}]);

Cv0[\[Beta]_, \[Xi]_] = \[Beta]^2 D[
    Log[Z1D[\[Beta], \[Xi], 1]], {\[Beta], 2}];

CSS[\[Beta]_, a_] = 
  Cv0[\[Beta], 18] + 
   a/2 \[Beta]^2 D[\[Beta]^2/Z1D[\[Beta], 18, 1] D[
        Z1D[\[Beta], 18, 1], {\[Beta], 2}] - 
      Log[Z1D[\[Beta], 18, 1]]^2, {\[Beta], 2}];
Quiet@Plot[{
   CSS[\[Beta], 0],
   CSS[\[Beta], 5/100],
   CSS[\[Beta], 1/10],
   CSS[\[Beta], 5/10]
   },
  {\[Beta], 0, 5},
  WorkingPrecision -> 50,
  PlotPoints -> 100,
  MaxRecursion -> 3,
  PlotLegends -> "Expressions"
  ]

I think that setting MaxRecursion artificially low (to 1) might be responsible for the gaps you see in your plot. I have it set to 3 here. J.M.'s suggestion of Exclusions -> None also works even with MaxRecursion -> 1.

Also, I don't know enough about the Mathematics of your functions to know for sure, but it sure looks like the result is complex beyond about 4.5 depending on the exact function.

Block[{$MaxExtraPrecision = 500},
 N[{CSS[9/2, 0], CSS[9/2, 5/100], CSS[9/2, 1/10], CSS[9/2, 5/10]}, 50]
]

$-79699.364150898969209802685927039062671148814585580$

$ 8.35688067895158097516109070604051 17173041706645372\times 10^7+12519.146845612096555038077041258063228049288586 i$

$1.6721731294318251847243161680673727340875456210533\times 10^8+25038.29369122419311007615408251612645609857717 i$

$8.3640536217251618823899729477739452329445740578499\times 10^8+125191.46845612096555038077041258063228049288586 i$

I'm using Block[{$MaxExtraPrecision = 500}, ...] to tell Mathematica that it can use a precision up to 550 while evaluating the numbers in order to ensure that it returns numbers with a precision of at least 50. Without the Block, I run into a warning that says $MaxExtraPrecision = 50 reach while evaluating. Fortunately, even without the Block the numbers are basically the same, but they're nowhere near a precision of 50. In fact, their precision is around 10.

For small values of $\beta$, the precision remains high. But as the value of $\beta$ increases, it seems like the precision diminishes rapidly unless you allow Mathematica to use extra precision. This doesn't seem to affect the plot, but I thought it was interesting so I wanted to mention it. It also leads me to believe that Mathematica isn't mistaken when it says those high $\beta$ functions are complex.