Maximum number of terms in the definition of the bilinear form?

Question

Hello!

I have been using FEniCS to solve the following linear system of PDEs:

Ax \partial_{x} \psi + Ay \partial_{y} \psi + \sigma \psi= Q

Where Ax and Ay are square matrices of size L.

I wrote the variational form using a mixed formulation which I coded as follows (the point of my question is not how to code the bilinear form so don't worry about understanding all the terms: the only important thing is that I do: a += ... a bunch of times):

# Define function spaces and mixed (product) space
DG = FunctionSpace(mesh, "DG", 1)
W = MixedFunctionSpace([DG]*L0)

# Define trial and test functions
psi = TrialFunctions(W)  # solution
v   = TestFunctions(W)   # test function
...
for i in range(0, L):
    for j in range(0, L):

        if abs(Ax[i][j]) > DOLFIN_EPS or abs(Absx[i][j]) > DOLFIN_EPS:
            a += -Ax[i][j]*psi[j]*v[i].dx(0)*dx
            a += (Ax[i][j]*avg(psi[j]) + 0.5*Absx[i][j]*jump(psi[j]))*jump(v[i])*abs(dot(ex,n("+")))*dS
            a += 0.5*(Ax[i][j]*psi[j]*(dot(ex,n))+ Absx[i][j]*psi[j]*abs(dot(ex,n)))*v[i]*ds

        if abs(Az[i][j]) > DOLFIN_EPS or abs(Absz[i][j]) > DOLFIN_EPS:
            a += -Az[i][j]*psi[j]*v[i].dx(1)*dx
            a += (Az[i][j]*avg(psi[j]) + 0.5*Absz[i][j]*jump(psi[j]))*jump(v[i])*abs(dot(ez,n("+")))*dS
            a += 0.5*(Az[i][j]*psi[j]*(dot(ez,n))+ Absz[i][j]*psi[j]*abs(dot(ez,n)))*v[i]*ds

a += (Sigma_t+1/(c*dt))*psi[i]*v[i]*dx
...

When the size of Ax and Ay is small (typically 16x16), it works great.
But when I increase their size to a larger number (typically 36x36), I get the following error:

RuntimeError: maximum recursion depth exceeded while calling a Python
object

So my question is: is there a maximum number of terms that can be added to the bilinear form? It really surprises me because those matrices above are pretty sparse and the number of terms I actually add is typically around 200 which does not seem to be unreasonable to me. Is it?

Thank you so much!
Vincent

timm · Answer 1 · Apr 30, 2014

Best answer

You can change the recursion limit in Python, I'd go easy on it … I imagine that it is possible to hang things up. On my system, the default value for the recursion limit is 1000.

import sys
sys.getrecursionlimit()
sys.setrecursionlimit(something_large)

However, I'm not so sure that the above will fix your problem … I tried to make a simple recursive program in Python and got a slightly different error:

def factorial(n):
  if n == 1:
    return 1
  else:
    return n*fact(n - 1)

factorial(1000)
RuntimeError: maximum recursion depth exceeded

answered Apr 30, 2014 by timm FEniCS User (2,100 points)
selected Apr 30, 2014 by V_L

Yes that works just fine! Thanks a lot for this quick answer!
I'm pretty new to Python and my google search hadn't been successful... :)

commented Apr 30, 2014 by V_L FEniCS User (4,440 points)

Maximum number of terms in the definition of the bilinear form?

1 Answer