### What's meaning of the error message, "The provided data did not satisfy the prerequisites"?

289
views
1
10 months ago by
I have two questions: 1.why I get a NaN?
2. After uncommenting, What's meaning of the error message, "The provided data did not satisfy the prerequisites"?

When comment
result1 = mysolver_1(mesh, f1)​

, I get the result and get a NaN. Or, I get the below error:
*** Error: Unable to compute matrix factorisation.
*** Reason: The provided data did not satisfy the prerequisites.
*** Where: This error was encountered inside EigenLUSolver.cpp.
*** Process: 0
***
*** DOLFIN version: 2017.1.0
*** Git changeset: unknown

Here is my simple source code. Forgive me it is also long, but you don't need to understand the whole code.

import pdb
from fenics import *
from mshr import *
import numpy as np
from pprint import pprint
from scipy import sparse

# form source terms
def all_Expression():
cell = 'tetrahedron'
import sympy as sp
x, y, z = sp.symbols('x[0] x[1] x[2]')
lamda = 1.0
p = x - 0.5
U = (x*(x-1) * y*(y-1) * z*(z-1) * (y-0.5)*(z-0.5))**3
p_x = p.diff(x, 1)
p_y = p.diff(y, 1)
p_z = p.diff(z, 1)
U_x = U.diff(x, 1)
U_y = U.diff(y, 1)
U_z = U.diff(z, 1)
Laplace_U = sum(U.diff(xi, 2) for xi in (x, y, z))
Laplace_Laplace_U = sum(Laplace_U.diff(xi, 2) for xi in (x, y, z))
#---------------source terms-----------------------
p_e = Expression(sp.printing.ccode(p), degree=1, cell=cell)
f1 = Expression((sp.printing.ccode(- p_x), \
sp.printing.ccode(- p_y), \
sp.printing.ccode(- p_z)), degree=1, cell=cell)
f2 = Expression(sp.printing.ccode(- Laplace_Laplace_U - lamda - p), degree=4, cell=cell)

return p_e, f1, f2

def mysolver_1(mesh, f):
parameters.linear_algebra_backend = 'Eigen'
h = mesh.hmax()

DG = FunctionSpace(mesh, 'DG', 0)
RT = FunctionSpace(mesh, 'RT', 1)

p = TrialFunction(DG)
v = TestFunction(RT)

bc = DirichletBC(RT, Constant([0.0, 0.0, 0.0]), DomainBoundary())

a = p * div(v) * dx
L = inner(f, v) *dx
a_ = p * dx

A, b = assemble_system(a, L, bc)
A_ = assemble(a_)

rows, columns, values = as_backend_type(A).data()
csrA = sparse.csr_matrix((values, columns, rows))
AA = sparse.vstack([csrA, A_.array().T])
b_ = list(b.array()); b_.append(0)
bb = np.array(b_)

from scipy.sparse.linalg import lsqr
pp = lsqr(AA, bb)[0]
p = Function(DG)
p.vector().set_local(pp)

plot(p, title='solution p')
return p, h

def mysolver_2(mesh, f, p):
h = mesh.hmax()

CG1 = FiniteElement('CG', mesh.ufl_cell(), 1)
R = FiniteElement('R', mesh.ufl_cell(), 0)
element = MixedElement([CG1, R])
V = FunctionSpace(mesh, element)

bc = DirichletBC(V.sub(0), Constant(0.0), DomainBoundary())

phi, lamda = TrialFunctions(V)
psi, mu = TestFunctions(V)

a = inner(grad(phi), grad(psi))*dx - lamda*psi*dx \
- phi*mu*dx
L = p*psi*dx + f*psi*dx

U = Function(V)
solve(a == L, U, bc)

phi, lamda = U.split()
plot(phi, title='phi')
plot(lamda, title='lamda')
return phi, lamda, h

if __name__ == '__main__':
p_e, f1, f2 = all_Expression()

domain = Box(Point(0,0,0), Point(1,1,1)) - Box(Point(0,0,0), Point(1,0.5,0.5))
mesh = generate_mesh(domain, 1)

pdb.set_trace()
for i in range(3):
#----------------------------------------error here----------------------------------------------------
result1 = mysolver_1(mesh, f1)

result2 = mysolver_2(mesh, f2, p_e)

mesh = refine(mesh)

interactive()


Community: FEniCS Project
written 10 months ago by Garth Wells
I am sorry for the inappropriate question.
Now I make the question clear.
Thanks a lot.
written 10 months ago by fanronghong

2
10 months ago by

The problem actually occurs for me in solver 2, rather than solver 1.

I believe that the reason is that the mesh you generate with "generate_mesh(domain, 1)" is so coarse that every node is on the boundary, so they are all constrained by DirichletBC. This results in a locking for the CG1-R mixed space and thus a singular system.