### Maxwell Equations Boundary Conditions

Hello,

since my previous post was a little to stumpy I will retry to describe my problem.

I want to solve Maxwells equations for a scattering problem. The underlying geometry is a small metallic (perfect conductor) square with surface $S_2$`S`_{2} which is surrounded by a larger square with surface $S_1$`S`_{1} which acts as an absorbing boundary. The inner of the larger square minus the smaller square is the area $\Omega$Ω . I want to calculate the scattered field $\vec{E}^{sc}$→`E`^{sc} knowing the incoming field $\vec{E}^{in}$→`E`^{in} . The functional for this problem is given by: $F\left(\vec{E}^{sc}\right)=F_1\left(\vec{E}^{sc}\right)+F_2\left(\vec{E}^{sc}\right)$`F`(→`E`^{sc})=`F`_{1}(→`E`^{sc})+`F`_{2}(→`E`^{sc}) with $F_1\left(\vec{E}^{sc}\right)=\int_{\Omega}\nabla\times\vec{E}^{sc}\cdot\nabla\times\vec{E}^{sc}dA+k_0^2\int_{\Omega}\vec{E}^{sc}\cdot\vec{E}^{sc}dA$`F`_{1}(→`E`^{sc})=∫_{Ω}∇×→`E`^{sc}·∇×→`E`^{sc}`d``A`+`k`_{0}^{2}∫_{Ω}→`E`^{sc}·→`E`^{sc}`d``A` and $F_2\left(\vec{E}^{sc}\right)=i\cdot k_0\cdot\int_{S_1}\vec{n}\times\vec{E}^{sc}\cdot\vec{n}\times\vec{E}^{sc}ds$`F`_{2}(→`E`^{sc})=`i`·`k`_{0}·∫_{S1}→`n`×→`E`^{sc}·→`n`×→`E`^{sc}`d``s` . The first part accounts for the integration of the Helmholtz equation and the second parts accounts for the absorbing boundary condition. The solution can now be found by minimizing $F\left(\vec{E}^{sc}\right)$`F`(→`E`^{sc}) under the condition that $\vec{n}\times\vec{E}^{sc}=-\vec{n}\times\vec{E}^{in}$→`n`×→`E`^{sc}=−→`n`×→`E`^{in} at $S_2$`S`_{2}. My problem is to implement the boundary condition at $S_2$`S`_{2}.

Any suggestions or ideas are welcome :) Thanks for your help!

```
from fenics import *
from mshr import *
import numpy as np
import scipy.sparse as sps
import matplotlib.pyplot as plt
from linear_algebra import convert_matrix_dolfin_scipy, convert_vector_dolfin_numpy, solve_scipy_system
"""
Define natural constants
"""
c0 = 299792458
mu0 = 4*np.pi*1e-7
eps0 = 1/(c0**2*mu0)
"""
Define Problem parameters
"""
# physical properties
eps_r = 1
mu_r = 1
freq = 1e9
lam_0 = c0/freq
k_0 = 2*np.pi*lam_0
E_0 = 1.
r_0 = np.array([0, -10.])
k = k_0 * np.array([1., 0.])
"""
Define parameters for solvers and backends
"""
parameters['linear_algebra_backend'] = 'Eigen'
"""
Import the needed geometry
"""
mesh = Mesh('geometry/square/square.xml')
domains = MeshFunction('size_t', mesh, 'geometry/square/square_physical_region.xml')
boundaries = MeshFunction('size_t', mesh, 'geometry/square/square_facet_region.xml')
"""
Define incoming wave as matrix expression
E = E_0*p*exp(-i*k*(r-r_0)), p = cross(e_r, e_z)
with (((Re(E_x), Re(E_y)), ((Im(E_x), Im(E_y)))
"""
E_in = Expression(((('E_0 * (x[1] - r_0_y) * cos(- k_x * (x[0] - r_0_x) - k_y * (x[1] - r_0_y)) / sqrt( pow(x[0] - r_0_x, 2) + pow(x[1] - r_0_y, 2))'),\
('E_0 * (x[1] - r_0_y) * sin(- k_x * (x[0] - r_0_x) - k_y * (x[1] - r_0_y)) / sqrt( pow(x[0] - r_0_x, 2) + pow(x[1] - r_0_y, 2))')),\
(('-E_0 * (x[0] - r_0_x) * cos(- k_x * (x[0] - r_0_x) - k_y * (x[1] - r_0_y)) / sqrt( pow(x[0] - r_0_x, 2) + pow(x[1] - r_0_y, 2))'),\
('-E_0 * (x[0] - r_0_x) * sin(- k_x * (x[0] - r_0_x) - k_y * (x[1] - r_0_y)) / sqrt( pow(x[0] - r_0_x, 2) + pow(x[1] - r_0_y, 2))'))),\
E_0 = E_0, k_x = k[0], k_y = k[1], r_0_x = r_0[0], r_0_y = r_0[1], degree = 1, domain = mesh.ufl_domain())
E_in_abs = Expression((('E_0 * (x[1] - r_0_y) / sqrt( pow(x[0] - r_0_x, 2) + pow(x[1] - r_0_y, 2))'),\
('-E_0 * (x[0] - r_0_x) / sqrt( pow(x[0] - r_0_x, 2) + pow(x[1] - r_0_y, 2))')),\
E_0 = E_0, r_0_x = r_0[0], r_0_y = r_0[1], degree = 1, domain = mesh.ufl_domain())
"""
Define function space for test and trial function
"""
# Function space
func_space = FunctionSpace(mesh, 'N1curl', 1)
# Trial and Test function
trial_func = TrialFunction(func_space)
test_func = TestFunction(func_space)
# Measures
dx = Measure('dx', domain=mesh, subdomain_data=domains)
ds = Measure('ds', domain=mesh, subdomain_data=boundaries)
# Normal vector to facet
face_norm = FacetNormal(mesh)
"""
Define the bilinear forms, assemble and convert to scipy matrix
"""
bilform_curl = dot(curl(trial_func), curl(test_func))*dx
bilform_dot = dot(trial_func, test_func)*dx
bilform_face = (trial_func[0]*face_norm[1]-trial_func[1]*face_norm[0])*(test_func[0]*face_norm[1]-test_func[1]*face_norm[0])*ds(1)
# Assemble the matrices
mat_curl = assemble(bilform_curl)
mat_dot = assemble(bilform_dot)
mat_face = assemble(bilform_face)
# Convert to scipy sparse
mat_curl = convert_matrix_dolfin_scipy(mat_curl)
mat_dot = convert_matrix_dolfin_scipy(mat_dot)
mat_face = convert_matrix_dolfin_scipy(mat_face)
# Add up matrix parts
mat = 1/mu_r*mat_curl - k_0**2*eps_r*mat_dot + 1j*k_0*mat_face
```

Thank you for the answer!

Since I am solving a 2D problem there is no natural way to express the cross product. If I use the 'DirichletBC' it expects a 2-vector not a value for the z-component.

### 1 Answer

Mathematical note: cross (exterior) product in 2D makes sense, but the result doesn't stay in the same exterior algebra (it is 0-form, not 1-form)...

Play with

`DirichletBC's`

method `get_boundary_values()`

(it returns global degree of freedom index and value of that dof) and you'll see it does what I said.
Hello,

Am I asking the question wrong? Is there any information missing? I really need some help to proceed with the simulation.

Thanks - Christopher