Implementing code with Dirichlet and Neumann conditions together on same problem

dolfinx
1

I try following this “Combining Dirichlet and Neumann conditions” from tutorial Combining Dirichlet and Neumann conditions — FEniCSx tutorial because my boundary conditions on left side is dirichlet x[0]=1 and top,bottom&right are neumann x[1] = 0 but it has error, Could anyone help me with this how to correct it, Thank you.

#mesh
mesh = create_unit_square(MPI.COMM_WORLD, 10, 10)
V = FunctionSpace(mesh, ("Lagrange", 1))
u = TrialFunction(V)
v = TestFunction(V)
a = ufl.inner(ufl.grad(u), ufl.grad(v)) * ufl.dx

#dirichletBC
def u_exact(x):
    return 1

def boundary_D(c):
    return np.isclose(x[0], 0)

dofs_D = locate_dofs_geometrical(V, boundary_D)
u_bc = Function(V)
u_bc.interpolate(u_exact)
bc = [dirichletbc(u_bc, dofs_D)]

#NeumannBC
x = SpatialCoordinate(mesh)
g = Constant(mesh, default_scalar_type(0))
f = Constant(mesh, default_scalar_type(-6))
L = ufl.inner(f,v)*ufl.dx - ufl.inner(g,v) *ufl.ds

from dolfinx.fem.petsc import LinearProblem
problem = fem.petsc.LinearProblem(a, L, bcs=bc, petsc_options={"ksp_type": "preonly", "pc_type": "lu"})
uh = problem.solve()
2

Please elaborate on what you mean by

What is the error?

3

Sorry, this is an error for #dirichletBC, I try to make the boundary on the left side to be dirichlet with value 1

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
File /usr/local/dolfinx-complex/lib/python3.10/dist-packages/dolfinx/fem/function.py:397, in Function.interpolate(self, u, cells, nmm_interpolation_data)
    395 try:
    396     # u is a Function or Expression (or pointer to one)
--> 397     _interpolate(u, cells)
    398 except TypeError:
    399     # u is callable

File /usr/lib/python3.10/functools.py:889, in singledispatch.<locals>.wrapper(*args, **kw)
    886     raise TypeError(f'{funcname} requires at least '
    887                     '1 positional argument')
--> 889 return dispatch(args[0].__class__)(*args, **kw)

File /usr/local/dolfinx-complex/lib/python3.10/dist-packages/dolfinx/fem/function.py:373, in Function.interpolate.<locals>._interpolate(u, cells)
    372 """Interpolate a cpp.fem.Function"""
--> 373 self._cpp_object.interpolate(u, cells, nmm_interpolation_data)

TypeError: interpolate(): incompatible function arguments. The following argument types are supported:
    1. (self: dolfinx.cpp.fem.Function_complex128, f: numpy.ndarray[numpy.complex128], cells: numpy.ndarray[numpy.int32]) -> None
    2. (self: dolfinx.cpp.fem.Function_complex128, u: dolfinx.cpp.fem.Function_complex128, cells: numpy.ndarray[numpy.int32], nmm_interpolation_data: Tuple[List[int], List[int], List[float], List[int]]) -> None
    3. (self: dolfinx.cpp.fem.Function_complex128, expr: dolfinx::fem::Expression<std::complex<double>, double>, cells: numpy.ndarray[numpy.int32]) -> None

Invoked with: <dolfinx.cpp.fem.Function_complex128 object at 0x7f08000ff3f0>, <function u_exact at 0x7f07c06788b0>, array([  0,   1,   2,   3,   4,   5,   6,   7,   8,   9,  10,  11,  12,
        13,  14,  15,  16,  17,  18,  19,  20,  21,  22,  23,  24,  25,
        26,  27,  28,  29,  30,  31,  32,  33,  34,  35,  36,  37,  38,
        39,  40,  41,  42,  43,  44,  45,  46,  47,  48,  49,  50,  51,
        52,  53,  54,  55,  56,  57,  58,  59,  60,  61,  62,  63,  64,
        65,  66,  67,  68,  69,  70,  71,  72,  73,  74,  75,  76,  77,
        78,  79,  80,  81,  82,  83,  84,  85,  86,  87,  88,  89,  90,
        91,  92,  93,  94,  95,  96,  97,  98,  99, 100, 101, 102, 103,
       104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116,
       117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 128, 129,
       130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142,
       143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, 154, 155,
       156, 157, 158, 159, 160, 161, 162, 163, 164, 165, 166, 167, 168,
       169, 170, 171, 172, 173, 174, 175, 176, 177, 178, 179, 180, 181,
       182, 183, 184, 185, 186, 187, 188, 189, 190, 191, 192, 193, 194,
       195, 196, 197, 198, 199], dtype=int32), ((), (), (), ())

Did you forget to `#include <pybind11/stl.h>`? Or <pybind11/complex.h>,
<pybind11/functional.h>, <pybind11/chrono.h>, etc. Some automatic
conversions are optional and require extra headers to be included
when compiling your pybind11 module.

During handling of the above exception, another exception occurred:

IndexError                                Traceback (most recent call last)
Cell In[54], line 12
     10 dofs_D = locate_dofs_geometrical(V, boundary_D)
     11 u_bc = Function(V)
---> 12 u_bc.interpolate(u_exact)
     13 bc = [dirichletbc(u_bc, dofs_D)]

File /usr/local/dolfinx-complex/lib/python3.10/dist-packages/dolfinx/fem/function.py:402, in Function.interpolate(self, u, cells, nmm_interpolation_data)
    400 assert callable(u)
    401 x = _cpp.fem.interpolation_coords(self._V.element, self._V.mesh.geometry, cells)
--> 402 self._cpp_object.interpolate(np.asarray(u(x), dtype=self.dtype), cells)

IndexError: invalid axis: 0 (ndim = 0)
4

This should be a vectorized function,
i.e.
It takes in a (3, num_points) array and returns an array of length num_points.

def u_exact(x):
    return np.ones(x.shape[1], dtype=np.float64)

Would work.

You could Also just change the array data directly
u_bc.x.array[:] = 1

5

It worked :laughing: thank you so much but i still dont understand much about neumann conditions as i want to make the 3 sides of boundary conditions except the left side to be a constant value(neumann conditions). i’ve only seen tutorials that specify top and bottom as neumann but never seen any of those similar to my problem.

from original post

#NeumannBC
x = SpatialCoordinate(mesh)
g = Constant(mesh, default_scalar_type(0))
f = Constant(mesh, default_scalar_type(-6))
L = ufl.inner(f,v)*ufl.dx - ufl.inner(g,v) *ufl.ds

there is no error occurred, is this correct?? i’m not sure.

6

There is no fundamental difference in specifying top/bottom or left/right. Simply change how you locate the boundary facets that are input to the Dirichlet-bc or dolfinx.mesh.meshtags-constructor (if using it as subdomain_data in ufl.ds).

Additionally, please note that ufl.ds with no subdomain_data and subdomain_id is an integration measure over the whole boundary.

However, as for every Dirichlet boundary condition specified and applied to the variational form, these boundary integrals are removed from the assemble matrix/vector. See for instance; Application of Dirichlet boundary conditions — FEniCS Tutorial @ Sorbonne for more details regarding Dirichlet boundary conditions.

1 Like