Problem with project in 3D

Hi everyone,

I have developed a model in 2D CFD model in FEniCS which works fine, now I am trying it on 3D geometries. However, in 3D I get an error in the following lines of script, which has to do with the project() function:

velocity_degree = 2
V_X2V= VectorFunctionSpace(mesh, 'CG', velocity_degree)

SM = Function(V_X2V)
SM0 = Constant((0,)*mesh.geometry().dim())

# if tstep = 1
    body_force = project(SM0, V_X2V)
# else
    body_force = project(SM, V_X2V)

If tstep is not 0, SM is a function that varies over the nodes. For both the project() lines I get the following error:

Error: Unable to successfully call PETSc function 'KSPSolve'
Reason: PETSc error code is: 76 (Error in external library)
Where: This error was encountered inside /build/dolfin-fyfCiA/dolfin-2017.2.0.post0/doolfin/la/PETScKrylovSolver.cpp

I already tried to change project() to interpolate(), but this only works for tstep=0. For body_force = interpolate(SM, V_X2V) I get:

TypeError: in method 'Function_interpolate', argument 2 of type 'dolfin::GenericFunction const &'

Does anyone know what I am doing wrong, and how to fix/workaround this?

Thanks in advance!

It could be that the solver is running out of memory in 3D. I believe the default is a direct solver, which I try to avoid in 3D. Try specifying an iterative solver, e.g.,

from dolfin import *

N = 10
mesh = UnitCubeMesh(N,N,N)
V = VectorFunctionSpace(mesh,"CG",1)
d = mesh.geometry().dim()

u = project(Constant(d*(1,)),V,
            solver_type="cg",preconditioner_type="jacobi")
1 Like