Solving Maxwell equations - FacetNormal - Index out of bounds error

7 months ago by


I am trying to solve a scattering problem in 2D. Therefore I want to solve Maxwells equations with FEM. At the beginning the calculation shall take place on an empty square with absorbing boundary conditions. The resulting bilinear form has a term of the type  $\int_{\left\{\partial\Omega\right\}}\left(\vec{n}\times\vec{E}\right)\left(\vec{n}\times\vec{T}\right)ds${∂Ω}(n×E)(n×T)ds with E the trialfunction and T the testfunction. Assembling this bilinear form results in an "index out of bound error". Below is a minimal example and the error message (FEniCS 2017.1.0 installed via anaconda).


from dolfin import *

mesh = UnitSquareMesh(10, 10)

# function space
func_space = FunctionSpace(mesh, 'N1curl', 1)
# Trial and Test function
trial_func = TrialFunction(func_space)
test_func = TestFunction(func_space)
# Normal vector to facet
face_norm = FacetNormal(mesh)

# bilinear form
bilform = dot(cross(trial_func, face_norm), cross(test_func, face_norm))*ds

A = PETScMatrix()
assemble(bilform, tensor=A)
Index out of bounds.
Traceback (most recent call last):
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/dolfin/compilemodules/", line 142, in jit
    result = ffc.jit(ufl_object, parameters=p)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ffc/", line 204, in jit
    module = jit_build(ufl_object, module_name, parameters)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ffc/", line 120, in jit_build
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/dijitso/", line 160, in jit
    header, source, dependencies = generate(jitable, name, signature, params["generator"])
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ffc/", line 66, in jit_generate
    prefix=module_name, parameters=parameters, jit=True)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ffc/", line 141, in compile_form
    prefix, parameters, jit)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ffc/", line 183, in compile_ufl_objects
    analysis = analyze_ufl_objects(ufl_objects, kind, parameters)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ffc/", line 94, in analyze_ufl_objects
    for form in forms)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ffc/", line 94, in <genexpr>
    for form in forms)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ffc/", line 174, in _analyze_form
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 247, in compute_form_data
    form = apply_algebra_lowering(form)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 185, in apply_algebra_lowering
    return map_integrand_dags(LowerCompoundAlgebra(), expr)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 58, in map_integrand_dags
    form, only_integral_type)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 39, in map_integrands
    for itg in form.integrals()]
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 39, in <listcomp>
    for itg in form.integrals()]
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 46, in map_integrands
    return itg.reconstruct(function(itg.integrand()))
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 57, in <lambda>
    return map_integrands(lambda expr: map_expr_dag(function, expr, compress),
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/corealg/", line 37, in map_expr_dag
    result, = map_expr_dags(function, [expression], compress=compress)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/corealg/", line 86, in map_expr_dags
    r = handlers[v._ufl_typecode_](v, *[vcache[u] for u in v.ufl_operands])
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 76, in cross
    return as_vector((c(1, 2), c(2, 0), c(0, 1)))
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/algorithms/", line 75, in c
    return Product(a[i], b[j]) - Product(a[j], b[i])
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/", line 449, in _getitem
    all_indices, slice_indices, repeated_indices = create_slice_indices(component, shape, self.ufl_free_indices)
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/", line 166, in create_slice_indices
    error("Index out of bounds.")
  File "/opt/anaconda3/envs/fenics/lib/python3.6/site-packages/ufl/", line 172, in error
    raise self._exception_type(self._format_raw(*message))
ufl.log.UFLException: Index out of bounds.
Community: FEniCS Project

2 Answers

7 months ago by
The cross product is not defined for two dimensional vectors. If you try the same on a 3D mesh, it should work.
A workaround for 2D is to manually define your desired quantities.
7 months ago by
Thank you for the fast reply. How do I know which linear algebra parts work in 2D? I have another term that is given as  $\int_{\Omega}\left(\nabla\times E\right)\left(\nabla\times T\right)dA$Ω(×E)(×T)dA which compiles just fine.
That's odd. According to the UFL manual that shouldn't work either. Cross product and curl are only defined for vectors of length three.
written 7 months ago by klunkean  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »