Applications of Summa.m

A package for Symbolic and Infinite Sums Manipulation

The main reason that prompted me to write the tiny package Summa was to allow my lazy self to check, without significant effort, the intermediate steps in the application of the power method to ordinary differential equations. Reindexing, collecting and simplifying sums may be quite straigthforward but it is nonetheless tedious when performed by hand. Incidentally, the manipulations required for that particular task extend to other fields, different from ODE solving and a few examples are given here. Despite the fact that the package is rather crude (it allows only for simple manipulations, without exerting any sort of control on user input - so beware of what you ask it to do!) there are several possible applications. The application of Summa.m to the solution of ordinary differential equations with the power method and with Frobenius' method is described in a separate notebook "Power Method with Summa.m".

A side effect of this package is that it is possible to use it to typeset formulae containing sums, preventing  any automatic and premature simplification by Mathematica; version 2 of Summa.m also allows the user to switch this behavior on and off at will, either globally or locally (i.e. limited to a single espression or part of an expresion). With version 2.3 it is also possibile to switch on and off the evaluation of When the expressions derived from subsequent manipulations of sums are needed, they can be obtained step by step as if 'by hand', thus avoiding direct typesetting of each intermediate formulae. But this is not merely a typesetting package: it can profitably be used to do some 'serious' symbolic manipulations on Sums as this notebook should now make clear.

In what follows, the adjectives 'symbolic', 'infinite' and 'numeric' are used to denote sums with one or both symbolic, infinite or numerical (integer, to be exact) endpoints.

Initialization

This loads the package (remember to place it into the path of Mathematica):

<<Summa.m

The symbols exported are given by ?Summa`*, but after the renaming occured in version 2.3 a more sensible list can be found with

?Summa`Sum*

The basic object defined in the package is Summa, a container that substitutes Sum and prevents evaluation. If you choose to disable the built-in Sum Evaluation to perform your symbolic manipulations, you won't even need to know a thing about it since it will be used transparently to the user (who will operate on Sum objects and will see output formatted as Sums.

?Summa

The other two main functions are devoted to enabling and disabling (either locally and globally) the original built-in evaluation associated with Sum. They are described in the next paragraph

Basic manipulation of finite symbolic sums

The following will show how to use the package Summa to find a closed form solution for the sum of an arithmetic progression. In this we will adopt the method devised by Gauss.

\(TraditionalForm\`\(	S = ∑\+\(k = 0\)\%n\((a + b\ k)\)\)\)

We start by writing our Sum with the standard order of summation: from 0 to n

Sa = Sum[a + b k, {k, 0, n}]

We can look at the sum in a different way, thanks to the commutative property. In fact we can begin to to add its terms from the end. We can thus change the index from k to n-k

and use the function SumReverse in order to state the equivalence between \(TraditionalForm
      \`\(c\_n + c\_\(n - 1\) +  ... \) + c\_1 + c\_0 = 
        \(c\_0 + c\_1 +  ... \) + c\_\(n - 1\) + c\_n\)

SumReIndex[Sa, k -> n - k]

and write a function SumReverse to exchange the sum's limits.We do it by usign Summa in the body of the procedure in order to be able to use it even if later on we should choose to enable ordinary Sum evaluation. We could have used Sum, if we wanted to, because every occurence of Sum, when Sum evaluation is disabled, is converted into a Summa object.

SumReverse[Summa[el_, {k_, ko_, k1_}]] := Summa[el, {k, k1, ko}]

Now we use our new function to state the equivalence between \(TraditionalForm
      \`\(c\_n + c\_\(n - 1\) +  ... \) + c\_1 + c\_0 = 
        \(c\_0 + c\_1 +  ... \) + c\_\(n - 1\) + c\_n\):

Sb = SumReverse[%%]

Since sums Sa and Sb are the same as sum S, if we add them together we get twice the sum we are trying to express in closed form. We can then see that even if we don't know how to express S in closed form, we can easily spot the value of the sum Sa + Sb:

SumSimplify[2  S == Sa + Sb]

We are almost done. We could let Mathematica act on the sum at the rhs to find its closed form, but it would not be fair; after all, if we wanted mathematica to do all the work for us we cold have evaluated the sum right at the first step. Since we want to push our manipulation as far as possible, we take the index-free term (2a+b n)  out of the sum, by means of BringOut

BringOut[%, 2  a + b n]

Now all we have to do is to count from 0 to n. That would be... let's see... (n+1). At this point we could even let Mathematica do this step for us, thus saving us a manual replacement of the two elementary sums

EnableSumEvaluation[%]

We can finally solve for S and find the closed form solution for our arithmetic progression problem:

Solve[%, S] // Simplify

That's it. Of course Mathematica knew this answer all along. But we have been able to perform every single step of the solution as if by hand, verifying our technique without having to carry on tedious substitutions with pencil and paper. Granted, the example is quite trivial but covering the boring steps of sum manipulation is indeed one of the aims of the Summa package.

Perturbation method for finite sums

The following examples are taken from "Concrete Mathematics", by Graham, Knuth, Patashnik. In chapter 2 the authors show how to apply the perturbation method to get a closed form for the sum of the first n powers of x and of element in the form k (x^k).

The perturbation method applied to the sum

(Subscript[S, n]) = (∑_(k = 0)^n) (Subscript[a, k])

consists in writing the sum (Subscript[S, n + 1])with (n+1) terms in two different ways: one that takes the last term out thus exposing the sum (Subscript[S, n])of the first n terms, and another one that takes the first term out. We then rearrange the latter sum's index in order to have the same starting point as the sum (Subscript[S, n])

\(TraditionalForm\`\(	
    S\_\(n + 1\) = 
      \(\(∑\+\(k = 0\)\%\(n + 1\)a\_k = 
            \(∑\+\(k = 0\)\%n a\_k + a\_\(n + 1\) = 
              S\_n + a\_\(n + 1\)\)\)
	S\_\(n + 1\) = 
        \(∑\+\(k = 0\)\%\(n + 1\)a\_k = 
          \(a\_0 + ∑\+\(k = 1\)\%\(n + 1\)a\_k = 
            a\_0 + ∑\+\(k = 0\)\%n a\_\(k + 1\)\)\)\)\)\)

And then we equate the two forms to obtain an equation

\(TraditionalForm\`\(	
    S\_n + a\_\(n + 1\) = a\_0 + ∑\+\(k = 0\)\%n a\_\(k + 1\)\)\)

that can be solved if we can rewrite the last sum as a function of initial sum (Subscript[S, n])(when powers of x are involved, multiplying or dividing by the correct powers of x can work the trick, as we shall soon see)

\(TraditionalForm\`\(	S\_n + a\_\(n + 1\) = f(S\_n) + a\_0\)\)

We thus have an equation in (Subscript[S, n])whose solution is the closed form we were looking for.

Geometric Progression

This example makes use of the function BringOut[sum,x] introduced in version 2.2 of the package. Here we show all the passages in one single pass. We start with (Subscript[S, n])

(Subscript[S, n])=(∑_(k = 0)^n) (x^k)

and we generate all the steps needed to find a closed form for this sum. Please note that on certain older versions Solve could choke on the subscripted variable.

eq = (Subscript[S, n]) + (x^(n + 1)) == Sum[(x^k), {k, 0, n + 1}]eq = SumSplit[eq, 1]eq = SumIndexAlign[eq, {k, 0}]eq = BringOut[eq, x]eq = eq /. Sum[x^k, {k, 0, n}] -> (Subscript[S, n])Solve[eq, (Subscript[S, n])]

Once we have this all done, we could intersperse the output with our comments. Instructions on the passages to be carried on has spared us a dull input session of bidimensional forms. The gain in time is better appreciated when treating sums whose terms depend on different values of the index.

Two methods to find a sum's closed form

Starting from our expression for the closed form of the above sum, we can derive a further relationship by differentiating both sides of the folowing equality:

eq = (∑_(k = 0)^n) (x^k) == ((1 - (x^(n + 1)))/(1 - x));

SumD[eq, x] // Simplify

We can align the exponent of the sum in the lhs without changing the endpoints. Using BringIn to multiply the addendum by x and premultiplying the whole sum by 1/x we get:

BringIn[%, x]

by multiplying both sides by x we get the closed form for the sum of k*(x^k). We use a pure function and Map to perform this task on both sides

x #& /@ %

We could have obtained the same expression by means of the perturbation method applied to the sum investigated. Following Graham, Knuth and Patashnik we start from

(Subscript[S, n])=(∑_(k = 0)^n) k (x^k)

and we perturb the sum (Subscript[S, n + 1])in the following way:

eq = 
      (Subscript[S, n]) + (n + 1) (x^(n + 1)) == 
        Sum[k (x^k), {k, 0, n + 1}]

The first term of the sum at the rhs is 0*x^0=0. By expliciting this and by realigning the index at 0 we get

\(eq = SumSplit[eq, 1]; 
eq = SumIndexAlign[eq, {k, 0}]\)

Now, to express the sum in the above equation in term of our initial (Subscript[S, n]), we expand the sum with ExpandAllSum (that iterates ExpandSum until everything has been 'unknotted')

eq = SumExpandAll[eq]

If we divide each addendum by x, premultiplying the sum by the same amount, we can get the exponent of the powers of x inside the sums aligned to just k.

eq = BringOut[eq, x]

Now all we have to do is to substitute for the known sums: the former totals to((1 - (x^(n + 1)))/(1 - x)), while the latter is the sought unknown (Subscript[S, n]).

eq = 
      eq /. {
Sum[ x^k, {k, 0, n}] -> 
            ((1 - (x^(n + 1)))/(1 - x)), 
          
Sum[k x^k, {k, 0, n}] -> (Subscript[S, n])}

We're done. We only have to solve for (Subscript[S, n])to get our answer:

Solve[eq, (Subscript[S, n])] // Simplify

Least square fit

Let's find the expression of the coefficients of the line

\(TraditionalForm\`\(	y\_fit = β0 + x\ β1\)\)

that fits a set of points {(Subscript[x, k]), (Subscript[y, k])} (k=1 to n) following the least square method. In order to do this we have to minimize the sum of square errors

\(TraditionalForm\`\(	
    SSE = ∑\+\(k = 1\)\%n\((y\_k - y\_fit)\)\^2\)\)

where (Subscript[y, fit]) is given by the equation above. Direct substitution of (Subscript[y, fit])in the formula for SSE yields:

\(Clear[x, y]; 

    sse = Sum[\((y\_k - \((β0 + β1\ x\_k)\))\)^2, {k, 1, n}]; 
    \)

All we have to do is to compute the partial derivatives with respect to the coefficients and to find out for which values the derivatives go to zero. Our first attempt at writing the first equation is

SumD[sse, β0] == 0 // SumExpand // BringOut

It's nice but not quite what we want. Let's get rid of the common factor 2 and simplify the sum of 1's from 1 to n by allowing Mathematica to locally operate on it as a Sum object (we must keep in mind that even though Sums are displayed, the internal representation consists of Summa objects). We restate our code in this way:

eq1 = 
      (SumD[sse, β0]/2 == 0 // SumExpand // BringOut) // 
        EnableSumEvaluation

The second equation is more easily obtained:

eq2 = SumD[sse, β1]/2 == 0 // SumExpand // BringOut

Solving for β0 and β1 yields the desired solution

coeff = Solve[{eq1, eq2}, {β0, β1}]

We can try to further simplify the expressions we obtained. Let's work on β1 first

\(b1 = β1 /. coeff[\([1]\)]; 
num = Numerator[b1]/n // Apart; 

    den = Denominator[\ b1]/n // Apart; \)

b1 = num/den

\(\(meanRule = {
	
        Sum[x\_k, {k, 1, n}] -> n\ \ 𝔼[x], 
	
        Sum[y\_k, {k, 1, n}] -> n\ 𝔼[y]}; \)\)

b1 = b1 /. meanRule // Simplify

We can do the same with β0

\(b0 = β0 /. coeff[\([1]\)]; 
num = Numerator[b0]/n // Apart; 

    den = Denominator[\ b0]/n // Apart; \)

b0 = num/den /. meanRule

But, once we know β1, it is better to find β0 in the following way

b0 = 
      β0 /. Solve[eq1, β0][[1]] /. meanRule // 
        Simplify

Solving Ordinary Differential Equations by series expansion

More on this is to be found in the separate notebook "The Power Method with Summa". 
Let's try to solve the following ODE with the power method, in a neighbor of x=0

eq = y''[x] - 2 (x + 3) y'[x] - 3  y[x] == 0;

Since the coefficients are polynomials, they are already in the form  of power series of x. We now substitute y[x] and its derivatives with power series of x=(x-0) :

\({ys = Sum[a[n] x^n, {n, 0, Infinity}], 
	y1s = DSum[ys, x], 
	
      y2s = DSum[y1s, x]}\)

We get

eqs = eq /. {y[x] -> ys, y'[x] -> y1s, y''[x] -> y2s}

eqs = SumExpandAll[eqs]

Now we write all the sums in terms of powers of x^n, by realigning the exponents.

eqs = SumExponentAlign[eqs, x, n -> n]

then we explicit the first terms (they turn out to be zero) in order to have all the sum starting from the same point. BreakSum will do this

SumSplit[eqs, 0]

After simplification we get a single sum

res = SumSimplify[%]

whose coefficient need to be zero in order for the equality to hold for every x. The equation is

eqRecur = SumCoefficient[res, x^n] == 0

And its solution will yield a recursion relation that will allow us to compute the coefficients of the original power series, once the first two are given.

Solve[eqRecur , a[n + 2]][[1]]

sol = % /. n -> n - 2 // Simplify

This will help in computing the coefficients

\(Clear[a, expr]; 
expr[n_] = a[n] /. sol; 

    a[n_] := \(a[n] = expr[n]\)\)

Once the first two are set, the subsequent ones are automatically computed and stored.

\(a[0] = 1; a[1] = 1; 
Table[a[n], {n, 0, 10}]\)

Here's a truncated form of a particular solution with a[0]=a[1]=1. We must remember to enable the built-in Sum evaluation when we compute the sum, lest we incur in an infinite recursion

yp = Sum[a[n] x^n, {n, 0, 10}] // EnableSumEvaluationPlot[Evaluate[N[yp]], {x, - .5,  .5}];

Two trivial typesetting applications

A step by step illustration of the perturbation method

As shown before, it is easy to write the simplifications needed to find the sum of the geometric progression by means of the perturbation method. Suppose we want to invert the lhs and rhs of every eqautions before adding them to a handout. We can do that by just inverting the starting equations and re-running the code for the simplification steps. By executing this cell

eq = Sum[(x^k), {k, 0, n + 1}] == (Subscript[S, n]) + (x^(n + 1))eq = SumSplit[eq, 1]eq = SumIndexAlign[eq, {k, 0}]eq = BringOut[eq, x]eq = eq /. Sum[x^k, {k, 0, n}] -> (Subscript[S, n])Solve[eq, (Subscript[S, n])]

we get all the output we need.We can now copy and paste the various steps in the text avoiding the manual rewriting of every equation. By converting every cell to text and then to traditional form, and by substituting the double equal of the Mathematica syntax we can get a pleasant looking set of equations that can be interspersed with text. The result is this:

We start by perturbing the sum (Subscript[S, n]) = (∑_(k = 0)^n) (x^k) , incrementing its endpoint. We can then write the perturbed form in terms of the original sum, whose value we wish to compute:

(∑_(k = 0)^(n + 1))(x^k)=(x^(n + 1)) + (Subscript[S, n])

What we need to do now is trying to find (Subscript[S, n])in the lhs of this equation. Let's take the first term out of the sum, to begin with

(∑_(k = 1)^(n + 1))(x^k)+1=(x^(n + 1)) + (Subscript[S, n])

If we now shift the index by 1, we find a form that resembles the initial sum.

(∑_(k = 0)^n) (x^(k + 1))+1=(x^(n + 1)) + (Subscript[S, n])

Factoring x out of the sum yields what we were looking for: the geometric sum has made its appearance in the left hand side of the equation

x (∑_(k = 0)^n) (x^k)+1=(x^(n + 1)) + (Subscript[S, n])

\(TraditionalForm\`\(	x\ S\_n + 1 = x\^\(n + 1\) + S\_n\)\)

It is then trivial to find its value by solving for (Subscript[S, n]).

\(TraditionalForm\`\(	S\_n = \(x\^\(n + 1\) - 1\)\/\(x - 1\)\)\)

This is a very trivial example and could have easily be input by handtyping every formula. But when expressions begin to be fairly complicated and nontrivial simplification occurs, such as in the demonstration application of the power series method to ordinary differential equations, the method can significantly reduce the time required and the probability of a typo.

Sums of three elements in compact form: vectors in Euclidean spaces (*)

When Summa.m is loaded, sums with finite numeric (integer, to be exact) endpoints are evaluated in the standard way by Mathematica. This could lead to the impossibility to manipulate a certain kind of sum. For example, consider the decomposition of a vector into its cartesian components in a three-dimensional Euclidean space. This is what we get if we input the following sum in a fresh Mathematica session (or immediately after loading Summa.m):

Sum[(Subscript[x, j]) (Subscript[𝕦, j]), {j, 1, 3}]

Which is, most of the time, a desirable behavior. But what if we want to perform some operations on the vector in its closed form? We should disable automatic evaluation of sums with integer endpoints. Option Extent->Numeric has been introduced in DisableSumEvaluation and EnableSumEvaluation to do just that. By disabling sum evaluation with Extent->Numeric (or with Extent->All) we prevent every form of  'sum takeover' by Mathematica's evaluator.

DisableSumEvaluation[Extent -> Numeric]

If we input the above sum now, we will see that it will be converted to a Summa object that won't be evaluated until the user explicitly enables the evaluation process. This means that it will be possible to manipulate the sum as if it were a symbolic or infinite sum. For example this is a timedependent vector r[t] that describes a curve in the three dimensional space with versors (Subscript[𝕦, j])

r = Sum[(Subscript[x, j])[t] (Subscript[𝕦, j]), {j, 1, 3}]

And this is the representation of the tangent vector.

SumD[r, t]

Enabling and disabling sum evaluation can be used to build traditional form expression with standard mma syntax.

\(expr = Sqrt[Sum[x\_j^2, {j, 1, 3}]]; 

    expr == EnableSumEvaluation[expr]\)

Again this example is outrageously trivial. But with fairly complicated expression, the user can save him or herself a considerable amount of retyping.

FourierCoefficients

What follows has shown to be a computation that can change depending on the version of Mathematica used, due to the different implementations of Integrate and Simplify  (and probably of Coefficient), I guess. If it won't work in your version, it's probably because these functions have been made too smart. This is the result obtained with Mathematica version 3.

The Fourier series of a function f[x] can be given in terms of Sin and Cos functions:

fourierSum[x_] = 
      (Subscript[a, 0]) + Sum[(Subscript[a, n]) Cos[n x] + (Subscript[b, n]) Sin[n x], {n, 1, ∞}]

To find the coefficients (Subscript[a, 0]), an and bn we rely upon the orthogonality property of  such  functions:  the scalar product with  the corresponding function of index  m  will set to zero every term except the one of index m. In this case the scalar product is an integral from -π to +π of the product of the functions. Using IntegrateSum instead of Integrate willl allow to thread the operation over the whole equation

Computing (Subscript[a, 0])

To find (Subscript[a, 0]) we compute the scalar product with Cos[0 x]=1, that is to say we just need to integrate over the interval -π__π

SumIntegrate[f[x] == fourierSum[x], {x, -π, π}] // Simplify

Since n is integer and not zero, we have that Sin[n π]/n=0 and thus the series in the right hand side is zero. Hence we obtain

eq1 = % /. Sin[n Pi] -> 0 // SumSimplify

Solve[eq1, (Subscript[a, 0])]

Computing (Subscript[a, m])(m≥1)

To find the coefficients am, with m≥1 we compute the scalar product with Cos[m x]

\(SumIntegrate[f[x] Cos[m\ x] == fourierSum[x] Cos[m\ x], 

        			{x, \(-π\), π}] // FullSimplify\)

The first term on the right hand side is  zero. The same is true for the terms in the series with n≠m, since they are sinuses of multiples of π divided by non zero quantities.We are left with just one term: the one with n=m which can be extracted with SumToTerm (that simply returns the general element inside Summa) . Since m+n is nonzero we have no problems in getting rid of the part with Sin[(n+m)π].

eq1 = 
      SumToTerm[%] /. Sin[m π] -> 0 /. Sin[(m + n) Pi] -> 0

We can't do the same for the Sin[(m-n)π]/(m-n) part. We have to take a limit:

eq2 = 
      eq1[[1]] == Limit[eq1[[2]] /. m - n -> z, z -> 0] /. n -> m

and then we can solve for (Subscript[a, m])

Solve[%, (Subscript[a, m])]

Computing (Subscript[b, m])(m≥1)

In the same way we can find the coefficients (Subscript[b, m])

\(SumIntegrate[f[x] Sin[m\ x] == fourierSum[x] Sin[m\ x], 

        			{x, \(-π\), π}] // FullSimplify\)

eq1 = SumToTerm[%] /. Sin[(m + n) Pi] -> 0

eq2 = 
      eq1[[1]] == Limit[eq1[[2]] /. m - n -> z, z -> 0] /. n -> m

Solve[%, (Subscript[b, m])]

Integrating by series, through differentiation

This is an interesting example that emerged from some notes I took years ago. It shows how to compute an integral by means of the series integration theorem, but the closed form actually found is obtained by means of the series differentiation theorem. Since we do not want Mathematica to mess with our integrals, producing its output, we define here an overriden version of integrate, that does essentialy nothing, and act as a wrapper. The code to do that in included in the appendix to this section.

Here is the integral to compute, for a generic value of the real parameter z:

\(TraditionalForm\`\(	
    F\ \((z)\) = 
      ∫\_\(-\(1\/2\)\)\%\(+
            \(1\/2\)\)\(e\^\(z\ x\)\/\(1 - x\)\) 
          dx\)\)

The integrand is continuos and bounded for every value of z and for x≠1. We can find a series equivalent of the integrand function by observing that when |x|<1 the function 1/(1-x) is the sum of the geometric series (∑_(k = 0)^∞)(x^k). We can then write our integral as

step1 = 
      Integrate[Exp[z x] Sum[x^k, {k, 0, ∞}], {x, -1/2, 1/2}]

We can bring the exponential factor inside the sum with index k. The functions of the resulting series are bounded over the interval of integration, and the convergence in uniform.

step2 = step1 // BringIn

It is then possible to integrate by series, thus reducing the integral of the series to the series of integrals of the functions. If we disable the builtin integration routine, we can go a little bit farther on our own. Exchanging the signs of integration and summation we get

\($EnableIntegrate = False; 

    step3 = step2\  /. 
        Integrate[Sum[el_, it_], d_] :> 
          Sum[Integrate[el, d], it]\)

We could have left Mathematica compute the integrals

(Subscript[ℐ, k])(z)=(∫_(-
                  (1/2))^((1/2)))(e^(x z)) (x^k) 
                dx

for us, but there is a neater form for the results that can be achieved by considering that

∂(((Subscript[ℐ, k])(z))/∂)z=(∂/∂)z 
                (∫_(-
                      (1/2))^((1/2)))(e^(x z)) (x^k) 
                    dx
=(∫_(-
                      (1/2))^((1/2)))(∂/∂)z 
                    ((e^(x z)) (x^k)) 
                    dx=(∫_(-(1/2))^((1/2)))x 
                        (e^(x z)) (x^k) 
                      dx=(∫_(-
                        (1/2))^((1/2)))(e^(x z)) 
                        (x^(k + 1)) dx
=(Subscript[ℐ, k + 1])(z)

We can iterate this relation k times to espress the k-th element in terms of the 0-th element

\(TraditionalForm\`\(	
    \(ℐ\_k\)(z) = 
      \(∂\(\(ℐ\_\(k - 1\)\)(z)\)\/∂z = 
        \(∂\^2\(\( ℐ\_\(k - 2\)\)(z)\)\/∂z
                \^2 = \(∂\^3\(\( ℐ\_\(k - 3\)\)(z)
                  \)\/∂z\^3 = 
            \(… = 
              \(∂\^\(k - 1\)\(\(ℐ\_1\)(z)
                      \)\/∂z\^\(k - 1\) = 
                ∂\^k\(\( ℐ\_0\)(z)\)\/∂z
                      \^k\)\)\)\)\)\)\)

Now all we need to compute is the value of the integral for the particular case k=0. The integral is so easy it could be done by hand. But we have kept Mathematica too much time under the carpet, so let's feed it the following:

\(lhs = Integrate[Exp[z\ x] x^k /. k -> 0, {x, \(-1\)/2, 1/2}]; 
    
$EnableIntegrate = True; 

    rhs = Integrate[Exp[z\ x] x^k /. k -> 0, {x, \(-1\)/2, 1/2}] // 
        FullSimplify; 
lhs == rhs\)

Our integral is thus

F (z)=(∫_(-(1/2))^(+
                  (1/2)))((e^(z x))/(1 - x)) 
                dx=(∑_(k = 0)^∞)(Subscript[ℐ, k])(z)=(∑_(k = 0)^∞)(∂^k)(( (Subscript[ℐ, 0])(z)
                                )/∂)(z^k)=(∑_(k = 0)^∞)((∂^k)/∂)(z^k)((sinh(z/2))/(z/2))

And it seems to be suffering from poor numerical stability near z=0.

\(F[z_] = Sum[2\ D[Sinh[z/2]/z, {z, k}], {k, 0, 50}]; 

    Plot[F[z], {z, 10, 15}, PlotRange -> All, Frame -> True]; \)

\(F[z_] = 
      Sum[\ D[\((E\^\(z/2\) - E\^\(\(-z\)/2\))\)/
            z, {z, k}], {k, 0, 50}]; 

    Plot[F[z], {z, 10, 15}, PlotRange -> All, Frame -> True]; \)

Appendix: Code to override Integrate

$EnableIntegrate = False; 
$Auto = False;

Unprotect[Integrate]; \(
\)Integrate[args__ /; $EnableIntegrate == False] := 
      
Block[
{$EnableIntegrate = True}, 
        
integrate[args]
]\(
\)Protect[Integrate];

Format[integrate[el_, it_]] := 
      
Block[
{$EnableIntegrate = True}, 
        
If[$Auto, 
Integrate[el, it], 
          
HoldForm[Integrate[el, it]]
]
          
]