# Error with HenrysBC and zero activation energy

**URL:** <https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113>\
**Category:** User Support\
**Tags:** festim2\
**Created:** [April 7, 2025, 11:06pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113 "2025-04-07T23:06:16Z")\
**Posts on this page:** 9\
**Page:** 1

<div class="post-metadata">

**Author:** ![VVKulagin](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/vvkulagin/32/51_2.png) [@VVKulagin](https://festim.discourse.group/u/VVKulagin)\
**Post date:** [April 7, 2025, 11:06pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/1 "2025-04-07T23:06:16Z")

</div>

Hi all!

I tried running a test problem using `HenrysBC` using a zero Henrys constant activation energy (E\_H=0). Here is MWE:

> **MWE**
>
> ```auto
> import festim as F
> import numpy as np
> 
> model = F.HydrogenTransportProblem()
> 
> vertices = np.linspace(0, 10, 50)
> 
> model.mesh = F.Mesh1D(vertices)
> 
> material = F.Material(D_0=1, E_D=0)
> 
> left_surface = F.SurfaceSubdomain1D(id=1, x=0)
> volume = F.VolumeSubdomain1D(id=2, borders=[0, 10], material=material)
> model.subdomains = [left_surface, volume]
> 
> H = F.Species("H")
> model.species = [H]
> 
> model.boundary_conditions = [
> F.HenrysBC(subdomain=left_surface, H_0=1, E_H=0, pressure=1, species=H)
> ]
> 
> model.temperature = 10
> 
> model.settings = F.Settings(atol=1e-10, rtol=1e-10, final_time=30)
> 
> model.settings.stepsize = F.Stepsize(0.05)
> 
> model.initialise()
> model.run()
> 
> ```

which raises the following error

> **Error message**
>
> ```auto
> Could not extract MPI communicator for Expression. Maybe you need to pass a communicator?
> 
> ---------------------------------------------------------------------------
> AttributeError Traceback (most recent call last)
> Cell In[12], line 29
> 25 model.settings = F.Settings(atol=1e-10, rtol=1e-10, final_time=30)
> 27 model.settings.stepsize = F.Stepsize(0.05)
> ---> 29 model.initialise()
> 30 model.run()
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/festim/hydrogen_transport_problem.py:288, in HydrogenTransportProblem.initialise(self)
> 285 self.create_implicit_species_value_fenics()
> 287 self.define_temperature()
> --> 288 self.define_boundary_conditions()
> 289 self.convert_source_input_values_to_fenics_objects()
> 290 self.convert_advection_term_to_fenics_objects()
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/festim/hydrogen_transport_problem.py:577, in HydrogenTransportProblem.define_boundary_conditions(self)
> 570 if isinstance(bc, boundary_conditions.ParticleFluxBC):
> 571 bc.create_value_fenics(
> 572 mesh=self.mesh.mesh,
> 573 temperature=self.temperature_fenics,
> 574 t=self.t,
> 575 )
> --> 577 super().define_boundary_conditions()
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/festim/problem.py:116, in ProblemBase.define_boundary_conditions(self)
> 114 for bc in self.boundary_conditions:
> 115 if isinstance(bc, F.DirichletBCBase):
> --> 116 form = self.create_dirichletbc_form(bc)
> 117 self.bc_forms.append(form)
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/festim/hydrogen_transport_problem.py:595, in HydrogenTransportProblem.create_dirichletbc_form(self, bc)
> 592 else:
> 593 function_space_value = bc.species.collapsed_function_space
> --> 595 bc.create_value(
> 596 temperature=self.temperature_fenics,
> 597 function_space=function_space_value,
> 598 t=self.t,
> 599 )
> 601 # get dofs
> 602 if self.multispecies and isinstance(bc.value_fenics, (fem.Function)):
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/festim/boundary_conditions/dirichlet_bc.py:246, in FixedConcentrationBC.create_value(self, function_space, temperature, t, K_S)
> 242 kwargs["T"] = temperature
> 244 # store the expression of the boundary condition
> 245 # to update the value_fenics later
> --> 246 self.bc_expr = fem.Expression(
> 247 self.value(**kwargs),
> 248 helpers.get_interpolation_points(function_space.element),
> 249 )
> 250 self.value_fenics.interpolate(self.bc_expr)
> 252 # if K_S is provided, divide the value by K_S (change of variable method)
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/dolfinx/fem/function.py:131, in Expression. __init__ (self, e, X, comm, form_compiler_options, jit_options, dtype)
> 129 if comm is None:
> 130 try:
> --> 131 mesh = ufl.domain.extract_unique_domain(e).ufl_cargo()
> 132 comm = mesh.comm
> 133 except AttributeError:
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/ufl/domain.py:242, in extract_unique_domain(expr)
> 240 def extract_unique_domain(expr):
> 241 """Return the single unique domain expression is defined on or throw an error."""
> --> 242 domains = extract_domains(expr)
> 243 if len(domains) == 1:
> 244 return domains[0]
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/ufl/domain.py:232, in extract_domains(expr)
> 230 """Return all domains expression is defined on."""
> 231 domainlist = []
> --> 232 for t in traverse_unique_terminals(expr):
> 233 domainlist.extend(t.ufl_domains())
> 234 return sorted(
> 235 join_domains(domainlist),
> 236 key=lambda D: (D.topological_dimension(), D.ufl_cell(), D.ufl_id()),
> 237 )
> 
> File ~/anaconda3/envs/vv-festim-report-env/lib/python3.13/site-packages/ufl/corealg/traversal.py:141, in traverse_unique_terminals(expr, visited)
> 139 """Traverse unique terminals."""
> 140 for op in unique_pre_traversal(expr, visited=visited):
> --> 141 if op._ufl_is_terminal_:
> 142 yield op
> 
> AttributeError: 'float' object has no attribute '_ufl_is_terminal_'
> 
> ```

---

<div class="post-metadata">

**Author:** ![remidm](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/remidm/32/4_2.png) [@remidm](https://festim.discourse.group/u/remidm)\
**Post date:** [April 7, 2025, 11:19pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/2 "2025-04-07T23:19:19Z")

</div>

What version of FESTIM is this using?

---

<div class="post-metadata">

**Author:** ![VVKulagin](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/vvkulagin/32/51_2.png) [@VVKulagin](https://festim.discourse.group/u/VVKulagin)\
**Post date:** [April 7, 2025, 11:34pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/3 "2025-04-07T23:34:34Z")

</div>

Right, sorry

It is the 2.0a2 version.

---

<div class="post-metadata">

**Author:** ![remidm](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/remidm/32/4_2.png) [@remidm](https://festim.discourse.group/u/remidm)\
**Post date:** [April 7, 2025, 11:49pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/4 "2025-04-07T23:49:39Z")

</div>

Do you have the error with SievertsBC?

---

<div class="post-metadata">

**Author:** ![VVKulagin](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/vvkulagin/32/51_2.png) [@VVKulagin](https://festim.discourse.group/u/VVKulagin)\
**Post date:** [April 7, 2025, 11:58pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/5 "2025-04-07T23:58:09Z")

</div>

Yes, the same problem.

---

<div class="post-metadata">

**Author:** ![remidm](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/remidm/32/4_2.png) [@remidm](https://festim.discourse.group/u/remidm)\
**Post date:** [April 8, 2025, 1:23pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/6 "2025-04-08T13:23:01Z")

</div>

I made a `FixedConcentrationBC` version giving the same error:

```python
import festim as F
import numpy as np
import ufl

model = F.HydrogenTransportProblem()

vertices = np.linspace(0, 10, 50)

model.mesh = F.Mesh1D(vertices)

material = F.Material(D_0=1, E_D=0)

left_surface = F.SurfaceSubdomain1D(id=1, x=0)
volume = F.VolumeSubdomain1D(id=2, borders=[0, 10], material=material)
model.subdomains = [left_surface, volume]

H = F.Species("H")
model.species = [H]

def arrhenius(T, x0, Ex):
    return x0 * ufl.exp(-Ex / F.k_B / T)

model.boundary_conditions = [
    F.FixedConcentrationBC(
        subdomain=left_surface, value=lambda T: arrhenius(T, 1, 0), species=H
    ),
]

model.temperature = 1000

model.settings = F.Settings(atol=1e-10, rtol=1e-10, final_time=30)

model.settings.stepsize = F.Stepsize(0.05)

model.initialise()
model.run()

```

---

<div class="post-metadata">

**Author:** ![remidm](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/remidm/32/4_2.png) [@remidm](https://festim.discourse.group/u/remidm)\
**Post date:** [April 8, 2025, 1:27pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/7 "2025-04-08T13:27:35Z")

</div>

And here’s a pure `dolfinx` version:

```python
import dolfinx
from mpi4py import MPI
import ufl

mesh = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 10, 10)
V = dolfinx.fem.functionspace(mesh, ("CG", 1))

T = dolfinx.fem.Function(V)

def arrhenius(T, x0, Ex):
    return x0 * ufl.exp(-Ex / T)

bc_expr = dolfinx.fem.Expression(
    arrhenius(T=T, x0=1, Ex=0), V.element.interpolation_points()
)

```

---

<div class="post-metadata">

**Author:** ![remidm](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/remidm/32/4_2.png) [@remidm](https://festim.discourse.group/u/remidm)\
**Post date:** [April 8, 2025, 2:03pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/8 "2025-04-08T14:03:32Z")

</div>

```python
bc_expr = dolfinx.fem.Expression(
    dolfinx.fem.Constant(mesh, 1.0),
    V.element.interpolation_points(),
)

```

This works though, so I guess you could use a `FixedConcentrationBC` instead and don’t bother with the arrhenius expression of K\_H

---

<div class="post-metadata">

**Author:** ![remidm](https://yyz2.discourse-cdn.com/free1/user_avatar/festim.discourse.group/remidm/32/4_2.png) [@remidm](https://festim.discourse.group/u/remidm)\
**Post date:** [April 8, 2025, 3:37pm UTC](https://festim.discourse.group/t/error-with-henrysbc-and-zero-activation-energy/113/9 "2025-04-08T15:37:53Z")

</div>

@VVKulagin the workaround for now is:

```python
solubility = 1.0
pressure = 1.0

model.boundary_conditions = [
    F.FixedConcentrationBC(
        subdomain=left_surface, value=solubility * pressure, species=H
    ),
]

```