 1. Introduction
 2. Basic equations and starting points
 3. The propagation of Xrays in a vacuum
 4. Oriented Gaussian beams for describing the Xray propagation through optical elements
 5. An estimate of the acceptable step h and the method error
 6. Technique for modeling the propagation and focusing of Xray waves in terms of oriented Gaussian beams
 7. Conclusions
 References
 1. Introduction
 2. Basic equations and starting points
 3. The propagation of Xrays in a vacuum
 4. Oriented Gaussian beams for describing the Xray propagation through optical elements
 5. An estimate of the acceptable step h and the method error
 6. Technique for modeling the propagation and focusing of Xray waves in terms of oriented Gaussian beams
 7. Conclusions
 References
research papers
Oriented Gaussian beams for highaccuracy computation with accuracy control of Xray propagation through a multilens system
^{a}Gdansk University of Technology, Faculty of Applied Physics and Mathematics, Gdansk, Poland, ^{b}Center for Functionalized Magnetic Materials (FunMagMa), Immanuel Kant Baltic Federal University, Kaliningrad, Russia, and ^{c}Immanuel Kant Baltic Federal University, Kaliningrad, Russia
^{*}Correspondence email: pawwojda@pg.edu.pl
A highly accurate method for calculating Xray propagation is developed. Within this approach, the propagating wave is represented as a superposition of oriented Gaussian beams. The direction of wave propagation in each Gaussian beam agrees with the local direction of propagation of the Xray wavefront. When calculating the propagation of Xray waves through lenses, the thin lens approximation is applied. In this approximation, the wave parameters change discontinuously when the wave passes through a lens; the corresponding explicit formulae are derived. The theory is applied to highly accurate calculation of the focusing of Xrays by a system of many beryllium lenses. Fine structure of the wave electric field on the focal plane is revealed and studied. The fine structure is formed due to the diffraction of waves at the edges of the lens apertures. Tools for controlling the calculation accuracy are proposed. The amplitude of the electric field on the focal plane and the focal spot width are shown to be very sensitive to the quality of the calculation, while the best focus position can be obtained even from simple calculations.
Keywords: Xray wave; Xray optics; Xray propagation; Helmholtz equation; Gaussian beams; wave propagation direction; lens; focusing; control of calculation accuracy; diffraction.
1. Introduction
At present, Xray optics are developing rapidly due to the possibilities of visualization of very small objects that they provide. To visualize the structure of a very small object, it is first necessary to focus the Xrays. Currently, commonly used focusing elements are bent crystals, multilayers or mirrors, capillaries, waveguides, refractive lenses and diffraction elements, such as Fresnel zone plates or Laue multilayer lenses (Snigirev et al., 1996; Roth et al., 2017).
In this paper, focusing in Xray optics using Xray refraction mechanisms will be considered. The first successes in this field were presented by Snigirev et al. in 1996, where a composite refractive lens was used (Snigirev et al., 1996). Since then, many experiments have been performed using focused Xrays. The use of lenses with a spherical or parabolic concave profile is one of the most popular methods of Xray focusing nowadays. Recently, researchers have also used diamond lenses for Xray focusing (Roth et al., 2017).
The materials used for Xray focusing have some common features. First of all, they are characterized by a complex n = 1 − δ + iβ, where δ is responsible for the refraction of Xrays in lenses and β for the attenuation of the The ideal material for making Xray lenses should have a δ value as large as possible and an β as small as possible. In addition, the material should allow the technological possibility of manufacturing lenses with a small curvature radius R. Knowing the material properties, it is possible to calculate with high accuracy the focal length f of an optical system consisting of N lenses with a given surface curvature R using the formula f = . In addition to the listed properties, the materials must be strong and inexpensive. Unfortunately, it is not easy to satisfy the latter two conditions. Beryllium is one of the most popular materials currently used for focusing Xrays; the cost of manufacturing a beryllium lens is about 5000 euros. Unfortunately, while the lens is in use the lens material is subject to oxidation, and after a while the lens has to be replaced. In order to effectively focus Xrays at relatively small distances and for achieving a large magnification, use of complex optical systems consisting of several dozen lenses is necessary. Therefore, the experiments are associated with significant costs. A way to reduce this cost would be to find or develop cheaper and more effective materials for lens manufacturing. Performing more extensive theoretical studies based on mathematical modeling is another way. This work is devoted to the development of a more advanced mathematical apparatus for theoretical studies and the application of advanced methods for calculating the Xray focusing by a system of many lenses.
At present, in theoretical studies in Xray optics, the most popular way to solve the paraxial wave equation is by fast Fourier transform (FFT) (Goodman, 1996). FFT is a very effective calculation method. However, the exact solution of the paraxial equation has the form of an infinite Fourier series, and this series has to be truncated when we move to the FFT, which leads to error. The magnitude of the error obviously depends on the circumstances, but the issue of accuracy and reliability of calculation has been poorly investigated. The methods for estimating the Fourier series truncating error are well developed for functions specified analytically, and the truncation error is expressed in terms of derivatives of the decomposed function. In calculations in Xray optics, we deal with functions digitized on a mesh, and classical methods of estimating the Fourier series truncation error are difficult to apply to functions specified in such a form. The Nyquist principle considers the necessary digitizing frequency of a function, but this is not an estimate of the Fourier series truncating error. In addition, studies show that various focal spot characteristics have significantly different sensitivities to the quality of the calculation (Kshevetskii et al., 2016). Therefore, considering the accuracy of calculating specific parameters instead of the calculation accuracy on the whole might be more appropriate. The function describing the electric field of the Xray wave passing through a system of several tens of lenses becomes rapidly oscillating and requires use of a very dense numerical mesh for its digitization and many Fourier series terms in calculations. In addition, it has been proved that, for a larger number of lenses, the faster the oscillating wavefield after the lenses, and the more Fourier series terms have to be taken into account for an acceptable description (Kshevetskii et al., 2016). These circumstances show the limited possibilities of Fourier methods and stimulate the development of new methods for solving the problems of Xray optics. New methods should provide high accuracy of calculations in complex cases and should be equipped with tools for controlling the accuracy and reliability of calculations (Kshevetskii et al., 2016; Kshevetskii & Wojda, 2015).
In this paper, a method for calculating the propagation and focusing of Xray waves based on the use of oriented Gaussian beams is presented. The idea of using Gaussian beams in optics is not new. The technique of using beams in optics has been given elsewhere (Kogelnik, 1965; Keller & Streifer, 1971; Deschamps, 1972; Deschamps et al., 1983). In later works (Chabory et al., 2005, 2010; Ghannoum et al., 2009), Gaussian beams (a Gaussian shooting algorithm) were used in radiolocation for finding the location of an object. Xray waves are characterized by an exceptionally high wave frequency; therefore, the method using Gaussian beams needs to be adapted to solve the problems of Xray optics.
Here, we set the incident wave in the form of a simple coherent wavefront, and we do not touch on the problem of real sources of Xray waves, partially coherent and nonstrictly monochromatic. However, Gaussian beams can also be used to study partially coherent radiation beams, as exemplified by the Gaussian Shell Model (Deschamps, 1972; Deschamps et al., 1983). There exists a numerical method based on the analytical method (Bahrdt, 2007), which resembles to some extent the method presented in this paper. It does not use Gaussian beams but allows simulation of electric field propagation `along geometrical rays'.
In Xray optics, the paraxial wave equation (Goodman, 1996; Kohn, 2012) is now widely used for calculating the electric wavefield. The paraxial equation is approximately derived from the Helmholtz equation (Levy, 2000; Ishimaru, 1991; Babich & Buldyrev, 1991a,b). Gaussian beams are exact solutions of the paraxial equation. Since the widths of the beams are much larger than the Xray wavelength, Gaussian beams can also be treated as almost exact solutions of the Helmholtz equation. As we shall see, practically all solutions of the problems of Xray optics can be constructed in the form of specially selected sums of Gaussian beams. We use oriented Gaussian beams, that is, Gaussian beams propagating at certain angles to the optical axis of the optical system. The condition that the propagation directions of the beams have to be orthogonal to the phase front of the simulated wave determines the inclination angles.
The developed approach combines properties of geometric optics (the Gaussian beams have clear directions of wave propagation and are localized in transverse directions) and of wave optics (the Gaussian beams are almost exact solutions of the Helmholtz wave equation). The method is well suited for solving problems with many lenses, and allows us to highly accurately calculate not only the focusing but also the finestructure effects on the focal plane that result from Xray diffraction at the edges of lens apertures.
Also, the paper proposes tools for controlling the calculation accuracy and shows the possibility of calculating various focusing parameters with the required accuracy. It is revealed that the halfwidth of the focal spot and, in particular, the amplitude of the electric field on the focal plane are very sensitive to the calculation quality. The position of the best focus can be obtained even with simple calculations.
2. Basic equations and starting points
The propagation of a monochromatic electromagnetic wave with frequency ω_{0} in a medium with complex n = 1 − δ + iβ is described by the Helmholtz equation (Levy, 2000; Ishimaru, 1991; Babich & Buldyrev, 1991a,b),
where E is the electric field, r = (x, y, z), k_{0} is the wavenumber, c is the speed of light, β is an and δ is a refractive decrement (β and δ are nonnegative).
More convenient for practical applications, a simplified wave equation can be derived from the Helmholtz equation (1) (Leontovich, 1944). Let us consider the case when the wave propagates along the xaxis, and the characteristic scales l_{y}, l_{z} of the wave along the axes y and z are much larger than the l_{x} of the wave along the xaxis: , . In this case, equation (1) can be greatly simplified. If we substitute
into equation (1) we obtain the equation for the function A(r). We suppose that A(r) varies only slowly with variable r. Neglecting the small term of a second derivative of the function A(r) with respect to the xvariable, we arrive at the paraxial equation
Equation (3) was first proposed by Leontovich in 1944 (Leontovich, 1944) for describing the propagation of a monochromatic electromagnetic wave within the paraxial approach. Equation (3) is often called the paraxial equation.
Equations (1) and (3) describe the propagation of monochromatic waves with frequency ω_{0}. Considering quasimonochromatic waves characterized by a set of frequencies ω close to ω_{0}, we can construct the solution of such a problem in the form of a superposition of waves with different ω. Thus, a solution of the general problem of propagation of nonmonochromatic waves still leads to equations (1) or (3).
3. The propagation of Xrays in a vacuum
Let us consider equation (3) for the case of a vacuum, that is, for complex n = 1, then equation (3) takes the following form,
This equation has important applications in optics in general (Goodman, 1996; Levy, 2000), where it describes the propagation of electromagnetic waves, and in acoustics (Babich & Buldyrev, 1991a,b). Equation (4) has an exact partial solution in the form of a Gaussian beam (Kogelnik, 1965),
Here σ is a real parameter, which controls the Gaussian beam width. At x = x_{0}, the solution (5) takes the form of a usual Gauss function,
We can consider (6) as a boundary condition given at the plane x = 0 to the equation (4). In this case, the formula (5) gives a solution to equation (4) that meets the boundary condition (6). Numerical examples of solutions of equation (4) with other boundary conditions can be found elsewhere (Kshevetskii et al., 2016; Kshevetskii & Wojda, 2015).
Applying the formula (5) to the formula (2), we obtain an approximate solution
of the Helmholtz equation (1). In our investigation, the condition (σk_{0})^{2} ≃ 10^{9} ≫ 1 is satisfied, and this means that the function (7) satisfies the Helmholtz equation (1) with very high accuracy; this is an almost exact solution of equation (1). Therefore, we shall often name (7) as a solution of the Helmholtz equation, for convenience's sake. We also shall call for convenience the solution (7) of the Helmholtz equation a Gaussian beam. The solution (7) of the Helmholtz equation (1) describes the wave propagation along the OXaxis; we assume that this axis coincides with the optical axis of the optical system.
Waves in Xray optics can propagate at small angles to the OXaxis. Therefore, it is desirable to generalize the particular solution (7) of the Helmholtz equation (1) to the case when the wave propagates at an arbitrary angle to the OXaxis. Let the vector indicate the direction of wave propagation and let the vector be directed along the OZaxis. We denote = , = . Then the generalization of the formula (7) to the case when the wave propagates along the direction indicated by the vector is written as follows,
Here r_{0} = (x_{0}, y_{0}, z_{0}), = . Since we are only interested in Gaussian beams propagating at small angles to the OXaxis, then with high accuracy k_{x} ≃ k_{0}. In the general case, however, k_{y} ≠ 0 and k_{z} ≠ 0. As before, σ is the parameter that determines the Gaussian beam width.
Since the formula (8) yields an almost exact particular solution of the Helmholtz equation (1) in the case (k_{0}σ)^{2} ≫ 1, there is a temptation to construct a general solution of equation (8) of the problem of Xray propagation in the form of a superposition of particular solutions of the kind (8), with appropriately chosen parameters.
We introduce in space a system of lines given by the equations (y = y_{m}, z = z_{l}), where y_{m} and z_{l} are real numbers and y_{m+1} = y_{m} + h, z_{l+1} = z_{l} + h, and h is the step of the mesh κ. Also, we denote each of the introduced lines by the pair (y_{m}, z_{l}).
Suppose we know the digitized electric wave field E on the plane x = x_{0}. That is, on the plane x = x_{0}, we have the mesh κ of points (y_{m}, z_{l}) where we know the values of the field E. We denote E_{0ml} = E(x_{0}, y_{m}, z_{l}). The problem is to find out the field E in the whole space. In principle, at each mesh point (y_{m}, z_{l}) we can specify its own direction of the Gaussian beam propagation determined from physical considerations. The general solution of the Xray wave propagation problem
Here = (x_{0},y_{m},z_{l}), is a wavevector of the (m,l)beam, C_{ml} is the factor that determines the (m,l)beam amplitude, and is a unit vector determining the propagation direction of waves in the (m,l)beam.
We have to reasonably determine the wavevectors , the multipliers C_{ml}, the unit vectors and the parameter σ of Gaussian beams. First of all, we define the wavevectors of the Gaussian beams used. We represent the Xray wavefield in the plane x = x_{0} in the neighborhood of the point = (x_{0},y_{m},z_{l}) in the form E = , where ϕ is a real local phase of the wave. Expanding the wave phase ϕ in a Taylor series in the neighborhood of the point , we obtain
We see that the definition of the local wavevector by the relations
is reasonable. From equation (4) it follows that k_{ml,x} can be represented as
Obviously, the unit vectors in (9) are directed along the vector, that is, = . It is also obvious that = , = . On the κ grid, the components of the vector can be approximated by the finite differences k_{ml,y} ≃ , k_{ml,z} ≃ , where h is the mesh step.
In our work, the Xrays propagate at very small angles to the OXaxis and for all points of the plane x = x_{0} the approximate relations ≃ 0 are valid because the vectors are almost orthogonal to the plane x = x_{0}; respectively, we obtain
We consider on the plane x = x_{0} an arbitrary closed star domain Ω containing the point of the grid and we find the coefficients C_{ml} from the condition
The relation (13) is assumed to be valid for any domain Ω whose diameter is much greater than both h and σ. From (13) then follows
From (14) it follows that C_{ml} ≃ E_{0ml}h^{2}/(2πσ^{2}), where E_{0ml} ≃ E(x_{0}, y_{m}, z_{l}). The lefthand side of the formula (14) is obtained by applying the usual quadrature formula to (13) for calculating the integral over the region Ω. When calculating the integrals on the right side of (14), assuming the beam width σ to be small compared with the diameter of the Ω region, we used the fact that the Gaussian functions decrease quickly with increase of [(y − y_{m})^{2} + (z − z_{l})^{2}], and we replaced the integrals over Ω of the Gaussian functions with integrals over the entire plane R^{2}.
In the derivation (14) we used the simplifying approximation . In a more general case, the wave can propagate at a small angle to the OXaxis. Understanding the principle of derivation of the formula for the factors C_{ml}, we below formulate the reasoning that does not use the assumption . We draw through the point a plane perpendicular to the vector . We introduce a special local coordinate system. Let the axis be along the direction indicated by the vector , the axis be along = and the axis be along = . We have taken into account here that the vector may be nonparallel to the optical axis of the system. We consider an arbitrary closed starshape domain in the plane containing the point of the plane . We assume that the domain is not very large; in this case, ≃ is fulfilled for all points (m,l) of the mesh inside .
We repeat the above derivation for Ω′ and obtain for C_{ml} a formula that coincides with the previous one,
but the new derivation of the formula (15) is more accurate and takes into account that the vectors may be nonparallel to the optical axis, but can form small angles with the optical axis.
In Fig. 1 for the planar case, an example of approximation of the function E(x_{0},y) = (E_{0} = 1.6 × 10^{7} V m^{−1}, α = 100 µm, k_{0} = ω_{0}/c, ω_{0} = 6π × 10^{18} Hz) on the line x = x_{0}, using a linear combination of the functions
is shown. The coefficients of the linear superposition are determined by the formula C_{m} = E_{0m}h/[(2π)^{1/2}σ], the components of the wavevector k_{m} are determined by the formula (11), r_{m} = (x_{0}, y_{m}), k_{m} = k_{m}, = . Some difference between these formulas and the above formulae (9) and (15) is explained by the fact that in this example we have only one transverse coordinate. We see that the system of functions (16) used allows us to approximate a given function very accurately on the line x = x_{0}. Since we apply a superposition of almost exact solutions of the Helmholtz equation, which satisfy the boundary condition with high accuracy, this method allows solving the problem with high accuracy.
4. Oriented Gaussian beams for describing the Xray propagation through optical elements
Now consider the propagation of Xrays through lenses. To begin with, the consideration of the propagation of Xrays through a single lens will suffice. Xrays penetrate any material and propagate in any material with very little reflection and absorption. This allows the boundary condition of the continuity of the function E to be used at the lens–air boundaries (Kohn, 2002; Lengeler et al., 1999). Materials used in Xray optics for lens manufacturing have such small parameters β and δ that equation (3) in the case of propagation of an Xray wave inside the lens can be represented in the following form,
The substitution
transforms equation (17) into a standard paraxial equation,
The method of solving equation (19), and also the Helmholtz equation (1), is described in Section 2. We can use the formulas (18), (5), (7), (10) and (11) for obtaining the solution to equation (1) to describe the propagation of an Xray wave through a lens.
The formulas (18) and (19) show that if an Xray wave encounters a lens with a local thickness d(y, z) then at the exit from the lens at the lens surface the Xray wave obtains an additional phase ϕ(y, z) = k_{0}(−δ + iβ)d(y, z), which is also equivalent to the fact that the function describing the electric field obtains the factor . The scheme explaining the local thickness of the lens is shown in Fig. 2.
Since the lens thickness is much less than the focal length, using the thin lens approximation is permissible and convenient. Within this approximation we take the lens center at the coordinate x_{C} of the lens on the OXaxis and assume that the lens lies on the plane x = x_{C}. Let the notations (x_{C} − 0) and (x_{C} + 0) call the limit on the left and the limit on the right to the point x_{C}; these notations will be used for writing the wave fallen on the thin lens from the left and the wave passed through the thin lens to the right. When an Xray wave passes through a thin lens, the wave phase (and, correspondingly, the electric field) changes abruptly. We suppose that a narrow Gaussian beam G(r, x_{C} − 0, y_{m}, z_{l}) described by the formula (7) falls onto the lens from the left. Also, let the width σ of a Gaussian beam be much smaller than the lens aperture. At the lens output at the lens surface, we obtain the wave
We have expanded the lens thickness in a Taylor series in the neighborhood of (y_{m}, z_{l}) and have taken into account linear terms of the series. Since the lens is considered as a thin one, the wave parameters changed a lot when the wave passed through the lens.
Now we consider the same question in terms of oriented Gaussian beams (8). The Xray wave , having r_{ml} = (x_{C}, y_{m, }, z_{l}) and the unit vector e_{1, ml} indicating the direction of wave propagation parallel to the OXaxis, fell on the lens from the left. At the output, after the lens, we obtained the wave . The parameters of this wave are determined by the formula (10), which leads to
Thus, the wave has changed its propagation direction and phase, but otherwise has remained unchanged. The addition to the wavevector in (21) has been generated by the factor before G(r, x_{C} + 0, y_{m}, z_{l}) in (20).
Now we generalize the obtained result to the case when the wave in the form of a superposition of narrow Gaussian beams falls from the left on the plane x = x_{C} containing a lens. Then after the lens we obtain the wave
This wave is a superposition of oriented Gaussian beams in which the parameters and are determined by the formula (21) and
In the approach used, the formula (22) describes very accurately the propagation of Xray waves in the air after the lens. High accuracy of the formulas (22), (23) and (21) is due to the lens symmetry, so that the lens center is clear. Also, the explicit formulas for very narrow Gaussian Xray beams make it possible to calculate explicitly the phase shifts with allowance for fine details.
In the derivation of the formulas (21) and (23), we assumed that the Xray wave propagates in the lens along the optical axis. This is a justified approximation, since Xray waves have a great penetrating power and practically do not change their propagation direction. Nevertheless, if we use a lot of lenses, about a hundred or more, the focal length of the optical system is relatively small. In this case, when calculating the local phase shifts, it is necessary to take into account the local slopes of the wave propagation direction with respect to the optical axis. Corresponding amendments can be taken into account; they will somewhat complicate the formulas. If we use only a few dozen lenses, the local slopes of the direction of wave propagation inside the lens are very small and practically do not affect the result.
We hope that the calculation method presented here will be convenient and adequate for studying the effects of lens defects (foreign inclusions, oxides, caverns) on focusing and imaging. A theoretical study of the influence of lens defects on focusing and imaging may yield more results than an experimental one, since we are unlikely to be able to obtain from the experiments the internal lens structure in detail. We believe that modeling the internal structure of lenses as random would be preferable. Obviously, the effect of one small defect in the lens will be small, but the effects of many defects in the lens can accumulate, especially if we apply a lot of lenses. Therefore, in calculations, high accuracy is required. When examining the impact of lens defects on the focusing and imaging, it is reasonable to not consider lenses as thin ones but to cut the lenses into layers. This will allow the introduction of a local wavevector inside each lens, and a better understanding of the wave behavior when the wave will meet a defect or multiple defects.
5. An estimate of the acceptable step h and the method error
When (k_{0}σ)^{2} ≃ 10^{9} ≫ 1, the Gaussian beams (8) satisfy well the Helmholtz equation (1). Therefore, we can consider them as practically exact solutions of the Helmholtz equation (1). A linear superposition of exact solutions (12) for x ≥ x_{0} gives an exact solution (9) of the Helmholtz equation (1). Thus, the error of the approximate solution is determined by the error in calculating and the error of approximating the boundary condition with a series (12). From the formula (10) it follows that the error in calculating the local wavevectors is O(h^{2}). This means that the error decreases in proportion to h^{2} when h decreases. As for the approximations (12)–(15), the analysis shows that the error in this case also decreases in proportion to h^{2}. Below we will apply the Runge rule (28) to check the calculation accuracy. The derivation of the Runge rule is given by Kshevetskii et al. (2016) and Kraus & Langer (2007). We have to substitute k = 2 into the Runge rule formula (28) because the approximation order is equal to two.
The problem of focusing Xrays is solved step by step. First we calculate the wavefield in the plane in front of the lenses. Then we compute the Xray propagation through the lenses. After that we calculate the wavefield after the lenses on the plane perpendicular to the optical axis. We finally use the wave electric field on the plane after the lenses for calculating the focusing. As shown in the work of Kshevetskii et al. (2016), in the case of many lenses the wavefield on the plane after the lenses is fastoscillating. The calculation quality is determined by how accurately we are able to describe this fastoscillating field.
To obtain an estimate of the acceptable step h, we assume that the phase of the wave incident on the lens system is approximately constant. The wave phase is changed significantly after the Xray wave has passed through the lenses. The formula
describes the dependence on the transverse coordinates of the thickness of each lens within the aperture. Here W_{sm} is the smallest distance between the parabolic surfaces of the lens and R is the curvature radius of the lens surfaces. According to the formulas (18) and (19), when the Xray wave has propagated through a system of N lenses, the Xray wave has the phase
dependent on the transverse coordinates. We have neglected the wave attenuation.
The approximation (9) assumes that the function (25) locally in some neighborhood of each point and at distances much larger the grid steps h can be approximated by a plane. Then we can introduce the local wavevector and the approximation (9) is effective. Obviously, the function (25) can be approximated in the neighborhood of each mesh point with a plane if within the aperture the conditions
are satisfied. The conditions (26) give a restriction on the step h of the mesh κ,
Actually we use in our calculations.
In the work of Kshevetskii et al. (2016), the question of a fundamental restriction to the step h was considered, and the condition was derived. Here a is the lens aperture. The condition (27) does not depend on the lens aperture and in practice it turns out to be substantially softer than the condition described by Kshevetskii et al. (2016). It allows the grid step h to be increased by about an order of magnitude. This possibility is explained by the fact that we use in (9) inclined Gaussian beams, which, due to their slopes, take into account fast oscillations on the OYZplane of the wave electric field. In the derivation of (27), the lenses are assumed to be ideal, having a parabolic shape. When calculating the Xray propagation through nonidealshaped lenses, the h step may have to be reduced for catching and investigating the effects of nonideal lenses.
A special issue concerns the influence of the edges of lens apertures on focusing and imaging. The derivatives of the n change abruptly at the edges of the lens apertures, and the wave phase φ also changes abruptly. Therefore, the condition (26) is violated at the edges of lens apertures, and even , tend to infinity. Let us consider all points (m, l) in the sum (9), which are in a thin ring of width ∊ containing the contribution of the edge of the lens aperture. The contribution of this thin ring into the total wavefield does not exceed E^{*}∊/σ^{2}. Here E^{*} is the maximum value of E in this ring. Since the width ∊ of the considered ring can be taken as arbitrarily small, the contribution of the edge of the lens aperture in the formula (9) to the general wave picture can be neglected. This conclusion can also be justified by the fact that the beam of Xray waves incident on the lens usually has an effective width that is smaller than the lens aperture size.
Nevertheless, the influence of the edges of lens apertures on the wave propagation and focusing exists. Formula (9) allows calculation with high accuracy of this influence on the wavefield in the focal plane; this will be demonstrated later.
5.1. Volume of calculations
The computational complexity (number of operations) of the presented method using Gaussian beams is estimated as M_{I}M_{F}. Here M_{I} is the number of points digitizing the electric field after the lenses. M_{F} is the number of points at which we calculate the electric field after the lenses at the distance we are interested in. After the Xrays have propagated through a system of several tens of lenses, the wavefield becomes fastoscillating, and a rather dense numerical grid is required for its digitization, which determines the number M_{I}. The number M_{F} is determined only by the desirable resolution of the picture of interest.
At present, Fourier methods are often used for calculations in Xray optics. It is known that the computational complexity of the FFT is estimated as Mlog_{2}(M), where M is the number of points used in the FFT method (Cooley et al., 1967). The FFT method is very fast. However, the FFT technology based on FFT is somewhat lacking in numerical robustness and stability (Chubar et al., 2017). Currently, work is underway to develop software that allows automatic control of the accuracy of calculations using FFT (Chubar et al., 2017). Therefore, the question of calculation accuracy control is an important part of modern research related to the computational problems of Xray optics. When the number of lenses increases, the frequency of oscillations of the wavefield in the plane behind the lenses perpendicular to optical axis increases fast. Therefore, the number of Fourier terms that have to be taken into account also increases fast. Standard methods for estimating the error of truncating Fourier series are applied to functions specified analytically, and they are expressed in terms of derivatives of the decomposed function. In Xray optics, we deal with functions digitized on a mesh on a plane, and applying standard formulas for estimating the truncating error is difficult. In addition, studies on the calculation accuracy of Xray focusing have revealed that the sensitivity of various parameters characterizing the focusing to the calculation quality is very different (Kshevetskii et al., 2016). Therefore, evaluating the accuracy of calculations is a serious problem and the calculation accuracy issue should be formulated for specific characteristics of focusing, and not as a general question.
Thus, the method developed in this paper has some advantages in cases where we need to guarantee the reliability and accuracy of calculations. At the same time, the developed method is quite universal and uncomplicated.
6. Technique for modeling the propagation and focusing of Xray waves in terms of oriented Gaussian beams
The propagation and focusing of Xray waves is investigated for the following conditions. We consider the situation that some source of Xray waves with wave frequency ω_{0} = 6π × 10^{18} Hz produces in the plane x = x_{0} (x_{0} = 0 m) a Gaussian distribution of electric wavefield intensity with a width = 100 µm. For beryllium lenses and given frequency of Xray waves, β = 3.1801 × 10^{−10} and δ = 2.2156 × 10^{−6}. Then the Xrays propagate further to a distance of 40 m, and then the waves propagate through a system of 30 parabolic lenses made of beryllium with a curvature radius of 50 µm and smallest width of 30 µm. After the Xrays have propagated through the lenses, the Xrays are focused, and we study them in the focal plane.
The general scheme of calculating the wave propagation is shown in Fig. 3. The problem is solved as follows. The specified wave electric field from a source of monochromatic Xray waves is digitized on the mesh on the plane x = x_{0}, and then the Xray propagation to the region x > x_{0} is calculated using formulas (9), (8), (10) and (15). The resulting wave solution is digitized again on the mesh on the plane x = x_{1}. The wave propagation through the first lens is calculated using the formulas (22), (21) and (23). The resulting wavefield after the first lens is digitized on the mesh on the plane x = x_{2}. The wave propagation through the second lens is calculated using formulas (22), (21) and (23). The obtained wavefield after the second lens is digitized on the mesh on the plane x = x_{3}. The formulas (22), (21) and (23) are used for calculating the wave propagation through the third lens, and so on. We repeat this until the wavefield is calculated after the last lens corresponding to the coordinate x_{30}. The resulting wavefield after the last lens is used to find the position of the best focus. It is digitized on the mesh on the plane spaced at a distance corresponding to the position of the best focus, which we denote x_{31}. Finally, we investigate the behavior of the digitized electric field on the focal plane x = x_{31}.
6.1. Results of calculations
The calculations are performed for a onedimensional case. We consider the case where the wave electric field of the Xray beam from the synchrotron far in front of the lenses in the plane x = x_{0} is described by the formula E(0,y) = (where E_{0} = 1.6 × 10^{7} V m^{−1} and α = 100 µm). Table 1 presents the calculation results obtained by approximating the wave electric field with the help of 1000, 2000 and 4000 oriented Gaussian beams, respectively. The widths σ = h of these Gaussian beams are equal to 1 µm, 0.5 µm and 0.25 µm, respectively. The E = A function takes the largest value at the distance of 0.36658 m after the lens system. Table 1 lists the values of E and the FWHM at the distance equal to the focal length (x = x_{31}).

Using the data in Table 1 and applying the Runge rule, we estimate the calculation errors of the results given in Table 1. The Runge rule (Kraus & Langer, 2007) allows to estimate computational errors of values independent of coordinates, such as the FWHM and the focal length. The electric field depends on the coordinates, but the maximum value of E in the focal spot is independent of the coordinates. Therefore, we can control the accuracy of calculating this physical quantity using the Runge rule, too. The error Q of computing the value of Z within the Runge rule can be estimated using the approximate equality
where Z_{h1} and Z_{h2} are approximate values of Z obtained with the space steps h_{1} and h_{2}, respectively, and . The derivation of the formula (28) is given by Kshevetskii et al. (2016).
According to the analysis performed in Section 5, the oriented Gaussian beam method is a method of the second order of accuracy with respect to the step h. Therefore, we have to take k = 2 in (28). Applying the Runge rule (28) to the results of computations for the case h_{1} = σ_{1} = 1 µm and and for the case h_{2} = σ_{2} = 0.5 µm, we have found out that the calculation error of the FWHM is approximately equal to 1.5% of its value. The error of computation of a maximal value of the wave electric field in the focal spot is 9% of its value. These results show that the amount of M = 1000 Gaussian beams with h = σ = 1 µm is not sufficient for a highprecision calculation of the maximum E in the focal spot. The results of the calculation for the case of h_{1} = σ_{1} = 0.5 µm and h_{2} = σ_{2} = 0.25 µm show that the FWHM calculation error is approximately 0.4%, and the error in calculating the maximum value of the wave electric field is 0.6%. Therefore, we see that the choice of h = σ = 0.5 µm gives highly accurate results.
Figs. 4 and 5 shows the behavior of the function E at the focal plane x = x_{31}. The curve shape in the center of Fig. 4 is close to the Gaussian function, and this is a socalled focal spot. Fig. 5 shows the fine structure of the function E, formed in some neighborhood of the focal spot. This fine structure consists of a set of Gaussianlike peaks. The maximum value of the function E in this fine structure is at least two orders of magnitude less than the maximum value of the function E in the focal spot. The values of the local maxims E in the fine structure decrease with increasing distance from the center. To make sure that this fine structure is not a computational artifact, but is a physical effect, in Fig. 5 the results of calculations with the ordinary and double steps h are pictured. The results of calculations with ordinary and double steps coincide with high accuracy. This indicates a physical nature of the calculated fine structure.
Comparing the results obtained with different widths α of the initial Gaussian distribution of the electric wavefield at x = x_{0}, we have found that the amplitudes of the Gaussianlike peaks of the fine structure in the focal plane x = x_{31} correlate with E at the edges of the lens apertures. This observation shows that the fine structure of the wave electric field in the focal plane has arisen as a result of Xray diffraction at the edges of the lens apertures.
7. Conclusions
A highprecision method for solving Xray optics problems, based on the use of superposition of oriented Gaussian beams, is developed. A Gaussian beam is an exact solution to the paraxial equation of optics, and for appropriate beam parameters (the wavelength is much smaller than the beam width) it can be interpreted as an almost exact solution of the Helmholtz equation. In an oriented Gaussian beam, the waves propagate at a given angle to the optical axis. Some relations that determine the optimal choice of Gaussian beam parameters and formulas for linear superposition coefficients are obtained. Also, the formulas that determine the variation of the parameters of oriented Gaussian beams and the superposition coefficients when the Xray wave propagates through the lens are derived. The propagation of a Gaussian beam through a lens leads to changes in the propagation direction of the wave and the complex amplitude of the electric field in the beam.
For the method presented, tools for controlling the accuracy and reliability of calculations are developed. It is shown that different characteristics of the focal spot have a significantly different sensitivity to the computation quality. Among the parameters analyzed (, FWHM, focus length), the maximum amplitude of the electric field in the focal spot is the most sensitive and exacting to the computation quality. A clear formula for the mesh step (the number of Gaussian beams used) required for highprecision calculations is derived.
As an application, the problem of focusing of Xrays by a system of 30 beryllium lenses in an Xray microscope was solved. The high accuracy of the presented method made possible calculating not only the Xrays focusing but also the fine structure of the wave electric field on the focal plane. The fine structure arises from the diffraction of Xrays at the edges of the lens apertures.
The developed highaccuracy modeling method is quite universal and is designed to solve various problems of Xray optics and, in particular, to study the effect of various defects in lenses (cavities, foreign inclusions, oxides) on focusing and imaging. The presented method allows solving problems with complex lens shapes or calculating the propagation of Xray waves through a sample with a complex internal structure, and with control of the accuracy of calculations.
References
Babich, V. M. & Buldyrev, V. S. (1991a). Asymptotic Methods in ShortWavelength Diffraction Theory. Oxford: Alpha Science International. Google Scholar
Babich, V. M. & Buldyrev, V. S. (1991b). ShortWavelength Diffraction Theory: Asymptotic Methods. Berlin: Springer Verlag. Google Scholar
Bahrdt, J. (2007). Phys. Rev. ST Accel. Beams, 10, 060701. CrossRef Google Scholar
Chabory, A., Sokoloff, J., Bolioli, S. & Combes, P. F. (2005). C. R. Phys. 6, 654–662. CrossRef CAS Google Scholar
Chabory, A., Sokoloff, J., Bolioli, S. & Elis, K. (2010). Proceedings of the 4th European Conference on Antennas and Propagation (EuCAP'2010), 12–16 April 2010, Barcelona, Spain, pp. 1–5. Google Scholar
Chubar, O., Rakitin, M., ChenWiegart, Y., Chu, Y. S., Fluerasu, A., Hidas, D. & Wiegartr, L. (2017). Proc. SPIE, 10388, 1038805. Google Scholar
Cooley, J., Lewis, P. & Welch, P. (1967). IEEE Trans. Audio Electroacoust. 15, 76–79. CrossRef Google Scholar
Deschamps, G. (1972). Proc. IEEE, 60, 1022–1035. CrossRef Google Scholar
Deschamps, J., Courjon, D. & Bulabois, J. (1983). J. Opt. Soc. Am. 73, 256. CrossRef Google Scholar
Ghannoum, I., Letrou, C. & Beauquet, G. (2009). Proceedings of RADAR 2009 – International Radar Conference `Surveillance for a Safer World', 12–16 October, 2009, Bordeaux, France. Google Scholar
Goodman, J. W. (1996). Introduction to Fourier Optics, 2nd ed. New York: McGrawHill. Google Scholar
Ishimaru, A. (1991). Electromagnetic Wave Propagation, Radiation and Scattering. Englewood Cliffs: Prentice Hall. Google Scholar
Keller, J. B. & Streifer, W. (1971). J. Opt. Soc. Am. 61, 40–43. CrossRef Google Scholar
Kogelnik, H. (1965). Appl. Opt. 4, 1562–1569. CrossRef Google Scholar
Kohn, V. G. (2002). J. Environ. Toxicol. Publ. Heal. Lett. 76, 701–704. Google Scholar
Kohn, V. G. (2012). J. Synchrotron Rad. 19, 84–92. Web of Science CrossRef CAS IUCr Journals Google Scholar
Kraus, J. & Langer, U. (2007). Editors. Lectures on Advanced Computational Methods in Mechanics. p. 164. Berlin: de Gruyter. Google Scholar
Kshevetskii, S., Wojda, P. & Maximov, V. (2016). J. Synchrotron Rad. 23, 1305–1314. CrossRef CAS IUCr Journals Google Scholar
Kshevetskii, S. P. & Wojda, P. (2015). Math. Appl. 43, 193–206. Google Scholar
Lengeler, B., Schroer, C., Tümmler, J., Benner, B., Richwin, M., Snigirev, A., Snigireva, I. & Drakopoulos, M. (1999). J. Synchrotron Rad. 6, 1153–1167. Web of Science CrossRef IUCr Journals Google Scholar
Leontovich, M. A. (1944). Izv. Akad. Nauk. SSSR Ser. Phys. 8, 16. Google Scholar
Levy, M. (2000). Parabolic Equation Methods for Electromagnetic Wave Propagation. Londo: The Institution of Electrical Engineers. Google Scholar
Roth, T., Alianelli, L., Lengeler, D., Snigirev, A. & Seiboth, F. (2017). MRS Bull. 42, 430–436. CrossRef Google Scholar
Snigirev, A., Kohn, V., Snigireva, I. & Lengeler, B. (1996). Nature, 384, 49–51. CrossRef CAS Web of Science Google Scholar
© International Union of Crystallography. Prior permission is not required to reproduce short quotations, tables and figures from this article, provided the original authors and source are cited. For more information, click here.