How to use the shape parameter of dolfinx.fem.FunctionSpace?

variational formulation
1

Hello,

It seems that dolfinx.fem.FunctionSpace{.python} can take an optional tuple{.python} as a shape{.verbatim} parameter. How is this chosen? In the unit test, it is:

gdim = mesh.geometry.dim
V = FunctionSpace(mesh, ("Lagrange", 1, (gdim,)))

which makes me think that I could always use this construct, but doing this (from the tutorial) fails:

from mpi4py import MPI
from dolfinx import mesh
from dolfinx.fem import FunctionSpace
from dolfinx import fem


domain = mesh.create_unit_square(
      MPI.COMM_WORLD, 8, 8, mesh.CellType.quadrilateral)

fe_type = "Lagrange"
fe_degree = 1
gdim = domain.geometry.dim
# Original
# V = FunctionSpace(domain, (fe_type, fe_degree))
V = FunctionSpace(domain, (fe_type, fe_degree, (gdim,)))

uD = fem.Function(V)
uD.interpolate(lambda x_in: 1 + x_in[0]**2 + 2 * x_in[1]**2)

Traceback (most recent call last):
  File "/usr/lib/python3.11/site-packages/dolfinx/fem/function.py", line 397, in interpolate
    _interpolate(u, cells)
  File "/usr/lib/python3.11/functools.py", line 909, in wrapper
    return dispatch(args[0].__class__)(*args, **kw)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/usr/lib/python3.11/site-packages/dolfinx/fem/function.py", line 373, in _interpolate
    self._cpp_object.interpolate(u, cells, nmm_interpolation_data)
TypeError: interpolate(): incompatible function arguments. The following argument types are supported:
    1. (self: dolfinx.cpp.fem.Function_float64, f: numpy.ndarray[numpy.float64], cells: numpy.ndarray[numpy.int32]) -> None
    2. (self: dolfinx.cpp.fem.Function_float64, u: dolfinx.cpp.fem.Function_float64, cells: numpy.ndarray[numpy.int32], nmm_interpolation_data: Tuple[List[int], List[int], List[float], List[int]]) -> None
    3. (self: dolfinx.cpp.fem.Function_float64, expr: dolfinx::fem::Expression<double, double>, cells: numpy.ndarray[numpy.int32]) -> None

Invoked with: <dolfinx.cpp.fem.Function_float64 object at 0x7f13febf3570>, <function <lambda> at 0x7f13fea002c0>, array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12, 13, 14, 15, 16,
       17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33,
       34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50,
       51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63], dtype=int32), ((), (), (), ())

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<string>", line 19, in __PYTHON_EL_eval
  File "<string>", line 18, in <module>
  File "/usr/lib/python3.11/site-packages/dolfinx/fem/function.py", line 402, in interpolate
    self._cpp_object.interpolate(np.asarray(u(x), dtype=self.dtype), cells)
RuntimeError: Interpolation data has the wrong shape/size.

Thank you kindly for your response.

My system

* FEniCSx software
dolfinx: 0.7.0.dev0_r27554.cfeffe0-1
basix: 0.7.0.dev0_r945.1117a8d-1
ufl: 2023.2.0.dev0_r3562.77ae57c-1
ffcx: 0.7.0.dev0_r7077.1d27238-1

* Dependencies
python: 3.11.5-2
python-numpy: 1.26.0-1
petsc: 3.19.5.1171.g37df9106526-1
hdf5-openmpi: 1.14.2-1
boost: 1.83.0-2
adios2: 2.8.3-5
scotch: 7.0.4-1
pybind11: 2.11.1-1
python-build: 1.0.1-1
python-cffi: 1.15.1-4
python-cppimport: 22.08.02.r6.g0849d17-1
python-installer: 0.7.0-3
python-mpi4py: 3.1.4-3
python-pytest: 7.4.2-1
python-scikit-build: 0.17.6-1
python-setuptools: 1:68.0.0-1
python-wheel: 0.40.0-3
xtensor: 0.24.0-2
xtensor-blas: 0.20.0-1

* Operating system
Arch Linux: 6.5.4-arch2-1
2

The shape parameter is our way of replacing VectorFunctionSpace and TensorFunctionSpace. If you Give the argument shape=(gdim,) you mage a function space with gdim components, where each component is from a first order Lagrange space.

For scalar valued function spaces, the shape is (), (which is the default).

1 Like
3

Thanks, I’m still trying to get my head around this (and it’s not necessarily related to DOLFINx):

As I understand (I can be wrong),

  1. both the (shown parts of the) tutorial and test file have 1st-degree polynomial Lagrange interpolation functions;
  2. in the tutorial, the domain is 2D;
  3. in the test file, the domain is 3D;
  4. a 2D domain requires a function space with 2 components (gdim = 2);
  5. a 3D domain requires a function space with 3 components (gdim = 3);
  6. mesh.create_unit_square creates a 2D domain;
  7. mesh.create_unit_cube (not shown here, but found in the test file) creates a 3D domain;
  8. the shape parameter of FunctionSpace is a tuple of the form (gdim,) (or () for scalars);
  9. when shape is (2,), FunctionSpace mages a space with 2 components (needed for 2D), and
  10. when shape is (3,), FunctionSpace mages a space with 3 components (needed for 3D).

Thus, what I am trying to understand is:

  1. Which of these is wrong?
  2. why does shape = () work for a 2D domain (as in the tutorial, if it is meant for scalars)? and
  3. why does shape = (2,) not work for a 2D domain (if shape = (3,) works for 3D)?

Thank you very much. I appreciate your time and patience (not to mention the amazing work to put FEniCSx together).

4

Imagine you are solving a diffusion equation. Then the concentration parameter (c is represented by as scalar value. Even if your computational grid is 2D or 3D, the Function that you want to represent is a scalar value, thus, dolfinx.fem.functionspace(mesh, ("Lagrange", 1)) create a fucntion space for a scalar valued function. I.e. if you create a function in this space c, c:\Omega\subset\mathbb{R}^N\mapsto \mathbb{R} where N is the geometrical dimension of the mesh (1, 2 or 3).

Let us then assume that we are either in 2D or 3D and want to represent a fluid velocity u. u has as many components as the geometrical dimension of the mesh, thus dolfinx.fem.functionspace(mesh, ("Lagrange", 1, (mesh.geometry.dim, )) creates a space for functions u, u:\Omega\subset\mathbb{R}^N\mapsto\mathbb{R}^N.

There are also cases when you want to represent functions in spaces that are not equal to the geometrical dimension. For instance assume that you have the concentration of 10 scalar species co-existing in space at any point. Then dolfinx.fem.functionspace(mesh, ("Lagrange", 1, (10, )) creates a
c:\Omega\subset\mathbb{R}^N\mapsto\mathbb{R}^{10}.

Of course you can also make tensor valued function spaces as well, using (N, M).

4 Likes
5

Ok. Thank you very much! That is totally clear (amazingly thorough and quick).

6

Thank you for this interesting topic and the amazing explanation! However, what is not clear to me is the potential difference between dolfinx.fem.functionspace(mesh, ("Lagrange", 1, (2, )) and MixedElement. Imagine I’m solving a diffusion-reaction problem using 1D elements, where I have the concentration of 2 scalar species co-existing in space at any point. The two scalar species interact with each other. If I use dolfinx.fem.functionspace(mesh, ("Lagrange", 1, (2, )), do I need to use something like a staggered scheme to ensure the equilibrium between the 2 species at each time step? In other words, is dolfinx.fem.functionspace(mesh, ("Lagrange", 1, (2, )) more like defining 2 FunctionSpaces or more like using a MixedElement? Thank you!

7

Mixed element is a generalization of the shape parameter, as you can have different elements for each component of you space.

8

Thank you so much for the prompt response! I see. In my case, I only need linear element for both components/species of the space. I had some difficulties in using ‘split’ and ‘assign’ of MixedElement correctly. So glad to know dolfinx.fem.functionspace(mesh, ("Lagrange", 1, (2, )) which appears to be suitable for my case. Many thanks!