From 6cfcc39c813480ac9614962063a6c8191a1223cf Mon Sep 17 00:00:00 2001 From: johannesring Date: Thu, 7 Sep 2023 15:40:24 +0200 Subject: [PATCH] =?UTF-8?q?Added=20'local=5Frhs'=20option=20to=20'local=5F?= =?UTF-8?q?project'=20and=20'calculate=5Fand=5Fprint=5Fflow=5Fproperties'?= =?UTF-8?q?=20functions=20as=20recommended=20by=20J=C3=B8rgen?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/fsipy/simulations/simulation_common.py | 15 +++++++++++---- 1 file changed, 11 insertions(+), 4 deletions(-) diff --git a/src/fsipy/simulations/simulation_common.py b/src/fsipy/simulations/simulation_common.py index 1e8a9110..e810c457 100644 --- a/src/fsipy/simulations/simulation_common.py +++ b/src/fsipy/simulations/simulation_common.py @@ -211,13 +211,15 @@ def peval(f: Function, x: Union[float, np.ndarray]) -> np.ndarray: return yglob -def local_project(f: Function, V: FunctionSpace) -> Function: +def local_project(f: Function, V: FunctionSpace, local_rhs: bool = False) -> Function: """ Project a given function 'f' onto a finite element function space 'V' in a reusable way. Args: f (Function): The function to be projected onto 'V'. V (FunctionSpace): The finite element function space. + local_rhs (bool, optional): If True, solve using a local right-hand side assembly. + If False (default), solve using a global right-hand side assembly. Returns: Function: The projected solution in 'V'. @@ -229,13 +231,16 @@ def local_project(f: Function, V: FunctionSpace) -> Function: solver = LocalSolver(a_proj, b_proj) solver.factorize() projected_u = Function(V) - solver.solve_local_rhs(projected_u) + if local_rhs: + solver.solve_local_rhs(projected_u) + else: + solver.solve_global_rhs(projected_u) return projected_u def calculate_and_print_flow_properties(dt: float, mesh: Mesh, v: Function, inlet_area: float, mu_f: float, - n: Expression, dsi: Measure) -> None: + n: Expression, dsi: Measure, local_rhs: bool = False) -> None: """ Calculate and print flow properties. @@ -247,13 +252,15 @@ def calculate_and_print_flow_properties(dt: float, mesh: Mesh, v: Function, inle mu_f (float): Fluid dynamic viscosity. n (dolfin.Expression): FacetNormal expression. dsi (dolfin.Measure): Measure for inlet boundary. + local_rhs (bool, optional): If True, solve using a local right-hand side assembly. + If False (default), solve using a global right-hand side assembly. Returns: None """ # Calculate the DG vector of velocity magnitudes DG = FunctionSpace(mesh, "DG", 0) - V_vector = local_project(sqrt(inner(v, v)), DG).vector().get_local() + V_vector = local_project(sqrt(inner(v, v)), DG, local_rhs).vector().get_local() h = mesh.hmin() # Calculate flow rate at the inlet