**NumericalMath****`****NSeries****`**

The Mathematica function Series finds a power series expansion of a function about a point by evaluating derivatives. However, sometimes you only want a numerical approximation to a series and the derivatives are difficult to evaluate. You could use ND to evaluate the derivatives, but that would involve repeating most of the work (i.e., evaluating the function itself at several points) for each derivative. If the function is analytic and can be evaluated in the complex plane, the solution is to use the NSeries.

Numerically finding series expansions.

NSeries works by evaluating the function at discrete points on a circle centered at the expansion point. It then uses InverseFourier to express the function as a sum of orthogonal basis functions, each of which is trivial to integrate around the circle. By Cauchy's integral formula it is seen that InverseFourier essentially gives us the required series coefficients; very little additional work is required and NSeries does it all for you. Because it is based on Cauchy's integral formula, care should be taken to ensure that no singularities of the function occur inside the circle of sample points.

Options for NSeries.

This loads the package.
In[1]:= **<< NumericalMath`NSeries`**

This gives an approximation to the series expansion of about the origin. Note that NSeries

is unable to recognize small numbers that should in fact be zero.
In[2]:= **NSeries[Exp[x], {x, 0, 3}]**

Out[2]=

Chop and Rationalize can often clean up the result.
In[3]:= **Rationalize[Chop[%]]**

Out[3]=

If the radius is too small cancellation becomes a problem.
In[4]:= **Rationalize[Chop[NSeries[Exp[x], {x, 0, 5},**

Radius -> 1/8]]]

Out[4]=

Increasing the WorkingPrecision will usually solve this problem.
In[5]:= **Rationalize[Chop[NSeries[Exp[x], {x, 0, 5},**

WorkingPrecision -> 40, Radius -> 1/8]]]

Out[5]=

If the radius is too large, the sample points become too scattered and it starts looking like poles exist.
In[6]:= **Rationalize[Chop[NSeries[Exp[x], {x, 0, 5},**

Radius -> 4]]]

Out[6]=

NSeries also recognizes poles.
In[7]:= **Chop[NSeries[Zeta[s], {s, 1, 5}]]**

Out[7]=