Re: [firedrake] imposing Dirichlet BCS
Hi Francis, Usually we use the boundary markers on the incoming mesh. Where are you getting your more general mesh from? Regards, David On Wed, 10 May 2017 at 15:37 Francis Poulin <fpoulin@uwaterloo.ca> wrote:
Hello,
We are trying to translate part of the mixed-poisson Fenics demo into Firedrake. I know that there is already a Firedrake demo that does this but for the special case of a square. I am hoping we can do this more generally.
Below are the lines of code from Fenics. Can someone suggest how this could be implemented in Firedrake?
Cheers, Francis
https://fenics.readthedocs.io/projects/dolfin/en/latest/demos/mixed-poisson/...
# Define function G such that G \cdot n = gclass BoundarySource(Expression): def __init__(self, mesh, **kwargs): self.mesh = mesh def eval_cell(self, values, x, ufc_cell): cell = Cell(self.mesh, ufc_cell.index) n = cell.normal(ufc_cell.local_facet) g = sin(5*x[0]) values[0] = g*n[0] values[1] = g*n[1] def value_shape(self): return (2,) G = BoundarySource(mesh, degree=2) # Define essential boundarydef boundary(x): return x[1] < DOLFIN_EPS or x[1] > 1.0 - DOLFIN_EPS
bc = DirichletBC(W.sub(0), G, boundary)
------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo
email: fpoulin@uwaterloo.ca Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637 <(519)%20888-4567>
-- Dr David Ham Department of Mathematics Imperial College London
Hello David, At the moment I am using UnitSquareMesh but at some point I would like to consider different geometries. What I tried is using n = FacetNormal(mesh) x, y = SpatialCoordinate(mesh) g = (sin(pi*x)*sin(pi*y))*n[0] + cos(pi*x)*cos(pi*y)*n[1] g1 = g*n[0] g2 = g*n[1] bcs = DirichletBC(W.sub(0), as_vector([g1, g2]), "on_boundary") but that failed miserably with the following, (firedrake) fpoulin@fpoulin-Gazelle:~/software/poulinrhebergen/codes/firedrake/stokes$ python stokes.py UFL:ERROR Integral of type cell cannot contain a ReferenceNormal. Traceback (most recent call last): File "stokes.py", line 96, in <module> bcs = DirichletBC(W.sub(0), as_vector([g1, g2]), "on_boundary") File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/bcs.py", line 57, in __init__ self.function_arg = g File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/bcs.py", line 125, in function_arg g = projection.project(g, self._function_space) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/projection.py", line 101, in project form_compiler_parameters=form_compiler_parameters) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/solving.py", line 124, in solve _solve_varproblem(*args, **kwargs) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem appctx=appctx) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/variational_solver.py", line 286, in __init__ super(LinearVariationalSolver, self).__init__(*args, **kwargs) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/variational_solver.py", line 156, in __init__ pre_function_callback=pre_f_callback) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/solving_utils.py", line 333, in __init__ form_compiler_parameters=fcp) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/assemble.py", line 143, in create_assembly_callable collect_loops=True) File "<decorator-gen-280>", line 2, in _assemble File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/utils.py", line 62, in wrapper return f(*args, **kwargs) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/assemble.py", line 192, in _assemble kernels = tsfc_interface.compile_form(f, "form", parameters=form_compiler_parameters, inverse=inverse) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 193, in compile_form number_map).kernels File "/home/fpoulin/software/firedrake/src/PyOP2/pyop2/caching.py", line 200, in __new__ obj = make_obj() File "/home/fpoulin/software/firedrake/src/PyOP2/pyop2/caching.py", line 190, in make_obj obj.__init__(*args, **kwargs) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 112, in __init__ tree = tsfc_compile_form(form, prefix=name, parameters=parameters) File "/home/fpoulin/software/firedrake/src/tsfc/tsfc/driver.py", line 46, in compile_form fd = ufl_utils.compute_form_data(form) File "/home/fpoulin/software/firedrake/src/tsfc/tsfc/ufl_utils.py", line 56, in compute_form_data do_estimate_degrees=do_estimate_degrees, File "/home/fpoulin/software/firedrake/src/ufl/ufl/algorithms/compute_form_data.py", line 382, in compute_form_data _check_facet_geometry(self.integral_data) File "/home/fpoulin/software/firedrake/src/ufl/ufl/algorithms/compute_form_data.py", line 168, in _check_facet_geometry error("Integral of type %s cannot contain a %s." % (it, cls.__name__)) File "/home/fpoulin/software/firedrake/src/ufl/ufl/log.py", line 172, in error raise self._exception_type(self._format_raw(*message)) ufl.log.UFLException: Integral of type cell cannot contain a ReferenceNormal. ------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo email: fpoulin@uwaterloo.ca Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637 ________________________________ From: firedrake-bounces@imperial.ac.uk [firedrake-bounces@imperial.ac.uk] on behalf of David Ham [David.Ham@imperial.ac.uk] Sent: Wednesday, May 10, 2017 10:41 AM To: firedrake Subject: Re: [firedrake] imposing Dirichlet BCS Hi Francis, Usually we use the boundary markers on the incoming mesh. Where are you getting your more general mesh from? Regards, David On Wed, 10 May 2017 at 15:37 Francis Poulin <fpoulin@uwaterloo.ca<mailto:fpoulin@uwaterloo.ca>> wrote: Hello, We are trying to translate part of the mixed-poisson Fenics demo into Firedrake. I know that there is already a Firedrake demo that does this but for the special case of a square. I am hoping we can do this more generally. Below are the lines of code from Fenics. Can someone suggest how this could be implemented in Firedrake? Cheers, Francis https://fenics.readthedocs.io/projects/dolfin/en/latest/demos/mixed-poisson/... # Define function G such that G \cdot n = g class BoundarySource(Expression): def __init__(self, mesh, **kwargs): self.mesh = mesh def eval_cell(self, values, x, ufc_cell): cell = Cell(self.mesh, ufc_cell.index) n = cell.normal(ufc_cell.local_facet) g = sin(5*x[0]) values[0] = g*n[0] values[1] = g*n[1] def value_shape(self): return (2,) G = BoundarySource(mesh, degree=2) # Define essential boundary def boundary(x): return x[1] < DOLFIN_EPS or x[1] > 1.0 - DOLFIN_EPS bc = DirichletBC(W.sub(0), G, boundary) ------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo email: fpoulin@uwaterloo.ca<mailto:fpoulin@uwaterloo.ca> Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637<tel:(519)%20888-4567> -- Dr David Ham Department of Mathematics Imperial College London
This is a duplicate of issue #981: "Using FacetNormal (and other things only defined on the boundary) in DirichletBCs currently not possible" https://github.com/firedrakeproject/firedrake/issues/981 On 10/05/17 15:46, Francis Poulin wrote:
Hello David,
At the moment I am using UnitSquareMesh but at some point I would like to consider different geometries.
What I tried is using
n = FacetNormal(mesh)
x, y = SpatialCoordinate(mesh) g = (sin(pi*x)*sin(pi*y))*n[0] + cos(pi*x)*cos(pi*y)*n[1] g1 = g*n[0] g2 = g*n[1]
bcs = DirichletBC(W.sub(0), as_vector([g1, g2]), "on_boundary")
but that failed miserably with the following,
(firedrake) fpoulin@fpoulin-Gazelle:~/software/poulinrhebergen/codes/firedrake/stokes$ python stokes.py UFL:ERROR Integral of type cell cannot contain a ReferenceNormal. Traceback (most recent call last): File "stokes.py", line 96, in <module> bcs = DirichletBC(W.sub(0), as_vector([g1, g2]), "on_boundary") File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/bcs.py", line 57, in __init__ self.function_arg = g File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/bcs.py", line 125, in function_arg g = projection.project(g, self._function_space) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/projection.py", line 101, in project form_compiler_parameters=form_compiler_parameters) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/solving.py", line 124, in solve _solve_varproblem(*args, **kwargs) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/solving.py", line 153, in _solve_varproblem appctx=appctx) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/variational_solver.py", line 286, in __init__ super(LinearVariationalSolver, self).__init__(*args, **kwargs) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/variational_solver.py", line 156, in __init__ pre_function_callback=pre_f_callback) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/solving_utils.py", line 333, in __init__ form_compiler_parameters=fcp) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/assemble.py", line 143, in create_assembly_callable collect_loops=True) File "<decorator-gen-280>", line 2, in _assemble File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/utils.py", line 62, in wrapper return f(*args, **kwargs) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/assemble.py", line 192, in _assemble kernels = tsfc_interface.compile_form(f, "form", parameters=form_compiler_parameters, inverse=inverse) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 193, in compile_form number_map).kernels File "/home/fpoulin/software/firedrake/src/PyOP2/pyop2/caching.py", line 200, in __new__ obj = make_obj() File "/home/fpoulin/software/firedrake/src/PyOP2/pyop2/caching.py", line 190, in make_obj obj.__init__(*args, **kwargs) File "/home/fpoulin/software/firedrake/src/firedrake/firedrake/tsfc_interface.py", line 112, in __init__ tree = tsfc_compile_form(form, prefix=name, parameters=parameters) File "/home/fpoulin/software/firedrake/src/tsfc/tsfc/driver.py", line 46, in compile_form fd = ufl_utils.compute_form_data(form) File "/home/fpoulin/software/firedrake/src/tsfc/tsfc/ufl_utils.py", line 56, in compute_form_data do_estimate_degrees=do_estimate_degrees, File "/home/fpoulin/software/firedrake/src/ufl/ufl/algorithms/compute_form_data.py", line 382, in compute_form_data _check_facet_geometry(self.integral_data) File "/home/fpoulin/software/firedrake/src/ufl/ufl/algorithms/compute_form_data.py", line 168, in _check_facet_geometry error("Integral of type %s cannot contain a %s." % (it, cls.__name__)) File "/home/fpoulin/software/firedrake/src/ufl/ufl/log.py", line 172, in error raise self._exception_type(self._format_raw(*message)) ufl.log.UFLException: Integral of type cell cannot contain a ReferenceNormal.
------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo
email: fpoulin@uwaterloo.ca Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637
------------------------------------------------------------------------ *From:* firedrake-bounces@imperial.ac.uk [firedrake-bounces@imperial.ac.uk] on behalf of David Ham [David.Ham@imperial.ac.uk] *Sent:* Wednesday, May 10, 2017 10:41 AM *To:* firedrake *Subject:* Re: [firedrake] imposing Dirichlet BCS
Hi Francis,
Usually we use the boundary markers on the incoming mesh. Where are you getting your more general mesh from?
Regards,
David
On Wed, 10 May 2017 at 15:37 Francis Poulin <fpoulin@uwaterloo.ca <mailto:fpoulin@uwaterloo.ca>> wrote:
Hello,
We are trying to translate part of the mixed-poisson Fenics demo into Firedrake. I know that there is already a Firedrake demo that does this but for the special case of a square. I am hoping we can do this more generally.
Below are the lines of code from Fenics. Can someone suggest how this could be implemented in Firedrake?
Cheers, Francis
https://fenics.readthedocs.io/projects/dolfin/en/latest/demos/mixed-poisson/...
# Define function G such that G \cdot n = g class BoundarySource(Expression): def __init__(self, mesh, **kwargs): self.mesh = mesh def eval_cell(self, values, x, ufc_cell): cell = Cell(self.mesh, ufc_cell.index) n = cell.normal(ufc_cell.local_facet) g = sin(5*x[0]) values[0] = g*n[0] values[1] = g*n[1] def value_shape(self): return (2,)
G = BoundarySource(mesh, degree=2)
# Define essential boundary def boundary(x): return x[1] < DOLFIN_EPS or x[1] > 1.0 - DOLFIN_EPS
bc = DirichletBC(W.sub(0), G, boundary)
------------------ Francis Poulin Associate Professor Department of Applied Mathematics University of Waterloo
email: fpoulin@uwaterloo.ca <mailto:fpoulin@uwaterloo.ca> Web: https://uwaterloo.ca/poulin-research-group/ Telephone: +1 519 888 4567 x32637 <tel:%28519%29%20888-4567>
-- Dr David Ham Department of Mathematics Imperial College London
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
participants (3)
- 
                
                David Ham
- 
                
                Francis Poulin
- 
                
                Miklós Homolya