1. Introduction

It is well known in physics that sometimes it takes much less effort to study an ideal infinite system than a finite one. This fact has deep roots and far-reaching consequences. One of them is that infinite-theory formalism is frequently used to study real systems which are finite ones.

Consider for example an infinite photonic crystal containing point sources placed periodically with the same period as the system. Obviously, the field generated by these sources should be also periodic and it is sufficient in this case to consider the fields in only one unit cell. The same should be true for the identical system of finite size when boundary effects can be neglected. The important question is: are the boundary effects really weak in finite photonic crystals? One can not answer on this question in advance. There are photonic systems where the boundary effects are important for any size of the system because of the strong interference between boundary and volume. Recently, the impact of the boundary on the local density of states (LDOS) and the frequency (Lamb) shift was investigated in the work [1] for two-dimensional photonic crystals. It was shown that the LDOS and the Lamb shift are sensitive functions of the shape of the cluster. In reality, the boundary effects can be screened due to multiple imperfections and loses existing copiously in photonic crystals [2]. Without the losses, when scattering by particles is strong, the boundaries should affect the fields inside the cluster. This is not the case for an ideal (without defects) infinite system where the particles interact but never scatter as a boundary does.

In this paper we calculate the local density of states (LDOS) in a photonic cluster by using two methods. The first method is the classical local perturbation method (see for example works [3] and [4]) ideally suited for finite structures made of small particles. The second one is a new method taking into account the periodicity of the scatterers in the cluster. In addition it also uses elements of the local perturbation method (LPM). In the following we will call the second method the LPM-Bloch method. The method is based on the assumption that the scatterers (or particles forming the scatterers) are small compared to the incident wavelength and that fields inside a cluster are almost periodic, implicitly implying what the boundary effects are weak. The advantage of the new method is that a much larger cluster can be calculated because the computer memory demand is much less then for the LPM.

In this paper we will compare the two methods for the case when the results of the LPM are exact (with the high numerical accuracy). For this purpose we study LDOS in the clusters made of particles which are small compared to the incident wavelength (see also [4] for comparison between LPM and exact results). We come to the surprising conclusion that the results of the LPM and the LPM-Bloch method do not converge when the size of the cluster grows.

It is interesting to note that in the work [5] the LDOS was studied for the finite photonic structures with the help of the periodic boundary conditions (PBC). It was shown that in band gap region the LDOS tends to the LDOS of an infinite system when the size of the cluster increases. From the physical point of view, the application of the PBC means that one considers an infinite system where the calculation accuracy grows with the density of k vectors (implicitly taken into account) and the density of the vectors is proportional to the cluster’s size. From a calculation with the PBC it is difficult to extract the properties of a realistic cluster.

2. Theory of the LPM-Bloch method

The local density of states ρ in a dielectric medium is proportional to the imaginary part of the trace of the Green’s tensor in polarization space and defined as (see, for example [6] and interesting discussions therein about the definition of the LDOS)

where ω is the circular frequency, c is the velocity of light in vacuum, ε ₀ is the permittivity of the host medium, and Ĝ(r,r′) is the Green’s tensor at observation point r with source placed at point r′. Im and Tr denote imaginary part and trace respectively. Note that the formula (1) holds in any dielectric medium without gain. The expression (1) shows that the LDOS can be easily computed when the Green’s tensor Ĝ(r,r′) is known. In the following we will calculate the Green’s tensor by using the LPM and by using the periodicity of the system.

The equation for the quasiperiodic Green’s tensor Ĝ(r,r′,k) in the finite periodic medium with periodically placed sources is defined by the following equation (see for example [7] where an infinite system was considered)

Here ⊗ defines tensor product, ε(r) is the permittivity of the medium containing scatterers, Î is the 3×3 unitary tensor in polarization space, and δ is the Dirac delta function. The vector r _n≡∑³ _i=1 m_i a _i, where m_i are integers and a _i are the basis vectors of the system having a Brillouin zone with volume V_BZ. Note that the vector k defines the quasiperiodicity or Bloch property of the tensor Ĝ(r,r′,k) such that in the case of the infinite system or a finite system with periodic boundary conditions

and the Green’s tensor Ĝ(r,r′) is found from the formula

Thus, solving the equation (2) relative Ĝ(r,r′,k) one can find the Green’s tensor Ĝ(r,r′) with the help of (4). We assume that the condition (3) is valid for finite size photonic structures. In this case the tensor Ĝ(r,r′,k) can be found with the help of the LPM as it is shown in Appendix. Below we present the final result for the Green’s tensor (the details of the solution of the Eq. (2) can be found in Appendix)

where $k_0^2=\frac{{\omega }^2}{c^2}{\epsilon }_0$ , $\hat{G}(0,\mathbf{r}\prime ,\mathbf{k})=\hat{G}{({\mathbf{r}}_1,\mathbf{r}\prime ,\mathbf{k})\mid }_{{\mathbf{r}}_1=0}$ , and f_n(q) is the Fourier transform of the function describing the shape of the scatterer (see also Appendix) and

is the quasiperiodic Green’s tensor of the free space. The formula (5) is the main result of this section. It explicitly shows that the total Green’s tensor is a sum of the free space term Ĝ ₀ and term due to the scattering by particles. The later is proportional to the contrast of the particles ε_n-ε ₀ and depends on the particle’s shape f_n. The interference between the scatterers as well as their resonance properties are taken into account by only one term Ĝ(0,r′,k) that should be calculated by solving the system of 9 linear equations obtained by substituting r by r _n into (5). In contrast, the LPM requires to solve the system of 9×N ₁×N ₂×N ₃ linear equations that is computer memory and time consuming. The formula (5) can be essentially simplified for identical (ε_n=ε ₁, f_n(q)=f ₁(q)) spherically symmetric (f_n(q)=f_n(|q|)) scatterers. In this case the second integral in (5) can be explicitly calculated with the help of the of method of residues. Finally we obtain

 

Fig. 1. The normalized scalar LDOS of the finite 1D periodic system versus size parameter ωd/2πc for observation point x=0.1 µm. The system consists of 21 layers centered at x=0. The permittivities of the layers are ε_sc=1.5 (Fig. 1(a)) and ε_sc=3 (Fig. 1(b)) and their width is L=0.0125 µm. The period of the system is d=1 µm and the permittivity of the host medium is ε ₀=1. The red curve shows the LPM result (“exact result”) and the blue dashed one represents the LPM-Bloch result calculated by using 201 subdivisions in k-space for numerical integration.

Download Full Size | PDF

where α=2π ²(ε ₁-ε ₀)f ₁(k ₀)/ε ₀ and R_n=|r-r _n| is the distance between the observation point and n-th scatterer. Note that f_n(k ₀)⋍V_n/8π ³, where V_n is volume of the n-th scatterer. Formula (7) gives the Green’s tensor of a medium periodically filled with N identical spherically symmetric particles which are small compared to the incident wavelength. Note that, when the scatterers are not sufficiently small (k ₀ L≥0.2, for example), they should be subdivided into smaller particles with size L_sub such that k ₀ L_sub≤0.1. It is interesting to note that the expression (7) describes the Green’s tensor of infinite photonic crystal when N→∞. Moreover, the formula explicitly shows that even small losses in the host medium (small positive imaginary part of k ₀) limit the number of interacting particles in infinite photonic crystal.

It should be emphasized here that the main assumption of the LPM-Bloch method is the existence of the quasiperiodicity (expression (3)) for a finite size photonic system. Obviously, the assumption is not correct for the lossless systems where light freely propagates from one boundary to another. Weak interference between scatterers (due to small contrast or large distance, for example) and small losses allow to treat finite system as infinite one because of the effective absence of the boundaries. In this case the expression (3) should be correct. When there is a deviation from the quasiperiodicity, the error of the LPM-Bloch method increases with the number of particles (and the size of the cluster as well) as can be seen from the formula (7). Additional discussion about influence of the quasiperiodicity on the error of the LPM-Bloch method can be found at the end of the Sec. 3.2.

 

Fig. 2. The normalized scalar LDOS of the finite 2D system versus size parameter ωd/2πc at the point x=0.01 µm, y=0. The system consists of 13×13 (Fig. 2(a)) and 21×21 (Fig. 2(b)) circles with radius L=0.01 µm and permittivity ε_sc=9. The period of the system is d=0.05 µm and the permittivity of the host medium is ε ₀=1. The scatterers are arranged into square lattice centered at x=y=0. The red curve shows the LPM result (“exact result”). The blue dashed and solid lines represent the LPM-Bloch results calculated by using 11 and 81 subdivisions in k-space (along each axis) for numerical integration.

Download Full Size | PDF

3. Numerical examples

In this section we present the numerically calculated LDOS for photonic clusters by using the LPM and the LPM-Bloch formalism. The error of the LPM-Bloch method, on one hand, increases with number of the scatterers and their contrast. On the other hand, the dimensionality of the system and the distance between the particles effectively decrease the influence of the boundaries (as well as interference between the scatterers) and as a result, the error of the LPM-Bloch method drops. In the case of the vector wave situation becomes more complicated due to the strong near field interference between the particles. When the distance between the particles grows, the interference sharply goes down and the accuracy of the LPM-Bloch calculation effectively increases.

3.1. Scalar case: 1-3D systems

Consider first the scalar wave propagation in simplest one-dimensional (1D) system formed by periodically placed layers. The Fig. 1(a) and 1(b) show the LDOS calculated for 1D systems with the help of the LPM and the LPM-Bloch formalism. The Fig. 1(a) shows the normalized LDOS versus size parameter ωd/2πc for the system consisting of 21 layers at observation point x=0.1 µm. The layers are centered at x=0 and have width L=0.0125 µm and permittivity ε_sc=1.5. The period of the system is d=1 µm and the permittivity of the host medium is ε ₀=1. The red curve shows the LPM result and the blue dashed one represents the LPM-Bloch result. One can see the very good agreement between the results until ωd/2πc=1.5. When we increase the permittivity, the width, or the number of layers the accuracy of the LPM-Bloch method decreases. As an example, the Fig. 1(b) shows the LDOS for the same system consisting of layers with ε_sc=3. Our calculations show that period does not affect the accuracy. This is due to infinite interference length between layers in 1D system. For 2D and 3D systems, in contrast, the interference between particles decreases with distance between them. The decreasing of the accuracy of the LPM-Bloch method with increasing the contrast or the number of layers indicates that lossless finite one-dimensional system has no LDOS quasiperiodicity.

 

Fig. 3. The normalized scalar LDOS of the finite 3D system versus size parameter ωd/2πc at the point x=0.01 µm, y=z=0. The system consists of 13×13×13 (Fig. 3(a)) and 21×21×21 (Fig. 3(b)) spheres with radius L=0.01 µm and permittivity ε_sc=9. The period of the system is d=0.05 µm and the permittivity of the host medium is ε ₀=1. The particles are arranged into simple cubic lattice centered at r=0. The red curve shows the LPM result (“exact result”). The blue dashed and solid lines represent the LPM-Bloch results calculated by using 11 and 81 subdivisions in k-space (along each axis) for numerical integration.

Download Full Size | PDF

Another interesting example is the two-dimensional (2D) system consisting of the small circles arranged into square lattice. The Fig. 2(a) and 2(b) show the LDOS calculated for the 2D systems with the help of the LPM and the LPM-Bloch formalism. The Fig. 2(a) shows the normalized LDOS versus size parameter ωd/2πc for the system consisting of 13×13 circles at observation point x=0.01 µm, y=0. The circles are centered at x=y=0 and have radius L=0.01 µm and permittivity ε_sc=9. The period of the system is d=0.05 µm and the permittivity of the host medium is ε ₀=1. The red curve shows calculations made with the help of the LPM and the blue ones (dashed and solid) represent the LPM-Bloch result. One can see relatively good agreement between the results despite the high index contrast. The Fig. 2(b) shows the LDOS for the same system consisting of 21×21 circles. One can see that the accuracy of the calculations deteriorates with the number of particles (as in 1D case), however less pronounced. The reason is that 2D systems are more open than 1D ones and the interference in 2D systems decreases as r ^-1/2 while in 1D systems it oscillates.

The last example is the three-dimensional (3D) system consisting of the small spheres arranged into cubic lattice. The Fig. 3(a) and 3(b) show the LDOS calculated for the 3D systems with the help of the LPM and the LPM-Bloch formalism. The Fig. 3(a) shows the normalized LDOS versus size parameter ωd/2πc for the system consisting of 13×13×13 spheres at observation point x=0.01 µm, y=z=0. The spheres are centered at r=0 and have radius L=0.01 µm and permittivity ε_sc ε_sc=9. The period of the system is d=0.05 µm and the permittivity of the host medium is ε ₀=1. The red curve shows the result calculated with the help of the LPM and the blue ones (dashed and solid) represent the LPM-Bloch result. One can see relatively good agreement between the results. The Fig. 3(b) shows the LDOS for the same system consisting of 21×21×21 spheres. One can see that there is no obvious loss of accuracy due to the number of the particles (for this period of the system). We would like to emphasize that we used typical values for the parameters and other values bring us to the same conclusions.

 

Fig. 4. The normalized vector LDOS of the finite 3D system versus size parameter ωd/2πc at the point x=0.01 µm, y=z=0. The system consists of 9×9×9 spheres with radius L=0.01 µm and permittivity ε_sc=9. The periods of the systems are d=0.025 µm (Fig. 4(a)) and d=0.05 µm (Fig. 4(b)) and the permittivity of the host medium is ε ₀=1. The particles are arranged into simple cubic lattice centered at the point r=0. The red curve shows the LPM result (“exact result”). The blue dashed and solid lines represent the LPMBloch results calculated by using 11 and 81 subdivisions in k-space (along each axis) for numerical integration.

Download Full Size | PDF

3.2. Vector case: 3D systems

The Fig. 4(a) and 4(b) show the vector LDOS calculated for 3D system with the help of the LPM and the LPM-Bloch formalism. It show the normalized LDOS versus size parameter ωd/2πc for the system consisting of 9×9×9 spheres at observation point x=0.01 µm, y=z=0. The spheres are centered at r=0 and have radius L=0.01 µm and permittivity ε_sc=9. The periods of the systems are d=0.025 µm and d=0.05 µm respectively. The permittivity of the host medium is ε ₀=1. The red curves show the result calculated with the help of the LPM and the blue ones (dashed and solid) represent the LPM-Bloch result. One can see relatively good agreement between the results. However, there is distinctive difference between the LPM and the LPM-Bloch curves. The difference does not vanish with increasing of the number of subdivisions in k-space and it becomes practically constant (offset to some extent) for low frequencies. When the period grows the difference decreases considerably. The Fig. 5 shows the relative error Δ defined as $\Delta =\frac{{\rho }_{\mathrm{LPM}}-{\rho }_{\mathrm{LPM}-\mathrm{Bloch}}}{{\rho }_{\mathrm{LPM}-\mathrm{Bloch}}}$ versus normalized period of the cluster d/2L at wavelength λ=14.5 µm (static regime). Solid line shows the relative error Δ while the dashed one shows cubically decreasing function ~(d/2L)^-3 for comparison. The Fig. 5 shows that the results predicted by the LPM and the LPM-Bloch formalism become practically identical (Δ=0.015) for values d/2L≥1.5. This indicates that the photonic cluster made of spheres can be treated as infinite one (at least for the LDOS calculations) when d≥3L.

The best way to understand the nature of the difference at low frequencies is to compare the Green’s tensors calculated with the help of the LPM and with the LPM-Bloch method. For our purposes it will be sufficient to consider the cluster made of identical spheres. In this case the Green’s tensor calculated with the help of the LPM takes the following form (the result can be easily obtained form the formula (11) in work [4])

 

Fig. 5. The relative error Δ versus normalized period of the cluster d/2L at wavelength λ=14.5 µm. The solid line shows Δ while the dashed one shows cubically decreasing function ~ (d/2L)^-3 for comparison. The LDOS was calculated at the point x=0.01 µm, y=z=0 for the system consisting of 9×9×9 spheres with radius L=0.01 µm and permittivity ε_sc ε_sc=9. The permittivity of the host medium is ε ₀=1. The particles are arranged into simple cubic lattice centered at the point r=0.

Download Full Size | PDF

where Ĝ(r _n,r′) is the Green’s tensor inside n-th scatterer. The Green’s tensor calculated with the help of the LPM-Bloch method is given by the formula (7). Comparing the formulas (7) and (8) we can see that the expressions will coincide when $\hat{G}({\mathbf{r}}_n,\mathbf{r}\prime )={\int }_{\mathrm{BZ}}\frac{\hat{G}(0,\mathbf{r}\prime ,\mathbf{k})}{V_{\mathrm{BZ}}}{e}^{i\mathbf{k}·{\mathbf{r}}_n}d\mathbf{k},$ , i.e. when the condition of quasiperiodicity (3) will be satisfied. For the static regime when the incident wavelength is much larger then the cluster’s size (k ₀ N ^1/3 d≪1) the tensor $\hat{G}(0,\mathbf{r}\prime ,\mathbf{k})={\hat{G}({\mathbf{r}}_1,\mathbf{r}\prime ,\mathbf{k})\mid }_{{\mathbf{r}}_1=0}$ , satisfies the following equation

while the tensor Ĝ(0,r′) can be obtained from the following system of equations

where one should use $\hat{G}(0,\mathbf{r}\prime )=\hat{G}{({\mathbf{r}}_1,\mathbf{r}\prime )\mid }_{{\mathbf{r}}_1=0}$ . Note that the sums in the expressions (9) and (10) are due to interference between the particles forming the cluster and in the absence of the interference the expressions are identical (after the integration in k space required for the formula (9)). The solution of the equation (10) can be found by substituting of the tensors Ĝ(r _n,r′) in this equation by the expression (10) for the points r _n. For the tensor Ĝ(0,r′) we have

where we wrote down only three first terms. Comparing the expression (11) with the expression (9) (after integration in k space) one can see that they differ by the term proportional to α ²/r ⁶ _jn. The factor α ² r ^-6 _jn is of the order of 1 for smallest distance between the particles (d=2L) when the error is maximal. When the period between the particles grows to d=3L, the factor α ² r ^-6 _jn decreases in 10 times that roughly corresponds to the situation presented in the Fig. 5. This analysis shows that the LPM-Bloch formalism does not take into account the second order interference effects which are very strong for short periods and exist only in finite size structures (for infinite structures the formalism is exact).

4. Conclusions

We have introduced a new approximate LPM-Bloch method calculating fields in the photonic clusters with weak interference between particles. We have demonstrated that the LPM-Bloch method is an accurate numerical tool for the calculation of the LDOS in such photonic structures.

The boundaries affect the LDOS inside a photonic cluster when interference between scatterers is strong. These boundary effects should be taken into account when the LDOS is calculated in low loss and high contrast photonic clusters.