Memory management with petsc4py

Linear Algebra
1

Hi everyone,

I am having a problem with freeing up memory in legacy fenics. As discussed here, with newer versions of petsc4py one should call the destroy method of the objects in order to free up memory ( garbage collection different between petsc4py 3.17.4 and 3.18.3 (#1309) · Issues · PETSc / petsc · GitLab ). However, when I assemble a matrix or vector into the fenics classes, access the underlying petsc4py object with .mat() or .vec() (optionally do some calculations, solve a linear system), and then try to free up the memory by calling A.mat().destroy(), the memory does not get freed and the object is still in memory.

I have the following MWE which shows the problem:

from fenics import *
from petsc4py import PETSc

mesh = UnitSquareMesh(8, 8)
comm = mesh.mpi_comm()

V = FunctionSpace(mesh, "CG", 1)
u = TrialFunction(V)
v = TestFunction(V)

a = dot(grad(u), grad(v)) * dx + u * v * dx
L = Constant(1.0) * v * dx

### Matrices
## First variant
A_fenics = assemble(a)
A_petsc = as_backend_type(A_fenics).mat()

as_backend_type(A_fenics).mat().destroy()
PETSc.garbage_cleanup(comm)
print(A_petsc.getSizes())  # expected: Segfault


## Second variant
A_fenics = PETScMatrix(comm)
A_petsc = A_fenics.mat()
assemble(a, tensor=A_fenics)

A_fenics.mat().destroy()
PETSc.garbage_cleanup(comm)
print(A_petsc.getSizes())  # expected: Segfault

### Vectors
## First variant
b_fenics = assemble(L)
b_petsc = as_backend_type(b_fenics).vec()

as_backend_type(b_fenics).vec().destroy()
PETSc.garbage_cleanup(comm)
print(b_petsc.getSize())  # expected: Segfault

## Second variant
b_fenics = PETScVector(comm)
b_petsc = b_fenics.vec()
assemble(L, tensor=b_fenics)

b_fenics.vec().destroy()
PETSc.garbage_cleanup(comm)
print(b_petsc.getSize())  # expected: Segfault

The memory only gets freed if I call A_petsc.destroy() or b_petsc.destroy(). This is against my understanding of the code, as .mat() or .vec() should return a pointer to the petsc4py object. And if I modify, e.g., the values of the objects, then I see that these changes are reflected in both A_petsc and A_fenics, or in b_petsc and b_fenics, respectively.

Does anyone have a clue to what is going on here and how the memory can be freed correctly?
For context: Of course I don’t want to generate a segfault in my program, but for this test it is an easy way to check whether the memory is still allocated or actually free.

Thanks a lot in advance,
Sebastian

2

I would suggest rolling back to an earlier version of PETSc, prior to the introduction of manual memory management

as legacy DOLFIN has not been maintained to keep up with modern PETSc changes. 3.17 came out in 2022 (Changes: 3.17 — PETSc v3.23.6-578-gf29e3a7d8e99 documentation), which is ~3 years after active development of Legacy DOLFIN stopped.

For instance .mat() returns a PETScBaseMatrix, not a PETSc4py matrix directly (Bitbucket).
I guess there is a magic cast for it through: Bitbucket
as alot of legacy FEniCS was prior to the fully fledged PETSc4py interface.

3

Thanks for your answer.

I am still not sure what is happening in the background, but I could resolve my actual issue with calling both the .destroy() methods of the solvers (and matrices / vectors) that were not needed anymore in addition to PETSc.garbage_cleanup() after the solves were done - even with newer versions of PETSc (I am using 3.20 currently).

So the memory is indeed freed up, even though the MWE above suggests otherwise. Maybe just the pointers still point to the object in memory so it can be accessed, but will be overwritten eventually?