# Singularity in gradient caused by Dirichlet boundary condition

I am looking for a mathematical explanation for the singularity caused by a Dirichlet boundary condition partially imposed at a boundary.

For instance $$\nabla^2u=0 ~ \text{in}~\Omega$$ where $$\Omega$$ is a rectangle of 10 by 4. Dirichlet boundary condition $$u=0$$ is imposed at $$y=0 ~\text{and}~ y = 4 ~\text{and}~ x < 2$$ and a flux is applied at the boundary $$x=10$$

Plotting the flux, you can see the singularity where the Dirichlet boundary condition ends

I plotted this with the following FEniCS code

from fenics import *

mesh = RectangleMesh(Point(0, 0), Point(10, 4), 400, 160, diagonal='crossed')

V = FunctionSpace(mesh, 'CG', 1)
u, v = TrialFunction(V), TestFunction(V)

def MyDirichlet(x, on_boundary):
return on_boundary and (x[0] < 2.0) and (x[0] > 0.0)

class MyNeumann(SubDomain):
def inside(self, x, on_boundary):
return on_boundary and near(x[0], 10.0)

subdomain = MeshFunction("size_t", mesh, mesh.topology().dim() - 1)
subdomain.set_all(0)
myneumann = MyNeumann()
myneumann.mark(subdomain, 1)
ds = Measure('ds', domain=mesh, subdomain_data=subdomain)

bc = DirichletBC(V, Constant(0.0), MyDirichlet)

L = Constant(1.0)*v*ds(1)

u_sol = Function(V)

solve(a==L, u_sol, bcs=bc)
File("sing_example.pvd") << u_sol
Vflux = VectorFunctionSpace(mesh, 'CG', 1)
File("sing_flux.pvd") << project(grad(u_sol), Vflux)


• Just out of curiosity, is there any particular reason to use this weird configuration of boundary conditions for Laplace equation? Jan 3, 2020 at 0:56
• No specific reason. Jan 3, 2020 at 3:10
• So, is this for research or right now just for learning and then moving on to solving a real problem in the future? Jan 3, 2020 at 3:37
• Both, the real problem is a convection diffusion problem where the advection field enters the right side of the rectangle and exists through the left side. The idea of placing the Dirichlet boundary conditions this way is to keep them away from the right side (the inlet). I have not thought of an alternative way to pose this real problem, but in any case I would like to know why this phenomenon happens. Jan 3, 2020 at 15:59
• As far as I know, you need to have at least some sort of boundary condition specified on every part of the boundary. Since you do not specify any conditions on $x=0$ or $(x>2)(x<10)$ I have no idea what Fenics will do. It probably depends on the details of solve.
– knl
Feb 1, 2020 at 21:15

As others have mentioned, the boundary-conditions are not clear. But such singularities in the gradient do occur in the neighborhood of a discontinuity in boundary conditions. Here's a mathematical explanation. Take the following prototypical boundary condition: $$u = 0$$ for $$y = 0, \ x > 0$$ and $$u = \pi$$ for $$y=0, \ x < 0$$. Then a solution of the Laplace equation for $$y\geq 0$$, which obeys the boundary-conditions and is bounded at infinity is simply $$u = \phi$$, where $$\phi$$ is the angle from the positive $$x$$ axis. The gradient is then $$\nabla u = \frac{1}{\rho}\mathbf{1}_\phi$$, where $$\rho$$ is the distance from the origin and $$\mathbf{1}_\phi$$ is a unit-vector in the $$\phi$$ direction, and it is singular near the origin. In this example, the flux is normal to the $$x$$ axis boundary, so either one of the Dirichlet conditions (but not both) may be readily replaced by a Neumann condition without altering the result.
To put it in more intuitive terms, the potential at $$y = 0^+$$ must jump discontinuously between two different values, therefore the gradient must be infinite near the jump.