Imposing boundary conditions for a system of equations

1

Hi everyone,

I am trying to impose the perfect conductor boundary conditions
n :heavy_multiplication_x: E = 0, n :black_small_square: H = 0, on the boundary,
on Maxwell’s equations
dE/dt - curl H = J,
dH/dt + curl E = 0 (parameters are set to 1 for simplicity)
using commands

bcD = [DirichletBC(W.sub(0), Constant((0.0, 0.0, 0.0)), DomainBoundary())]
bcB = [DirichletBC(W.sub(1), Constant((0.0, 0.0, 0.0)), DomainBoundary())]
bc = [bcD, bcB].

But I received these complaints from the compiler:

File “EM2.py”, line 58, in
solve(a == L, U, bc)
File “/usr/lib/python3/dist-packages/dolfin/fem/solving.py”, line 220, in solve
_solve_varproblem(*args, **kwargs)
File “/usr/lib/python3/dist-packages/dolfin/fem/solving.py”, line 235, in _solve_varproblem
= _extract_args(*args, **kwargs)
File “/usr/lib/python3/dist-packages/dolfin/fem/solving.py”, line 339, in _extract_args
bcs = _extract_bcs(args[2])
File “/usr/lib/python3/dist-packages/dolfin/fem/solving.py”, line 400, in _extract_bcs
raise RuntimeError(“solve variational problem. Unable to extract boundary condition arguments”)
RuntimeError: solve variational problem. Unable to extract boundary condition arguments
(FEniCS version 2019.1.0, Ubuntu 18.04.4 LTS)

I also want to add other variables to couple with the Maxwell’s equations, they will have their boundary conditions.

I have copied below a simple testing program. When I change the solve command to solve (a == L, U, bcD) it works. However I still want to specify both boundary conditions.

# Maxwell's Equations
# dE/dt - curl H = J
# dH/dt + curl E = 0

from dolfin import *
from ufl import nabla_div

# Create mesh
mesh = UnitCubeMesh(10, 10, 10)

# Build function spaces
D1 = FiniteElement("N1curl", mesh.ufl_cell(), 1)
B1 = FiniteElement("RT", mesh.ufl_cell(), 1)
DB = D1 * B1
W = FunctionSpace(mesh, DB)

# Boundaries
bcD = [DirichletBC(W.sub(0), Constant((0.0, 0.0, 0.0)), DomainBoundary())]
bcB = [DirichletBC(W.sub(1), Constant((0.0, 0.0, 0.0)), DomainBoundary())]
bc = [bcD, bcB]

# Initial conditions
U_0 = Expression(('sin(pi*x[1])*sin(pi*x[2])', 'sin(pi*x[2])*sin(pi*x[0])', '0', '-pi*sin(pi*x[0])*cos(pi*x[2])', \
    'pi*sin(pi*x[1])*cos(pi*x[2])', 'pi*sin(pi*x[2])*(cos(pi*x[0])-cos(pi*x[1]))'), degree = 2)
U_n = project(U_0, W)
E_n, H_n = split(U_n)

# Define variational problem
(E, H) = TrialFunctions(W)
(D, B) = TestFunctions(W)

E_ex = Expression(('sin(pi*x[1])*sin(pi*x[2])*exp(-t)', 'sin(pi*x[2])*sin(pi*x[0])*exp(-t)', '0'), degree = 2, t = 0)
H_ex = Expression(('-pi*sin(pi*x[0])*cos(pi*x[2])*exp(-t)', 'pi*sin(pi*x[1])*cos(pi*x[2])*exp(-t)', \
    'pi*sin(pi*x[2])*(cos(pi*x[0])-cos(pi*x[1]))*exp(-t)'), degree = 2, t = 0)
J = Expression(('-(2*pi*pi+1)*sin(pi*x[1])*sin(pi*x[2])*exp(-t)', \
    '-(2*pi*pi+1)*sin(pi*x[2])*sin(pi*x[0])*exp(-t)', '0'), degree = 2, t = 0)

dt = 0.1
a = inner(E, D)*dx + inner(H, B)*dx - dt*inner(H, curl(D))*dx + dt*inner(curl(E), B)*dx
L = dt*inner(J, D)*dx + inner(E_n, D)*dx + inner(H_n, B)*dx

# U = (E, H)
U = Function(W)
t = 0
for n in range(10):

    # Update current time
    E_ex.t = 0.1*n
    H_ex.t = 0.1*n
    J.t = 0.1*n

    # Solve variational problem
    solve(a == L, U, bc)

    # Update previous solution
    U_n.assign(U)

Thank you for considering my question. I am very grateful.

Yu

1 Like
2

The main issue here seems to be that you are wrapping the boundary conditions in multiple lists. It should be:

bcD = DirichletBC(W.sub(0), Constant((0.0, 0.0, 0.0)), DomainBoundary())
bcB = DirichletBC(W.sub(1), Constant((0.0, 0.0, 0.0)), DomainBoundary())
bc = [bcD, bcB]
1 Like
3

Thank you so much, dokken, you’ve offered a great help.

Yu

4

Hi everyone,

Sorry for replying to this old post. I would like to ask if anyone could give me some hints and/or references about why imposing perfect electric conductor (PEC) condition is equivalent, in fenics, to impose the field to be zero on the boundary, as done by @Andromath here:

Thank you very much.

5

By the nature of the Nedelec element enforcing homogeneous BCs on the boundary ensures \vec{n} \times \vec{u} = 0 . Cf. Peter Monk’s book.

1 Like