How to use same vector-valued mathematical expressions with both ufl and numpy?

[August]

Hi all,

I’m using this kind of technique to avoid having to write an expression both in ufl (eg for generating the exact right hand side f) and numpy (eg for pointwise evaluations):

import dolfinx
import ufl
from mpi4py import MPI
import numpy as np

mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 5, 5)

def u_exact(m):
    return lambda x: m.sin(m.pi * x[0]) * m.sin(m.pi * x[1])

x = ufl.SpatialCoordinate(mesh)
f = ufl.div(ufl.grad(u_exact(ufl)(x)))
pts = np.array([[0.5, 0.1], [0.1, 0.9]])
u0 = u_exact(np)(pts)

What I cannot get to work is how to do this if u_exact is a vector-valued function. I’m trying to be clever with mixing np.stack and ufl.as_vector but without success. Is there an example somewhere I can look at?

Thanks!

[francesco-ballarin]

For legacy dolfin there used to be GitHub - MiroK/ulfy: Traslate UFL to SymPy to Expression (ufl ↔ sympy, but then from sympy one could use sympy.lambdify to numpy). I am not aware of similar projects for dolfinx.

[bleyerj]

The fem.Expression mechanism would be a solution if your points are defined via the same coordinates on the reference element. I am not sure if something is similar with points defined by their global coordinates

[bleyerj]

Maybe this ?

Evaluating Finite Element Function

The points have the wrong shape (single point, the function expects multiple points): from mpi4py import MPI import numpy as np from dolfinx import mesh from dolfinx import fem from dolfinx.fem import FunctionSpace from dolfinx import default_scalar_type import ufl from dolfinx.fem.petsc import LinearProblem from dolfinx import geometry domain = mesh.create_unit_square(MPI.COMM_WORLD, 8, 8, mesh.CellType.quadrilateral) V = FunctionSpace(domain, ("Lagrange", 1)) uD = fem.Function(V)…

[dokken]

For point evaluation, you need to specify what cell the point is in, as specified in: Implementation — FEniCSx tutorial

One can of course go through dolfinx.fem.Expression for more advanced code, such as shown in: Evaluation of UFL expression at mesh points - #2 by dokken and subsequent posts.

Other similar resources:

Point Sources (Redux)

Following my post above, with a minor adaptations, a point source can be implemented as follows: # Create a point source for Poisson problem # Author: Jørgen S. Dokken # SPDX-License-Identifier: MIT from mpi4py import MPI from petsc4py import PETSc import dolfinx import dolfinx.fem.petsc import numpy as np import ufl def compute_cell_contributions(V, points): # Determine what process owns a point and what cells it lies within mesh = V.mesh _, _, owning_points, cells = dolfinx.cp…

Evaluate_basis_all in dolfinx

The DOLFINx way to do this would be to: Figure out which cell the point x is in (using dolfinx.compute_colliding_cells, as done in: Implementation — FEniCSx tutorial) Pull these coordinates back to the reference element, using mesh.geometry.cmap.pull_back (dolfinx/fem.cpp at 48ece0087588f5f274a36aec223abb129fa59e8e · FEniCS/dolfinx · GitHub) Use basix.FiniteElement(element_of choice).tabulate(1, X_ref) to tabulate at the reference points. An example of the last step can be found in: asimov-cu…