Error interpolating on sub of vector function

Hello,

I believe this is a simple question but have not found a solution. I have a fem.Function with a vector function space and I am trying to interpolate a constant value onto one of the subspaces. In this attempt, I get a long error message (pasted below) and the kernel keeps running, producing the error messages, until I manually stop it. The MWE is

import numpy
from dolfinx import mesh, fem, nls, plot, cpp
import ufl
from mpi4py import MPI
from petsc4py import PETSc
import pyvista

import gmsh
from dolfinx.io import gmshio
gmsh.initialize()

dx_min = 0.1
dx_max = 0.1
cyl = gmsh.model.occ.addDisk(0, 0, 0, 1, 1)
gmsh.model.occ.synchronize()
gdim = 2
gmsh.model.addPhysicalGroup(gdim, [cyl], 1)
gmsh.option.setNumber("Mesh.CharacteristicLengthMin",dx_min)
gmsh.option.setNumber("Mesh.CharacteristicLengthMax",dx_max)
gmsh.model.mesh.generate(gdim)

# convert to fenics format
gmsh_model_rank = 0
mesh_comm = MPI.COMM_WORLD
domain, cell_markers, facet_markers = gmshio.model_to_mesh(gmsh.model, mesh_comm, gmsh_model_rank, gdim=gdim)

# vector function space
W = fem.VectorFunctionSpace(domain, ("CG",1))

# create function and set values of sub
u = fem.Function(W)
u.sub(0).interpolate(fem.Constant(domain, PETSc.ScalarType(1)))

Error messages:

Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. WARNING:UFL:Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation. Couldn't map 'c_0' to a float, returning ufl object without evaluation.

Thanks!

If you want to interpolate a constant you could either:

# Map from a sub function (cheap, but only works for constants)
Wx, map_to_W = W.sub(0).collapse()
ux = dolfinx.fem.Function(V)
ux.x.set(1)
u.x.array[map_to_W] = ux.x.array
# Map with an expression (not cheap for constants,but more general)
c = dolfinx.fem.Constant(domain, PETSc.ScalarType(1))
sub_expr = dolfinx.fem.Expression(c, W.element.interpolation_points())
u.sub(0).interpolate(sub_expr)

Thank you for the quick reply! I used the map with expression approach and get more errors. The full code is:

import numpy
from dolfinx import mesh, fem, nls, plot, cpp
import ufl
from mpi4py import MPI
from petsc4py import PETSc
import pyvista

import gmsh
from dolfinx.io import gmshio
gmsh.initialize()

dx_min = 0.1
dx_max = 0.1
cyl = gmsh.model.occ.addDisk(0, 0, 0, 1, 1)
gmsh.model.occ.synchronize()
gdim = 2
gmsh.model.addPhysicalGroup(gdim, [cyl], 1)
gmsh.option.setNumber("Mesh.CharacteristicLengthMin",dx_min)
gmsh.option.setNumber("Mesh.CharacteristicLengthMax",dx_max)
gmsh.model.mesh.generate(gdim)

# convert to fenics format
gmsh_model_rank = 0
mesh_comm = MPI.COMM_WORLD
domain, cell_markers, facet_markers = gmshio.model_to_mesh(gmsh.model, mesh_comm, gmsh_model_rank, gdim=gdim)

# vector function space
W = fem.VectorFunctionSpace(domain, ("CG",1))

# create function and set values of sub
u = fem.Function(W)
c = fem.Constant(domain, PETSc.ScalarType(1))
sub_expr = fem.Expression(c, W.element.interpolation_points())
u.sub(0).interpolate(sub_expr)

and the Error messages are:

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
Cell In[1], line 33
     31 u = fem.Function(W)
     32 c = fem.Constant(domain, PETSc.ScalarType(1))
---> 33 sub_expr = fem.Expression(c, W.element.interpolation_points())
     34 u.sub(0).interpolate(sub_expr)

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/dolfinx/fem/function.py:117, in Expression.__init__(self, ufl_expression, X, form_compiler_params, jit_params, dtype)
    114 else:
    115     raise RuntimeError(f"Unsupported scalar type {dtype} for Expression.")
--> 117 self._ufcx_expression, module, self._code = jit.ffcx_jit(mesh.comm, (ufl_expression, _X),
    118                                                          form_compiler_params=form_compiler_params,
    119                                                          jit_params=jit_params)
    120 self._ufl_expression = ufl_expression
    122 # Prepare coefficients data. For every coefficient in form take
    123 # its C++ object.

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/dolfinx/jit.py:56, in mpi_jit_decorator.<locals>.mpi_jit(comm, *args, **kwargs)
     51 @functools.wraps(local_jit)
     52 def mpi_jit(comm, *args, **kwargs):
     53 
     54     # Just call JIT compiler when running in serial
     55     if comm.size == 1:
---> 56         return local_jit(*args, **kwargs)
     58     # Default status (0 == ok, 1 == fail)
     59     status = 0

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/dolfinx/jit.py:211, in ffcx_jit(ufl_object, form_compiler_params, jit_params)
    208     r = ffcx.codegeneration.jit.compile_coordinate_maps(
    209         [ufl_object], parameters=p_ffcx, **p_jit)
    210 elif isinstance(ufl_object, tuple) and isinstance(ufl_object[0], ufl.core.expr.Expr):
--> 211     r = ffcx.codegeneration.jit.compile_expressions([ufl_object], parameters=p_ffcx, **p_jit)
    212 else:
    213     raise TypeError(type(ufl_object))

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/ffcx/codegeneration/jit.py:193, in compile_expressions(expressions, parameters, cache_dir, timeout, cffi_extra_compile_args, cffi_verbose, cffi_debug, cffi_libraries)
    182 """Compile a list of UFL expressions into UFC Python objects.
    183 
    184 Parameters
   (...)
    188 
    189 """
    190 p = ffcx.parameters.get_parameters(parameters)
    192 module_name = 'libffcx_expressions_' + \
--> 193     ffcx.naming.compute_signature(expressions, _compute_parameter_signature(p)
    194                                   + str(cffi_extra_compile_args) + str(cffi_debug))
    195 expr_names = [ffcx.naming.expression_name(expression, module_name) for expression in expressions]
    197 if cache_dir is not None:

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/ffcx/naming.py:62, in compute_signature(ufl_objects, tag)
     59 rn.update(dict((d, i) for i, d in enumerate(domains)))
     61 # Hash on UFL signature and points
---> 62 signature = ufl.algorithms.signature.compute_expression_signature(expr, rn)
     63 object_signature += signature
     64 object_signature += repr(points)

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/ufl/algorithms/signature.py:118, in compute_expression_signature(expr, renumbering)
    114 def compute_expression_signature(expr, renumbering):  # FIXME: Fix callers
    115     # FIXME: Rewrite in terms of compute_form_signature?
    116 
    117     # Build hashdata for all terminals first
--> 118     terminal_hashdata = compute_terminal_hashdata([expr], renumbering)
    120     # Build hashdata for full expression
    121     expression_hashdata = compute_expression_hashdata(expr, terminal_hashdata)

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/ufl/algorithms/signature.py:64, in compute_terminal_hashdata(expressions, renumbering)
     61     data = expr._ufl_signature_data_(renumbering)
     63 elif isinstance(expr, Constant):
---> 64     data = expr._ufl_signature_data_(renumbering)
     66 elif isinstance(expr, Argument):
     67     data = expr._ufl_signature_data_(renumbering)

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/ufl/constant.py:72, in Constant._ufl_signature_data_(self, renumbering)
     69 def _ufl_signature_data_(self, renumbering):
     70     "Signature data for constant depends on renumbering"
     71     return "Constant({}, {}, {})".format(
---> 72         self._ufl_domain._ufl_signature_data_(renumbering), repr(self._ufl_shape),
     73         repr(renumbering[self]))

File /usr/local/anaconda3/envs/fenicsx-env/lib/python3.10/site-packages/ufl/domain.py:121, in Mesh._ufl_signature_data_(self, renumbering)
    120 def _ufl_signature_data_(self, renumbering):
--> 121     return ("Mesh", renumbering[self], self._ufl_coordinate_element)

KeyError: Mesh(VectorElement(FiniteElement('Lagrange', triangle, 1, variant='equispaced'), dim=2, variant='equispaced'), 0)

I cannot reproduce this with dolfinx v0.6.0. This has been fixed in `dolfinx.fem.Expression` containing only a `dolfinx.fem.Constant` fails · Issue #2342 · FEniCS/dolfinx · GitHub

Ah, yes. I installed FEniCS with conda following the instructions https://github.com/FEniCS/dolfinx#conda and it appears I have dolfinx version 0.5.2. I’ll get that updated.

Thank you so much!

The latest version in conda-forge is v0.6.0: Fenics Dolfinx :: Anaconda.org

Thanks. I see now it was update 19 days ago. I installed 21 days ago =). The code works no, thank you so much!