Facing an error after running the code

[A123-bansal123]

I am getting the following error after running my code, and I have tried different examples, but I still get the same error. Please help me in resolving this error.

[sghelichkhani]

Have you tried export OMP_NUM_THREADS = 1 before running your code, as suggested by the Warning message?

[A123-bansal123]

yes

[dham]

Can we have a minimal failing example please? I.e. a short piece of Python code that exhibits the error, plus whatever command you are using to run it.

It looks somewhat like you are using threads at the Python level (which is different from the OMP_NUM_THREADS thing). That isn't going to work. If you want to run Firedrake in parallel then you need to use MPI.

[A123-bansal123]

I remove some terms and try to solve simple Navier-Stokes equation code but still it shows some error:

error :

File "NS2.py", line 31, in <module>
    solve(F == 0, up, bcs = bcs)
  File "PETSc/Log.pyx", line 115, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "PETSc/Log.pyx", line 116, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "/home/aparna/software/firedrake/src/firedrake/firedrake/adjoint/solving.py", line 53, in wrapper
    output = solve(*args, **kwargs)
  File "/home/aparna/software/firedrake/src/firedrake/firedrake/solving.py", line 129, in solve
    _solve_varproblem(*args, **kwargs)
  File "/home/aparna/software/firedrake/src/firedrake/firedrake/solving.py", line 168, in _solve_varproblem
    problem = vs.NonlinearVariationalProblem(eq.lhs, u, bcs, J, Jp,
  File "PETSc/Log.pyx", line 115, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "PETSc/Log.pyx", line 116, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "/home/aparna/software/firedrake/src/pyadjoint/pyadjoint/tape.py", line 99, in wrapper
    return function(*args, **kwargs)
  File "/home/aparna/software/firedrake/src/firedrake/firedrake/adjoint/variational_solver.py", line 27, in wrapper
    self._ad_adj_F = adjoint(dFdu)
  File "PETSc/Log.pyx", line 115, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "PETSc/Log.pyx", line 116, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func
  File "/home/aparna/software/firedrake/src/firedrake/firedrake/ufl_expr.py", line 250, in adjoint
    return ufl.adjoint(form, reordered_arguments)
  File "/home/aparna/software/firedrake/src/ufl/ufl/formoperators.py", line 144, in adjoint
    return compute_form_adjoint(form, reordered_arguments)
  File "/home/aparna/software/firedrake/src/ufl/ufl/algorithms/formtransformations.py", line 467, in compute_form_adjoint
    raise ValueError("Expecting bilinear form.")
ValueError: Expecting bilinear form.

[dham]

This is an unrelated issue to your first error. The problem here is the third to last line where you redefine up. This means that the variable you are solving for (the new up) doesn't appear in the Form in your equation. If you remove the third to last line the issue goes away.

As an aside, using a Lagrange multiplier to deal with the pressure singularity is usually a bad idea. It's much better to tell the solver about the nullspace. See: https://www.firedrakeproject.org/solving-interface.html#solving-singular-systems

[A123-bansal123]

Sir, I tried different solutions, but they didn't converge.

Secondly, I want to ask if there's an alternate way to write 'f' in terms of 'u_exact' and 'p_exact.' In FEniCS, we used SymPy functions. Can we use similar functions in Firedrake?

error:

Error norm based on H_1 norm of u
0.125 -0.008220155056495093
0.0625 -0.0010646453333461113
0.03125 -9.314984231669773e-05
0.015625 -7.083432977015917e-06

Error norm based on L_2 norm of p
0.125 0.0691566819597118
0.0625 0.004536742939314852
0.03125 0.0003210832206884258
0.015625 2.5450133532849243e-05

Error norm based on L_2 norm of u
0.125 -0.019223186690572276
0.0625 -0.0015972680133422314
0.03125 -0.0001184367165971343
0.015625 -8.209070626590626e-06

[A123-bansal123]

Agree with you, but it is just a code for the simple Navier-Stokes equation, and I don't know where I made a mistake in formulation, as this formulation is well-known and commonly used. As this is a nonlinear part, there might be a mistake in solve(). Please provide some direction if possible.

I think that expressions are not supported in Firedrake, so instead, we use interpolate. Could you provide an example of generating expressions by passing them to Python's eval()?

Thank you for you response and help.

[dham]

Well we don't have a separate class for expressions, because you can just use UFL expressions. You don't actually need to interpolate.

I'm not sure what you mean by a mistake in solve. If you mean that you suspect a bug in solve, that's possible but exceptionally unlikely because the testing of solve is very extensive. I would need to see a specific example where you demonstrate solve getting the answer wrong (you would easily be able to show this by demonstrating that the residual is nonzero).

I don't have a sympy example to hand but what I had in mind was simply using sympy, and running eval on the string representation of the result. I have been reminded that a more robust thing to do is to use the .subs method on the sympy expression to replace sympy variables with the corresponding Firedrake symbols.

[A123-bansal123]

okay Thankyou so much for your response

I also try this example, this one is somehow easier than previous one and here error doesn't decreases as the mesh becomes finer and finer. I don't understand why this happen.
In this example when nx=4, error = 0.5250143495193417
nx=8, error = 0.2946536881739494
nx=16, error = 0.28247202758183304
nx=32, error = 0.2820600180290477
nx=64, error = 0.2820427874418266

from firedrake import *
nx=16
M = UnitSquareMesh(nx, nx)
# Define Function Spaces
V = VectorFunctionSpace(M, "Lagrange", 2)
W = FunctionSpace(M, "Lagrange", 1)
LM = FunctionSpace(M, "Real" , 0)
Z = V * W * LM
# Define Test and trial function
up = Function(Z)
u, p, rho = split(up)
v, q, lamda = TestFunctions(Z)
x, y = SpatialCoordinate(M)
# define f
f = Function(V)
source = as_vector((400*x*y*y*y*y*y*y+100*x*x*x*x*x*x*x-100*x*x*x*y*y*y*y+120*x*y,1200*x*x*y*y*y*y*y-100*x*x*x*x*y*y*y+100*y*y*y*y*y*y*y+60*x*x-60*y*y))
f.interpolate(source)
#define exact solution u
u_e = Function(V)
velocity = as_vector(((20*x*y*y*y),(5*x*x*x*x-5*y*y*y*y)))
u_e.interpolate(velocity)
# define exact pressure p
p_e = Function(W)
p_e.interpolate((60*x*x*y-20*y*y*y-5))
# define weak formulation
F = (inner(grad(u), grad(v))*dx +  inner(dot(grad(u), u), v) * dx  - dot(p, div(v))*dx - dot(div(u), q)*dx + rho*q*dx + lamda*p*dx - inner(f, v) * dx)
# define boundary condition
bcs = [DirichletBC(Z.sub(0), u_e, (1, 2, 3, 4))]
# solve the linear system with direct solver 
solve(F == 0, up, bcs=bcs)
u, p, rho  = up.subfunctions
 # define error for p and u
erH3=sqrt(assemble(inner(grad(u_e-u),grad(u_e-u))*dx)) 
print(erH3)

[dham]

When the error is fixed but not converging, this often means that it's converged to something other than what you are expecting. I would plot the solution, the analytical solution and the difference between them and see if that sheds light on the issue. See: https://www.firedrakeproject.org/visualisation.html#plotting-with-matplotlib

[A123-bansal123]

Thank you for helping me , the issue am facing is resolved.