Error with NSum : it returns NSum::nsnum: Summand (or its derivative) f[n] is not numerical at point n=17

Question

Consider the following example (I had a lot of trouble to find a minimal working example, I think it is compactified enough now).

Omega0 = 1.
t=2
nAvg = 10.

Omegan[n_] := Omega0*Sqrt[n + 1]

f[n_] := Piecewise[{{Cos[Omegan[n]*t/2]^2*
     Abs[Exp[-nAvg/2]*Sqrt[nAvg]^n/Sqrt[Factorial[n]]], 
    0 <= n <= 20}}, 0]

NSum[f[n], {n, 1, 100}]

If you run this short script, it should return you :

NSum::nsnum: Summand (or its derivative) [...] big message [...] is not numerical at point n = 16

This problem I am facing occurs only with some specific function. It occurs with this complicated looking function I gave you but if you try simpler one the script may just work correctly.

---

My questions :

First: I would like to understand why I have this error.

Second: How to solve it ?

Extra question

Is it really more efficient to use NSum[] than N[Sum[]]. Because I have read (I don't remember where) that when Mathematica sees N[Sum[]] he understands that the sum has to be done numerically (instead of trying symbolic method THEN approximating numerically).

---

Extra infos :

I have already seen Is this a bug of NSum?

With some functions it solves the problem to add NSumTerms->number, with some other it doesn't. The thing is I would like to be able to face this problem "in general" so for this I need to understand what is happening (I read the documentation and I don't).

In short : how to do numeric summation in general with mathematica ? In my specific case I have functions that may be piecewise defined. In all generality my function can be a product/sum of piecewise functions so it is not obvious at first view to know the boundary of the sum without looking more carefully, which I would like to avoid.

Answer 1

First, I explored the meaning of the error message:

NSum[f[n], {n, 1, 100}]

NSum::nsnum: Summand (or its derivative) [Piecewise] <<1>> is not numerical at point n = 16.

The error message suggests evaluating f and its derivative at 16:

f[16]
f'[16]
(*
  0.0454989
  -0.0802852 + 0.0523825 Derivative[1][Re][16]
*)

Aha! The derivative of f has a Re'[16] in it? It must come from some nonanalytic complex function. The only one in the definition of f[] is Abs[]. Let's try the relatively new RealAbs[] in its place.

f[n_] := Piecewise[{{Cos[Omegan[n]*t/2]^2*
     RealAbs[Exp[-nAvg/2]*Sqrt[nAvg]^n/Sqrt[Factorial[n]]], 
    0 <= n <= 20}}, 0]

NSum[f[n], {n, 1, 100}]
(*  2.69803  *)

Success!

---

Reply to comments:

Most methods of Sum[] seek a symbolic solution. The methods of NSum use numerical methods to predict the sum without adding up all the terms. Such numerical methods are slower than adding up a small number of terms (less than, say, 10^5 or so depending on the summand). They also tend to rely on the summand being analytic (representable by power series), so NSum seems ill-suited for Piecewise[] summands whose support is bounded by a small interval.

Sum[] with Method -> "Procedural" simply adds up all the terms. It is similar to Total@Table[N@f[n], {n, 1, 100}] in both method and speed. (First suggested by @Mariusz.)

Examples: The summands Piecewise[{{(-1.)^n/n, n < 70}}, 0.] produces a Real number in all cases; Piecewise[{{(-1)^n/n, n < 70}}] produces an Integer when n >= 70 and otherwise either a Real in NSum[] because it is numericized or a Rational number normally. NSum[] uses extrapolation to get -Log[2.], which is incorrect but only by a little; however if we shorten the domain to n < 60, the sum is unacceptable. One can inspect the difference in timings for oneself.

nn = 10^6;
summand = Piecewise[{{(-1)^n/N@n, n < 70}}];
nres = NSum[summand, {n, nn}] // AbsoluteTiming
Sum[N@summand, {n, nn}, Method -> "Procedural"] // AbsoluteTiming
Total@Table[N@summand, {n, nn}] // AbsoluteTiming
(*
{0.126167, -0.693147}
{0.016046, -0.700341}
{0.019534, -0.700341}
*)

nn = 10^6;
summand = Piecewise[{{(-1)^n/n, n < 70}}, 0.];
nres = NSum[summand, {n, nn}] // AbsoluteTiming
Sum[N@summand, {n, nn}, Method -> "Procedural"] // AbsoluteTiming
Total@Table[N@summand, {n, nn}] // AbsoluteTiming
(*
{0.015743, -0.693147}
{0.015114, -0.700341}
{0.016413, -0.700341}
*)

nn = 10^7;
summand = Piecewise[{{(-1)^n/n, n < 70}}, 0.];
nres = NSum[summand, {n, nn}] // AbsoluteTiming
Sum[N@summand, {n, nn}, Method -> "Procedural"] // AbsoluteTiming
Total@Table[N@summand, {n, nn}] // AbsoluteTiming
(*
{0.015895, -0.693147}
{0.149556, -0.700341}
{0.168185, -0.700341}
*)

nn = 10^7;
summand = Piecewise[{{(-1)^n/n, n < 60}}, 0.];
nres = NSum[summand, {n, nn}] // AbsoluteTiming
Sum[N@summand, {n, nn}, Method -> "Procedural"] // AbsoluteTiming
Total@Table[N@summand, {n, nn}] // AbsoluteTiming
(*
{0.167812, -0.693175 - 0.00530469 I}
{0.14409, -0.70155}
{0.16297, -0.70155}
*)