Memory usage using Assemble

1

Hello everyone,

I am a beginner with FEniCS and I’ve tried to use it for couple months to develop a topology optimization based framework for an electromagnetic problem. I currently have a version of my code which works fine, however it seems that I have a problem of memory leak in the function where I compute the FEA and my objective function, at each iteration step. Unfortunately, it keeps my algorithm from reaching convergence…

As pointed out in this discussion : Is there memory leak in FEniCS and how to fix it? - #3 by bmf It seems that the problem comes from the fact that assemble create a new tensor each time it is called, and I call it several times since I am using finite differences to perform my sensitivity analysis. By looking at the assemble() documentation I also saw that I should use the tensor argument to avoid creating a new one each time it is called, like this :

A = PETScMatrix()
assemble(a, tensor=A)

However, doing so makes my jupyter kernel crash without any error. Here is a minimal example which reproduces the error, where bcs are some DC boudary conditions :

f = interpolate(Constant(1.0e-5), V)  # the source term for the PDE for elec
u = TrialFunction(V)
v = TestFunction(V)
s = Function(V)
i = 0
Q = PETScMatrix()
while i <10000 :
    L = f*v*dx # ~= 0
    F = inner(grad(v),-grad(u))*dx
    solve(F == L, s, bcs)
    Q = assemble(s*dx, tensor = Q)
    i +=1

Can anyone help me with this problem ?
Thanks in advance,
Pierre

2

s * dx is a rank 0 form and will be assembled as a scalar. Attempting to assemble it into a matrix is likely causing the error.

1 Like
3

Thank you for your answer, indeed I didn’t even thought about that :sweat_smile:, however there is a small mistake in my MWE and the result of the assemble should be a scalar for example Q should be define as :
Q = assemble(s*s*dx, tensor = Q)

But I have no idea on how to define a scalar valued tensor, could you help me on this please ?

4

Performance will greatly improve by using preallocation of vectors and matrices.

However, scalars are cheap, there’s little to no reason to preallocate a double precision float for your form. I’d just use Q = assemble(s*s*dx) for rank 0 forms.

You can probably eek out a little more performance by preallocating the UFL form though. Symbolic algebra is expensive. E.g.

# ... do stuff
interesting_form = s * s * dx
# ... do more stuff
for t in time_loop:
    # do time loop things ...
    Q = assemble(interesting_form)

Also, a general comment on performance: solve() is very expensive if you’re calling it every step in your loop. You’re creating and destroying all the necessary datastructures with each call. Consider LinearVariationalProblem or setting up your linear algebra system before the loop.

1 Like
5

In fact, the problem was coming from the fact that I was unintentionally using the assemble() version from dolfin-adjoint and not the one from dolfin, it was therefore keeping a copy at each iteration (and there are a lot of iterations :grin:). Now everything seems to be working fine without an excessive memory usage.

Thank you for your help ! And I will have a look at LinearVariationalProblem it could indeed help fasten my convergence which is kinda long, thanks for the tips !

1 Like