Physics-Informed Neural Systems for the Simulation of EUV Electromagnetic Wave Diffraction from a Lithography Mask

P H Y S I C S - I N F O R M E D N E U R A L S Y S T E M S F O R T H E S I M U L A T I O N O F E U V E L E C T R O M A G N E T I C W A V E D I FF R A C T I O N F RO M A L I T H O G R A P H Y M A S K V asiliy A. Es’kin ∗ , Egor V . Ivanov Department of Radiophysics, Uni versity of Nizhn y Novgorod 23 Gagarin A v e., Nizhny No vgorod 603022, Russia vasiliy.eskin@gmail.com , iev90078@gmail.com March 18, 2026 A B S T R A C T Physics-informed neural networks (PINNs) and neural operators (NOs) for solving the problem of diffraction of Extreme Ultra violet (EUV) electromagnetic wa ves from contemporary lithography masks are presented. A novel hybrid W av eguide Neural Operator (WGNO) is introduced, based on a wa ve guide method with its most computationally expensi ve components replaced by a neural network. T o e valuate performance, the accuracy and inference time of PINNs and NOs are compared against modern numerical solvers for a series of problems with known exact solutions. The emphasis is placed on in vestigation of solution accurac y by considered artificial neural systems for 13.5 nm and 11.2 nm wav elengths. Numerical experiments on realistic 2D and 3D masks demonstrate that PINNs and neural operators achie ve competiti v e accuracy and significantly reduced prediction times, with the proposed WGNO architecture reaching state-of-the-art performance. The presented neural operator has pronounced generalizing properties, meaning that for unseen problem parameters it delivers a solution accuracy close to that for parameters seen in the training dataset. These results provide a highly ef ficient solution for accelerating the design and optimization workflo ws of next-generation lithography masks. Key words Deep Learning · Physics-informed Neural Netw orks · Neural Operator · Predicti ve modeling · Computational physics · Extreme Ultraviolet Lithography 1 Introduction Extreme ultraviolet lithography is a critical component in the de velopment of semiconductor manuf acturing processes. This technology enables the creation of smaller and more adv anced semiconductor chips by utilizing shorter wa velengths of light to etch intricate designs onto silicon wafers. Through this method, manufacturers are able to continue scaling down the size of transistors in accordance with Moore’ s La w . Integrated circuits are made using EUV lithography , a process analogous to photographic printing, in which the patterns that will become layers of an integrated circuit are exposed on a semiconductor w afer , one layer at a time. The light pattern on the wafer is formed due to the reflection of EUV electromagnetic beam from the mask (see Figure 1). Modern integrated circuits require more than 80 different masks in their production. Due to diffraction and interference phenomena, the patterns on the mask do not match the light pattern on the wafer (see Figure 1 (a)). T o obtain the desired pattern on the wafer , a multi-stage technology of optical proximity correction (OPC) of mask is used (Figure 1 (b)). One of the stages of optical approximation correction is calculating electromagnetic fields in the area of the mask location. The diffracted electromagnetic wa ves can be calculated rigorously by using electromagnetic (EM) simulators. Howe ver , these calculations require massi ve computational resources. ∗ Corresponding author: V asiliy Alekseevich Es’kin (v asiliy .eskin@gmail.com) A preprint T o accelerate the EM simulations, various approximation models ha ve been proposed, including the “domain decompo- sition method” and the “M3D filter” [1]. These models break down a mask pattern into tw o-, one-, and zero-dimensional patterns. In these models, the electromagnetic (EM) field of a mask pattern is calculated by combining the EM fields from 2D, 1D, and 0D components. These models are currently employed in numerous EUV lithography simulation tools [1]. Howe ver , these models are still very resource-intensiv e and do not take into account the nonlocality of electromagnetic interaction. Recently , many attempts hav e been made to simulate the 3D effects of masks using deep neural networks (DNNs) such as con v olutional neural networks (CNNs) [2], generati ve a dversarial network (GAN) or U-Net [3, 4]. Howe ver , these neural network approaches are based on supervised learning, require fairly e xtensi ve datasets, hav e a significant training time, and often do not demonstrate the necessary accuracy of the solution and the de gree of its generalization. In this paper , we present physics-informed neural networks [5 – 7] and neural operators [8, 9] for solving the problem of diffraction of extreme ultraviolet electromagnetic wa ves from 3D masks. The training of these neural systems is performed in an unsupervised manner , lev eraging the gov erning physical equations directly in the training process. W e ev aluate the performance and accuracy of the neural network systems on benchmark problems and on realistic 2D and 3D mask models for current industrial lithography systems (13.5 nm) [10] and of promising ones (11.2 nm) [11, 12]. The paper is structured as follo ws. In section 2, physical problem is formulated and basic equations are gi ven. Section 3 includes the description of solution methods, which are the finite element method, wa v eguide method, methods based on physics-informed neural networks and neural operators. Section 4 presents the e valuation of these methods on problems with exact solutions, as well as numerical experiments applying these methods to simulate electromagnetic wa ve dif fraction from realistic masks. Finally , in Section 5 concluding remarks are given. Figure 1: EUV lithography with (a) and without (b) optical proximity correction of mask. 2 F ormulation of the problem and basic equations Consider a layered mask of thickness D , consisting of J layers, located in free space and enclosed in an interv al [ − D , 0] along the z -axis of a Cartesian coordinate system ( x, y , z ), as shown in Fig. 2 (a). The mask structure is periodic along the x and y axes with periods L x and L y , respectiv ely . Each j th layer of the mask is filled with a medium, which is uniform in the z direction and has a dielectric permitti vity of the j -th layer is ε j = ε j ( x, y ) . The dielectric permittivities of the media were obtained from experimental data published in [13, 14]. W e consider the dif fraction of an electromagnetic monochromatic plane w av e with angular frequenc y ω from a mask (see Fig. 2(a)). The magnetic fields in the incident w av e are giv en, with exp( iω t ) time dependence dropped, by H ( i ) = H 0 exp[ − i ( k 0; x x + k 0; y y − k 0; z z )] , where H 0 is the magnetic field amplitude, k 0; x , k 0; y , and k 0; z are the components of the wav e vector k 0 in free space ( k 0 =  k 2 0; x + k 2 0; y + k 2 0; z  1 / 2 , k 0 = ω /c , where c is the speed of light 2 A preprint Figure 2: Geometry of the problem in free space), and the superscript ( i ) denotes the incident w av e. T o simplify further consideration of the problem, we assume that the field of incident wa ve is periodic along the x and y axes with periods L x and L y , respectiv ely . W e assume that the media of the mask are time-independent and nonmagnetic (the magnetic permeability µ = 1 ). Maxwell equations for a monochromatic field in the Gaussian system of units are written as follows: rot E = − ik 0 H , rot H = ik 0 ε E , div ( ε E ) = 0 , div H = 0 . (1) It can be shown that the vector components E and H of electromagnetic field in a homogeneous media along the z axis within each j th layer are expressed in terms of magnetic components H x and H y which, in turn, satisfy the follo wing system of equations: ∆ H x + k 2 0 ε j H x + 1 ε j ∂ ε j ∂ y  ∂ H y ∂ x − ∂ H x ∂ y  =0 , ∆ H y + k 2 0 ε j H y + 1 ε j ∂ ε j ∂ x  ∂ H x ∂ y − ∂ H y ∂ x  =0 . (2) The problem is solved within a computational domain subject to boundary conditions. Periodicity conditions are satisfied on the lateral boundaries H ( − L x / 2 , z ) = H ( L x / 2 , z ) , H ( − L y / 2 , z ) = H ( L y / 2 , z ) . (3) At the top ( z = z max ) and bottom ( z = z min ) of the domain, Sommerfeld radiation conditions are applied to ensure that scattered wa ves propagate outwards without unphysical reflections. These are formulated for every spatial F ourier 3 A preprint harmonic of the field as follow: L x / 2 Z − L x / 2 L y / 2 Z − L y / 2  ∂ H ⊥ ∂ z ± ik z ; mn H ⊥  e i ( κ x mx + κ y ny ) dxdy = 0 , (4) where “ + ” and “ − ” correspond to z = z max and z = z min , respectively , k z ; mn = ( k 2 0 − κ 2 x m 2 − κ 2 y n 2 ) 1 / 2 (branch Re ( k z ; mn ) ≥ 0 and Im ( k z ; mn ) ≤ 0 is taken), κ x = 2 π /L x , κ y = 2 π /L y , m and n ∈ Z . 3 Solution Methods 3.1 Finite Element Method (FEM) The Finite Element Method is a powerful and general numerical technique for solving partial differential equations. The core principle of FEM inv olv es discretizing the computational domain into a mesh of smaller , simpler subdomains called finite elements. W ithin each element, the unknown physical field is approximated by a set of basis functions. Assembling the equations for all elements leads to a large global system of linear equations, which is then solved to obtain the approximate solution across the entire domain. While robust, FEM can be computationally expensi ve for problems like EUV diffraction, which require very fine meshes to resolve the small wa velength features. In our comparativ e studies, we utilize the open-source solver FreeFem++ [15]. 3.2 W av eguide Method The W a ve guide (WG) method serves as our high-fidelity reference solv er [2, 16, 17]. It leverages the layered structure of the mask. W ithin each layer , the field is represented as a superposition of wa ve guide modes. The solution we will be found in the following form [1, 16, 17] H x ( x, y , z ) = Z ( z ) h x ( x, y ) , H y ( x, y , z ) = Z ( z ) h y ( x, y ) . (5) W e assume that the z dependence can be described by Z ( z ) = exp ( − ik z z ) . (6) The equations for the components h x and h y are ∆ ⊥ h x + k 2 0 εh x − k 2 z h x + 1 ε ∂ ε ∂ y  ∂ h y ∂ x − ∂ h x ∂ y  = 0 , ∆ ⊥ h y + k 2 0 εh y − k 2 z h y + 1 ε ∂ ε ∂ x  ∂ h x ∂ y − ∂ h y ∂ x  = 0 . (7) In this subsection the subscript of layer j is dropped for brevity . The functions h x and h y will be found in the form of a Fourier series: h x ( x, y ) = ∞ X m,n = −∞ B mn ψ mn , h y ( x, y ) = ∞ X m,n = −∞ C mn ψ mn , (8) ψ mn = exp ( − iκ x mx − iκ y ny ) . 4 A preprint From the system of Equations (7) we hav e the follo wing equations − X m,n ( κ 2 x m 2 + κ 2 y n 2 ) B mn ψ mn + k 2 0 X l,p ε lp ψ lp X m,n B mn ψ mn − X l,p ˆ ε ( y ) lp ψ lp X m,n iκ x mC mn ψ mn − X m,n iκ y nB mn ψ mn ! = k 2 z X m,n B mn ψ mn , − X m,n ( κ 2 x m 2 + κ 2 y n 2 ) C mn ψ mn + k 2 0 X l,p ε lp ψ lp X m,n C mn ψ mn + X l,p ˆ ε ( x ) lp ψ lp X m,n iκ x mC mn ψ mn − X m,n iκ y nB mn ψ mn ! = k 2 z X m,n C mn ψ mn . (9) Here ε mn = 1 L x L y L x / 2 Z − L x / 2 L y / 2 Z − L y / 2 ε ( x, y ) ψ mn dxdy , (10) ˆ ε ( x ) mn = 1 L x L y L x / 2 Z − L x / 2 L y / 2 Z − L y / 2 1 ε ( x, y ) ∂ ε ( x, y ) ∂ x 1 ψ mn dxdy , ˆ ε ( y ) mn = 1 L x L y L x / 2 Z − L x / 2 L y / 2 Z − L y / 2 1 ε ( x, y ) ∂ ε ( x, y ) ∂ y 1 ψ mn dxdy . Simplifying the equations, we obtain: − X m,n ( κ 2 x m 2 + κ 2 y n 2 ) B mn ψ mn + k 2 0 X l,p X m,n ε lp B mn ψ ( m + l )( n + p ) − i X l,p X m,n ˆ ε ( y ) lp ( κ x mC mn − κ y nB mn ) ψ ( m + l )( n + p ) = k 2 z X m,n B mn ψ mn , − X m,n ( κ 2 x m 2 + κ 2 y n 2 ) C mn ψ mn + k 2 0 X l,p X m,n ε lp C mn ψ ( m + l )( n + p ) + i X l,p X m,n ˆ ε ( x ) lp ( κ x mC mn − κ y nB mn ) ψ ( m + l )( n + p ) = k 2 z X m,n C mn ψ mn . (11) Then, we obtain the following equations: − ( κ 2 x m 2 + κ 2 y n 2 ) B mn + k 2 0 X l,p ε ( m − l )( n − p ) B lp − i X l,p ˆ ε ( y ) ( m − l )( n − p ) ( κ x lC lp − κ y pB lp ) = k 2 z B mn , − ( κ 2 x m 2 + κ 2 y n 2 ) C mn + k 2 0 X l,p ε ( m − l )( n − p ) C lp + i X l,p ˆ ε ( x ) ( m − l )( n − p ) ( κ x lC lp − κ y pB lp ) = k 2 z C mn . (12) In the system of equations (12) the values k 2 z must be read as eigen values, B mn and C mn are elements of eigen vectors of the matrix on the left side of these equations. The system of Equations (12) is transformed into following large algebraic eigen v alue problem for each layer j (see [17]) M ( j ) layer  B ( j ) C ( j )  =  k ( j ) z  2  B ( j ) C ( j )  , (13) 5 A preprint where M ( j ) layer is the matrix deri ved from the left-hand side of Equations (12), B ( j ) and C ( j ) are the eigen v ectors formed by B mn and C mn . The eigen v alues  k ( j ) z  2 giv e the propagation constants. The total field in each j th layer is a linear combination of modes " H ( j ) x H ( j ) y # = 2 N X p =1 k ( j ) z ; p h A ( j ) p ;1 e ik ( j ) z ; p z + A ( j ) p ;2 e − ik ( j ) z ; p z i × N x ,N y X m,n = − N x , − N y " B ( j ) p,mn C ( j ) p,mn # ψ mn , (14) where N x and N y are maximum numbers of the Fourier harmonics along the x and y axes, N = (2 N x + 1)(2 N y + 1) . The reflected field ( z > 0 ) is written as: " H ( r ) x H ( r ) y # = N x ,N y X m,n = − N x , − N y k z ; mn " A ( r ) x ; mn A ( r ) y ; mn # ψ mn e − ik z ; mn z . (15) The transmitted field ( z < − D ) is written as: " H ( t ) x H ( t ) y # = N x ,N y X m,n = − N x , − N y k z ; mn " A ( t ) x ; mn A ( t ) y ; mn # ψ mn e ik z ; mn z . (16) Here A ( j ) p ;1 , A ( j ) p ;2 , A ( r ) x ; mn , A ( r ) y ; mn , A ( t ) x ; mn and A ( t ) y ; mn are the unkno wn coef ficients. By satisfying continuity conditions of the tangential field components at each layer interface (at z = z min; j and z = z max; j for j th layer) H ( j − 1) x    z = z max; j = H ( j ) x    z = z max; j , H ( j − 1) y    z = z max; j = H ( j ) y    z = z max; j , E ( j − 1) x    z = z max; j = E ( j ) x    z = z max; j , E ( j − 1) y    z = z max; j = E ( j ) y    z = z max ; j , H ( j +1) x    z = z min; j = H ( j ) x    z = z min; j , H ( j +1) y    z = z min; j = H ( j ) y    z = z min; j , E ( j +1) x    z = z min; j = E ( j ) x    z = z min; j , E ( j +1) y    z = z min; j = E ( j ) y    z = z min; j , a global system of linear equations is formed: ˆ MA = R , (17) where ˆ M is matrix of system of equations obtained from the boundary conditions, A is the vector of unkno wn mode amplitudes (noted abo ve) and R is determined by the incident field. Solving this large linear system is the most computationally expensi ve part of the WG method. Figure 3 shows a visual summary of the standard W ave guide method. Here, the ’Solver’ block is the primary computational bottleneck. 3.3 Physics-Inf ormed Neural Networks According to PINN approach [5], which relies on the universal approximation theorem [18, 19], we represent the latent variables H x and H y with a deep neural network u θ ( x ) , where θ denote all trainable parameters (e.g., weights and biases) of the network. W e use multilayer perceptron (MLP) as deep neural network. Finding optimal parameters is an optimization problem, which requires the definition of a loss function such that its minimum gi ves the solution of the PDEs (2). The physics-informed model is trained by minimizing the composite loss function which consists of the local residuals of the differential equation o ver the problem domain and its boundary as sho wn below: L ( θ ) = λ bc L bc ( θ ) + λ r L r ( θ ) , (18) 6 A preprint Figure 3: Schematic diagram of the wa ve guide method. Here “F” is forward F ourier transform, “Eigen” is calculating eigen v alues and eigenv ectors of system (13) for every layer , “Field calculation” is calculating field with Equations (14) – (16). where L bc ( θ ) = 1 N bc N bc X i =1    B [ u θ ]  x ( bc ) i     2 , (19) L r ( θ ) = 1 N r N r X i =1    R [ u θ ]  x ( r ) i     2 , (20) R [ u 1 ] := ∆ u 1 + k 2 0 εu 1 + 1 ε ∂ ε ∂ y ·  ∂ u 2 ∂ x − ∂ u 1 ∂ y  , (21) R [ u 2 ] := ∆ u 2 + k 2 0 εu 2 + 1 ε ∂ ε ∂ x ·  ∂ u 1 ∂ y − ∂ u 2 ∂ x  . (22) Here L bc is a boundary loss term that corresponds to the boundary conditions, and L r is a residual loss term that corresponds to non-zero residuals of the governing PDEs (see details in [5, 6]), λ bc and λ r are hyperparameters, which allo w for separate tuning of the learning rate for each of the loss terms in order to improv e the con vergence of the model (see the procedure for choosing these parameters and examples [5 – 7]; in our calculations λ bc = N x /L x and λ r = 1 ), B is a boundary operator corresponding to boundary conditions: u ( x, y , z ) = u ( x + L x , y , z ) , u ( x, y , z ) = u ( x, y + L y , z ) , L x / 2 Z − L x / 2 L y / 2 Z − L y / 2  ∂ u ∂ z + ik z ; mn u  1 ψ mn dxdy        z = z max = 0 , L x / 2 Z − L x / 2 L y / 2 Z − L y / 2  ∂ u ∂ z − ik z ; mn u  1 ψ mn dxdy        z = z min = 0 , (23) where the first tw o conditions are the periodicity condition of the solution, the last two are the Sommerfeld conditions at the upper ( z max ) and lower ( z min ) boundaries of the domain under consideration. n x ( bc ) i o N bc i =1 and n x ( r ) i o N r i =1 are sets of points corresponding to boundary condition domain and PDE domain, respectiv ely . These points can be the vertices of a fixed mesh or can be randomly sampled at each iteration of a gradient descent algorithm. All required gradients w .r .t. input variables ( x ) and network parameters ( θ ) can be efficiently computed via automatic dif ferentiation [20] with algorithmic accuracy , which is defined by the accurac y of computation system. 7 A preprint The optimization problem can be defined as follows: θ ∗ = arg min θ L ( θ ) , (24) where θ ∗ are optimal parameters of the neural netw ork which minimize the discrepanc y between the e xact unkno wn solution u and the approximate one u θ θ θ ∗ . A schematic diagram of the general PINN method, which is used throughout this paper , is shown Figure 4. This neural network learning method is unsupervised, meaning that training is carried out without labeled datasets. Figure 4: Schematic diagram of the general PINN method. 3.4 W av eguide Neural Operator The development and use of neural networks primarily focus on training mappings between finite-dimensional Euclidean spaces or finite sets. This approach has been extended to neural operators that learn mappings between functional spaces [8, 9, 21, 22]. For partial differential equations, neural operators directly learn a mapping from any giv en functional parametric dependence to a solution, thereby learning an entire class of PDE solutions, as opposed to classical methods which solve one specific instance of an equation. The process of solving the giv en problem can be considered as operator G : A → U between function spaces A and U . In the case of the mask problem, A and U are an infinite function space of functions of incident field and an infinite function space of functions of reflected, transmitted and in-mask fields, respectiv ely . The WG method is solution operator G † of the giv en dif fraction problem: G † : A → U . (25) W e aim to approximate G † by constructing a parametric map G θ θ θ : A → U , θ θ θ ∈ Θ , (26) for some finite-dimensional parameter space Θ by choosing θ θ θ † ∈ Θ so that G θ θ θ † ≈ G † . Let’ s construct G θ θ θ from G † by replacing the most resource-intensiv e part of G † with a deep neural network. Specifically , we train a multi-layer perceptron (MLP) to learn the mapping from the right-hand side v ector R and the mask parameters (Fourier coef ficients of permitti vity) to the solution vector A of the linear system (17): A θ = MLP  R , { ε ( j ) mn } J j =1 ; θ  , (27) where ε ( j ) mn is giv en by Equation (10). The optimization problem can be defined as follow: θ † = arg min θ L ( θ ) , (28) L ( θ ) = ∥ ˆ MA θ − R ∥ 2 2 . In this formulation of the problem, the neural operator G θ θ θ is mesh–independent, and its training takes place in a latent (Fourier) space, inherits the physical structure of the WG method, and directly targets the primary computational bottleneck. A schematic diagram of the WGNO is shown Figure 5. 8 A preprint Figure 5: Schematic diagram of the wave guide neural operator . Here “F” is forward Fourier transform, “Eigen” is calculating eigen v alues and eigen vectors of system (13) for e very layer , “Field calculation” is calculating field with Equations (14) – (16). 4 Numerical Experiments W e e v aluated the performance of the proposed methods on se veral test problems. The TE polarization of the incident wa ves w as used for all considered cases. The neural networks were implemented in PyT orch [23] (version 2.6 under CUD A 12.4) and the training was carried out on a node with GPU Nvidia GeForce GTX 1660 T i and CPU Intel Core i5-9300H. For the PINNs, we used a MLP with 2 hidden layers and 128 neurons each with hyperbolic tangent activ ation function ( tanh ). T o sa ve computational resources and to accelerate con vergence for problems with periodic boundary conditions along x (in 2D case) it was constructed a Fourier feature embedding of input coordinates in the follo wing form v ( x ) = [1 , cos( k x x ) , sin( k x x ) , cos(2 k x x ) , sin(2 k x x ) , ..., cos( mk x x ) , sin( mk x x )] , (29) where k x = 2 π /L x and m is some non-negati ve integer . W ith such embedding any trained approximation of problem solution by PINN exactly satisfies periodic conditions (see details in [6, 24, 25]). For WGNO, the MLP had 2 hidden layers. 4.1 V alidation on T est Pr oblems T o ev aluate the accuracy of the solutions provided by numerical solvers (FreeFem++, WG) and solvers based on artificial neural networks, we hav e tested the methods on three 2D problems with kno wn analytical solutions: (1) plane wa ve propagation in a homogeneous medium, (2) reflection from a single interface (Figure 6 (a)), and (3) reflection from a dielectric layer (Figure 6 (b)). As the main performance metrics for the considered solvers, we e v aluated the relati ve L 2 error (see [6] and appendix A) and the time to provide a solution (the time to calculate a solution for numerical solver or the inference time for neural networks). In addition, training times are marked separately for neural networks. W e assume the media are homogeneous in the y direction. An electromagnetic wa ve of TE polarization propag ates in free space with wave vector k = k ( i ) x x 0 − k ( i ) z z 0 ( k ( i ) x and k ( i ) z are positi ve values, and k ( i ) x = k ( i ) z ). The medium belo w interface and in the dielectric layer have permitti vity ε = 4 . The dielectric layer thickness is π /k ( i ) z . The coordinates ( x, z ) of the points in the computational domain satisfy the relations: x ∈ [ − π /k ( i ) x , π /k ( i ) x ] and z ∈ [ − π /k ( i ) z , π /k ( i ) z ] ( z ∈ [ − 2 π /k ( i ) z , π /k ( i ) z ] for the dielectric layer). For the numerical solv er we divided e very axis into 100 interv als, as a result we had about 10 4 collocation points. N x equals 10 for WG and WGNO methods at all test experiments. For the PINN and the WGNO approach we divided the segment along e very axis into 100 intervals. In our experiments for the PINN, we used tw o-stage learning which consisted of 1000 epochs of optimization with the Adam optimizer with learning rate 10 − 3 and 5 epochs of optimization with the LBFGS optimizer . For WGNO we used 1000 epochs of optimization with the Adam optimizer with learning rate 10 − 3 and 10 − 5 . W e trained randomly initialized neural networks using the Glorot scheme [26], repeating the training 7 times with dif ferent random seeds. The results are summarized in T able 1. The WG method achieves machine precision, as expected. FreeFem++ provides sufficient computational accuracy for applied (engineering) needs. The PINN achie ves reasonable accuracy , b ut it is less accurate than FEM in the case of the third problem (dielectric layer). The proposed WGNO demonstrates high performance, with errors sev eral orders of magnitude smaller than the PINN, highlighting the effecti veness of the hybrid 9 A preprint Figure 6: Problems with exact solutions: (a) is problem of electromagnetic wa ves reflection from a single interface, (b) is problem of electromagnetic wa ves reflection from a dielectric layer . physics-based approach. Comparisons of the electric field distributions of the exact solution and the solution gi ven by solvers are gi ven in Appendix B. The fields are normalized to the electric field amplitude of the incident w av e. T able 1: Performance comparison on 2D test problems. Inference times are 1 . 2 × 10 − 3 s and 1 . 7 × 10 − 4 s for the PINN and the WGNO, respectiv ely . Problem W a ve guide FreeFem++ Rel. L 2 Error T ime (s) Rel. L 2 Error T ime (s) 1. W a ve propagation 2 . 814 × 10 − 15 1 . 0 × 10 − 3 2 . 323 × 10 − 3 3.467 2. W a ve scattering on interf ace 2 . 730 × 10 − 15 2 . 0 × 10 − 3 3 . 818 × 10 − 3 3.385 3. Reflection from layer 6 . 110 × 10 − 15 6 . 0 × 10 − 3 8 . 085 × 10 − 3 5.033 Problem PINN WGNO Rel. L 2 Error T raining (s) Rel. L 2 Error T raining (s) 1. W a ve propagation 1 . 694 × 10 − 4 ± 5 . 3 × 10 − 5 855 4 . 032 × 10 − 8 ± 1 . 6 × 10 − 8 3 . 7 2. W a ve scattering on interf ace 1 . 596 × 10 − 3 ± 7 . 9 × 10 − 4 1532 1 . 785 × 10 − 7 ± 2 . 8 × 10 − 7 3 . 7 3. Reflection from layer 5 . 507 × 10 − 2 ± 4 . 4 × 10 − 3 1221 4 . 709 × 10 − 5 ± 5 . 3 × 10 − 5 11 4.2 2D Lithography Mask Simulation W e simulated a realistic 2D EUV mask (see Fig. 2) for wa velengths ( λ ) of 13.5 and 11.2 nm. These wav elengths correspond to the operating wa velengths of current industrial lithography systems [10] and of promising ones [11, 12]. 4.2.1 W av elength 13.5 nm The mask included an absorber consisting of a T aBO layer and a T aBN layer (10 nm for PINN, and 60 nm for WGNO), a thin Ru layer , a Mo/Si layer for PINN and a 31-layer Mo/Si mirror for the WGNO (see Figure 2 (b)). The parameters 10 A preprint T able 2: Parameters of the mask for wa velength λ =13 . 5 nm which used to calculations. Method Properties Media Layer T aBO T aBN Ru Mo Si ε 0 . 857 − i 0 . 079 0 . 861 − i 0 . 071 0 . 785 − i 0 . 030 0 . 853 − i 0 . 012 0 . 998 − i 0 . 004 PINN thickness 10 nm 10 nm 2 nm 3 nm 4 nm WGNO thickness 10 nm 60 nm 2 nm 3 nm 4 nm of these masks are given in T able 2. The absorbers had a hole with a width L x / 2 so the permittivities of the absorber in layer j ( j is 1 or 2) are described by: ε j ( x ) = 1 2  tanh  x + a d  − tanh  x − a d  (1 − ε ) + ε. (30) Here ε is dielectric permitti vity of layer , 2 a is length of hole ( a = L x / 4 ) and d is the thickness of the transition region at the edge of the hole ( d = λ/ 10 ). Note that we obtained the permittivities of the media from experimental data [13, 14].The angle of incidence of the wa ve is 6 ◦ , which corresponds to the angle of incidence of the wa ves on the mask in modern lithographic systems. As can be seen from Figure 7, which shows the dependence of the reflection coef ficient (po wer of reflected wa v es normalized to the po wer of incidence wav es) on the angle of incidence of TE EUV wa ve for Mo/Si mirror of 60 layers, this angle lies on the plateau of the maximum of this dependence. Figure 7: Dependence of reflectivity of TE EUV w av e (wa velength 13.5 nm) on incident angle for Mo/Si mirror of 60 layers. Thicknesses of Mo and Si layers are giv en in T able 2. Maximum of reflectivity is 0.75 at incident angel 8.86 ◦ . Figures 8 and 9 show the comparisons of the calculated E y field, normalized to the amplitude of the incident wa ve, predicted by PINN and WGNO, respectively . It follo ws from the demonstrated results that while the PINN captures the general wa ve pattern, its absolute error is significant. In contrast, the WGNO result is visually indistinguishable from the reference WG solution, with a very small absolute error (see Fig. 9). 4.2.2 W av elength 11.2 nm For the w av elength of 11.2 nm, the top two layers of the mask follo w the same configuration as described above: an absorber consisting of a T aBO layer (10 nm) and a T aBN layer (10 nm for PINN, 60 nm for WGNO). Since a standardized solution for mirrors at this wavelength has not yet been established and related research is still ongoing [27 – 29], we considered sev eral potential mirror designs for further experiments. Figure 10 shows the reflecti vity as a function of the angle of incidence for Mo/Be, Ru/Be, and Ru/Be/Sr mirrors, each consisting of 60 layers with thicknesses specified in T able 3. As can be seen, the Ru/Be/Sr mirror exhibits the highest reflectivity . Consequently , this structure was 11 A preprint Figure 8: (a) is real part of k in the media normalized to k 0 . Comparison of reference WG solution (b) and PINN solution (c) for a 2D mask at 13.5 nm, (d) is absolute error . Figure 9: (a) is real part of k in the media normalized to k 0 . Comparison of reference WG solution (b) and WGNO solution (c) for a 2D mask at 13.5 nm, (d) is absolute error . chosen as the mask substrate (see Figure 2 (c)). For the PINN-based solver , a single layer of this structure was utilized, whereas for the WGNO, the mirror consisted of a 31-layer Ru/Be/Sr . The dielectric permittivities of T aBO and T aBN are ε = 0 . 906 − i 0 . 064 and ε = 0 . 909 − i 0 . 060 respecti vely . As in the previous mask case, the absorbers had a hole with width L x / 2 so the permittivities of the absorbers are described by Eqn. (29). The angle of incidence of the wa ve is 6 ◦ . Figures 11 and 12 show the comparisons of the calculated E y field, normalized to the amplitude of the incident wa ve, predicted by PINN and WGNO, respectiv ely , for the wa velength 11.2 nm. The quantitati ve results for the both wavelengths are summarized in T able 4. The PINN struggles with this comple x problem, yielding errors of several percent after long training times. The WGNO, howe ver , achieves e xcellent accuracy with a short training time. T able 3: Parameters of layers of mirror for wa velength λ =11 . 2 nm. Mirror Properties Media Layer Mo Ru Be Sr ε 0 . 91 − i 0 . 009 0 . 872 − i 0 . 012 1 . 025 − i 0 . 003 0 . 986 − i 0 . 002 Mo/Be thickness 2.22 nm – 3.5 nm – Ru/Be thickness – 2.01 nm 3.72 nm – Ru/Be/Sr thickness – 1.7 nm 2.7 nm 1.34 nm 12 A preprint Figure 10: Dependences of reflectivity of TE EUV w av e (wa velength 11.2 nm) on incident angle for Mo/Be, Ru/Be and Ru/Be/Sr mirrors of 60 layers. Thicknesses of layers are giv en in T able 3. Maximum of reflectivities are: 0.77 at incident angel 10.71 ◦ for Mo/Be, 0.79 at incident angel 11.06 ◦ for Ru/Be, and 0.81 at incident angel 10.96 ◦ for Ru/Be/Sr . Figure 11: (a) is real part of k in the media normalized to k 0 . Comparison of reference WG solution (b) and PINN solution (c) for a 2D mask at 11.2 nm, (d) is absolute error . T able 4: Performance on 2D lithograph y mask simulation. Inference times are 1 . 2 × 10 − 3 s and 1 . 8 × 10 − 4 s for PINN and WGNO, respectiv ely . λ (nm) Method Rel. L 2 Error T raining (s) 13.5 PINN 3 . 1 × 10 − 2 1003 WGNO 3 . 8 × 10 − 6 4 . 3 11.2 PINN 3 . 3 × 10 − 2 1252 WGNO 3 . 5 × 10 − 7 7 . 7 4.3 3D Lithography Mask Simulation Finally , we extended our analysis to the full, challenging 3D mask problem. The mask structure and parameters of the incident field are as in pre vious case for the WGNO. The nonuniform permitti vities of the absorber in layer j are 13 A preprint Figure 12: (a) is real part of k in the media normalized to k 0 . Comparison of reference WG solution (b) and WGNO solution (c) for a 2D mask at 11.2 nm, (d) is absolute error . described by: ε j ( x, y ) = 1 4  tanh  x + a d  − tanh  x − a d  ×  tanh  y + b d  − tanh  y − b d  (1 − ε )+ ε, (31) where ε is permittivity of absorber of j th layer , a = L x / 4 , b = L y / 4 , d = λ/ 10 . Figures 13 and 14 show cross-sections of the field calculated by our WGNO compared to the reference solution at 13.5 nm and 11.2 nm, respectiv ely . The agreement is excellent, demonstrating the capability of our method to handle the full 3D problem with high fidelity . The performance metrics for the 3D case, sho wn in T able 5, are ev en more compelling. The WGNO maintains high accuracy while achie ving a speedup of over 200 times compared to the rigorous WG solver . The entire training process for the 3D case took less than 20 seconds, demonstrating the remarkable effi ciency and scalability of our proposed operator . T able 5: Performance of WGNO on 3D lithography mask simulation. Training time less 20 s. W a velength Rel. L 2 Error Inference (s) T raining time (s) 13.5 nm 3 . 1 × 10 − 7 2 . 24 × 10 − 4 19.4 11.2 nm 2 . 8 × 10 − 6 2 . 19 × 10 − 4 18.9 4.3.1 Generalization For the e v aluation of generalization properties of WGNO we use the measure of generalization µ (see in [7]) which is giv en in appendix C. W e in vestigate these properties of the described neural network using the 2D example 4.2.1 (the wa velength is 13.5 nm). W e used the following parameter of the considered problem a ∈ [0 . 99 L x / 4 , L x / 4] , which is half the length of the hole in the absorbers. Consider ho w the generalization parameter of the solutions provided by WGNO depends on the number of partitions N q of the area of change of the parameter a . The calculation results are presented in the T able 6. As can be seen from this table, the measure of generalizing properties µ increases with increasing size of N q , and this measure of generalization quickly reaches values close to 1. 5 Conclusion In this work, we have sho wn that physics-informed neural networks and neural operators can achieve high accuracy for complex dif fraction problems of EUV lithography simulation. W e determined that while training times can v ary , the inference time is extremely small. For relati vely thin masks, PINNs can pro vide solutions accurate enough for initial 14 A preprint Figure 13: Comparison of reference WG solution ((a),(d)) and WGNO solution ((b),(e)) for a 3D mask at 13.5 nm, (c) and (f) are absolute errors. (a)–(c) are in cross-section y = 0 , (d)–(f) are in the cross-section x = 0 . Figure 14: Comparison of reference WG solution ((a),(d)) and WGNO solution ((b),(e)) for a 3D mask at 11.2 nm, (c) and (f) are absolute errors. (a)–(c) are in cross-section y = 0 , (d)–(f) are in the cross-section x = 0 . calculations. Most significantly , we have established that our proposed W ave guide neural operator is a state-of-the-art solver for this application. Using the example of solving a 2D problem, it is shown that the presented WGNO has 15 A preprint T able 6: Dependence of measure of generalization µ of WGNO on the number of partitions N q of the area of change of the parameter a . N q Measure of generalization µ 5 0.94225 10 0.98473 20 0.99607 50 0.99922 100 0.99978 pronounced generalizing properties: for problem parameters unseen during training, it provides a solution accuracy close to that of the parameters present in the training dataset. It provides highly accurate solutions for a full 3D mask, surpassing modern approaches in both accuracy and inference time, all while requiring only a very short training period. This lev el of performance makes the WGNO a highly promising tool for accelerating the OPC design c ycle in semiconductor manufacturing. Furthermore, we are confident that the approaches we have de veloped are not limited to lithography and can be applied to study the properties and optimization of other comple x electromagnetic structures, such as metamaterials and composite materials. 16 A preprint A Relative L 2 error T o e v aluate the accuracy of the approximate solution obtained with the help of considered methods, the values of the solution of problems (2)–(4) calculated by the numerical solver or predicted by the neural network at gi ven points are compared with the exact solution or the values calculated on the basis of classical high-precision numerical methods. As a measure of accuracy , the relative total L 2 error of prediction is taken, which can be e xpressed with the follo wing relation ϵ error = ( 1 N e N e X i =1 [ u θ θ θ ( x i , y i ) − u ( x i , y i )] 2 ) 1 / 2 × ( 1 N e N e X i =1 [ u ( x i , y i )] 2 ) − 1 / 2 , (32) where { x i , y i } N e i =1 is the set of ev aluation points taken from the domain, u θ θ θ and u are the predicted and reference solutions, respectiv ely , N e is the number of ev aluation points. B Comparisons of the electric field distributions of the exact solution and the solutions given by the solvers Figures 15 – 18 present a comparison between the exact solutions of the test problems and those obtained using FreeMem++, the wa ve guide method, PINN, and WGNO. The selected test cases include: plane wa ve propagation in a homogeneous medium, reflection from a single interface (Figure 6 (a)), and reflection from a dielectric layer (Figure 6 (b)). C Measure of generalization pr operties T o e valuate the generalization properties of the approximate solution obtained using the physics-informed methods, two sets of solutions predicted by the neural network on homogeneous grids of this parameter are required for the problem parameter q . The first set of solution values is calculated for the grid of the parameter q in steps of ∆ q for those values of q at which the training was carried out. W e denote this set of parameters as { q i } N q i =1 . The second set of solution values is calculated for the calculated grid of the parameter q in steps of ∆ q , of fset by ∆ q / 2 , i.e. at points as far aw ay as possible from the nearest points of q where the training was carried out, b ut not out of the training domain. W e will denote this set of parameters as  q i +1 / 2  N q − 1 i =1 . As a measure of generalization properties, the ratio of the relative total L 2 error of prediction of set  q i +1 / 2  N q − 1 i =1 to the ratio of the relativ e total L 2 error of prediction of set { q i } N q i =1 is taken, which can be expressed with the follo wing relation µ = ϵ 1 [ u θ , u ] ϵ 1 / 2 [ u θ , u ] . (33) Here ϵ 1 / 2 [ u θ , u ] =    N q − 1 X j =1 N e X i =1  u θ ( x i , q j +1 / 2 ) − u ( x i , q j +1 / 2 )  2    1 / 2 ×    N q − 1 X j =1 N e X i =1  u ( x i , q j +1 / 2 )  2    − 1 / 2 , ϵ 1 [ u θ , u ] =    N q X j =1 N e X i =1 [ u θ ( x i , q j ) − u ( x i , q j )] 2    1 / 2 ×    N q X j =1 N e X i =1 [ u ( x i , q j )] 2    − 1 / 2 , where { x i } N e i =1 is the set of e v aluation points taken in the computational domain, { q i } N q i =1 and  q i +1 / 2  N q − 1 i =1 are sets taken from domain [ q min , q max ] , u θ and u are the predicted and reference solutions, respectiv ely , the second argument of the functions u θ and u denotes the parameter of the dif ferential equation for which these functions were calculated. Thus, the closer the value of the measure of generalization properties µ is to one, the neural network outputs solutions to the problem at the parameters of the problem at which it was not trained, closer in accurac y to its solutions to the problem at the parameters at which it was trained. 17 A preprint Figure 15: FreeFem++ solver . The left panels show the e xact solutions, the middle panels sho w the solutions calculated using the gi ven method, and the right panels show the absolute dif ference between the exact and calculated solutions. The relativ e error ϵ error are giv en in the table 1. Panels (a), (b) and (c) correspond the problems of wa ve propagation in a homogeneous medium, reflection of an electromagnetic wav es from a interface between two media and reflection from a dielectric layer , respecti vely . References [1] H. T anabe, A. Jinguji, and A. T akahashi, “W eakly guiding approximation of a three-dimensional wa veguide model for extreme ultraviolet lithography simulation, ” J. Opt. Soc. Am. A , vol. 41, no. 8, pp. 1491–1499, Aug 2024. [Online]. A v ailable: https://opg.optica.org/josaa/abstract.cfm?URI=josaa- 41- 8- 1491 [2] ——, “ Accelerating extreme ultraviolet lithography simulation with weakly guiding approximation and source position dependent transmission cross coefficient formula, ” J ournal of Micro/Nanopatterning , Materials, and Metr ology , vol. 23, no. 1, pp. 014 201–014 201, 2024. [3] V . Medv ede v , A. Erdmann, and A. Roßk opf, “3d mask simulation and lithographic imaging using physics-informed neural networks, ” 2023. [Online]. A vailable: https://publica.fraunhofer .de/handle/publica/450903 [4] V . M. Medv ede v , A. Erdmann, and A. Roßk opf, “3d euv mask simulator based on physics- informed neural networks: effects of polarization and illumination, ” 2024. [Online]. A vailable: https://publica.fraunhofer .de/handle/publica/477957 18 A preprint Figure 16: The same is shown in Figure 15 for the w av eguide method. [5] M. Raissi, P . Perdikaris, and G. E. Karniadakis, “Physics-informed neural networks: A deep learning framew ork for solving forward and inv erse problems in volving nonlinear partial dif ferential equations, ” Journal of Computational Physics , vol. 378, pp. 686–707, 2019. [Online]. A vailable: https://linkinghub .else vier .com/retriev e/pii/S0021999118307125 [6] V . A. Es’kin, D. V . Davydo v , E. D. Egorova, A. O. Malkhanov , M. A. Akhukov , and M. E. Smorkalov , “About Modifications of the Loss Function for the Causal Training of Physics-Informed Neural Networks, ” Doklady Mathematics , vol. 110, no. S1, pp. S172–S192, dec 2024. [Online]. A vailable: https://link.springer .com/10.1134/S106456242460194X [7] V . A. Es’kin, A. O. Malkhanov , and M. E. Smorkalov , “ Are two hidden layers still enough for the physics-informed neural networks?” Journal of Computational Physics , vol. 537, p. 114085, 2025. [Online]. A v ailable: https://www .sciencedirect.com/science/article/pii/S0021999125003687 [8] L. Lu, P . Jin, G. Pang, Z. Zhang, and G. E. Karniadakis, “Learning nonlinear operators via deeponet based on the univ ersal approximation theorem of operators, ” Natur e machine intelligence , v ol. 3, no. 3, pp. 218–229, 2021. [9] Z. Li, N. K ov achki, K. Azizzadenesheli, B. Liu, K. Bhattacharya, A. Stuart, and A. Anandkumar , “Fourier neural operator for parametric partial differential equations, ” 2021. [Online]. A v ailable: https://arxiv .org/abs/2010.08895 [10] A. Y akshin, E. Louis, P . Görts, E. Maas, and F . Bijkerk, “Determination of the layered structure in mo/si multilayers by grazing incidence x-ray reflectometry , ” Physica B: Condensed Matter , vol. 283, no. 1-3, pp. 143–148, 2000. 19 A preprint Figure 17: Approach based on the PINN. The left panels sho w the exact solutions, the middle panels show the calculated solutions using the giv en method, right panels are absolute difference between exact solutions and calculates solutions. The relativ e error ϵ error are giv en in the table 1. Panels (a), (b) and (c) correspond the problems of wa ve propagation in a homogeneous medium, reflection of an electromagnetic wav es from a interface between two media and reflection from a dielectric layer , respecti vely . [11] N. I. Chkhalo, “New Concept for the De velopment of High-Performance X-ray Lithography, ” Russian Micr oelectronics , v ol. 53, no. 5, pp. 397–407, oct 2024. [Online]. A v ailable: https://link.springer .com/10.1134/ S1063739724600511 [12] I. Abramov , S. Golubev , E. Gospodchikov , A. Shalashov , A. Perekalov , A. Nechay , and N. Chkhalo, “Laser discharge in a high-pressure jet of heavy noble gas: Expansion of emitting v olume promises an efficient source of euv light for lithography , ” Phys. Rev . Appl. , vol. 23, p. 024004, Feb 2025. [Online]. A vailable: https://link.aps.org/doi/10.1103/PhysRe vApplied.23.024004 [13] B. Henke, E. Gullikson, and J. Davis, “X-ray interactions: Photoabsorption, scattering, transmission, and reflection at e = 50-30,000 e v , z = 1-92, ” Atomic Data and Nuclear Data T ables , vol. 54, no. 2, pp. 181–342, 1993. [Online]. A v ailable: https://www .sciencedirect.com/science/article/pii/S0092640X83710132 [14] Center for X-Ray Optics. Lawrence Berkeley National Laboratory. (1993–2025) Index of refraction. [Online]. A v ailable: https://henke.lbl.go v/optical_constants/getdb2.html 20 A preprint Figure 18: The same is shown in Figure 17 for the w av eguide neural operator . [15] F . Hecht, “Ne w dev elopment in freefem++, ” Journal of numerical mathematics , v ol. 20, no. 3-4, pp. 251–266, 2012. [16] C.-M. Y uan and A. J. Strojwas, “Modeling optical microscope images of integrated-circuit structures, ” Journal of the Optical Society of America A , vol. 8, no. 5, pp. 778–790, 1991. [17] K. D. Lucas, H. T anabe, and A. J. Strojw as, “Efficient and rigorous three-dimensional model f or optical lithography simulation, ” J ournal of the Optical Society of America A , vol. 13, no. 11, pp. 2187–2199, 1996. [18] K. Hornik, M. Stinchcombe, and H. White, “Multilayer feedforward networks are uni versal approximators, ” Neural Networks , v ol. 2, no. 5, pp. 359–366, 1989. [Online]. A vailable: https://www .sciencedirect.com/science/ article/pii/0893608089900208 [19] G. V . Cybenko, “ Approximation by superpositions of a sigmoidal function, ” Mathematics of Contr ol, Signals and Systems , vol. 2, pp. 303–314, 1989. [Online]. A vailable: https://api.semanticscholar .or g/CorpusID:3958369 [20] A. Griewank and A. W alther , Evaluating Derivatives , 2nd ed. Society for Industrial and Applied Mathematics, 2008. [Online]. A v ailable: https://epubs.siam.org/doi/abs/10.1137/1.9780898717761 [21] Z. Li, N. K ov achki, K. Azizzadenesheli, B. Liu, K. Bhattacharya, A. Stuart, and A. Anandkumar , “Neural operator: Graph kernel network for partial differential equations, ” 2020. [Online]. A v ailable: https://arxiv .org/abs/2003.03485 21 A preprint [22] H. Y ang, Z. Li, K. Sastry , S. Mukhopadhyay , M. Kilgard, A. Anandkumar , B. Khailany , V . Singh, and H. Ren, “Generic lithography modeling with dual-band optics-inspired neural networks, ” 2022. [Online]. A vailable: https://arxiv .org/abs/2203.08616 [23] A. Paszke, S. Gross, F . Massa, A. Lerer, J. Bradb ury , G. Chanan, T . Killeen, Z. Lin, N. Gimelshein, L. Antiga, A. Desmaison, A. Köpf, E. Y ang, Z. DeV ito, M. Raison, A. T ejani, S. Chilamkurthy , B. Steiner , L. Fang, J. Bai, and S. Chintala, PyT or ch: An Imperative Style, High-P erformance Deep Learning Library . Red Hook, NY , USA: Curran Associates Inc., 2019. [24] S. W ang, S. Sankaran, and P . Perdikaris, “Respecting causality for training ph ysics-informed neural networks, ” Computer Methods in Applied Mechanics and Engineering , vol. 421, p. 116813, 2024. [Online]. A vailable: https://www .sciencedirect.com/science/article/pii/S0045782524000690 [25] S. Dong and N. Ni, “A method for representing periodic functions and enforcing exactly periodic boundary conditions with deep neural netw orks, ” J ournal of Computational Physics , v ol. 435, p. 110242, 2021. [Online]. A v ailable: https://www .sciencedirect.com/science/article/pii/S0021999121001376 [26] X. Glorot and Y . Bengio, “Understanding the difficulty of training deep feedforward neural networks, ” in Pr oceedings of the Thirteenth International Conference on Artificial Intellig ence and Statistics , ser . Proceedings of Machine Learning Research, Y . W . T eh and M. Titterington, Eds., v ol. 9. Chia Laguna Resort, Sardinia, Italy: PMLR, 2010, pp. 249–256. [Online]. A v ailable: https://proceedings.mlr .press/v9/glorot10a.html [27] N. I. Chkhalo, S. A. Gusev , A. N. Nechay , D. E. Pariev , V . N. Polkovnikov , N. N. Salashchenko, F . Schäfers, M. G. Sertsu, A. Sokolov , M. V . Svechniko v , and D. A. T atarsky , “High-reflection mo/be/si multilayers for euv lithography , ” Opt. Lett. , v ol. 42, no. 24, pp. 5070–5073, Dec 2017. [Online]. A v ailable: https://opg.optica.org/ol/abstract.cfm?URI=ol- 42- 24- 5070 [28] V . N. Polko vnikov , N. N. Salashchenko, M. V . Svechniko v , and N. I. Chkhalo, “Beryllium-based multilayer x-ray optics, ” Physics-Uspekhi , vol. 63, no. 1, p. 83, jan 2020. [Online]. A vailable: https://doi.org/10.3367/UFNe.2019.05.038623 [29] R. A. Shaposhniko v , V . N. Polkovnik ov , N. N. Salashchenk o, N. I. Chkhalo, and S. Y . Zue v , “Highly reflecti ve ru/sr multilayer mirrors for wav elengths 9–12 nm, ” Opt. Lett. , vol. 47, no. 17, pp. 4351–4354, Sep 2022. [Online]. A v ailable: https://opg.optica.org/ol/abstract.cfm?URI=ol- 47- 17- 4351 22