Effective Medium Theory for Multi-Component Materials Based on Iterative Method

Abstract

For biomedical applications in the terahertz band, composites such as macromolecule compounds, biotissues and phantoms are studied. A description of dielectric properties of composite materials using mathematical models has its own fundamental and technological importance. In this work, we present an iterative effective medium theory for multi-component materials. The model has good performance in describing composite materials with more than two components. The theory is evaluated by comparing with the complex permittivity of three different composite materials. A comparison with other commonly used models is given in the form of relative errors.

1. Introduction

Terahertz (THz) radiation is the electromagnetic waves locating in the frequency range of 0.1–10 THz. With its non-ionizing and non-destructive feature due to the low photon energy, THz radiation has attracted attention from researchers for many years, especially in the field of biomedicine, such as molecule [1,2,3], protein [4], cells [5] and cancer diagnosis [6,7,8,9,10]. One important part of studying biomedical samples is to understand their dielectric properties. For this purpose, one of the most using investigation method is THz time-domain spectroscopy (TDS). By recording the waveform of THz radiation which propagates through a sample, the dielectric properties of the sample, such as relative complex permittivity, can be easily measured and extracted.

However, the dielectric properties of a biomedical sample may vary along with the changes of sample compositions and structures. To have a better understanding of the relation between composition and dielectric properties, numerical descriptions of composite materials are needed. The task is usually solved by mathematically modeling the dependence of the dielectric properties on the volume concentration of the components and the structure of the composites. Theoretically, the modeling can be done by considering every particle in the sample individually, but the complexity is enormous. An alternative way is to investigate the effective dielectric properties of the sample macroscopically by using effective medium theory (EMT). Various mathematical models in EMT have been developed with the consideration of the structures [11] and the volume fractions of components [12].

It should be noted that the shape and the size of the components are required for some of the EMT models. For example, the Bruggeman (BG) model has different mathematical expression depending on the geometry of the particles [12]. On the other hand, the Landau–Lifshitz–Looyenga (LLL) model does not consider the shape of particles inside of a composite, but the model has a restriction on component sizes. The LLL model yields accurate enough estimations when the particles is larger than 50 $\mathsf{\mu }$m [13]. Moreover, most EMT models in practice are used to describe two-component composite. It is necessary to develop a new concept of modeling that is more suitable for multi-component composites.

In this paper, an iterative method is proposed for modeling the relative complex permittivity of materials with more than three-components. The concept of the proposed method is given in detail. In total, three different composites were studied. The results of modeling using the iterative method were compared with other popular EMT models for validation.

2. Basics

Let us imagine a composite with fine particles filled in a host and the electric field covers larger space than the scale of the inhomogeneity, the composite can be considered as isotropic and characterized by effective complex permittivity ${\widehat{\epsilon }}_{eff}$. With the average field strength and induction of the electric field over the entire space, the Maxwell equation can be written as

$$\widehat{\mathbf{D}}={\widehat{\epsilon }}_{eff}\widehat{\mathbf{E}}.$$

(1)

The local value of the electric field strength and complex permittivity can be written as:

$$\widehat{\mathbf{E}}=\langle \widehat{\mathbf{E}}\rangle +\delta \widehat{\mathbf{E}},$$

(2)

$$\widehat{\epsilon }=\langle \widehat{\epsilon }\rangle +\delta \widehat{\epsilon },$$

(3)

where $\widehat{\epsilon }={\widehat{\epsilon }}_{eff}$, and $\langle \widehat{\epsilon }\rangle =\frac{1}{V}{\int }_V\widehat{\epsilon }\left(r\right)\mathit{dV}$ is the average value of the complex permittivity over the entire space. By substituting Equations (2) and (3) into Equation (1) and solving the equations with $div\overrightarrow{D}=0$, we can write

$$\langle \widehat{\epsilon }\rangle div(\delta \widehat{\mathbf{E}})=-\langle \widehat{\mathbf{E}}\rangle \nabla \left(\delta \widehat{\epsilon }\right).$$

(4)

In the article [14], a detailed derivation of the solution of this equation is given, which has the form:

$${\widehat{\epsilon }}_{eff}=\langle \widehat{\epsilon }\rangle -\frac{{\langle \delta \widehat{\epsilon }\rangle }^2}{3\langle \widehat{\epsilon }\rangle },$$

(5)

where $\delta \widehat{\epsilon }$—in the general case, is the difference between the complex permittivities of particles and the host [14]. By raising the right side of Equation (5) to the $\beta$ power and carrying out the expansion in the second order, we can get

$${\left({\widehat{\epsilon }}_{eff}\right)}^{\beta }={\langle \widehat{\epsilon }\rangle }^{\beta }\left(1+\frac{(\beta -1)\beta {\langle \delta \widehat{\epsilon }\rangle }^2}{2{\langle \widehat{\epsilon }\rangle }^2}\right).$$

(6)

It is worth to be noted that, with a small $\delta \widehat{\epsilon }$, this formula can be converted to the LLL model [15], the complex refractive index model [16] and the model for layered structures [11] by defining $\beta =1/3,1/2$ and $\pm 1$, respectively. If we take the contribution of the particle volume fraction into account, then using the Taylor series expansion, with small $\delta \widehat{\epsilon }$, Equation (6) can be extended for the N-component mixture:

$${\left({\widehat{\epsilon }}_{eff}\right)}^{\beta }=\sum _{i=1}^N{\phi }_i{\left({\widehat{\epsilon }}_i\right)}^{\beta }.$$

(7)

where ${\widehat{\epsilon }}_{eff}$, ${\widehat{\epsilon }}_i$—are the effective complex permittivity of the composite and its components, respectively. ${\phi }_i$—the volume fractions of the components.

3. Derivation

The main disadvantage of Equation (7) is that it poorly approximates anisotropic composite. Considering the shape of particles and their anisotropy complicates the task of calculating dielectric properties. Fortunately, this problem can be solved using the Maxwell–Garnett and BG models [12].

It is known that when an external electric field ${\widehat{\mathbf{E}}}_{ex}$ is applied to the dielectric, it becomes polarized. The total complex value of the electric field strength $\widehat{\mathbf{E}}$ in the mixture of anisotropic particles:

$$\widehat{\mathbf{E}}=\langle {\widehat{\mathbf{E}}}_d\rangle +{\widehat{\mathbf{E}}}_{ex},$$

(8)

where $\langle {\widehat{\mathbf{E}}}_d\rangle =\frac{1}{V}{\int }_V{\widehat{\mathbf{E}}}_d\left(r\right)\mathit{dV}$ is the average complex value of the dipole field over the entire space. From Equation (8) we may obtain the expression to determine the total dipole moment $\overrightarrow{d_{tot}}$ [12]:

$$\overrightarrow{d_{tot}}=\frac{V}{4\pi }\left({\widehat{\epsilon }}_{eff}-{\widehat{\epsilon }}_h\right)\left(1-\frac{4\pi {\alpha }_p\varrho }{3{\widehat{\epsilon }}_h}\right){\widehat{\mathbf{E}}}_{ex},$$

(9)

where ${\widehat{\epsilon }}_h$ is the complex permittivity of the host material, $\varrho$ is the number of particles in a unit volume, and ${\alpha }_p$ is the polarizability, which is calculated as

$${\alpha }_p=\frac{x_1{x}_2{x}_3}{3}\frac{{\widehat{\epsilon }}_h({\widehat{\epsilon }}_i-{\widehat{\epsilon }}_h)}{{\widehat{\epsilon }}_h+{\nu }_p({\widehat{\epsilon }}_i-{\widehat{\epsilon }}_h)},$$

(10)

where $x_1,x_2,x_3$ are semiaxes of ellipsoid, ${\nu }_p$ is the depolarization factor, which in the general case is a tensor [14].

The value of this parameter for various particle forms is given in Figure 1.

From Equations (9) and (10) we get:

$$\frac{({\widehat{\epsilon }}_{\left(MG\right)p}-{\widehat{\epsilon }}_h)}{{\widehat{\epsilon }}_h+{\nu }_p({\widehat{\epsilon }}_{\left(MG\right)p}-{\widehat{\epsilon }}_h)}=\phi \frac{({\widehat{\epsilon }}_i-{\widehat{\epsilon }}_h)}{{\widehat{\epsilon }}_h+{\nu }_p({\widehat{\epsilon }}_i-{\widehat{\epsilon }}_h)},$$

(11)

where $\phi =\frac{4\pi \varrho }{3}{x}_1{x}_2{x}_3$ is the volume fraction of particles and ${\widehat{\epsilon }}_{\left(MG\right)p}={\widehat{\epsilon }}_{eff}$ is the effective complex permittivity of a two-component composite getting from the Maxwell–Garnett model [12].

BG model was later developed for high concentrations of particles [12]. Equation (11) is used for a N-component mixture:

$$\sum _{i=1}^N{\phi }_i\frac{{\widehat{\epsilon }}_i-{\widehat{\epsilon }}_{\left(Br\right)p}}{{\widehat{\epsilon }}_{\left(Br\right)p}+{\nu }_p({\widehat{\epsilon }}_i-{\widehat{\epsilon }}_{\left(Br\right)p})}=0,$$

(12)

where ${\widehat{\epsilon }}_{\left(Br\right)p}$ is the effective complex permittivity of the composite getting from the BG model.

Even though the above models are widely used (especially the BG model), in practice the simulation result of these models are not accurate for multi-component composites. Homogenization may be used to solve this problem. The multi-component composites can be homogenized if we consider the object as “stacking dolls”, in which the N-component mixture is considered iteratively. That is, for a multi-component mixture, we first consider two of the components. These two chosen components form a effective component and its effective complex permittivity is calculated. Then, the third component and the previous effective component form a new effective component, whose effective complex permittivity may also be calculated accordingly. Furthermore, afterwards the forth, fifth component and so on. The choice of component in each step is defined by the process of sample preparation. This approach can be expressed as

$${\widehat{\epsilon }}_{eff}={\widehat{\epsilon }}_{h_N}=f_{N-1}\left({\widehat{\epsilon }}_{p_{N-1}},f_{N-2}\left({\widehat{\epsilon }}_{p_{N-2}},{\widehat{\epsilon }}_{h_{N-2}}\right)\right),$$

(13)

$${\widehat{\epsilon }}_{h_{N-1}}=f_{N-2}\left({\widehat{\epsilon }}_{p_{N-2}},f_{N-3}\left({\widehat{\epsilon }}_{p_{N-3}},{\widehat{\epsilon }}_{h_{N-3}}\right)\right),$$

(14)

$$⋮$$

$${\widehat{\epsilon }}_{h_3}=f_2\left({\widehat{\epsilon }}_{p_2},{\widehat{\epsilon }}_{h_2}\right)=f_2\left({\widehat{\epsilon }}_{p_2},f_1\left({\widehat{\epsilon }}_{p_1},{\widehat{\epsilon }}_{h_1}\right)\right),$$

(15)

$${\widehat{\epsilon }}_{h_2}=f_1\left({\widehat{\epsilon }}_{p_1},{\widehat{\epsilon }}_{h_1}\right),$$

(16)

where N—serial number of components. For example, for a three-component medium, the number of iterations is equal to 3.

4. Validation

To confirm the validity of this iterative method, we investigated three different types of composite. The experimental complex permittivities of the composites and each component were obtained from published articles or from experiments. The simulated complex permittivities of the composites were then compared with the experimental ones and the relative error between them was calculated by the following formula:

$$\Delta \left(f\right)=\left|\frac{{\widehat{\epsilon }}_{exp}\left(f\right)-{\widehat{\epsilon }}_{theory}\left(f\right)}{{\widehat{\epsilon }}_{exp}\left(f\right)}\right|100\%$$

(17)

where $\Delta \left(f\right)$—the relative deviation, ${\widehat{\epsilon }}_{exp}\left(f\right)$, ${\widehat{\epsilon }}_{theory}\left(f\right)$—the experimental and simulated effective complex permittivity, respectively.

4.1. Titannium Dioxide, Silver and Low Density Polyethylene

The first sample is a three-component material titanium dioxide (TiO${}_2$) (%), silver (Ag) and low density polyethylene (LDPE) with a thickness of 1.37 mm [17]. The complex permittivities for composite particles were taken from [18,19,20]. The authors of this work [17] show that TiO${}_2$ nanoparticles are spherical in shape. For mathematical modeling, the BG model for spherical particles was used: in the first iteration, the effective complex permittivity for a host LDPE ${\widehat{\epsilon }}_h$ with silver ${\widehat{\epsilon }}_2$ was calculated. In the second iteration, the complex permittivity from the first iteration was substituted to the BG model as the new host, and TiO${}_2$ was considered as the doped component and its complex permittivity ${\widehat{\epsilon }}_3$ was also substituted. The new three-component model is written as:

$${\phi }_1\frac{{\widehat{\epsilon }}_h-{\widehat{\epsilon }}_{{Br}_1}}{{\widehat{\epsilon }}_{{Br}_1}+{\nu }_p({\widehat{\epsilon }}_h-{\widehat{\epsilon }}_{{Br}_1})}+{\phi }_2\frac{{\widehat{\epsilon }}_2-{\widehat{\epsilon }}_{{Br}_1}}{{\widehat{\epsilon }}_{{Br}_1}+{\nu }_p({\widehat{\epsilon }}_2-{\widehat{\epsilon }}_{{Br}_1})}=0$$

(18)

$$(1-{\phi }_3)\frac{{\widehat{\epsilon }}_{{Br}_1}-{\widehat{\epsilon }}_{{Br}_2}}{{\widehat{\epsilon }}_{{Br}_2}+{\nu }_p({\widehat{\epsilon }}_{{Br}_1}-{\widehat{\epsilon }}_{{Br}_2})}+{\phi }_3\frac{{\widehat{\epsilon }}_3-{\widehat{\epsilon }}_{{Br}_2}}{{\widehat{\epsilon }}_{{Br}_2}+{\nu }_p({\widehat{\epsilon }}_3-{\widehat{\epsilon }}_{{Br}_2})}=0$$

(19)

where $v_p=1/3$, ${\phi }_2$, ${\phi }_3$, and ${\phi }_1=1-{\phi }_2-{\phi }_3$ are the volume fraction of Ag, TiO${}_2$, and LDPE (host), respectively. First the complex value of ${\widehat{\epsilon }}_{B{r}_1}$ was obtained, then it was substituted into Equation (19). ${\widehat{\epsilon }}_{B{r}_2}$ is the estimated complex permittivity of the sample.

The comparison of the experimental and simulated effective complex permittivity of TiO${}_2$Ag/LDPE and the relative errors of the EMT models are shown in Figure 2. The iterative method gives an average of 0.3% and 25% error when simulating the real and imaginary parts of the effective complex permittivity, respectively. On the other hand, the LLL model has very big errors. The BG model has similar performance as the iterative method, the BG model still has slightly higher errors as 26% error when simulating the imaginary parts of the effective complex permittivity.

4.2. Gelatin, Water, Glycerin, Starch, Bentonite

The second sample is a five-component biopolymer containing gelatin, water, glycerin, starch and 9% bentonite with a thickness of 0.13 mm [21]. Similar as in Section 4.1, BG model with iterative method was used. For this, the effective complex permittivities of gelatin–water, gelatin–water–glycerol, gelatin–water–glycerol–starch, and gelatin–water–glycerin–starch–bentonite were calculated in stages. The results were compared with the standard LLL model Equation (7) and the BG model Equation (12) for the five-component mixture, where the depolarization coefficient ${\nu }_p$ equals to 1/3. The iteration is shown schematically in Figure 3.

The modeling results are shown in the Figure 4. Even though the dispersion properties of this material contain a resonance peak, which is not considered in the simulation [21], the results of numerical calculations correspond to experimental curves (Figure 4). As noted earlier, EMT models are applicable only for homogeneous media. However, this five-component composite is anisotropic as mentioned in the paper [21]. In this regard, big errors occurred in the modeling of the imaginary part of complex permittivity. However, an average deviation of 4% was obtained for the real part of complex permittivity, which is relatively low.

4.3. Gelatin, Water and Oil

The third sample is a three-component phantom made of gelatin (10.9%), water (50.1%) and oil (39%) [22]. The dielectric properties of gelatin and oil were taken from [23,24], respectively. The BG model with iterative method was applied with the following sequence: at the first stage, the complex permittivity of gelatin-water system was calculated, then at the second stage, the complex permittivity of the phantom was calculated by taking the complex permittivity of oil into account. In both stages, two-component BG model of spherical particles was used. The results were compared with the complex permittivity of the phantom from the article [22].

Figure 5 shows the experimental and estimated complex permittivity of phantom and the relative errors of the EMT models. The LLL model and the iterative method have the similar performance, while BG model has much worse performance.

5. Conclusions

In this article, we present a novel concept of applying EMT models based on iterative method. The validation was carried out by comparing the experimental and simulated data of three different composite materials. The results show that the iterative method has better performance comparing with other standard EMT models. In different models, the proposed iterative method is applicable in the THz frequency range.