Symbolic conditional sum

Is there a way to calculate, or at least state properly, a symbolic sum which iterates over a sequence conditionally?

For example, if a(n) is a sequence (of reals) and c(n) is a "condition" (boolean function defined on the integers), then I wish to compute, or at least state, the sum over all a(n)'s for which c(n).

Formally, I'd like to do something like this:

n = Symbol('n', integer=True, positive=True)
a = 1 / n**2
c = Eq(n%3, 1)
## s = Sum(a, (n, 0, oo), condition=c)

So s should be the following sum:

1/1 + 1/16 + 1/49 + 1/100 + 1/169 + ...

Of course, in this case I can define s manually, but I wonder whether I can do it automatically somehow.

Edit:

By manually I mean

Sum(1/(3*n+1)**2, (n, 0, oo))

This can be evaluated (and it is about 1.12173301393634).

An attempt I've made (which failed) is as follows:

Sum(Eq(n%3, 1) * (1/n**2), (n, 0, oo))

Trying to evaluate this using .evalf() resulted

AttributeError: 'BooleanFalse' object has no attribute '_eval_evalf'

Edit (2):

Here's another attempt, which yields a wrong result:

p = Piecewise((1/(n**2), Eq(n%3, 1)), (0, True))
Sum(p, (n, 0, oo)).evalf()
## returns 1.00000000000000

Either I've done something wrong, or this is a sympy bug.

Solution

Using a Piecewise is definitely the correct way to do this. The wrong result is a bug in SymPy, which you should report at https://github.com/sympy/sympy/issues/new. Your second method won't work because in SymPy booleans are not implicitly integers (True and False are not 1 and 0).