FD-BPM-based rotational symmetry light field modeling method under cylindrical coordinate system

Abstract

The invention belongs to the technical field of optical calculation and simulation, and particularly relates to an FD-BPM-based rotational symmetry light field modeling method under a cylindrical coordinate system. According to the method, the model of the optical element is created by taking Matlab as a platform, and the method comprises the steps of setting the morphology and the refractive index of the optical element and the wavelength, the amplitude and the wavefront of an incident light field. And calculating a near-field light field and a light field at any position in a far field of the rotationally symmetrical optical element based on a beam propagation method and a finite difference method with a boundary condition. According to the method, various optical information modulated by the optical element can be obtained, including intensity and phase information of any position in a far field and a near field, so that the transmittance, the focusing efficiency, the focal depth, the focal spot size and the like can be calculated. According to the method, the accuracy and the calculation efficiency of optical simulation are greatly improved, the research and development period can be remarkably shortened in design and optimization of large-aperture and short-wavelength focusing optical elements, and the experiment cost is reduced.

Description

Technical Field

[0001] The present invention belongs to the technical field of optical computing and simulation, and particularly relates to a method for modeling a rotationally symmetric light field in a cylindrical coordinate system. Background Art

[0002] In various imaging systems, optical elements such as focusing and diffraction elements are crucial components. Among them, optical performance simulation and analysis calculation are the primary and key steps in the design and research and development of optical elements, and play an important guiding role in the design optimization, processing preparation, and performance evaluation of optical elements. With the improvement of the processing and design levels of optical devices, the previous simple calculation method of thin grating approximation is no longer applicable, especially for sub-wavelength devices, devices with large aspect ratios in topography, and other devices that need to consider waveguide effects and volume effects.

[0003] In recent years, the construction pace of synchrotron radiation light sources in China has accelerated, and technologies such as X-ray microscopic imaging and X-ray probe microarea analysis have developed vigorously, providing powerful characterization tools and richer research methods for frontier disciplines and basic disciplines such as biology, medicine, archaeology, chemistry, and physics. X-ray focusing imaging devices represented by zone plates are gradually developing towards nanoscale resolution (feature size). Considering the characteristics of short X-ray wavelengths and small feature sizes of X-ray optical devices, although traditional simulation methods such as FDTD are still feasible in near-field calculations, the limitation of computing power makes them inadequate in far-field calculations. In addition, in the design of X-ray optical devices, most of the current mainstream X-ray zone plate devices, such as Kinoform topography zone plates, Kinoform-Fresnel composite zone plates, zone plate multiplication rectangular zone plates, stepped or trapezoidal zone plates, and even free-form rotationally symmetric zone plates, are designed based on the Kirz theory to complete the topography design under the thin grating approximation. For thick lenses, soft X-rays, and X-ray optical devices in high-resolution cases, the Kirz theory is not perfect and is difficult to play an effective guiding role in design.

[0004] In addition, sometimes the wavelength is much smaller than the size of the component. Common optical simulation software such as FDTD and Comsol often have problems such as large computational volume, low accuracy, and long calculation time during simulation. The Takagi-Taupin description and the beam propagation method (BPM) are common methods for studying short-wavelength optical devices, but the Takagi-Taupin equation is an analytical expression under ideal conditions and cannot analyze the focusing situation of actual optical elements. Most of today's BPM methods are also based on Fourier transform or generalized Fourier transform and cannot well solve the problem of its inherent periodic boundary. A small number of BPMs using the finite difference method do not reasonably optimize and model the three-dimensional rotationally symmetric light field, resulting in the inability to calculate the three-dimensional light field or low calculation efficiency and huge computing power consumption.

[0005] Therefore, it is necessary to establish an optical calculation method for optical elements that conforms to the actual boundary, can introduce the actual topography of the optical elements, has a high calculation dimension, high precision, and fast calculation speed. Since general two-dimensional optical elements such as zone plates, circular Kinoform optical elements, etc. all have the characteristic of rotational symmetry, the present invention proposes to discretize the Helmholtz equation in the cylindrical coordinate system by the finite difference method, and transform the tridiagonal coefficient matrix linear equations into recurrence formulas for solution, reasonably introduce the central symmetry boundary condition and the transparent boundary condition, supplement the missing definite solution equations in the finite difference method, and obtain the finite difference beam propagation method that is fast, accurate, and does not require solving large-scale matrices. Summary of the Invention

[0006] The purpose of the present invention is to propose a modeling method for a rotationally symmetric light field in the cylindrical coordinate system that can introduce the actual topography of the optical element, has a high calculation dimension, high precision, fast calculation speed, and the result conforms to the actual boundary, so as to solve the problems encountered by other methods mentioned in the above background technology.

[0007] The rotationally symmetric light field modeling method in the cylindrical coordinate system provided by the present invention is based on the finite difference (FD) and the beam propagation method (BPM), as shown in Figure 17 shown, the specific steps are as follows:

[0008] (1) Set the incident light field and the optical element parameters; among them:

[0009] The incident light field parameters include the wavelength λ, amplitude, and waveform; the incident light field is used as the initial wave, represented by φ 0 (r), where r represents the radial coordinate in the cylindrical coordinate system, and is used to calculate the subsequent near field and far field;

[0010] The optical element parameters include the optical element thickness t, radius R, focal length f, topography function, and the complex refractive index of the material used for the optical element (including the real part representing the phase shift and the imaginary part representing the absorption); the optical element parameters are used to generate the optical element model to be calculated, and are used to calculate the subsequent near field, the light field φ t (r) on the exit surface, far field, and the light field φ f (r) on the focal plane;

[0011] (2) Use the slowly varying amplitude approximation: E(r,z) = φ(r,z)exp(-jn 0 kz), and rewrite the Helmholtz equation in the cylindrical coordinate system:

[0012]

[0013] into:

[0014]

[0015] where j is the imaginary unit, k is the wave vector modulus, n is the refractive index of the optical device material, and n 0 is the refractive index of the immersion environment of the optical device, generally considered to be the refractive index of vacuum or air, i.e., n 0 = 1;

[0016] (3) Set the number of sampling points along the radial direction, near-field propagation direction, and far-field propagation direction of the optical element:

[0017] r i = iΔr, z m = mΔz, (0 ≤ i ≤ I, 0 ≤ m ≤ M)

[0018] where Δr and Δz are the minimum sampling lengths respectively, and I and M are the maximum number of sampling points in the radial and propagation directions;

[0019] (4) Given φ i (r), v i+1 (r) is derived by the following steps:

[0020] (4.1) Using the Crank-Nicholson method, rewrite Equation (2) as:

[0021]

[0022] where:

[0023] represents the optical field at the point (iΔr, mΔz), and the coefficients are respectively:

[0024]

[0025]

[0026] where α is the absorption coefficient of the material for light in this wavelength band;

[0027] (4.2) Using the recurrence relation, assume:

[0028]

[0029] where the coefficients are:

[0030]

[0031] where,

[0032] The initial value of the recurrence formula is:

[0033]

[0034] From this, it can be deduced

[0035] (5) Calculate The specific process is as follows:

[0036] Let \(i = 0\) in equation (3), then we have:

[0037]

[0038] Considering the rotational symmetry of the optical field, we have:

[0039] \(\varphi(r,z)=\varphi(-r,z)\).

[0040] Then equation (13) can be written as:

[0041]

[0042] Combined with the coefficients listed in step (4), we have:

[0043]

[0044] (6) Calculate The specific process is as follows:

[0045] Apply the transparent boundary condition:

[0046]

[0047] Where:

[0048]

[0049] (7) Generate a high-precision optical field model, including:

[0050] Optical field intensity in the near field and far field: Intensity = E 2 ;

[0051] Phase distribution:

[0052] Transmittance: That is, the ratio of the light intensity Inten out emitted from the optical element to the light intensity Inten in of the incident light;

[0053] Focusing efficiency: That is, the sum of the optical field intensity distribution by radius at the focal plane divided by the light intensity of the incident light;

[0054] Depth of focus: That is, the ratio of the wavelength to the square of the numerical aperture;

[0055] Focal spot size: Generally characterized by measuring the full width at half maximum of the focal spot in the calculation results.

[0056] In the present invention:

[0057] The wavelength of the light field is visible light, X-ray, ultraviolet or infrared.

[0058] The material of the optical element is metal, inorganic material, organic material or composite material.

[0059] The optical element has a rotationally symmetric morphology, including a zone plate, a Kinofom optical element, a Fresnel lens, a convex lens or a concave lens, or an optical element composed of several combinations thereof.

[0060] Compared with the prior art, the beneficial effects of the method of the present invention are as follows:

[0061] First, the present invention uses the Crank-Nicholson finite difference method to discretize the Helmholtz equation in the cylindrical coordinate system. Considering the rotational symmetry of the optical device, the calculation dimension is reduced and the calculation efficiency is improved.

[0062] Second, the recurrence formula method is used to solve the problem of solving the linear equation system with a tridiagonal coefficient matrix, which greatly reduces the demand for memory space during the operation and speeds up the operation speed.

[0063] Third, the central symmetry boundary condition and the transparent boundary condition are used to supplement the missing definite solution conditions of the linear equation system, and the result conforms to the optical phenomenon at the real physical boundary.

[0064] Fourth, since the Crank-Nicholson finite difference method adopted is unconditionally convergent, in theory, this method can calculate the light field at infinity, that is, this method can uniformly calculate the near field and the far field.

[0065] Fifth, compared with the existing BPM method based on discrete Hankel transform (an optical modeling and calculation method based on Hankel transform and beam propagation method

CN113281900A

[0066] Figure 1 (a) is a schematic diagram of the light field calculated by the present invention.

[0067] Figure 1 (b) is a schematic diagram of the morphology of the rectangular Fresnel HSQ zone plate calculated in Example 1.

[0068] Figure 1 (c) is a topographical schematic diagram of the zone-multiplied rectangular Fresnel HSQ zone plate calculated in Example 2.

[0069] Figure 2 It is a distribution diagram of the far-field electric field intensity of the rectangular Fresnel HSQ zone plate in Example 1.

[0070] Figure 3 It is a normalized distribution diagram of the focused spot intensity of the rectangular Fresnel HSQ zone plate in Example 1.

[0071] Figure 4 It is a distribution diagram of the cross-sectional intensity of the focused spot of the rectangular Fresnel HSQ zone plate in Example 1.

[0072] Figure 5 It is a distribution diagram of the light field intensity inside the optical element of the zone-multiplied rectangular Fresnel zone plate in Example 2.

[0073] Figure 6 It is a distribution diagram of the far-field light field intensity of the zone-multiplied rectangular Fresnel zone plate in Example 2.

[0074] Figure 7 It is a normalized distribution diagram of the focused spot intensity of the zone-multiplied rectangular Fresnel zone plate in Example 2.

[0075] Figure 8 It is a distribution diagram of the cross-sectional intensity of the focused spot of the zone-multiplied rectangular Fresnel zone plate in Example 2.

[0076] Figure 9 It is a topographical schematic diagram of the cross-section along the radial direction of the kinoform-Fresnel composite zone plate in Example 3.

[0077] Figure 10 It is a distribution diagram of the intensity of the light field emerging from the optical device of the kinoform-Fresnel composite zone plate in Example 3.

[0078] Figure 11 It is a distribution diagram of the far-field light field intensity of the optical device of the kinoform-Fresnel composite zone plate in Example 3.

[0079] Figure 12 It is a normalized distribution diagram of the focused spot intensity of the optical element of the kinoform-Fresnel composite zone plate in Example 3.

[0080] Figure 13 It is a distribution diagram of the cross-sectional intensity of the focused spot of the optical element of the kinoform-Fresnel composite zone plate in Example 3.

[0081] Figure 14It is the distribution diagram of the far-field light field intensity of the zone-multiplied zone plate with a beamstop in Example 4.

[0082] Figure 15 It is the normalized distribution diagram of the focused spot intensity of the zone-multiplied zone plate with a beamstop in Example 4.

[0083] Figure 16 It is the distribution diagram of the cross-sectional intensity of the focused spot of the zone-multiplied zone plate with a beamstop in Example 4.

[0084] Figure 17 It is the flowchart of the present invention. Detailed implementation manners

[0085] The present invention will be further described below in conjunction with the accompanying drawings and embodiments, but the present invention is not limited to the examples. Any simple change to the calculation parameters in the embodiments belongs to the protection scope of the present invention.

[0086] Example 1: Calculate the focusing process of a rectangular Fresnel zone plate optical element at an energy of 100 eV. The specific steps are as follows:

[0087] (1) Set the energy of the incident light field to 100 eV, and the plane wave with unit amplitude is incident vertically. The material of the rectangular Fresnel zone plate optical element is HSQ, with a diameter of 100 μm, a thickness of 380 nm, the width of the outermost ring is 100 nm, the focal length is about 800 μm, and the morphology is an ideal rectangular zone distribution, as Figure 1 (b) shown.

[0088] (2) Use the beam propagation method and the finite difference method to calculate the light field in the zone plate. The sampling points along the radial direction of the optical element are 2048, the sampling points in the near-field propagation direction are 1024, and the sampling points in the far-field propagation direction are 1024.

[0089] (3) Continue to iterate and calculate the far-field light intensity distribution. The calculation result is as Figure 2 shown. It can be seen that there is a clear focusing light path and multiple secondary focal points, and there is slight direct light on the light path.

[0090] (4) Generate a high-precision light field model, and calculate that the focusing efficiency is 11%, the resolution is 105 nm, and the calculated focal spot results are as Figure 3 and Figure 4 shown. The waveform is relatively ideal, and the full width at half maximum is in line with the theoretical value.

[0091] Example 2: Calculate the focusing process of a zone-multiplied rectangular Fresnel zone plate optical element at an energy of 5.5 keV. The specific steps are as follows:

[0092] (1) Set the energy of the incident light field to 5.5 keV and the plane wave with unit amplitude is incident perpendicularly. The material of the rectangular Fresnel zone plate optical element template is HSQ, the material grown by ALD used in the zone multiplication technology is Pt, the diameter is 100 μm, the thickness is 900 nm, the width of the outermost ring is 30 nm, the focal length is about 13.6 mm, and the morphology is an ideal rectangular zone distribution, as Figure 1 (c) shown.

[0093] (2) Use the beam propagation method and the finite difference method to calculate the light field inside the zone plate. The number of sampling points along the radial direction of the optical element is 2048, the number of sampling points in the near-field propagation direction is 1024, and the number of sampling points in the far-field propagation direction is 1024. The calculation result of the light field inside the optical element is as Figure 5 shown.

[0094] (3) Continue the iteration to calculate the far-field light intensity distribution. Set multiple focal order selection apertures (OSA) in the far field, and the obtained calculation result is as Figure 6 shown. It can be seen that there is a clear focusing optical path and multiple secondary focal points, and there is slight direct light on the optical path.

[0095] (4) Generate a high-precision light field model. The calculated focusing efficiency is 13.25%, the resolution is 31.5 nm, and the calculated focal spot results are as Figure 7 and Figure 8 shown. The waveform is relatively ideal, and the full width at half maximum is consistent with the theoretical value.

[0096] Example 3: Calculate the focusing process of the Kinoform-Fresnel composite zone plate optical element at an energy of 500 eV. The specific steps are as follows:

[0097] (1) Set the energy of the incident light field to 500 eV and the plane wave with unit amplitude is incident perpendicularly. The material of the composite zone plate template is HSQ, the material grown by ALD used in the zone multiplication technology is Pt, the diameter is 10 μm, the thickness is 160 nm, the width of the outermost ring is 10 nm, the focal length is about 40 um, and the cross-sectional morphology along the radius is as Figure 9 shown.

[0098] (2) Use the beam propagation method and the finite difference method to calculate the light field inside the zone plate. The number of sampling points along the radial direction of the optical element is 2048, the number of sampling points in the near-field propagation direction is 1024, and the number of sampling points in the far-field propagation direction is 1024. The light field intensity distribution exiting from the optical device is as Figure 10 shown.

[0099] (3) Continue the iteration to calculate the far-field light intensity distribution. Set multiple focal order selection apertures (OSA) in the far field, and the obtained calculation result is as Figure 11 shown. It can be seen that the direct light of the part with the middle Kinoform morphology is relatively strong, and there are multiple secondary focal points on the optical path.

[0100] (4) Generate a high-precision optical field model, and the calculated focusing efficiency is 4%, and the resolution is 11.3 nm. The calculated focal spot results are as Figure 12 and Figure 13 shown. The waveform is relatively ideal, and the full width at half maximum (FWHM) is in line with the theoretical value.

[0101] Example 4: Calculate the focusing process of a zone-multiplying zone plate optical element with a central through-beamstop at an energy of 500 eV. The specific steps are as follows:

[0102] (1) Set the incident optical field energy to 500 eV, and a plane wave with a unit amplitude is incident vertically. The material of the compound zone plate template is HSQ, the material used in the zone-multiplying technology grown by ALD is Pt, the diameter is 10 μm, the thickness is 160 nm, the width of the outermost ring is 10 nm, the diameter of the beamstop is 3 μm, and the focal length is about 40 μm.

[0103] (2) Use the beam propagation method and the finite difference method to calculate the optical field inside the zone plate. The sampling points along the radial direction of the optical element are 2048, the sampling points in the near-field propagation direction are 1024, and the sampling points in the far-field propagation direction are 1024.

[0104] (3) Continue to iterate, calculate the far-field light intensity distribution, and set multiple focus order selection apertures (OSAs) in the far field. The calculated results are as Figure 14 shown. It can be seen that the through-light in the center of the optical path is blocked by the beamstop, and there are multiple secondary foci in the optical path.

[0105] (4) Generate a high-precision optical field model, and the calculated focusing efficiency is 3.7%, and the resolution is 11.8 nm. The calculated focal spot results are as Figure 15 and Figure 16 shown. The waveform is relatively ideal, and the full width at half maximum (FWHM) is in line with the theoretical value.

Claims

1. A method for modeling a rotationally symmetric light field in a cylindrical coordinate system, characterized in that: It is based on finite difference (FD) and beam propagation method (BPM). The specific steps are as follows: (1) Setting the incident light field and optical element parameters; wherein: The incident light field parameters include wavelength λ, amplitude, and waveform; the incident light field is used as the initial wave and is represented by φ0(r), where r represents the coordinate along the radial direction in the cylindrical coordinate system and is used to calculate the subsequent near field and far field; The optical element parameters include the optical element thickness t, radius R, focal length f, morphology function, and complex refractive index of the material used for the optical element (including the real part representing the phase shift and the imaginary part representing the absorption); the optical element parameters are used to generate the optical element model to be calculated, and are used to calculate the subsequent near field and exit surface light field φ t (r), far field, focal plane light field φ f (r); (2) Using the slowly varying amplitude approximation: E(r,z) = φ(r,z)exp(-jn0kz), the Helmholtz equation in the cylindrical coordinate system is: Rewritten as: Wherein, j is the imaginary unit, k is the wave vector modulus, n is the refractive index of the optical device material, and n0 is the refractive index of the environment in which the optical device is immersed; (3) Set the number of sampling points along the radial direction, near-field propagation direction, and far-field propagation direction of the optical element: r i =iΔr,z m =mΔz,(0≤i≤I,0≤m≤M) Among them, Δr and Δz are the minimum sampling lengths, I and M are the maximum number of sampling points in the radial and propagation directions, respectively; (4) Given φ i (r), φ i+1 (r) is derived from the following steps: (4.1) Using the Crank-Nicholson method, equation (2) can be written as: in: Represents the light field at point (iΔr, mΔz), the coefficients are: Among them, α is the absorption coefficient of the material for light in this band; (4.2) Using the recursive relation, we can conclude that: The coefficients are: in, The initial value of the recursive formula is: From this we can infer (5) Calculation The specific process is: Let i = 0 in formula (3), then: Considering the rotational symmetry of the light field, we have: φ(r,z)=φ(-r,z). So formula (13) can be written as: Combining the coefficients listed in step 4), we have: (6) Calculation The specific process is: Apply a transparency boundary condition: in: k p =|κ r |+ik im (18) (7) Generate a high-precision light field model, including: Light field intensity in near field and far field: Intensity = E 2 ; Phase distribution: Transmittance: That is, the light intensity emitted from the optical element Inten out The intensity of the incident light in The ratio of Focusing efficiency: That is, the light field intensity at the focal plane is summed according to the radius distribution and divided by the intensity of the incident light; Depth of Focus: That is, the ratio of wavelength to the square of numerical aperture; Focus spot size: It is usually characterized by measuring the half-height width of the focus spot in the calculation results.

2. The method according to claim 1, characterized in that The wavelength of the light field is visible light, X-ray, ultraviolet or infrared.

3. The method according to claim 1, characterized in that The optical element material is metal, inorganic material, organic material or composite material.

4. The method according to claim 1, characterized in that: The optical element has a rotationally symmetrical shape, and includes a zone plate, a Kinofom optical element, a Fresnel lens, a convex lens or a concave lens, or an optical element in which several of the elements are combined.

Images