Analytical Methods for Causality Evaluation of Photonic Materials.

We comprehensively review several general methods and analytical tools used for causality evaluation of photonic materials. Our objective is to call to mind and then formulate, on a mathematically rigorous basis, a set of theorems which can answer the question whether a considered material model is causal or not. For this purpose, a set of various distributional theorems presented in literature is collected as the distributional version of the Titchmarsh theorem, allowing for evaluation of causality in complicated electromagnetic systems. Furthermore, we correct the existing material models with the use of distribution theory in order to obtain their causal formulations. In addition to the well-known Kramers-Krönig (K-K) relations, we overview four further methods which can be used to assess causality of given dispersion relations, when calculations of integrals involved in the K-K relations are challenging or even impossible. Depending on the given problem, optimal approaches allowing us to prove either the causality or lack thereof are pointed out. These methodologies should be useful for scientists and engineers analyzing causality problems in electrodynamics and optics, particularly with regard to photonic materials, when the involved mathematical distributions have to be invoked.

1. Introduction

The time-domain response of electromagnetic system has to be causal [1,2,3,4], i.e., a system cannot respond before the excitation starts. Hence, causality restricts the characteristics of any system, not only in the time domain but also in the frequency domain. Therefore, it is often stated in literature that real and imaginary parts of various complex parameters and characteristics (e.g., susceptibilities and frequency responses) are related by the Kramers–Krönig (K–K) integral relations [4]. The K–K relations are also valid between the modulus logarithm and the argument of the frequency-domain response (e.g., between the phase delay and the attenuation of a solution to a wave-propagation problem). However, it occurs that integrations in the K–K relations are not straightforward, and sometimes even impossible to perform. It stems from singular, highly oscillatory, or even diverging integrands, and infinite integration limits in the Hilbert transformation [5]. Therefore, the integrals cannot be evaluated in many important practical cases because they are divergent or undefined [6]. The problems concerning causality and the K–K relations are studied not only in electrodynamics [7,8,9,10,11] but also in acoustics [12,13,14,15,16], solid mechanics [17,18] and the control theory [19].

Although the subject of causality and dispersion of electromagnetic waves in photonic materials (e.g., dielectrics) was initiated in 1927 by the seminal papers of Kramers and Krönig [20,21], it has remained active in literature up to now [22]. Furthermore, classical books, devoted solely to this subject, are published [4]. As one can notice, although the book [4] was published in 1972, it has remained the initial point for many investigations referring to the subject, and constitutes the basic reference in this research area. On the other hand, a discussion related to origins of the Titchmarsh theorem recently appeared in literature [23]. This theorem establishes the equivalence between causality and the K–K relations. It occurs that it is a compilation of two important theorems, i.e., the Paley–Wiener theorem and the Marcel Riesz theorem. Therefore, the Authors of [23] have perceived the need for studying the subject with the use of rigorous mathematical tools. This approach, although often longer than the one presented in physics literature, has the advantage of being more precise.

Let us present some recent results related to causality in electrodynamics and optics. In [17], the K–K relations, which are valid for a general class of linear homogeneous or inhomogeneous media, are derived. That is, the proof of the K–K relations proceeds without a priori knowledge of wave velocity in the medium supporting the wave, when the frequency goes towards infinity. In [10], the K–K relations are derived for the effective index of modes propagating in optical waveguides. When material dispersion and absorption can be neglected within the frequency range of interest, the evanescent modes introduce an effective loss term in the K–K relations, meaning that these relations are valid even if material absorption is negligible within the frequency range of interest. In [11], a novel approach to the standard one is proposed for the derivation of the K–K relations for linear optical properties. That is, this approach is not based on contour integration and the Cauchy integral formula. Although this derivation still employs analytic behavior of the property under consideration, it employs only elementary properties of the Hilbert transformation to obtain the second formula of the K–K pair, from the Herglotz representation of the optical property as a Herglotz function. In [24], linear-response laws and causality within the time and frequency domains are analyzed in electrodynamics. It is demonstrated that one can violate causality in the frequency domain by making a vanishing-absorption approximation. Our contribution to this subject relies on generalization of the K–K relations for evaluation of causality in power-law media. In general, for square-integrable functions of the frequency, the validity of classical K–K relations is equivalent to causality in the time domain [1,4,17]. Things get complicated when the K–K relations are verified between the modulus logarithm and the argument. Then, the considered function does not belong to the class of square-integrable functions, and one can employ classical K–K relations with subtractions, but their satisfaction does not mean that the originally considered function is causal [13,17]. That is, the dispersion relation can be formulated for the considered system. It is based on the assumption that the subtracted logarithm of the response is causal, but the considered response function does not have to be causal. Therefore, we address this issue in [25], where the K–K relations are generalized towards non square-integrable functions. That is, the K–K relations with one and two subtractions are formulated. Their validity and satisfaction of additional assumptions imply causality of the considered function. Then, the formulated theory is applied to the analysis of electromagnetic media characterized by power-law frequency dispersion [25] and fractional-order (FO) models [26]. However, in our opinion, it might still be unclear for community members how to apply the developed mathematical tools to causality evaluation.

Therefore, we have decided to set in order the theory of causality and prepare the review of various analytical techniques allowing for evaluating causality of photonic materials. It includes the classical methods and models known from electrodynamics handbooks [7], as well as recently published results on causality evaluation [25,26]. Furthermore, we also consider various causality tests, which allow for establishing minimal limits for losses, for metamaterials [27,28]. Therefore, one can easily find their own approach to causality evaluation when a new model is formulated in classical electrodynamics. However, it is worth noticing that the procedure of canonical quantization of macroscopic electromagnetism [29] requires satisfaction of the K–K relations by dielectric functions of a linear, inhomogeneous, magnetodielectric medium. Hence, our review should hopefully also be useful in quantum-electrodynamics research. The paper begins with a short introduction of the notation, definitions, and basic mathematical theorems used throughout the paper. In Section 3, the problem of causality is presented, with various mathematical models applicable in electrodynamics and optics. However, we do not focus on physical motivations for these models, but we rather formulate mathematical approaches allowing for causality evaluation. Analytical methods for causality evaluation are presented in Section 4, whereas a set of examples is presented in Section 5. In this paper, the most important models of photonic materials are analyzed in terms of causality with the use of precise mathematical tools. Furthermore, the set of various distributional theorems presented in literature is collected as the distributional version of the Titchmarsh theorem, allowing us to evaluate causality of complicated electromagnetic systems on a mathematically rigorous basis. Hopefully, researchers interested in evaluating causality of novel models can use our paper as a template for their own mathematical derivations and proofs.

2. Review of Background Mathematics

2.1. Basic Notations

In this paper, a standard engineering notation is employed, which denotes an imaginary unit as . For the complex number , we denote its real part as , and its imaginary part as . The (open) right half-plane is denoted as . We denote the space of all functions which are holomorphic in the right half-plane and bounded by a polynomial as . To be more precise: means that the function G is holomorphic in the right half-plane, and there exist and a constant such that for all for all satisfying .

Then, we define the Fourier transformation for the absolutely integrable function , using the formulation applied in electrical sciences whereas the inverse Fourier transformation is given by

In order to change the settings of the Fourier transformation to the one widely applied in mathematics and physics, one should replace the imaginary unit j with in (2) and (3) (where ) [8]. Various literature sources which we refer to sometimes apply the mathematical i-convention; hence we convert those results to the engineering j-convention.

Among various elementary functions, we often refer to the Heaviside step function , defined as for and for . Then we refer to the signum function defined as for and for . One should be aware that these functions appear as regular tempered distributions and, as such, should be identified up to the Lebesque measure zero sets. Hence, in both cases, it actually does not matter how the functions are defined for .

In this paper, we often refer to the concept of Hölder continuous functions. We call the function a Hölder continuous one (with an exponent ), if there exists such a constant that for all . If the condition (4) is satisfied on every interval (with a constant possibly depending on M), we can say that the function is locally Hölder continuous. We should mention here that the case corresponds to Lipschitz continuous functions. One can also notice that the class of locally Hölder continuous functions covers all the functions with a continuous derivative, but it is essentially larger (with the function as an example of Hölder continuous function with the exponent ).

In this paper, the Hilbert transformation [5] of the function is applied. It is defined as

The classical domain of this definition is the space for any . The definition (5) can be formulated for a wide class of distributions, which is discussed below.

2.2. Fractional Calculus

In modeling of dielectrics, one can find FO derivatives, hence the notation is fixed below. Various approaches to fractional calculus are considered in literature, i.e., Riemann-Liouville, Caputo, Grünwald–Letnikov and Marchaud, to mention just a few. For appropriate definitions, refer to classical monographs [30,31,32]. In this paper, when the fractional derivative of order is applied, its Marchaud definition is used (refer to [31] Sections 5.4 and 5.5) for ( and ). In (6), the function f is assumed to be sufficiently smooth, e.g., with bounded by a function which is not growing too quickly in .

The main reasons behind the usage of the Marchaud derivative in our considerations are as follows [33,34]:

The Marchaud derivative of the order for the function f exists, if is bounded and locally Hölder with an appropriate exponent.

For the Marchaud derivative of the order of the exponential function (where is fixed), one obtains

The Marchaud and Grünwald–Letnikov derivatives coincide for a very broad class of functions.

The Marchaud derivative satisfies the semigroup property for all the f functions for which this definition coincides with the Grünwald–Letnikov definition, i.e., where .

It is demonstrated in [33,34] that, in order to obtain the equivalence between the results in the time and frequency domains, the FO derivative modeling electromagnetic systems should be representable in the phasor domain (i.e., satisfy (7)) and satisfy the semigroup property. From this point of view, we considered in [33,34] the following definitions of FO derivatives applied for the electromagnetic modeling: Riemann–Liouville, Caputo, Liouville–Caputo, Liouville, Marchaud, Grünwald–Letnikov, Caputo–Fabrizio, Atangana–Baleanu, Atangana–Koca–Caputo and the conformable derivative. Out of these most popular approaches, only the Grünwald–Letnikov and Marchaud definitions (which are actually equivalent for a wide class of functions) satisfy the semigroup property and are naturally representable in the phasor domain. Therefore, we employ the Marchaud derivative in our research focusing on causality evaluation of photonic materials.

2.3. Distribution Theory

The distribution theory is applied in our investigations; therefore the notation is fixed below. It is a very formal mathematical theory, which is widely applied in many branches of science and engineering. However, the very formal formulation of this theory, e.g., proposed in purely mathematical books [35,36], is not straightforwardly applicable in applied physics and engineering. Even though alternative attitudes allow for different views on distributions, they are based on the same foundations. Therefore, we refer a reader to literature sources [37,38,39], which provide different perspectives on the distribution theory.

The support of the continuous function is the set [37]

From now on, denotes the space of test functions of the class with compact support, endowed with appropriate topology (i.e., with the formal definition of all convergent sequences within the space of test functions). The topology is given by an appropriate family of seminorms, as described in Section 1.2.6 of [38]. The space dual to , i.e., the space of distributions, is denoted as . The linear continuous functional f on is denoted using the dual-pair notation , where and . The space of Schwartz functions, i.e., rapidly decreasing functions, is denoted as . Its dual space, i.e., the space of tempered distributions, is denoted as . The Fourier transformation is defined for all tempered distributions by the formula [38] Section 3.1.4 for all . The support of the distribution f is defined as the set being a complement of the largest open set U, on which f vanishes [38] Section 1.3.1

Let denote the space of all functions. Then its dual is the space of distributions with compact support.

We should also refer to some spaces of test functions (and related spaces of distributions), which are useful when discussing the Hilbert transformation. Let us take and such that . Then, the space of test functions consists of all such that all the derivatives for all with the topology defined by the family of seminorms inherited from Sobolev spaces (for details, one is referred to (1.57) in [35]). Let denote the space dual to . The special case is with being dual to . The derivative of the order k of the function (or distribution) f is denoted as . For the space , one can formulate the following theorem:

([35] Theorem 1.26). The distribution u belongs to

One can find a detailed discussion of properties of the spaces in Section 10.2 in [5].

The relation between different subspaces of the space of distribution can be summarized as (see [5] Section 10.2) where .

Let us define the distribution , which is needed in the sequel, as follows:

As shown in [5] Section 10.7, the distribution belongs to for all .

In [35] Lemma 1.8, the space is defined as the set of all tempered distributions such that . Following this definition, is defined as the space of such tempered distributions that , and the Fourier transformation is a one-to-one mapping between and .

If one considers two distributions , then its convolution is not always well-defined (refer to the discussion in Section 10.6 of King’s book [5]). However, when it is possible to define the convolution of distributions f and g, it is defined as the distribution h such that for the test function . There are some cases, when the convolution of distributions is well-defined (refer to the end of Sections 10.6 and 10.7 in [5]), that is:

if and , then

if and , then

if and , then

if , then exists

if , and then exists, and , where .

The Hilbert transformation of the distribution is defined as the distribution (see [5] Equation (10.83)) satisfying the formula where . For the function , the Hilbert transformation can be written as the convolution which is defined in the distributional sense (see [5] Equation (10.102) and the following discussion ending at Equation (10.121)).

We should bear in mind some important properties of the distributional Hilbert transformation. That is, the following properties are valid for :

For any distribution , one obtains (cf. Section 10.9 in [5])

Let be a function. Then one obtains (see (4.30) in [5])

For any , the Hilbert transformation commutes with the distributional derivative (see (10.173) and the entire Section 10.10 in [5])

The Hilbert transformation of the Dirac delta is given by (see (10.85) in [5])

The distribution belongs to , and one obtains (see (10.87) in [5])

A similar result to the one given above can be stated for the translated distribution . That is, one obtains

Although the above property is quite obvious, we were not able to find it given directly in literature. Therefore, a short proof of this property is given. Let us take the test function and check (we use the definition (14), the variable change property (17) for the transformation of function, and additionally the change of variables ) which completes the proof of (21).

We also need to remember topology in the space of tempered distributions: One says that the sequence of tempered distributions converges to a tempered distribution if for any Schwartz function [37].

3. Causality in Electrodynamics and Optics

Let us consider Maxwell’s equations in the time domain where and denote, respectively, the electric- and magnetic-field intensities, and denote, respectively, the electric- and magnetic-flux densities, and denote, respectively, the current and charge densities. In this paper, we focus on dielectric properties of isotropic electromagnetic media, i.e., photonic materials, whose properties are described in the frequency domain. It stems from the fact that the number of models of dielectric properties of media is much larger than for magnetic properties. For such mathematical models, it is crucial to approximate physical characteristics by using causal formulas. However, the presented methods and tools can be extended towards magnetic characteristics. Furthermore, the obtained solutions of Maxwell’s equations should also be causal. Finally, having the frequency response of media in the wave-propagation (or wave-guiding) problem, one can evaluate causality using the theorems presented below, supported by methods of their usage.

3.1. Basic Definitions

The function (generally ) or the distribution is called causal if its support . The Fourier transform is called a causal transform if . In other words, one can assume that is causal if it can be represented as , where is the Heaviside step function and for [40]. In practical terms, one should also assume that is the function whose Laplace transform has a non-degenerate region of convergence.

In terms of physics, causality means that the effect does not precede the cause. Hence, if one considers the electromagnetic system whose time-domain function describes the system response to the Dirac delta excitation, then causality means that for . It also means that response of the system depends only on excitation values from the past. To sum up, the mathematical definitions of causality formulated above closely follow physical understanding of this term.

3.2. Dielectric Models

Transformation of (23)–(26) into the frequency domain gives where tilde denotes phasor representation of the physical quantity, i.e., , and denotes the angular frequency. To solve this set of equations, one additionally needs constitutive relations between and as well as between and . Hence, one can write

Then one can also write that where and are, respectively, the electric and magnetic polarizations, and are, respectively, the vacuum permittivity and permeability, and are, respectively, the electric and magnetic susceptibilities of the medium. As it has already been mentioned, our considerations are focused on isotropic dielectric media. Therefore, below we consider general functions, which are mainly related to dielectrics, but can also represent, in an obvious way, magnetic models of media (i.e., magnetic susceptibility).

In electrodynamics, it is required that is a causal transform [7], because all dielectrics are not able to polarize instantaneously in response to an applied field. Alternatively, one can require that the function is a causal transform.

Let us formulate the first considered dielectric model, in which the transform is a real constant. That is where is the constant relative permittivity. Then let us assume the ohmic conduction, i.e., the current flow whose density depends on the electric-field intensity as follows:

In (36), denotes the electrical conductivity. Therefore, assuming that , (30) can be written as where includes ohmic losses within the complex permittivity . This result can be generalized towards where are functions of frequency and describes all losses in the electric field. If the losses in a dielectric material stem from ohmic conduction, then is proportional to . Otherwise, if losses stem from the bound charge and dipole relaxation phenomena, then is not proportional to . Analogously, one can write for the permeability where are functions of frequency.

Let us formulate a few popular models describing the response of dielectrics with the use of permittivity:

Djordjevic-Sarkar relationship for lossy dielectrics [6,41] where . The first term in (41) is the relative permittivity at very high frequencies, the second term is the broadband logarithmic term, and the third term comes from conductivity. In (41), , , are model parameters. This formula gives a simple closed-form expression, which approximates the measured permittivity of the popular FR-4 substrate used for manufacturing printed circuit boards.

Westerlund relationship for FO capacitors [42,43] where . It allows for formulating the constitutive relation (31) in the time domain with the use of FO derivative as

Although this model does not explain the nature of internal processes in dielectrics, it reproduces and predicts their behavior much better than any other theory (according to the Authors of [42]). Therefore, it is referred to as an ‘engineering’ model of dielectrics. Furthermore, this model allows for obtaining the electrical characteristics of FO capacitors, (refer to [43,44]).

Power-law relationship for porous media [45] where , , A is a constant, and is close to 1.0 in a low frequency region, and is within the range of 0–0.5 in a high frequency region. The model (44) describes the permittivity of porous media such as wet soils and sedimentary rocks, which has been observed to be considerably different than in the case of water and parent minerals.

3.3. Dielectric Relaxation

Let us formulate several popular dielectric models based on electric susceptibility in the complex domain ( is the relaxation time, and where ):

Debye [46,47]

This model is frequently used to describe simple dielectric characteristics of electromagnetic media arising from bipolar relaxation. The Formula (45) is characterized by a single relaxation time, which is capable of handling materials with high-water content. However, experimental studies show that the relaxation behavior of a wide range of dielectrics strongly differs from the Debye relaxation formula. Therefore, a number of phenomena such as broadness, asymmetry and excess in dielectric dispersion has motivated the development of empirical response functions described below, such as Cole-Cole, Cole-Davidson, Havriliak-Negami, and Raicu [46].

Lorentz [7,47] where is the frequency of a pole pair (the undamped resonant frequency of the medium), and is the damping coefficient. The model is based on the classical theory of light-matter interaction, and describes the frequency-dependant polarization due to bound charges. That is, bindings between electrons and nucleus in atom are treated similarly to those of the mass-spring harmonic-oscillator system. It is worth mentioning that any function obeying the K–K relations can be approximated as a superposition of Lorentzian functions, to any precision [48]. Therefore, the Lorentzian function (46) can be considered as a general building block for implementing causal susceptibilities of various materials, e.g., metamaterials.

Lorentz in high-frequency limit [7] where and is the plasma frequency of medium. This model results from (46), assuming that the frequency is far above the highest resonant frequency in the medium.

Lorentz in high-frequency limit with static magnetic induction [7] where and . In (48), is the frequency of precession of a charged particle in magnetic field. This model is the extension of (47), which involves an interaction between static magnetic field and tenuous electronic plasma of uniform density, when transverse waves propagate parallel to the direction of .

Lorentz in FO generalization [49] where is the termed bulk plasma frequency associated with electrons, and is the damping coefficient. This model extends the classical Lorentz model (46) for a dielectric material with the use of FO derivatives, but it is formulated in the frequency domain.

Drude [47] where and . In (50), is the Drude pole frequency and is the inverse of the pole relaxation time. This model results from the application of kinetic theory to electrons in solids for optical frequencies. It can be obtained from the aforementioned Lorenz model (i.e., harmonic-oscillator model) when the restoration force is removed (i.e., free electrons are assumed which are not bound to a particular nucleus).

Cole-Cole [46,50,51,52]

This model has been developed as an empirical extension of the Debye model (45), which can be obtained for .

Cole-Davidson [46,53]

This model has been developed as an empirical extension of the Debye model (45), which can be obtained for .

Havriliak-Negami [46,54,55]

This model extends Cole-Cole (51) and Cole-Davidson (52) models.

Raicu [46,56] where is the relative dielectric increment in the Raicu model. This model extends (53) by including the additional parameter .

Universal dielectric response [57,58,59] where and is a positive constant. In general, this model is valid for , where is the loss-peak frequency. It describes the observed behavior of dielectric properties demonstrated by solid-state systems. That is, it involves power-law scaling of dielectric properties with frequency, which is widely observed in nature.

A particular dielectric model can be a weighted sum of the several susceptibility characteristics presented above. For instance, the characteristic can be a weighted sum of several Lorentz functions (46), defined for different pole pairs and damping coefficients [7]. In such a case, if each of the components is causal, then the considered model is causal as well.

3.4. Frequency Response

If one obtains a causal solution to Maxwell’s Equations (23)–(26), then each of the components of the vector field , , , can depend only on the previous excitation values , . For instance, the relation between the excitation and the response of electromagnetic system (, ) in the frequency domain can be written as follows [8,60]:

In (56) and (57), , denote dyadic Green’s functions of electric-electric and magnetic-electric type, respectively. Because the considered electromagnetic-field systems are linear and time-invariant, it is possible to write the relation between a single excitation component () and a single output component () in the frequency domain. Then, one obtains where denotes frequency response of the system.

For instance, let us consider one-dimensional (1-D) propagation of a monochromatic plane wave along the z direction in a medium described by material parameters and . Then one can write the following Helmholtz equation for and :

In (59), is the square of complex-valued wavenumber. In optics, the refractive index is mainly used to describe an electromagnetic medium, which is a dimensionless number describing how fast the light travels through the medium. That is, , where is the velocity of light in the medium and c is the velocity of light in the vacuum. Hence, one can also write that . If the refractive index is a complex number for a given angular frequency, its real part indicates the phase velocity, whereas its imaginary part describes the attenuation of electromagnetic waves in the medium. Let us consider the signalling problem [61,62,63], where the electric field in the homogeneous Equation (59) is excited by a source at the spatial-domain boundary. Hence one obtains whose general solution is given by in (61), and denote, respectively, complex amplitudes of waves propagating in the and directions, and and are complex roots of . Considering wave propagation (or guiding) in the direction only, the propagation constant depends on the choice of and functions, and is selected as the one with a positive real part. Hence one obtains where is equal to . Such a solution is physically equivalent to impinging of the plane wave on the half-space constituting a medium described by and . Then the wave is transferred into the medium and its time-domain waveform can be obtained as described in [25]. Alternatively, one can consider (62) as a general solution of the wave-guiding problem for, e.g., optical waveguide. Assuming the fixed length of the wave-propagation distance , and taking and , the Formula (62) can be considered as the relation (58) where

Such a function is usually required to be a causal transform in electrodynamics. Furthermore, one can require that is relativistically causal. That is, not only the inverse Fourier transform of is equal to zero for , but the inverse Fourier transform of is also equal to zero for .

Let us consider FO models in electrodynamics, which start with constitutive relations as follows [43,44,63]:

Equation (43) is repeated here as (65) for the sake of completeness. When these relations are applied to Maxwell’s Equations (23)–(26), then one obtains the FO Maxwell’s equations in the following form:

In order to analyze the 1-D propagation of a monochromatic plane wave along the z-direction, one can apply the phasor-domain representation and arrive at the special version of the Helmholtz Equation (59) with

With the additional assumption of the lack of current and charge sources within the considered space, and assuming that there is no power dissipation due to Joule’s heating, i.e., the current density is related to the electric-field intensity by the classical Ohm law () with the conductivity , the frequency response (i.e., the transfer function in the frequency domain) is given by (cf. [63]) where .

For , and , one can write (67)–(70) in a compact form where is the Riemann-Silberstein (RS) vector in the time-fractional electrodynamics [64] and . Then, one can write the diffusion-wave equation in time domain for space without sources

Assuming the plane-wave, spherical and cylindrical symmetries of solutions to (77), one obtains, respectively, the following transfer functions describing the 1-D wave propagation [64]:

In (80), the function is the Bessel function of the first kind of zero order, i.e.,

4. Methods and Analytical Tools for Causality Evaluation

Let us consider a complex-valued transfer function or distribution, which is the Fourier transform of a certain time-domain function or distribution. It is worth noticing that we do not assume a priori that is real-valued.

Approaches to causality evaluation are presented in Figure 1. Having , one can evaluate causality by way of applying the Paley–Wiener theorem, calculating the inverse Fourier transformation, finding a holomorphic extension to the right half-plane, or checking various forms of the K–K relations. One should notice that, for the considered function , not every approach can easily be applied to prove either causality or lack thereof.

Figure 1

Approaches to causality evaluation.

4.1. Paley–Wiener Theorem

Let us start our considerations from the Paley–Wiener theorem, which allows for characterisation of the modulus of the complex-valued function in terms of causality.

(Paley–Wiener, [65] Theorem XII). Let

One should notice that the Paley–Wiener theorem does not state that the complex-valued function is a causal transform. It states that, for the modulus satisfying (82), the causal transform exists with the same modulus. It also states that if does not satisfy (82), then is surely not a causal transform. This theorem is a valuable tool, which can be used to prove that the transfer function is not a causal transform.

4.2. Calculation of Inverse Fourier Transformation

The simplest approach to causality evaluation relies on calculating the inverse Fourier transformation of the function . Then is not causal if its support is not contained in (if is a continuous function, then it is enough to show that it is not equal to zero for a certain ). Alternatively, is causal. The method seems to be very simple, but it may be really difficult to calculate the inverse Fourier transformation. In some cases, one knows the exact formula for the (inverse) Fourier transform for a given function or distribution. On the other hand, in numerous cases the exact formula is unknown.

In some cases, one can easily prove lack of causality by referring to the properties of the inverse Fourier transform. For instance, if the function is an function, then the inverse Fourier transform is a continuous function. Hence, it is enough to find a single point such that , e.g., to show that . This idea may not be applied directly when is an function or a distribution which is not represented by an function. With the definition of the Fourier transformation extended to the above-mentioned domains, one identifies the result up the sets of measure zero. In this case, without continuity of the result, showing that the inverse Fourier transform is non-zero at a single point does not prove lack of causality.

4.3. Holomorphic Extensions and K–K Relations

The classical perspective on causality is provided by the Titchmarsh theorem, which works for functions (see Theorem 1.6.1 in [4]):

If a square-integrable function

The inverse Fourier transform

One should notice that the relations (83) and (84) hold in the sense of elements of the space, i.e., the equalities hold for almost all .

The relations (83) and (84) are also referred to as the K–K relations or the dispersion relations. This theorem delivers two aforementioned approaches to prove causality, i.e, searching for an appropriate holomorphic extension of to the right-half plane, and proving the validity of the K–K relations (83) and (84).

If the function is the Fourier transform of the real-valued function , then it is hermitian, i.e., it has an even real part and an odd imaginary part. Therefore, the K–K relations (83) and (84) can be formulated for almost all by the following integrals on :

It should be mentioned that the growth assumption in the case (ii) of Theorem 3 (i.e., that the maps belong to for all , and that all the norms are uniformly bounded by the same constant ) is vitally important. The sole existence of a holomorphic extension of the function may not be sufficient for its causality. The case of is a good example. This function naturally extends to the holomorphic function , where , while the inverse Fourier transform is not a causal function. The uniform boundedness of the norm is the violated condition. It is because , thus, for , the norm can be arbitrarily large.

The general distributional version of the K–K relations, given as Theorem 5 in [25] (following Theorem 3.10 in [35]), is presented below. This version is a generalization of the well-known K–K relations with subtractions (as described in Section 1.7 of [4]). The procedure of subtractions works for such a Fourier transform , which is not necessarily in . However, when divided by some polynomial of degree k, it belongs to . The generalization towards distributions is sometimes required because the division itself can introduce a singularity, resulting in a function which is not locally square integrable (see the discussion in Section 4.2 in [25]). In the distributional version, it is not required that the division result, i.e., , is in , but it can be a distribution belonging to a class for some , which can be broader than the functional space (see Example 1 in [25]).

Let us assume that or

One should notice that (88) can be written with reference to the distributional version of the Hilbert transformation as

The other important detail hidden in this theorem is related to the growth condition on the holomorphic extension. One should remember that the condition means that is holomorphic in , and that its growth is controlled by some polynomial as the condition (1) states.

A slightly less general version of Theorem 4, which can be used in the case of distribution F being represented as a locally integrable function, is formulated below.

(Theorem 6, [25]). Let us assume that for where

In some cases, the K–K relations can be evaluated for the transfer function logarithm. This attitude can result in sufficient conditions for causality.

(Theorem 8, [25]). Let us assume that the function

the function

Then

The ‘holomorphic-extension’ approach can be generalized towards functions with a polynomial growth. The first theorem can be found in [35] as Theorem 2.7.

(Theorem 2.7, [35]). If

The theorem formulated above gives rise to practical sufficient conditions for causality of the transform .

Let us assume that there exists a function

f is locally integrable, with a growth in

the functions

Then there exists

This is a direct consequence of Theorem 7. Let us take any and notice that

The functions are definitely continuous (as sections of holomorphic function), and hence locally integrable. Moreover, each of the integrals and exists. Because for any and almost all , as well as the function is (pointwise, almost everywhere) convergent to the function , we can refer to the Dominated Convergence Theorem (see, e.g., (2.206) in [5]) and state that

It means that the tempered distributions converge to G in topology. Thanks to Theorem 7, we know that is a causal distribution. □

If the condition (ii) of the above theorem holds for

Furthermore, a more general version (with growth restrictions which are not necessarily polynomial) can be formulated as follows:

(see [66] Theorem 3.8). Let us suppose that

for each

there exists

Then there exists such a distribution

The set of distributional theorems presented in this section can be collected as a distributional version of the Titchmarsh theorem. It provides the conditions for the tempered distribution to be supported in , due to the properties of its Fourier transform . Let us mention that each of the tempered distributions can be represented as , where (see the discussion following Definition 3.2 in [35]).

Let us assume that the tempered distribution

for the distribution

the distribution F is the boundary value in the

the following relation is satisfied for a certain polynomial

The equivalence is stated in Theorem 7 above (i.e., Theorem 2.7. in [35]). As one can see, (iii) is equivalent to (87), which is equivalent to (i) by Theorem 4. □

As it has been mentioned in

Then one can calculate the integral for

(Dielectric model with constant permittivity (35)). This model implies that electric susceptibility is a causal distribution.

(Dielectric model with ohmic losses (38)). This model is not valid for

This directly implies that is a causal distribution.

(Debye relaxation model (45)). For this relaxation model, one can calculate the inverse Fourier transformation as follows [

For t < 0,

(Lorentz model (46)). For this model, one can calculate the inverse Fourier transformation [ for

For t < 0,

(Lorentz in high-frequency limit-model (47)). This model requires some explanation, which is important from the formal mathematical perspective. The function principal part formalism (see, e.g., the discussion in Section 10.1 in [

It might be shown that (in the distributional sense) is given by

For t < 0, one obtains

Hence the model (

One is also referred to Example 20 given below, where causality of this model is discussed from the perspective of the K–K relations.

(Lorentz in high-frequency limit with static magnetic induction-model (48)). One can notice that the Formula ( where

Then we use the following properties of the Fourier transformation:

Hence one obtains

Finally, the inverse Fourier transformation of (

However, (

Then one obtains in the time domain

Although (

One is also referred to Example 21 given below, where causality of this model is discussed from the perspective of the K–K relations.

(Drude model (50)). This formula is undefined for

For such a distributional extension, one can find the inverse Fourier transform (well known in literature [

Hence, for

(Djordjevic-Sarkar relationship for lossy dielectrics (41)). This model is undefined for

Using the time-domain representation of the function (

It is a causal distribution.

(Cole-Cole model (51)). For

Then one can calculate

Let us denote one can notice that, for

One can notice that for all

This expression can easily be estimated when

The last estimate shows that all the functions

For where for a certain positive constant

Let us notice that

One can also notice that by a constant. It proves that

(Cole-Davidson model (52)). One can easily notice that, for

In the case of

Therefore one can state that

(Havriliak-Negami model (53)). This model can naturally be extended towards the function because

(Raicu model (54)). This model has natural holomorphic extension to

It leads to where (one can take the exponent

Hence it is estimated by a constant in any half-plane

Let us first observe that, if

Let us now assume that at least one of the constants

Hence, all the assumptions of Theorem 8 are satisfied and the transform

(Lorentz model in FO generalisation (49)). This model has a natural holomorphic extension given by

The open question is whether this extension is well-defined over with an unknown

Hence, no matter if

Hence the extension ( for any

5. Causality Evaluations

5.3. Applications of Holomorphic Extensions

Now, let us take a look at the transfer functions which are derived for the formulation of time-fractional electrodynamics based on the RS vector.

(Plane-wave and spherical symmetries of solutions to diffusion-wave equation formulated based on RS vector-models (78) and (79)). For the functions

The assumptions of the Titchmarsh Theorem 3 are satisfied in the point (ii) for both functions; hence it directly proves their causality.

(Cylindrical symmetry of solution to diffusion-wave equation formulated based on RS vector-model (80)). The function valid for

The behavior in infinity described by (

As one can see, due to the estimate ( and because of the trivial estimate one can focus on the estimate of

On the other hand, by the Hankel’s Asymptotic Expansions (see [

(Westerlund relationship for FO capacitors [42,43] and power-law relationship for porous media [45]-models (42) and (44)). Both Equations ( for certain constants is a causal transform (subtracting 1 does not influence causality, since it is a transform of the Dirac delta, i.e., a causal distribution). One can see that

Let us refer to Theorem 5 for

Let us calculate the Hilbert transforms and

Let us refer to a certain integral formula taken from (3.241.3) in [

Taking

Similarly, taking

Now, we take the Formula (

It means that and the real parts on both sides of (

It means that which proves that the imaginary parts of (

Let

Let us observe that, by (

Hence the relation ( one can see that (referring to (

However, it is well-known that assigning an integrable function its primitive. It is natural to consider this operator as a causal operator, as it converts causal integrable functions to causal locally integrable functions. The problem is that multiplication by the function where

However, this formula does not work for all the tempered distributions. For instance, let us take and

Let us fix for

The left-hand side of the above inequality is obviously

Let us now calculate the convolution for the test function

It means that

Hence, both sides of (

5.4. Applications of K–K Relations

Having reviewed two abstract cases, let us look at some of the models from the perspective of the K–K relations.

(Lorentz in high-frequency limit-model (47)). Because we treat the function

It means that

In these derivations, the Formulas ( is a causal transform.

(Lorentz in high-frequency limit with static magnetic induction–model (48)). As before, one can see that the Formula (

Let us check the K–K relations for this function. Since

It implies that results in a causal transform.

5.5. How to Prove Lack of Causality?

Obviously, one can prove that an appropriate holomorphic extension of the transform does not exist, but it does not look like an easy task. One should return to the K–K relations instead, and prove that one of the equalities (83) or (84) (or in the case of hermitian transforms (85) or (86)) is not satisfied in . Because the equalities are given between the elements of space, it means that violation of any of the conditions (83) and (84) for the transform in a single point does not prove that the transform is not causal in general. Equations (83) and (84) are in sense; hence such equalities are valid almost everywhere. Fortunately, in certain cases, it occurs that the relations (83) and (84) are valid for all .

We now provide the theorem stated by Wood in 1929 ([71] Theorem I, see also [5] Section 3.4.1):

(Wood, [71] Theorem I). Assuming that

the integrals

f is locally Hölder continuous of the order

Then the function g is also locally Hölder continuous with an exponent α and holds for all

Hence, for the Hölder continuous functions or , with a behavior in as required by Theorem 11, violation of any of the relations (83) or (84) in a single point proves that they do not hold in . This attitude is used in the following example:

In [

Hence one obtains

As it is shown in Lemma 2 in [ and

Because both functions are locally Hölder continuous, the K–K relations hold true in

In the case of distributions, Theorem 11 can also be helpful. Let us first consider the following example:

Let us return to the case of where the sign depends on the parity of belongs to and

In order to state that

As one can see, the function for all

This example can be a good starting point for some general observations in the context of Theorem 4, providing an easily verifiable necessary condition for causality.

Let us assume that

Then, if the real (respectively imaginary) part of G is a locally Hölder continuous

5.6. Causality Tests for Refractive Index

One can formulate the causality problem for media described by the refractive index . From the point of view of Maxwell’s equations, the velocity of electromagnetic wave is described by and . Hence one obtains that , refer to (59). Let us assume that the dielectric-relaxation function satisfies the assumptions of the Titchmarsh Theorem 3. Therefore, among others, it belongs to and is causal. Then has a holomorphic extension to the right half-plane, and the permittivity also has a holomorphic extension to the complex right half-plane (the same considerations are applicable to the permeability ). One should notice that this condition formulated for settings implies that the considered function is holomorphic in the upper half-plane. We should mention that the existence of holomorphic extension is not a sufficient condition for causality, i.e., the behavior in ∞ is important as well. In general, assuming that and belong to does not imply that the product . Hence, one has to be careful when deciding whether is an function. Similarly, if the holomorphic extensions and satisfy the assumptions (ii) of Theorem 3, it does not mean that these assumptions are satisfied by the product . This means that we may not form conclusions about causality based only on the existence of holomorphic extension. Some assumptions concerning the behavior of the product for a fixed and are needed as well. As proposed by Stockman [27], the K–K relations can be written for the complex refractive index as

Although the assumption does not imply that exists, in practical terms it is natural to assume that . Due to this assumption regarding the electric susceptibility , as well as the magnetic one, one obtains that when [74]. The usual physical justification behind this assumption is that an incident, oscillating, electromagnetic field entering any medium stimulates the charges in that medium to oscillate (light-matter interaction). However, for very high frequencies of the incident field () the charges of the medium cannot respond, because they have a finite mass, hence their inertia. As a result, for those very high frequencies, it is as if the field ‘sees’ a vacuum (), because it effectively does not interact with the medium at all.

From this additional assumption (i.e., , ), as well as the assumption that both functions and are bounded (it happens, e.g., when they are continuous), one can conclude that the product belongs to . In this case, one knows that belongs to . It allows one to say that the relations (183) and (184) imply causality of the transform by the classical Titchmarsh Theorem 3. From the formal perspective, the assumptions related to limits of the functions , as and their boundedness can be relaxed (e.g., by saying that one of these functions is any function, and the other one is a bounded function), but it is far less natural and does not allow to conclude that when .

A separate discussion is needed when the susceptibilities and are not functions, e.g., when they are tempered distributions. In this case, one has to be careful when defining the product , as the product of two distributions is not necessarily well-defined. One can surely define as a tempered distribution when both and are represented by locally integrable functions, which are bounded as by some polynomial. In this case, the K–K relations can be checked in the distributional sense with the use of Theorem 4, but it requires dividing the function or by an appropriate power of . Let us also notice that the K–K relations for and can be different. In the case considered above, when belongs to , there is no need to divide by any positive power of (i.e., one may take in the Formula (87) or (88)). On the other hand, is not the function in , and in order to verify any of the Equations (87) or (88), one should take .

In literature one can also find the K–K relations formulated for

However, as it can be noticed in [27], the function can be not holomorphic in the right half-plane , even if is holomorphic. For instance, when approaches zero, then the derivative of its square root does not exist. It clearly shows that causality of is not necessarily equivalent to the case of causality. Apart from the problem with being holomorphic in , there is also the problem of its behavior as . If , then it does not necessarily mean that . Moreover, if , then it is not necessarily true that . In order to draw any causality conclusions from the relations (185) and (186) formulated for , one should know that belongs to . Despite these issues, the K–K relation of the type (185) is successfully used for the interpretation of experimental results in [75]. That is, the effect of a femtosecond-laser-induced electronic band-gap shift on the refractive index is explicitly studied with the use of the K–K relations. Clearly, from this relation, a change in the absorption described by curve in turn affects . It is worth mentioning that the K–K relations (185) and (186) formulated for are valid for a single mode propagation in waveguides (e.g., optical), for which one can write where .

Let us now assume that is causal (e.g., it is an function which satisfies (185) and (186)). Then, if is such a distribution that the convolution exists, and if it is a tempered distribution (the convolution of two causal distributions is causal as well), then is also causal. However, the opposite theorem is false in general. Nevertheless, if one knows that the function is holomorphic in and it does not achieve 0, then is holomorphic as well.

Still, one can prove that the passivity of media implies that is holomorphic in the right half-plane [76]. The problem of choice between and for applications in optics is debated in [77]. That is, the consideration of propagation of optical pulses with the use of complex index of refraction is inconvenient in general. Therefore, when calculating the wave vector, one can take either or as a velocity of wave propagation. Whereas the difference between the two is negligible for small losses, it is significant in other cases. The analysis of pulse propagation demonstrates that the use of results in a wave vector different than that actually exhibited by the propagating pulse. On the other hand, the definition always correctly calculates the wave vector of pulse, hence it is preferred in optical investigations. Moreover, for negative refraction media, when the sign of permittivity or permeability changes as a frequency function, one should notice that the derivative of does not exist when approaches zero. One should take all these issues into consideration when either the K–K relations (183) and (184) or (185) and (186) are applied.

It is worth noticing that the K–K relations are also applicable—to a certain extent—in nonlinear optics. That is, the response function in the time domain should also be causal for nonlinear media. In general, nonlinear complex susceptibilities can have poles not only on a half of the complex-frequency plane; however, there are cases when it happens. In such a case, when holomorphic properties are available on a half of the complex-frequency plane, the standard K–K relations (83) and (84) can be useful [78]. The review of nonlinear K–K relations in optics and photonics is presented in [79], where it is shown that the nonlinear dispersion relations have a common form that can be understood in terms of the linear K–K relations (83) and (84) applied to a new electromagnetic system consisting of the material and the perturbation of its parameters. As noticed in [80], the nonlinear K–K relations are useful in optics and photonics to predict that an enhancement in the nonlinear optical absorption for a specific wavelength usually leads to a decrease in the nonlinear optical refraction associated with a considered material.

Equation (183) can be used to derive the condition of negative refraction with no (or low) loss at the observation frequency [27], i.e.,

It is derived assuming that, at and near the observation frequency , the material is transparent (e.g., the losses are compensated by gain), which mathematically implies that and . Furthermore, it is assumed that the negative refractive index is implied by the condition , where and denote, respectively, the phase and group velocity. Due to these limitations, the condition (187) is replaced in [28] by the condition which does not require that and at the observation frequency. The condition (187), obtained from the K–K relations, implies that compensation of the optical losses or significant reduction, by any means (material or structural) of the imaginary part of permittivity and permeability, can also change the real parts of these quantities in such a way that the negative refraction disappears [27]. Concerning (183) and (187), care should be exercised for the case when both the real part of permittivity and the real part of permeability are negative–simultaneously and within the same frequency region, because in that case, although the product is positive, one should nonetheless select the root with negative sign of the real part of and ensure that, in the absence of a gain mechanism the medium does remain passive [81]. Further, one can notice that (187) stipulates that it is impossible to have a loss-free or amplifying (, ) medium with a negative refractive index (, , ) for all the frequencies–or else (187) would not hold true. However, (187) does not preclude the possibility that such a medium exists within a finite-bandwidth frequency region, outside which the product in the numerator in the integral of (187) could make a sufficiently ‘negative’ contribution, so that the overall (187) still holds true. Indeed, such lossless negative-refractive-index media have been reported in the past, both experimentally [82] and numerically [83,84,85].

5.7. Summary

The described analytical methods for causality evaluation of dielectric models of photonic materials are summarised in Table 1.

Table 1

Summary of causality evaluation methods for photonic materials (IFT—inverse Fourier transformation, KKR—K–K relations, HE—holomorfic extension, PWT—Paley-Wiener theorem).

Model	Equation	Method	Example
dielectric with constant permittivity	(35)	IFT	2
dielectric with ohmic losses	(38)	IFT	3
Djordjevic-Sarkar for lossy dielectric [6,41]	(41)	IFT	9
Westerlund [42,43]	(42)	KKR	17
power-law for porous media [45]	(44)	KKR	17
generalized power-law [25]	(65) and (66)	PWT	1
Debye relaxation [46,47]	(45)	IFT	4
Lorentz [7,47]	(46)	IFT	5
Lorentz in high-frequency limit [7]	(47)	IFT, KKR	6, 20
Lorentz with static magnetic induction [7]	(48)	IFT	7
Lorentz in FO generalization [49]	(49)	HE	14
Drude [47]	(50)	IFT	8
Cole-Cole [46,50,51,52]	(51)	HE	10
Cole-Davidson [46,53]	(52)	HE	11
Havriliak-Negami [46,54,55]	(53)	HE	12
Raicu [46,56]	(54)	HE	13

6. Conclusions

In this article, a comprehensive analysis of mathematical techniques for causality evaluation of photonic materials is presented. It includes not only the approaches valid for the functions, i.e., those for which the Titchmarsh theorem can be useful, but also the functions to which the distribution theory and the FO calculus have to be applied. We present a set of theorems applicable for causality evaluations, as well as specific examples showing how to use this mathematical machinery. Furthermore, the set of various distributional theorems presented in literature is collected as the distributional version of the Titchmarsh theorem, allowing us to evaluate causality of complicated electromagnetic systems on a mathematically rigorous basis. In addition to the well-known K–K relations, we have also outlined four further methodologies, namely application of the Paley–Wiener theorem, calculation of the inverse Fourier transformation, identification of holomorphic extensions to the right half-plane, and check of the K–K relations for the natural logarithm of a system’s frequency response. The collection of these methodologies—otherwise scattered in a wide range of pertinent literature—may prove useful for scientists and engineers investigating causality problems in electrodynamics and optics.