# Combined Modeling of the Optical Anisotropy of Porous Thin Films

^{1}

^{2}

^{*}

Next Article in Journal

Previous Article in Journal

Previous Article in Special Issue

Previous Article in Special Issue

Research Computing Center, M. V. Lomonosov Moscow State University, Leninskie Gory, 1, 119991 Moscow, Russia

Moscow Center of Fundamental and Applied Mathematics, Leninskie Gory, 1, 119234 Moscow, Russia

Author to whom correspondence should be addressed.

Received: 14 May 2020 / Revised: 23 May 2020 / Accepted: 25 May 2020 / Published: 28 May 2020

(This article belongs to the Special Issue Multiscale Modeling of Emerging Two-Dimensional (2D) Materials and Devices)

In this article, a combined approach for studying the optical anisotropy of porous thin films obtained by the glancing angle deposition is presented. This approach combines modeling on the atomistic and continuum levels. First, thin films clusters are obtained using the full-atomistic molecular dynamics simulation of the deposition process. Then, these clusters are represented as a medium with anisotropic pores, the shapes parameters of which are determined using the Monte Carlo based method. The difference in the main components of the refractive index is calculated in the framework of the anisotropic Bruggeman effective medium theory. The presented approach is tested and validated by comparing the analytical and simulation results for the model problems, and then is applied to silicon dioxide thin films. It is found that the maximum difference between the main components of the refractive index is 0.035 in a film deposited at an angle of 80°. The simulation results agree with the experimental data reported in the literature.

The deposition in a vacuum is the most widely used technique for thin film production [1]. The structural and optical properties of the deposited films depend on the deposition condition, in particular, on the angle between the direction of atoms flow and normal to the substrate surface [2]. Deposition at large angles, when atoms move almost parallel to the substrate, leads to the formation of glancing angle deposited (GLAD) films consisting of the separate nanostructures with different shapes and dimensions [3]. Due to the low reflectance and anisotropy of the refractive index, GLAD films are widely used in optical instruments [4,5,6,7,8].

To describe the anisotropic properties of the deposited films, the anisotropic Bruggeman effective medium approach is usually used [6,9,10,11,12,13,14]. The Bruggeman formalism was initially developed for the calculation of the effective dielectric constant of a medium with randomly distributed spherical inclusions [15]. The values of the dielectric constant of medium and inclusions are different. In Reference [16], the approach was generalized for the case of the ellipsoidal inclusions. The anisotropy of a medium with inclusions was taken into account through depolarization factors depending on the ellipsoids shape parameters.

The voids between separated nanostructures in GLAD films can be considered as inclusions of different shapes and dimensions in the films’ matter. Within the framework of the Bruggeman’s formalism, these inclusions are represented by a set of identical ellipsoids oriented in the same direction. Since the experimental determination of the shape parameters of ellipsoids is still a challenge, these parameters are used to fit experimental results within effective medium models [6,10,11,12,13,14]. On the other hand, shape parameters can be calculated through the Cartezian coordinates of atoms in the clusters obtained by the simulation of the thin film deposition process. Currently, the different methods of the atomistic level, including the ballistic approach [17,18], molecular dynamics (MD) [19,20], and kinetic Monte Carlo (kMC) [14] method were applied to study the deposition process, structural and mechanical properties of GLAD films. In Reference [14], shape parameters were defined for the titanium dioxide GLAD films using the structure tensor, calculated by the integral of density gradient [21]. The main components of this tensor correspond to the averaged inverse shape parameters. The planar birefringence, important in GLAD films applications for polarization optic [22,23], was calculated using the shape parameters and was compared with the experimental data from Reference [14].

In the present work, a combined approach for modeling of the optical anisotropy of porous thin films is presented. First, the MD-based procedure developed earlier in References [24,25,26] is used to simulate the high-energy deposition of silicon dioxide thin films. Further, the averaged pore shape parameters determining the anisotropy properties of porous thin films are calculated using the Monte Carlo based method. The method is quite general and is characterized by high computational efficiency. It can be applied to any porous atomistic clusters, obtained using MD, kMC or ballistic approach and consisting of up to 10^{7}–10^{8} atoms. Finally, the values of the difference between the main components of the refractive index ∆n are calculated in the framework of the anisotropic Bruggeman effective medium approach. The approach is validated by comparing the analytical and simulation results for model problems. The obtained results agree with experimental data for silicon dioxide thin films.

This section presents a combined approach to calculating ∆n. We start from the description of the MD procedure that is used in the simulation of thin film growth. Then a brief overview of the anisotropic Bruggeman effective medium theory as applied to the porous atomistic clusters is given. In the final subsection, an original Monte Carlo based approach to the calculation of the averaged shape parameters is described.

The simulation of thin film growth is the most time-consuming step in the ∆n calculation. Nevertheless, the use of parallel computations allows one to increase the dimensions of clusters in the full-atomistic MD simulation up to several tens of nanometers. An increase in cluster dimensions is especially significant in the case of GLAD films consisting of separate nanoscale structures.

The geometry of the deposition process is shown in Figure 1. The simulation area includes parts of the substrate, film and gas phase above growing film. An increase in the deposition angle α leads to an increase in the porosity of the film and anisotropy of the structure. In particular, deposition at α = 80° leads to the formation of an anisotropic structure with separated slanted columns (Figure 1). For this reason, the dependence of ∆n on α is investigated in a wide range of α from 0° to 80°. In this work, we use clusters that were simulated in our previous works [25,26], except for the cluster deposited at an angle α = 50°.

All silicon dioxide film clusters were obtained using the MD-based step by step method, as described in Reference [27]. The energy of interatomic interactions was calculated in the frame of the DESIL force field [27]:
where q_{i(j)} charge of i(j)-th atom, q_{O} = −0.5q_{Si} = −0.65e, A_{ij} and B_{ij}, parameters of Lennard-Jones potential for the van der Waals interaction, r_{ij} − interatomic distance, A_{SiO} = 4.6 × 10^{−8} kJ·nm^{12}/mol, A_{SiSi} = A_{OO} = 1.5×10^{−6} kJ·nm^{12}/mol, B_{SiO} = 4.2 × 10^{−3} kJ·nm^{6}/mol, B_{SiSi} = B_{OO} = 5 × 10^{−5} kJ·nm^{6}/mol. MD simulation is performed using the GROMACS package [28]. The leap-frog algorithm is used for the integration of particles motion equations [29].

U = q_{i}q_{j}/r_{ij} + A_{ij}/r_{ij}^{12} − B_{ij}/r_{ij}^{6}

Horizontal dimensions of the substrate are 20 nm and 60 nm, the vertical thickness of the substrate is taken equal to 6 nm (Figure 1). At each injection step, silicon and oxygen atoms with the stoichiometric proportion of 1:2 are inserted randomly at the top of the simulation box. Then injected atoms move toward the substrate and form new sub-layers on the surface of the growing film. The maximum value of the deposited film thicknesses is about 45 nm, which corresponds to about 4000 injection steps. The initial values of the silicon and oxygen atoms kinetic energies are E(Si) = 10 eV and E(O) = 0.1 eV, which corresponds to the high-energy deposition processes like ion beam sputtering (IBS) [3]. At each injection step, the NVT (constant number of particles, volume and temperature) ensemble is used. The vertical dimension of the simulation box is increased by 0.01 nm after each injection step in order to compensate for the growth of film thickness. The horizontal dimensions of the box remain unchanged. The Berendsen thermostat [30] is applied to keep the simulation cluster temperature, T = 300 K, constant. The duration of one injection step is 6 ps, and the time step of MD modeling is 0.5 fs.

The simulation was carried out using the equipment of the shared research facilities of HPC computing resources at Lomonosov Moscow State University [31].

In the framework of the anisotropic Bruggeman effective medium approach, the pores in the structure are approximated by a set of ellipsoids [9]. These ellipsoids are randomly distributed over the film volume and have the same shape parameters a_{x}, a_{y}, a_{z} and the same orientations (Figure 2). The original method for the calculation of shape parameters is described in Section 2.3.

After determining the values of a_{x}, a_{y} and a_{z}, the depolarization factors L_{x}_{(y,z)} are calculated as follows [9,16]:
where x^{’} is the integration parameter. These factors depend only on the ratios of a_{x}, a_{y} and a_{z}. Indeed, using the substitution x^{’} = a_{x}^{2}t, we obtain for L_{x}:
where t is the dimensionless variable. Similar equations can be obtained for others depolarization factors.

$${L}_{x\left(y,z\right)}=\frac{{a}_{x}{a}_{y}{a}_{z}}{2}\underset{0}{\overset{\infty}{{\displaystyle \int}}}\frac{d{x}^{\prime}}{\left({x}^{\prime}+{a}_{x}^{2}\left({a}_{y}^{2},{a}_{z}^{2}\right)\right)\sqrt{\left({x}^{\prime}+{a}_{x}^{2}\right)\left({x}^{\prime}+{a}_{y}^{2}\right)\left({x}^{\prime}+{a}_{z}^{2}\right)}}$$

$${L}_{x}=\frac{1}{2}\underset{0}{\overset{\infty}{{\displaystyle \int}}}\frac{dt}{\left(1+t\right)\sqrt{\left(1+t\right)\left(1+t{a}_{x}^{2}/{a}_{y}^{2}\right)\left(1+t{a}_{x}^{2}/{a}_{z}^{2}\right)}}$$

The main components of the dielectric constant tensor ε_{x(y,z)} of the medium consisting of k materials are calculated as follows [9,16]:
where f_{i} = V_{i}/V and ε_{i} are the volume fraction and dielectric constant of the i-th material. The components of the refractive index n_{x}_{(y,z)} are calculated using the Maxwell relation between refractive index and dielectric constant n^{2} = ε.

$$\sum}_{i=1}^{k}{f}_{i}\frac{{\epsilon}_{i}-{\epsilon}_{x\left(y,z\right)}}{{\epsilon}_{x\left(y,z\right)}+{L}_{x\left(y,z\right)}\left({\epsilon}_{i}-{\epsilon}_{x\left(y,z\right)}\right)}=0$$

For the case of porous films Equation (4) is simplified as follows:
where ε_{df} = 2.22 [32] is the dielectric constant of dense silicon dioxide film, f_{1} is the free volume fraction of film. Value of f_{1} is calculated using the film density dependence on the deposition angle ρ(α):

$${f}_{1}\frac{1-{\epsilon}_{x\left(y,z\right)}}{{\epsilon}_{x\left(y,z\right)}+{L}_{x\left(y,z\right)}\left(1-{\epsilon}_{x\left(y,z\right)}\right)}+\left(1-{f}_{1}\right)\frac{{\epsilon}_{df}-{\epsilon}_{x\left(y,z\right)}}{{\epsilon}_{x\left(y,z\right)}+{L}_{x\left(y,z\right)}\left({\epsilon}_{df}-{\epsilon}_{x\left(y,z\right)}\right)}=0$$

f_{1}(α) = 1 − ρ(α)/ρ(0)

In Equation (6), it is assumed that a normally deposited film (α = 0) has no voids and the density of the dense fraction in films, deposited under different α, is equal to ρ(0).

A sample of the growing film structure is shown in the right part of Figure 3. The elongated pores are formed between columns and are highlighted in dark blue. It is seen that the pores are slanted approximately in the same direction.

In the framework of the anisotropic Bruggeman effective medium approach, it is necessary to calculate the averaged pores dimensions a_{x}, a_{y} and a_{z} along the axes of the coordinate system. Furthermore, these parameters are considered as ellipsoids shape parameters in calculating the depolarization factors, see Equations (2) and (3). The coordinate system can be defined by the plane of incidence (x,z) and the plane (x,y) slanted to the substrate plane at an angle β [9]. The Y axis is directed parallel to the substrate plane and perpendicular to the incidence plane, due to the symmetry of the deposition process geometry (Figure 3). The Z axis should be oriented along the direction in which the pores are elongated. This direction is given by the angle β_{m} at which the ratio a_{z}/a_{x} reaches its maximum.

- The procedure for calculating the a
_{x(y,z)}and n_{x(y,z)}values consists of the following steps: - The initial value of β angle is chosen.
- Set of N points with random coordinates (x
_{i}, y_{i}, z_{i}) inside the cluster is generated. For each of the N_{in}points located inside the pores, the pore dimensions a_{xi}, a_{yi}, and a_{zi}are calculated, as shown in Figure 3. The periodic boundary conditions in the horizontal plane are used (see dotted rectangles in the top right part of Figure 3). In the vertical direction, the dimensions of the open pores are limited by the upper boundary of the film. - The averaged pores dimensions a
_{x(y,z)}, are calculated as follows:$${a}_{x\left(y,z\right)}\left({N}_{in}\right)=\frac{{{\displaystyle \sum}}_{i=1}^{{N}_{in}}{a}_{x\left(y,z\right)i}}{{N}_{in}}$$ - Step 2 and 3 are repeated with growing N
_{in}values. The procedure is finished when a_{x(y,z)}(N_{in}) are varying insignificantly with a further increase in N_{in}. The ratio a_{z}/a_{x}is calculated. - The values of L
_{x}_{(y,z)}and n_{x(y,z)}are calculated using Equations (3) and (4). The difference of the main components of the refractive index ∆n is calculated. - The value of β is changed, and steps 2–5 are repeated.
- The a
_{x(y,z)}(N_{in}) values calculated using β, corresponding to the maximum of a_{z}/a_{x}ratio, are the shape parameters.

Let us prove that Equation (7) gives the a_{x(y,z)} values averaged over the pores volume. Indeed, the volume element dV(a_{z}), including all of the points giving the same value of a_{z}, is defined in Figure 3:

dV(a_{z}) = a_{z}(x,y)dxdy

In the frame of the Monte Carlo method [33] and in the N_{in} →∞ limit, the next condition is satisfied:
where N(a_{z}) is the number of points, giving a_{z} value. Thus, the summation in (7) can be replaced by integration:

$$\frac{N\left({a}_{z}\right)}{{N}_{in}}=\frac{dV\left({a}_{z}\right)}{{V}_{p}}$$

$${a}_{z}=li{m}_{{N}_{in}\to \infty}\frac{{{\displaystyle \sum}}_{i=1}^{{N}_{in}}{a}_{zi}}{{N}_{in}}=\frac{{\displaystyle \iint}{a}_{z}\left(x,y\right)dV\left({a}_{z}\right)}{{V}_{p}}=\frac{{\displaystyle \iint}{a}_{z}^{2}\left(x,y\right)dxdy}{{V}_{p}}$$

The values of a_{y(z)} are calculated in a similar way. For the case of ellipsoidal pores of different dimensions, oriented in the same direction, Equation (10) gives values of a_{x(y,z)} which are one and a half times higher than the averaged values of the ellipsoids semi-axes (Appendix A). However, this difference does not matter, since only ratios of a_{x}, a_{y} and a_{z} are required in Equation (3).

The pore boundary search algorithm takes into account that every atom in the clusters is considered as a hard sphere with a certain radius R. The values of the atomic radii are discussed in the text above Table 1, in Section 3. Thus, the pore surface consists of fragments of spherical surfaces centered on the atoms. To calculate the a_{x(y,z)I} values that are required in Equation (7), the next algorithm is used:

- The coordinate system is centered in the probe point, which is randomly chosen.
- The distance from each atom to the X axis is equal to ${r}_{A}=\sqrt{{z}_{A}^{2}+{y}_{A}^{2}}$.
- If r
_{A}is less than atomic radius R, the distance from the probe point to the point of intersection of the X axis with the sphere of radius R, centered on the atom, is calculated as follows (Figure 4):$${x}_{2A}={x}_{A}-\sqrt{{R}^{2}-{r}_{A}^{2}}$$ - The minimum value of x
_{2A}is equal to x_{2i}. - The x
_{1i}value is calculated similarly. Then the a_{xi}value is calculated as a_{xi}= x_{2i}– x_{1i}. - The y
_{1(2)i}and z_{1(2)i}values are calculated in a similar way. The values of a_{x}, a_{y}and a_{z}are calculated according to Equation (7).

To validate the method for calculating shape parameters, two model geometries of the pores are considered: A spherical pore and a set of ellipsoidal pores of various dimensions (right side of Figure 5). Pores are created in the atomistic cluster of silicon dioxide film obtained by high-energy deposition. Horizontal dimensions of the cluster are 20 nm and 60 nm, the vertical dimension is 30 nm, and the directions of the coordinate axes can be seen in Figure 3. Atoms whose centers are inside the pores are removed.

The diameter of the spherical pore D is taken equal to 6 nm. For a continuous medium, the integral on the right-hand side of Equation (10) can be calculated analytically. This gives a_{x} = a_{y} = a_{z} = 3D/4 = 4.5 nm. The results of numerical simulation are shown in the plots in the upper part of Figure 5. Values of a_{x(y,z)} converge to approximately 4 nm with an increase in the number of probe points. The small difference between the results of analytical and numerical solutions is due to the difference between the pore surface and the spherical surface. Indeed, the pore surface is formed by the fragments of spheres centered on the atoms closest to the center of the pore (Figure 5, right part). Since these fragments are partially located inside a sphere of diameter D, the calculated shape parameters are smaller than the analytical parameters for a continuous medium.

The shape parameters are also calculated for the three ellipsoidal pores. The principal axes of all ellipsoids are directed along the coordinate system axes (Figure 3), the value of β is equal to 30°. Ellipsoids are placed inside a cluster so as to avoid their overlapping. To diverse the shape of the ellipsoids, values of their semi-axes are taken equal to a_{x}_{1;2;3} = 8; 7; 4 nm, a_{y}_{1;2;3} = 4; 5; 5 nm, a_{z}_{1;2;3} = 2; 3; 6 nm.

The dependences of shape parameters ratios on the number of probe points are shown in the lower part of Figure 5. These ratios converge fast to the theoretical values with an increase in the number of probe points. Some differences between the numerical and analytical results are explained in the same way as for the spherical pore.

To apply the developed approach to the silicon dioxide films, the values of the atomic radii of silicon and oxygen atoms should be determined. As follows from Equation (5), the components of the dielectric constant and refractive index tensors depend on the free volume fraction of the film. The free volume is defined as volume outside spheres centered on the atoms of the cluster. Therefore, the free volume is fully determined by the Cartesian coordinates of the atoms and values of the atomic radii. On the other hand, the dependence of free volume fraction on the deposition angle is determined by Equation (6). Thus, the value of the atomic radii should be chosen so as to reproduce this dependence.

Earlier in the porosity calculation [34] the van der Waals radii of oxygen and silicon atoms were taken equal to R_{vdW}(O) = 0.152 nm [35] and R_{vdW}(Si) = 0.21 nm [36]. However, these radii leads to f_{1}(α) values significantly higher than those obtained using Equation (6), see the bottom row in Table 1. This is due to the large number of small pores appearing between adjacent atoms if atomic spheres do not cover the entire volume of the film. Thus, the values of the atomic radii should be increased in comparison with the values of van der Waals radii. It was found that radii equal to 0.21 nm for oxygen and 0.27 nm silicon allow one to reproduce f_{1}(α) dependence obtained using Equation (6) (third and fourth rows, Table 1).

The convergence of the pore shapes parameters in the thin film clusters with an increase in the number of probe points N is shown in Figure 6. A value of N is approximately ten thousand, which is enough to calculate the values of a_{x(y,z)}. This calculation takes about ten minutes for clusters consisting of 10^{6} atoms. Since the algorithm has linear scaling with the number of atoms, it takes about one hour for clusters with 10^{7} atoms. As can be seen in Figure 6, an approximate estimation of ratios of a_{x(y,z)} values can be done much faster. In any case the calculation of the of a_{x(y,z)} values requires much less time than the atomistic simulation of the deposition process using MD or kMC methods. Besides, the algorithm can be easily parallelized since the probe points in different parts of the clusters can be chosen independently.

Dependencies of the shape parameters and depolarization factors on the orientation of the coordinate system axes are shown in Figure 7. Due to the symmetry of the deposited structure, the differences in these values achieve maximum when the Z axis is directed almost parallel to the slanted columns (see the right part of Figure 7). This orientation of axes ensures the maximum value of a_{z} close to 4 nm (see the left part of Figure 7). On the contrary, with this axes’ orientation, the value of a_{x} reaches its minimum, since the X axis is directed almost perpendicular to the elongated pores. Value of a_{y} changes insignificantly with an increase of rotation angle β.

The dependencies of the depolarization factors on the angle β are explained in the same way. In accordance with the theoretical results [9,16], the sum of all factors is constant and equal to 1. At a rotation angle β_{m}, the factor L_{x} is maximum, since the ratio of a_{x}/a_{z} is minimal, see Equation (3). The column tilt angle γ can be approximately calculated as 90° − βm (see the right part of Figure 7), which gives about 50° for γ. The tangential rule [17] predicts γ = 53° for the deposition angle α = 70°, which is close to our value of 50°.

The results of the refractive index calculation are shown in Figure 8 and in Table 2 and are compared with experimental data known from the literature. For α = 50° the ∆n value is less than 0.003, i.e., is at the level of anisotropy of the normally deposited SiO_{2} films [32]. The dependences of the refractive index on the rotation angle β for the deposition angles α = 60°, 70°, 80° are similar. The difference ∆n = n_{z} − n_{x} reaches its maximum at a certain value of the rotation angle β = β_{m}. The value of β_{m} increases with increasing of α. This is explained by an increase of the tilt angle of the separated columns in GLAD films with the increase in the deposition angle [3]. Since the pores are elongated along these columns, the angle between the substrate plane and Z semi-axis of the ellipsoids also increases.

As expected, the ∆n value monotonically increases with increasing the deposition angle, due to the increase in the film porosity. The difference between the n_{z} and n_{x} values reaches a maximum value of 0.035 when α is equal to 80°, and β_{m} is approximately 55°.

The calculated values of ∆n are less than experimental ones for all of the deposition angles (Table 2). This can be explained as follows. As it is mentioned in Section 2.2, dimensions of the open pores in the vertical direction are limited by the upper boundary of the film. This leads to the underestimation of the a_{z}/a_{x} ratio calculated for the coordinate system orientation corresponding to the rotation angle β_{m} (Figure 7). To take into account, the contribution of the large elongated pores to the a_{z}/a_{x} ratio and ∆n, the simulation with clusters of larger dimensions should be performed.

The calculated values of the refractive index n = (n_{x} + n_{y} + n_{z})/3 are slightly higher than the experimental ones (Table 2). A possible reason for this is an underestimation of the porosity of the clusters deposited in our MD simulation runs. This also leads to the underestimation of the calculated ∆n values compared to the experimental data, see the paragraph above.

The error in the refractive index in the experiments cited in Table 2 is about 0.005 [6]. When modeling, various sets of random numbers were used in the Monte Carlo method. This led to a difference in the refractive index values less than 0.002. The convergence of the calculated values with an increase in the number of probe points was ensured in our simulation experiments, see plots in Figure 6.

To summarize, the reached level of the simulation accuracy is close to that reported in Reference [14]. To increase the accuracy of the developed approach, it is necessary to increase the dimensions of simulation clusters. This requires the use of more high-performance computational recourses. It should be emphasized that the model does not have any adjustable parameters used to reproduce the experimental results.

In this study, the combined modeling of the optical anisotropy of porous silicon dioxide thin films was carried out. Atomistic clusters representing the structure of the films are obtained using the full-atomistic MD simulation of the deposition process. The pores in films structure are approximated by a set of ellipsoids having the same shape parameters and the same orientation. These parameters are calculated in the framework of the Monte Carlo based method, and are further used in the anisotropic Bruggeman effective medium theory to determine the difference of the main components of the refractive index ∆n.

The ∆n values are calculated for the high-energy deposited clusters. Deposition angles vary from 0° to 80°, which correspond to the normal and glancing angle deposition, respectively. Normally deposited film is considered as a homogeneous structure without void volume. The increase in the deposition angle leads to the decrease of the film density, due to pore formation. The maximum value of the difference of the refractive index main components is equal to 0.035. While the reached level of the simulation accuracy is close to that reported earlier [14], the calculated values of ∆n are still somewhat lower than the experimental ones. Possible reasons for this are discussed, and it is pointed out that even larger atomistic clusters should be considered to increase the accuracy of simulation experiments.

The developed method of combined modeling of the anisotropic properties of atomistic clusters can be applied to other oxide films. For the calculation of ∆n values, information is required about the deposition conditions—including the energy of the incoming atoms, the value of the deposition angle, and the substrate temperature.

Conceptualization, A.V.T., F.V.G, V.B.S.; methodology, F.V.G; software, F.V.G. and V.B.S.; validation, F.V.G.; formal analysis, F.V.G; investigation, F.V.G.; resources, A.V.T. and V.B.S.; data curation, A.V.T. and V.B.S.; writing—original draft preparation, F.V.G; writing—review and editing, A.V.T. and V.B.S; visualization, F.V.G.; supervision, A.V.T. and V.B.S.; project administration, V.B.S.; funding acquisition, A.V.T. and V.B.S. All authors have read and agreed to the published version of the manuscript.

The work was supported by the Russian Science Foundation (Grant No. 19-11-00053).

The authors declare no conflict of interest.

Consider a medium with N ellipsoidal pores that have different shape parameters and are oriented in the same direction. For this case, the semi-axis value averaged over the pores volume is calculated as follows:
where V_{i} = (4π/3)a_{i}b_{i}c_{i} is the volume of the i-th ellipsoid, a_{i}, b_{i} c_{i} are the semi-axes of i-th ellipsoid.

$$\langle a\rangle =\frac{{{\displaystyle \sum}}_{i=1}^{N}{V}_{i}{a}_{i}}{{{\displaystyle \sum}}_{i=1}^{N}{V}_{i}}$$

On the other hand, the averaged pore dimensions parameter a_{x} is calculated as follows (see Equation (10):
where integration in each term is performed over the volume of the i-th pore. The boundary of the i-th ellipsoid is defined by the equation:

$${a}_{x}=li{m}_{{N}_{in}\to \infty}\frac{{{\displaystyle \sum}}_{i=1}^{{N}_{in}}{a}_{xi}}{{N}_{in}}=\frac{{{\displaystyle \sum}}_{i=1}^{N}{\displaystyle \iint}{a}_{xi}^{2}\left(y,z\right)dydz}{{{\displaystyle \sum}}_{i=1}^{N}{V}_{i}}$$

$${\left(\frac{x}{{a}_{i}}\right)}^{2}+{\left(\frac{y}{{b}_{i}}\right)}^{2}+{\left(\frac{z}{{c}_{i}}\right)}^{2}=1$$

The a_{xi} value is determined using Equation (A3):

$${a}_{xi}\left(y,z\right)=2{a}_{i}\sqrt{1-{\left(\frac{y}{{b}_{i}}\right)}^{2}-{\left(\frac{z}{{c}_{i}}\right)}^{2}}$$

So:
where yi = b_{i}(1 − (z/c_{i})^{2})^{1/2}. Substituting Equation (A5) into Equation (A2), we obtain:

$$\iint}{a}_{xi}^{2}\left(y,z\right)dydz=\underset{-{c}_{i}}{\overset{{c}_{i}}{{\displaystyle \int}}}\left(\underset{-yi}{\overset{yi}{{\displaystyle \int}}}{a}_{xi}^{2}\left(y,z\right)dy\right)dz=\frac{3{a}_{i}}{2}{V}_{i$$

$${a}_{x}=\frac{3{{\displaystyle \sum}}_{i=1}^{N}{V}_{i}{a}_{i}}{2{{\displaystyle \sum}}_{i=1}^{N}{V}_{i}}=\frac{3}{2}\langle a\rangle $$

The a_{y(z)} values can be calculated in the same way. Thus, the ratios of the shape parameters are equal to the ratios of the averaged principal semi-axes:

$$\frac{{a}_{x}}{{a}_{y}}=\frac{\langle a\rangle}{\langle b\rangle};\frac{{a}_{x}}{{a}_{z}}=\frac{\langle a\rangle}{\langle c\rangle}$$

- Piegari, A.; Flory, F. Optical Thin Films and Coatings; Woodhead Publishing: Cambridge, UK, 2013. [Google Scholar]
- Robbie, K.; Brett, M.J.; Lakhtakia, A. Chiral sculptured thin films. Nature
**1996**, 384, 616. [Google Scholar] [CrossRef] - Hawkeye, M.M.; Brett, M.J. Glancing angle deposition: Fabrication, properties, and applications of micro- and nanostructured thin films. J. Vac. Sci. Technol.
**2007**, 25, 1317–1335. [Google Scholar] [CrossRef] - Zhao, Y.P.; Ye, D.X.; Wang, P.I.; Wang, G.C.; Lu, T.M. Fabrication of Si nanocolumns and Si square spirals on self-assembled monolayer colloid substrates. Int. J. Nanosci.
**2002**, 1, 87. [Google Scholar] [CrossRef] - Woo, S.-H.; Park, Y.J.; Chang, D.H.; Sobahan, K.M.A.; Hwangbo, C.K. Wideband Antireflection Coatings of Porous MgF2 Films by Using Glancing Angle Depositio. J. Korean Phys. Soc.
**2007**, 51, 1501–1506. [Google Scholar] [CrossRef] - Tkachenko, V.; Marino, A.; Otón, E.; Bennis, N.; Otón, J.M. Morphology of SiO
_{2}films as a key factor in alignment of liquid crystals with negative dielectric anisotropy. Beilstein J. Nanotechnol.**2016**, 7, 1743–1748. [Google Scholar] [CrossRef] [PubMed] - Trottier-Lapointe, W.; Zabeida, O.; Schmitt, T.; Martinu, L. Ultralow refractive index optical films with enhanced mechanical performance obtained by hybrid glancing angle deposition. Appl. Opt.
**2016**, 55, 8796–8805. [Google Scholar] [CrossRef] - Barranco, A.; Borras, A.; Gonzalez-Elipe, A.R.; Palmero, A. Perspectives on oblique angle deposition of thin films: From fundamentals to devices. Prog. Mater. Sci.
**2016**, 76, 59–153. [Google Scholar] [CrossRef] - Schmidt, D.; Schubert, M. Anisotropic Bruggeman effective medium approaches for slanted columnar thin films. J. Appl. Phys.
**2013**, 114, 083510. [Google Scholar] [CrossRef] - Liang, D.; Sekora, D.; Rice, C.; Schubert, E.; Schubert, M. Optical anisotropy of porous polymer film with inverse slanted nanocolumnar structure revealed via generalized spectroscopic ellipsometry. Appl. Phys. Lett.
**2015**, 107, 071908. [Google Scholar] [CrossRef] - Beydaghyan, G.; Buzea, C.; Cui, Y.; Elliott, C.; Robbie, K. Ex situ ellipsometric investigation of nanocolumns inclination angle of obliquely evaporated silicon thin films. Appl. Phys. Lett.
**2005**, 87, 153103. [Google Scholar] [CrossRef] - Nerbø, I.S.; Roy, S.L.; Foldyna, M.; Kildemo, M.; Sønderg˚ard, E. Characterization of inclined GaSb nanopillars by Mueller matrix ellipsometry. J. Appl. Phys.
**2010**, 108, 014307. [Google Scholar] [CrossRef] - Schmidt, D.; Schubert, E.; Schubert, M. Optical properties of cobalt slanted columnar thin films passivated by atomic layer deposition. Appl. Phys. Lett.
**2012**, 100, 011912. [Google Scholar] [CrossRef] - Badorreck, H.; Steinecke, M.; Jensen, L.; Ristau, D.; Jupé, M.; Müller, J.; Tonneau, R.; Moskovkin, P.; Lucas, S.; Pflug, A.; et al. Correlation of structural and optical properties using virtual materials analysis. Opt. Express
**2019**, 27, 22209–22225. [Google Scholar] [CrossRef] [PubMed] - Bruggeman, D.A.G. Berechnung verschiedener physikalischer Konstanten von heterogenen Substanzen. I. Dielektrizitätskonstanten und Leitfähigkeiten der Mischkörper aus isotropen Substanzen. Ann. Physik (Leipzig)
**1935**, 416, 636–664. [Google Scholar] [CrossRef] - Polder, D.; van Santen, J.H. The effective permeability of mixtures of solids. Physica
**1946**, 12, 257. [Google Scholar] [CrossRef] - Smy, T.; Vick, D.; Brett, M.J.; Dew, S.K.; Wu, A.T.; Sit, J.C.; Harris, K.D. Three-dimensional simulation of film microstructure produced by glancing angle deposition. J. Vac. Sci. Technol. A
**2000**, 18, 2507–2512. [Google Scholar] [CrossRef] - Grüner, C.; Liedtke, S.; Bauer, J.; Mayr, S.G.; Rauschenbach, B. Morphology of Thin Films Formed by Oblique Physical Vapor Deposition. ACS Appl. Nano Mat.
**2018**, 1, 1370–1376. [Google Scholar] [CrossRef] - Luo, Y.; Lin, M.; Zhou, N.; Huang, H.; Tsai, C.-T.; Zhou, L. Molecular dynamics simulation study of the microstructure of a-Si:H thin film grown by oblique-angle deposition. Phys. B Condens. Matter
**2018**, 545, 80–85. [Google Scholar] [CrossRef] - Joe, M.; Moon, M.-W.; Oh, J.; Lee, K.-H.; Lee, K.-R. Molecular dynamics simulation study of the growth of a rough amorphous carbon film by the grazing incidence of energetic carbon atoms. Carbon
**2012**, 50, 404–410. [Google Scholar] [CrossRef] - Arseneau, S. Junction Analysis: Representing Junctions through Asymmetric Tensor Diffusion. Ph.D. Thesis, McGill University, Montréal, QC, Canada, 2006. [Google Scholar]
- Oliver, J.B.; MacNally, S.; Smith, C.; Hoffman, B.N.; Spaulding, J.; Foster, J.; Papernov, S.; Kessler, T.J. Fabrication of a glancing-angle-deposited distributed polarization rotator for ultraviolet applications. In Advances in Optical Thin Films VI; International Society for Optics and Photonics: Boulder, CO, USA, 2018; Volume 10691, p. 106911C. [Google Scholar] [CrossRef]
- MacNally, S.; Smith, C.; Spaulding, J.; Foster, J.; Oliver, J.B. Glancing-angle-deposited silica films for ultraviolet wave plates. Appl. Opt.
**2020**, 59, A155–A161. [Google Scholar] [CrossRef] - Grigoriev, F.V.; Sulimov, A.V.; Katkova, E.V.; Kochikov, I.V.; Kondakova, O.A.; Sulimov, V.B.; Tikhonravov, A.V. Computational Experiments on Atomistic Modeling of Thin Film Deposition. Appl. Opt.
**2017**, 56, C87–C90. [Google Scholar] [CrossRef] - Grigoriev, F.V.; Sulimov, V.B.; Tikhonravov, A.V. Atomistic simulation of the glancing angle deposition of SiO
_{2}thin films. J. N. Cr. Sol.**2019**, 512, 98–102. [Google Scholar] [CrossRef] - Grigoriev, F.V.; Sulimov, V.B.; Tikhonravov, A.V. Structure of Highly Porous Silicon Dioxide Thin Film: Results of Atomistic Simulation. Coatings
**2019**, 9, 568. [Google Scholar] [CrossRef] - Grigoriev, F.V.; Sulimov, A.V.; Kochikov, I.V.; Kondakova, O.A.; Sulimov, V.B.; Tikhonravov, A.V. Supercomputer Modeling of the Ion Beam Sputtering Process: Full-Atomistic Level; SPIE—The International Society for Optical Engineering: Jena, Germany, 2015; Volume 9627, p. 962708. [Google Scholar] [CrossRef]
- Abraham, M.J.; Murtola, T.; Schulz, R.; Páll, S.; Smith, J.C.; Hess, B.; Lindahl, E. GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX
**2015**, 1, 19–25. [Google Scholar] [CrossRef] - Hockney, R.W.; Goel, S.P.; Eastwood, J. Quiet High Resolution Computer Models of a Plasma. J. Comp. Phys.
**1974**, 14, 148–158. [Google Scholar] [CrossRef] - Berendsen, H.J.C.; Postma, J.P.M.; van Gunsteren, W.F.; DiNola, A.; Haak, J.R. Molecular-Dynamics with Coupling to an External Bath. J. Chem. Phys.
**1984**, 81, 3684–3690. [Google Scholar] [CrossRef] - Voevodin, V.V.; Antonov, A.S.; Nikitenko, D.A.; Shvets, P.A.; Sobolev, S.I.; Sidorov, I.Y.; Stefanov, K.S.; Voevodin, V.V.; Zhumatiy, S.A. Supercomputer Lomonosov-2: Large Scale, Deep Monitoring and Fine Analytics for the User Community. Supercomput. Front. Innov.
**2019**, 6, 4–11. [Google Scholar] [CrossRef] - Jiang, Y.; Liu, H.; Wang, L.; Liu, D.; Jiang, C.; Cheng, X.; Yang, Y.; Ji, Y. Optical and interfacial layer properties of SiO2 films deposited on different substrate. Appl. Opt.
**2014**, 53, A83–A87. [Google Scholar] [CrossRef] - Metropolis, N.; Ulam, S. The Monte Carlo Method. J. Am. Stat. Assoc.
**1949**, 44, 335–341. [Google Scholar] [CrossRef] - Grigoriev, F.V.; Sulimov, V.B.; Tikhonravov, A.V. Dependence of the thin films porosity on the deposition conditions: Results of the molecular dynamics simulation. Mosc. Univ. Phys. Bulletin.
**2019**, 74, 171–175. [Google Scholar] [CrossRef] - Bondi, A. van der Waals Volumes and Radii. J. Phys. Chem.
**1964**, 68, 441. [Google Scholar] [CrossRef] - Allinger, N.L. Calculation of Molecular Structure and Energy by Force-Field Methods. Adv. Phys. Org. Chem.
**1976**, 13, 1–82. [Google Scholar] [CrossRef] - Yang, S.; Zhang, Y. Spectroscopic ellipsometry investigations of porous SiO2 films prepared by glancing angle deposition. Surf. Interface Anal.
**2013**, 45, 1690–1694. [Google Scholar] [CrossRef] - Melninkaitis, A.; Grinevičiūtė, L.; Abromavičius, G.; Mažulė, L.; Smalakys, L.; Pupka, E.; Ščiuka, M.; Buzelis, R.; Kičas, S.; Tolenis, T. Next-generation allsilica coatings for UV applications. In Laser-Induced Damage in Optical Materials; International Society for Optics and Photonics: Boulder, CO, USA, 2017; Volume 10447. [Google Scholar] [CrossRef]
- Lee, J.; Seong, T.-Y.; Amano, H. Oblique-Angle Deposited SiO2/Al Omnidirectional Reflector for Enhancing the Performance of AlGaN-Based Ultraviolet Light Emitting Diode. ECS J. Solid State Sci. Technol.
**2020**, 9, 026005. [Google Scholar] [CrossRef]

A | 0° | 40° | 50° | 60° | 70° | 80° |
---|---|---|---|---|---|---|

ρ(α) | 2.40 | 2.35 | 2.24 | 1.92 | 1.75 | 1.36 |

f_{1}(α) Equation (6) | 0 | 0.020 | 0.067 | 0.200 | 0.271 | 0.433 |

f_{1}(α)^{1} | 0.023 | 0.032 | 0.080 | 0.195 | 0.269 | 0.436 |

f_{1}(α)^{2} | 0.205 | 0.208 | 0.261 | 0.366 | 0.429 | 0.640 |

A | 60° | 70° | 80° |
---|---|---|---|

∆n (calc.) | 0.015 | 0.025 | 0.035 |

∆n (exp.) | 0.025^{1} | 0.04^{1} | 0.05^{1} |

n (calc). | 1.391 | 1.353 | 1.270 |

n (exp). | 1.42^{1}1.31 ^{3}1.33 ^{4} | 1.34^{1}1.325 ^{2}1.21 ^{3} | 1.22^{1}1.20 ^{2}1.16 ^{3}1.25 ^{4} |

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).