Error after repeatedly create solver

Hi, I’m trying to run an elliptic equation solver with random coefficients for many times (to sample the coefficient). But I got the following error. I saw the same error in this post, but I still do not know the solution. Error setting PETSc options repeatedly.

Attached a mini example of what I have. It get to this error after about 170 times.

Traceback (most recent call last):
File “/home/shared/mini.py”, line 57, in
solver.solve(sample)
File “/home/shared/mini.py”, line 47, in solve
problem = fem.petsc.LinearProblem(a, L, bcs=self.bcs, petsc_options={“ksp_type”: “preonly”, “pc_type”: “lu”, “pc_factor_mat_solver_type”: “mumps”})
File “/usr/lib/petsc/lib/python3/dist-packages/dolfinx/fem/petsc.py”, line 589, in init
opts[k] = v
File “PETSc/Options.pyx”, line 23, in petsc4py.PETSc.Options.setitem
File “PETSc/Options.pyx”, line 91, in petsc4py.PETSc.Options.setValue
petsc4py.PETSc.Error: error code 55

Primary job terminated normally, but 1 process returned
a non-zero exit code. Per user-direction, the job has been aborted.


mpirun detected that one or more processes exited with non-zero status, thus causing
the job to be terminated. The first process to do so was:

Process name: [[50127,1],0]
Exit code: 1

google drive with the mini example. mini_example - Google Drive

A very simple rewrite of the problem avoids this error and speeds it up massively:

import numpy
import ufl

from mpi4py import MPI
from dolfinx import mesh
from dolfinx import fem
from dolfinx.fem import FunctionSpace
from petsc4py.PETSc import ScalarType


class EllipticRandomCoefSolver:
    def __init__(self, comm):
        self.comm = comm

        self.domain = mesh.create_unit_square(
            comm, 64, 64, mesh.CellType.quadrilateral)
        self.tdim = self.domain.topology.dim

        self.V = FunctionSpace(self.domain, ("CG", 1))
        dofs_left = fem.locate_dofs_geometrical(
            self.V, self.dirichlet_boundary_left)
        dofs_right = fem.locate_dofs_geometrical(
            self.V, self.dirichlet_boundary_right)

        bc_left = fem.dirichletbc(ScalarType(0.0), dofs_left, self.V)
        bc_right = fem.dirichletbc(ScalarType(1.0), dofs_right, self.V)
        self.bcs = [bc_left, bc_right]
        u = ufl.TrialFunction(self.V)
        v = ufl.TestFunction(self.V)
        x = ufl.SpatialCoordinate(self.domain)

        f = ufl.cos(2*numpy.pi*x[0])*ufl.sin(2*numpy.pi*x[1])
        g = fem.Constant(self.domain, ScalarType(0))
        self.sample = fem.Constant(self.domain, ScalarType(0))
        k = self.sample*ufl.cos(2*numpy.pi*x[0])*ufl.sin(2*numpy.pi*x[1])+2.0

        a = ufl.dot(k*ufl.grad(u), ufl.grad(v)) * ufl.dx
        L = f * v * ufl.dx - g * v * ufl.ds
        self.problem = fem.petsc.LinearProblem(a, L, bcs=self.bcs, petsc_options={
            "ksp_type": "preonly", "pc_type": "lu", "pc_factor_mat_solver_type": "mumps"})

    @staticmethod
    def dirichlet_boundary_left(x):
        return numpy.isclose(x[0], 0)

    @staticmethod
    def dirichlet_boundary_right(x):
        return numpy.isclose(x[0], 1)

    def solve(self, sample):
        self.sample.value = sample
        self.uh = self.problem.solve()


if __name__ == "__main__":
    comm = MPI.COMM_SELF
    solver = EllipticRandomCoefSolver(comm)
    samples = numpy.random.uniform(0, 1, 1000)
    for idx, sample in enumerate(samples):
        print(idx)
        solver.solve(sample)

In general, you should avoid compiling forms inside of loops.

Thank you dokken, problem solved :slight_smile: