Unable to compute the curl of a vector field

I see, yes.
I tried 2 approaches, still unable to get it working.
First approach is by implementing the function curl_2d(), as a copy-paste of the example you provide.

Jecurl_flux_calculator = Expression(curl_2d(J_vector), Jcurl.element.interpolation_points())

but I get the same kind of errors:

  File "/root/shared/debug_discourse/curl_MWE.py", line 188, in <module>
    Jecurl.interpolate(Jecurl_flux_calculator)
  File "/usr/local/dolfinx-real/lib/python3.10/dist-packages/dolfinx/fem/function.py", line 397, in interpolate
    _interpolate(u, cells)
  File "/usr/lib/python3.10/functools.py", line 889, in wrapper
    return dispatch(args[0].__class__)(*args, **kw)
  File "/usr/local/dolfinx-real/lib/python3.10/dist-packages/dolfinx/fem/function.py", line 388, in _
    self._cpp_object.interpolate(expr._cpp_object, cells)
RuntimeError: Function value size not equal to Expression value size

The 2nd approach is to use Nate’s code:

Jecurl_flux_calculator = Expression(as_vector((J_vector[1].dx(0)-J_vector[0].dx(1))), Jcurl.element.interpolation_points() )

yielding a different error:

Traceback (most recent call last):
  File "/root/shared/debug_discourse/curl_MWE.py", line 186, in <module>
    Jecurl_flux_calculator = Expression(as_vector((J_vector[1].dx(0)-J_vector[0].dx(1))), Jcurl.element.interpolation_points() )
  File "/usr/local/lib/python3.10/dist-packages/ufl/tensors.py", line 312, in as_vector
    raise ValueError("Expecting rank 1 tensor.")
ValueError: Expecting rank 1 tensor.

Consider the following code:


import ufl    
def curl_2d(a: fem.Function):
    return ufl.as_vector((0, 0, a[1].dx(0) - a[0].dx(1)))

Jcurl = functionspace(domain, ("DG", deg-1, (3,)))
Jecurl = Function(Jcurl)
Jecurl_flux_calculator = Expression(curl_2d(J_vector), Jcurl.element.interpolation_points())
Jecurl.interpolate(Jecurl_flux_calculator)

I have chosen the shape of Jcurl to be (3,), as you have defined curl_2D` as a 3D vector.
You could alternatively call:


Jcurl = functionspace(domain, ("DG", deg-1))
Jecurl = Function(Jcurl)
Jecurl_flux_calculator = Expression(curl_2d(J_vector), Jcurl.element.interpolation_points())
Jecurl.interpolate(Jecurl_flux_calculator)

where you define Jcurl to be a scalar and let curl_2D return a scalar value.

1 Like