I get the following type error when trying to define dirichletbc.
The code is a modified version of this example Setting multiple Dirichlet, Neumann, and Robin conditions
The mathematics:
Equation and problem definition
For a domain \Omega \subset \mathbb{R}^n with boundary \partial
\Omega = \Gamma_{D} \cup \Gamma_{N} \cup \Gamma_{M}, the Poisson equation with
particular boundary conditions reads:
\nabla^{2}u =0 \quad {\rm in} \ \Omega, \\ u = 1 \quad {\rm on} \ \Gamma_{D}, \\ \nabla u \cdot n = 0 \quad {\rm on} \ \Gamma_{N}, \\ \nabla u \cdot n = h \quad {\rm on} \ \Gamma_{M}
a(u, v) := \int_{\Omega} \nabla u \cdot \nabla v \, {\rm d} x + \int_{\Gamma_{M}} r u v \, {\rm d} x \\
L(v) := \int_{\Omega} 0 v \, {\rm d} x + \int_{\Gamma_{M}} h v \, {\rm d} s + \int_{\Gamma_{N}} 0 v \, {\rm d} s.
- \Omega = [0,2] \times [0,1] (a rectangle)
- \Gamma_{D} = \{(x, 0) \subset \partial \Omega\}
- \Gamma_{M} = \{(2, y) \cup(x, 1) \subset \partial \Omega\}
- \Gamma_{N} = \{(0, y) \subset \partial \Omega\}
- h = 0.1*(u-0) = r*(u-0)
The full code:
import numpy as np
import pyvista
from dolfinx.fem import (Constant, Function, FunctionSpace, assemble_scalar,
dirichletbc, form, locate_dofs_topological)
from dolfinx.fem.petsc import LinearProblem
from dolfinx.mesh import create_unit_square, locate_entities, meshtags, create_rectangle, CellType
from mpi4py import MPI
from petsc4py.PETSc import ScalarType
from ufl import (FacetNormal, Measure, SpatialCoordinate, TestFunction, TrialFunction,
div, dot, dx, grad, inner, lhs, rhs)
from dolfinx.io import XDMFFile
from dolfinx.plot import create_vtk_mesh
mesh = create_rectangle(comm=MPI.COMM_WORLD,
points=((0.0, 0.0), (2.0, 1.0)), n=(32, 16),
cell_type=CellType.triangle)
u_ex =Constant(mesh, ScalarType(1)) #dirichlet
x = SpatialCoordinate(mesh)
# Define physical parameters and boundary condtions
s = Constant(mesh, ScalarType(0))
f = Constant(mesh, ScalarType(0))
n = FacetNormal(mesh)
g = Constant(mesh, ScalarType(0))
kappa = Constant(mesh, ScalarType(1))
r = Constant(mesh, ScalarType(0.1))
# Define function space and standard part of variational form
V = FunctionSpace(mesh, ("CG", 1))
u, v = TrialFunction(V), TestFunction(V)
F = kappa * inner(grad(u), grad(v)) * dx - inner(f, v) * dx
boundaries = [(1, lambda x: np.isclose(x[1], 0)),
(2, lambda x: np.logical_or(np.isclose(x[0], 2),np.isclose(x[1], 1))),
(3, lambda x: np.isclose(x[0], 0))]
facet_indices, facet_markers = [], []
fdim = mesh.topology.dim - 1
for (marker, locator) in boundaries:
facets = locate_entities(mesh, fdim, locator)
facet_indices.append(facets)
facet_markers.append(np.full_like(facets, marker))
facet_indices = np.hstack(facet_indices).astype(np.int32)
facet_markers = np.hstack(facet_markers).astype(np.int32)
sorted_facets = np.argsort(facet_indices)
facet_tag = meshtags(mesh, fdim, facet_indices[sorted_facets], facet_markers[sorted_facets])
ds = Measure("ds", domain=mesh, subdomain_data=facet_tag)
# use a simple class to define bc
class BoundaryCondition():
def __init__(self, tipo, marker, values):
self._type = tipo
if tipo == "Dirichlet":
facets = facet_tag.find(marker)
dofs = locate_dofs_topological(V, fdim, facets)
self._bc = dirichletbc(value=Constant(mesh, ScalarType(1)), dofs=dofs)
elif tipo == "Neumann":
self._bc = inner(values, v) * ds(marker)
elif tipo == "Robin":
self._bc = values[0] * inner(u-values[1], v)* ds(marker)
else:
raise TypeError("Unknown boundary condition: {0:s}".format(type))
@property
def bc(self):
return self._bc
@property
def type(self):
return self._type
# Define the Dirichlet condition
boundary_conditions = [BoundaryCondition("Dirichlet", 1, u_ex),
BoundaryCondition("Robin", 2, (r, s)),
BoundaryCondition("Neumann", 3, g)]
Thank you!