Applicability of Effective Medium Approximations to Modelling of Mesocrystal Optical Properties

Rigorous superposition T-matrix method is used to compute light interaction with mesocrystalline structures. The results are used to validate the applicability of effective medium theories for computing the effective optical constants of mesocrystal structures composed of optically isotropic materials. It is demonstrated that the Maxwell-Garnett theory can fit the rigorous simulation results with an average accuracy of 2%. The thus obtained refractive indexes can be used with any electromagnetic simulation software to represent the response of mesocrystals composed of optically small primary particles arranged into a cubic type lattice structures.


Introduction
Mesoscopically structured crystals, or mesocrystals, are highly ordered nanoparticle superstructures, showing organisation of nanocrystals in a crystallographic register [1,2].Mesocrystals composed of optical materials can significantly extend the range of material properties available for optical devises and applications.Optical properties of mesocrystals can be fine tuned by changing the size of nanocrystallites, the mesocrystal's constitutive parts, their density and internal arrangement [3].
In general, optical structures can be divided into three main categories with respect to the wavelength of the electromagnetic radiation-with dimensions much larger, comparable or much smaller than the wavelength of the electromagnetic wave.Depending on the mesocrystal's nanocrystallites size, it should fit in one of the categories.Theoretical descriptions of the three mentioned types of mesocrystals can be different.Structures of the first type, with nanocrystallites much larger than the incident light's wavelength, can be well described by the ray optics.Mesocrystals with building blocks comparable to the wavelength are the subject of the photonic crystal research, a subject that was intensively studied over the last three decades [4].Mesocrystals with nanocrystallite sizes significantly smaller than the wavelength of interacting light, should in principle be classified as metamaterials [5], because the light does not feel individual nanocrystals and only experience the collective response of the structure.Conventional metamaterials are composed of artificial atoms with strong optical response, which is often due to resonances in particles or structures made of metals or high dielectric constant materials [6].Usually, interactions between metamaterial-atoms spoil the strong individual resonances.To keep the interaction low, the inter element distances are made relatively large, which in return leads to relatively low volume fraction of the meta-atoms in conventional metamaterials.In mesocrystals, to the contrary, the material volume fraction can reach the values close to the 70% [7,8].The standard effective medium theories like averaging, Maxwell-Garnet [9] or Bruggeman [10] approximations are believed to adequately describe particulate systems with low particle volume fractions [11].Thus, the aim of this paper is to study how well the above mentioned effective medium theories can describe the light interaction with mesocrystals.The discussion here is limited to mesocrystals with homogeneous optical response, which implies that the primary particle material has an isotropic refractive index and the primary particles are arranged in a lattice with cubic symmetry.For applications relying on birefringent mesocrystals, which are possible in principle but were not discussed so far [2], an extension of the proposed approach accounting for the anisotropy of the primary particles or the mesoscopic lattice would be required.
The results of a rigorous electromagnetic simulation tool should be compared to those obtained from the effective medium approximations.Ideally, one would simulate the light interaction from a homogeneous medium onto a macroscopic slab of mesocrystalline material.Unfortunately, the modern computers are not powerful enough to simulate such an object.The problem is in the fine discretization needed to mesh the nanocrystallites and the number of nanocrystallites to be taken in consideration, which very quickly fills up the memory of any computer.Another possibility is to simulate the response of a finite size mesocrystalline structure, like, for example, a mesocrystal of spherical shape, and compare the results with the Lorentz-Mie simulations combined with effective medium theories.This approach was applied in [12][13][14][15] to study the accuracy of the effective medium description of spherical objects composed of random particulate media.Simulations presented in this paper are closely related to one of the prospective applications of mesocrystals-the functional pigments [3], where the optical response of the pigment particles can be tuned by controlling the mesocrystal's density and shape, and the reported results can be used to optimize the properties of such systems.
The variety of available electromagnetic simulation methods [16] like Separation of Variables (SVM), Finite Difference Time Domain (FDTD), Finite Element (FEM), Discrete Dipole Approximation (DDA), Null-Field (NFM) to name a few, is caused by the diversity of materials, system topologies, applications and interrogations.The relative simplicity and straightforwardness of object representation made the FDTD and FEM approaches very popular in the recent time.Nevertheless, the drawback of this simplicity and straightforwardness are the necessity to deal with reflections from the computational domain boundaries and the need to discretize not only the structures but also the free space, which in the case of loose structures like particle agglomerates results in numerical loads unmanageable by modern computers.Alternatives, lacking the mentioned problems with the computation window boundaries and free space discretization, make use of the Transition Matrix (T-matrix) approach, which summarizes the optical response of an object in form of a matrix, relating the field expansions of the incident and scattered fields, see Section 4.4 for more details.The superposition T-matrix method [17] is applied in this paper to simulate the light scattering on spherical mesocrystalline ZnO nanoparticles.The T-matrix simulation results are then compared with the Lorenz-Mie light scattering simulations [18] for effective dielectric constants representing the mesocrystal.Details on the Lorenz-Mie and T-matrix simulations as well as effective medium theories used to compute the dielectric constants of mesocrystals are presented in the Materials and Methods section of the paper.

Light Scattering on Mesocrystalline Particles
In our study, mesocrystals are represented by dense aggregates of spherical nanocrystals.Figure 1a shows the geometry of a mesocrystalline particle composed of N p = 675 identical spheres with diameter d = 15 nm arranged in a face-centered-cubic (fcc) lattice with the lattice constant a = 15 nm.The geometry was defined by selecting only the whole particles fitting within a sphere of 160 nm diameter centered on a lattice node.Optical properties of the mesocrystalline particles of Figure 1a are determined, in the first place, by the properties of constituent nanocrystals.However, as the neighbour nanocrystals are in a close contact, the interaction between the nanocrystallites influence the overall performance of the aggregated particle.Light scattering on such an object can be computed using the superposition T-matrix approach [17], described in Section 4.4.The method takes the positions and dimensions of the nanospheres and the dielectric constants of the matrix and nanocrystals materials, which are water and ZnO in this case, as input parameters.The output of the method is a set of the field expansion coefficients (which are too complex and bulky to be presented here) allowing us to compute electromagnetic fields at any point in space.The field expansion coefficients are usually used to compute extinction cross section of the simulated objects, which can be directly compared to light transmission measurements of diluted particle suspensions [3,8].The extinction cross section spectrum obtained with the superposition T-matrix method for the ZnO mesocrystalline particle of Figure 1a in water is shown in  as the neighbour nanocrystals are in a close contact, the interaction between the nanocrystallites ¿ influence the overall performance of the aggregated particle.Light scattering on such object can be computed using the superposition T-matrix approach [17], described in Subsection 4.4.The method takes the positions and dimensions of the nanospheres and the dielectric constants of the matrix and nanocrystals materials, which are water and ZnO in this case, as input parameters.The output of the method is a set of the field expansion coefficients (which are too complex and bulky to be presented here) allowing to compute electromagnetic fields at any point in space.The field expansion coefficients are usually used to compute extinction cross section of the simulated objects, which can ¼ be directly compared to light transmission measurements of diluted particle suspensions [3,8].The The superposition T-matrix approach requires significant investments in equipment and also takes a certain amount of computation time.For one single wavelength point, the simulations described above require approximately one hour of computation time, which results in approximately one day for the spectrum presented in Figure 2a.In contrast, a Lorenz-Mie simulation provides an extinction cross section spectrum of a spherical particle within a fraction of a second.To apply the Lorenz-Mie theory to the mesocrystalline sphere of Figure 1a, it should be replaced by an equivalent volume sphere of Figure 1b and the particulate medium composed of the nano-spheres should be replaced by an effective dielectric constant.Popular homogenization approaches, the recipes for replacing complex inhomogeneous media with effective parameters, are presented in Section 4.2.In addition to the dielectric constants representing the constituent materials, the effective medium expressions (4)-( 6) depend on the material filling fraction, f , which is the ratio of the volume occupied by the material to the total volume of the structure.Touching spheres arranged in the fcc-lattice occupy f = 0.74 of the structure volume, the densest packing of equal spheres.The diameter of the effective sphere, D, having the same filling fraction and amount of the sphere's materials is given by the expression which is obtained directly from the equality V s f = N p V p , with V s and V p representing the volumes of the effective sphere and the individual spherical nanocrystal.
The diameter of the effective or representing sphere for the discussed mesocrystal particle composed of 675 d = 15 nm nano-spheres is D = 145.5 nm. Figure 2a shows the extinction cross section spectra computed with the Lorenz-Mie theory for the effective sphere filled by 74% with ZnO for the effective dielectric constants computed by the volume averaging (4), Maxwell-Garnett ( 5) and Bruggeman (6) approximations.The corresponding extinction cross section are presented by the dotted, solid and dashed lines respectively.The values of the effective dielectric constants real and imaginary parts are shown in Figure 6a,b respectively.
It should be mentioned that, while fitting experimental results by simulated extinction cross sections, spectra are usually normalized by the extinction's maximum value or the leftmost point of the spectrum.This normalization would make the similar spectra of Figure 2a even closer to each other.The particle size dispersion in experimental samples makes it almost impossible to judge which effective medium theory is better suited to present the mesocrystalline medium.In our case, only simulation results are compared, and the aim is to perfectly match the approximate Lorenz-Mie simulations with the rigorous T-matrix extinction cross section results.The comparison is presented in Figure 2b, where the deviation of the approximate extinction cross section is plotted as a function of the wavelength, with C ext and C tm ext being the approximate and the T-matrix extinction cross sections.The best performing approach is the Maxwell-Garnett approximation providing the accuracy on the order of 2% over the entire visible spectrum.The next one is the Bruggeman effective medium approximation with the average deviation of about 5% and the worst case of 15% deviation.The volume fraction averaging is 15%, with the worst deviation at around 30%.

Fitting the T-Matrix Results with Effective Medium Theories
The filling fraction f opt exactly fitting the Lorenz-Mie effective medium simulation to the values obtained with the T-matrix method was determined for the spectral region from 300 to 700 nm, in order to study the potential of the three effective medium approximations.The task was performed by a numerical optimization procedure, where the values of f were changed in the range from 0 to 1 to set the deviation (2) to zero.For every new value of f the effective sphere diameter and the effective dielectric constant were computed anew and the extinction cross section was then obtained from the Lorenz-Mie simulation.As it can be seen from Figure 3, for the UV part of the spectrum, the volume averaging approach (4) can not provide the exact fit of the T-matrix results, while in the visible the parameter f opt takes values from 0 to 0.3.Using the Bruggeman approximation (6), the fit f opt can be found for any wavelength in the range of interest, although the values of the best fitting filling fraction f opt differ from wavelength to wavelength with the overall spread from 0.5 to 0.7.And only for the Maxwell-Garnett approximation (5) the best fitting filling fraction f opt is almost constant over the computed spectral region, which is exactly the behaviour that is expected from the volume filling fraction.The factual filling fraction f = 0.74 is outside of the optimum range of f opt , which is between 0.75 and 0.77.This shift of f opt is probably due to the boundary effects, which can be explained by the fact, that particles close to the surface respond differently to the electromagnetic excitation compared to those deep inside the medium, because of the missing neighbour particles.For larger structures with a larger amount of particles the volume of the border layer will become small compared to the entire structure volume, and the range of optimum f opt will shift to the factual value.It is difficult, however, to verify this hypothesis, because the discussed structure is the biggest one that can be handled by the computer at my disposal.

Size and Structure Variations
Here the stability of the Maxwell-Garnett approximation with respect to the size and internal composition of mesocrystalline particles is studied.The two different nano-sphere's diameters are 15 and 20 nm, the more drastic changes of the diameters are not practical, because for the chosen mesocrystalline dimensions of 150 nm somewhat smaller constituent nanocrystals will result in a numerically unmanageable geometry, while for a slightly larger nanocrystallite diameter the mesocrystalline particle will degenerate to a collection of several nanocrystallites.The geometric parameters of the structures considered in this subsection are listed in Table 1, and the resulting geometries are shown in Figure 4.  (a-f) Spherical particle geometries, described by parameters listed in Table 1.
The superposition T-matrix approach and the Lorenz-Mie theory were used to compute the extinction cross sections of the listed structures.The filling fraction f was assumed to have the theoretical value of 0.74.The obtained extinction cross sections are shown in Figure 5a, where the vertical axis is in the logarithmic scale to fit all the six structures on one single plot.Structures (b) and (e) of Figure 4  Letters (a-f) in the legends refer to the structure geometries of Table 1 and Figure 4.
The extinction cross sections obtained by the T-matrix approach are graphically indistinguishable from the corresponding Lorenz-Mie results relying on the Maxwell-Garnett effective medium approximation.In order to visualize the discrepancy between the methods, deviations computed by (2) are shown in Figure 5b.The average deviation from the reference T-matrix value is below 3% for all the six geometries.The largest discrepancies are observed for the small structures where, as it was already mentioned, the border layer occupies a significant part of the structure.Nevertheless, the deviation is within 5% even for the problematic wavelength region near 400 nm.For larger structures-(b), (e) and (c), (f)-the deviations gradually decrease to the value of 2%.According to Figure 2b, the volume averaging and Bruggeman approximations applied to the easy structure (c) are less accurate than the Maxwell-Garnett results from the difficult structure (d).

Discussion
The possibility to describe optical properties of mesocrystalline particles using three effective medium approximations-volume average, Maxwell-Garnett and Bruggeman-was investigated.The Maxwell-Garnet effective medium approximation was able to match the extinction cross section obtained by the T-matrix method with the accuracy of around 2%.The accuracy of the effective medium approximation depends on the size of the primary particles and the homogeneity of their distribution in the medium [12][13][14][15].In mesocrystals homogeneity is an intrinsic property provided by the high ordering of the primary particles, it should be noted that in particulate random media with high fraction of the inclusion material lack of homogeneity reduce the accuracy of the approach [15].The appropriateness of the primary particle size can be estimated by computing the size parameter [12] x = πn m d/λ, where n m is the matrix refractive index, d is the primary particle diameter and λ is the wavelength of light.The effective medium concept works best for primary particles with size parameter lower than 0.3, as x reaches 0.5 noticeable deviations of the exact solution and the effective medium approximation can be expected.Dealing with spectra, the size parameter should be estimated for the smallest wavelength, in our case for d = 20 nm primary particles in water, the size parameter computed for λ = 300 nm is x = 0.28.The filling fraction optimization demonstrated that optimum values of the filling fraction f opt for the Maxwell-Garnett approximation is almost constant over the simulated wavelength range, whereas for the other two methods f opt changes with the wavelength, making the initial meaning of the volume fraction concept inapplicable.The fact that optimum values f opt for the Maxwell-Garnett approximation are higher than the factual filling fraction of 0.74 means that using the value of 0.76 instead will result in values closer to the reference T-matrix values, driving the deviation δ c to sub-percent values.This possibility was not pursued because the purpose of this paper was not to perfectly fit a particular mesocrystalline structure by a Lorenz-Mie scattering simulation, but rather to develop a recipe giving a set of dielectric constants suitable to adequately describe the optical properties of dense mesocrystalline media.The recipe should rely on the parameters of the mesocrystals (like, for example, grain size and material filling fraction) that can be obtained from the experiment.It was demonstrated on hand of rigorous T-matrix simulations for structures (c) and (f) of Figure 4 (see Figure 5a) that the nanocrystallite's size is not that important for the mesocrystal optical performance.The factual filling fraction, f , however, is not more than 2% away from the optimized value f opt .Thus, the most accurate effective dielectric constants of mesocrystalline media can be obtained from the Maxwell-Garnett Equation ( 5).The input parameters-filling fraction and the dielectric constants of the nanocrystals and matrix-can be determined experimentally or from the literature.The thus obtained effective dielectric constants can be used in combination with any available electromagnetic solver.However, it should be noted that for structures containing small size features with lateral dimensions of several mesocrystal lattice constants, deviation of the simulation results from the real structure behaviour can be observed, because in such structures nanocrystallite's properties might prevail over the collective effects inside mesocrystals.

Material Optical Properties
ZnO is a material which is broadly used in various industries and applications like, for example, cosmetics, photo voltaics and optical random lasers [19].ZnO can also be grown in mesocrystalline form, as it is demonstrated in [7].Mesocrystalline ZnO particles can be synthesized from a zinc precursor salt dissolved in dimethylformamide (DMF) with addition of polyvinylpyrrolidone (PVP) under vigorous stirring.The reported particles were on average 160 nm large, and the characteristic dimension of the obtained nanocrystals according the the X-ray diffraction analysis was 15 nm, the values used to define the geometry of Figure 1a.Depending on the chemical's relative concentrations, temperature and stirring conditions the synthesis results in the formation of particles with different sizes of constituent ZnO nanocrystals as well as with different particle shapes ranging form spheroids to rod or bone-shaped particles [3,8].
The material parameter which is important for the electromagnetic modelling of optical effects is called dielectric constant or permittivity, .Or alternatively, the index of refraction n = √ is used to characterize the material's optical properties.ZnO is a birefringent material characterized by two refractive index values-ordinary n o and extraordinary n e .According to Park and Schneider [20] the birefringence (∆n = n e − n o ) of ZnO is less than 0.02 in the entire visible range.As the corresponding values of n e and n o are close to 2, ZnO can be considered optically isotropic.ZnO birefringence increases for the UV radiation [20] making the isotropy assumption less accurate.Nevertheless, the considered wavelength range is extended into the UV to 300 nm to evaluate the performance of the averaging procedures for isotropic materials with strong dispersion and absorption.In this study, the ZnO refractive index set reported by Bergström et al. [21] is used, the corresponding real and imaginary parts of the dielectric constant ZnO are shown in Figure 6 by dotted lines with circles.
For electromagnetic simulations of particles embedded into a material other than air the dielectric constant of the surrounding medium-the matrix-should also be taken into account.Here, the ZnO mesocrystalline particles are assumed to be dispersed in water.The refractive index of water n H 2 O is assumed to be constant over the entire wavelength range from 300 nm to 700 nm and equal to 1.33, which results in the permittivity value of H 2 O = 1.77.The voids inside the mesocrystals are filled with the matrix material, such that properties of the two materials-ZnO and water-are participating in the electromagnetic simulations of the mesocrystalline structures.

Effective Medium Approximations
Optical properties of a mixture of two optically different materials are different from those of the ingredients.One of the most simple and popular dielectric constant averaging approaches is the volume fraction weighted averaging The method is often used to compute dielectric constants of suspensions and solutions of materials with similar dielectric constants.
The Maxwell-Garnett effective medium approximation [9] was deduced for colloidal gold particles, but is still applicable to a variety of other materials [11].The effective dielectric constant in this approximation is given by the expression Finally, Bruggeman [10] suggested a different effective medium mixing rule which, contrary the effective medium approximation of Maxwell-Garnett, is symmetric with respect to the host medium and inclusion material.This symmetry is often seen as a facility to extend the applicability of the Bruggeman approximation to a broader range of filling fraction values.Expression (6) constitutes a quadratic equation in bg eff , the analytic expressions of the solutions, however, are quite bulky and will not be listed in this paper.The required values can be easily computed numerically.It should also be noted, that quadratic equations usually have two solutions, and in case the physical solution representing a non-amplifying medium with positive imaginary part ( > 0) of the dielectric constant should be selected.
The average dielectric constant spectra computed using (4)-( 6) for ZnO mesocrystal in water with the filling fraction f = 0.74 are presented in Figure 6 by the dotted, solid and dashed lines respectively.

Lorenz-Mie Scattering
Lorenz-Mie theory [18] provides the full-vectorial solution of the macroscopic Maxwell's equations for spherical objects.One of the best descriptions of the method can be found in the book of Bohren and Huffman [22].The results presented in this paper were obtained with an own implementation of the method in MATLAB R using the recipes from [23,24] to improve the numerical stability of the code for systems with large refractive index differences.

T-Matrix Method
The Transition Matrix (T-matrix) approach summarizes the optical response of an object in the form of a matrix, relating the field expansions of the incident and scattered fields.The approach can be seen as a generalization of the Lorenz-Mie theory to objects other than spheres [25].With the T-matrix computed, the system behavior under different excitation conditions can be evaluated at minimum computation costs.The T-matrix of a particle ensemble can be composed from the T-matrixes of the components [17].Initially, the T-matrix approach was developed in the framework of the Null-Field method [25][26][27][28] or Extended Boundary Condition method [29], and is often not distinguished from them.However, the recent developments suggest that the T-matrix formalism can also be deduced from other simulation techniques [16].Over the years, analytic extensions were developed allowing us to compute scattering from axially symmetric particles [30], superellipsoids [31], rotated particles [32], particle ensembles [17,33], orientation averaged scattering [30], particles close or deposited on a surface [34] for particles with known T-matrix.The progress in the field is reflected in the Comprehensive T-matrix reference database [35][36][37][38][39][40], books [27,28,32,41], and collections of simulation programs [42,43].The results presented in the paper under the name "T-matrix" rely on the superposition T-matrix method [17,33], which was programmed in MATLAB R following [27,32] for defining off diagonal elements of the matrix and using the Lorenz-Mie code described before in Section 4.3 for the diagonal blocks.The code was used to simulate the response of composite particles [44], particulate media [45] and photonic structures [46].To illustrate the package performance, a 250 nanocrystals ensemble accounting for 3 orders of vector spherical harmonics for constituent particles requires per wavelength 40 min of computation time at an eight CPU system (2×IntelXeon E5462 @2.80 GHz, 32 GB DDR2 RAM).The size of the computed T-matrix files on the hard drive is approximately 0.5 GB in this case.The choice of the platform was determined in the first place by the convenience of visualization and post processing of the T-matrix simulation results in MATLAB R , most probably, by the cost of performance which can be achieved with lower level programming languages [42,43].
Figure 2a by open circles connected by the dotted line.The material parameters used in the simulation are described in Section 4.1.(a) (b)

½Figure 2 .
Figure 2. Extinction cross sections of mesocrystalline particle computed using the superposition T-matrix method and the Lorenz-Mie theory (a); Deviation of the Lorenz-Mie results from the exact superposition T-matrix extinction cross section (b).

Figure 3 .
Figure 3. Optimum filling fraction f to match the T-matrix extinction spectrum.

Figure 5 .
Figure 5. Extinction cross sections (a) and relative deviations of extinction cross sections (b) of mesocrystalline particles computed with the Maxwell-Garnet approximation and the T-matrix method.Letters (a-f) in the legends refer to the structure geometries of Table1and Figure4.

Figure 6 .
Figure 6.Real (a) and imaginary (b) parts of the dielectric permittivities.The dotted line with empty circles represent the values of pure ZnO (orientation averaged).The dotted, solid and dashed lines represent the effective permittivities of mesocrystal structure with f = 0.74 volume fraction of ZnO immersed in water with H 2 O = 1.33 2 computed by the averaging, Maxwell-Garnett and Bruggeman approximations.