Statistical Analysis on the Structural Size of Simulated Thin Film Growth with Molecular Dynamics for Glancing Angle Incidence Deposition

: For the purpose of a deeper understanding of thin ﬁlm growth, in the last two decades several groups developed models for simulation on the atomistic scale. Models using molecular dynamics as their simulation method already give results comparable to experiments, however statistical analysis of the simulations themselves are lacking so far, reasoned by the limits imposed by the computational power and parallelization that can only be used in lateral dimensions. With advancements of software and hardware, an increase in simulation speed by a factor of up to 10 can be reached. This allows either larger structures and/or more throughput of the simulations. The paper analyses the signiﬁcance of increasing the structure size in lateral dimensions and also the repetition of simulations to gain more insights into the statistical ﬂuctuation contained in the simulations and how well the coincidence with the experiment is. For that, glancing angle incidence deposition (GLAD) coatings are taken as an example. The results give important insights regarding the used interaction potential, the structure size and resulting important differences for the density, surface morphology, roughness and anisotropy. While larger structures naturally can reproduce the real world in more detail, the results show which structure sizes are needed for these aspects without wasting computational resources.


Introduction
The simulation and modeling of processes is an important area far beyond thin film growth. However, the computation of the deposition of thin films is a challenge, because for a concise model several orders of magnitude in space and time have to be covered by a multi-scale model. As described in [1][2][3], this involves several models from direct simulation Monte Carlo (DSMC), to atomistic models like molecular dynamics (MD) or kinetic Monte Carlo (kMC), to quantum mechanical density functional theory (DFT) calculations. Within the framework of the theoretical description, the complete virtualization of the coating process place enormous demands on both the computational technology and the numerical procedures. In addition, the computational effort and numerical precision have to be carefully balanced for each individual sub-area of the virtual coating. Atomistic models are employed to explore the nanostructure properties of the thin films [1,4,5], while glancing angle deposition (GLAD) is of high interest due to its concise nanostructure [6][7][8]. Results that correspond well between experiment and simulation have already been obtained [6,7] taking into account the size of the simulated structures.
The size is limited due to computational resources, and although parallelized algorithms are used, this allows mainly the increase of the lateral structure dimensions. The growth height is more restricted by the computational resources, because the deposition of atoms can be parallelized to an only very limited extent. Therefore, it is important to choose a structural size, which describes the properties of the material well, while not exhausting resources more than needed. Despite advances in computer technology, software and algorithms, the thickness of structures calculated by classical molecular simulation are still about one order of magnitude away from the structures of experimentally produced single layers, while the factor for multi-layers considering optical thicknesses is about two to three orders of magnitude. The technological advances will not be able to close this gap completely at least for multi layers in the foreseeable future. Therefore, it is crucial to evaluate the significance of the virtual structures as accurately as possible.
The intention of this manuscript is to conclude which structural sizes are useful to consider for virtual material investigations by a quantitative statistical analysis of several simulations and different structure sizes. Also, the impact of the chosen interaction potential is investigated. Only in this way can the maximum benefit by correlation of the theoretical results with the experiments be derived. The following presentation makes a substantial step towards improved knowledge of the significance of the simulations.

Materials and Methods
A major challenge for the simulation is to describe correctly the nanostructure, which occurs especially in thermal coating processes. In particular, the application of GLAD coatings relies on the proper use of these structural properties, making this coating process ideally suited for evaluation of the simulation results. For the investigation the thin film growth of TiO 2 and SiO 2 is analyzed, while the recent publications from Badorreck et al. and Grineviciute et al. [6,7] serve as a basis. The growth of the thin films is investigated by means of atom deposition within the framework of classical molecular dynamics.
Badorreck et al. [6] investigated GLAD coatings of TiO 2 with regard to the birefringence properties. However, the former performed simulations with deposition under 0 • , 30 • , 50 • and 70 • exhibit a small size in lateral dimension of 7 × 8 nm 2 ; therefore, here the results will be compared to structures grown on a substrate area of 30 × 20 nm 2 . In addition, the impact of using different interaction potentials will be analyzed here, namely the potentials from Matsui et al. [9] and Zhang et al. [10]. The latter has the advantage, that it can be combined with potentials for other materials, e.g., SiO 2 , from the same publication, as was performed in [5]. The substrate temperature is assumed to be at 300 K. More sophisticated potentials, e.g., modified embedded atom method (MEAM) [11] like from Lee et al. [12] are not considered here, because the computational effort is at least more than one order of magnitude higher compared to the potentials from Matsui and Zhang, which are of the Buckingham type [13], and therefore are not feasible for the investigation.
Grineviciute et al. [7] analyzed also birefringent properties, here of SiO 2 serial bideposition coatings, with the emphasis on the effect through continued coating onto these porous structures, which is essential for the application in all silica GLAD coatings. Here, the first layer grown by serial bi-deposition is performed for two different substrate areas of 16 × 7 nm 2 and 30 × 20 nm 2 , while the growth of the second layer is repeatedly simulated for the one with 16 × 7 nm 2 . The interaction potential from Zhang et al. [10] is used here. The substrate temperature is also assumed to be at 300 K.
The algorithm for the deposition simulation consists of placing several atoms at the top surface of the structure. For each iteration, the starting positions are computed by a randomly positioned grid, where each position on the grid has its own random position variation. Simultaneously, it is ensured that the distance between each position is larger than the interaction cutoff, which is sufficient for low deposition energies, exceeding thermic energies not significantly. With the direction of their velocity vectors, the atoms are directly moved towards the surface keeping the limits of the interaction cutoff of 1.2 nm, so that the computationally expensive part of the MD is minimized. An illustration of this state from different viewpoints is shown in Figure 1 with the atoms to be deposited viewed in green, traveling into the direction of the atoms marked in blue as indicated by the pink arrow. One MD deposition iteration is performed within Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) [14] in two stages, the first within a NVE (constant number of particles, volume and energy), the second within a NVT (constant number of particles, volume and temperature) ensemble. The integration time step is 1 fs, while the number of integration steps is determined automatically for each deposition iteration depending on how fast the atoms reach the surface and usually lies in the range of 10 4 for the NVE step and 50 for the NVT step with a damping parameter of 0.01. Long-range coulomb interactions are calculated by Wolf's summation [15]. A more detailed description of the algorithm can be found in [6]. The computational effort used for the following simulations ranges, depending on the structure size, from one month to about six months on a system capable of 10 TFlops.  [14] in two stages, the first within a NVE (constant number of particles, volume and energy), the second within a NVT (constant number of particles, volume and temperature) ensemble. The integration time step is 1 fs, while the number of integration steps is determined automatically for each deposition iteration depending on how fast the atoms reach the surface and usually lies in the range of 10 4 for the NVE step and 50 for the NVT step with a damping parameter of 0.01. Longrange coulomb interactions are calculated by Wolf's summation [15]. A more detailed description of the algorithm can be found in [6]. The computational effort used for the following simulations ranges, depending on the structure size, from one month to about six months on a system capable of 10 TFlops.
(a) (b) Figure 1. (a) Illustration from the top of a deposition step with a glancing angle of incidence of 50° directed to the top-left. Atoms to be deposited (green) are placed at the top and moved towards the structure keeping the limits of the interaction cutoff. The blue marked atoms indicate the direction of the depositing atoms and would be hit without any interaction, also indicated by the pink arrow. (b) 3D view of the same illustration with the dimension of the simulation box. Graphics generated with NGL [16,17].
The analysis of the finally grown thin films is then performed for the density behavior and the surface roughness. Further their anisotropy is analyzed as described in Ref. [6].

TiO2
GLAD coatings are manufactured under very specific conditions. The particle energy must be kept low and the coating temperature is low so that the structures can develop in a suitable manner. In addition, the angle of incidence must be set in such a way that the shadowing effect leads to the tailored growth of columns. Therefore, we will compare the results for structures grown for different deposition angles of 0°, 30°, 50° and 70° for the two different substrate areas of 7 × 8 nm 2 and 30 × 20 nm 2 . Comparing the smaller with the larger substrate the surface area increases by a factor of 10. Due to the computational effort, the smaller structures are simulated with the potential from Matsui et al. [9], while for the larger ones the potential from Zhang et al. [10] is used. To allow still a comparison between the structures despite the different used potentials, the structure for 0° and 7 × 8 nm 2 substrate size is simulated also with the potential from Zhang.
The structures in Figure 2 and especially the structures grown under 70° show that the larger structures exhibit a higher variation of the microstructure than can be modeled by the smaller ones. While the small structure for 70° deposition angle effectively returns one slanted column, which is periodically continued at the lateral borders and actually can be considered as a wall, the larger structure exhibits several individual and complete columns. The analysis of the finally grown thin films is then performed for the density behavior and the surface roughness. Further their anisotropy is analyzed as described in Ref. [6].

TiO 2
GLAD coatings are manufactured under very specific conditions. The particle energy must be kept low and the coating temperature is low so that the structures can develop in a suitable manner. In addition, the angle of incidence must be set in such a way that the shadowing effect leads to the tailored growth of columns. Therefore, we will compare the results for structures grown for different deposition angles of 0 • , 30 • , 50 • and 70 • for the two different substrate areas of 7 × 8 nm 2 and 30 × 20 nm 2 . Comparing the smaller with the larger substrate the surface area increases by a factor of 10. Due to the computational effort, the smaller structures are simulated with the potential from Matsui et al. [9], while for the larger ones the potential from Zhang et al. [10] is used. To allow still a comparison between the structures despite the different used potentials, the structure for 0 • and 7 × 8 nm 2 substrate size is simulated also with the potential from Zhang.
The structures in Figure 2 and especially the structures grown under 70 • show that the larger structures exhibit a higher variation of the microstructure than can be modeled by the smaller ones. While the small structure for 70 • deposition angle effectively returns one slanted column, which is periodically continued at the lateral borders and actually can be considered as a wall, the larger structure exhibits several individual and complete columns.
A simple qualitative consideration of the structures shows for the large structures that the direction of growth seems to change slightly with an increasing layer thickness resulting in curved columns. Such details cannot be seen in the smaller structures. This already indicates a qualitative difference. A simple qualitative consideration of the structures shows for the large structures that the direction of growth seems to change slightly with an increasing layer thickness resulting in curved columns. Such details cannot be seen in the smaller structures. This already indicates a qualitative difference.
The investigation is performed for two different interaction potentials, to also get an idea of the impact of the potential itself onto the structural properties. For a better understanding of the comparison between the thin film growth for different substrate sizes, first the impact of both used interaction potentials from Matsui et al. [9] and Zhang et al. [10] is explored. In Figure 3, a comparison of the density is shown for the small substrate size of 7 × 8 nm 2 , performed four times with the Matsui potential and one time with the Zhang potential. For the thin film growth using the Matsui potential the mean of the density profiles with about 3 g/cm 3 (excluding the region of the substrate and initial growth zone) is plotted with the standard deviation as shaded region. Compared to the growth with the Zhang potential, that results in a significantly lower density of about 2.25 g/cm 3 over the full height. The difference is much larger than the standard deviation of about ±0.25 g/cm 3 shown for the Matsui potential. While here no standard deviation for the simulation with the Zhang potential can be given, the analysis of the following 30 × 20 nm 2 structures divided into several 7 × 8 nm 2 lateral-sized parts with regard to the density profile indicates a very similar standard deviation as shown for the Matsui potential. Compared with experimental results [6], the structures grown using the Zhang potential obviously underestimate the density. The cross check for high energetic grown structures at 10 eV shows also a difference in density for both potentials of 4.05 g/cm 3 and 3.25 g/cm 3 for the Matsui and Zhang potential, respectively. Although the structures in Figure 2 are grown with different potentials, we are analyzing the results without correction factors.
The investigation is performed for two different interaction potentials, to also get an idea of the impact of the potential itself onto the structural properties. For a better understanding of the comparison between the thin film growth for different substrate sizes, first the impact of both used interaction potentials from Matsui et al. [9] and Zhang et al. [10] is explored. In Figure 3, a comparison of the density is shown for the small substrate size of 7 × 8 nm 2 , performed four times with the Matsui potential and one time with the Zhang potential. For the thin film growth using the Matsui potential the mean of the density profiles with about 3 g/cm 3 (excluding the region of the substrate and initial growth zone) is plotted with the standard deviation as shaded region. Compared to the growth with the Zhang potential, that results in a significantly lower density of about 2.25 g/cm 3 over the full height. The difference is much larger than the standard deviation of about ±0.25 g/cm 3 shown for the Matsui potential. While here no standard deviation for the simulation with the Zhang potential can be given, the analysis of the following 30 × 20 nm 2 structures divided into several 7 × 8 nm 2 lateral-sized parts with regard to the density profile indicates a very similar standard deviation as shown for the Matsui potential. Compared with experimental results [6], the structures grown using the Zhang potential obviously underestimate the density. The cross check for high energetic grown structures at 10 eV shows also a difference in density for both potentials of 4.05 g/cm 3 and 3.25 g/cm 3 for the Matsui and Zhang potential, respectively. Although the structures in Figure 2 are grown with different potentials, we are analyzing the results without correction factors.
Looking at the density profiles in Figure 4 for the structures from Figure 2, two interesting aspects are revealed. On the one hand, the density profiles for the small structures are very noisy. As already shown in Figure 3, a standard deviation of ±0.25 g/cm 3 can be expected, therefore a variation of up to ±0.5 g/cm 3 is not uncommon. Here, the increase of the structure size results in a considerable improvement, with a reduction of the density variation to about ±0.1 g/cm 3 . Therefore, the dependency of the density on the deposition angle can be considered with much more significance. While in principle also the small structures show this dependency, the fluctuation makes it especially hard to determine a difference between 0 • and 30 • . The large structures show here very clearly that the density for 30 • is a bit less than for a 0 • deposition angle.  Looking at the density profiles in Figure 4 for the structures from Figure 2, two interesting aspects are revealed. On the one hand, the density profiles for the small structures are very noisy. As already shown in Figure 3, a standard deviation of ±0.25 g/cm 3 can be expected, therefore a variation of up to ±0.5 g/cm 3 is not uncommon. Here, the increase of the structure size results in a considerable improvement, with a reduction of the density variation to about ±0.1 g/cm 3 . Therefore, the dependency of the density on the deposition angle can be considered with much more significance. While in principle also the small structures show this dependency, the fluctuation makes it especially hard to determine a difference between 0° and 30°. The large structures show here very clearly that the density for 30° is a bit less than for a 0° deposition angle.
On the other hand, the larger structures show a significantly lower overall density than the small structures. With the results shown in Figure 3, a significant amount of the difference of 0.75 g/cm 3 can be attributed to applying different potentials. However, the reduction of the density is still larger with a remaining additional drop of about 0.5 g/cm 3 . While the density roughly varies between 2 and 3 g/cm 3 for the small ones, for the larger ones this is reduced to 1 to 1.75 g/cm 3 . In sum that is about a factor of 2. Compared to the experimental results presented in [6], the density for the structures grown with the Zhang potential are significantly too low. Consequently, a modification of the potential seems to be necessary for its application in the future.    . Density profiles in growth direction z of TiO2 for an incidence angle of 0° and a substrate area of 7 × 8 nm 2 . The comparison between using different interaction potentials is shown. For the Matsui potential [9] the standard deviation shown as shaded region is computed from four independent simulations. With the Zhang potential [10] a lower density is obtained.
Looking at the density profiles in Figure 4 for the structures from Figure 2, two interesting aspects are revealed. On the one hand, the density profiles for the small structures are very noisy. As already shown in Figure 3, a standard deviation of ±0.25 g/cm 3 can be expected, therefore a variation of up to ±0.5 g/cm 3 is not uncommon. Here, the increase of the structure size results in a considerable improvement, with a reduction of the density variation to about ±0.1 g/cm 3 . Therefore, the dependency of the density on the deposition angle can be considered with much more significance. While in principle also the small structures show this dependency, the fluctuation makes it especially hard to determine a difference between 0° and 30°. The large structures show here very clearly that the density for 30° is a bit less than for a 0° deposition angle.
On the other hand, the larger structures show a significantly lower overall density than the small structures. With the results shown in Figure 3, a significant amount of the difference of 0.75 g/cm 3 can be attributed to applying different potentials. However, the reduction of the density is still larger with a remaining additional drop of about 0.5 g/cm 3 . While the density roughly varies between 2 and 3 g/cm 3 for the small ones, for the larger ones this is reduced to 1 to 1.75 g/cm 3 . In sum that is about a factor of 2. Compared to the experimental results presented in [6], the density for the structures grown with the Zhang potential are significantly too low. Consequently, a modification of the potential seems to be necessary for its application in the future.   On the other hand, the larger structures show a significantly lower overall density than the small structures. With the results shown in Figure 3, a significant amount of the difference of 0.75 g/cm 3 can be attributed to applying different potentials. However, the reduction of the density is still larger with a remaining additional drop of about 0.5 g/cm 3 . While the density roughly varies between 2 and 3 g/cm 3 for the small ones, for the larger ones this is reduced to 1 to 1.75 g/cm 3 . In sum that is about a factor of 2. Compared to the experimental results presented in [6], the density for the structures grown with the Zhang potential are significantly too low. Consequently, a modification of the potential seems to be necessary for its application in the future.
Regarding the surface morphology in Figures 5 and 6 for the small and larger structures respectively, the morphology shows that the growth for 0 • starts with many columnar structures (mean height of 5 nm), while for 70 • fewer columns are recognized. This can be considered a direct result of the shadowing effect. Also, thicker structures reduce the number of columns, which is caused by the coalescence or overgrowth of the columns. The comparison of the morphology and roughness of the small structures with the larger ones show fewer columns and much less roughness for the small ones.
Regarding the surface morphology in Figures 5 and 6 for the small and larger structures respectively, the morphology shows that the growth for 0° starts with many columnar structures (mean height of 5 nm), while for 70° fewer columns are recognized. This can be considered a direct result of the shadowing effect. Also, thicker structures reduce the number of columns, which is caused by the coalescence or overgrowth of the columns. The comparison of the morphology and roughness of the small structures with the larger ones show fewer columns and much less roughness for the small ones. mean height/nm 0° 30° 50° 70° 5 9 12 Figure 5. Images of the surface morphology for a substrate area of 7 × 8 nm 2 and different mean heights. The scale is shifted for different heights, but its range is kept the same to make all images comparable. mean height/nm 0° 30° 50° 70° 5 9 12 Figure 6. Images of the surface morphology for a substrate area of 30 × 20 nm 2 for different mean heights. The scale is shifted for different heights, but its range is kept the same to make all images comparable.
In Figure 7, the structures from Figure 2 are analyzed with regard to the root mean square (rms) roughness [18]. This is calculated by Welch's method [19] as implemented in NumPy/SciPy [20,21]. The maximum roughness for 7 × 8 nm 2 is reached after 5 to 10 nm while with a larger angle of incidence the maximum is reached faster. When the maximum is reached, the roughness in the following is characterized by a large fluctuation around the maximum with no significant difference between the incidence angles. For 30 × 20 nm 2 the roughness is much higher already for a height of 10 nm and a plateau cannot be identified yet. In addition, the incidence angle of the deposited atoms shows a significant impact on the roughness, the value for 70° compared to 0° at a height of 12 nm is about twice Figure 5. Images of the surface morphology for a substrate area of 7 × 8 nm 2 and different mean heights. The scale is shifted for different heights, but its range is kept the same to make all images comparable.
Regarding the surface morphology in Figures 5 and 6 for the small and larger structures respectively, the morphology shows that the growth for 0° starts with many columnar structures (mean height of 5 nm), while for 70° fewer columns are recognized. This can be considered a direct result of the shadowing effect. Also, thicker structures reduce the number of columns, which is caused by the coalescence or overgrowth of the columns. The comparison of the morphology and roughness of the small structures with the larger ones show fewer columns and much less roughness for the small ones. mean height/nm 0° 30° 50° 70° 5 9 12 Figure 5. Images of the surface morphology for a substrate area of 7 × 8 nm 2 and different mean heights. The scale is shifted for different heights, but its range is kept the same to make all images comparable. mean height/nm 0° 30° 50° 70° 5 9 12 Figure 6. Images of the surface morphology for a substrate area of 30 × 20 nm 2 for different mean heights. The scale is shifted for different heights, but its range is kept the same to make all images comparable.
In Figure 7, the structures from Figure 2 are analyzed with regard to the root mean square (rms) roughness [18]. This is calculated by Welch's method [19] as implemented in NumPy/SciPy [20,21]. The maximum roughness for 7 × 8 nm 2 is reached after 5 to 10 nm while with a larger angle of incidence the maximum is reached faster. When the maximum is reached, the roughness in the following is characterized by a large fluctuation around the maximum with no significant difference between the incidence angles. For 30 × 20 nm 2 the roughness is much higher already for a height of 10 nm and a plateau cannot be identified yet. In addition, the incidence angle of the deposited atoms shows a significant impact on the roughness, the value for 70° compared to 0° at a height of 12 nm is about twice Figure 6. Images of the surface morphology for a substrate area of 30 × 20 nm 2 for different mean heights. The scale is shifted for different heights, but its range is kept the same to make all images comparable.
In Figure 7, the structures from Figure 2 are analyzed with regard to the root mean square (rms) roughness [18]. This is calculated by Welch's method [19] as implemented in NumPy/SciPy [20,21]. The maximum roughness for 7 × 8 nm 2 is reached after 5 to 10 nm while with a larger angle of incidence the maximum is reached faster. When the maximum is reached, the roughness in the following is characterized by a large fluctuation around the maximum with no significant difference between the incidence angles. For 30 × 20 nm 2 the roughness is much higher already for a height of 10 nm and a plateau cannot be identified yet. In addition, the incidence angle of the deposited atoms shows a significant impact on the roughness, the value for 70 • compared to 0 • at a height of 12 nm is about twice as large. This indicates that the lateral size of 7 × 8 nm 2 cannot image the full surface morphology and roughness. For a representative image of the roughness, structures have to cover all characteristic features completely. If the structure and, therefore, the surface area, is smaller, the relevance of the analysis is reduced.
The analysis of the anisotropy is presented in Table 1. While the fill factors are corresponding to the density distribution of the simulations shown in Figure 4, the birefringence show no similar behavior comparing small and large structures. For the substrate size 7 × 8 nm 2 , the birefringence steadily increases with larger deposition angles, and this is not reproduced for 30 × 20 nm 2 . The result of the birefringence ∆n for 0 • angle of incidence and 30 × 20 nm 2 is larger than for 7 × 8 nm 2 , although here a value closer to 0 is expected due to smaller statistical fluctuation. The structure for 30 • deposition angle shows a larger value of ∆n, while for 50 • and 70 • the values are between 0 • and 30 • . This can be a result of the extremely porous structure, which prevents a pronounced anisotropy.  The analysis of the anisotropy is presented in Table 1. While the fill factors are corresponding to the density distribution of the simulations shown in Figure 4, the birefringence show no similar behavior comparing small and large structures. For the substrate size 7 × 8 nm 2 , the birefringence steadily increases with larger deposition angles, and this is not reproduced for 30 × 20 nm 2 . The result of the birefringence ∆ for 0° angle of incidence and 30 × 20 nm 2 is larger than for 7 × 8 nm 2 , although here a value closer to 0 is expected due to smaller statistical fluctuation. The structure for 30° deposition angle shows a larger value of ∆ , while for 50° and 70° the values are between 0° and 30°. This can be a result of the extremely porous structure, which prevents a pronounced anisotropy.
In sum the analysis shows, that regarding especially the Zhang potential for TiO2, there is a need for further investigation and optimization of the potential parameters, due to the produced densities, which are much too low. However, the different structure sizes are still comparable qualitatively. For reproducing features like several columns, changed growth during deposition and imaging the surface roughness to a much better extent, the use of a substrate size of at least 30 × 20 nm 2 can be recommended.

SiO2
For SiO2, an initial porous layer, grown by serial bi-deposition, is continued with a second layer. For the second layer deposition configurations for different incidence angles of 0°, 30° and 50° with a rotating substrate are chosen. Each simulation is performed three times, to evaluate the variation of the continued growth caused by the randomized deposition algorithm. The structures are shown in Figure 8. While the lower thirds of the structures are nearly identical, the upper two thirds exhibit differences. In sum the analysis shows, that regarding especially the Zhang potential for TiO 2 , there is a need for further investigation and optimization of the potential parameters, due to the produced densities, which are much too low. However, the different structure sizes are still comparable qualitatively. For reproducing features like several columns, changed growth during deposition and imaging the surface roughness to a much better extent, the use of a substrate size of at least 30 × 20 nm 2 can be recommended.

SiO 2
For SiO 2 , an initial porous layer, grown by serial bi-deposition, is continued with a second layer. For the second layer deposition configurations for different incidence angles of 0 • , 30 • and 50 • with a rotating substrate are chosen. Each simulation is performed three times, to evaluate the variation of the continued growth caused by the randomized deposition algorithm. The structures are shown in Figure 8. While the lower thirds of the structures are nearly identical, the upper two thirds exhibit differences.
When looking at the density profiles of the structures in Figure 9, an improvement to the results in previous works from the authors in collaboration with the Vilnius university in [7] can be presented. The plot shows the mean and standard deviation of the density profiles evaluated for the three repeated simulations. Here, the impact of the second layer deposition onto the first layer can be better distinguished for the 30 • and 50 • cases in the growth height region between 15 and 20 nm. First layer 0° 30° Figure 8. SiO2 structures, with the same initial first zig-zag layer (left), followed by a layer with a deposition angle of 0°, 30° and 50° with rotating substrate. Graphics generated with NGL [16,17] When looking at the density profiles of the structures in Figure 9, an improvement to the results in previous works from the authors in collaboration with the Vilnius university in [7] can be presented. The plot shows the mean and standard deviation of the density profiles evaluated for the three repeated simulations. Here, the impact of the second layer deposition onto the first layer can be better distinguished for the 30° and 50° cases in the growth height region between 15 and 20 nm. The standard deviation shows that the impact and its variation is minimal for 50°, while for 30° the impact is larger and with a higher fluctuation. Due to this variation the difference for 30° and 50° was hard to determine in [7]. For the height region of the second layer (25 to 50 nm) a slightly more dense structure is obtained for 0° compared with 50°, while the density for 30° is fluctuating between both. This is also an improvement on [7], where no significant difference between all three deposition angles can be quantified.
The surface morphology and roughness are presented in Figures 10 and 11 and show a large increase of the surface roughness for the first layer (zig zag), which is followed by a decrease with the continuation of the second layers. For 0° incidence, the largest decrease in surface roughness can be reached, and for 50° the value is constantly larger, while for 30° the roughness is fluctuating between the values for 0° and 50°.
The evaluation considering the anisotropic properties of the structures is also important for obtaining an idea of the random variations. Table 2 give the results for the first layer region between 7 and 20 nm (which is impacted by the second layer) and the second layer region between 25 and 47 nm, respectively. The fill factor is determined within limits of the standard deviation between 0.6% and 0.9% for 0° and 30°. Interestingly, for 50° the impact on the first layer has a very small variation of 0.05%, which is also caused by the fact that the impact itself is small, while the variation of the second layer is significantly First layer 0° 30° Figure 8. SiO2 structures, with the same initial first zig-zag layer (left), followed by a layer with a deposition angle of 0°, 30° and 50° with rotating substrate. Graphics generated with NGL [16,17] When looking at the density profiles of the structures in Figure 9, an improvement to the results in previous works from the authors in collaboration with the Vilnius university in [7] can be presented. The plot shows the mean and standard deviation of the density profiles evaluated for the three repeated simulations. Here, the impact of the second layer deposition onto the first layer can be better distinguished for the 30° and 50° cases in the growth height region between 15 and 20 nm. The standard deviation shows that the impact and its variation is minimal for 50°, while for 30° the impact is larger and with a higher fluctuation. Due to this variation the difference for 30° and 50° was hard to determine in [7]. For the height region of the second layer (25 to 50 nm) a slightly more dense structure is obtained for 0° compared with 50°, while the density for 30° is fluctuating between both. This is also an improvement on [7], where no significant difference between all three deposition angles can be quantified.
The surface morphology and roughness are presented in Figures 10 and 11 and show a large increase of the surface roughness for the first layer (zig zag), which is followed by a decrease with the continuation of the second layers. For 0° incidence, the largest decrease in surface roughness can be reached, and for 50° the value is constantly larger, while for 30° the roughness is fluctuating between the values for 0° and 50°.
The evaluation considering the anisotropic properties of the structures is also important for obtaining an idea of the random variations. Table 2 give the results for the first layer region between 7 and 20 nm (which is impacted by the second layer) and the second layer region between 25 and 47 nm, respectively. The fill factor is determined within limits of the standard deviation between 0.6% and 0.9% for 0° and 30°. Interestingly, for 50° the impact on the first layer has a very small variation of 0.05%, which is also caused by the fact that the impact itself is small, while the variation of the second layer is significantly The standard deviation shows that the impact and its variation is minimal for 50 • , while for 30 • the impact is larger and with a higher fluctuation. Due to this variation the difference for 30 • and 50 • was hard to determine in [7]. For the height region of the second layer (25 to 50 nm) a slightly more dense structure is obtained for 0 • compared with 50 • , while the density for 30 • is fluctuating between both. This is also an improvement on [7], where no significant difference between all three deposition angles can be quantified.
The surface morphology and roughness are presented in Figures 10 and 11 and show a large increase of the surface roughness for the first layer (zig zag), which is followed by a decrease with the continuation of the second layers. For 0 • incidence, the largest decrease in surface roughness can be reached, and for 50 • the value is constantly larger, while for 30 • the roughness is fluctuating between the values for 0 • and 50 • .
The evaluation considering the anisotropic properties of the structures is also important for obtaining an idea of the random variations. Table 2 give the results for the first layer region between 7 and 20 nm (which is impacted by the second layer) and the second layer region between 25 and 47 nm, respectively. The fill factor is determined within limits of the standard deviation between 0.6% and 0.9% for 0 • and 30 • . Interestingly, for 50 • the impact on the first layer has a very small variation of 0.05%, which is also caused by the fact that the impact itself is small, while the variation of the second layer is significantly larger with 1.59%. The variation of the birefringence is very small for the first layer, however here only the impact of the second layer on to the first one is obtained. For the second layer, the variation is slightly larger, while the influence from the first layer is still significant at least for 0 • . In sum, the standard deviation is much lower than the birefringence values computed for the first layers and in the range of the changes that describe a decreased birefringence for smaller incidence angles. The number of m = 3 simulations for each angle have to be considered as a compromise between computational effort and statistical significance. Assuming a normal distribution, for calculation of a 95% confidence interval, the multiplicative factor with regard to the Student's t distribution with m-1 degrees of freedom results in being 4.30, while the factor of 1.96 would be reached for m = ∞. Therefore, the statistical significance could be doubled roughly, but with a much higher computational effort. larger with 1.59 %. The variation of the birefringence is very small for the first layer, however here only the impact of the second layer on to the first one is obtained. For the second layer, the variation is slightly larger, while the influence from the first layer is still significant at least for 0°. In sum, the standard deviation is much lower than the birefringence values computed for the first layers and in the range of the changes that describe a decreased birefringence for smaller incidence angles. The number of m = 3 simulations for each angle have to be considered as a compromise between computational effort and statistical significance. Assuming a normal distribution, for calculation of a 95 % confidence interval, the multiplicative factor with regard to the Student's t distribution with m-1 degrees of freedom results in being 4.30, while the factor of 1.96 would be reached for m = ∞. Therefore, the statistical significance could be doubled roughly, but with a much higher computational effort.     Figure 10. Surface morphology of the zig-zag structure from serial bi-deposition (first row) for 4, 7 and 14 nm mean height, and continued growth with glancing angle deposition under 0°, 30° and 50° with rotating substrate for 21, 32 and 40 nm mean height. Figure 11. Standard deviation σ of the height during growth of the first layer (zig zag) and the continued growth of the second layer with different coating conditions on the right including also the fluctuation (standard deviation) as shaded region. Figure 11. Standard deviation σ of the height during growth of the first layer (zig zag) and the continued growth of the second layer with different coating conditions on the right including also the fluctuation (standard deviation) as shaded region.

Discussion and Conclusions
The results show several aspects that have to be considered for performing and evaluation of the simulated thin film growth. The analysis of TiO 2 reveal, that the chosen interaction potential has a significant impact onto the structural properties, most notably the density and at least as a secondary effect also on the surface roughness and anisotropy. Here, the potential from Zhang et al. produces densities that are much too low compared to experimental values. While showing this already for the small structures, the effect is increased for the larger structures. This corresponds partly to the findings from Grigoriev et al. [8] (p. 7), who investigated SiO 2 thin film growth with high energetic Si, and observed a small decrease of the density for larger structures. However, our results for TiO 2 and the decrease in density for larger structures are more pronounced.
Nevertheless, the increase of the lateral structure size from 7 × 8 nm 2 to 30 × 20 nm 2 results in a significantly better description, which means much less fluctuation, of the density profiles in the growth direction and the development of the surface rms roughness. Furthermore, qualitatively the description of the rms roughness is enhanced, so that much higher values are reached and show a dependence on different angles of incidence. The small structures of 7 × 8 nm 2 cannot reproduce that. The analysis of the anisotropy is not clear. Here, further investigations have to be conducted into how large the influence of the different chosen potentials really is and which effects can be solely attributed to the structure size.
The analysis for repeated simulations of SiO 2 GLAD coatings gives results that are more reliable and also an estimate of the significance of the simulated structural properties and their difference. Therefore, repetitions of the simulations can help to gain more insights, while e.g., the columnar structure is probably not a complete representation of the real world, as learned from the analysis regarding TiO 2 structures. The size of the dimension in y with 7 nm does not allow the buildup of independent columns, instead forms of walls due to the periodic boundary conditions are grown. The density of the SiO 2 structures grown with the potential from Zhang et al. [10] results in experimentally reasonable values.
The presented data allows an estimation of the significance of structural properties depending on the structure size. While for a detailed dependence more data are needed, the results already give an important idea of how large the structures have to be for the respective property. Furthermore, the importance of the choice of the potential is emphasized, which needs to be thoroughly tested and compared to experimental findings, especially regarding the density.