NEAR-FIELD IMAGING OF IMPEDANCE GRATING SURFACES
By
Yuqi Hong

A DISSERTATION
Submitted to
Michigan State University
in partial fulfillment of the requirements
for the degree of
Applied Mathematics – Doctor of Philosophy
2014

ABSTRACT
By
Yuqi Hong
This dissertation investigates the inverse rough surface scattering problem in near-field optical imaging, the purpose of which is to reconstruct a scattering surface with a resolution
beyond the diffraction limit. The surface is assumed to be a small and smooth deformation
of a plane surface. Based on a transformed field expansion, the boundary value problem with
complex scattering surface is converted into a successive sequence of a two-point boundary
value problems in the frequency domain, whereby an analytic solution for the direct scattering problem is derived from the method of integrated solution. By neglecting the high
order terms in the power series expansion, the nonlinear inverse problem is linearized and
an explicit inversion formula is obtained. A spectral cut-off regularization is adopted to suppress the exponential growth of the noise in the evanescent wave components, which carry
high spatial frequency of the scattering surface and contribute to the super resolution in the
near-field regime. The method works for sound soft, sound hard, and impedance surfaces,
and requires only a single illumination at a fixed frequency and is realized efficiently by
the fast Fourier transform. Numerical results show that the method is simple, stable, and
effective to reconstruct scattering surfaces with sub-wavelength resolution.

To Dad and Mom

iii

ACKNOWLEDGMENTS
First and foremost, I want to thank my parents, for loving me and taking care of me from the
start of my life, for always supporting me on my dream path, and for being a wonderful family.
If I look back through twenty-seven years of my life, my parents have always been there,
giving me the best of everything and inspiring me. Without their love and encouragement,
I would have never pursued a PhD in mathematics, let alone finished it. Without their
guidance and patience, my academic path would have never have gone so far and smooth.
My greatest academic thanks goes to Professor Gang Bao, who has my thesis advisor for
four years, for his guidance and support, for his kindness and patience during my graduate
study at Michigan State University. I would like to thank him for giving me this interesting
thesis topic and encouraging me all the time. In these four years, from the qualifying exam
to dissertation defense, from thesis topic to job suggestions, from research field to trivial
paperworks, he has spent so much time and effort on me that I would have never reached
this point without him. In addition, his generous support of research assistantship gave me
the boost to go further on my academic path without distraction.
The man with the second greatest impact on my research is Professor Peijun Li. His guidance
and tutelage energized my research and set me on a path to achievement. I also learned the
subtleties of research from working with him. I would like to express my gratitude to my
comprehensive exam and defense committee professors, Gang Bao, Yingda Cheng, Di Liu,
Jianliang Qian, Moxun Tang, Baisheng Yan, Zhengfang Zhou, for devoting time and concern
to me. I also thank the Department of Mathematics for giving me so much valuable teaching
experience and so many professors with whom to interact.
I also want to thank Justin Droba, for taking care of me and loving me, for teaching me
physics and coding, for always being there for me and being a sweet husband, and for
everything he has done for me.
iv

I also want to thank Dennis Droba and Diane Droba, for being a sweet father-in-law and
mother-in-law, for taking care of me, and for so many times of picking me up from Detroit
airport and driving me the hour and a half back to Lansing and cooking for me.
My journey through graduate school would not have been nearly as enjoyable without the
good company that continuously surrounded me, sharing in the good moments and providing
advice and emotional support. I would especially like to acknowledge my friends, Xianfeng
Hu, Jun Lai, Yuliang Wang, Zixuan Wang, Qi Tang, Eric Wolf, Hongli Gao, all of whom
were always there when I needed them.
Finally, I acknowledge the alumni students and postdocs of Prof. Bao’s group not already
mentioned: Guanghui Hu, Xianliang Hu, Xiang Xu, Junshan Lin, Zhengfu Xu, Xinghui
Zhong. I hope our paths will cross again one day.

v

TABLE OF CONTENTS
LIST OF TABLES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii
LIST OF FIGURES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . viii
Chapter 1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1

Chapter 2 Background Knowledge
2.1 Maxwell’s Theory . . . . . . . .
2.2 Notations . . . . . . . . . . . . .
2.3 Wave Propogation . . . . . . . .
2.4 The Vector Helmholtz Equation
2.5 Fundamental Polarizations . . .
2.6 Impedance Boundary Condition
2.7 Analysis . . . . . . . . . . . . . .

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

6
6
9
10
13
16
18
19

Chapter 3 Forward Problem . .
3.1 Model problem . . . . . . . .
3.2 Transformed field expansion
3.3 Fourier series expansion . . .

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

22
22
34
38

.
.
.
.

.
.
.
.

Chapter 4 Inverse Problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
4.1 Reconstruction formula . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
Chapter 5 Numerical Experiments . . . .
5.1 Numerical Example 1: Low Frequency
5.2 Numerical Example 2: High Frequency
5.3 Summary and Conclusion . . . . . . . .

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

51
52
54
55

APPENDIX . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
REFERENCES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64

vi

LIST OF TABLES
Here
Table 2.1

Quantities of Maxwell’s Equations . . . . . . . . . . . . . . . . . . . . .

7

Table 2.2

Fundamental Electromagnetic Constants . . . . . . . . . . . . . . . . .

8

vii

LIST OF FIGURES
Here
Figure 1.1

Placement of Boundary Conditions . . . . . . . . . . . . . . . . . . . . .

4

Figure 3.1

Problem Geometry . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

Figure 5.1

Numerical Experiment #1: ε = 0.01 . . . . . . . . . . . . . . . . . . . . . 52

Figure 5.2

Numerical Experiment #1: h = 0.1λ . . . . . . . . . . . . . . . . . . . . 53

Figure 5.3

Numerical Experiment #2: ε = 0.01 . . . . . . . . . . . . . . . . . . . . . 54

Figure 5.4

Numerical Experiment #2: h = 0.1λ . . . . . . . . . . . . . . . . . . . . 55

viii

Chapter

1

Introduction
The aim of this dissertation is to present an effective mathematical model and accurate
computational algorithm for the grating surface problem that arises in near-field imaging with
general impedance boundary condition. The impedance boundary condition has become an
important concept in the modeling of wave propagation. It can be used for obstacles that are
partially or totally penetrable for time-harmonic scattering problems in electromagnetism.
The famous sound hard and sound soft boundary conditions are special cases of the general
impedance boundary condition.
A grating surface is a periodic structure with a non-local perturbation such that the new
surface lies within a small distance of the original plane. The perturbation is assumed to
be smooth and very small. The input information is measured at a constant distance from
above the grating surface orders in magnitude smaller than the wavelength of the incident
radiation. This is the defining feature of near-field optics, which focuses on features of light
or waves occuring in regions smaller than the wavelength. The hallmark of near-field optics is
that it breaks the diffraction limit and is able to obtain images at subwavelength resolution.
It also examines evanescent wave modes, which are crucial in near-field optical microscopes.
Near-field optics also has broad application in biology, chemistry, and materials science.
This dissertation limits discussion to time-harmonic waves field incident normally incident
1

upon on a grating surface from above. A homogeneous medium fills the space above the
grating surface. The method we will present also works for general incidence, but the computation is so much more complicated that it obscures the elegance of the method. Consequently, we limit presentation here to normal incidence only. The grating surface problem is
just one kind of scattering surface problem. Scattering is a physical process in which some
forms of radiation, such as light, sound, or moving particles, are redirected from their original trajectories by localized non-uniformities in the medium through which they are passing.
For example, the color of the sky is formed by dispersive scattering of sunlight by particles and molecules in the atmosphere, with stronger scattering at short (blue) wavelengths
leading to the sky’s distinctive hue. The scattering surface problem is so important that
it has been applied broadly across modern science and technology. For instance, detection
of the interior of the earth is possible by measuring the travel time of waves produced by
seismographs. The technique is invaluable to oil companies, who need to determine the best
location to drill. Beyond applications in the energy industry, medical imaging uses X-rays,
ultrasound waves, and other electromagnetic waves to help diagnose suspected tumors and
other malignant growths. Scattering problems also model animal echolocation, or the bio
sonar of bats, dolphins, and toothed whales. This idea can also be used to help blind humans
find their way around.
Paired with every inverse problems is a forward problem. For scattering, the forward scattering problem is to determine how radiation or particles behave with preknowledge of the
scatterer. That is, the goal to find the total field given complete information about the
grating surface and incident wave. Conversely, inverse scattering seeks to determine the
characteristics of an object, such as shape or internal structure, solely from data concerning
the scattering of incoming radiation, waves, or particles. Mathematically, the inverse scattering problem is to reconstruct the function that defines the scattering surface given the
total field and incident wave.

2

The scattering surface is called an unbounded or infinite rough surface when the perturbation is non-local and not confined to a compact region. The problem has been studied and discussed by a large number of researchers, including Voronovich[40], Chew and
Warnick[41], also Chandler-Wilde and Monk[18], DeSanto and Martin[24], Li, Wu and
Zheng[32], Milder[36], DeSanto and Wombell[42],Chandler-Wilde and Zhang[19], Bao and
Lin[14]. Bao and Lin[14] considered an inverse scattering problem with a local perturbation
on a perfect electric conductor (PEC) is considered. Since the perturbation is local, the
Sommerfeld radiation condition can be imposed in their model. However, a subsequent work
by Bao and Li[10] examined inverse scattering upon a PEC perturbed nonlocally. In this
scenario, the usual Sommerfeld radiation condition is inapplicable. So instead, the transparent boundary condition is imposed. Two related papers examine matters arising in nonlocal
PEC scatters: Bao and Li[9] analyzed the convergence and error estimate; and Li and Shen
[30]proved well-posedness and uniqueness of the weak solution of the corresponding forward
problem using a variational approach. Furthermore, the near-field imaging of infinite rough
surface scattering problem in dielectric media was tackled by Bao and Li[11]. Other related
works in near-field optics include Carney and Schotland in [17] and [16], Courjon[22], Bainier
and Courjon[23].
A cavity wall is a local perturbation of a planar surface in which the perturbation lies below
the unperturbed portion of the plane (cf. [26]). Many beautiful results in different offshoots
of the cavity problems have been established. Among them are Li[28], Bao, Gao and Li[29],
Li and Wood[31], Li,Wu and Zheng[33], Bao, Gao, Lin and Zhang[8], Ammari, Bao, and
Wood[2].
If the scattering surface is periodic structure, it is called a grating surface. Bao, Cowsar,
and Masters[5] and Bao, Dobson and Cox[6], Petit[39], Bao and Friedman[7] have established
excellent studies on mathematical analysis and numerical approaches for diffraction gratings.
Kirsch[27], Bao[3] and Ammari[1] published uniqueness theorems for grating problems. Chen

3

and Wu[20] and Bao, Chen and Wu[4] crafted computational strategies based on the finite
element method (FEM). Bao, Li and Lv [12] and Bao, Li and Wu[13] published an elegant
work about inverse diffraction grating problem. The inverse grating surface problem with
PEC boundary condition has recently been studied in work of Cheng, Li and Wang[21]. This
thesis focuses on grating problem with general impedance boundary condition, which is a
more accurate approximation to the exact electromagnetic properties of the material.
Excepting the introduction, this dissertation is comprised of four chapters. The first introduces mathematics and physics background, and the second establishes the mathematical
model and develops a solver for the forward problem. Chapter 4 derives an inversion formula for reconstruction of the grating surface function. The final chapter presents numerical
experiments that highlight the accuracy and efficiency of Chapter 3’s method.
As will be shown in Section 2.4, the mathematical model of the scattering problem is a
two-dimensional Helmholtz equation with two different types of boundary conditions:
1. A transparent boundary condition on the inverse problem’s measurement surface
2. A universally applicable impedence boundary condition on the grating surface
The placement of these boundary conditions is depicted visually below in Figure 1.1.

Γ

S

Figure 1.1 Placement of Boundary Conditions. On the measurement surface Γ,
shown in red, a transparent boundary condition is placed; on the grating surface S (in
blue), an impedence condition.

4

The neccessity of the transparent boundary condition arises from the invalidity of the usual
Sommerfeld radiation condition because of the nonlocal perturbation under consideration.
To deal with the boundary value problem on a curved domain, a transformed field expansion
is applied and the domain is transformed into a rectangular region. A formal power series
expansion converts a boundary value problem concerning a complex scattering surface into
a much simpler progression of two-point boundary value problems. Works by Nicholls and
Reitich [38, 37], Malcolm and Nicholls[35], Bruno and Reitich[15] are the primary references
for the transformed field expansion and boundary perturbation method. The boundary
perturbation method and boundary integral operators are applied to an inverse problem with
periodic scattering surface in Malcolm and Nicholls[34]. An explicit solution of the forward
problem is obtained by the integration solution method. The uniqueness and wellposedness
of the solution of this problem is proven.
Chapter 4 of this thesis derives an inversion formula for the grating surface. Discarding
the high order terms in the power series creates a linearization of the nonlinear inverse
problem. By calculating just the first two terms in the forward problem, an elegant yet
accurate expression for grating surface is achieved. The exponential growth in the inversion
formula is cured by a spectral cut-off regularization that reflects the ill-posedness of the
inverse problem.
Chapter 5 of this thesis presents numerical experiments that demonstrate the accurate stability and efficiency of this method for reconstructing the grating surface with subwavelength
resolution. When the deformation parameter is fixed, higher accuracy can be achieved by
taking measurements closer to the surface. This reinforces the purpose of near-field imaging:
measurements taken too far away lead to amplification of the exponential in the reconstruction formula and loss of valuable information contained in evanescent wave mode.

5

Chapter

2

Background Knowledge
2.1

Maxwell’s Theory

Maxwell’s equations are precise mathematical expressions of the laws of nature. There are
many different ways using different notations to display the same mathematical formulation,
just like using different language to express a same idea. Here we express Maxwell equation
in the following mathematical form:
∂D
+ J,
∂t
∂B
,
∇×E=−
∂t

∇×H=

(2.1)
(2.2)

∇ ⋅ D = ρ,

(2.3)

∇ ⋅ B = 0,

(2.4)

where E, B, H, D, J and ρ are all functions of space x and t.
Equation (2.1) is Ampere’s law or the generalized Ampere circuit law. Equation (2.2) is
Faraday’s law or Faraday’s magnetic induction law. Equation (2.3) is Coulomb’s law or
Gauss’ law for electric fields. Equation (2.4) is Gauss’ law or Gauss’ law for magnetic fields.
We refer to E as electric fields, D as electric displacement, and H and B as magnetic fields.

6

Symbol

Meaning

Units

H(x, t)

Magnetic field

amperes/m

E(x, t)

Electric field

volts/m

B(x, t)

Magnetic flux density

webers/m2

D(x, t)

Electric displacement

coulombs/m2

ρ(x, t)

Electric charge density

coulombs/m2

J(x, t)

Electric current density

amperes/m2

Table 2.1 Quantities of Maxwell’s Equations. This table lists the physical quantities
appearing in (2.1)–(2.4) with their meanings and their units.
Ampere’s original formulation of (2.1) lacked the electric displacement contribution

∂D
∂t ,

a

term added by Maxwell in his pioneering 1861 work “On Physical Lines of Force.” While
leading to correct results within magnetostatics, Ampere’s law without Maxwell’s correction
fails when time-varying charges are present in several ways:
1. Capacitors, which result in open circuits for direct current sources, can only be employed in alternating circuit systems. When the original Ampere’s Law is applied
to such a system, the result is a discontinuous current profile. Maxwell’s correction
restores the proper behavior.
2. Taking the divergence of both sides of (2.1), using that ∇ ⋅ (∇ × A) = 0, and calling
upon (2.3) results in the continuity equation

∇ ⋅ J(x, t) = −

∂
ρ(x, t)
∂t

(2.5)

which states that if current flows out of a region (i.e., ∇ ⋅ J > 0), then the charge
density in that region must correspondingly decrease. This amounts to a conservation
of charge. In Ampere’s original formulation,

∇ ⋅ J(x, t) ≡ 0
7

which implies that the outward flux in space of current always vanishes. Obviously,
this is correct only if charge does not vary in time—e.g., in direct current systems.
3. Faraday’s Law (2.2) states that a time-varying magnetic field induces a time-varying
electric field. With the displacement term in (2.1), the corrected Ampere’s Law, in conjunction with a constitutive relation between D and E, reciprocates this phenomenon:
time-varying electric fields induce time-varying magnetic fields. This relationship between the time-varying electric and magnetic fields constitutes the foundation of electromagnetic wave theory and led Maxwell to prediction of electromagnetic waves.
In linear materials lacking magnetic dipole moments, D and E are related by a dielectric
permittivity ε, and the magnetic fields B and H are related by a magnetic permeability µ:

D = εE

(2.6)

B = µH

(2.7)

We call (2.6) and (2.7) the constitutive relations. In free space (a vacuum), µ = µ0 and ε = ε0 .
The numerical values of these fundamental constants are given below in Table 2.2.

Quantity

Value

Units

ε0

8.85 × 10−12

farad/meter

µ0

4π × 10−7

henry/meter

Table 2.2 Fundamental Electromagnetic Constants.

Inside a material medium, the permittivity ε is determined by the electrical properties of
the medium and the permeability µ is decided by the magnetic properties of the medium.
In particular, for nonmagnetic materials, µ = µ0 , but there is no physical material for which
ε = ε0 . Both ε and µ are more properly rank-2 tensors (3x3 matrices) with components that

8

zero out according to the material’s crystal symmetry. For example, in uniaxial crystals,
⎛ ε (x)
⎜ xx
⎜
ε(x) = ⎜
0
⎜
⎜
⎜
0
⎝

0
εyy (x)
0

⎞
⎟
⎟
⎟
0
⎟
⎟
⎟
εzz (x) ⎠
0

Only in homogeneous and isotropic materials will ε remain a scalar constant and in such
materials the electric field E is parallel to D and the magnetic fields B and H are parallel† .
In just isotropic materials, ε will be a constant tensor. Furthermore, because all materials
display dispersion, ε also depends on frequency ω.
In this thesis, however, discussion will be limited to homogeneous and isotropic materials so
that ε ≡ ε. Because all incident waves are assumed to be monochromatic, we may ignore the
effects of dispersion and set ε to its value appropriate for the input frequency.

2.2

Notations

A vector A has a magnitude and a direction, which can be represented graphically by
a straight-line element of length proportional to the magnitude of A and with an arrow
pointing the direction of A. In a Cartesian coordinate system, which is also the rectangular
coordinate system, A can be written in terms of the three Cartesian components Ax , Ay , Az :
A = xˆAx + yˆAy + zˆAz

(2.8)

where Ax , Ay , Az are the projections of A onto the x, y, z axes. Denote the unit vector in
the direction of x, y, z axes as xˆ, yˆ, zˆ. The del operator ∇ is a vector differential operator,
†

More strictly, H is the magnetic field and B is the magnetic flux density. In nonmagnetic materials,
however, the two are proportional so the distinction is not important to us. We refer to both as “magnetic
field” as is common in the literature.

9

written as

∇ = xˆ

∂
∂
∂
+ yˆ + zˆ
∂x
∂y
∂z

(2.9)

When operating on a scalar function φ(x, y, z), the result is a vector

∇φ = xˆ

∂φ
∂φ
∂φ
φ + yˆ φ + zˆ
∂x
∂y
∂z

(2.10)

Since xˆ ⋅ xˆ = 1, yˆ ⋅ yˆ = 1, zˆ ⋅ zˆ = 1,we have
∇2 = ∇ ⋅ ∇ =

∂2
∂2
∂2
+
+
∂x2 ∂y 2 ∂z 2

(2.11)

The curl of a vector field H is a vector defined as

∇ × H = (ˆ
x

∂
∂
∂
+ yˆ + zˆ ) × H =
∂x
∂y
∂z

xˆ

yˆ

zˆ

∂
∂x

∂
∂y

∂
∂z

Hx Hy Hz
= xˆ (

∂Hy ∂Hx
∂Hx ∂Hz
∂Hz ∂Hy
−
) + yˆ (
−
) + zˆ (
−
)
∂y
∂z
∂z
∂x
∂x
∂y

The divergence of a vector function is a scalar, defined as
∂
∂
∂
+ yˆ + zˆ ) ⋅ (ˆ
xDx + yˆDy + zˆDz )
∂x
∂y
∂z
∂Dx ∂Dy ∂Dz
=
+
+
∂x
∂y
∂z

∇ ⋅ D = (ˆ
x

2.3

Wave Propogation

In the physics of wave propagation, a plane wave (also spelled planewave) is a constantfrequency wave whose wavefronts (surfaces of constant phase) are infinite parallel planes of

10

constant peak-to-peak amplitude normal to the phase velocity vector.
The wavenumber (also wave number) is the spatial frequency of a wave, either in cycles per
unit distance or radians per unit distance. It can be envisioned as the number of waves that
exist over a specified distance—analogous to frequency being the number of cycles or radians
per unit time. In terms of angular frequency, the wavenumber κ is given by

κ=

2π ω
=
λ
v

(2.12)

where λ is the wavelength, ω is angular frequency of the wave.
Theorem 2.1 (Dispersion Relation). The wavenumber κ is given by

κ2 = µ0 ε0 ω 2 = ω 2 /c2
where µ0 ε0 =

1
c2

and c = 3 × 108 m/sec is the speed of light in free space.

The dispersion relation provides an important connection between the spatial frequency κ
and the temporal frequency ω.
Definition 2.1. Fields for which the time variation is sinusoidal or periodic are called
time-harmonic fields. Time-harmonic waves function can be written in a form of A(x, t) =
A(x)eiωt , where ω is the wave frequency.
Plane wave is one kind of time-harmonic fields and is called a “monochromatic” wave because it consists of a single frequency (color) of radiation. An incident wave may also be
“polychromatic,” or composed of a finite number of monochromatic waves. More generally,
any field u(x) may be decomposed into an infinite combination of monochromatic waves by
means of the Fourier transform:
uˆ(ξ) = ∫ u(t)e−2πiξt dt
R

11

By evaluating the above expression at a particular frequency ω = 2πξ, we can extract the
component of u at that frequency.
Theorem 2.2 (Gauss’s Law). The area integral of the electric field over any closed surface
is equal to the net charge enclosed in the surface divided by the permittivity of space.

∮S E dx =

Qenclosed
ε0

In other words, the total of the electric flux out of a closed surface is equal to the charge
enclosed divided by the permittivity. It can be written in the differential form

∇⋅E=

ρ
ε0

The electric flux through an area is defined as the electric field multiplied by the area of
the surface projected in a plane perpendicular to the field. Gauss’s Law is a general law
applying to any closed surface. It is an important tool since it permits the assessment of the
amount of enclosed charge by mapping the field on a surface outside the charge distribution.
For geometries of sufficient symmetry, it simplifies the calculation of the electric field.
Remark 2.1. First taking the divergence and then integrating both sides of the constitutive
relation (2.6) over a surface region S, we have in free space (ε = ε0 )
∮S ∇ ⋅ D dx = ε0 ∮S ∇ ⋅ E dx
∮∂S D ⋅ n dx = Qenclosed
with the equality on the right side coming from the Divergence Theorem; the left side, from
Gauss’ Law. Therefore, Gauss’ Law gives meaning to the electric displacement: D accounts
for the effects of charge contained within materials.

12

∎

2.4

The Vector Helmholtz Equation

We are interested in Maxwell equations in regions devoid of source—that is in regions where
J = 0 and ρ = 0. This doesn’t mean there is no source anywhere in all space, instead, source
exist outside the regions of interest in order to produce fields in these regions. Employing
the constitutive relations (2.6) and (2.7) and this source-free assumption in (2.1)–(2.4), the
Maxwell Equations become
∂
E
∂t
∂
∇ × E = −µ0 H
∂t

∇ × H = ε0

(2.13)
(2.14)

∇⋅E=0

(2.15)

∇⋅H=0

(2.16)

In terms of scalar partial differential equations, Ampere’s Law (2.13) takes the form
∂
∂
∂
Hz − Hy = ε0 Ex
∂y
∂z
∂t
∂
∂
∂
Hx −
Hz = ε0 Ey
∂z
∂x
∂t
∂
∂
∂
Hy − Hx = ε0 Ez
∂x
∂y
∂t

(2.17)
(2.18)
(2.19)

Faraday’s Law (2.14) becomes
∂
∂
∂
Ez − Ey = −µ0 Hx
∂y
∂z
∂t
∂
∂
∂
Ex −
Ez = −µ0 Hy
∂z
∂x
∂t
∂
∂
∂
Ey − Ex = −µ0 Hz
∂x
∂y
∂t

13

(2.20)
(2.21)
(2.22)

Finally, the two the Gauss’ Laws for electric, (2.15), and magnetic, (2.16), fields are
∂
∂
∂
Ex + Ey + Ez = 0
∂x
∂y
∂z
∂
∂
∂
Hx + Hy + Hz = 0
∂x
∂y
∂z

(2.23)
(2.24)

Taking the derivative with respect to time of (2.17) and plugging in (2.21) (2.22), we get

µ 0 ε0

∂2
∂ ∂
∂
∂ ∂
∂
E
=
−
(
E
−
E
)
+
(
E
−
Ez )
x
y
x
x
∂t2
∂y ∂x
∂y
∂z ∂z
∂x
∂2
∂2
∂2
= ( 2 + 2 + 2 ) Ex
∂x
∂y
∂z

(2.25)
(2.26)

Finally, we plug in (2.23). Similarly, we can derive the following equations:
∂2
∂2
∂2
∂2
+
+
−
µ
ε
) Ex = 0
0 0
∂x2 ∂y 2 ∂z 2
∂t2
∂2
∂2
∂2
∂2
( 2 + 2 + 2 − µ0 ε0 2 ) Ey = 0
∂x
∂y
∂z
∂t
2
2
2
∂
∂
∂
∂2
( 2 + 2 + 2 − µ0 ε0 2 ) Ez = 0
∂x
∂y
∂z
∂t

(

(2.27)
(2.28)
(2.29)

which can be collected into the vector expression

∇2 E − µ0 ε0

∂2
E=0
∂t2

(2.30)

This is the vector Helmholtz wave equation. If the wave is time-harmonic, then by means of
the Fourier transform,

∂2
∂t2 E(x, t)

= (iω)2 E(ω, t), the Helmholtz equation can be written in

frequency space†

∇2 E + µ0 ε0 ω 2 E = 0
†

(2.31)

We implicitly always work in frequency space (after transforming in time only), so we do not make any
notational distinctions between the field and its transformed version aside from the change of symbol from
t to ω.

14

Use the dispersion relation, κ2 = ω 2 µ0 ε0 , we have
(∇2 + κ2 )E = 0

(2.32)

There is also another way to derive Helmholtz equation from Maxwell equation in free space.
For all time-harmonic wave in free space,
∂B
= −iωB
∂t
∂D
∇×H=
= iωD = iωεE
∂t
∇×E=−

which implies

−

We multiply by

1
µ

1
∇×E=B
iω

on both sides,

−

1
1
∇×E= B
iωµ
µ

and then take the curl. Using (2.4) and noting that

−

B
µ

= H, we have

1
∇ × (∇ × E) = ∇ × H = iωεE
iωµ
⇒ ∇ × (∇ × E) − ω 2 εµE = 0

Using the vector identity

∇ × ∇ × A = ∇(∇ ⋅ A) − ∇2 A

15

and calling upon the Gauss Law ∇ ⋅ E = 0 for charge-free spaces, we have
∇2 E + ω 2 εµE = 0

(2.33)

Due to the dispersion relation κ2 = ω 2 µ0 ε0 , the equation becomes
∇2 E + κ2 E = 0

(2.34)

This is precisely what we had before.

2.5

Fundamental Polarizations

Diffractive optics is an emerging technology with many applications. Some of the important
applications include the design and fabrication of optical elements such as corrective lenses,
anti-reflective interfaces, beam splitters, and sensors. A particularly important application
is replacing conventional lenses by diffractive gratings which are designed and fabricated by
interference fringes on holographic plates.
The basic electromagnetic theory of gratings has been studied extensively since Rayleighs
time in the early 20th century. Recent advance has been greatly accelerated due to several
new approaches and numerical methods including differential methods, integral methods,
analytical continuation, variational method, and others.
Diffractive optical elements, as opposed to the traditional optical lenses, have many advantages. They are light, small, and inexpensive. Often diffractive structures exhibit certain
periodicity. There are two classes of grating structures: linear gratings (one-dimensional
gratings) and crossed gratings (biperiodic or two-dimensional gratings).
The scalar approach is valid only as long as the wavelength of the incident light is small
compared with characteristic transverse dimension of the structure. However, in many ap16

plications the structure dimensions are of the order of the wavelength of the incident light.
In that case, the vectorial electromagnetic character of light becomes of interest. We are
then speaking of the vector theory. In the vector model, polarization properties specified by
the electric and magnetic fields in accordance with Maxwells equations must be studied.
Suppose that a grating is illuminated by a plane wave which lies in the (x, y) plane. The
electromagnetic fields are assumed to be independent of z. We consider the following two
fundamental cases of polarizations: transverse electric (TE) polarization and transverse magnetic (TM) polarization.
In TE polarization, the electric field is parallel to the grooves or is pointed to the z direction,
i.e.,
⎛
0
⎜
⎜
E=⎜
0
⎜
⎜
⎜
⎝ u(x, y)

⎞
⎟
⎟
⎟
⎟
⎟
⎟
⎠

(2.35)

It follows from the Maxwell equation for the electric field

∇ × (∇ × E) − κ2 E = 0

that we may derive the two-dimensional scalar Helmholtz equation

∆u + κ2 u = 0.

In TM polarization, the magnetic field is parallel to the grooves or is pointed to the z
direction, i.e.,
⎛
0
⎜
⎜
H=⎜
0
⎜
⎜
⎜
⎝ u(x, y)

17

⎞
⎟
⎟
⎟
⎟
⎟
⎟
⎠

It follows from the Maxwell equation for the magnetic field
∇ × (κ−2 ∇ × H) − H = 0

that we may derive
∇ ⋅ (κ−2 ∇u) + u = 0.
Because κ is independent of x, this is the same equation as for TE polarization.

2.6

Impedance Boundary Condition

The impedance boundary condition is one which relates the tangential components of the
electric and magnetic fields via an impedance factor, which is a function of the properties
of the surface and maybe the incident field upon it. The approximate boundary conditions
can be used to simulate the material properties of a surface, and is very useful in simplifying
the analytical and numerical solution of scattering problems. A simple kind of approximate
boundary condition is impedance boundary condition, which is given by

n × (∇ × E) −iη n × (n × E) = 0,
normal comp.
of ∇×E∝H

tangential comp.
of E

where n is the unit outward normal on the surface and where η is a function of the material
properties called impedance parameter.
The unit outward normal is n = (n1 , n2 , 0) where n21 + n22 = 1. Using the form of E in the TE
polarization (2.35), a simple calculation yields
⎛
0
⎜
⎜
n × (∇ × E) = ⎜
0
⎜
⎜
⎜
∂u
∂u
⎝ −n1 ∂x − n2 ∂y
18

⎞ ⎛ 0
⎟ ⎜
⎟ ⎜
⎟=⎜ 0
⎟ ⎜
⎟ ⎜
⎟ ⎜
⎠ ⎝ −∂n u

⎞
⎟
⎟
⎟
⎟
⎟
⎟
⎠

where
∂n u = n ⋅ ∇u
Next,
⎛ n ⎞ ⎛ nu
⎜ 1 ⎟ ⎜ 2
⎜
⎟ ⎜
⎟ ⎜
n × (n × E) = ⎜
⎜ n2 ⎟ × ⎜ −n1 u
⎜
⎟ ⎜
⎜
⎟ ⎜
⎝ 0 ⎠ ⎝ 0

⎞ ⎛ 0 ⎞
⎟ ⎜
⎟
⎟ ⎜
⎟
⎟=⎜ 0 ⎟
⎟ ⎜
⎟
⎟ ⎜
⎟
⎟ ⎜
⎟
⎠ ⎝ −u ⎠

Combining above identities, we obtain the impedance boundary condition for the TE polarization
∂n u − iηu = 0.

2.7

Analysis

Denote α = (α1 , α2 , ..., αN )T ∈ ZN
+ , where Z+ is the set of non-negative integers, we set
N

∣α∣1 = ∑ ∣αi ∣ and for u ∈ C ∣α∣1 (Ω), define
i=1

∂ ∣α∣1 u
∂ αu
=
∂xα ∂x1 α1 , ..., ∂xN αN
Definition 2.2. An open, connected set Ω ⊂ RN ,N = 1, 2, 3 will be referred to as a domain.
The fundamental Sobolev spaces are denoted W s,p (Ω), where s ∈ Z+ , 1 ≤ p < ∞ and Ω is a
domain in RN . These spaces are defined by
W s,p (Ω) = {u ∈ Lp (Ω) ∣ ∂ α u ∈ Lp (Ω) for all ∣α∣1 ≤ s}

Each of these spaces is equipped with the norm
1

⎛
⎞p
∣∣ u ∣∣W s,p (Ω) = ∑ ∫ ∣∂ α u(x)∣p dV (x)
⎝∣α∣1 ≤s Ω
⎠

19

(2.36)

An important case occurs when p = 2. Denote the Sobolev Spaces
H s (Ω) = {u ∶ Dα u ∈ L2 (Ω) for all ∣α∣ ≤ s}

and the periodic functional space

Hp1 (Ω) = {u ∈ H 1 (Ω) ∶ u(0, y) = u(2π, y)},
Hp1 (Ω) is a subspace of H 1 (Ω). To define spaces on the boundary, we first need Parseval’s
Theorem:
Theorem 2.3 (Parseval). Let f, g ∈ L2 [0, 2π] with f (0) = f (2π) and g(0) = g(2π) (i.e., f
and g are 2π-periodic) with Fourier series
f (x) = ∑ f (n) einx
n∈Z

g(x) = ∑ g (n) einx
n∈Z

Then we have
2π

∫0

f (x)g(x) dx = ∑ f (n) g (n)
n∈Z

This theorem allows us to define the trace functional space for periodic functions u defined
on boundary Γ with Fourier coefficient u(n) :
2

H s (Γ) = {u ∈ L2 (Γ) ∶ ∑ (1 + n2 )s ∣u(n) ∣ < ∞},
n∈Z

which is a Banach space equipped with the norm
1/2
2

∣∣ u ∣∣s,Γ = [2π ∑ (1 + n2 )s ∣u(n) ∣ ]
n∈Z

The dual space associated with H s (Γ) is the space H −s (Γ) with respect to the scalar product
20

in L2 (Γ) defined by
⟨u, v⟩Γ = ∫ uvdx = 2π ∑ u(n) v (n)
Γ

n∈Z

The equivalence of Fourier series coefficients and the standard L2 inner product is also
established by Parseval’s Theorem.
Theorem 2.4 (Trace Theorem). Assume Ω is bounded Lipschitz domain and

1
p

< s ≤ 1,

¯ has a unique continuous extension as a linear
the mapping γT (u) = u∣∂Ω defined on C ∞ (Ω)
operator from W s,p (Ω) onto W s−1/p,p (∂Ω). Moreover,
W01,p (Ω) = {u ∈ W 1,p (Ω) ∣ γT (u) = 0}.

Theorem 2.5 (Lax-Milgram). Let H be a Hilbert space equipped with inner product ⟨ ⋅ , ⋅ ⟩H
and norm ∣∣ ⋅ ∣∣H . Assume that B ∶ H × H → R is a bilinear mapping, for which there exist
constants α, β > 0 such that
∣B[u, v]∣ ≤ α ∣∣ u ∣∣H ∣∣ v ∣∣H

(2.37)

for all u, v ∈ H and
2

β ∣∣ u ∣∣H ≤ B[u, u]

(2.38)

for all u ∈ H. Let f ∶ H → R be a bounded linear functional on H. Then there exists a unique
element u ∈ H such that
B[u, v] = f (v)
for all v ∈ H.
Remark 2.2. Condition (2.37) says that the bilinear form is bounded. Condition (2.38)
∎

says that it is coercive.

21

Chapter

3

Forward Problem
3.1

Model problem

In this chapter, we shall introduce a mathematical model, define some notations for the
grating problem by a surface, and present a boundary value problem for the surface grating
model. Let us first describe the problem geometry, depicted visually below in Figure 3.1.

Γ

Ω
y

S
x

Λ
Figure 3.1 Problem Geometry. The grating surface S (in red) is assumed to be
periodic with period Λ. The region Ω is the void between the scattering surface and upper
boundary Γ. In the inverse problem of the next chapter, Γ will be the measurement surface.

The grating surface S is assumed to be invariant in the z direction and periodic in the x

22

direction with period Λ. We assume that the grating surface can be represented by the curve

S = {(x, y) ∈ R2 ∶ y = f (x), 0 < x < Λ},

where f is a periodic function with period Λ. Furthermore, we assume the grating surface
is small and smooth deformation of a plane surface, therefore f takes the form

f (x) = εg(x),

g ∈ C 2 (R),

where ε is sufficiently and is called the surface deformation parameter. This deformation is
not restricted to the case of a local perturbation of a plane surface. Rather, it is considered
to be a general perturbation of a planar surface, with nonlocal perturbation permitted.
Let h > max0<x<Λ f (x). In the next chapter, this will become the measurement distance, but
this term has no meaning for the forward problem. Define the upper surface

Γ = {(x, y) ∈ R2 ∶ y = h, 0 < x < Λ},

For now, it should just be considered a plane close to but not touching grating; eventually, a
transparent boundary condition will be placed here imposed. The last region to be described
is the space bounded below by the grating surface S and bounded above by the line Γ

Ω = {(x, y) ∈ R2 ∶ f (x) < y < h, 0 < x < Λ}.

The domain Ω is where the Helmholtz equation that characterizes propagation of the scattered field will be formulated.
Let Ω be filled with a homogeneous medium with a constant wavenumber κ > 0. An incoming
plane wave uinc = ei(αx−βy) is incident upon S from the top, where θ ∈ (−π/2, π/2) is the

23

incident angle and

α = κ sin θ

β = κ cos θ

For simplicity, this thesis shall only consider normal incidence since our method requires
only a single incident wave. Hence, α = 0, β = κ and the incident wave is reduced to
uinc = e−iκy .

It is easily verified that the incident wave satisfies the Helmholtz equation in the whole space

(∆ + κ2 )uinc = 0 in R2

(3.1)

As derived in Section 2.5, the diffraction of a time-harmonic electromagnetic wave in the TE
polarization can be modeled by the two-dimensional Helmholtz equation:

(∆ + κ2 )u = 0 in Ω,

(3.2)

where u is the total field consisting the incident field uinc and the diffracted field ud , i.e.,

u = uinc + ud .

(3.3)

Adding (3.1) and (3.2) and making use of the relation (3.3), it can been seen that the
diffracted field satisfies the same Helmholtz equation as uinc and u:

(∆ + κ2 )ud = 0 in Ω,

24

(3.4)

For a grating with the impedance boundary condition, the total field u satisfies

(∂n − iη)u = 0 on S,

(3.5)

where n = (n1 , n2 ) is the unit normal vector on S, given explicitly as
n1 = √

f ′ (x)
1 + [f ′ (x)]2

n2 = − √

,

1
1 + [f ′ (x)]2

,

(3.6)

η ∈ C is the impedance coefficient and is assumed to be known. To ensure existence of a
unique solution for the direct problem, we assume that Re η > 0 and Im η ≥ 0. When η = 0, the
impedance condition becomes the sound hard condition. The impedance boundary condition
will reduce to sound soft condition—perhaps better known as a perfect electric conductor
(PEC)—when η = ∞ . We are going to deal with the general impedance boundary condition
and present a general solution for any η. Our result will be applicable to specific cases, like
sound soft boundary condition and sound hard boundary condition.
It is known that for any given periodic function u(x) with period Λ, it has a Fourier series
expansion
u(x) = ∑ un ein( Λ )x ,
2π

n∈Z

where the Fourier coefficients are given by

un =

Λ
1
2π
u(x)e−in( Λ )x dx.
∫
Λ 0

Motivated by uniqueness, we shall seek the periodic solution of u, which implies ud is also
periodic with same period. Therefore, we may write the Fourier series of ud in x:

ud (x, y) = ∑ udn (y)ein(2π/Λ)x ,
n∈Z

25

(3.7)

When (3.7) is plugged into (3.4), an ordinary differential equation of ud is derived,
2πn 2 d
d2 udn
2
+ [κ − (
) ] un (y) = 0
dy 2
Λ
The diffracted wave is required to be a bounded outgoing wave in Ω. Mandating this, the
solution for each n can be readily obtained as

udn (y) = An eiβn y

(3.8)

Upon substitution into the Fourier series (3.7), we see that ud is a (infinite) combination of
plane waves:
ud (x, y) = ∑ An eiαn x+iβn y ,

(3.9)

n∈Z

valid for y > max0<x<Λ f (x). An is a complex number and the factors in the exponent are
given explicitly by
⎧
⎪
2
2 1/2
⎪
⎪
⎪ (κ − αn )
βn = ⎨
⎪
⎪
2
2 1/2
⎪
⎪
⎩ i(αn − κ )

2πn
αn =
Λ

for κ > ∣αn ∣,
for κ < ∣αn ∣.

Now we will introduce a transparent boundary condition on Γ. Taking the partial derivative
of ud with respect to y and then evaluating at y = h yields
∂y ud (x, h) = iβn ∑ An eiαn x+iβn h

(3.10)

n∈Z

Now given any u on Γ, define a boundary operator T which maps the Dirichlet data u(x, h)
to the Neumann data ∂y u(x, h):
(T u)(x) = ∑ iβn u(n) eiαn x .
n∈Z

26

For ud on Γ, we have
(T ud )(x) = ∑ iβn An eiαn x+iβn y

(3.11)

n∈Z

The transparent boundary condition on Γ is

∂y ud = T ud

(3.12)

∂y (u − uinc ) = T (u − uinc )

(3.13)

which is

Equivalently, a transparent boundary condition for the total field which bridges Dirichlet
data and Neumann data is derived:

(∂y − T )u = ψ

on Γ,

(3.14)

where ψ = ∂y uinc − T uinc = −2iκe−iκh .
To summarize, the surface scattering model can be reduced to the following boundary value
problem for the total field:

⎧
⎪
⎪
⎪
(∆ + κ2 )u = 0 in Ω,
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎨(∂n − iη)u = 0 on S,
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
(∂ − T )u = ψ on Γ,
⎪
⎩ y

(3.15)

This thesis considers two problems: the direct diffractive grating problem and the inverse
diffractive grating problem. The direct problem is to determine the total field u, with
knowledge of the incident field and the grating function f . The aim of this dissertation is
to study the inverse diffractive grating problem, which is to reconstruct the grating surface
function f from the measurement of the total field u on Γ when given the incident field uinc .
In particular, we are interested in the inverse problem in near-field regime, in which the
measurement distance h is much smaller than the wavelength λ = 2π/κ of the incident wave,
27

in the presence of noise.
We point out that u∣Γ should belong to H 1/2 . The trace theorem guarantees that if the
variational problem has a weak solution in H 1 , then the boundary information will belong
to H 1/2 . Since u∣Γ ∈ H 1/2 ⊂ L2 (Γ), the Fourier series of u will converge in the sense of L2 .
The forward problem has a C ∞ solution by elliptic regularity, since the boundary here
Lipschitz continuous (but not smooth because of the corners).
Before proceeding further, we prove a result about the boundary operator T , below.
Lemma 3.1. The linear boundary operator T : H 1/2 (Γ) → H −1/2 (Γ) is continuous.
Proof. For any u, v ∈ H 1/2 (Γ), we have
⟨T u, v⟩Γ = i2π ∑ βn (1 + n2 )−1/2 (1 + n2 )1/4 u(n) (1 + n2 )1/4 v (n)
n∈Z

Define
1

∣βn ∣
κ2 − n2 2
G(n) ≜
=
∣
∣
1 + n2
(1 + n2 )1/2
Claim. G(n) is bounded on Z.
Proof. Extend G to all of R by

1

κ2 − x2 2
G(x) = ∣
∣
1 + x2
and compute

√
√
⎧
2
⎪
⎛
x
1
+
x
κ2 − x2 ⎞
⎪
⎪
⎪
−
+
if ∣x∣ < ∣κ∣
⎪
⎪
⎪
1 + x2 ⎠
⎪ 1 + x2 ⎝ κ2 − x2
′
√
√
G (x) = ⎨
2
⎪
⎛
x
1
+
x
x2 − κ2 ⎞
⎪
⎪
⎪
+
if ∣x∣ ≥ ∣κ∣
⎪
⎪
⎪
1 + x2 ⎝ x2 − κ2
1 + x2 ⎠
⎪
⎩

Now observe that for a fixed κ, G′ has the following properties:
1. G′ (x) < 0 for x < −κ
28

(3.16)

2. G′ (x) > 0 for −κ < x < 0
3. G′ (x) < 0 for 0 < x < κ
4. G′ (x) > 0 for x > κ
Accordingly, G(x) is increasing on [−κ, 0] ∪ [κ, ∞] and decreasing on (−∞, −κ) ∪ (0, κ).
Consequently,

G(x) ≤ max{G(0), G(∞), G(−∞)}
= max{κ, 1}
Let C ≜ max{κ, 1} + 1 so that for all n ∈ Z, ∣G(n)∣ ≤ C. This proves the claim.

∎

Now, applying the result of the previous claim to (3.16), we have
∣⟨T u, v⟩Γ ∣ ≤ C ∣2π ∑ (1 + n2 )1/4 u(n) (1 + n2 )1/4 v (n) ∣
n∈Z

Then, by the Cauchy-Schwarz inequality, the above becomes

∣⟨T u, v⟩Γ ∣ ≤ 2πC

√
∑

2
(1 + n2 )1/2 ∣u(n) ∣

n∈Z

√

∑ (1 + n2 )1/2 ∣v (n) ∣

2

n∈Z

≤ 2πC ∣∣ u ∣∣1/2,Γ ∣∣ v ∣∣1/2,Γ
This shows that the L2 (Γ) duality pairing of T u with H 1/2 (Γ) functions forms a bounded
∗

linear functional on H 1/2 (Γ). Therefore, we conclude that T u ∈ (H 1/2 (Γ)) = H −1/2 (Γ).
Further, we have the bound

∣∣ T u ∣∣−1/2,Γ ≤

∣⟨T u, v⟩Γ ∣
v∈H 1/2 (Γ) ∣∣ v ∣∣1/2,Γ
sup

≤ C ∣∣ u ∣∣1/2,Γ

29

which establishes that T is continuous from H 1/2 (Γ) to H −1/2 (Γ).

∎

Lemma 3.2. The estimates

Re⟨T u, u⟩Γ ≤ 0

and Im⟨T u, u⟩Γ ≥ 0

hold for any u ∈ H 1/2 (Γ).
Proof. By the definition of the transparent boundary operator T , we have for any u ∈ H 1/2 (Γ)
that
2

⟨T u, u⟩Γ = i2π ∑ βn ∣u(n) ∣ .
n∈Z

Taking the real part gives
1/2

Re⟨T u, u⟩Γ = −2π ∑ ( ∣αn ∣ − κ2 )

2

∣u(n) ∣ ≤ 0

∣αn ∣>κ

and taking the imaginary part gives
1/2

2
Im⟨T u, u⟩Γ = 2π ∑ (κ2 − ∣αn ∣ )

2

∣u(n) ∣ ≥ 0,

∣αn ∣<κ

∎

which yields the proof.

We are going to focus on the variational formulation for the direct diffractive grating problem.
Define the periodic functional space

Hp1 (Ω) = {u ∈ H 1 (Ω) ∶ u(0, y) = u(2π, y)}

(3.17)

which is a subspace of H 1 (Ω).
Now we multiply the Helmholtz equation by the complex conjugate of a test function v ∈
Hp1 (Ω):
∆uv + κ2 uv = 0
30

Then we integrate over Ω and integrate by parts, obtaining
∂u
2
∫Ω ∇u ⋅ ∇v dx − ∫∂Ω ∂n v dS = κ ∫Ω uv dx
Applying the boundary conditions on Γ and S, we get
∂u
∂u
2
∫Ω ∇u ⋅ ∇v dx − ∫S ∂n v dS − ∫Γ ∂n v dS = κ ∫Ω uv dx
which implies

∫Ω ∇u ⋅ ∇v dx − iη ∫S uv dS − ⟨T u, v⟩Γ − ⟨ρ, v⟩Γ = κ ∫Ω uv dx
2

1
Therefore we derive the variational formulation for the forward problem: find u ∈ HS,p
(Ω),

such that
aΩ (u, v) = ⟨ρ, v⟩Γ for all v ∈ Hp1 (Ω)

(3.18)

where
aΩ (u, v) = ∫ ∇u ⋅ ∇v dx − iη ∫ uv dS − ⟨T u, v⟩Γ − κ2 ∫ uv dx
S

Ω

(3.19)

Ω

Lemma 3.3. The direct scattering problem has at most one solution.
Proof. Because the problem is linear, if u1 and u2 are both solutions for this direct scattering
problem, then so is u1 − u2 . It therefore suffices to show that u = 0 in Ω if ψ = 0 (no source
term). It follows from Green’s theorem that we have

0 = ∫ (u∆¯
u − u¯∆u) = ∫
Ω

∂Ω

(u∂n u¯ − u¯∂n u)

= ∫ (u∂n u¯ − u¯∂n u) + ∫ (u∂n u¯ − u¯∂n u)
S

Γ

= −2iRe(η) ∫ ∣u∣2 − 2iIm⟨T u, u⟩Γ ,
S

which yields that u = 0 on S since Re(η) > 0 and Im⟨T u, u⟩Γ ≥ 0. The impedance boundary
31

condition on S gives that ∂n u = 0 on S, which implies that u ≡ 0 in Ω.
Theorem 3.1. The direct scattering problem has a unique weak solution in Hp1 (Ω).
Proof. Decompose the bilinear form into aΩ = a1 − κ2 a2 , where
v − ⟨T u, v⟩Γ
v − iη ∫ u¯
a1 (u, v) = ∫ ∇u ⋅ ∇¯
S

Ω

and
a2 (u, v) = ∫ u¯
v.
Ω

Simple calculations yield

Rea1 (u, u) = ∫ ∣∇u∣2 + Im(η) ∫ ∣u∣2 − Re⟨T u, u⟩Γ ≥ ∫ ∣∇u∣2
Ω

S

Ω

and
Ima1 (u, u) = −Im(η) ∫ ∣u∣2 − Im⟨T u, u⟩Γ ≤ −Im(η) ∫ ∣u∣2 .
S

S

We conclude that a1 is coercive from
∣a1 (u, u)∣ ≥ C(∣Rea1 (u, u)∣ + ∣Ima1 (u, u)∣)
2

≥ C (∫ ∣∇u∣2 + ∫ ∣u∣2 ) = C ∣∣ u ∣∣1,Ω ,
Ω

S

where the last inequality may be obtained from the Poincare inequality.
Define an operator K ∶ L2 (Ω) → H 1 (Ω) by
a1 (Ku, v) = a2 (u, v) for all v ∈ Hp1 (Ω).

32

∎

Using the Lax-Milgram Theorem, we have

∣∣ Ku ∣∣1,Ω ≤ C ∣∣ u ∣∣0,Ω .
Thus K is bounded from L2 (Ω) to H 1 (Ω). Because H 1 (Ω) is compactly embedded into
L2 (Ω), K ∶ L2 (Ω) → L2 (Ω) is a compact operator.
Define a function ξ ∈ L2 (Ω) by requiring ξ ∈ H 1 (Ω) and satisfying
a1 (ξ, v) = ⟨ρ, v⟩Γ

for all v ∈ Hp1 (Ω).

It follows from the Lax-Milgram theorem again that

∣∣ ξ ∣∣1,Ω ≤ C ∣∣ ρ ∣∣0,Ω .

Using the operator K, we can see that the direct scattering problem is equivalent to finding
u ∈ L2 (Ω) such that
(I − κ2 K)u = ξ.
It follows from the Fredholm alternative and the uniqueness result that we have

∣∣ u ∣∣0,Ω ≤ C ∣∣ ξ ∣∣0,Ω ,

which gives
u = ξ + κ2 Ku.
The proof is completed by noting that K is a bounded operator from L2 (Ω) to H 1 (Ω).

33

∎

3.2

Transformed field expansion

In this section, we shall introduce the transformed field expansion to convert the domain to
a rectangular region to facilitate derivation of the analytic solution for the direct problem.
The expansion of solution is given as a power series of the parameter ε and is crucial in the
inversion formula.
Consider the changes of variables:

x˜ = x,

y˜ = h (

y−f
),
h−f

which maps the domain Ω to the rectangle

D = {(˜
x, y˜) ∈ R2 ∶ 0 < y˜ < h, 0 < x˜ < Λ}.
Notice that the grating surface y = f (x) is mapped into the plane surface y˜ = 0 while the
measurement plane y = h is mapped to itself as y˜ = h under this change of variables. Now
the grating problem can be restated in this transformed coordinate. First, we prove the
following.
Theorem 3.2. Let w(˜
x, y˜) = u(x, y). Then w is a periodic function with period Λ.
Proof. For arbitrary x ∈ R, y ∈ [0, h], we have w(˜
x+Λ, y˜) = u(x+Λ, y), and u(x+Λ, y) = u(x, y)
since u is periodic, therefore w(˜
x + Λ, y˜) = u(x, y) where u(x, y) is w(˜
x, y˜) by definition. We
can conclude w(˜
x + Λ, y˜) = w(˜
x, y˜), which implies w is periodic with Λ.
Now we are going to focus on the scattering problem for w. By the chain rule,
∂w ∂w ∂ x˜ ∂w ∂ y˜
=
+
∂x ∂ x˜ ∂x ∂ y˜ ∂x
∂w ∂w
−f ′
=
+
h(y − h)
∂ x˜ ∂ y˜
(h − f )2

ux =

34

∎

=

∂w ∂w ′ h − y
+
hf
∂ x˜ ∂ y˜
(h − f )2

and
∂w ∂w ∂ x˜ ∂w ∂ y˜
=
+
∂y ∂ x˜ ∂y ∂ y˜ ∂y
∂w h
=
∂ y˜ h − f

uy =

So we have
∂x = ∂ x˜ − f ′ (
∂y = (

h − y˜
) ∂ y˜,
h−f

h
) ∂ y˜.
h−f

Therefore,
′

∂w ∂w ′ h − y
uxx = (
+
hf
)
∂ x˜ ∂ y˜
(h − f )2 x
∂ 2w
∂ 2w
∂w
h−y
h−y ′
(hf ′
)
=
+
hf ′
+
2
∂ x˜∂x ∂ y˜∂x
(h − f )
∂ y˜
(h − f )2 x
2

′

∂ 2w
∂ 2w
∂w
∂ 2w
h−y
h−y
′ h−y
+
2
hf
+
(hf ′
) +
(hf ′
)
=
2
2
2
2
∂ x˜
∂ x˜∂ y˜
(h − f )
∂ y˜
(h − f )
∂y
(h − f )2 x
where the derivative in the last term
′

h−y
h(h − y)f ′′ (h − f )2 + hf ′ (h − y)2(h − f )f ′
(hf
)
=
(h − f )2 x
(h − f )4
h(h − y)
= (f ′′ (h − f ) + 2f ′2 )
(h − f )3
′

Note that
y−f
h−f
h(h − y) h(h − f + f − y) h h−f − h h−f
h − y˜
=
=
=
(h − f )3
(h − f )3
(h − f )2
(h − f )2

35

So
′

h−y
h − y˜
) = (f ′′ (h − f ) + 2f ′2 )
(hf
2
(h − f ) x
(h − f )2
′

We may differentiate again and obtain the needed second derivatives:
2

∂ 2w
∂ 2w
h − y˜ ∂w
∂ 2w
h−y
′ h−y
+
2
hf
+
(hf ′
) + [(f ′′ (h − f ) + 2f ′2 )
]
uxx =
2
2
2
2
∂ x˜
∂ x˜∂ y˜
(h − f )
∂ y˜
(h − f )
(h − f )2 ∂y
=

∂ 2w
˜ ∂ 2w
˜ 2 ∂ 2w
h − y˜ ∂w
′h−y
′h−y
+
2f
+
(f
)
+ [(f ′′ (h − f ) + 2f ′2 )
]
2
∂ x˜
h − f ∂ x˜∂ y˜
h − f ∂ y˜2
(h − f )2 ∂y

∂ ∂w h
(
)
∂y ∂ y˜ h − f
h2 ∂ 2 w
=
(h − f )2 ∂ y˜2

uyy =

Since uxx + uyy + κ2 u = 0, it can be derived that w, upon dropping the tilde, satisfies the
following equation in D:

(c1 ∂xx + c2 ∂yy + c3 ∂xy + c4 ∂y + c1 κ2 ) w = 0,

(3.20)

where

c1 = (h − f )2

c3 = −2f ′ (h − y)(h − f )

c2 = [f ′ (h − y)]2 + h2

c4 = −(h − y)[f ′′ (h − f ) + 2(f ′ )2 ].

Combing (3.5) and (3.6), the impedance boundary condition becomes
f
f ′ (1 − ) ∂x w − (1 + (f ′ )2 ) ∂y w
h
f
1/2
− iη (1 + (f ′ )2 ) (1 − ) w = 0 on y = 0.
h

36

(3.21)

The transparent boundary condition (3.14) reduces to
f
∂y w = (1 − ) (T w + ρ)
h

on y = h.

(3.22)

By the small perturbation assumption, we consider a formal expansion of w in a power series
of ε:
∞

w(x, y; ε) = ∑ wk (x, y) εk .

(3.23)

k=0

Substituting f = εg into each of c1 , c2 , c3 , and c4 and inserting the power series expansion
(3.23) into (3.20)–(3.22), we derive the recursion equation for wk :
(∆ + κ2 )wk = vk

in D,

(3.24)

where
1
vk = [2g∂xx + 2g ′ (h − y)∂xy + g ′′ (h − y)∂y + 2κ2 g]wk−1
h
1
− 2 [g 2 ∂xx + (g ′ )2 (h − y)2 ∂yy + 2gg ′ (h − y)∂xy − (2(g ′ )2 − gg ′′ ) (h − y)∂y + κ2 g 2 ]wk−2 .
h
The impedance boundary condition (3.21) becomes

(∂y + iη)wk = ϕk

where

on y = 0,

(3.25)

g
g
ϕk = [g ′ ∂x − iη ( )] wk−1 − [g ′ ( ) ∂x + (g ′ )2 ∂y ]wk−2
h
h
⌊ k2 ⌋

1
g
− iη ∑ ( 2 )(g ′ )2l × [wk−2l − ( ) wk−2l−1 ] .
l
h

l=1

Here ⌊l⌋ is the largest integer not greater than m and the generalized binomial coefficient
1
2

( )=
l

1 1
2(2

− 1)⋯( 12 − l + 1)
.
l!

37

The transparent boundary condition (3.22) reduces to

∂y wk − T wk = ψk

on y = h,

(3.26)

where
ψ0 = ψ,
g
ψ1 = − ( ) (T w0 + ψ),
h
g
ψk = − ( ) T wk−1 ,
h

⎫
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎬
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
k ≥ 2. ⎪
⎪
⎭

(3.27)

It is understood that wk = 0 for all k < 0 in the above recurrence relations. For each wk ,
the equation (3.24) and the boundary conditions (3.25) and (3.26) require only the previous
solutions wk−1 , wk−2 , . . . , w0 . Hence the transformed diffraction problem (3.24)–(3.27) can be
solved recursively starting from k = 0.

3.3

Fourier series expansion

For fixed k, we shall derive an analytic solution of wk for the problem (3.24)–(3.27), considering vk , ϕk , ψk as known functions. Since wk , vk , ϕk , ψk are periodic functions in x with
period Λ, they have the following Fourier series expansions
(n)

wk (x, y) = ∑ wk (y)eiαn x ,
n∈Z
(n)

vk (x, y) = ∑ vk (y)eiαn x ,
n∈Z
(n)

ϕk (x, y) = ∑ ϕk (y)eiαn x ,
n∈Z
(n)

ψk (x) = ∑ ψk eiαn x .
n∈Z

38

Substituting these expansions into (3.24), (3.25), and (3.26), we obtain a two-point boundary
(n)

value problem for wk :
⎫
⎪
⎪
⎪
0 < y < h,⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
(n)
⎪
⎪
dwk
(n)
(n)
⎬
+ iηwk = ϕk , y = 0,
⎪
dy
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
(n)
⎪
⎪
dwk
⎪
(n)
(n)
⎪
⎪
− iβn wk = ψk , y = h.
⎪
⎪
dy
⎭
(n)

d2 w k
(n)
(n)
+ βn2 wk = vk ,
2
dy

(3.28)

By taking the Fourier series, the two-dimensional grating problem for the Helmholtz equation
reduces to a one-dimensional two-point boundary value problem in the y direction, which
can be solved analytically. Using Theorem B.1 in Appendix, we obtain an explicit solution
of (3.28).
Theorem 3.3. The two-point boundary value problem (3.28) has a unique solution, which
is given by
(n)

(n)

(n)

(n)

(n)

h

(n)

(n)

wk (y) = K1 (y)ϕk − K2 (y)ψk + ∫ K3 (y, z)vk (z) dz
0
where
(n)

K1 (y) =
(n)

K2 (y) =

eiβn y
i(βn + η)
eiβn (h−y) (n)
K (y)
i(βn + η)

⎧
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
(n)
K3 (y, z) = ⎨
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎩

eiβn (y−z) (n)
K (z), z < y,
i(βn + η)
eiβn (z−y) (n)
K (y), z > y,
i(βn + η)

39

(3.29)

and
K (n) (t) =

(βn + η) (βn − η) 2iβn t
+
e
.
2βn
2βn

It follows from Theorem 3.3 that the power series (3.23) along with the solution representation (3.29) gives an analytic solution of the direct grating problem (3.20)–(3.22).
Remark 3.1. Theorems 3.1 and 3.3 together give an explicit analytic solution of the grating
(n)

surface problem. From the definition of Kj , it is clear that when ∣αn ∣ < κ, the solution
component at n is a propagating wave. When ∣αn ∣ > κ, the solution component becomes an
∎

evanescent wave, which decays exponentially in the y direction.

Remark 3.2. Notice that Re η ≠ 0 and Im η ≠ 0 while βn is either a real number or a pure
imaginary number. Hence we have η ≠ βn for all n ∈ Z.

∎

Remark 3.3. When the impedance boundary condition is reduced to sound soft condition
and f ∈ C 2 , the convergence of power series for w has been proved in [9].
(n)

Based on (3.29), we shall derive explicit expressions for solutions w0
an important role in the reconstruction formula.

40

∎
(n)

and w1 , which play

Chapter

4

Inverse Problem
4.1

Reconstruction formula

The previous chapter presented an analytic solution of forward problem based on a transformed field expansion and power series. We turn attention to the inverse problem. Suppose
that the total field is measured on the boundary Γ with noise. Assume that the noisy data
takes the form
uδ (x, h) = u(x, h) + O(δ),
where u(x, h) is the exact data and δ is the noise level.
Since the the information is collected on y = h, it is natural to evaluate the power series
(3.23) at y = h. Noting w(x, h) = u(x, h) and wδ (x, h) = uδ (x, h), we have
wδ (x, h) = w0 (x, h) + εw1 (x, h) + O(ε2 ) + O(δ).

(4.1)

Rearranging (4.1), and dropping O(ε2 ) and O(δ) yield
εw1 (x, h) = wδ (x, h) − w0 (x, h),

41

(4.2)

which is the linearization of the nonlinear inverse problem and enables us to find an explicit
reconstruction formula for the linearized inverse problem.
Based on (4.2), we next shall derive the analytical solution for the leading term w0 and then
deduce an equation relating w1 . From there, we can obtain the scattering surface function
f for an explicit inversion formula.
Recalling (3.2), (3.26), and (3.27), we have

v0 = 0,

ψ0 = −2iκe−iκh ,

ϕ0 = 0,

whose Fourier coefficients are
(n)

v0

= 0,

(n)

(n)

ϕ0 = 0,

ψ0

= −2iκe−iκh δ0n .

Here δkn is the Kronecker delta. It follows from (3.29) that
(n)

(n)

(n)

w0 (y) = −K2 ψ0

(n)

= 2iκe−iκh K2 (y)δ0n ,

(4.3)

which gives
(n)

(0)

w0 (x, y) = ∑ w0 (y)eiαn x = 2iκe−iκh K2 (y)
n∈Z

= e−iκy + (

κ + η iκy
)e .
κ−η

(4.4)

κ+η
Denote uref (x, y) = ( κ−η
) eiκy . It is clear that uref represents an outgoing reflected field, that

is, w0 (x, y) = uinc + uref . It can also be verified that w0 satisfies the impedance condition on
the flat plane y = 0.
Remark 4.1. Physically, the leading term w0 shows how incident field uinc interacts with
the impedance plane surface. It provides insight into how the incident wave is reflected by
42

the grating surface without any deformation (i.e., ε = 0), or when S is a smooth flat plane.
Mathematically, it satisfies the Helmholtz equation,

(∆ + κ2 )w0 = 0 in D,

with the impedance boundary condition

(∂y − iη)w0 = 0 on y = 0,

and transparent boundary condition

∂y w0 − T w0 = ψ

on y = h,

where ψ = −2iκe−iκh . It can be verified that the solution w0 consists of the incident uinc and
∎

the reflected field uref .
Now we are going to derive w1 using (3.2) and the equations of Remark 4.1,
1 ′′
[g (x)(h − y)∂y + 2κ2 g]w0
h
2κ2 −iκy
κ + η iκy
iκ(h − y) −iκy
κ + η iκy ′′
=
[e
+(
) e ] g(x) −
[e
−(
) e ] g (x)
h
κ−η
h
κ−η

v1 (x, y) =

and its Fourier coefficient
(n)

v1 (y) =

2κ2 −iκy
κ + η iκy
iκαn2 (h − y) −iκy
κ + η iκy
[e
+(
) e ] gn +
[e
−(
) e ] gn ,
h
κ−η
h
κ−η

(4.5)

where gn is the Fourier coefficient of g (note the difference in notation from the usual g (n) !).
Following (3.26) and (4.4) yield
g
ϕ1 (x) = [g ′ ∂x − iη( )]w0 (x, 0)
h

43

=−

iηg(x)
2iκη
w0 (x, 0) = − [
] g(x)
h
h(κ − η)

which has Fourier coefficients
(n)

ϕ1 = − [

2iκη
] gn .
h(κ − η)

(4.6)

Similarly, we have from (3.27) and (4.4) that
g
g
ψ1 (x) = − ( ) (T w0 + ψ0 ) = − ( ) (T w0 + ψ)
h
h
g
iκ −iκh
κ + η iκh
= − ( ) ∂y w0 (x, h) = [e
−(
) e ] g(x)
h
h
κ−η
which has Fourier coefficients
(n)

ψ1

=

iκ −iκh
κ + η iκh
[e
−(
) e ] gn .
h
κ−η

(4.7)

Letting k = 1 and evaluating at y = h in (3.29), we get
(n)

(n)

(n)

(n)

h

(n)

w1 (h) = K1 (h)ϕ1 − K2 (h)ψ1 + ∫

0

(n)

(n)

K3 (h, z)v1 (z)dz.

(4.8)

where
ieiβn h
,
η − βn
i
(n)
K2 (h) =
K (n) (h)
η − βn
i
i η + βn 2iβn h
=−
+
e
2βn 2βn η − βn
eiβn (h−z) (n)
(n)
K3 (h, z) =
K (z)
i(βn − η)
eiβn (h−z) eiβn (h+z) βn + η
=
+
2iβn
2iβn βn − η
(n)

K1 (h) =

(n)

(n)

(n)

(n)

h

(n)

(n)

By letting M1 = K1 (h)ϕ1 , M2 = K2 (h)ψ1 , M3 = ∫0 K3 (h, z)v1 (z)dz, (4.8) could

44

(n)

be rewritten as w1 (h) = M1 − M2 + M3 . We continue simplifying M1 ,M2 , and M3 :
(n)

(n)

M1 = K1 (h)ϕ1
=

2κηeiβn h
gn ,
(η − βn )(κ − η)h

(n)

(n)

M2 = K2 (h)ψ1
=

κgn
[(η − βn )(κ − η)e−iκh + (η + βn )e2iβn h (κ + η)eiκh
2βn h(κ − η)(η − βn )
− (η − βn )(κ + η)eiκh − (η + βn )e2iβn h (κ − η)e−iκh ],

M3 = ∫

h

0

=

(n)

(n)

K3 (h, z)v1 (z)dz

h
h
κeiβn h
αn2
κ2 gn eiβn h
Kn+ (z; η) dz +
(h − z)Kn− (z; η) dz
∫
∫
ihβn (βn − η)(κ − η) 0
(βn − η)(κ − η) 2hβn 0

I+

I−

≜ N1

≜ N2

where
Kn± (z; η) = [(βn − η)e−iβn z + (βn + η)eiβn z ] × [e−iκz (κ − η) ± eiκz (κ + η)]
Performing the integration yields

I+ =

(βn + η)(κ + η)
(βn + η)(κ − η)
(1 − ei(βn +κ)h ) +
(1 − ei(βn −κ)h )
(βn + κ)
(βn − κ)
(βn − η)(κ − η)
(βn − η)(κ + η)
(1 − e−i(βn +κ)h ) −
(1 − e−i(βn −κ)h )
−
(βn + κ)
(βn − κ)

I− = 2i(κ2 − η 2 ) +

1
A
2hβn

45

with
(βn + η)(κ − η)(βn + κ)
(βn + η)(κ + η)(βn − κ)
(1 − e+i(βn +κ)h ) −
(1 − e+i(βn −κ)h )
(βn + κ)
(βn − κ)
(βn − η)(κ − η)(βn − κ)
(βn − η)(κ + η)(βn + κ)
(1 − e−i(βn +κ)h ) +
(1 − e−i(βn −κ)h )
−
(βn + κ)
(βn − κ)

A=

Combining N1 and N2 , we get
M3 =

2iκ(κ2 − η 2 )eiβn h
κeiβn h
gn +
×
(η − βn )(κ − η)
2hβn (βn − η)(κ − η)
gn [ (βn + η)(κ + η)(1 − e+i(βn +κ)h ) − (βn + η)(κ − η)(1 − e+i(βn −κ)h )
− (βn − η)(κ − η)(1 − e−i(βn +κ)h ) + (βn − η)(κ + η)(1 − e−i(βn −κ)h )]

Subtracting M2 from M3 leads to
M3 − M2 =

2ηκeiβn h
2iκ(κ2 − η 2 )eiβn h
gn +
gn
(η − βn )(κ − η)
h(βn − η)(κ − η)

Finally, we calculate M1 − M2 + M3 .
(n)

w1 (h) = M1 − M2 + M3
2κηeiβn h
2iκ(κ2 − η 2 )eiβn h
2ηκeiβn h
gn −
gn +
gn
h(η − βn )(κ − η)
(βn − η)(κ − η)
h(βn − η)(κ − η)
−2iκ(κ + η)eiβn h
=
gn
(βn − η)
=

Therefore, an elegant equation relating Fourier coefficient of w1 and g is obtained:
(n)

w1 (h) =

2iκ(η + κ) iβn h
e gn .
(η − βn )

(4.9)

Remark 4.2. The impedance boundary condition (3.5) reduces to the perfectly conducting

46

boundary condition u = 0 on S as η → ∞, and (4.9) becomes
(n)

w1 (h) = 2iκeiβn h gn ,
∎

which is consistent with the result in [10].

We recall that f = εg so that fn = εgn , where fn is the Fourier coefficient of f . Plugging
(4.3) and (4.9) into (4.2), we deduce that

fn =

i(βn − η)
(n)
(n)
[w (h) − w0 (h)] e−iβn h ,
2κ(κ + η) δ

(n)

(4.10)
(n)

where wδ (h) is the Fourier coefficient of the noisy data wδ (x, h) and w0 (h) is the Fourier
coefficient of w0 (x, h) given as
(n)

w0 (h) = [e−iκh + (

κ + η iκh
) e ] δ0n .
κ−η

Remark 4.3. The inversion formula (4.10) includes cases for sound soft and sound hard
boundary conditions. For sound hard boundary conditions, i.e., η = 0, the inversion formula
reduces to
fn =

iβn
(n)
(n)
[w (h) − w0 (h)] e−iβn h ,
2κ2 δ

For the sound soft boundary condition, η = ∞, the inversion formula reduces to
fn = −

i
(n)
(n)
[w (h) − w0 (h)] e−iβn h .
2κ δ
∎

Remark 4.4. The inversion formula (4.10) accounts for both propagation wave modes and

47

evanescent wave modes. More explicitly, we have
√
⎧
√
i( κ2 − αn2 − η) (n)
⎪
(n)
2
2
⎪
⎪
[wδ (h) − w0 (h)] e−i κ −αn h
⎪
⎪
⎪ √2κ(κ + η)
fn = ⎨
√
⎪
αn2 − κ2 + iη (n)
⎪
(n)
⎪
α2n −κ2 h
⎪
(h)]
e
(h)
−
w
−
[w
⎪
0
δ
⎪
2κ(κ
+
η)
⎩

for ∣αn ∣ < κ
(4.11)
for ∣αn ∣ > κ

The low frequency modes of the scattering surface function f come from propagating (bounded)
waves, while the evanescent modes generate the high frequency modes. The factor on the
right hand side of (4.10), e−iβn h , reflects the ill-posedness of this problem. In other word, the
measurement noise will be exponentially amplified for high frequency modes ∣αn ∣ > κ.
It follows from the definition of βn and (4.10) that it is a well-posed problem to seek reconstruction of the Fourier coefficients fn for which ∣αn ∣ < κ, since small variations in the
measured data will not be amplified and lead to large errors in the reconstruction. However,
the resolution of the reconstruction f is limited by the given wavenumber κ. In contrast,
it is severely ill-posed to reconstruct those Fourier coefficients fn with ∣αn ∣ > κ, since small
variations in the data will be exponentially amplified and lead to huge errors in the reconstruction. However, these ill-behaving modes contribute to the super resolution of the
reconstructed function f .
To obtain a stable and super-resolved reconstruction, we may adopt a regularization to
suppress the exponential growth. One remedy is to make h as small as possible, i.e., measure
the data at the distance which is as close as possible to the grating surface, or bring the
scanning tip close to the scattering surface. This is exactly the idea of near-field optics. The
other remedy is to adopt a commonly used regularization, such as cut-off regularization or
Tikhonov regularization[25]. Following [14], we consider the spectral cut-off regularization.
Define the signal-to-noise ratio (SNR) by
SNR = min{ε−2 , δ −1 }.

48

For fixed h, the cut-off frequency ω is chosen in such a way that
e(ω

2 −κ2 )1/2 h

= SNR,

which implies that the spatial frequency will be cutoff for those below the noise level and
the surface deformation parameter. More explicitly, we have
log SNR 2
ω
= [1 + (
)]
κ
κh

1/2

,

(4.12)

which indicates ω > κ as long as SNR > 0 and super resolution may be achieved.
Taking into account the frequency cut-off, we may have a regularized reconstruction formulation:

fn =

i(βn − η)
(n)
(n)
[wδ (h) − w0 (h)] e−iβn h χn (αn ),
2κ(κ + η)

(4.13)

where the characteristic function
⎧
⎪
⎪
⎪
⎪
⎪ 1 for ∣αn ∣ ≤ ω,
χn (αn ) = ⎨
⎪
⎪
⎪
0 for ∣αn ∣ > ω.
⎪
⎪
⎩
Once fn have been computed, the grating surface function can be approximated by
f (x) ≈ ∑ fn eiαn x .

(4.14)

∣αn ∣≤ω

Hence, only two fast Fourier transforms are needed to reconstruct the grating surface func(n)

tion: one for the data to obtain wδ (h) and another to obtain the approximated f .
Remark 4.5. When the impedance boundary condition reduces to sound hard boundary

49

condition, i.e. η = 0, the regularized inversion formula (4.13) is reduced to
fn =

iβn
(n)
(n)
[wδ (h) − w0 (h)] e−iβn h χn (αn ),
2
2κ

(4.15)

For the sound soft boundary condition, i.e.,η = ∞, the regularized inversion formula (4.13)
is reduced to

fn = −

i
(n)
(n)
[wδ (h) − w0 (h)] e−iβn h χn (αn )
2κ

(4.16)

∎

This completes solution to the inverse problem.

50

Chapter

5

Numerical Experiments
We have already derived an elegant reconstruction formula for the inverse problem. In this
chapter, we are going to present two numerical examples to illustrate the effectiveness of
our method for algorithmic implementation and dependence of resolution on parameters
measurement distance h, surface deformation parameter ε, and the noise level δ.
The near-field data u(x, h) is created by solving the direct grating problem using the finite
element method with a perfectly matched layer. The wavenumber is taken as κ = 2π which
corresponds to the wavelength λ = 1. The grating period is Λ = 1, i.e., Λ = λ, and the
impedance constant η = 2π(1 + i), which is merely illustrative and stands for a complex
impedance. To get the data we need on boundary Γ, we use synthetic data: first, we solve
the forward problem by FEM then get the solution, then use these data as information to
reconstruct the scattering function f (x). We present two numerical examples. The first
example only contains low Fourier modes and the second example contains high frequency
modes which necessitate cut-off regularization will be applied.
In all the figures, the plots are rescaled with respect to the wavelength λ. Due to the
unstructured triangular meshes, the wave field data u(x, h) is not equally spaced with respect
to x. We construct a curve u(x, h) by using the natural cubic spline interpolation formula
based on the computed discrete data u(x, h). The curve u(x, h) is evaluated at equally
51

spaced points xj and used as our synthetic scattering data. For the stability test, some
relative random noise is added to the data, i.e., the near-field measurement is updated with

uδ (x, h) = u(x, h)(1 + δ rand),
where rand represents a normally distributed random number in [−1, 1]. In the examples,
approximately 1% noise is added to the data.

5.1

Numerical Example 1: Low Frequency

The exact grating surface is given by f (x) = εg(x), where
g(x) = 0.5 sin(4πx) + 0.5 sin(6πx).

This is a simple example since the grating profile function g only contains a couple of low
Fourier modes. We examine the dependence of reconstruction results on the parameters h
and ε. Figures 5.1 and 5.2 display the reconstructed surfaces against the exact surfaces.

(a)

(b)

Figure 5.1 Numerical Experiment #1: ε = 0.01. Reconstruction of a grating with
finite Fourier modes: the exact surface (solid line) is plotted against the reconstructed
surface (dashed line) using different h and the same ε = 0.01: (a) h = 0.1λ; (b), h = 0.2λ.

52

Figure 5.1 above compares the reconstruction using two different measurement distances.
The examples of this figure fix ε = 0.01. (a) and (b) plot the reconstructed surfaces by
using h = 0.1λ and h = 0.2λ, respectively. It is clear that smaller h of (a) gives a better
reconstruction than what is seen in (b). In particular, the fine features of the grating surface
are completely recovered and sub-wavelength resolution is achieved when using h = 0.1λ. A
larger cut-off frequency ω may be used in the reconstruction when the measurement distance
h is smaller, allowing for inclusion of higher-frequency modes, and hence, more detail.
The next comparison fixes h = 0.1λ and consider reconstructions for two different ε. The
result is displayed in Figure 5.2 below.

(a)

(b)

Figure 5.2 Numerical Experiment #1: h = 0.1λ. Reconstruction of a grating with
finite Fourier modes: the exact surface (solid line) is plotted against the reconstructed
surface (dashed line) using different ε and the same h = 0.1λ: (a) ε = 0.01; (b), ε = 0.05.

Figure 5.2 shows the reconstruction surfaces from using ε = 0.01 and ε = 0.05 in (a) and (b),
respectively. Our method yields better results with the smaller ε of (a) than with larger ε
of (b). Because we use a linearized model (i.e., truncate all terms of ε2 or higher), smaller ε
gives a more accurate approximation of the original nonlinear model problem.

53

5.2

Numerical Example 2: High Frequency

This time, the exact grating surface is given by f (x) = εg(x), where
g(x) = 0.5esin(4πx) + 0.4esin(6πx) − 1.5.

This is a harder example as the grating profile function g contains much higher Fourier
modes. It is expected that larger cut-off frequency ω is desirable and thus even smaller h
is necessary in order to capture all the surface features. As in the preceding example, we
present comparisons of different ε and h values.

(a)

(b)

Figure 5.3 Numerical Experiment #2: ε = 0.01. Reconstruction of a grating with
infinite Fourier modes: the exact surface (solid line) is plotted against the reconstructed
surface (dashed line) using different h and the same ε = 0.01: (a) h = 0.05λ; (b), h = 0.1λ.

Figure 5.3 shows the reconstructed surfaces for h = 0.05λ and h = 0.1λ, respectively, both
with ε = 0.01. For the larger h value in (b), the reconstruction is unable to capture the fine
features in the middle part of the exact grating surface. Recall that this choice of h was
sufficient for the grating surface in Example 1 (the wavelength is the same in both examples).
In contrast, all the detailed features are almost reconstructed in (a) when using h = 0.05λ
because larger ω is indeed allowed to be taken.
54

(a)

(b)

Figure 5.4 Numerical Experiment #2: h = 0.1λ. Reconstruction of a grating with
infinite Fourier modes: the exact surface (solid line) is plotted against the reconstructed
surface (dashed line) using different ε and the same h = 0.1λ: (a) ε = 0.01; (b), ε = 0.05.

Our final comparison, above in Figure 5.4, compares the results by using different ε. Taking
the same h = 0.1λ, (b) and (c) show the reconstructed surfaces using ε = 0.01 and ε = 0.05,
respectively. Again, it can be seen from the subplots of Figure 5.4 that smaller ε yields better
reconstruction than larger ε due to the approximation error in the linearization procedure.
Nonlinear corrections may be needed in order to improve the reconstructions.

5.3

Summary and Conclusion

We have presented a simple, stable, and effective computational method for the imaging
of impedance grating surfaces and achieved sub-wavelength resolution. The grating surface model was assumed to be a small and smooth deformation of a planar surface. Using
the transformed field expansion, the complex grating surface problem was converted into
a forward-marching sequence of two-point boundary value problems. Thanks to the power
series expansion and integration method, we were able deduce an analytic solution for the
direct problem. By dropping higher order terms in power series, we linearized the nonlinear inverse problem and obtained an explicit reconstruction formula for both propagation
55

and evanescent wave modes. To suppress the exponential growth, an appropriate cut-off
frequency was chosen from SNR analysis after proper consideration of the deformation parameter, noise level, and the measurement distance. Fast Fourier transforms contributed to
implementing the reconstruction formula.
In conclusion, we have shown that the method presented works for general impedance grating surfaces. We presented two examples, one of which has finite Fourier modes and another
one has infinite Fourier modes, and investigated how the parameters influence the reconstructions. The numerical results show that super resolution may be achieved by using small
measurement distance, which is exactly the principle of near-field optical imaging.

56

APPENDIX
blank space

57

APPENDIX
A.1

Integration solution method

For self-containedness, the integrated solution method is briefly introduced to solve a twopoint boundary value problem. We refer to Zhang [43] for the details of the integrated
solutions of ordinary differential equation system and two-point boundary value problems.
Consider the two-point boundary value problem
⎧
⎪
⎪
u′ (y) + M (y)u(y) = f (y),
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎨
A0 u(y)∣y=0 = r0 ,
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
B1 u(y)∣y=h = s1 ,
⎪
⎩

(A.1)

where f (y) ∈ Cm are m-dimensional vector fields, r0 ∈ Cm1 and s1 ∈ Cm2 are given m1 - and
m2 -dimensional vector fields, respectively, M (y) ∈ Cm×m is an m×m matrix, and A0 ∈ Cm1 ×m
and B1 ∈ Cm2 ×m are full rank matrices with m1 + m2 = m, i.e., rankA0 = m1 and rankB1 = m2 .
Let Φ(y) be the fundamental matrix of the system
⎧
⎪
⎪
cΦ′ (y) + M (y)Φ(y) = 0,
⎪
⎪
⎪
⎨
⎪
⎪
⎪
⎪
Φ(0) = Im ,
⎪
⎩

(A.2)

where Im is the m × m identity matrix.
Theorem A.1. The two-point boundary value problem (A.1) has a unique solution if and
only if

⎡
⎤
⎢ A0 ⎥
⎥
⎢
⎥ ≠ 0.
det ⎢⎢
⎥
⎢
⎥
⎢B1 Φ(h)⎥
⎣
⎦

58

(A.3)

Let the pair of functions {A(y), r(y)} and {B(y), s(y)} be the integrated solutions of the
problems (A.1), then there exist D0 (A, y) ∈ Cm1 ×m1 and D1 (B, y) ∈ Cm2 ×m2 such that
⎧
⎪
⎪
′
⎪
⎪
⎪ A = AM + D0 A, A(0) = A0 ,
⎨
⎪
⎪
⎪
r′ = Af + D0 r, r(0) = r0 ,
⎪
⎪
⎩

(A.4)

⎧
⎪
⎪
′
⎪
⎪
⎪ B = BM + D1 B, B(h) = B1 ,
⎨
⎪
⎪
⎪
s′ = Bf + D1 s, s(h) = s1 .
⎪
⎪
⎩

(A.5)

Theorem A.2. If the two-point boundary value problem (A.1) has a unique solution, then
the matrix

⎡
⎤
⎢ A(y) ⎥
⎢
⎥
⎢
⎥ ∈ Cm×m
⎢
⎥
⎢
⎥
⎢B(y)⎥
⎣
⎦

is nonsingular.

Theorem A.3. The two-point boundary value problem (A.1) is equivalent to the linear
system

A.2

⎡
⎤
⎡
⎤
⎢ A(y) ⎥
⎢r(y)⎥
⎢
⎥
⎢
⎥
⎢
⎥ u(y) = ⎢
⎥.
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢B(y)⎥
⎢s(y)⎥
⎣
⎦
⎣
⎦

(A.6)

A two-point boundary value problem

Consider a two-point boundary value problem
⎧
⎪
⎪
u′′ + β 2 u = v, 0 < y < h,
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎨ u′ − iηu = ϕ, y = 0,
⎪
⎪
⎪
⎪
⎪
⎪
⎪
′
⎪
⎪
⎩ u − iβu = ψ, y = h,

59

(B.1)

where β ∈ C, η ∈ C, and η ≠ β.
Let u1 = u and u2 = u′ . The two-point boundary value problem (B.1) can be equivalently
formulated into
⎧
⎪
⎪
u′ + M u = v,
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎨ A0 u(0) = ϕ,
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎩ B1 u(h) = ψ,
where

⎡ ⎤
⎢u ⎥
⎢ 1⎥
u = ⎢⎢ ⎥⎥ ,
⎢u ⎥
⎢ 2⎥
⎣ ⎦

⎡ ⎤
⎢0⎥
⎢ ⎥
v = ⎢⎢ ⎥⎥ ,
⎢v ⎥
⎢ ⎥
⎣ ⎦

(B.2)

⎡
⎤
⎢ 0 −1⎥
⎢
⎥
⎥,
M = ⎢⎢
⎥
⎢β 2 0 ⎥
⎢
⎥
⎣
⎦

and
A0 = [−iη

B1 = [−iβ

1],

1].

Theorem A.4. The two-point boundary value problem (B.1) has a unique solution given by

u(y) = K1 (y)ϕ − K2 (y)ψ + ∫

h

0

K3 (y, z)v(z)dz,

where
K1 (y) =

eiβy
,
i(β − η)

⎧
⎪
⎪
⎪
⎪
⎪
⎪
⎪
K3 (y, z) = ⎨
⎪
⎪
⎪
⎪
⎪
⎪
⎪
⎩

K2 (y) =

eiβ(h−y)
K(y),
i(β − η)

eiβ(y−z)
K(z), z < y,
i(β − η)
eiβ(z−y)
K(y), z > y.
i(β − η)

Here
K(t) =

(β − η) (β + η) 2iβt
+
e .
2β
2β

60

(B.3)

Proof. It can be verified that there exists a non-singular matrix P such that
P −1 M P = N,

where

⎡
⎤
⎢−iβ 0 ⎥
⎢
⎥
⎥,
N = ⎢⎢
⎥
⎢
⎥
⎢ 0 iβ ⎥
⎣
⎦

⎡
⎤
⎢1
1 ⎥⎥
⎢
⎥,
P = ⎢⎢
⎥
⎢
⎥
⎢iβ −iβ ⎥
⎣
⎦

1
P −1 =
2iβ

⎡
⎤
⎢iβ 1 ⎥
⎢
⎥
⎢
⎥.
⎢
⎥
⎢
⎥
⎢iβ −1⎥
⎣
⎦

A simple calculation yields that the fundamental matrix of (A.2) is
⎡ iβy
⎤
⎢e
⎥
⎢
⎥ −1
⎢
⎥P ,
Φ(y) = P ⎢
⎥
⎢
⎥
−iβy
⎢
e ⎥⎦
⎣
which gives

⎡
⎤
⎢ A0 ⎥
−iη
1
⎢
⎥
⎥=
det ⎢⎢
⎥
⎢
⎥
⎢B1 Φ(h)⎥
−iβe−iβh e−iβh
⎣
⎦

= i(β − η)e−iβh ≠ 0.

It follows from Theorem A.1 that the two-point boundary value problem (B.2) and thus
(B.1) has a unique solution.
Let {A(y), r(y)} and {B(y), s(y)} be the integrated solutions of the problems (B.2). Taking
D0 = iβ,

D1 = −iβ,

we obtain from (A.4) that the integrated solutions satisfy
⎧
⎪
⎪
′
⎪
⎪
⎪ A = AM + iβA, A(0) = A0 ,
⎨
⎪
⎪
⎪
r′ = Av + iβr, r(0) = ϕ,
⎪
⎪
⎩

61

(B.4)

and
⎧
⎪
⎪
′
⎪
⎪
⎪ B = BM − iβB, B(h) = B1 ,
⎨
⎪
⎪
⎪
s′ = Bv − iβs, s(h) = ψ.
⎪
⎪
⎩

(B.5)

Upon solving the initial value problems (B.4) and (B.5), we obtain the integrated solutions
for A and B:
A = [A1 (y) A2 (y)],

B = [B1 (y) B2 (y)],

(B.6)

where

A1 (y) =

i(β − η) i(β + η) 2iβy
−
e ,
2
2

A2 (y) =

(β − η) (β + η) 2iβy
+
e ,
2β
2β

and
B1 (y) = −iβ,

B2 (y) = 1.

Once A and B are available, we may solve the initial value problems (B.4) and (B.5) and
obtain the integrated solutions for r and s:

r(y) = eiβy ϕ + ∫

y
0

eiβ(y−z) A2 (z)v(z)dz,

s(y) = eiβ(h−y) ψ − ∫

h

eiβ(z−y) v(z)dz.

(B.7)
(B.8)

y

It follows from Theorem A.3 that the two-point boundary value problem (B.2) is equivalent
to the linear system

⎡
⎤⎡ ⎤ ⎡ ⎤
⎢A1 A2 ⎥ ⎢ u ⎥ ⎢r⎥
⎢
⎥⎢ ⎥ ⎢ ⎥
⎢
⎥⎢ ⎥ = ⎢ ⎥.
⎢
⎥⎢ ⎥ ⎢ ⎥
⎢
⎥⎢ ⎥ ⎢ ⎥
⎢B1 B2 ⎥ ⎢u′ ⎥ ⎢s⎥
⎣
⎦⎣ ⎦ ⎣ ⎦

62

An application of Cramer’s rule yields

u=

rB2 − sA2
.
A1 B2 − B1 A2

(B.9)

A simple calculation yields
A1 B2 − B1 A2 = i(β − η).
Substituting (B.6)–(B.8) into (B.9), we deduce (B.3).

63

∎

REFERENCES
blank space

64

REFERENCES
[1] H. Ammari. Uniqueness theorems for an inverse problem in a doubly periodic structure.
Inverse Problems, 11(4):823, 1995.
[2] H. Ammari, G. Bao, and A. Wood. Analysis of the electromagnetic scattering from a
cavity. Japan Journal of Industrial and Applied Mathematics, 19(2):301–310, 2002.
[3] G. Bao. A unique theorem for an inverse problem in periodic diffractive optics. Inverse
Problems, 10:335–340, 1994.
[4] G. Bao, Z. Chen, and H. Wu. adaptive finite element method for diffration gratings.
Journal of the Optical Society of America A, 22:1106–1114, June 2005.
[5] G. Bao, L. Cowsar, and W. Masters. Mathematical modeling in optical science, volume 22 of Frontiers in Applied Mathematics. SIAM, 2001.
[6] G. Bao, D. C. Dobson, and J. A. Cox. Mathematical studies in rigorous grating theory.
Journal of the Optical Society of America A, 12(5):1029–1042, May 1995.
[7] G. Bao and A. Friedman. inverse problems for scattering by periodic structure, volume
132 of Archive for Rational Mechanics and Analysis. Springer, 1995.
[8] G. Bao, J. Gao, J. Lin, and W. Zhang. Mode matching for the electromagnetic scattering
from three-dimensional large cavities. IEEE transactions on antennas and propagation,
60(4):2004–2010, April 2012.
[9] G. Bao and P. Li. Convergence analysis in near-field imaging.
[10] G. Bao and P. Li. Near-field imaging of infinite rough surfaces. SIAM Journal on
Applied Mathematics, 73(6):2162–2187, 2013.
[11] G. Bao and P. Li. Near-field imaging of infinite rough surfaces in dielectric media. SIAM
Journal on Imaging Sciences, 7(2):867–899, 2014.
[12] G. Bao, P. Li, and J. Lv. Numerical solution of an inverse diffraction grating problem
from phaseless data. Journal of the Optical Society of America A, 30(3):293–299, Mar
2013.
[13] G. Bao, P. Li, and H. Wu. A computational inverse diffraction grating problem. Journal
of the Optical Society of America A, 29(4):394–399, Apr 2012.
[14] G. Bao and J. Lin. Near-field imaging of the surface displacement on an infinite ground
plane. Inverse Problem and Imaging, 7:377–396, May 2013.
[15] O. P. Bruno and F. Reitich. Numerical solution of diffraction problems: a method of
variation of boundaries. Journal of the Optical Society of America A, 10(6):1168–1175,
Jun 1993.
65

[16] P. S. Carney and J. C. Schotland. Inverse scattering for near-field microscopy. Applied
Physics Letters, 77(18):2798–2800, 2000.
[17] P. S. Carney and J. C. Schotland. Near-field tomography. In Gunther Uhlmann, editor,
Inside Out Inverse Problems and Applications, pages 133–168. Cambridge University
Press, 2003.
[18] S. N. Chandler-Wilde and P. Monk. Existence, uniqueness and variational methods
for scattering by unbounded rough surfaces. SIAM Journal on Mathematical Analysis,
37(2):598618, 2005.
[19] S. N. Chandler-Wilde and B. Zhang. A uniqueness result for scattering by infinite rough
surfaces. SIAM Journal of Applied Mathematics, 58:17741790, 1998.
[20] Z. Chen and H. Wu. An adaptive finite element method with perfectly matched absorbing layers for the wave scattering by periodic structures. SIAM Journal of Numerical
Analysis, 41:799–826, 2003.
[21] T. Cheng, P. Li, and Y. Wang. Near-field imaging of perfectly conducting grating
surfaces. preprint.
[22] D. Courjon. Near-field Microscopy and Near-Field Optics. Imperial College Press, 2003.
[23] D. Courjon and C. Bainier. Near field microscopy and near field optics. Reports on
Progress in Physics, 57(10):989, 1994.
[24] J. A. DeSanto and P.A. Martin. On the derivation of boundary integral equations
for scattering by an infinite one-dimensional rough surface. Journal of the Acoustical
Society of America, 102:67–77, 1997.
[25] H. Engl, M. Hanke, and A. Neubauer. Regularization of Inverse Problems. Kluwer
Academic Publishers, 1996.
[26] J. M. Jin. The Finite Element Method in Electromagnetics. Wiley & Sons, 2002.
[27] A. Kirsch. Uniqueness theorems in inverse scattering theory for periodic structures.
Inverse Problems, 10(1):145, 1994.
[28] P. Li. An inverse cavity problem for maxwell’s equations. Journal of Differential Equations, 252:3209–3225, 2012.
[29] P. Li, G. Bao, and J. Gao. Analysis of direct and inverse cavity scattering problems.
Numerical Mathematics: Theory, Methods and Applications, 4:419–442, 2011.
[30] P. Li and J. Shen. Analysis of the scattering by an unbounded rough surface. Mathematical Methods in the Applied Sciences, 35(18):2166–2184, 2012.
[31] P. Li and A. Wood. A two-dimensional helmhotlz equation solution for the multiple
cavity scattering problem. Journal of Computational Physics, 240:100–120, May 2013.

66

[32] P. Li, H. Wu, and W. Zheng. Electromagnetic scattering by unbounded rough surfaces.
SIAM Journal on Mathematical Analysis, 43(3):1205–1231, 2011.
[33] P. Li, H. Wu, and W. Zheng. An overfilled cavity problem for maxwell’s equations.
Mathematical Methods in the Applied Sciences, 35(16):1951–1979, 2012.
[34] A. Malcolm and D. P. Nicholls. A boundary perturbation method for recovering interface
shapes in layered media. Inverse Problems, 27(9):095009, 2011.
[35] A. Malcolm and D. P. Nicholls. A field expansions method for scattering by periodic
multilayered media. The Journal of the Acoustical Society of America, 129(4):1783–
1793, 2011.
[36] D. M. Milder. An improved formalism for wave scattering from rough surfaces. Journal
of the Acoustical Society of America, 89:529–541, 1991.
[37] D. P. Nicholls and F. Reitich. Shape deformations in rough-surface scattering: cancellations, conditioning, and convergence. Journal of the Optical Society of America A,
21(4):590–605, Apr 2004.
[38] D. P. Nicholls and F. Reitich. Shape deformations in rough surface scattering: Improved
algorithms. Journal of the Optical Society of America A, 21:606–621, 2004.
[39] R. Petit. Electromagnetic Theory of Gratings. Springer-Verlag, 1980.
[40] A.G. Voronovich. Wave scattering from rough surfaces. Springer series on wave phenomena. Springer-Verlag, 1994.
[41] K. F Warnick and W. C. Chew. Numerical simulation methods for rough surface scattering. Waves in Random Media, 11(1):R1–R30, 2001.
[42] R. J. Wombell and J. A. DeSanto. The reconstruction of shallow rough-surface profiles
from scattered field data. Inverse Problems, 7(1):L7, 1991.
[43] G. Q. Zhang. Integrated solutions of ordinary differential equation system and two-point
boundary value problems. Journal of Computational Mathematics, 3:245–254, 1981.

67