Boundary Conditions

The BoundaryConditions class collects all boundary constraints for a PDE problem. findiff supports Dirichlet, Neumann and Robin (mixed) boundary conditions in any number of dimensions.

Setup

import numpy as np
from findiff import Diff, PDE, BoundaryConditions

# 1D grid
x = np.linspace(0, 1, 100)
dx = x[1] - x[0]

# 2D grid
x = np.linspace(0, 1, 50)
y = np.linspace(0, 1, 50)
dx, dy = x[1] - x[0], y[1] - y[0]
X, Y = np.meshgrid(x, y, indexing='ij')

Creating a BoundaryConditions Object

Pass the grid shape to the constructor:

# 1D
bc = BoundaryConditions((100,))

# 2D
bc = BoundaryConditions((50, 50))

Dirichlet Boundary Conditions

Dirichlet conditions fix the function value at the boundary. Assign a scalar or array directly:

# 1D: u(0) = 0, u(1) = 1
bc = BoundaryConditions((100,))
bc[0] = 0
bc[-1] = 1

# 2D: fix values on edges
bc = BoundaryConditions((50, 50))
bc[0, :] = 0          # u = 0 at x = 0
bc[-1, :] = 1         # u = 1 at x = 1
bc[:, 0] = Y[:, 0]    # u = y at y = 0
bc[:, -1] = Y[:, -1]  # u = y at y = 1

Boundary slicing uses standard NumPy indexing. The assigned value can be:

  • A scalar (applied to every point on that boundary)

  • A 1D array matching the boundary size

Neumann Boundary Conditions

Neumann conditions fix a derivative at the boundary. Pass a 2-tuple of (derivative_operator, value):

bc = BoundaryConditions((100,))
bc[0] = 0                         # Dirichlet: u(0) = 0
bc[-1] = Diff(0, dx), 0           # Neumann: u'(1) = 0

In 2D, use the normal derivative operator for each edge:

bc = BoundaryConditions((50, 50))
bc[0, :] = Diff(0, dx), 0       # du/dx = 0 at x = 0
bc[-1, :] = 1.0                  # Dirichlet at x = 1
bc[:, 0] = 300.0                 # Dirichlet at y = 0
bc[1:-1, -1] = Diff(1, dy), 0   # du/dy = 0 at y = 1

Note

When specifying Neumann conditions on a partial slice like bc[1:-1, -1], make sure the corners are covered by other conditions so the system is not under-determined.

Robin (Mixed) Boundary Conditions

Robin conditions combine Dirichlet and Neumann:

\[\alpha \, u + \beta \, \frac{\partial u}{\partial n} = g\]

Specify them as a 4-tuple (alpha, diff_op, beta, g):

bc = BoundaryConditions((100,))
bc[0] = 1.0                              # Dirichlet: u(0) = 1
bc[-1] = (1, Diff(0, dx), 1, 3)          # Robin: u + u' = 3 at x = 1

For 2D problems:

bc = BoundaryConditions((50, 50))
bc[-1, :] = (1, Diff(0, dx), 0.5, g_edge)   # alpha*u + beta*du/dx = g
bc[0, :] = u_left
bc[:, 0] = u_bottom
bc[:, -1] = u_top

Alternatively, construct the Robin operator yourself using Identity() and pass it as a Neumann-style 2-tuple:

from findiff import Identity

robin_op = alpha * Identity() + beta * Diff(0, dx)
bc[-1] = robin_op, g

Summary of Syntax

Type

Syntax

Example

Dirichlet

bc[slice] = value

bc[0] = 0

Neumann

bc[slice] = (Diff(...), value)

bc[-1] = Diff(0, dx), 0

Robin

bc[slice] = (alpha, Diff(...), beta, g)

bc[-1] = (1, Diff(0, dx), 1, 3)