Update time dependent bcs

festim/boundary_conditions/dirichlet_bc.py

@@ -24,6 +24,8 @@ class DirichletBC:
             fenics format
         bc_expr (fem.Expression): the expression of the boundary condition that is used to
             update the value_fenics
+        time_dependent (bool): True if the value of the bc is time dependent
+        temperature_dependent (bool): True if the value of the bc is temperature dependent
 
     Usage:
         >>> from festim import DirichletBC
@@ -57,6 +59,30 @@ def value_fenics(self, value):
             )
         self._value_fenics = value
 
+    @property
+    def time_dependent(self):
+        if self.value is None:
+            return False
+        if isinstance(self.value, fem.Constant):
+            return False
+        if callable(self.value):
+            arguments = self.value.__code__.co_varnames
+            return "t" in arguments
+        else:
+            return False
+
+    @property
+    def temperature_dependent(self):
+        if self.value is None:
+            return False
+        if isinstance(self.value, fem.Constant):
+            return False
+        if callable(self.value):
+            arguments = self.value.__code__.co_varnames
+            return "T" in arguments
+        else:
+            return False
+
     def define_surface_subdomain_dofs(self, facet_meshtags, mesh, function_space):
         """Defines the facets and the degrees of freedom of the boundary
         condition

festim/hydrogen_transport_problem.py

@@ -570,7 +570,10 @@ def update_time_dependent_values(self):
                 self.temperature_fenics.interpolate(self.temperature_expr)
 
         for bc in self.boundary_conditions:
-            bc.update(t=t)
+            if bc.time_dependent:
+                bc.update(t=t)
+            elif self.temperature_time_dependent and bc.temperature_dependent:
+                bc.update(t=t)
 
     def post_processing(self):
         """Post processes the model"""

test/test_dirichlet_bc.py

@@ -404,3 +404,48 @@ def test_integration_with_a_multispecies_HTransportProblem(value_A, value_B):
     if isinstance(expected_value, Conditional):
         expected_value = float(expected_value)
     assert np.isclose(computed_value, expected_value)
+
+
+@pytest.mark.parametrize(
+    "input, expected_value",
+    [
+        (1.0, False),
+        (None, False),
+        (fem.Constant(mesh, 1.0), False),
+        (lambda t: t, True),
+        (lambda t: 1.0 + t, True),
+        (lambda x: 1.0 + x[0], False),
+        (lambda x, t: 1.0 + x[0] + t, True),
+        (lambda x, t, T: 1.0 + x[0] + t + T, True),
+        (lambda x, t: ufl.conditional(ufl.lt(t, 1.0), 100.0 + x[0], 0.0), True),
+    ],
+)
+def test_bc_time_dependent_attribute(input, expected_value):
+    """Test that the time_dependent attribute is correctly set"""
+    subdomain = F.SurfaceSubdomain1D(1, x=0)
+    species = F.Species("test")
+    my_bc = F.DirichletBC(subdomain, input, species)
+
+    assert my_bc.time_dependent is expected_value
+
+
+@pytest.mark.parametrize(
+    "input, expected_value",
+    [
+        (1.0, False),
+        (None, False),
+        (fem.Constant(mesh, 1.0), False),
+        (lambda T: T, True),
+        (lambda t: 1.0 + t, False),
+        (lambda x, T: 1.0 + x[0] + T, True),
+        (lambda x, t, T: 1.0 + x[0] + t + T, True),
+        (lambda x, t: ufl.conditional(ufl.lt(t, 1.0), 100.0 + x[0], 0.0), False),
+    ],
+)
+def test_bc_temperature_dependent_attribute(input, expected_value):
+    """Test that the temperature_dependent attribute is correctly set"""
+    subdomain = F.SurfaceSubdomain1D(1, x=0)
+    species = F.Species("test")
+    my_bc = F.DirichletBC(subdomain, input, species)
+
+    assert my_bc.temperature_dependent is expected_value

test/test_h_transport_problem.py

@@ -202,10 +202,8 @@ def test_update_time_dependent_values_temperature(T_function, expected_values):
         # TEST
         if isinstance(my_model.temperature_fenics, fem.Constant):
             computed_value = float(my_model.temperature_fenics)
-            print(computed_value)
         else:
             computed_value = my_model.temperature_fenics.vector.array[-1]
-            print(computed_value)
         assert np.isclose(computed_value, expected_values[i])
 
 
@@ -423,3 +421,58 @@ def test_post_processing_update_D_global():
 
     # TEST
     assert value_t_1 != value_t_2
+
+
+@pytest.mark.parametrize(
+    "temperature_value, bc_value, expected_values",
+    [
+        (5, 1.0, [1.0, 1.0, 1.0]),
+        (lambda t: t + 1, lambda T: T, [2.0, 3.0, 4.0]),
+        (lambda x: 1 + x[0], lambda T: 2 * T, [10.0, 10.0, 10.0]),
+        (lambda x, t: t + x[0], lambda T: 0.5 * T, [2.5, 3.0, 3.5]),
+        (
+            lambda x, t: ufl.conditional(ufl.lt(t, 2), 3 + x[0], 0.0),
+            lambda T: T,
+            [7.0, 0.0, 0.0],
+        ),
+    ],
+)
+def test_update_time_dependent_bcs_with_time_dependent_temperature(
+    temperature_value, bc_value, expected_values
+):
+    """Test that temperature dependent bcs are updated at each time step when the
+    temperature is time dependent, and match an expected value"""
+
+    # BUILD
+    H = F.Species("H")
+    volume_subdomain = F.VolumeSubdomain1D(id=1, borders=[0, 4], material=dummy_mat)
+    surface_subdomain = F.SurfaceSubdomain1D(id=1, x=4)
+    my_model = F.HydrogenTransportProblem(
+        mesh=test_mesh, species=[H], subdomains=[volume_subdomain, surface_subdomain]
+    )
+    my_model.t = fem.Constant(my_model.mesh.mesh, 0.0)
+    dt = fem.Constant(test_mesh.mesh, 1.0)
+
+    my_model.temperature = temperature_value
+    my_bc = F.DirichletBC(species=H, value=bc_value, subdomain=surface_subdomain)
+    my_model.boundary_conditions = [my_bc]
+
+    my_model.define_temperature()
+    my_model.define_function_spaces()
+    my_model.assign_functions_to_species()
+    my_model.define_markers_and_measures()
+    my_model.create_dirichletbc_form(my_bc)
+
+    for i in range(3):
+        # RUN
+        my_model.t.value += dt.value
+        my_model.update_time_dependent_values()
+
+        # TEST
+        if isinstance(my_model.boundary_conditions[0].value_fenics, fem.Constant):
+            computed_value = float(my_model.boundary_conditions[0].value_fenics)
+        else:
+            computed_value = my_model.boundary_conditions[0].value_fenics.vector.array[
+                -1
+            ]
+        assert np.isclose(computed_value, expected_values[i])