Where to find 'project'-function in dolfinx?

Hi,

I am trying to project values on particular points onto a function space so that I can compare it to error values that live on that function space. Fenics provides a function called ‘project’ that accomplishes this, but now I am working with fenicsx and I can’t find the function anywhere (I looked at dolfinx and several of its submodules). Does it still exist and if so, where can I find it?

1 Like

project is just syntactic sugar. The projection f_h of a function f onto the space V is found by computing the problem: find f_h \in V such that

(f_h, v) = (f, v) \quad \forall v \in V.

This is very easy to set up. See for example @michalhabera 's implementation in dolfiny: dolfiny/projection.py at master · michalhabera/dolfiny · GitHub.

The problem with the project syntax is the uninformed user will assume it will handle many esoteric cases. For example, incomplete data in f, singular components in V, obscure functions f, reusing the LU factorisation for multiple projections onto the same space, and so on… Each of these user dependent cases should be handled by whatever the user’s numerical problem requires, and not a catch-all project function.

4 Likes

Thanks for the explanation @nate. I am working with @bouncer on doing a rather easy calculation for which previously we used the project function from FEniCS. Let me rephrase the question:

We have a function with a reference solution and a function with the calculated solution:

scalar_element = FiniteElement("P", mesh.ufl_cell(), 2)
V = FunctionSpace(mesh, scalar_element)
u_n = Function(V, name="Temperature")
u_ref = Function(V, name="reference")

over the course of the simulation, u_n gets resolved. At some point we want to calculate the relative error between the reference and the calculated solution in the l2 norm. Previously we did this by:

error_normalized = (u_ref - u_approx) / u_ref
error_pointwise = project(abs(error_normalized), v)
error_total = sqrt(assemble(inner(error_pointwise, error_pointwise) * dx))

for which we needed the project functionality. Is there a simpler way to calculate the total error in dolfinx? I am guessing what we are trying to do is not the most efficient way.

See for instance: Error control: Computing convergence rates — FEniCSx tutorial
or:
Test problem 1: Channel flow (Poiseuille flow) — FEniCSx tutorial
ref

import numpy as np
def u_exact(x):
    values = np.zeros((2, x.shape[1]), dtype=PETSc.ScalarType)
    values[0] = 4*x[1]*(1.0 - x[1])
    return values
u_ex = Function(V)
u_ex.interpolate(u_exact)

L2_error = form(dot(u_ - u_ex, u_ - u_ex)*dx)
# do stuff
# ....
for i in range(num_steps):
    # Solve problem
    #....

    # Compute error at current time-step
    error_L2 = np.sqrt(mesh.comm.allreduce(assemble_scalar(L2_error), op=MPI.SUM))

I would also like to note that you do not need the project functionality for computing the error “pointwise”, i.e.

could be handled with dolfinx.fem.Expression and interpolate, see for instance: Problem with grad and interpolate - #2 by dokken

I am also not sure why you do not just call

error_normalized = (u_ref - u_approx) / u_ref
error_total = sqrt(assemble(inner(error_normalized, error_normalized) * dx))

in DOLFIN, or the equivalent version of it in DOLFINx.

1 Like

Thank you for the extremely helpful answer! We now have a crisp and clear error computation code:

error_pointwise = form((u_approx - u_ref) ** 2 * dx)
error_total = np.sqrt(mesh.comm.allreduce(assemble_scalar(error_pointwise), MPI.SUM))

which works perfectly!

1 Like