Virtual source for Lommel–Gauss beams

We introduce a group of virtual sources for generating Lommel–Gauss beams based on beam superposition to analyze nonparaxial propagation. We typically derive the paraxial approximation and integral representations of the nth-order Lommel–Gauss beams. The first three orders of the nonparaxial corrections for the on-axis field of the Lommel–Gauss beams are analytically obtained. The on-axis intensity distribution of the corresponding nonparaxial corrections is also provided.


Introduction
Beam propagation is governed by the Helmholtz wave equation, which can be separately solved in orthogonal coordinate systems. Four special solutions of the wave equation, namely, the exponential, Bessel, Mathieu, and parabolic functions, can be obtained in the Cartesian, cylindrical, elliptical, and parabolic coordinates, respectively [1]. The function forms of the four special solutions are independent of the propagation distance. Accordingly, these forms correspond to four types of important nondiffracting beams: ideal plane (Cosh), Bessel [2], Mathieu [3], and parabolic [4] beams, respectively. Notably, none of the Cosh, Bessel, Mathieu, and parabolic functions are integrally square. Therefore, ideal nondiffracting beams with an infinite extent and energy are not physically achievable. In fact, the Cosh-Gauss [1,5], Bessel-Gauss [1,6], Mathieu-Gauss [1,7], and parabolic-Gauss [1] beams, which have finite energy, can be achieved experimentally. These beams can be regarded as quasi nondiffracting beams because they can propagate over an extensive range without significant diffraction [1]. Particularly, Bessel beams are the first nondiffracting beam with practical value amongst the four types of beams, and are extensively applied in laser manipulation [8], medical imaging [9], and other fields. The classical Bessel beam is a centrally symmetric concentric circle. Recently, Kovalev and Zhao [10,11] theoretically introduced and generated a new type of beam called the Lommel beam, which can be mathematically expressed in terms of the highorder Lommel functions with two variables. Lommel beams are a linear superposition of Bessel beams with identical axial projections of the wave vector; the transverse intensity distribution of the former has a reflective symmetry with respect to the Cartesian coordinate axes. Particularly, the intensity pattern across the beam's section can be tuned continuously through a simple adjustment of the beam parameters [11]. Furthermore, another advantage of Lommel beams lies in the continuous change of its orbital angular momentum, whereas Bessel beams exhibit a discrete change. These distinct properties will provide immense potential for optical trapping, the rotation of microparticles [11], and multiplexing communication [12]. In order to use Lommel beams in experiments, the propagation properties of the Lommel beams should be analyzed. Similar to the Bessel-Gauss beams, only Lommel-Gauss beams [12,13] with finite energy actually exist. Hence, this study focuses on the Lommel-Gauss beams.
By constructing virtual source points [5][6][7][14][15][16], researchers have studied the nonparaxial propagation characteristics of several kinds of nondiffracting beams with important application value, such as the Bessel-Gauss [6], Hermite-Gauss [15], and Laguerre-Gauss [16] beams. To research the propagation characteristics of the Lommel-Gauss beams, we introduce a group of virtual sources for generating the nth-order Lommel-Gauss beams, thereby establishing the inhomogeneous Helmholtz equation. Thereafter, we derive the paraxial approximation and nonparaxial integral representations of the Lommel-Gauss beams by using the Fourier-Bessel transform pair and Weber integral formula. Lastly, we obtain the first three orders of nonparaxial corrections for the on-axis field of the Lommel-Gauss beams. Since the research results provide a theoretical foundation for accurately revealing the propagation nature of Lommel-Gauss beams, they lay a theoretical foundation for the application of the Lommel-Gauss beams.

Theory
The Lommel beams can be described through a series of Bessel functions of orders n+2p (p=0, 1, 2, K) with cylindrical coordinates as follows [10,11,13,17]: ]is the propagation factor, k 2p l = / is the wave number of the monochromatic light with wavelength l, b is the beam's scaling factor, c is the asymmetry parameter, n is the integer parameter defining the order of the Lommel beams, and J x n p 2 The reason why the Lommel-Gauss beams expressed in equation (1) can be expanded into the summation of the infinite terms of the Bessel functions of different order numbers in the cylindrical coordinate system can be understood as the principle of the independence and superposition of beam propagation. In the initial plane of z 0 = , the Lommel-Gauss beams can be expressed as follows [13]: where 0 w is the beam waist width of the Gaussian envelope in z=0.
In the physical space z 0 > , constructing a scalar Lommel-Gauss beam that propagates along the z-axis is preferred. We can assume that E z , , where k is the wave number, and d (·) is the Dirac delta function. A differential equation for the radial spectrum of U z , n p 2 r + ( ) is determined from equation (4) using the Fourier-Bessel transform pair [6]: ) is the radial spectrum of U z , ) into equation (5), we obtain the following expression: ). Equation (7) is a complex expression. Hence, we have to deal with it approximately. Under the restriction of k 2 2 h  , z can be expressed as a series expansion of a small amount of . We retain the leading term for the amplitude factor and the first two terms for the phase factor in equation (7). In this approximation, equation (7) transforms as follows:  (2) and (7), we know that the scale length of the variation of r is 0 w and the significant range of the variation of h is from 0 to a quantity of the order of 1 0 w / . Given that the beam waist 0 w is generally larger than the wavelength k 2p l = / , the condition of k 2 2 h  is especially easy to satisfy.
According to the Weber integral formula [18,19] The parameters S ex , ex r , and z ex can be determined by comparing equation (10) with equation (2). Thus, we can readily obtain the following expressions: By substituting equations (12) to (14) into equation (10), the paraxial approximation to U z , n p ,2 r ( ) can be obtained. Therefore, the paraxial approximation to E z , , n P , r j ( )of the Lommel-Gauss beams is written as follows: The additional subscript P represents the paraxial approximation.
By substituting equations (12) to (14) into equation (7), we also obtain the precise integral expression for U z , Thereafter, the on-axis field for the nth-order Lommel-Gauss beams is obtained using equations (3) The distribution of the field of 0 r = should be analyzed because researchers in the majority of cases mainly focus on the on-axis intensity distribution when using various types of beams in scientific experimentation.
For equation (17), we perform the series expansion of in the powers of h. When k 2 2 h  , the product of both series' terms up to the order k m 0, 1, 2, ) is retained, and the mth-order nonparaxial corrections are obtained [7]. The product of both series' terms up to the order k 0 6 w -( ) is determined to obtain the first three nonparaxial corrections of the nth-order Lommel-Gauss beams on the axis for m=3. Equation (17) can be presented as follows: ) is a confluent hypergeometric function, and the function G is a gamma function [18]. Thus, we can derive the following equation from equations (18) The analytical expressions for the first-order and secondorder nonparaxial corrections and the analytical expression of the zero-order paraxial approximation can be obtained using the same method.
As an example, we use the nonparaxial corrections for the on-axis field of the Lommel-Gauss beams using the typical parameters 632.8 nm, 80 m 1 l b = = -, and 60 m 0 w m = to calculate the on-axis intensity distribution of the Lommel-Gauss beams for n=0. The calculation results are shown in figure 1.
In figure 1, the normalized intensity between the nonparaxial corrections and paraxial approximation are not