How to Solve Poisson’s Equation with Robin Boundary Conditions

Poisson’s equation with Robin boundary conditions combines an interior source equation with a boundary law involving both the solution and its normal derivative. This type of condition models exchange with the surroundings: for example, heat transfer through a surface or diffusion through a permeable boundary.

The main points are the sign convention, the weak formulation, and the conditions guaranteeing uniqueness. The discussion below uses the standard convention

$$-\Delta u=f \quad \text{in } \Omega,$$

with outward unit normal $n$ and Robin condition

$$\frac{\partial u}{\partial n}+\alpha u=g \quad \text{on } \partial\Omega.$$

Here $\partial u/\partial n=\nabla u\cdot n$, $\alpha$ is a prescribed boundary coefficient, and $g$ is prescribed boundary data.

Deriving the weak formulation

Let $v$ be a sufficiently smooth test function. Multiply the PDE by $v$ and integrate over $\Omega$:

$$\int_\Omega (-\Delta u)v\,dx=\int_\Omega fv\,dx.$$

Green’s identity gives

$$\int_\Omega \nabla u\cdot\nabla v\,dx-\int_{\partial\Omega}\frac{\partial u}{\partial n}v\,ds =\int_\Omega fv\,dx.$$

The Robin condition says $\partial u/\partial n=g-\alpha u$. Substituting this expression produces

$$\int_\Omega \nabla u\cdot\nabla v\,dx +\int_{\partial\Omega}\alpha uv\,ds =\int_\Omega fv\,dx+\int_{\partial\Omega}gv\,ds.$$

Thus the weak problem is: find $u$ in an appropriate Sobolev space, usually $H^1(\Omega)$, such that this identity holds for every test function $v\in H^1(\Omega)$.

The corresponding bilinear form and linear functional are

$$a(u,v)=\int_\Omega \nabla u\cdot\nabla v\,dx +\int_{\partial\Omega}\alpha uv\,ds,$$

and

$$L(v)=\int_\Omega fv\,dx+\int_{\partial\Omega}gv\,ds.$$

This formulation is the starting point for finite-element methods.

Existence and uniqueness

For a bounded, connected Lipschitz domain, suppose $\alpha\in L^\infty(\partial\Omega)$ is nonnegative and satisfies $\alpha\geq\alpha_0>0$ on a subset of $\partial\Omega$ having positive surface measure. Then the bilinear form is coercive on $H^1(\Omega)$. For a disconnected domain, this positivity condition is needed on a subset of positive measure in each connected component. In practical terms, the gradient term controls variations of $u$, while the boundary term controls the otherwise-undetermined constant component.

Under standard boundedness assumptions—for example, $f\in L^2(\Omega)$, $g\in L^2(\partial\Omega)$, and $\alpha\in L^\infty(\partial\Omega)$—the Lax–Milgram theorem then gives a unique weak solution.

The special cases explain why the assumptions matter:

• Dirichlet condition: prescribing $u$ on the whole boundary removes constant freedom directly.
• Pure Neumann condition: if $\alpha=0$, adding a constant to a solution produces another solution. Solvability requires $$\int_\Omega f\,dx+\int_{\partial\Omega}g\,ds=0,$$ and the solution is unique only up to an additive constant.
• Positive Robin condition: if $\alpha>0$ on the whole boundary, constants are no longer invisible, so no Neumann-type compatibility condition is needed.
• Negative or sign-changing $\alpha$: the boundary term can destroy coercivity. A solution may fail to be unique, or the problem may become solvable only for data satisfying an additional orthogonality condition. This is an eigenvalue issue, not something that can be decided from the boundary formula alone.

If a different sign convention is used, such as $\Delta u=f$, every integration-by-parts formula should be rederived rather than copied. Likewise, some applications write the condition as $k\,\partial u/\partial n+h u=r$; dividing by $k$ is valid only where $k\neq0$.

Worked one-dimensional example

Consider the interval $0<x<1$ and solve

$$-u''(x)=1,$$

with homogeneous Robin conditions

$$-u'(0)+u(0)=0, \qquad u'(1)+u(1)=0.$$

These are the correct outward-normal conditions. At the left endpoint the outward normal points toward negative $x$, so $\partial u/\partial n=-u'(0)$; at the right endpoint it is $u'(1)$.

Integrating the differential equation twice gives

$$u(x)=-\frac{x^2}{2}+Cx+D.$$

Its derivative is $u'(x)=-x+C$. Apply the left boundary condition:

$$-u'(0)+u(0)=-C+D=0,$$

so $D=C$. Apply the right condition:

$$u'(1)+u(1)=(-1+C)+\left(-\frac12+C+D\right)=0.$$

Using $D=C$ gives

$$-\frac32+3C=0, \qquad C=D=\frac12.$$

Therefore

$$\boxed{u(x)=-\frac{x^2}{2}+\frac{x}{2}+\frac12}.$$

A quick check confirms the result: $-u''=1$, while $-u'(0)+u(0)=0$ and $u'(1)+u(1)=0$.

How the finite-element method uses Robin data

Suppose $u_h$ is an approximation in a finite-dimensional basis. The weak equation becomes a linear system

$$(K+R)\mathbf{U}=\mathbf{F}+\mathbf{G},$$

where $K$ comes from the volume integrals of gradient products and $R$ comes from the boundary integrals involving $\alpha$. The source vector $\mathbf F$ represents $f$, and $\mathbf G$ represents $g$.

Robin data therefore do not usually replace rows of the stiffness matrix in the way essential Dirichlet data do. Instead, they contribute naturally to the matrix and right-hand side through boundary integrals. In a finite-difference method, the boundary equation must be discretized separately, often using a one-sided approximation or a ghost point.

Common mistakes

1. Using $u'$ instead of the outward normal derivative. In one dimension, the derivative has opposite signs at the two endpoints.
2. Dropping the boundary integral. Integration by parts creates a boundary term; the Robin condition converts it into the $\alpha uv$ and $gv$ terms.
3. Applying the Neumann compatibility condition automatically. It is required for pure Neumann problems, but not generally when a positive Robin term is present.
4. Assuming any $\alpha$ guarantees uniqueness. Negative or sign-changing coefficients can lead to noncoercive or resonant problems.
5. Mixing PDE sign conventions. The weak form above is for $-\Delta u=f$. Starting with $\Delta u=f$ changes the signs.