Skip to content

Commit

Permalink
Merge branch 'fenicsx' into surface-reactions
Browse files Browse the repository at this point in the history
  • Loading branch information
RemDelaporteMathurin committed May 1, 2024
2 parents b1f2332 + 0b02ffa commit 936738d
Show file tree
Hide file tree
Showing 6 changed files with 46 additions and 12 deletions.
Empty file.
7 changes: 6 additions & 1 deletion festim/boundary_conditions/flux_bc.py
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,7 @@ def value_fenics(self, value):
value, (fem.Function, fem.Constant, np.ndarray, ufl.core.expr.Expr)
):
raise TypeError(
f"Value must be a dolfinx.fem.Function, dolfinx.fem.Constant, or a np.ndarray not {type(value)}"
f"Value must be a dolfinx.fem.Function, dolfinx.fem.Constant, np.ndarray or ufl.core.expr.Expr not {type(value)}"
)
self._value_fenics = value

Expand Down Expand Up @@ -141,6 +141,8 @@ class ParticleFluxBC(FluxBCBase):
is applied
value (float, fem.Constant, callable): the value of the particle flux
species (festim.Species): the species to which the flux is applied
species_dependent_value (dict): a dictionary containing the species that the value. Example: {"name": species}
where "name" is the variable name in the callable value and species is a festim.Species object.
Attributes:
subdomain (festim.SurfaceSubdomain): the surface subdomain where the particle flux
Expand All @@ -151,6 +153,8 @@ class ParticleFluxBC(FluxBCBase):
fenics format
bc_expr (fem.Expression): the expression of the particle flux that is used to
update the value_fenics
species_dependent_value (dict): a dictionary containing the species that the value. Example: {"name": species}
where "name" is the variable name in the callable value and species is a festim.Species object.
Usage:
Expand All @@ -160,6 +164,7 @@ class ParticleFluxBC(FluxBCBase):
>>> ParticleFluxBC(subdomain=my_subdomain, value=lambda t: 1 + t, species="H")
>>> ParticleFluxBC(subdomain=my_subdomain, value=lambda T: 1 + T, species="H")
>>> ParticleFluxBC(subdomain=my_subdomain, value=lambda x, t: 1 + x[0] + t, species="H")
>>> ParticleFluxBC(subdomain=my_subdomain, value=lambda c1: 2 * c1**2, species="H", species_dependent_value={"c1": species1})
"""

def __init__(self, subdomain, value, species, species_dependent_value={}):
Expand Down
12 changes: 10 additions & 2 deletions festim/exports/surface_flux.py
Original file line number Diff line number Diff line change
Expand Up @@ -25,13 +25,21 @@ def __init__(
) -> None:
super().__init__(field, surface, filename)

def compute(self, n, ds):
def compute(self, ds):
"""Computes the value of the surface flux at the surface
Args:
n (ufl.geometry.FacetNormal): normal vector to the surface
ds (ufl.Measure): surface measure of the model
"""

# obtain mesh normal from field
# if case multispecies, solution is an index, use sub_function_space
if isinstance(self.field.solution, ufl.indexed.Indexed):
mesh = self.field.sub_function_space.mesh
else:
mesh = self.field.solution.function_space.mesh
n = ufl.FacetNormal(mesh)

self.value = fem.assemble_scalar(
fem.form(
-self.D
Expand Down
10 changes: 4 additions & 6 deletions festim/hydrogen_transport_problem.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
import ufl
from mpi4py import MPI
import numpy as np
import tqdm.autonotebook
import tqdm.auto as tqdm
import festim as F

from dolfinx.mesh import meshtags
Expand Down Expand Up @@ -747,13 +747,14 @@ def run(self):

if self.settings.transient:
# Solve transient
self.progress = tqdm.autonotebook.tqdm(
self.progress = tqdm.tqdm(
desc="Solving H transport problem",
total=self.settings.final_time,
unit_scale=True,
)
while self.t.value < self.settings.final_time:
self.iterate()
self.progress.refresh() # refresh progress bar to show 100%
else:
# Solve steady-state
self.solver.solve(self.u)
Expand Down Expand Up @@ -820,10 +821,7 @@ def post_processing(self):
for export in self.exports:
# TODO if export type derived quantity
if isinstance(export, F.SurfaceQuantity):
export.compute(
self.mesh.n,
self.ds,
)
export.compute(self.ds)
# update export data
export.t.append(float(self.t))

Expand Down
27 changes: 25 additions & 2 deletions test/test_flux_bc.py
Original file line number Diff line number Diff line change
Expand Up @@ -88,14 +88,37 @@ def test_create_value_fenics_value(value, expected_value):
assert np.isclose(bc.value_fenics.value, expected_value)


def test_create_value_fenics_dependent_conc():
"""Test that the value_fenics of ParticleFluxBC is set correctly when the value is dependent on the concentration"""
# BUILD
left = F.SurfaceSubdomain1D(1, x=0)
my_species = F.Species("test")
my_species.solution = F.as_fenics_constant(12, mesh)
T = F.as_fenics_constant(1, mesh)
t = F.as_fenics_constant(0, mesh)
bc = F.ParticleFluxBC(
subdomain=left,
value=lambda c: 1.0 + c,
species=my_species,
species_dependent_value={"c": my_species},
)

# RUN
bc.create_value_fenics(mesh, T, t)

# TEST
assert isinstance(bc.value_fenics, ufl.core.expr.Expr)
assert bc.value_fenics == 1.0 + my_species.solution


def test_value_fenics_setter_error():
left = F.SurfaceSubdomain1D(1, x=0)
my_species = F.Species("test")
bc = F.ParticleFluxBC(subdomain=left, value=1.0, species=my_species)

with pytest.raises(
TypeError,
match="Value must be a dolfinx.fem.Function, dolfinx.fem.Constant, or a np.ndarray not <class 'str'>",
match="Value must be a dolfinx.fem.Function, dolfinx.fem.Constant, np.ndarray or ufl.core.expr.Expr not <class 'str'>",
):
bc.value_fenics = "coucou"

Expand Down Expand Up @@ -237,7 +260,7 @@ def test_create_value_fenics_type_HeatFluxBC(value, expected_type):
(lambda t: 50.0 if t < 1 else 0.0, 50),
],
)
def test_create_value_fenics_value(value, expected_value):
def test_create_value_fenics_value_HeatFluxBC(value, expected_value):
"""Test that"""
# BUILD
left = F.SurfaceSubdomain1D(1, x=0)
Expand Down
2 changes: 1 addition & 1 deletion test/test_surface_quantity.py
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,7 @@ def test_surface_flux_export_compute():
my_export.D = D

# RUN
my_export.compute(n=my_mesh.n, ds=ds)
my_export.compute(ds=ds)

# TEST
# flux = -D grad(c)_ \cdot n = -D dc/dx = -D * 4 * x
Expand Down

0 comments on commit 936738d

Please sign in to comment.