1. INTRODUCTION

The two-part tutorial on Maxwell Garnett approximation was published in JOSA A in 2016 [1,2]. Given the vast literature on the subject and the variety of applications in which Maxwell Garnett approximation continues to play an important role, the exposition in [1,2] was limited to a subset of topics that, in the author’s opinion, were conceptually important. However, one omission was particularly unfortunate: we did not cover in any detail disordered media. Although the elementary considerations presented in Section 2 of [1] apply equally to periodic and disordered composites, there is a class of problems that requires a more nuanced approach.

Namely, consider a medium consisting of a homogeneous host and many electrically small, non-spherical, randomly oriented inclusions. The dipole polarizabilities of the inclusions depend on their shapes and are, in general, tensorial. Therefore, we expect the shape of inclusions to remain important even in the low-concentration limit. Yet, since the composite is isotropic on average, its effective permittivity must be isotropic (scalar). In the preceding parts of the tutorial, we did not provide any formulas that are applicable to this case, and we will fill this gap below.

This third part of the tutorial is organized as follows. In Section 2, we discuss the homogenization formulas that were derived by us previously and explain why these results are inadequate to the problem formulated above. In Section 3, we briefly review the relevant literature. In Section 4, we summarize all assumptions that must be made to derive the general Maxwell Garnett-type approximation. Derivation of the main result is presented in Section 5. Here we obtain a formula, which allows for many special cases and encompasses the results derived in the first to parts of the tutorial as well as the case we are specifically interested in here. In Section 6, we derive various special cases of the above general result. Here we discuss in detail two-component mixtures (media consisting of only two materials) and multicomponent mixtures in which inclusions are of the same shape and randomly oriented in space. Finally, Section 7 contains a discussion.

As in the preceding two parts of the tutorial, the Gaussian system of units is used. Tensors are denoted by an overhead hat as in $\hat \alpha$. Principal values and axes of a tensor will be labeled by the variable $p$, where $p = 1,2,3$. The unit tensor is denoted by $\hat I$. Some quantities are tensorial in general but can be scalar under some conditions. For example, the effective permittivity ${\hat \epsilon _{\texttt{MG}}}$ is in general a tensor. However, if the effective medium is isotropic, we can write ${\hat \epsilon _{\texttt{MG}}} = \hat I{\epsilon _{\texttt{MG}}}$. In such cases, we say that the effective permittivity is reduced to a scalar and write ${\epsilon _{\texttt{MG}}}$ without a hat in place of ${\hat \epsilon _{\texttt{MG}}}$.

2. SURVEY OF FORMULAS

It is instructive to briefly recount the effective medium formulas that were already derived in the previous parts of this tutorial [1,2], discuss whether any of these formulas are good candidates or can be generalized easily to capture the physical properties of random composites of non-spherical inclusions, and then preview the more general result that will be derived below.

We can start the survey with the anisotropic mixing formula for ellipsoidal inclusions {Eq. (34) in [1]}. In a slightly rearranged form, this equation reads

Alternatively, we can use the isotropic mixing formula {Eq. (18) in [1]}, viz.,

Yet another option is to start from the expression containing the dipole polarizability of the inclusions $\alpha$ explicitly. In what follows, we will frequently use the dimensionless specific polarizability $\beta$ defined as

However, direct generalization of Eq. (5) as described above is not the right approach. In what follows, we will show that the correct mixing formula applicable to the problem at hand is

In what follows, we will show that Eq. (7) is a special case of a more general formula, which is stated in Eq. (23) below. Equations (1), (3), and (5) are also special cases of this general result. Moreover, Eq. (23) encompasses many other special cases including isotropic and anisotropic media and multicomponent mixtures as is discussed in more detail in Section 6.

3. HISTORICAL OVERVIEW

Accounting for non-sphericity of inclusions in Maxwell Garnett-type effective medium theories goes back to the work of Wiener (1912) [3]. In particular, Wiener has derived a formula equivalent to Eq. (1) for two-component mixtures of similarly oriented ellipsoids, and considered two extreme cases: $\nu = 0$ and $\nu = 1$. As can be easily verified, Eq. (1) yields the arithmetic average of $ {\epsilon _{\texttt i}}$ and $ {\epsilon _{\texttt h}}$ for $\nu = 0$ and the harmonic average for $\nu = 1$. These two extreme cases are known as Wiener bounds and play an important role in the theory of electromagnetic homogenization of composites even beyond the low concentration limit. We have illustrated the Wiener bounds in the first two parts of this tutorial [1,2] with numerical examples involving Maxwell Garnett, Bruggeman, and more rigorous homogenization theories.

However, Wiener’s treatise is lengthy, considers many different subjects, and to the best of our knowledge has not been published in an English translation. In 1953, Bragg and Pippard re-derived Wiener’s result (with a reference to the original 1912 work) to explain birefringence in hemoglobin crystals [4]. Wiener’s theory was again re-derived in 1989 and extended to the case when inclusions can contain interstitial voids filled with liquid (essentially, structured or multi-layer inclusions) by Oldenbour and Ruiz [5]. Below, we will follow the derivation these references conceptually but in a substantially more general setting.

In a further development, Stroud [6] and Stroud and Pan [7] applied the mean-field approximation to derive effective parameters of composite materials. The mean-field approach is to start from the rigorous integral equations for the electric field or the potential and derive much simpler algebraic equations by factoring certain correlators. In other words, a function inside an integral transform is replaced by an a priori unknown constant (the mean field) so that the integral can be evaluated; then the constant is found from the ensuing algebraic equations. Typically, the mean-field approximation results in Bruggeman-type symmetric mixing formulas. Some of the physical effects considered in [6,7] include coupling to a static magnetic field and radiative corrections (extending the theory beyond quasistatics; see also [8 –10]).

References [6,7] were limited to spherical isotropic inclusions. However, the formalism of these references is quite general, and it was extended to describe an effect that is of interest to us (wherein a medium contains anisotropic inclusions but is isotropic on average) in 1997 by Levy and Stroud [11]. In this reference, the inclusions were still assumed to be spherical, but the material of inclusions could be anisotropic with arbitrary orientations of the principal axes. The theory was further generalized to include anisotropic host and ellipsoidal inclusions by Levy and Cherkaev in 2013 [12]. The results derived below largely overlap with those of Levy and Cherkaev, although there are some differences. In particular, we do not consider intrinsic anisotropy of the materials while Levy and Cherkaev did not consider multicomponent mixtures or non-ellipsoidal shapes.

Some recent applications of Maxwell Garnett theory to randomly distributed non-spherical particles include calculations of conductivity of mineral aggregates [13] and investigation of optical spectra of composite cholesteric elastomers doped with metallic ellipsoids [14].

4. ASSUMPTIONS

A. Linearity and Frequency Domain

If the electromagnetic response of the medium is linear, as we assume here, it is possible to consider different electromagnetic frequencies independently. We can start by focusing on a single working frequency $\omega$. Any time- and position-dependent, strictly monochromatic field can be written in the form ${\textbf{F}}({\textbf{r}},t) = {\rm{Re}}[{{\textbf{F}}_\omega}({\textbf{r}}){e^{- {\rm{i}}\omega t}}]$, where the factor ${e^{- {\rm{i}}\omega t}}$ is known as the phasor. Once we substitute this ansatz into Maxwell’s equations, the phasor cancels out, and we are left with equations containing only spatial derivatives and the (generally, complex) time-independent fields ${{\textbf{F}}_\omega}({\textbf{r}})$. Non-monochromatic fields can be introduced by using the superposition principle, that is, by simply adding monochromatic solutions at different frequencies together to satisfy some initial conditions for the field. Below we assume that all fields oscillate at a working frequency $\omega$ and say that we are working in the frequency domain as opposed to the time domain wherein time dependence is considered explicitly.

In what follows, we will omit the argument $\omega$ in the notations. Thus, we write $\epsilon$ instead of ${\epsilon _\omega}$ for the dielectric permittivity, ${\textbf{E}}({\textbf{r}})$ instead of ${{\textbf{E}}_\omega}({\textbf{r}})$ for the electric field, etc. We will keep in mind, however, that the permittivities of the host and inclusions can depend on $\omega$ and the same is true for the effective permittivity of the composite. Moreover, the effective permittivity can exhibit optical resonances—the phenomenon wherein a physical quantity such as absorption changes rapidly or has sharp spectral features near one or more special (resonant) frequencies [15].

We should also keep in mind that, for any given geometry of the composite, the frequency that we can consider is not arbitrary but bounded from above by the geometrical conditions, which are discussed next.

B. Geometry

The material sample of a composite medium can be of arbitrary shape and does not need to be small compared to the free-space wavelength. However, we require that the statistical properties of the composite be spatially uniform. That is, each physically small volume of the composite must have the same statistical properties as the composite as a whole.

The physically small volume is defined as follows. Consider a region inside the sample, ${\mathbb{V}}$, of sufficiently “simple” shape [16] and volume $V$ and containing $N \gg 1$ randomly distributed and oriented, relatively small, non-overlapping inclusions. We assume that ${\mathbb{V}}$ is entirely contained inside a large sample whose size and shape are not important. We denote the region occupied by $n$th inclusion and its volume by ${\Omega _n}$ and ${v_n}$, respectively, where $n = 1,2, \ldots ,N$. The permittivity of the host medium at the working frequency is denoted by $ {\epsilon _{\texttt h}}$ and the permittivity of $n$th inclusion by ${\epsilon _n}$. In the case when all inclusions have the same permittivity (i.e., are made of the same material), we use the notation $ {\epsilon _{\texttt i}}$ for the inclusion permittivity.

The volume fraction of inclusions $f$ is given by

As any standard effective medium theory, Maxwell Garnett approximation is based on the quasistatic approximation. A necessary condition for this approximation to be applicable is that the characteristic size of each inclusion $a$ be sufficiently small so that $ka \ll 1$, where $k = \omega /c$ is the free-space wavenumber. However, in the case of disordered media, this condition is not sufficient; we will also require that $kL \ll 1$, where $L \gg a$ is the characteristic size of the physically small volume ${\mathbb{V}}$. That is, we assume that the quasistatic approximation holds in ${\mathbb{V}}$.

The assumptions that ${\mathbb{V}}$ is large enough to be statistically representative of the entire sample but small enough for the quasistatic approximation to be accurate are somewhat contradictory. However, these conditions can hold simultaneously at sufficiently low frequencies. If the working frequency is too high, not only Maxwell Garnett approximation can break, but the medium may become not electromagnetically homogeneous in principle. Such media cannot be characterized by spatially uniform effective parameters.

The physically small volume and various physical scales discussed above are illustrated in Fig. 1.

 

Fig. 1. Schematic illustration of the geometry for which Maxwell Garnett approximation is applicable. The physically small volume occupies region ${\mathbb{V}}$ of volume $V$. Inclusions occupy regions ${\Omega _n} \in {\mathbb{V}}$, $n = 1,2, \ldots ,N$, and the volume of $n$th inclusion is ${v_n}$. The characteristic size of ${\mathbb{V}}$, $L$ is much larger than the characteristic size of one inclusion, $a$, and much smaller than the free-space wavelength, $\lambda$. Inclusions are distributed sparsely so that the volume fraction of inclusions, $f$, is much less than unity. Dielectric permittivity at the working frequency is $ {\epsilon _{\texttt h}}$ in the host and ${\epsilon _n}$ in $n$th inclusion. When all inclusions have the same permittivity, we use the notation $ {\epsilon _{\texttt i}}$ for the permittivity of inclusions.

Download Full Size | PDF

C. Fundamental Assumption of Maxwell Garnett Theory

Even if the above conditions hold with arbitrary precision, we cannot solve the electromagnetic problem in a complex medium analytically; additional assumptions are needed to make progress. We now state the fundamental assumption of the Maxwell Garnett theory, which will allow us to derive the effective permittivity of the medium analytically. Unlike the conditions discussed above, this assumption is not about the medium or the electromagnetic frequency but rather about the electric field that is induced in the medium by external sources.

Namely, we assume that inclusions, when polarized by an externally applied field, do not produce any secondary fields (in addition to the external field that is already present in the host material) outside of their interior. This approximation was already stated in the first part of this tutorial [1], but now we formulate it more generally. In [1], we considered geometries for which the field generated by a polarized inclusion integrates to zero in certain regions outside of the inclusion. For example, the field of a polarized sphere integrates to zero in any co-centric spherical shell, which excludes the region occupied by the inclusion itself. For ellipsoidal inclusions, the integration region is an ellipsoidal volume. The rationale used in [1] was that, since the field in the host material integrates to zero, and since it is produced by many inclusions that are distributed uniformly over the entire composite, the secondary field outside of the inclusions cancels out at any point in the host. However, inside the inclusions, the depolarizing field does not integrate to zero and, therefore, must be taken into account.

There is, however, a problem with the above argument. For example, consider the simplest case of electric field produced by a uniformly polarized sphere. The difficulty is that, even in the quasistatic regime, this field is not integrable. The zero integral can be obtained by selecting a special integration region (a finite spherical shell). Choosing a different region would yield a different integration result. For this reason, we have encountered an ambiguity when deriving the anisotropic Maxwell Garnett mixing formula in [1]: the result seemed to depend on the choice of integration region. Of course, once we include the far-zone field into consideration (which makes sense at non-zero frequencies), the zero integration result becomes even more questionable.

One way to deal with the above difficulty is to acknowledge that the fundamental assumption of Maxwell Garnett theory is heuristic and cannot be rigorously justified mathematically. On one hand, this point of view is not very satisfying since it provides no proof that Maxwell Garnett mixing formulas are accurate to second order in $f$. We can be sure that these formulas are accurate to first order in $f$, but the main thrust of Maxwell Garnett-type theories is to obtain corrections to the $O(f)$ result. On the other hand, once we realize what the nature of approximation is, we can apply it easily to a very general class of problems, which is indeed what we will do in the remainder of this tutorial.

Thus, the main approximation of Maxwell Garnett theory is to disregard the electric field produced by polarized inclusions inside the host material. Simply stated, we assume that the field in the host is the externally applied field. The above point is not fully appreciated in the literature. A commonly encountered point of view is that the local field corrections that arise in Maxwell Garnett theory are due to electromagnetic interaction of inclusions. However, this is not so; the local field corrections simply account for the difference between the external and the average fields in a composite [17].

We recognize that the assumption stated above can hold to various degrees of accuracy. The “best” shapes of inclusions that respect this assumption are spheres or ellipsoids or, more generally, convex particles with a high degree of symmetry such as cubes or cylinders. Non-convex inclusions or inclusions without a center of symmetry can pose a problem. For example, a non-convex inclusion with the shape of two touching or slightly overlapping spheres can concentrate the electric field in the gap, which is occupied by the host medium. This concentration of the field outside of the region of inclusions contradicts the basic assumption. For this reason, Maxwell Garnett approximation can break down if inclusions cluster or otherwise come into close contact with each other: two touching inclusions are equivalent to one strongly non-convex inclusion, and this is problematic for the theory.

Ultimately, the problem of homogenization of random mixtures with an arbitrary geometry is very complicated, and verification of the Maxwell Garnett corrections beyond the first order in $f$ is still an open question. We will, however, proceed with calculations under the assumption that the field produced by polarized inclusions in the host material is zero.

5. DERIVATION

The derivation presented below is conceptually similar to the work of Wiener [3], Bragg and Pippard [4], and Oldenbourg and Ruiz [5], but substantially more general. Most importantly, we will allow for random orientations of inclusions and, thus, capture the somewhat counterintuitive case when the medium is isotropic on average but its effective permittivity depends on the shape of non-spherical inclusions. We will also allow for general or varying shapes of inclusions and for multicomponent mixtures.

Perhaps, the simplest and most general way to derive Maxwell Garnett approximation (the “first avenue” according to Sihvola [18]) is based on the following equation:

Let us start with computing $\langle {\textbf{E}}({\textbf{r}})\rangle$. As discussed above, the fundamental assumption of Maxwell Garnett approximation is that the electric field produced by polarized inclusions is averaged out to zero in the host medium and, therefore, can be completely disregarded. It follows from this assumption that

To establish a connection with more familiar quantities, we will relate ${\hat \kappa _n}$ to the dipole polarizability of $n$th inclusion ${\hat \alpha _n}$ and, further, to the dimensionless specific dipole polarizability ${\hat \beta _n}$. The polarizability ${\hat \alpha _n}$ relates the excess dipole moment [19] of $n$th inclusion, ${{\textbf{d}}_n}$, to the external field ${{\textbf{E}}_{\texttt{ext}}}$. The excess dipole moments are defined as follows. Following [1] (Section 2.D), we split the total polarization field ${\textbf{P}}({\textbf{r}})$ into two contributions as

The conventional dipole polarizability of an inclusion in a homogeneous host is introduced as the linear coefficient between the excess dipole moment and the external field, viz.,

6. SPECIAL CASES

A. Two Component Media

In the case of a two-component media, all inclusions have the same permittivity $ {\epsilon _{\texttt i}}$, and all spectral parameters ${s_n}$ are equal to $s$, where $s$ is given by Eq. (2). Under the circumstances, Eq. (23) simplifies to

Let us assume that $\langle \hat \beta \rangle$ is diagonalizable in the above sense either because all ${\hat \beta _n}$ are real or due to some symmetry. Then we can find a set of mutually orthogonal, real-valued vectors of unit length ${{\textbf{u}}_p}$ such that

By using Eq. (25), we can simplify Eq. (24) as

Next, we illustrate how these results work in the case of ellipsoids.

1. Similarly Oriented Ellipsoids

If all ellipsoids are similarly shaped and oriented in space, we have

2. Randomly Oriented Ellipsoids

If the ellipsoids are of the same shape but randomly oriented, then

(a) For spheres (${\nu _1} = {\nu _2} = {\nu _3} = 1/3$):

(b) For randomly oriented disks (${\nu _1} = {\nu _2} = 0$, ${\nu _3} = 1$):

(c) For randomly oriented needles (${\nu _1} = 0$, ${\nu _2} = {\nu _3} = 1/2$):

Equation (32a) is equivalent to Eq. (3). It is widely known and has been stated by us in various forms in [1]. However, Equations (32b) and (32c) are less known and were not derived in the first two parts of this tutorial. These results answer the question posed in Section 1 (the form of Maxwell Garnett mixing formula for randomly distributed and oriented non-spherical particles) for the extreme shapes of needles and disks. It can be seen that the result is not the same as the conventional isotropic formula [Eq. (32a)]; the three expressions in Eq. (32) differ even to first order in $f$. This conclusion is not surprising but not widely appreciated.

B. Isotropic Multicomponent Media with Inclusions of Fixed Shape

In general multicomponent media, each inclusion is characterized by shape, orientation, and permittivity. The averages in Eq. (23) must be computed with respect to the joint probability of these three variables, which can be a rather tedious procedure. However, if some of the variables are statistically independent of the others or fixed, the averaging is simplified.

Consider inclusions of a fixed shape whose permittivities can take discrete values ${\epsilon _m}$ ($m = 1,2,...,M$) with the probabilities ${P_m}$. We further assume that the inclusions are randomly distributed in space and randomly rotated with respect to a laboratory frame. Moreover, the spatial orientation and the permittivity of an inclusion are not correlated. Under these assumptions, the weight factors ${v_n}/\langle v\rangle$ in Eq. (21) are all equal to unity, and we can write

With account of Eq. (34), the general Maxwell Garnett formula [Eq. (23)] takes the form

Rotational averaging can be performed according to the general formula $b(s) = \frac{1}{3}{\rm{Tr}}\hat \beta (s)$, from which we find

1. Randomly Oriented Ellipsoids

All ellipsoids have three non-zero terms in the summation over $k$ with ${w_1} = {w_2} = {w_3} = 1/3$, but the depolarization factors ${\nu _1},{\nu _2}$ and ${\nu _3}$ depend on the aspect ratios. As in the previous subsection, we limit consideration to the three special cases of spheres, thin circular disks and thin needles.

(a) Spheres (${\nu _1} = {\nu _2} = {\nu _3} = 1/3$):

(b) Randomly oriented disks (${\nu _1} = {\nu _2} = 0$, ${\nu _3} = 1$):

(c) Randomly oriented needles (${\nu _1} = 0$, ${\nu _2} = {\nu _3} = 1/2$):

7. DISCUSSION

Random complex media are very complicated from the theoretical standpoint, and detailed analytical or numerical description of such objects beyond various mean-field approximations is a formidable task. In the three-part tutorial, we discussed the physical setting in which electromagnetic properties of a composite can be approximated by a homogeneous effective medium. The Maxwell Garnett approximation is perhaps the most widely used approach to deriving the effective permittivity of such media. The approximation is, however, applicable only in the limit of low volume fraction of inclusions, $f$. Whereas we can almost always be sure that the Maxwell Garnett theory is accurate to first order in $f$ (assuming the physical conditions described in Section 4 were met), the accuracy of higher-order corrections cannot be easily ascertained. Still, we can hope that, in many cases, especially with simple-shape and non-metallic inclusions, the Maxwell Garnett approximation is accurate at least to second order in $f$.

If the Maxwell Garnett approximation fails, the effective medium description can still be applicable, but one needs to use other, more general approaches to deriving the effective parameters. Generically, the problem is known as electromagnetic homogenization. In the second part of this tutorial [2], we have shown how effective medium parameters can be computed in periodic two-component media with arbitrary volume fractions. In principle, one can extend this approach to random media by generating a quasi-random cell (a computer model of the physically small volume) and replicating it periodically in all directions. The resultant medium would still be periodic with a quasi-random cell, and the methods of [2] can be applied to such cases. However, in any realistic setting, this formulation is intractable due to high computational complexity.

Of course, one should be conscious that complex media exhibit many physical phenomena that cannot be captured, even in principle, by an effective medium theory. Such phenomena include strong multiple scattering, speckles, and localization of light. In the physical regimes wherein these phenomena can be observed, Maxwell Garnett theory is definitely inapplicable, and drastically different theoretical approaches must be used. Therefore, before using the Maxwell Garnett theory, one should decide whether it is applicable. Material in Section 4 was designed to help with this determination, but it is by no means exhaustive. Physically, the question is always whether multiple scattering plays a significant role in a composite at the working frequency. If it does, Maxwell Garnett theory is inapplicable.

The mechanisms by which multiple scattering leads to breakdown of Maxwell Garnett approximation have been extensively studied in the literature. One observation is that multiple scattering in random composites results in strong fluctuations of the electric field on the scale of one inter-inclusion distance [15]. Maxwell Garnett, as well as other mean-field approximations that rely on factorization of field correlators (such as Bruggeman’s approximation or the theories of Stroud and Pan [6,7]), cannot capture such phenomena. As soon as the fast-fluctuating component of the field becomes comparable to the mean field, the respective approximation loses accuracy. Another way to look at the role of multiple scattering is the following. Maxwell Garnett mixing formulas contain the volume fraction of inclusions but not their specific locations. On the other hand, multiple scattering clearly depends on higher-order density correlation functions in a composite. One can view $f = \langle \theta ({\textbf{r}})\rangle$ as the first-order correlation function, where $\theta ({\textbf{r}})$ is unity inside an inclusion and zero in the host. But second-order scattering depends on the two-point correlation function ${\rho _2}({\textbf{R}}) = \langle \theta ({\textbf{r}})\theta ({\textbf{r}} + {\textbf{R}})\rangle$. There can be many composites with the same $f$ but different ${\rho _2}({\textbf{R}})$. According to Maxwell Garnett approximation, all such media have the same effective properties. But in practice the multiply scattered light will be reflected and transmitted differently by samples with same $f$ but different ${\rho _2}({\textbf{R}})$. Therefore, Maxwell Garnet approximation cannot be accurate if second-order scattering is non-negligible. As the strength of scattering increases, higher-order correlators also come to the fore, and complex resonant phenomena occur. The connection between the higher-order density correlation functions, multiple scattering, and breakdown of Maxwell Garnett approximation was pointed out and investigated in [22].

We, thus, conclude that second-order scattering between inclusions is not captured correctly by Maxwell Garnett approximation. In this respect, it would be useful to clear one common misconception. Namely, it is often stated that the local field correction [essentially, the denominator or the inverse in the main result Eq. (23) above] appears due to the electromagnetic interaction of inclusions. For example, this point of view was expressed in [12]. Although this misconception does not affect any results and can, therefore, be seen as a matter of interpretation, we think that it is useful to recognize that there is no interaction of inclusions in Maxwell Garnett theory at all [17]. The local field to which the inclusions are subjected is the same as the external (applied) field, and the latter would not be changed if we remove the inclusions. The only regions of space that are affected by the inclusions (according to the approximation) are those that the inclusions occupy. It is the contribution of these regions to the average electric field and displacement that is important in Maxwell Garnett theory.