From 44db013c021f58fbdc3a2525ab4d7ba50da47f49 Mon Sep 17 00:00:00 2001
From: karamoko <karamoko.samassa@inria.fr>
Date: Mon, 16 Jan 2023 17:06:42 +0100
Subject: [PATCH] fixing sympy expression creation from a string
---
src/bvpy/domains/abstract.py | 1 +
src/bvpy/ibvp.py | 40 +++++++++----------
src/bvpy/solvers/linear.py | 2 +-
src/bvpy/solvers/time_integration.py | 20 ++++++----
src/bvpy/utils/pre_processing.py | 8 ++--
.../bvpy_tutorial_5_reaction_diffusion.ipynb | 7 ++--
6 files changed, 43 insertions(+), 35 deletions(-)
diff --git a/src/bvpy/domains/abstract.py b/src/bvpy/domains/abstract.py
index 32917b9..2cea58d 100644
--- a/src/bvpy/domains/abstract.py
+++ b/src/bvpy/domains/abstract.py
@@ -403,6 +403,7 @@ class AbstractDomain(ABC):
self.model.addPhysicalGroup(dim, [self.volumes[0]])
self.model.mesh.generate(dim)
+
self.mesh, cell_markers, facet_markers = gmshio.model_to_mesh(self.model, MPI.COMM_SELF, 0)
# self.bdata = _generate_default_boundary_data(self.mesh)
self.ds = ufl.Measure('ds', domain=self.mesh)
diff --git a/src/bvpy/ibvp.py b/src/bvpy/ibvp.py
index 672889b..1260a39 100644
--- a/src/bvpy/ibvp.py
+++ b/src/bvpy/ibvp.py
@@ -75,7 +75,7 @@ class IBVP(BVP):
self.time = Time()
self.scheme = scheme(vform)
- self.previous_solution = self.solution.copy(deepcopy=True)
+ self.previous_solution = self.solution.copy()
self.scheme.set_previous_solution(self.previous_solution)
self._progressbar = ProgressBar()
self.solution_steps = None
@@ -112,31 +112,31 @@ class IBVP(BVP):
if store_steps:
self.solution_steps = OrderedDict()
self.solution_steps[self.time.current] =\
- self.solution.copy(deepcopy=True)
+ self.solution.copy()
if filename is not None:
# file = fe.XDMFFile(MPI_COMM_WORLD, filename)
- file = io.XDMFFile(MPI_COMM_WORLD, filename)
- file.parameters['functions_share_mesh'] = True
- file.write(self.solution, self.time.current)
- logger.setLevel(logging.WARNING)
- self._progressbar.set_range(t_i, t_f)
-
- while t_f-self.time.current > dt:
- self.step(dt, **kwargs)
- file.write(self.solution, self.time.current)
- if store_steps:
- self.solution_steps[self.time.current] =\
- self.solution.copy(deepcopy=True)
+ with io.XDMFFile(MPI_COMM_WORLD, filename, "w") as file:
+ file.write_mesh(self.solution.function_space.mesh)
+ # file.parameters['functions_share_mesh'] = True
+ file.write_function(self.solution, self.time.current)
+ logger.setLevel(logging.WARNING)
+ self._progressbar.set_range(t_i, t_f)
+
+ while t_f-self.time.current > dt:
+ self.step(dt, **kwargs)
+ file.write_function(self.solution, self.time.current)
+ if store_steps:
+ self.solution_steps[self.time.current] =\
+ self.solution.copy()
+ self._progressbar.update(self.time.current)
+ self._progressbar.show()
+
+ self.step(t_f-self.time.current, **kwargs)
+ file.write_function(self.solution, self.time.current)
self._progressbar.update(self.time.current)
self._progressbar.show()
- self.step(t_f-self.time.current, **kwargs)
- file.write(self.solution, self.time.current)
- self._progressbar.update(self.time.current)
- self._progressbar.show()
- file.close()
-
else:
logger.setLevel(logging.WARNING)
self._progressbar.set_range(t_i, t_f)
diff --git a/src/bvpy/solvers/linear.py b/src/bvpy/solvers/linear.py
index 746e90a..01ce774 100644
--- a/src/bvpy/solvers/linear.py
+++ b/src/bvpy/solvers/linear.py
@@ -96,7 +96,7 @@ class LinearSolver():
logger.info(' Preconditioner: '
+ str(self.parameters['preconditioner']))
logger.info(' Size : '
- + str(self._u.function_space().dim()))
+ + str(self._u.function_space.mesh.geometry.dim))
logger.info(' Solving LinearProblem ...')
# A = fe.assemble(self._a)
diff --git a/src/bvpy/solvers/time_integration.py b/src/bvpy/solvers/time_integration.py
index b232607..67565ef 100644
--- a/src/bvpy/solvers/time_integration.py
+++ b/src/bvpy/solvers/time_integration.py
@@ -22,6 +22,7 @@
# import fenics as fe
from ufl.classes import Zero
import ufl
+from dolfinx import fem
class Time():
"""Container class for Time.
@@ -104,19 +105,22 @@ class ImplicitEuler(FirstOrderScheme):
def apply(self, u, v, sol):
lhs, rhs = self._vform.get_form(u, v, sol)
if isinstance(rhs, Zero):
- dudt = ufl.dot((sol-self.previous_solution)
- / ufl.Constant(self.dt), v) * self.dx
+ # dudt = ufl.dot((sol-self.previous_solution)
+ # / ufl.Constant(self.dt), v) * self.dx
+ dudt = ufl.dot((sol-self.previous_solution)/self.dt, v) * self.dx
lhs += dudt
else:
- dudt = ufl.dot((u-self.previous_solution)
- / ufl.Constant(self.dt), v) * self.dx
+ # dudt = ufl.dot((u-self.previous_solution)
+ # / ufl.Constant(self.dt), v) * self.dx
+ dudt = ufl.dot((u-self.previous_solution)/self.dt, v) * self.dx
- dudt_lhs, dudt_rhs = ufl.lhs(dudt), ufl.rhs(dudt)
-
- lhs += dudt_lhs
- rhs += dudt_rhs
+ dudt_lhs, dudt_rhs = fem.form(ufl.lhs(dudt)), fem.form(ufl.rhs(dudt))
+ # lhs += dudt_lhs
+ # rhs += dudt_rhs
+ lhs = dudt_lhs
+ rhs = dudt_rhs
return lhs, rhs
diff --git a/src/bvpy/utils/pre_processing.py b/src/bvpy/utils/pre_processing.py
index 2208fa9..a9b1f24 100644
--- a/src/bvpy/utils/pre_processing.py
+++ b/src/bvpy/utils/pre_processing.py
@@ -24,7 +24,8 @@ import logging
# import fenics as fe
import numpy as np
-from sympy import Function, Symbol, lambdify, symbols, sympify
+import sympy
+from sympy import Function, Symbol
from sympy.parsing.sympy_parser import parse_expr
from types import FunctionType
@@ -45,8 +46,9 @@ def _expression_from_string(source, functionspace, degree, **kwargs):
for src in source:
expr.append(str(parse_expr(src, local_dict=_XYZ, evaluate=False)))
else:
- # expr = str(parse_expr(source, local_dict=_XYZ, evaluate=False))
- expr = sympify(source)
+ expr = parse_expr(source, local_dict=_XYZ, evaluate=False)
+ x = Symbol('x')
+ expr = sympy.lambdify(x, expr)
#TODO: is this part useful anymore ?
diff --git a/tutorials/bvpy_tutorial_5_reaction_diffusion.ipynb b/tutorials/bvpy_tutorial_5_reaction_diffusion.ipynb
index 6d49e2f..37d7ced 100644
--- a/tutorials/bvpy_tutorial_5_reaction_diffusion.ipynb
+++ b/tutorials/bvpy_tutorial_5_reaction_diffusion.ipynb
@@ -121,8 +121,9 @@
"outputs": [],
"source": [
"from bvpy.utils.pre_processing import create_expression\n",
+ "import sympy\n",
"\n",
- "init = create_expression(\"1+tanh(-x/.1)\", degree=4)"
+ "init = create_expression(\"1+tanh(-x/.1)\", degree=4)\n"
]
},
{
@@ -537,7 +538,7 @@
],
"metadata": {
"kernelspec": {
- "display_name": "Python 3.10.6 ('bvpyx-dev')",
+ "display_name": "Python 3 (ipykernel)",
"language": "python",
"name": "python3"
},
@@ -615,7 +616,7 @@
},
"vscode": {
"interpreter": {
- "hash": "4ca261fcaa8fbe534f2a6fc835c771b3b71066f7c237a4d8268a0ea3f65d5ab3"
+ "hash": "433263323d5545258b42c6f78b5e3d676fbc72f662b59d217ea7a5b51f0f4213"
}
}
},
--
GitLab