# Maxwell's equations UPML in FDTD with inhomogeneous media

I'm looking at matching the UPML (uniaxial perfectly matched layer) defined in Taflove&Hagness' Computational electrodynamics to an inhomogeneous media (inhomogeneous w.r.t. both $$\varepsilon$$ and $$\sigma$$, that extends into the PML), in FDTD with Yee's scheme (finite-difference time domain scheme, staggered in space and time by half a step, see https://en.wikipedia.org/wiki/Finite-difference_time-domain_method).

This UPML involves defining the additional quantities:

$$s_w = \kappa_w + \frac{\sigma_w}{i\omega\varepsilon_0}$$

where $$\kappa_w=1$$ usually, and $$\sigma_w$$ is an anisotropic conductivity that only changes along the $$w-$$axis, $$w\in\lbrace x,y,z\rbrace$$.

This is used in the following way, to obtain a "perfect match" to the interior domain, i.e. reflectionless transmission, such that the wave goes into the PML and starts dampening in amplitude due to the nonzero conductivity:

They (Taflove&Hagness) go on to detail how $$\sigma_x$$ is defined within the UPML (when "matched" to vacuum):

$$\sigma_x(x) = \left(\frac{x}{d}\right)^m \sigma_{x,\text{max}}$$

for a polynomial grading of order $$m$$ and a PML of thickness $$d$$. This of course assumes the grading starts with $$0$$, or $$\sigma_\text{min}=0$$. For $$\sigma\neq 0$$ inside the domain, one could of course instead define it as:

$$\sigma_x(x) = \sigma_{x,\text{min}} + \left(\frac{x}{d}\right)^m (\sigma_{x,\text{max}}-\sigma_{x,\text{min}})$$

Then, on pg. 304, it is mentioned:

Later in the book, on page 331, it is mentioned that the maximum sigma in the PML should be scaled by the $$\textbf{square root}$$ of the relative permittivity.

So my questions are:

1) Given a computational half-space (i.e. for $$x\leq0$$) with $$\varepsilon_\text{rel}=10$$ and $$\sigma$$=5, and the other half vacuum ($$\varepsilon_\text{rel}=1,\sigma=0$$), both extending into the UPML:

What would the sigma (say $$\sigma_y$$) grading look like inside the UPML for polynomial grading of order $$m=3$$? What are the $$\sigma_{min}$$ and $$\sigma_{max}$$ for each half space?

2) What about if $$\varepsilon=\varepsilon_0$$ in both halves, but $$\sigma$$=5 in one and $$\sigma=0$$ in the other? Is there any produced surface current inside the PML due to discontinuity? Why not? How is $$\sigma$$ graded polynomially? (i.e. what are $$\sigma_{min}$$ and $$\sigma_{max}$$?)

3) If we have only varying $$\varepsilon$$ in the domain (outside of the UPML), and constant $$\sigma$$, how does that match at the UPML? What is the polynomial grading of $$\sigma$$? (i.e. what are $$\sigma_{min}$$ and $$\sigma_{max}$$?)

All the information I'm trying to provide here is presented in the very same chapter of the cited book. I encourage you to go through it again.

In a general manner, remember that the PML is not at all a physical medium, it is just a numerical artifact to solve a pure numerical problem (the undesired reflection of the electromagnetic waves inside the domain); therefore, it can, in principle, have any value for the $$\sigma_w$$ and related parameters. Also note that, whenever there is a discontinuity in conductivity of (real physical) materials, there will be reflections, even if their wave impedance is the same. On the other hand, we want the PML to be as small as possible, to decrease the computational burden of the simulations. Therefore, it is desired that the PML starts with a relatively low conductivity (to minimize initial reflections) and also, at the same time, to atenuate the incoming wave as fast as possible. So, basically, we want a monotonically growing function to describe the $$\sigma_w$$ behaviour. Of course, this can be done in infinite ways. Extensive numerical experiments we carried out to find a very good expression for this function (see the references of the cited chapter), providing the polynomial gradings you showed. The order of the polynomial is defined by how good you need your PML to be.

As for the questions, the answer below should cover all three.

As you have pointed, $$\sigma_y$$ must be independent of $$x$$ and $$z$$. Also, as this values are arbitrary (because the PML is not physical), the simplest choice is to set $$\epsilon = \epsilon_0$$ for the entire PML domain. The value of $$\sigma_{y,max}$$ for a polynomial of order 3 would the be given by $$\sigma_{y,max}=-\frac{(3+1)\ln(R(0))}{2\eta d},$$ where $$R(0)$$ is the desired reflection error and $$d$$ is the PML thickness. Again, this equation is provided in the book same chapter, but I can't refer its number because apparently I am using a different edition than yours.

The formulation described above actually can be shown to be material independent. Therefore, in principle you could use for any physical system, including the ones described in 2 and 3.

The modified version of the $$\sigma_w$$ presented last is for applications involving UPML for nonlinear/dispersive materials, where, of course, the behaviour of permitivities is much more complex. In this case, I would recommend you implementing the Convolutional-PML, which is much more robust and efficient for this kind of simulations.

• Thanks for the answer, there are a few things unclear though. What is the minimal value of $\sigma_y$ at the beginning of the PML, let's say for example 1 that I presented above? If it's independent of the transversal directions, it will be "grading" the PML starting with the same $\sigma_{y,min}$, which means at both the mediums you have a discontinuity in $\sigma$ in the normal direction, going from the medium into the PML, which means you get reflections. Is this intended? – Daniel Jun 6 '19 at 7:42
• Furthermore, formula (7.78) I present is not used for a dispersive or nonlinear material, but simply a half-aluminum half-air example. (Section 7.9.3) Examples with more complex mediums start at section 7.11. Also, how is the impedance still matched if the $\varepsilon$ is discontinuous across the medium-PML interface? – Daniel Jun 6 '19 at 9:09
• As far as I have understood, both permitity and conductivity are included in the form of the PML tensor component proposed. In principle, I think this discontinuity in the PML conductivity could lead to some kind of unphysical charge/current, like the permitivity; nevertheless, the section 7.5.4 (named "Inhomogeneous media") explicitly state that this normalized description indeed leads to a material indepedent formulation of the PML (I have here the third edition of the book, which seems to be different than yours. In fact, my version does not even have the section 7.9.3 you pointed). – Bruno Anghinoni Jun 13 '19 at 22:31
• To be honest, I have not simulated physical lossy media, so I can not assure you, but I think if you describe your PML parameters based on the right media, for example, and use then for all PML domain, it should work. In this scenario, once you have the maximum conductivity with the equation I wrote, the minimum valuem (i.e. the first value) would be just $(1/d)^n$ times the maximum conductivity, assuming a n-th order polynomial grading. – Bruno Anghinoni Jun 13 '19 at 22:58