Question: I have a fairly large number (around 1000) of step functions, each with only two intervals. I'd like to sum them up and then find the maximum value. What is the best way to do this? I've tried out sympy, with code as follows:

```
from sympy import Piecewise, piecewise_fold, evalf 
from sympy.abc import x
from sympy.plotting import *
import numpy as np

S = 20

t = np.random.random(20)

sum_piecewise = None

for s in range(S):
    p = Piecewise((np.random.random(), x<t[s]), (np.random.random(), x>=t[s]))
    if not sum_piecewise:
        sum_piecewise = p 
    else:
        sum_piecewise += p

print sum_piecewise.evalf(0.2)

```

However, this outputs a large symbolic expression and not an actual value, which is what I want.

Answer: 
As it appears that you consider numerical functions, it is better (in terms of performance) to work with Numpy. Here's one approach:

```
import numpy as np
import matplotlib.pyplot as plt

np.random.seed(10)
S = 20 # number of piecewise functions

# generate S function parameters. 
# For example, the k-th function is defined as equal to 
# p_values[k,0] when t<t_values[k] and equal to
# p_values[k,1] when t>= t_values[k]
t_values = np.random.random(S)
p_values = np.random.random((S,2))

# define a piecewise function given the function's parameters
def p_func(t, t0, p0):
    return np.piecewise(t, [t < t0, t >= t0], p0)

# define a function that sums a set of piecewise functions corresponding to
# parameter arrays t_values and p_values
def p_sum(t, t_values, p_values):
    return np.sum([p_func(t, t0, p0) for t0, p0 in zip(t_values,p_values)])

```

Here is the plot of the sum of functions:

```
t_range = np.linspace(0,1,1000)
plt.plot(t_range, [p_sum(tt,t_values,p_values) for tt in t_range])

```

[![enter image description here](https://i.stack.imgur.com/XMTrf.png)](https://i.stack.imgur.com/XMTrf.png)

Clearly, in order to find the maximum, it suffices to consider only the `S` time instants contained in `t_values`. For this example,

```
np.max([p_sum(tt,t_values,p_values) for tt in t_values])

```

> 
> `11.945901591934897`
> 
> 
>
What about using [substitution](http://docs.sympy.org/latest/tutorial/basic_operations.html#substitution)? Try changing `sum_piecewise.evalf(0.2)` by `sum_piecewise.subs(x, 0.2)`