Summa: 
Symbolic sums manipulation

Avoiding Mathematica's evaluator for symbolic and infinite Sums

The purpose of this notebook is to show how the functionalities of Summa.m were implemented. The package defines an abstract object Summa (intended as a version of Sum that does not evaluate) and several procedures to manipulate expression containg such object. Since version 2 of the package, the user is allowed to disable and enable the ordinary Sum evaluation process either locally or globally, by automatically converting  every Sum into a Summa object or vice versa. By disabling Sum evaluation, every Sum is automatically converted into a Summa object before being fed to the evaluator: this approach allows the user to apply  manipulation routines developed for Summas directly to built-in Sums. Most of this notebook was written for an earlier version of Summa.m; even though additions were made as long as new functionalities were added to the package,  some of the functions herein described can have a slightly different implementation in the latest available version of the package.A new alternate naming has been added in order to have most functions acting on Sums begin with Sum* (notable exception are BringIn and BringOut).

Overriding Sum

The main change introduced with version 2 of this package is the ability to address Summa objects with exactly the same input syntax used for built-in Sums. In the rest of this notebook we will address, somewhat loosely, this behavior by "enabling and disabling the Sum evaluation process". The built-in evaluation process is not actualy being inhibithed: it just doesn't have a chance to see Sums anymore unless the user sets a specific environment variable that suppresses the conversions of Sum into Summa. Please refer to the notebook "Basic Overriding" for further information on this subject. A rather extensive notebook on the subject, along with a package that automates the generation of the code neeeded to override a built-in procedure or object is currently under development.

Overriding Sum

Let's see how the concepts expressed in the "basic overriding" notebook apply to Sum: first thing first we redefine Sum in terms of the user defined object Summa. We give our 'enabling' flag a unique name: $EnableSumEvaluation

\($EnableSumEvaluation = False; 
Protect[$EnableSumEvaluation]; \)

\(Unprotect[Sum]; 

    Sum[args__ /; $EnableSumEvaluation == False] := 
      
Block[
        
{$EnableSumEvaluation = True}, 
        
Summa[args]
]\)Protect[Sum];

Since we haven't written any code for Summa yet, calls to Sum will be diverted to an unevaluated Summa that simply acts as a wrapper: Note that this definitiom  will inhibit multiple sums evaluation too. At present, though, the package does not allow multiple dimensional sum manipulation, except for very simple nested sums.

{Sum[(Subscript[x, k]), {k, 2, ∞}], Sum[(Subscript[x, k]), {k, 2, 7}], 
      Sum[(Subscript[x, k]), {k, 5}], Sum[(Subscript[a, j]) ((Subscript[x, j])^k), {j, 3}, {k, 0, 2}]}

Later on we will write a code that leads to an unevaluated wrapper when the endpoints are either infinite or symbolic, but calls to the the built-in function Sum when both endpoints are integers. This is possibile because the call to Summa, when we input a Sum object, take place inside a Block with $EnableSumEvaluation=True. 
When $EnableSumEvaluation is set to False every instance of Sum is automatically converted to a Summa object. Hence, we can write the code for our symbolic maipulation functions just in terms of Summa. For example, the folowing function, that relies on replacement rules to bring a common factor out of the wrapper can be applied to Summa object, of course, but will work even when it is applied to a Sum :

\(outrule1 = {
        
Summa[f_ + g_, it : {k_Symbol, __}] :> 
          Summa[f, it] + Summa[g, it], 
        
Summa[c_\ f_, it : {k_Symbol, __}] :> 
          c\ Summa[f, it] /; FreeQ[c, k]
}; 

    outrule2 = 
      Summa[c_, it : {k_Symbol, __}] :> 
        c\ Summa[1, it] /; FreeQ[c, k]; 

    outrules = Flatten[{outrule1, outrule2}]; 

    BringOut[s_] := s //. outrules\)

BringOut[Summa[a k (Subscript[x, k]), {k, 2, ∞}]]

BringOut[Sum[a k (Subscript[x, k]), {k, 2, ∞}]]

We will delve into these functions with more details later on. For the time being, let's focus on the conversion between Summa and Sum

A step further: overriding nested sums

Nested sums pose a problem.When the Block structure is applied to the outer sum, those that appen to be inside it are evaluated as if they were ordinary Mathematica sums; in fact they fall under the scope of Block with $EnableSumEvaluation=True.Here's an example

Sum[Sum[a[j] a[k], {j, 1, n}], {k, 1, n}]

We can try to fix this by using a second Block that enclosed the arguments of every sum.This is the new version of the code

\(Unprotect[Sum]; 
Clear[Sum]; 

    Sum[args__] /; $EnableSumEvaluation == False := 
      Block[{$EnableSumEvaluation = True}, 
        
Summa[
          
Sequence@@
            Block[
{$EnableSumEvaluation = False}, 
              
tmp[args]]]
]\)Protect[Sum];

It seems to work

Sum[Sum[Sum[el, it1], it2], it3]

We can work some simplification, too

Sum[Sum[a[j] a[k], {j, 1, n}], {k, 1, n}]

% // BringOut

Enabling and disabling Sum evaluation

If we set the global variable $enablebuiltin to True, we can access to the builtin form of Sum, and still be able to access to the user defined code for Summa

\(Unprotect[$EnableSumEvaluation]; 
$EnableSumEvaluation = True; 

    Protect[$EnableSumEvaluation]; \)

Summa objects work as usual

Summa[a/k^2, {k, 1, ∞}]

And so do Sums, according to the predefined Mathematica behaviour

Sum[a/k^2, {k, 1, ∞}]

We can write a couple of procedures that unprotect and protect the environment variable for us.

EnableSumEvaluation[] := 
      
(
        
Unprotect[$EnableSumEvaluation]; 
        
$EnableSumEvaluation = True; 
        
Protect[$EnableSumEvaluation]; 
)
        
DisableSumEvaluation[] := 
      
(
        
Unprotect[$EnableSumEvaluation]; 
        
$EnableSumEvaluation = False; 
        
Protect[$EnableSumEvaluation]; 
)

Moreover, we can make this process localized to a single expression, regardless of the value of the environment variable $EnableSumEvaluation. For example, to locally enable Sum evaluation, we can let a procedure substitute every occurence of Summa with Sum when the builtìin evaluation is turned on and replacing the surviving Sums with Summa prior to outputting the result in an environment where the built-in sum evaluaton is disabled. We must give this procedure the attibute Holdfirst to avoid the evaluation of the expression passed

Attributes[EnableSumEvaluation] = {HoldFirst}; EnableSumEvaluation[expr_] := 
      
(
EnableSumEvaluation[]; 
        
First[{expr /. Summa -> Sum, 
              DisableSumEvaluation[]}] /. 
          Sum -> Summa
)
Attributes[DisableSumEvaluation] = {HoldFirst}; DisableSumEvaluation[expr_] := 
      
(
DisableSumEvaluation[]; 
        
First[{expr /. Sum -> Summa, 
              EnableSumEvaluation[]}] /. Summa -> Sum
)

Here we use First to give back a result after the last eleement of the list has been evaluated. 
The procedures implemented in the package are slightly different and have options to make them more flexible (the can print a warning message and leave the built-in evaluation process on or off once the expression has been returned). See the notebook 'Usage of Summa' for more information on these procedures and when to use them. 
Before we continue we must override ordinary sum evaluation with our own procedure:

DisableSumEvaluation[]

Enabling and disabling Sum enhancements (*)

To be added in the next major version.

Defining the Abstract Object Summa

Summa[el,{x,n1,n2}]

Definition

We definine a Summa object in such a way that it will not be evaluated unless both the start and end values of the index are numerical (integers, to be precise). In the latter case it will be evaluated as a built-in Sum. This is the idea:

Summa[el_, {k_Symbol, k0_Integer, k1_Integer}] := 
      Sum[el, {k, k0, k1}]

And this is the implementation that works for multidimensional sums as well, regardless of the form used for the iterators:

Summa[el_, 
        it : ({_Symbol, _Integer, _Integer} | {_Symbol, _Integer})
               .. ] := Sum[el, it]

This is what we get at this point for symbolic (or infinite) and finite sums

{f1 = Summa[(Subscript[a, k]) x^k, {k, 0, ∞}], 
      Summa[(Subscript[a, k]) x^k, {k, ∞}], Summa[x^k, {k, n}]}

{s1 = Summa[(Subscript[a, k]) x^k, {k, 0, 4}], Summa[(Subscript[a, k]) x^k, {k, 4}], 
      Summa[x^k, {k, 2}]}

Oh, and what about sums with multiple indices? In order to be evaluated as Mathematica Sums, every iterator should have integer endpoints:

{Sum[(Subscript[a, k]) ((Subscript[x, k])^j), {k, ∞}, {j, 1, 2}], 
      Sum[(Subscript[a, k]) ((Subscript[x, k])^j), {k, 3}, {j, 1, n}]}{Sum[(Subscript[a, k]) ((Subscript[x, k])^j), {k, 3}, {j, 1, 2}]}

We can add a line to convert the compact iterator form into its full-blown equivalent. The following line will take care of any sequence of iterators, replacing those with only two elements with their expanded counterparts. A condition has to be fullfilled in order to avoid an infinite recursion (we are now defining a Summa in terms of itself, so we need something to tell the transformed version from the untouched one).

Summa[el_, it : ({_, _} | {_, _, _}) .. ] /; 
        EvenQ[Times@@(Map[Length, {it}, -1])] := 
      Summa[el, 
        ReleaseHold[Hold[it] /. ({k_, k1_} :> {k, 1, k1})]]

{Sum[(Subscript[x, k]), {k, n}], Sum[(Subscript[a, k]) ((Subscript[x, k])^j), {k, 3}, {j, ∞}]}

This definiton will force us to write the call to Summa in a form that won't rely on iterators espressed as {k,ks} (even if ks means 'one or more elements'), because that form will be immediately reverted to {k,1,ks}. We will use the more general form for the iterator and, if necessary, we will specify the nature of the first element.

it : {_, _, _}it : ({_, _} | {_, _, _})it : ({_, _} | {_, _, _}) .. it : ({k_, _} | {k_, _, _}) .. it : ({k_Symbol, _} | {k_Symbol, _, _}) ..

A step further: a switch for sums with numerical endpoints

Since we might want to manipulate sums with numerical endpoints in their compact form, without letting mma evaluate them, it may prove useful subordinate the trasformation of Summas into Sums to the value of a boolean environmental variable. We set the value of this variable to True, by default

$EnableNumericSum = True;

And add a condition to the definition of Summa

Summa[el_, 
          it : ({_Symbol, _Integer, _Integer} | {_Symbol, _Integer})
                 .. ] /; $EnableNumericSum == True := Sum[el, it]

Then we rig the procedures EnableSumEvaluation and DisableSumEvaluation with some extra code to set the right value of this variable. The code is not shown here, since there is nothing special to it. More on this can be found in the Usage notebook.

Output Format(s)

We might not like the way our redefined sums are displayed. We are accustomed to see sums as Sum objects, and on version 3 or later two dimensional formatting represents a useful visual aid it would be a pity not to use. It seems a good idea to define for Summa  a format similar to that used by Sum. Since the function Format needs to call the built-in version of Sum, we must enclose such a call inside a Block that sets $EnableSumEvaluation to True.

\(Format[Summa[el_, it : {k_Symbol, __}]] := 
      
Block[{$EnableSumEvaluation = True}, 
		
        Sum[el, it]
]\)

When $EnableSumEvaluation is False,this definition will apply to Sums too,since in that case every Sum is automatically transformed into a Summa object before it can be evaluated,and then formatted accordingly.

{Sum[(Subscript[x, k]), {k, 2, ∞}], Summa[(Subscript[x, k]), {k, 2, ∞}], 
      Summa[(Subscript[x, k]), {k, ∞}]}

InputForm[%]

There is a shortcoming for this method, though: since the formatting results from Mathematica's  built-in evaluation of  Sum (while the sum's true form is a Summa object), every convergent sum known to Mathematica will be represented as the value toward which the sum actually converges. For example, the following sum is known to converge to e, hence it will be displayed in this way:

Summa[1/k!, {k, 0, ∞}]

Even if its internal representation is that of an unevaluated Summa object

InputForm[%]

To avoid such idiosyncracies we can wrap HoldForm around Sum. Better yet, we can let the user choose whether he wants to see the unevaluated or evaluated form of the sums he is working with (so that he can see if it is the case to apply the predefined built-in evaluation to certain parts of the expression). Switching from one behavior to the other one is made by means of the private variable $Auto,  that can be set by the exported funtion SetSumFormat. Since we want our procedure to be as general as possibile, we also allow for multiple iterators. We don't have to worry about mixed iterators (such as {k, a,b} and {k, b}) since that will be taken Our code becomes:

\($Auto = True; 

    Format[Summa[el_, it : {_Symbol, __} .. ]] := 
      
Block[
        
{$EnableSumEvaluation = True}, 
        
If[$Auto, 
Sum[el, it], 
          
HoldForm[Sum[el, it]]
]
          
]\)

When $Auto is set to False, we can see the unevaluated form of our sums

$Auto = False; {Summa[(Subscript[x, k]), {k, 2, ∞}], 
      Summa[1/k!, {k, 0, ∞}], 
      Summa[1/k!, {k, ∞}], Sum[(Subscript[a, k]) (x^k), {k, 0, 2}]}

While when it is set to True we see the result of Mathematica's evaluation (and sometimes even an error message warning about the non convergence of certain sums):

$Auto = True; {Summa[(Subscript[x, k]), {k, 2, ∞}], 
      Summa[1/k!, {k, 0, ∞}], 
      Summa[1/k!, {k, ∞}], Sum[(Subscript[a, k]) (x^k), {k, 0, 2}]}

We must keep in mind, though, that the internal representation of the sums is still a Summa object

InputForm[%]

OutputForm[%%]

If we need to convert Sum in Summa or viceversa we must resort to EnableBuiltinEvaluation and DisableBuiltinEvaluation. Our formatting can take care of multidimensional sums as well.

\($Auto = False; 
    
{Summa[\(a\_k\) x\_k\^j, {k, 2, ∞}, {j, 0, 2}], 
        Summa[j, {j, k}, {k, 0, 2}], 
        Summa[\(a\_k\) x\_k\^j, {k, 2}, {j, 0, 2}]} // ColumnForm\)

\($Auto = True; 
    
{Summa[\(a\_k\) x\_k\^j, {k, 2, ∞}, {j, 0, 2}], 
        Summa[j, {j, k}, {k, 0, 2}], 
        Summa[\(a\_k\) x\_k\^j, {k, 2}, {j, 0, 2}]} // ColumnForm\)

Same caveat as before:

InputForm[%]

SetSumFormat, defined in the package, allows to specify which format to adopt by default

\(SetSumFormat[Automatic] := \($Auto = True\); 

    SetSumFormat[Inhibited] := \($Auto = False\); \)

\(SetSumFormat[Automatic]; 
Summa[x^k, {k, 0, ∞}]\)

\(SetSumFormat[Inhibited]; 
Summa[x^k, {k, 0, ∞}]\)

Other Sum formats (*)

There are plans to add other sum format type such as

\(\(	∑\+\(a ≤ k ≤ b\)p\_k
	
        \(∑\+\(k ∈ S\)p\_k\)\)\)

Since some of this formats would require the rewriting of most of the associated manipulation procedures, only a formatting option is to be added. The idea is to allow the conversion of expression from one form into the other (so some informartion wuold be stored and hidden inside the Summa structure) and use it for typesetting purpose sonly.
But this is is still an unimplemented feature.

SetSumFormat[ftype_] := 
      
Switch[ftype, 
        

Automatic, 
        
Format[Summa[el_, {k_Symbol, k0_, k1_}]] := 
          
Block[{$EnableSumEvaluation = True}, 
            
Sum[el, {k, k0, k1}]
], 
        

Inhibited, 
        
Format[Summa[el_, {k_Symbol, k0_, k1_}]] := 
          
Block[{$EnableSumEvaluation = True}, 
            
HoldForm[Sum[el, {k, k0, k1}]]
              
], 
        

Inequality, 
        
Format[Summa[el_, {k_Symbol, k0_, k1_}]] := 
          
Block[{$EnableSumEvaluation = True}, 
            
HoldForm[Sum[el, {k, k0, k1}]]
              
], 
        

Set, 
        
Format[Summa[el_, {k_Symbol, k0_, k1_}]] := 
          
Block[{$EnableSumEvaluation = True}, 
            
HoldForm[Sum[el, {k, k0, k1}]]
              
], 
        

Characteristic, 
        
Format[Summa[el_, {k_Symbol, k0_, k1_}]] := 
          
Block[{$EnableSumEvaluation = True}, 
            
HoldForm[Sum[el, {k, k0, k1}]]
              
]
]

Extracting sum terms

GetElement[expr]

The following is  a command that allows to extract coefficients from a Summa wrapper. It's a very crude procedure and attention should be paid in order not to mix coefficients from different sums; moreover the meaning of the coeffienct extraction procedures must be clear in the mind of the user since these procedures can yield meaningless equalities. GetElement simply deletes the Summa[_, {k,a,b}] structure:

\(getElement = Summa[el_, it_] :> el; 

    GetElement[s_] := s /. getElement\)

And it does so regardless of the meaning of the expression passed to it. Here is an equation with a sound mathematical meaning:

espr = 
      2 - 3  x + 5  x^n + Summa[a[k] x^k, {k, 0, ∞}] == 0

And here is the (nonsensical but as we will se later, useful) result of GetElement.

GetElement[espr]

The user must be aware of what he did and consider this result with due care, since this equality does not follow from the above equation. It is nonetheless useful, as it will be shown later, to recover the expression of the general coefficient a[k].

Extracting coefficients from the powers of x

Use with care: it won't work with expressions that have no summa. It does not penetrate equations, either.

GetCoefficient[Summa]
GetCoefficient[expr,x,0]

GetCoefficient is more an hack than a real procedure: it somehow manages to extract the coefficients of specified powers of a given variable x from simple polynomials and single sums inside an expression. It has been slightly modified since version 1.x of this package and should be used with the utmost care in order to avoid meaningless results.

The following code simply applies GetCoefficient to the element inside the single sum passed, either alone, or as a part of a bigger expression. Please note that the code is designed to cope with only one sum at the time (it will not apply properly to nested sum, either). Hence, it's up to the user to isolate the sum before applying GetCoefficient.

GetCoefficient[Summa[el_, it_], x_, k___] := 
      GetCoefficient[el, x, k]

GetCoefficient[s_  /; Length[Cases[s, _Summa, ∞]] == 1 , 
        x_, k___] := 
      GetCoefficient[Cases[s, _Summa, ∞][[1]], x, k]

Now we define the 'real' code for GetCoefficient. On version 4 and higher Coefficient is nearly perfect for our purposes, and we can resort to it with just a slight make-up. It is mostly a matter of renaming, here, to provide an interface that could be used unaltered on previous versions too. We add a small rule that set sums with zero elements to zero in order to avoid results in the form Sum[0,{k,0,∞}].

\(GetCoefficient[el_, x_] := Coefficient[el /. Summa[0, _] -> 0, x]; 
    
GetCoefficient[el_, x_, n_] := 
      Coefficient[el /. Summa[0, _] -> 0, x, n]; \)

But, rule apart, why do we have to bother about renaming Coefficient in the first time? The answer is that on version 3 and below it can produce results that need further manipulation. For these versions of Mathematica, when the exponent is set to zero or is symbolic, Coefficient gives results we don't actually need for our specific purposes. Moreover, it can lead to spurious results when the expression and the the power of x are both symbolic. The undesired behaviour is exemplified by the following results

Coefficient[2 - 3  x + 5  x^n, x, 0]

Coefficient[2 - 3  x + 5  x^n, x, n]

Coefficient[4 + c x^(k - 3), x^k, 1]

(on version 4 and higher - I have yet to check the results in version 5, though - the above commands result in 2, 5, 0, respectively). At any rate, even if version 3 results can be explained by a 'superficial' parsing of the expressions, they give place to an undesired behaviour that does not fit our purposes.  Moreover, even in version 4 the behaviour of Coefficient, though symbolically correct, is not fully consistent with what we specifically need; for example we would like to get 2 as an answer to:

Coefficient[2 - 3  x + 5  x^n, x^n, 0]

We thus need some extra code to take care of the undesired behaviour of Coefficient. 
	- Enters the hack -
We can for instance check that the result does not contain the specified variable (please note that this definition actually replaces the one, with the same signature, previously introduced.

GetCoefficient[el_, x_, n_] := 
      
Block[
        
{res = Coefficient[el /. Summa[0, _] -> 0, x^n]}, 
        
If[FreeQ[res, x], res, 0]
];

Since this rule applied to n=0 lead to an invalid call to Coefficient (x^0 evaluates to 1 and Coefficient won't be passed a symbolic variable) we have to add one more rule to take care of this.  Incidentally this line of code fixes the behaviour that can give back a function of x when the exponent is set to 0 even on later versions of Mathematica :

GetCoefficient[el_, x_, 0] := 
      el /. {x^_ -> 0, x -> 0} /. Summa[0, _] -> 0;

The main rule works fine as long as x is a symbol. When the user specifies a power of x, we must split it in order to separate the symbol from the exponent (we could make the signature of these procedure more rigorous by specifying a particular pattern matching scheme, but we are assuming that the user knows what she's doing). We use a default patern to pack in a single lina a call to GetCoefficient[el,x^k], too.

GetCoefficient[el_, x_^k_, n_:  1] := GetCoefficient[el, x, k*n];

The actual code in the package

We now have a code that can consistently work on different versions of Mathematica. We can add a conditional statement that checks the version number so that the 'hack' will be evaluated only on the versions that really need it.

\(\( (*common\ code*) 
      
GetCoefficient[Summa[el_, it_], x_, k___] := 
      \(GetCoefficient[el, x, k]

          GetCoefficient[
            s_\  /; Length[Cases[s, _Summa, ∞]] == 1\ , x_, 
            k___] := 
        GetCoefficient[\(Cases[s, _Summa, ∞]\)[\([1]\)], x, k]

          
\(GetCoefficient[el_, x_] := 
            Coefficient[el /. Summa[0, _] -> 0, x]; 
          \)
\(GetCoefficient[el_, x_^k_, n_:  1] := 
            GetCoefficient[el, x, k*\ n]; 
          \)
\(GetCoefficient[el_, x_, 0] := 
            \(el /. {x^_ -> 0, x -> 0}\) /. Summa[0, _] -> 0; 
          \)\)\)\)

If[$VersionNumber ≥ 3, 
      

  
        
GetCoefficient[el_, x_, n_] := 
        Coefficient[el /. Summa[0, _] -> 0, x, n], 
      

  
        
GetCoefficient[el_, x_, n_] := 
        
Block[
          
{
            res = Coefficient[el /. Summa[0, _] -> 0, x^n]}, 
          
If[FreeQ[res, x], res, 0]
]
          
]

Testing GetCoefficient

We can test it, now.

GetCoefficient[2 - 3  x + 5  x^n, x, 0]

GetCoefficient[2 - 3  x + 5  x^n, x, 1]

GetCoefficient[2 - 3  x + 5  x^n, x^n, 0]

GetCoefficient[2 - 3  x + 5  x^n, x^n, 1]

GetCoefficient[2 - 3  x + 5  x^n, x, n]

GetCoefficient[2 - 3  x + 5  x^n - 7  x^(2  n), x^n, 2]

GetCoefficient[espr, x^k]

Basic Operations on Summa Objects

Almost every procedure defined in the package Summa.m is inherently rule-based and can be applied to every Summa object contained in the expression passed to it. This was done to ease the manipulation of equations, avoiding the burned of a mapping on both sides of the equality sign. The procedures are designed for onedimensional sums only. Multidimensional sums can be treated only when they can be explicitly expressed as nested onedimensional sums.

Linear combination and simplification

BringOut and ExpandSum (SumExpand)

BringOut[expr]

The following two rules bring everything that does not explicitly depend upon the index out of the sum

\(\(outrule1 = {

        Summa[f_ + g_, it : {k_Symbol, __}] :> 
          Summa[f, it] + Summa[g, it], 

        Summa[c_\ f_, it : {k_Symbol, __}] :> 
          c\ Summa[f, it] /; FreeQ[c, k]
}; \)\)

outrule2 = 
      Summa[c_, it : {k_Symbol, __}] :>  c Summa[1, it] /; FreeQ[c, k];

BringOut[s_] := s //. Flatten[{outrule1, outrule2}]

NEW! in Summa version 2.2

outrule3[x_] = 
      Summa[el_, it : {k_Symbol, __}] :>  x Summa[el/x, it] /; FreeQ[x, k];

BringOut[s_, x_] := s /. outrule3[x]

Here's an example

f1 = Summa[(x - 1)^2  x x^k (Subscript[a, k]), {k, 0, ∞}]

f2 = BringOut[f1]

BringOut[Summa[f[x], {k, 0, ∞}]]

ExpandSum[expr]

ExpandSum expands what is inside every sum and  makes a separate sum for each term. It is usually used in connection with BringIn (see later) to perform simplifications and reindexing

ExpandSum[s_] := ExpandAll[BringOut[s]]

f3 = ExpandSum[f2]

BringIn and SimplifySum (SumSimplify)

BringIn[espr]

BringIn[expr] brings the eventual factor pre-multiplying a sum inside the sum symbol. In doing so it performs an automatic expansion of what is inside.

inrule1 =  
      c_ Summa[el1_, it : {k_Symbol, __}] :> 
        Summa[If[$AutoExpand, c el1 // Expand, c el1], it] /; FreeQ[c, k];

To avoid this default expansion, the internal variable $AutoExpand has to be set to False. This is a variable that's not accessible to the user and will be modified only through the options of the calling procedures.

$AutoExpand = True;

\(\(BringIn[s_] := s /. inrule1
\)\)

NEW! in version 2.2

inrule3[x_] = 
      Summa[el_, it : {k_Symbol, __}] :>  
        1/x Summa[ x el, {k, lim}] /; FreeQ[x, k];

BringIn[s_, x_] := s /. inrule2[x]

Here's an example

ExpandSum[f1]

BringIn[%]

SimplifySum[espr]

SimplifySum[espr] reduce all the sums with the same  indexes (that is with the same index name, start and final values) to one single sum.

\(\(inrule2 = 
      Summa[el1_, it_] + Summa[el2_, it_] :> 
		
        Summa[If[$AutoSimplify, 
				el1 + el2 // Simplify, 
				
            el1 + el2], it]; \)\)

The variable $AutoSimplify is not accessibible to the user and will be modified by means of Options of the functions calling the rule.

$AutoSimplify = True;

\(\(SimplifySum[s_] := \(s /. inrule1\) //. inrule2
\)\)

This is how it works:

expr = ExpandSum[f1 + x f2] == 0

SimplifySum[expr]

SimplifySum[espr]

Sometimes it may be useful to operate a simplification of each sum separately. For this task we can use the procedure

\(\(SimplifyEachSum[s_] := 
      s /. \ Summa[el_, it_] :> Summa[Simplify[el], it]
\)\)

SimplifyEachSum[BringIn[expr]]

Reindexing

BreakSum (SumSplit)

BreakSum[espr,n]

BreakSum[espr,n0] splits all the sums contained in espr in order to get one sum that  starts from the value n0. The sums with numerical starting and ending points are automatically expanded (as is for Sum).

\(\(splitrule[n_] = 
      Summa[el_, {k_, ko_, k1_}] :> 
	
        Summa[el, {k, ko, n - 1}] + Summa[el, {k, n, k1}]; \)\)

\(\(BreakSum[s_, n_Integer] := s /. splitrule[n]
\)\)

Here is a sample sum and its equivaled decomposition with one part whose index starts form 0

(f0 = Summa[(Subscript[a, k]) x^k, {k, -∞, ∞}]) 
      == (f01 = BreakSum[f0, 0])

Here we create a sum whose index ranges from 2 to the highest value

(f1 = Summa[(x - 1)^2  x x^k (Subscript[a, k]), {k, 0, ∞}]) 
      == (f12 = BreakSum[f1, 2])

As is the case for the other procedures, BreakSum too can be applied to composite expressions, such as equations

BreakSum[f0 + f1 == 0, 2]

%  //  SimplifySum

ReIndexSum (SumReIndex)

ReIndexSum[ espr, k-> ( k+offset) ]

ReIndexSum[espr,k->new[k]] changes the index as indicated in new[k], where new[k] is a 1:1 function of k over the set of possible values permitted (a permutation of the indices k). The transformation affects every sum contained in espr that has index k. It is not allowed to change the index name (a simple replacement can do that) and new[k] must explicitly depend on k. Linear functions are most commonly used to perform this task. 
The new version allows for a more flexible manipulation of the indices. When the endpoints are finite, we can use solve to allow for more general functions new[k].

\(\(ruleIndex[k_, newk_] /; \(! FreeQ[newk, k]\) = 
      Summa[el_, {k, ko_, k1_}] :> 
	
        Summa[
el /. k -> newk, 
          
{k, k /. First[Solve[newk == ko, k]], 
            k /. First[Solve[newk == k1, k]]}
]; 
    \)\)

It is possible to check for the existence of the solution (Semantica seems to be a neat way to do that, but its beyond the scope of this notebook).

Care must be exerted in order to adapt the function to the finite set of indices. Also, Solve can cope well with infinities, but it could be advisable to give it a little help in order to avoid to specifiy that a quantity has to be real and positive. For example, suppose we perform the following linear transformation of the index

Summa[(Subscript[a, k]), {k, 0, ∞}] /. ruleIndex[k, a k + b]

The sign of the coefficient a is important in making the correct choice of Infinity. We could nonetheless fix up a solving function that not only checks for the invertibility of the index transformation, but automatically assumes that every symbolic parameter contained therein be positive. Maybe this will be added in a future realease. For the time being we rig up ReIndexSum in the usual rule-based fashion.

ReIndexSum[s_, k_ -> j_] := 
      s /. ruleIndex[k, j]

Let's see a couple of examples: this shifts the index of the following sum by two positions to the right

f1 + x f2 // BringInReIndexSum[%, k -> k - 2]

AlignIndex (SumIndexAlign)

AlignIndex[ espr, {k,k0} ]

We can define a rule that will align the indexes of every summa object contained in expr. This is strictly a change of variable, and should not be confused with BreakSum, that explicitates part of the sums in order to have them start from (or end to) the same value of the index.

\(\(indexShift[{k_, o_}] = 
      Summa[el_, {k, ko_, k1_}] :> 
	
        Summa[el /. k -> k - o + ko, {k, o, k1 - ko + o}]; \)\)

\(\(AlignIndex[s_, {k_, o_}] := s /. indexShift[{k, o}]
\)\)

As an example consider

f = 
      Summa[(Subscript[a, k]), {k, 2, ∞}] - 
        5  Summa[(Subscript[b, k]), {k, -1, ∞}]

AlignIndex[f, {k, 0}]

We can overload this procedure in order to handle the sums with indices ranging down to minus infinity. A 'fake' option will suite our needs.

\(\(indexDownShift[{k_, o_}] = 
      Summa[el_, {k, ko_, k1_}] :> 
	
        Summa[el /. k -> k - o + k1, {k, ko - k1 + o, o}]; \)\)

\(\(AlignIndex[s_, {k_, o_}, Direction -> Down] := 
      s /. indexDownShift[{k, o}]

	\)\)

Let's test it

f = 
      Summa[(Subscript[a, k]), {k, -∞, 2}] - 
        5  Summa[(Subscript[b, k]), {k, -∞, -1}]

AlignIndex[f, {k, 0}, Direction -> Down]

AlignExponent (SumExponentAlign)

AlignExponent[espr, x, k->k0]

AlignExponent allows to shift the indexes of all the series contained in expr in such a way that all the exponents of x are the same. This line makes it easier to treat non-numerical basis (but it has to be applied after due simplification)

GetExponent[espr_, x_] := 
      Flatten[Cases[{espr}, x^n_  -> n, ∞]][[1]]

Here's how it differs from the built-in function Exponent

\(espr = \((x - 1)\)^n + 3  x^5; 
    
{{"Expression", "Exponent", "GetExponent"}, {3  x^5, 
          Exponent[espr, x], GetExponent[espr, x]}, {\((x - 1)\)^n, 
          Exponent[espr, x - 1], GetExponent[espr, x - 1]}} // TableForm\)

The rule that does most of the re-aligning work is:

\(\(powerShift[x_, k_Symbol -> kin_] = 
      Summa[el_, {k, k0_, k1_}] :> 
	
        Summa[el /. k -> kin - GetExponent[el, x] + k, {k, 
            k0 - kin + GetExponent[el, x], k1 - kin + GetExponent[el, x]}]; 
    \)\)

The procedure to employ has an option, Explode, that can be used to disable the default simplification

\(Options[AlignExponent] = {Explode -> True}; 

    AlignExponent[s_, x_, k_ -> kin_, opts___] := 
	
      Module[
		{opt1, res}, 
		
        opt1 = \(Explode /. {opts}\) /. Options[AlignExponent]; 
		
        res = If[opt1, 
				BringIn[ExpandSum[s]], 

              				$Explode = False; BringIn[s]
				] /. 
            powerShift[x, k -> kin]; 
		$Explode = True; 
		res
	]
        \)

Some examples follow

f = Summa[(Subscript[a, k]) x^(k - 2), {k, 0, ∞}]

AlignExponent[f, x, k -> k + 1]

f + f1

AlignExponent[f + f1, x, k -> k]

f = Summa[(Subscript[a, k]) x^(k - 1/2), {k, 1, ∞}]

AlignExponent[f, x, k -> k - 3/2]

Here is an example in which it is not devisable to perform an automatic simplification :

s1 = Summa[(Subscript[a, k]) (x - xo)^(k - 2), {k, 0, ∞}]

AlignExponent[s1, x - xo, k -> k + 1]

The fact is that the expansion performed inside the sum destroys the (x-xo)^(k-2) term, hence the matching is not correct. Setting Explode to False will avoid this problem

AlignExponent[s1, x - xo, k -> k + 1, Explode -> False]

Term by Term Differentiation and Integration

Basic term by term Differentiation

The following procedures are intended for using in the 'by hand' resolution of ordinary differential equations by means of the power series and Frobenius methods. The differentiation of the sums in these instances is performed at the beginning of the manipulation, and it is usually done on each series separately. Yet the procedures allow for  somewhat a higher flexibility.

DSum

DSum[espr,x]
DSum[espr,{x,n}]

Option: Inside (True/False)

The basic rules to do this contemplate the basic properties of the differentiation operation:  we resort to a dummy head "d" to perform our simplifications in subsequent steps. Since version 2.1 the restriction that the variable with respect to which differentiate be a symbol has been removed. This will allow to differentiate with respect to indexed variables such as x[j] or (Subscript[x, j]). Since no error control is performed, it is up to the user not to submit nonsensical variables to the procedure.
First thing first we implement the basic linearity relationship

Clear[f, f1, f2, c, x, d, diffRule0, diffRule1, diffRule2, diffRule]

\(\(diffRule0[x_] = 

      	{
	d[c_\ f_] :> c\ d[f]\  /; \ FreeQ[c, x], 
	
        d[f1_ + f2_] :> d[f1] + d[f2]
	}; \)\)

After that we can put other more specific rules. The product rule is the only one necessary to make the command suitable for the application of the power method to ordinary differential equations. Other rules can be added, though, most notably the rules for the derivative of the quotient and of the power of a function.

\(\(diffRule1[x_] = 

      	{
	d[f1_\ f2_] :> d[f1] f2 + f1\ d[f2], 
	
        d[f_^c_] :> c\ f^\((c - 1)\)\ d[f]\  /; \ FreeQ[c, x], 
	
        d[f1_/f2_] :> \((d[f1]\ f2 - f\ 1  d[f2])\)/f2^2, 
	
        d[f1_[f2_]] :> \ d[f2] D[f1[dummy], dummy] /. dummy -> f2
	}; 
    \)\)

diffRule[x_] = Flatten[{diffRule0[x], diffRule1[x]}];

Lastly, we need a rule to differentiate inside the sums (it is assumed that the user knows when this is allowed), substituting for d when all of the previous properties have been applied. Note that the expressions that do not contain Summa object are differentiated by Mathematica in the usual way, thus applying every conceivable rule. This means that the rules given in diffRule2 are limiting only in the case of expressions containing sums.

\(\(diffRule2[x_] = 

      	{
	d[Summa[el_, it_]] :> Summa[D[el, x], it], 
	
        d[f_] :> D[f, x] /; \ FreeQ[f, Summa]
	}; \)\)

The procedure DSum, collecting all the steps, is given an option so that the user can choose whether to  bring inside  or not all  the factors premultiplying sums. Note that diffRule is applied recursively, while diffRule2 is not.

Clear[DSum]\(Options[DSum] = {Inside -> True}; 

    DSum[s_, x_, opts___] := 
	
      Module[
		{opt1}, 
		
        opt1 = \(Inside /. {opts}\) /. Options[DSum]; 
		
        \(d[If[opt1, BringIn[s], s]] //. diffRule[x]\) /. diffRule2[x]
	]
        \)

This overloads the procedure in order to handle equations directly

DSum[lhs_ == rhs_, x_, opts___] := 
      DSum[lhs, x, opts] == DSum[rhs, x, opts]

And this allows for multiple differentiation:

DSum[s_, {x_, n_}, opts___] := Nest[DSum[#, x, opts]&, s, n]

Now we can test it, bearing in mind that it was devised for basic manipulation only

\(Table[DSum[Summa[\(a\_k\) x^k, {k, 0, ∞}], {x, n}], 

        	{n, 0, 3}] // TableForm\)

Obviously, the first n terms are zero, as can be seen by

\(Table[
        BreakSum[
		
          DSum[Summa[\(a\_k\) x^k, {k, 0, ∞}], {x, n}], n], 

        	{n, 0, 3}] // TableForm\)

The following shows an application of the Inside option: setting Inside to False causes the product to be differentiated according to the rule D[f g]=D[f] g+f D[g]

s = x^2  Summa[(Subscript[a, k]) x^k, {k, 0, ∞}]

DSum[ s, x]

DSum[ s, x, Inside -> False]

We can recover the previous result by means of a suitable simplification

SimplifySum[BringIn[%]]

DSum can be applied to expression not containing sums... (in this case the Mathematica D operator will take care of the differentiation)

\(Clear[f, g, x]; 
DSum[f[x], x]\)

DSum[f[x] g[x], x]

...and to mixed expressions, even in the form of equations

expr = 
      Summa[(Subscript[a, k]) x^(k - 2) (x - Sin[x]), {k, 5, ∞}] + 
        2  (Subscript[a, 1 ])/x^3

DSum[expr, x]

expr1 = f[x] (g[x] + Summa[h[x, k], {k, 0, ∞}]) == 0

DSum[expr1, x]

expr2 = f[x] + g[x] Summa[h[x, k], {k, 0, ∞}] == 0

DSum[expr2, x]

espr = 
      Summa[((Subscript[a, k - 5]) + 2  (Subscript[a, k]) - (Subscript[a, k - 3]) - (Subscript[b, k - 4]) + 
                (Subscript[b, k - 3])) (-x^(k - 3)), {k, 1, 
            ∞}] - 2  (Subscript[a, 0])/(x^3) + (Subscript[a, 0]) - 2  (Subscript[a, 2])/x - 
        2  (Subscript[a, 1])/x + (Subscript[a, 1]) x - 2  (Subscript[a, 4]) x + x (Subscript[b, 0]) - (Subscript[b, 0]) - x (Subscript[b, 1]) - 
        2  (Subscript[a, 3])

DSum[espr, x]

DEachSum

DEachSum[expr]

In certain instances it is better to resort to a differentiation of what is inside the sums, by using the procedures DEachSum. It  operates on  sums only, by crudely applying the differentiation operator inside Summa.

DEachSum[s_, x_] := s /. Summa[el_, it_] :> Summa[D[el, x], it]\(DEachSum[s_, {x_, n_}] := 
      s /. Summa[el_, it_] :> Summa[D[el, {x, n}], it]
		\)

Basic term by term Integration

IntegrateSum[espr,x]
IntegrateSum[espr,{x,a,b}]

Option: Inside (True/False)

Following the definitions given to implement DSum, we can introduce the following rule

\(\(intRule1[x_Symbol | {x_Symbol, _, _}] = 

      	{
	in[f1_ + f2_] :> in[f1] + in[f2], 
	
        in[c_*f_] :> c\ in[f]\  /; \ FreeQ[c, x]\ , 
	
        in[f1_\ f2_] :> in[BringIn[Expand[f1\ f2]]]\ 
	}; \)\)

A second rule operate the actual integration, by moving the operator inside Summa objects and letting Mathematica integrate the non-Summa parts.

\(\(intRule2[x_] = 

      	{
	in[Summa[el_, it_]] :> Summa[Integrate[el, x], it], 
	
        in[f_] :> Integrate[f, x]
	}; \)\)

We can see that intRule1 does not perform any integration

espr = 
      3  x^2 - Sin[x] + 3  x y Summa[(Subscript[b, n]) x^n, {n, 2, ∞}]in[espr]  //. intRule1[x]

Since that is a job for intRule2.

(in[espr]  //. intRule1[x]) /. intRule2[x]

This is another expression, integrated

in[Summa[(Subscript[c, k]) x^k, {k, 0, ∞}]] //. intRule1[x] /. 
      intRule2[x]

Both rules can be put into a single procedure, allowing an option for taking outside the integral sign everything that does not depend on the variable of integration. As before, intRule1 is applied repeatedly while intRule2 is applied just once.

\(Options[IntegrateSum] = {Inside -> True}; 

    IntegrateSum[s_, x_, opts___] := 
	
      Module[
		{opt1}, 
		
        opt1 = \(Inside /. {opts}\) /. Options[IntegrateSum]; 
		
        If[opt1, 
			\(in[BringIn[s]]\  //. intRule1[x]\) /. intRule2[x], 
          
			BringOut[\(in[s] //. intRule1[x]\) /. intRule2[x]]
		]

        	]\)

Definite integration is allowed

\(f = Summa[\(a\_k\) x^\((2  k)\), {k, 0, ∞}]; 
    
{f, 
	IntegrateSum[f, x], 
	IntegrateSum[f, {x, 0, 1}], 
	
        IntegrateSum[f, {x, a, b}]} // ColumnForm\)

And, in version 3, we can get conditional integrals

IntegrateSum[espr, {x, 0, 1}]

% /. Re[n] > -2 -> True

Prior to integrating, factors premultiplying the sum are automatically brought inside

s = x y z^2 Summa[(Subscript[a, k]) z^k, {k, 0, Infinity}]

IntegrateSum[s,z]

Setting the Option Inside->False will prevent such a behaviour

IntegrateSum[s,z,Inside->False]

To force threading over equations we need only add this line of code:

\(\(IntegrateSum[lhs_ == rhs_, x_, opts___] := 
\ \ 
      IntegrateSum[lhs, x, opts] == IntegrateSum[rhs, x, opts]	
\)\)

Clear[f]IntegrateSum[
      f[x] Exp[I m x] == 
        Summa[(Subscript[c, n]) Exp[-I n x], {n, 0, ∞}] Exp[I m x], 
      x]

Other examples follow

Clear[f] s = f[x] Summa[(Subscript[a, k]) x^k, {k, 0, Infinity}]

IntegrateSum[s,x]

s = f[x] Summa[(Subscript[a, k]) g[x], {k, 0, ∞}]

IntegrateSum[s, x]

BringOut[%]

IntegrateSum[s, x, Inside -> False]

Selectively Applying and Mapping Functions (*)

These are untested and undocumented procedures that will be added to a later version of the package.

SumPosition

SumPosition

We might also want to map a function inside specific sums in a given expression. Suppose we want to shift the index of the second sum only in

expr = 
      Summa[(Subscript[a, k]) x^k, {k, 1, ∞}] + 3  Cos[x] - 2 == 
        Sqrt[3  x^2] Summa[Cos[k x], {k, 1, ∞}] - 
          2  Summa[(Subscript[b, k]) x^(k - 2), {k, 2, ∞}]

withoout having to change the index in the other sums appearing in expr. We need a procedure that can do what SumMap does but only at specific positions.

We have to redefine the meaning of position in expr to mean something simpler than what it really means to Mathematica. Position k means 'the k-th Sum appearing in expr". In this particular case, the Sum we want to act upon is that in 'summa position' 1. Note that we have to refer to the output given by Mathematica since the position in the fullform of the expression passed to it can be different from that input by the user. This is a side effect of the canonical reordering of functions, such as Plus, with attributes like Orderless. As we can see, the position of the Summa we want to zero in is

Position[expr, Summa[Cos[k x], {k, 1, ∞}]]

So, how can we identify its 'summa position' as 1? One way is to find the positions of every occurence of Summa in expr

Position[expr, _Summa]

Then we use the position of the element Position[expr,_Summa] in this table as the position for our sum. We need to look up what is the actual position associated with the ordinal number 1 in a table like this

{Range[Length[%]], %} // Transpose

hence the position identifying the k-th Sum is

SumPosition[expr_, k_] := Part[Position[expr, Summa[_, _]], k]

We can identify the position of the first sum only or of the first and third sum appearing in expr. Please note that the actual position is referred to the fullform of the expression, so it may necessary to have a look at FullForm[expr] to identify the exact place where we want to map our function.

SumPosition[expr, 1]SumPosition[expr, {1, 3}]

SumMap

SumMap

Sometimes it may be useful to apply a function to the element of a sum without having to worry about the Sum(ma) wrapper., The procedures SumMap and SumApply implement a simple mapping scheme to do just that. The difference with respect to Map is that without a level specification (or with a specification of 0) SumMap will wrap the function around the element of every Summa object contained in the expression passed to it.

SumMap[f_, expr_] := 
      expr /. Summa[el_, it_] :> Summa[f@el, it]SumMap[f_, expr_, 0] := SumMap[f, expr]SumMap[f_, expr_, n_] := 
      expr /. Summa[el_, it_] :> Summa[Map[f, el, n], it]

It may be useful to perform different kind of transfrormations on the elements of the sum. For example, we can map a TrigToExp transformation inside every Sum

expr = 
      Summa[(Subscript[a, k]) x^k, {k, 1, ∞}] + 3  Cos[x] - 
        Sqrt[3  x^2] Summa[Cos[k x], {k, 1, ∞}]

SumMap[TrigToExp, expr]

Collect, Apart, and other simplifications functions can be wrapped on the elements, mapped and applied inside every sum.

SumMapAt

MapAt[f,expr,SumPosition[expr,n]]
SumMapAt

We might also want to map a function inside specific sums in a given expression. Suppose we want to shift the index of the second sum only in

expr = 
      Summa[(Subscript[a, k]) x^k, {k, 1, ∞}] + 3  Cos[x] - 2 == 
        Sqrt[3  x^2] Summa[Cos[k x], {k, 1, ∞}] - 
          2  Summa[(Subscript[b, k]) x^(k - 2), {k, 2, ∞}]

withoout having to change the index in the other sums appearing in expr. We need a procedure that can do what SumMap does but only at specific positions. This is a task that in ordinary Mathematica world is accomplished by the function MapAt

?MapAt

All we have to do now is to map our functions at the specific positions given, remembering that we need to map them inside the s

MapAt[f, expr, SumPosition[expr, {1, 3}]]

SumMapAt[f_, expr_, pos_] := MapAt[f, expr, SumPosition[expr, pos]]SumMapAt[f_, expr_, pos_, lev_] := 
      MapAt[Map[f, expr, lev], expr, SumPosition[expr, pos]]

SumMapAt[f, expr, 2]

x + s - 2

SumApply e SumApplyAt

SumApply

SumApply[f_, expr_, n_:  1] := 
      expr /. Summa[el_, it_] :> Summa[Apply[f, el, n], it]

Further Enhancements (*)

Adding new names and options to the procedures

Version 2.3 adds a new naming for Summa.m's procedures. Most of them now start with the string Sum. The older names are retained because the author is unbelievalby lazy and wanted to make sure that older notebooks would work without having to change any name. The following shows a correspondence between new and old names. (The actual package has the equalities reversed since the naming was changed in the body of the package)

\(GetExponent; 
SumToTerm = GetElement; 

    SumCoefficient = GetCoefficient; 
    
SumExpand = ExpandSum; 
SumExpandAll = ExpandAllSum; 
    
SumSimplify = SimplifySum; 
SumFullSimplify = FullSimplifySum; 

    SumSimplifyEach = SimplifyEachSum; 

    SumFullSimplifyEach = FullSimplifyEachSum; 
BringIn; 

    BringOut; 
SumSplit = BreakSum; 
    
SumReindex = ReindexSum; 
SumIndexAlign = AlignIndex; 
    
SumExponentAlign = AlignExponent; 
SumD = DSum; 

    SumDEach = DEachSum; 
    \  (*\ to\ be\ embedded\ in\ SumD\ with\ option\ SeparateSums -> 
        True*) 
SumIntegrate = IntegrateSum; \)

Moreover, in order to give to the manipulation procedures the same options of the built-in procedures they are extending to sum manipulation, part of the code of the package had to be modified. Since most procedurs in Summa.m are rule based, it has been necessary to modify all of the rules that called the procedure with options. Here it is shown how to fix just two procedure: BringIn and SumSimplify.

\(inrule1[opts___] = 
      c_\ Summa[el1_, it : {k_, _, _}] :> 
        Summa[If[$Explode, Expand[c\ el1, FilterOptions[Expand, opts]], 
              c\ el1], it] /; FreeQ[c, k]; 

    inrule2[opts___] = 
      Summa[el1_, it_] + Summa[el2_, it_] :> 
        Summa[If[$Semplifica, 
            Simplify[el1 + el2, FilterOptions[Simplify, opts]], el1 + el2], 
          it]; 
inrule3[x_] = 
      Summa[el_, {k_, lim__}] :> 
        \((1/x)\)*Summa[x*el, {k, lim}] /; FreeQ[x, k]; 


    BringIn[s_, opts___] := s /. inrule1[opts]\)\(BringIn[s_, x_] := s /. inrule3[x]
\)SumSimplify[s_, opts___] := 
      s /. inrule1[opts] //. inrule2[opts] /. 
        Summa[0, _] -> 0

As we can see, not every rule need to be given options. Only those calling procedures that needs options.BringIn will accept the same options of Expand. In the particular case of inrule2 called by SumSimplify there are two built-in procedures carrying options: one is Simplify, and that was to be expected; the other one is Expand. By giving options to Expand too, and by selecting only valid options for each procedure through the use of FilterOptions, we can specify both kind of options in our call to SumSimplify.

Version 2.3 of the package does not allow the user to specify default values for the options associated to any particular Summa.m procedure. It is not difficult to code such an extension though. But since version 4 will need to associate such options to the built-in procedures themselves, the options are left 'read-only'.

Options[BringIn] = Options[Expand]Options[SumSimplify] = Union[Options[Simplify], Options[Expand]]

Other procedures, such as SumIntegrate, have options on their own. In this case they must be added to those of the built-in procedures they are internally using. FilterOptions is needed in order to avoid message errors

\(intRule1[x_Symbol | {x_Symbol, _, _}] = {
        in[f1_ + f2_] :> in[f1] + in[f2], 
        in[c_\ f_] :> c\ in[f] /; FreeQ[c, x], 
        in[f1_\ f2_] :> in[BringIn[Expand[f1\ f2]]]}; 

    intRule2[x_, opts___] = {
        in[Summa[el_, it_]] :> Summa[Integrate[el, x, opts], it], 
        in[f_] :> Integrate[f, x, opts]\  (*
          \  /; FreeQ[f, Summa]*) }; 


    SumIntegrate[s_, x_, opts___] := 
      Module[
{opt1, opt2}, 
        
opt1 = 
          \(Inside /. {opts}\) /. Options[SumIntegrate]; 
        
opt2 = FilterOptions[Integrate, opts]; 
        
If[opt1, 
          
\(in[BringIn[s]] //. intRule1[x]\) /. 
            intRule2[x, opt2], 
          
BringOut[
            \(in[s] //. intRule1[x]\) /. intRule2[x, opt2]]
            
]
]\)

Options[SumIntegrate] = 
      Union[{Inside -> True}, Options[Integrate]];

There's a little note to be made, regarding adding options: the automatic reordering of definitions performed by Mathematica is prevented by the specification of the pattern for options. If we simpy add opts___, function calls that are more specific are automatically placed before the more general ones. See for example:

Clear[f]f[expr_, x_, opts___] := FOO[expr, x]f[lhs_ == rhs_, x_, opts___] := 
      Equal[foo[lhs, x], bar[rhs, x]]

?f

f[Sin[t] == Cos[y], x]

But if we specify that the options have to be Rules, then Mathematica won't sort the definitions anymore, leaving them in the order the user entered them:

Clear[f]f[expr_, x_, opts___Rule] := FOO[expr, x]f[lhs_ == rhs_, x_, opts___Rule] := 
      Equal[foo[lhs, x], bar[rhs, x]]

?f

f[Sin[t] == Cos[y], x]

Of course there is any easy way out: specifying the right order in the package

Clear[f]f[lhs_ == rhs_, x_, opts___Rule] := 
      Equal[foo[lhs, x], bar[rhs, x]]f[expr_, x_, opts___Rule] := foo[expr, x]

?f

f[Sin[t] == Cos[y], x]

Since in the actual package the type Rule is specificed for the added options, this 'manual' reordering was necessary for the differentiation and integration procedures.

Enhancing built-in procedures (*)

To be added in a later, major, version.