Statistical information enhanced robust design method of optical thin film

Jianglin Dai; Gong Zhang; Shigeng Song; Zeyu Deng; Zihan Wei

doi:10.1364/OE.471998

1. Introduction

Since the error between the theoretical and actual values of the optical film structure para-meters is inevitable, except for the spectrum properties, the sensitivity of spectral has been considered as well in the design of optical thin films. However, the sensitivity of spectral coefficients to errors is not fully considered in the conventional optimization algorithms [1–5], which a small film thickness error could also result in great variations in production yield. A robust film system with high production yield should be designed for the practical fabrication.

At present, the single object optimization method [6–9] has been commonly applied in the robust design of optical thin films, which has only one merit function F in single object (even if this function is the sum of several sections). In the early robust design, the corresponding calculation of the first-order partial derivative of the film [9–11] acted as an additional item of the merit function. However, the accuracy of the first order approximation is significantly limited as the true parameter error cannot be considered significantly small [12], and the merit function is commonly located near the local extreme value, where the gradient is close to or equal to zero. Among the robust design methods proposed recently [13–18], Monte Carlo method [12,13] adopts the statistical experimental method of sampling simulation and selects the sample average of merit function error as the mathematical expectation. It is capable of achieving effective results for various error distributions, unfortunately, the time taken to get an accurate result required for simulation is often intolerable. The analytical approximation method [15,16] is based on the second-order Taylor expansion of tolerance type merit function. The mathematical expectation of the merit function error is calculated with the fast calculation model of the first and second-order partial derivatives of the film system. It exhibits the advantages of fast calculation speed and high precision. It is applicable to the system which the film parameter errors produced in the fabrication process complies with the zero mean normal distribution. (e.g., quartz crystal oscillator monitoring system). The limitation is that it does not consider the standard deviation of the merit function, therefore it cannot fully represent the degree of error dispersion. In 2018, Shang Qi Kuang et al. proposed a multi-objective optimization robust design method. He constructed two merit functions for broad-band extreme ultraviolet (EUV) thin films [19] to calculate the ideal performance and robustness respectively. With this method, the effectiveness of the method was verified in elevating the fabrication yield of EUV mirrors through two types of design examples of broadband Mo / Si multilayers. The above-mentioned robust design methods are mainly based on the simulation or calculation of the mathematical expectation. However, the effect of the error standard deviation on the practical production is not well considered in these methods.

In this paper we propose new merit functions for robust design from the perspective of fabrication yield, calculating the mathematical expectation and standard deviation of errors, as an attempt to improve the yield about 5%∼15% for different layers. A multi-objective optimization strategy was found to avoid the pre-allocation of weights and the conflict between the original merit function and the additional items due to single objective optimization.

2. Merit function of robust design

According to the automatic design method of optical thin films, the merit function refers to a quantitative index to determine the matching degree between the theoretical spectral coefficient and the target value, acting as the objective function in the process of optimization. The merit functions applied previously for conventional design are presented below [9]:

(1)$$F = {[\frac{1}{N}{\sum\limits_{k = 1}^N {{\omega _k}(Q({\lambda _k}) - {Q_0}({\lambda _k}))} ^n}]^{\frac{1}{n}}}.$$

Under the differences in the range of values of different spectral coefficients (e.g., $\mathrm{0\ < T\ < 1}$, $\mathrm{0\ < \varphi < 2\pi }$) and their different dimensions, the mentioned merit function can only be exploited to measure one spectral coefficient in many cases. These problems are solved with the following function

(2)$$F = {[\frac{1}{N}\sum\limits_{k = 1}^N {(\frac{{{\omega _k}(Q({\lambda _k}) - {Q_0}({\lambda _k}))}}{\Delta }} )^n}{]^{\frac{1}{n}}}.$$

where $\mathrm{\Delta }$ denotes the maximum permissible departure between $Q({\lambda _k})$ and ${Q_0}({\lambda _k})$. The presented merit function is transformed into the ratio function, unifying the order of magnitude of the merit function and achieving dimensionless. When used, the square tolerance merit function with $\textrm{n = 2}$ is generally selected.

In previous studies, the merit functions employed in robust design can fall to three types [9,15,18].

The merit functions are listed in Table 1, where m denotes the total number of layers; NN represents the total number of Monte Carlo simulations. The first type of merit function calculates the weighted sum of the first partial derivatives of different coatings by conventional merit function as an additional term of the merit function; the second type of merit function obtains the mathematical expectation of the merit function in the presence of error disturbance through Monte Carlo simulation. The third type of merit function is based on the second-order Taylor expansion of the conventional merit function. Under the assumption that the film parameter errors are independent and satisfy the zero mean normal distribution, the analytical solution of the mathematical expectation is taken as the merit function.

The accuracy of approximation of the first-order expansion employed in the first type of merit function is significantly limited, which cannot be applied to most robust designs. Though the accuracy of approximation of the second and third types of merit functions has been significantly improved, and each has its own advantages in the scope of application of the model and the calculation time, the key elements of the two types of merit functions are the mathematical expectation to determine the error of the merit function, but it can't reveal the dispersion degree of the error well. For the mentioned reason, building a merit function with complete statistical information to calculate the mathematical expectation and standard deviation can enhance the robustness of the design film system more effectively.

Given the spectral coefficient error attributed to the structural parameter error of the film,

(3)$$\Delta {Q_\lambda } = {Q_\lambda } - {Q^{\prime}_\lambda },$$

where ${Q_\lambda }$ denotes the practical spectral coefficient, ${Q^{\prime}_\lambda }$ represents the theoretical value. Based on the second order Taylor expansion [20], it yields:

(4)$$\begin{array}{l} \Delta {Q_\lambda } \approx \sum\limits_{i = 1}^m {\frac{{\partial {Q_\lambda }}}{{\partial {d_i}}}} \delta {d_i} + \sum\limits_{i = 1}^m {\frac{{\partial {Q_\lambda }}}{{\partial {n_i}}}} \delta {n_i}\\ + \frac{1}{2}\sum\limits_{i,j = 1}^m {\{ \frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {d_j}}}\delta {d_i}\delta {d_j} + \frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}\partial {n_j}}}\delta {n_i}\delta {n_j} + 2\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {n_j}}}\delta {d_i}\delta {n_j}\} } , \end{array}$$

where $\delta d$ denotes the thickness error vector; $\delta n$ represents the refractive index error vector. Due to the Hessian matrix being a symmetric matrix, thus

(5)$$\begin{array}{l} \sum\limits_{i,j = 1}^m {\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {d_j}}}\delta {d_i}\delta {d_j}} = 2\sum\limits_{i > j}^m {\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {d_j}}}\delta {d_i}\delta {d_j}} + \sum\limits_{i = 1}^m {\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}^2}}\delta {d_i}^2} \\ \sum\limits_{i,j = 1}^m {\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}\partial {n_j}}}\delta {n_i}\delta {n_j}} = 2\sum\limits_{i > j}^m {\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}\partial {n_j}}}\delta {n_i}\delta {n_j}} + \sum\limits_{i = 1}^m {\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}^2}}\delta {n_i}^2} , \end{array}$$

$$$$

For both sides of the Eq. (4), the mathematical expectation and standard deviation are obtained

(6)$$\begin{array}{l} {M_{\Delta {Q_\lambda }}} \approx \sum\limits_{i = 1}^m {\frac{{\partial {Q_\lambda }}}{{\partial {d_i}}}} E(\delta {d_i}) + \sum\limits_{i = 1}^m {\frac{{\partial {Q_\lambda }}}{{\partial {n_i}}}} E(\delta {n_i}) + \\ \frac{1}{2}\sum\limits_{i,j = 1}^m {\{ \frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {d_j}}}E(\delta {d_i}\delta {d_j}) + \frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}\partial {n_j}}}E(\delta {n_i}\delta {n_j}) + 2\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {n_j}}}E(\delta {d_i}\delta {n_j})\} } , \end{array}$$

(7)$$\begin{array}{l} {S^2}_{\Delta {Q_\lambda }} \approx \sum\limits_{i = 1}^m {\{ {{(\frac{{\partial {Q_\lambda }}}{{\partial {d_i}}})}^2}{\sigma ^2}(\delta {d_i}) + {{(\frac{{\partial {Q_\lambda }}}{{\partial {n_i}}})}^2}{\sigma ^2}(\delta {n_i})} \} + \sum\limits_{i > j}^m {\{ {{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {d_j}}})}^2}{\sigma ^2}(\delta {d_i}\delta {d_j}) + {{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}\partial {n_j}}})}^2}{\sigma ^2}(\delta {n_i}\delta {n_j})} \} \textrm{ + }\\ \textrm{ + }\frac{1}{4}\sum\limits_i^m {\{ {{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}^2}})}^2}{\sigma ^2}(\delta {d_i}^2) + {{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}^2}})}^2}{\sigma ^2}(\delta {n_i}^2)\} } + \sum\limits_{i,j = 1}^m {{{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {n_j}}})}^2}{\sigma ^2}(\delta {d_i}\delta {n_j})} , \end{array}$$

where E denotes the mathematical expectation, $\sigma $ represents the standard deviation. Now, F, M and S functions are obtained.

For the practical industrial production, in addition to the expected value of spectral coefficient that determine the actual degree of the error dispersion, there is also the standard deviation of spectral coefficient, which owns great significance for fabrication. When the errors of the film parameters are independent of each other and comply with the normal distribution of the zero mathematical expectation, most of the expectation in Eq. (6) reaches to 0, the expectation of squared error is equated with the standard deviation. the higher order moments ${\sigma ^2}(\delta {d_i}^2),{\sigma ^2}(\delta {n_i}^2)$ appear in the Eq. (7) and the statistic ${(\frac{{\delta {d_i} - \mu }}{\sigma })^2}$ follows the chi-square distribution ${\chi ^2}(1)$, hence the higher order moments can be expressed as ${\sigma ^2}(\delta {d_i}^2) = 2{\sigma ^4}(\delta {d_i}),{\sigma ^2}(\delta {n_i}^2) = 2{\sigma ^4}(\delta {n_i})$. Therefore, Eq. (6) and Eq. (7) can be given by

(8)$${M_{\Delta {Q_\lambda }}} = \frac{1}{2}\sum\limits_{i,j = 1}^m {\{ \frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}^2}}\sigma _{d,i}^2 + \frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}^2}}\sigma _{n,i}^2} \} ,$$

(9)$$\begin{array}{l} {S^2}_{\Delta {Q_\lambda }} \approx \sum\limits_{i = 1}^m {\{ (\frac{{\partial {Q_\lambda }}}{{\partial {d_i}}}} {)^2}\sigma _{d,i}^2 + {(\frac{{\partial {Q_\lambda }}}{{\partial {n_i}}})^2}\sigma _{n,i}^2\} + \sum\limits_{i > j}^m {{{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {d_j}}})}^2}\sigma _{d,i}^2\sigma _{d,j}^2} \\ \frac{1}{2}\sum\limits_{i = 1}^m {\{ {{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}^2}})}^2}\sigma _{d,i}^4 + {{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {n_i}^2}})}^2}\sigma _{n,i}^4} \} + \sum\limits_{i,j = 1}^m {{{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {n_j}}})}^2}\sigma _{d,i}^2\sigma _{n,j}^2} \\ \approx \sum\limits_{i = 1}^m {\{ (\frac{{\partial {Q_\lambda }}}{{\partial {d_i}}}} {)^2}\sigma _{d,i}^2 + {(\frac{{\partial {Q_\lambda }}}{{\partial {n_i}}})^2}\sigma _{n,i}^2\} \\ + \frac{1}{2}\sum\limits_{i,j = 1}^m {\{ {{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {d_j}}})}^2}\sigma _{d,i}^2\sigma _{d,j}^2 + } {(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {d_j}}})^2}\sigma _{d,i}^2\sigma _{d,j}^2 + 2{(\frac{{{\partial ^2}{Q_\lambda }}}{{\partial {d_i}\partial {n_j}}})^2}\sigma _{d,i}^2\sigma _{n,j}^2\} . \end{array}$$

Here the spectral coefficient error conforms to the normal distribution. In accordance with the statistical principle, the actual spectrum will be located in the standard deviation envelope with 68.26% probability and in the triple standard deviation envelope with 99.74% probability.

For other distribution of film parameter error, if some specific statistical parameters are known, the mathematical expectation and standard deviation of spectral coefficient increment can be determined by the mentioned formula as well by using Eq. (7), and the statistical distribution parameters under the corresponding distribution form can be subsequently obtained. For independent error distribution, it requires mathematical expectation and standard deviation; for correlated error distribution, it requires other high-order statistics.

According to the expectation and standard deviation of spectral coefficients, we propose two novel merit functions

(10)$$\begin{array}{l} {f_1} = \frac{1}{N}\sum\limits_{k = 1}^N {{{(\frac{{Q({\lambda _k}) + {M_{\Delta {Q_{{\lambda _{_k}}}}}} - {Q_0}({\lambda _k})}}{\Delta })}^2}} ,\\ {f_2} = \frac{1}{N}\sum\limits_{i = 1}^N {{{(\frac{{Q({\lambda _k}) + {M_{\Delta {Q_{{\lambda _{_k}}}}}} \pm a \bullet {S_{\Delta {Q_{{\lambda _k}}}}} - {Q_0}({\lambda _k})}}{\Delta })}^2}} . \end{array}$$

In Eq. (10), the coefficients a can be selected according to different statistical probabilities, which is corresponding to the error envelope of different probabilities. For normally distributed error, the coefficient is supposed to be set to 3, corresponding to the envelope with 99.74% probability. For other distributions, the coefficient a corresponding to a specific probability P is supposed to be determined according to the function $P = \int_{ - a}^a {\phi (x )} dx$, where $\phi (x )$ denotes the probability density function. The signs can be selected given the type of the designed film system. For instance, if the design target is the transmission spectrum T of antireflection film, we hope the lower limit of T could be higher, therefore the sigh is supposed to be selected as “-”, and that of high reflection film is “+”. For the case of a beamsplitter with a requirement of T = 50%, both equations using “+” and “-” signs respectively are supposed to be used. Since this design method is based on M.S functions, we name it M-S Robust method. With the fast analytical calculation method of the first and second order partial derivatives of the film parameters [21], the first-order partial derivative matrix and Hesse matrix of the presented formula can be calculated efficiently and accurately. Furthermore, Hesse matrix mostly calculates the continuous product of matrix

(11)$$A = {M_m}\ldots {M_{k + 1}}{D_k}{M_{k - 1}}\ldots {M_{j + 1}}{D_j}{M_{j - 1}}\ldots {M_1},$$

where ${M_i}$ denotes the matrix of the i-layer film; ${D_i}$ represents the first-order partial derivative matrix of the i-layer film.

Matrix A can be rewritten as

(12)$$A = O{D_k}P{D_j}Q.$$

where

(13)$$\begin{array}{l} O = {M_m}{M_{m - 1}}\ldots {M_{k + 1}},\\ P = {M_{k - 1}}{M_{k - 2}}\ldots {M_{j + 1}},\\ Q = {M_{j - 1}}{M_{j - 2}}\ldots {M_1}. \end{array}$$

Nested double cycle [21] is adopted to implement the calculation program: construct a four-dimensional matrix B, its elements $B(:,:,i,j) = {M_i}{M_{i - 1}}\ldots {M_j}$ are adopted to record the continuous products of the film matrix, which can be reused in the calculation to implement the fast calculation of Hesse matrix.

In the calculation of spectral coefficients of optical thin films, according to the matrix calculating method, the main part of the algorithm is a multiplication of two-dimensional matrices. An optical film with m layers requires m times multiplication in the calculation of theoretical spectrum. For matrix B, calculation of each element requires only one-time multiplication, however, for the entire B matrix, calculation requires m(m+1)/2 times multiplications, this is due to the symmetry of Hesse matrix. When calculating the first-order and second-order partial derivatives using the B matrix, matrix multiplication calculation of 2 m and 2 m(m+1) times is required respectively. Therefore, the calculation time of ${f_1}$ and ${f_2}$ is about 3 m times that of theoretical spectrum, which is suitable for most calculation. Since the elements with different indexes i are independent of each other when calculating matrix B and Hessian matrix, parallel computing can be used to accelerate the algorithm.

3. Multi-objective optimization strategy

Multi-objective optimization strategy in the present study avoids the pre-allocation of weights and the conflict between the original merit function and the additional items [21] due to single objective optimization. The multi-objective optimization strategy proposed in the present study first uses Non-Dominated Sorting Genetic Algorithm (NSGA-II) for initial optimization, and stops when reaching the preset maximum generation number or satisfies the terminating condition. Then, we transform the problem into an equivalent multi-objective programming problem, and the minimum maximum method based on Sequential quadratic programming (SQP) is used for local refined optimization.

Fig. 1. Flow diagram of NSGA-II algorithm

Download Full Size | PDF

Fig. 2. The theoretical transmission spectrum T of the (a) visual band and (b) infrared band antireflection coating at normal incidence, the mathematical expectation of transmission spectra T + M△T and standard deviation envelope T + M△T ± S△T under the standard deviation of film geometric thickness error of 1%.

Download Full Size | PDF

NSGA-II [22,23] is a multi-objective optimization algorithm based on the genetic algorithm. Figure 1 shows the flow diagram of NSGA-II algorithm. It introduces elite retention strategy and fast non-dominated sorting strategy to NSGA, which can obtain Pareto frontier with multiple solutions:

(1) Initialization. The initial population A₀ with n individuals is created, and the offspring population B₀ is generated by the genetic algorithm (crossover and mutation), and the two populations and population are merged into P_t;
(2) Fast non-dominated sorting of A. The crowding degree of the individuals in each non-dominated set is calculated. According to the non-dominated relationship and the crowding degree of the individuals, the tournament is conducted among the individuals. The winner is selected to form a novel population, and the next cycle process acts as the parent population;
(3) Repeat step (2) until having reached the preset maximum generation number or satisfied the terminating condition.

SQP [24,25] acts as a significantly effective method to process small and medium-sized nonlinear programming. It transforms the problem to be optimized into a series of relatively simple quadratic programming subproblems and solves them. The advantage of this method is that it has local superlinear convergence. To implement multi-objective optimization, the maximum value of the two objective merit functions is optimized with the maximum minimum method to minimize the maximum value.

(14)$$F(x) = \min \mathop {ma}\limits_i x({f_i}(x)).$$

4. Designs

In our study, we independently write the program of merit function for robust design and use gamultiobj and fminimax functions involved in MATLAB optimization toolbox to implement the above algorithm flow. The function gamultiobj is based on NSGA-II, it reserves the non-dominated set partially for the next generation in the elite retention and the default reservation number reaches 0.35. The function fminimax is optimized by the minimum maximum method by complying with SQP.

4.1 Design of visible and infrared dual band anti-reflection coatings

This design was one of the topics of the 1988 optical film design competition [26]. Its requirements are elucidated as follows: when designing antireflection coating in 420 - 680 nm and 10500 - 10700 nm, refractive index for substrate is 2.6 for the 400-700 nm band and 2.4 for the 10000-11000 nm band. The refractive index of high and low refractive index materials are respectively 2.32, 1.52 and 2.22, 1.43 in two wavebands. The thickness of the film should be at least 11 nm, the number of layers should not exceed 90 and the thickness error of the film should be 1%, while the refractive index error is not considered.

Suyong Wu et al. had presented the results of robust design based on second-order Taylor expansion calculation [27]. The design result refers to D = {S | 54.13H 28.04L 11.33H 146.44L 11.00H 54.37L 11.00H 166.27L 11.00H 55.90L 11.36H 155.68L 11.00H 58.11L 11.00H 153.77L 14.20H 54.68L 11.00H 143.46L 35.23H 22.66L 40.69H 110.51L 15.26H 32.49L 61.13H 85.51L | A}.

In the robust design, for rigorous comparison, the same film layer number (28) is employed as their design, the sampling points number of 420 - 680 nm and 10500 - 10700 nm reached 27 and 11, the spectral tolerance takes up 0.1%, the lower limit of film thickness is 11 nm, the upper limit is 1.5D, where D denotes the thickness of Wu’s design; the population size is 200, the crossover probability is 0.9, and the mutation probability reaches 0.1. In this scale problem, the computation time of each generation is about 0.79s, and the total time is about 6.5 minutes. The result of our robust design is {S | 53.69H 27.23L 11.39H 148.63L 11.00H 53.59L 11.00H 167.47L 11.00H 55.32L 11.24H 156.06L 11.00H 57.35L 11.06H 154.48L 13.95H 53.50L 11.00H 146.13L 34.99H 22.97L 39.78H 112.21L 15.24H 31.57L 61.59H 85.64L | A}. The main spectral performance parameters of the two results are listed in Table. 2.

4.2 Design of robust cut-off filter

We design a cut-off filter with 100% transmittance in 400-650 nm band and 0% transmittance in 700-900 nm band at normal incidence. The refractive index of incident medium and substrate is 1.0 and 1.52, respectively. The high and low refractive index materials are Ti₃O₅ and SiO₂, respectively. Considering material dispersion, the fitted Cauchy dispersion coefficient is listed in Table 3. It is assumed that the geometric thickness errors of all coatings comply with the normal distribution of 1% standard deviation and the error in refractive index is not considered.

Table 1. Merit functions employed in robust design

View Table | View all tables in this article

Table 2. Spectral performance of visible and infrared dual band anti-reflection coatings (%)

View Table | View all tables in this article

Table 3. Cauchy dispersion coefficient of Ti3O5 and SiO2

View Table | View all tables in this article

In the design, the target points of 400 - 650 nm and 700 - 900 nm are 26 and 21, respectively and the spectral accuracy coefficient is 1%. With the conventional design method, the initial structure of the film system refers to a stack of 21 layers with quarter wavelength L and h, and the reference wavelength reaches 700 nm. Combined with needle method optimization, the design result of 27 layers is achieved: D = {S | 135.66L 8.77H 30.36L 92.10H 10.74L 3.42H 142.13L 86.15H 143.07L 82.47H 140.30L 81.16H 139.13L 80.63H 138.70L 80.53H 138.79L 80.82H 139.55L 81.73H 141.42L 84.28H 146.33L 20.67H 2.20L 63.59H 74.30L | A}. In the robust design, the same film layer number (27) is employed, the lower limit of film thickness is 0.5D, the upper limit is 1.5D, where D denotes the thickness of conventional design. The population size is 200 and the crossover probability is 0.9, the mutation probability reaches 0.1. The calculation time of the respective generation population is about 0.74s, and the total time is about 6 minutes. For the significantly small thickness of the needle layer, the thickness standard deviation is set to 0.5 nm. The structure obtained by the robust design method is {S | 114.36L 8.15H 34.04L 82.65H 2.19L 13.23H 151.23L 85.55H 143.09L 82.33H 140.43L 81.07H 139.41L 80.60H 139.15L 80.60H 139.27L 81.05H 140.16L 81.98H 142.51L 84.76H 148.68L 29.49H 1.73L 56.04H 75.65L | A}. Table 4 shows the theoretical reflectivity and transmittance of the two structures. In order to verify the effectiveness of robust design, we utilized Monte Carlo simulation method to simulate the two structures for 20000 times and counted the fabrication yield of each band. The standard deviation of film thickness is set to 1%. Figure 3 shows the theoretical transmission spectrum. Figure 4 and Fig. 5 show the quartiles of the simulation results and the spectral envelope of 95% of the samples. Table 5 shows the fabrication yield of each band.

Fig. 3. (a) Transmission band. (b) Stop band theoretical transmission spectrum of conventional and robust designs of cut-off filter in 400-650 nm and 700-900 nm band.

Download Full Size | PDF

Fig. 4. (a) Transmission band. (b) Stop band result of 20000 Monte-Carlo simulations: the quartiles and the spectral envelope of 95% of the samples of conventional design of cut-off filter in 400-650 nm and 700-900 nm band.

Download Full Size | PDF

Fig. 5. (a) Transmission band. (b) Stop band result of 20000 Monte-Carlo simulations: the quartiles and the spectral envelope of 95% of the samples of M-S Robust design of cut-off filter in 400-650 nm and 700-900 nm band.

Download Full Size | PDF

Table 4. Theoretical transmittance and reflection of cut-off filter

View Table | View all tables in this article

Table 5. Cut-off filter fabrication yield of each band

View Table | View all tables in this article

4.3 Design of robust infrared dual-band filter

We designed an infrared dual-band filter with 100% transmittance in 3160-3460nm and 7570-7770nm band and 0% transmittance in 3700-6000nm band at normal incidence The refractive index of incident medium and 1.0. The low refractive index material is ZnS, the high refractive index material and structure are both Ge. The fitted Cauchy dispersion coefficient is listed in Table 6. It is assumed that the geometric thickness errors of all coatings also comply with the normal distribution of 1% standard deviation and the error of refractive index is not considered.

Table 6. Cauchy dispersion coefficient of Ge and ZnS

View Table | View all tables in this article

In the design, each band has 50 target points and the spectral accuracy coefficient is 1%. With the conventional design method, the result of 26 layers is achieved according to Liu’s [28] research: D={S |697.77L 42.43H 61.70L 141.78H 247.46L 420.83H 185.35L 265.13H 541.27L 261.92H 484.77L 248.16H 509.98L 278.93H 530.78L 425.26H 737.93L 230.69H 713.79L 465.56H 107.81L 105.06H 103.13L 214.23H 353.73L 2.7H | A}. In the robust design, the same film layer number (26) was employed, he lower limit of film thickness is 0.5D, the upper limit is 1.5D, where D denotes the thickness of conventional design. The population size is 200 and the crossover probability is 0.85, the mutation probability reaches 0.15. The calculation time of the respective generation is nearly 0.71s, and the total time is about 6 minutes. The structure obtained by the robust design method is {S | 691.50L 50.53H 118.70L 86.79H 341.63L 370.06H 171.24L 327.94H 492.8L 247.23H 501.63L 257.06H 514.82L 272.31H 569.08L 407.93H 666.10L 276.32H 656.93L 421.03H 94.52L 106.88H 168.63L 198.29H 272.5L 32.31H | A}. Table 7 shows the theoretical reflectivity and transmittance of the two structures. In effort to confirm the effectiveness of robust design, we also utilized Monte Carlo simulation method to simulate the two structures for 20000 times and counted the fabrication yield of each band. The standard deviation of film thickness is also set to 1%. Figure 6 and Fig. 7 show the quartiles of the simulation results and the spectral envelope of 95% of the samples. Table 8 shows the fabrication yield of each band.

Fig. 6. Results of 20000 Monte Carlo simulations: the quartiles and the spectral envelope of 95% of the samples of conventional design of infrared dual-band filter in (a) 3160-3460 nm, (b) 3700-6000 nm, (c) 7570-7770 nm band.

Download Full Size | PDF

Fig. 7. Results of 20000 Monte Carlo simulations: the quartiles and the spectral envelope of 95% of the samples of M-S Robust design of infrared dual-band filter in (a) 3160-3460 nm, (b) 3700-6000 nm, (c) 7570-7770 nm band.

Download Full Size | PDF

Table 7. Infrared duan-band filter fabrication yield of each band

View Table | View all tables in this article

Table 8. Infrared dual-band filter fabrication yield of each band

View Table | View all tables in this article

Table 9. Calculation time of robust designs (population size = 200)

View Table | View all tables in this article

Note, in the above design, the variable range in NSGA-II algorithm was set according to the results of traditional design, i.e. using the traditional design as the reference point. This results in similarities between the two designs. However, the design obtained by the robust design method demonstrated improved performance as compared to traditional designs. And the calculation times of robust designs are listed in Table 9.

5. Analysis

In the design of dual band antireflection coatings, Fig. 2 illustrates the theoretical spectrum of the film system in Table 1, the mathematical expectation of the spectrum in the presence of error and the envelope of standard deviation (the actual spectrum falls in it with 68.26% probability) under the assumption of the relative error of the film geometric thickness of the normal distribution type with the standard deviation of 1%. As indicated from the comparison of the data in Table 1 and Fig. 2, the mathematical expectation and standard deviation envelope of spectrum in 420 - 550 nm band are slightly worse than those of Wu, whereas the range of 550 - 680 nm is better than that of Wu. The average value of the whole wavelength range is higher than that of Wu and the average value is elevated by 0.0011% and 0.0018%, respectively. In the range of 10500 - 10700 nm, it is better than that of Wu and the average value increased by 0.0031% and 0.0034%, respectively. The spectral sensitivity of the two bands is lower than that of Wu, demonstrating that the results achieved with the proposed method are more robust.

In the design of cut-off filter, Fig. 4 and Fig. 5 show the spectrums of the two design results. It can be seen that the spectral area of interest: cut-off band, is very similar between the two design methods. In particular, the robust design method achieves calculation results in much reduced time with high fabrication yield. However, the ripples seen for robust design method is more pronounced compared to conventional design method, this is a comprise to achieve higher fabrication yield and more efficient calculation while maintaining maximum optimization near cut-off band. The data in Table 5 shows that in the 20000 Monte Carlo simulations, the fabrication yield of robust design results has increased by about 5% on average in the transmission band and about 2% in the stop band.

In the design of infrared dual-band filter, the spectral performance of the conventional design is similar to that of the robust design, however the robust design has better performance in the Monte Carlo simulations. The data in Table 8 shows that in 20000 Monte Carlo simulations, the fabrication yield of robust design results has increased by about 9.5% on average in the transmission band, especially in 3160-3460nm band, the fabrication yield that the transmission of each wavelength is better than 98% is improved by 12%. Figure 6 and Fig. 7 show the data of each wavelength in detail. In the range of 3400-3460nm, the quartile of transmittance increased by 2%. In the range of 7570-7770nm, the quartile of reflectivity decreased by 0.3%. Although the fabrication yield of cut-off band decreased by about 1%, from 99.99% to 98.68%, the fabrication yield of other bands increased significantly. This shows that robust design method has a good application in the design of infrared dual-band filter.

The above analysis shows that the robust design method proposed in this paper can effectively improve the fabrication yield of various thin films. However, the improvement effect is limited for films whose spectrum has a band that is significantly sensitive to errors, such as cut-off filters. For this kind of thin films, the weight setting at different wavelengths is necessary to be reconsidered.

6. Conclusion

In accordance with the principle of statistics, two types of novel assessment functions including the complete discrete information of spectral coefficients are built. A multi-objective optimization strategy based on NSGA-II and SQP is proposed to simultaneously optimize the mathematical expectation of spectral coefficient error and the standard deviation envelope. In addition, the robustness of the design results is verified by performing the design experiments of visible infrared dual band antireflection film, cut-off filter and infrared dual-band filter. As demonstrated by the results, the dual band antireflection coating developed with this robust design method exhibits stronger robustness and the spectral error envelope fluctuation of the cut-off filter is smaller. Furthermore, the method is easy to implement by software and capable of effectively improving the production yield of the film system with low calculation time and space cost.

Funding

111 Project (D21009); National Natural Science Foundation of China (11973040).

Disclosures

The authors declare no conflicts of interest.

Data Availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

References

1. A. V. Tikhonravov and M. K. Trubetskov, “Modern design tools and a new paradigm in optical coating design,” Appl. Opt. 51(30), 7319–7332 (2012). [CrossRef]

2. A. Tikhonravov and M. Trubetskov, “Modern status and prospects of the development of methods of designing multilayer optical coatings,” J. Opt. Technol. 74(12), 845–850 (2007). [CrossRef]

3. A. V. Tikhonravov, M. K. Trubetskov, T. V. Amotchkina, and A. Thelen, “Optical coating design algorithm based on the equivalent layers theory,” Appl. Opt. 45(7), 1530–1538 (2006). [CrossRef]

4. A. V. Tikhonravov, M. K. Trubetskov, and G. W. DeBell, “Application of the needle optimization technique to the design of optical coatings,” Appl. Opt. 35(28), 5493–5508 (1996). [CrossRef]

5. L. Li and J. Dobrowolski, “Computation speeds of different optical thin-film synthesis methods,” Appl. Opt. 31(19), 3790–3799 (1992). [CrossRef]

6. M. Trubetskov, “Design of Multilayer Coatings Using Deep Search Methods,” in Optical Interference Coatings(Optical Society of America2019), p. TC. 2.

7. A. Deliwala, “Automated Design of Optical Thin Films via Statistical Inference and Parallelized Computation,” in Optical Interference Coatings(Optical Society of America2019), p. TD. 5.

8. J. R. Birge, “Thin-film design refinement via efficient parallel combinatorial search,” in Optical Interference Coatings(Optica Publishing Group2016), p. TC. 2.

9. J. Dobrowolski, F. Ho, A. Belkind, and V. Koss, “Merit functions for more effective thin film calculations,” Appl. Opt. 28(14), 2824–2831 (1989). [CrossRef]

10. A. V. Tikhonravov, M. K. Trubetskov, and T. V. Amotchkina, “Theoretical notes on one magic reflectance value,” in Optical Interference Coatings(Optica Publishing Group2007), p. WB3.

11. A. V. Tikhonravov, M. K. Trubetskov, T. V. Amotchkina, and A. A. Tikhonravov, “Application of advanced optimization concepts to the design of high quality optical coatings,” in 19th Congress of the International Commission for Optics: Optics for the Quality of Life(SPIE2003), pp. 1061–1062.

12. H. A. Macleod and H. A. Macleod, Thin-film optical filters (CRC press, 2010).

13. Y.-J. Jen, M.-J. Lin, and Z.-H. Yu, “Optimized Angular Insensitive Filter by Admittance Tracing Method,” in Optical Interference Coatings(Optical Society of America2019), p. TC. 7.

14. T. Amotchkina, U. Brauneck, A. Tikhonravov, and M. Trubetskov, “Sensitivity-directed refinement for designing broadband blocking filters,” Opt. Express 23(5), 5565–5570 (2015). [CrossRef]

15. S. Wu, X. Long, and K. Yang, “Novel robust design method of multilayer optical coatings,” Int. Scholarly Res. Not. 2012 (2012).

16. W. Suyong, L. XingWu, and Y. Kaiyong, “Application of thin film errors sensitivity control concept to robust design of multilayer optical coatings,” High power laser&article beams 24(10), 2391–2399 (2012). [CrossRef]

17. V. Pervak, M. Trubetskov, and A. Tikhonravov, “Robust synthesis of dispersive mirrors,” Opt. Express 19(3), 2371–2380 (2011). [CrossRef]

18. J.-F. Tang and Q. Zheng, “Automatic design of optical thin-film systems—merit function and numerical optimization method,” J. Opt. Soc. Am. 72(11), 1522–1528 (1982). [CrossRef]

19. S.-Q. Kuang, X.-P. Gong, and H.-G. Yang, “Robust design of broadband EUV multilayer using multi-objective evolutionary algorithm,” Opt. Commun. 410, 805–810 (2018). [CrossRef]

20. S. A. Furman and A. V. Tikhonravov, Basics of optics of multilayer systems (Atlantica Séguier Frontieres, 1992).

21. S. Wu, X. Long, and K. Yang, “Accurate calculation and Matlab based fast realization of merit function's Hesse matrix for the design of multilayer optical coating,” Optoelectron. Lett. 5(5), 359–363 (2009). [CrossRef]

22. M. T. Jensen, “Reducing the run-time complexity of multiobjective EAs: The NSGA-II and other algorithms,” IEEE T. Evolu. Comput 7(5), 503–515 (2003). [CrossRef]

23. K. Deb, S. Agrawal, A. Pratap, and T. Meyarivan, “A fast elitist non-dominated sorting genetic algorithm for multi-objective optimization: NSGA-II,” in International conference on parallel problem solving from nature(Springer2000), pp. 849–858.

24. P. E. Gill, V. Kungurtsev, and D. P. Robinson, “A stabilized SQP method: superlinear convergence,” Mathematical Programming 163(1-2), 369–410 (2017). [CrossRef]

25. Z. Zhu, “A simple feasible SQP algorithm for inequality constrained optimization,” Applied mathematics and computation 182(2), 987–998 (2006). [CrossRef]

26. P. Baumeister, “Evaluation of the solutions for two design problems presented at the 1998 Optical Interference Coatings Conference,” Appl. Opt. 39(13), 2230–2234 (2000). [CrossRef]

27. W. Suyong and L. Xingwu, “Research on Some Key Technologies for Robust Design, Parameters Characterization and Reverse Engineering of Thin Film Optical Coatings,” Doctoral Thesis (National University of Defense Technology, Changsha Hunan P.R.China,2011), Page 68

28. L. Wenqi and F. Xiuhua, “Research and Fabrication of Infrared Dual-Band Filter in Detection System of Methane,” Master’s Thesis (Changchun University of Science and Technology,2019), Page 20–26

	Merit function	Key elements	Advantages
1	$F F = F + \sum_{i = 1}^{m} ω_{i} \| \frac{\partial F}{\partial d} \|$	First order partial derivative of spectral coefficient	Both approximate and analytical calculations can be realized quickly
2	$F F = \frac{1}{N N} \sum_{i = 1}^{N} ω_{i} F (d + Δ d)$	Monte Carlo simulation	It’s suitable for robust design under various error distributions
3	$F F = F + \frac{1}{2} \sum_{i = 1}^{m} (\frac{\partial^{2} F}{\partial {d_{i}}^{2}} {σ_{d i}}^{2} + \frac{\partial^{2} F}{\partial {n_{i}}^{2}} {σ_{n i}}^{2})$	Second order Taylor expansion of merit function	The analytical calculation has high precision and high speed

	420-680nm			10500-10700nm
	${\bar{R}}_{t h e o r y}$	${\bar{R}}_{e r r o r}$	${\bar{R}}_{e r r o r} - S_{Δ R}$	${\bar{R}}_{t h e o r y}$	${\bar{R}}_{e r r o r}$	${\bar{R}}_{e r r o r} - S_{Δ R}$
Wu’s desgin	0.0806	0.1048	0.1661	0.0382	0.0386	0.0411
M-S Robust design	0.0799	0.1037	0.1643	0.0351	0.0355	0.0377

$n (λ) = A + \frac{B}{λ^{2}} + \frac{C}{λ^{4}}$	A	B	C
SiO₂	1.426	9.377e-3	−4.489e-4
Ti₃O₅	2.196	5.4834e-2	3.401e-5

	400-650nm			700-900nm
	$\bar{T}$	$T_{min}$	$T_{max}$	$\bar{R}$	$R_{min}$	$R_{max}$
Conventional	99.940%	99.885%	99.980%	99.709%	98.074%	99.948%
M-S Robust design	99.927%	99.783%	99.999%	99.709%	98.052%	99.948%

		Conventional	M-S Robust
400-650nm	All points T > 98%	86.81%	91.01%
400-650nm	All points T > 99%	46.66%	51.21%
700-900nm	All points T < 2.5%	90.92%	91.7%

Statistical information enhanced robust design method of optical thin film

Abstract

1. Introduction

2. Merit function of robust design

3. Multi-objective optimization strategy

4. Designs

4.1 Design of visible and infrared dual band anti-reflection coatings

4.2 Design of robust cut-off filter

4.3 Design of robust infrared dual-band filter

5. Analysis

6. Conclusion

Funding

Disclosures

Data Availability

References

Data Availability

Cited By

Figures (7)

Tables (9)

Equations (15)

Optics Express

		Conventional	M-S Robust
3160-3460nm	All points T > 98%	62.70%	76.50%
3700-6000nm	All points T < 3%	99.99%	98.68%
7570-7770nm	All points T > 99%	88.89%	93.75%

$n (λ) = A + \frac{B}{λ^{2}} + \frac{C}{λ^{4}}$	A	B	C
Ge	4.229	0.0750	0.0337
ZnS	2.183	0	0

$n (λ) = A + \frac{B}{λ^{2}} + \frac{C}{λ^{4}}$	A	B	C
Ge	4.229	0.0750	0.0337
ZnS	2.183	0	0