Many laser applications require specific irradiance distributions to ensure optimal performance. Geometric optical design methods based on numerical calculation of two plano-aspheric lenses have been thoroughly studied in the past. In this work, we present an alternative new design approach based on functional differential equations that allows direct calculation of the rotational symmetric lens profiles described by two-point Taylor polynomials. The formalism is used to design a Gaussian to flat-top irradiance beam shaping system but also to generate a more complex dark-hollow Gaussian (donut-like) irradiance distribution with zero intensity in the on-axis region. The presented ray tracing results confirm the high accuracy of both calculated solutions and emphasize the potential of this design approach for refractive beam shaping applications.
© 2014 Optical Society of America
Many commercial, medical and scientific applications of the laser have been developed since its invention. Materials processing applications, like cutting or welding; and medical applications, like corneal surgery are just a few examples. Many of these applications require a specific beam irradiance distribution to ensure an optimal performance . There are several optical technologies that can be used for laser beam shaping, including simple apertures, diffractive optical elements , lenses  and mirrors . The right choice depends on the optical design problem. Refractive laser beam shapers are commonly used due to their high optical efficiency and simple structure .
Often, it is possible to apply geometrical methods to shape a laser beam profile. This common design approach is based on the ray mapping between the input and the output plane which takes into account conservation of energy, a constant optical path length of the geometric wave-front between the input and output planes, and the ray trace equations. Geometric ray mapping-based designs with two plano-aspheric lenses have been thoroughly studied and improved [6–10]. Figure 1 shows a typical Galilean configuration (without internal focusing) of a refractive laser beam shaper consisting of two plano-aspheric lenses transforming a Gaussian input into a flattop output intensity distribution with a gradual roll-off from the uniform to null region. The indicated geometric ray mapping then describes where an initial parallel ray with radial distance r is refracted at the first aspherical surface to be redirected in parallel direction with radial distance R at the second aspherical surface.
Even though analytic expressions of ray mapping functions do exist for some specific cases, the surface profiles of the lenses are still calculated numerically. The main idea of most common geometrical optics design approaches is to describe the laser beam shaping design problem by mathematical equations and then use numerical techniques to solve the particular equation. Solutions have been found in the form of an integral equation, by Kreuzer, ordinary differential equation, by Rhodes and Shealy, and Monge-Ampère equation, by Oliker. Usually, a set of discrete data points is obtained using numerical calculation techniques and the final profiles of the lenses are constructed by either polynomial or spline fitting.
In this work, we present a novel different design approach that describes the entire optical design problem based on functional differential equations. This design approach is related to an analytic design method for optimal imaging that enables the coupling of three ray-bundles with only two free-form lens surfaces in two and three dimensions [11, 12]. However, the original design approach worked so far only for spherical and plane wave-fronts. Therefore, further extensions are necessary and key to properly incorporate the ray mapping function.
Starting from the example of a Gaussian to flat-top irradiance beam shaping system in Fig. 1, a set of two functional differential equations are derived in Sec. 2. The mathematical problem is then converted into a set of linear algebraic equations that allows calculating analytic symbolic expressions for the lens profiles f and g described by two-point Taylor series expansions. To further demonstrate the versatility of this new approach, a second design problem is solved in Sec. 3 to generate a more complex dark-hollow Gaussian (donut-like) irradiance profile with zero intensity in the on-axis region. Here, the mapping function is given by a fractional Taylor series which has direct impact on the surface representation as well. The presented ray tracing evaluation confirm the high accuracy of both calculated solutions and emphasize the potential of this design approach for refractive beam shaping applications.
2. Analytic solution to convert a Gaussian into a flat-top irradiance distribution
If the lenses are rotationally symmetric and described by surfaces of revolution, this means that the problem can be examined in a plane. The input irradiance distribution is given by the Gaussian function13], the output intensity is chosen as a flattened Lorentzian (FL) distribution Fig. 1. Where q denotes the beam shape parameter and RFL is the beam width parameter of the flattened Lorentzian . The input encircled energy is calculated as a function of r as Fig. 2(a) for ε = 1, w0 = 2.366mm, RFL = 3.25mm and q = 16.
Now that the output coordinate value R can be calculated for every input coordinate r, it is possible to derive an integral representation for the sag of each aspherical lens that can be solved numerically following the work of Kreuzer, Shealy and Hoffnagle .
Instead of following this well-known design approach, an alternative solution strategy will be explained hereafter. As both incoming and outgoing geometric wave-fronts are plane (i.e. all rays are parallel), the two plane surfaces in Fig. 1 have no optical power and can be neglected by assuming that the in- and output rays are immersed in a medium with refractive index n2. This is shown in Fig. 2(b). Consider an arbitrary light ray path. The ray is emitted from a point v⃗0 on the plane incident wave-front and intersects the first surface at a point p⃗ = (r, f(r)), where it is refracted towards (R(r), g(R(r))) on the second surface and refracted again towards the plane outgoing wave-front (parallel to the z-axis). The path of the light rays are governed by Fermat’s principle which states that the optical path length (OPL) between two points is an extremum along the light ray containing these points. The OPL between points v⃗0, p⃗ and (R(r), g(R(r))) is given byEqs. (7) and (8) form a pair of functional differential equations for two unknown functions f(r) and g(R(r)).
Suppose (f, g) is an analytic and smooth solution to the functional differential equations (7) and (8), Taylor’s theorem implies that the functions must be infinitely differentiable and have a power-series representation. In our previous work , a Taylor series approach has been successfully used to solve a related set of functional differential equations by solving equationsFig. 2(a). This circumstance is directly linked to the fact that the evaluation of Eqs. (9) contains derivatives of the mapping function up to nth order, evaluated at that point r1. If the mapping function is now reconstructed from its (exact) derivatives at r1, the Taylor series of this mapping function would not converge for all values r (see Fig. 2(a)) and the same applies for the lens profiles f and g. Therefore, as alternative, a Taylor series expansion about two points r1 and r2 can be defined as 14]. The two-point Taylor series approximation of the mapping function converges well for all values r and is added to the graph in Fig. 2(a). The selection of points r1 and r2 depends strongly on the shape of the considered mapping function. To find a suitable pair of points, it has been found sufficient to verify that the two-point Taylor series approximation of the mapping function converges well for all values r. This was always the case even for manually selected r1 values close to the optical axis and r2 values close to the maximum radius. There is no strict mathematical proof that the same convergence also applies for the lens profiles. However, all evaluated designs demonstrated this behavior so far. The results presented in this work will provide strong evidence that the solutions for the lens profiles are analytic and smooth as well. Therefore, the surface profile functions can be expressed by two-point Taylor series expansions Fig. 2(b). For the first ray, z1 can be freely chosen without loss of generality, whereas the value z3 fixes the distance between the two lens surfaces. For the second ray, the position z2 can not be determined analytically and represents a single parameter that links the two local solutions at r1 and r2. This fact will be explained in more detail at the beginning of Sec. 2.1. The position z4 then follows from a constant OPL. All four slope values m1, m2, m3 and m4 at points ri (i = 1..4) follow from the ray mapping relationship and Snell’s law.
The initial conditionsFig. 2(b) then satisfy the conditional equations Di = 0 for i = 1, 2 and provide general solutions for the initial two-point Taylor series coefficients depending upon the single variable z2. In ascending order, it is now possible to calculate (2k + 1)th order Taylor polynomials of f(r) and g(r) by calculating the derivatives and evaluating equations Eqs. (13) solve Eqs. (14).
For every higher order n ≥ 1, the set of equations (14) results in four linear algebraic equations for two-point Taylor series coefficients pi, qi, ui and vi. By sorting and combining all terms, the equations can be expressed as a compact matrix equation for particular two-point Taylor series coefficients for arbitrary n ≥ 1. The elements of the linear algebraic equations consist of mathematical expressions only dependent on previously calculated Taylor series coefficients. Therefore, all two-point Taylor series coefficients can be calculated as exact symbolic expressions in ascending order - step by step. So far, no approximations have been made. The general solution scheme for n ≥ 1 allows to calculate successive Taylor series coefficients of f(r) and g(r) up to a very high order. The inevitable truncation of the infinite sum of terms will be the only approximation made. The radii of convergence for the Taylor expansions f(r) and g(r) are very important, as they indicate the maximum aperture that can be achieved for a given set of initial values. For the examples considered in this work, the radii of convergence are larger than the radial coordinate 0 < r < rmax of the lenses.
2.1. Evaluation of the analytic solution
The parameters of the in- and output intensity distributions are w0 = 2.366mm, RFL = 3.25mm and q = 16 as before. The two points of the Taylor series are r1 = 0.2mm and r2 = 3.6mm. The chosen distance between z1 = 0mm and z3 = 16mm defines the lens separation. For a chosen refractive index n2 = 1.5, this initial ray path construction determines the OPL. After applying the constant OPL condition for point r2, there still exists an infinite number of solutions that solve the design problem locally in the neighborhood of r2. In order to link the two local Taylor expansion solutions at r1 and r2, the correct relative positioning in z-direction of the two solution parts has to be identified. It is very important to note that both local solutions at the respective points are exact and calculated symbolically, only the relative position of the two parts of the overall solution has to be adjusted numerically. To find the appropriate value of z2, the implemented calculation of the Taylor coefficients for the lens profiles is defined as a function in Mathematica. After calculation of the Taylor series coefficients up to a certain order k, 10 equidistant rays between 0 mm and 4 mm (corresponds to ≈ 99.7% encircled input energy) are traced and evaluated. The spatial deviation of each ray is then defined by the mapping error R′i − R(ri). Where R′ i is the calculated intersection of the ith ray with the surface profile g and R(ri) is the ideal position defined by the mapping function. As a second figure of merit, the variation in optical path length is also calculated for all rays. As one possible overall merit function, m(z2) is given by the sum of the two root mean square values of the spatial deviation and variation in OPL of the 10 rays over the full entrance aperture.
Figure 3 shows the evaluation of the merit function m(z2) as a function of variable z2 for values between −0.5 and 0 mm and for a step size of 0.0025 mm. A logarithmic plot is chosen to better cover the wide range of values.
The curve in parameter space is characterized by a convex curvature which makes it easy and fast to find the optimal relative position of z2 using a local optimization in Mathematica.
To quantitatively evaluate the accuracy of the derived solution, the spatial and angular deviation is calculated for 200 rays between 0 and 4 mm along r-direction. This evaluation is performed for 7th (k=3), 13th (k=6) and 19th (k=9) order Taylor polynomials for f and g, respectively. Figure 4(a) shows the variation in optical path length (OPL) and Fig. 4(b) shows the spatial mapping error for 3 calculated solutions of different order. The root mean square values of the variation in OPL are 1.647μm (k=3), 0.617μm (k=6) and 0.058μm (k=9). The root mean square values of the mapping error are 0.559 × 10−2 (k=3), 0.061 × 10−2 (k=6) and 0.009 × 10−2 (k=9), respectively. For both graphs, the variation in OPL and mapping error reduce very well with an increasing order underlining the convergence and achievable high accuracy of the derived solution.
Furthermore, the maximum collimation error after refraction at the second plane surface is evaluated for the 3 different solutions. The maximum collimation error is 2.167mrad (k=3), 1.124mrad (k=6) and 0.265mrad (k=9). The achieved collimation quality for k=9 is clearly better than the maximum collimation error of 0.8mrad in . Further improvements, if necessary, can be achieved by increasing the order of the Taylor polynomials and/or by adjusting the positions of points r1 and r2, respectively.
3. Analytic solution for Gaussian to dark-hollow Gaussian refractive beam shapers
The existence of a closed-form expression for the mapping function was important for the design approach in Sec. 2 to work. The following example will demonstrate that this is not a necessary condition. The optical design problem of a refractive beam shaper that transforms a Gaussian into a dark-hollow Gaussian (DHG) intensity profile with zero intensity at the optical axis does not allow to calculate a closed-form expression for the mapping function, as will be shown in the following. Focused small ring patterns are of interest for some industrial applications such as laser cladding technology. Figure 5 shows the cross-section of the correspondent setup together with the in- and output intensity profiles.
The input irradiance distribution is again given by the Gaussian function as defined in Eq. (1). The output intensity distribution is given byEq. (15) reduces to the fundamental beam with beam waist w1 .
As done in the previous Section, it is possible to determine the relationship between output radial distance R and input coordinate r such that the encircled energy B(R) at the output is equal to the encircled energy A(r) at the input. The input encircled energy was given in Eq. (3) and the output encircled energy isFig. 6(a); for parameters w0 = 2.366mm, w1 = 2mm and n = 1.
Interestingly, the inverse mapping function R(r) is then described by a Taylor series expansionFig. 6(a)). When trying to use the same solution scheme for the approximated mapping function R(r) as in the previous section, the problem appears that even for a two-point Taylor series representation of the surfaces the profiles will still not converge for all values of r. This problem is directly linked to the fractional exponents of the mapping function R(r), which will become more clear later.
To solve this design problem, it is better to solve its inverse problem: where a dark-hollow Gaussian is transformed into a Gaussian beam. This inverse optical design problem is identical due to the symmetry of the optical design where both in- and outgoing geometric wave-fronts are plane. The inverted design is shown in Fig. 6(b), where the reversed directions of the rays are indicated. As done previously, this time from right to left, the OPL between points v⃗0, p⃗ = (R, g(R)) and (r(R), f(r(R))) is given byEqs. (20) and (21) at points R1 and R2 in ascending order Eq. (21) follows that it is also possible to calculate the second surface profile f(r(R)) from the integral equation Eq. (18) as argument in the first derivative g′(R), the unknown function f(R) can be exactly calculated via symbolic integration. The function f(R) is then described by a polynomial function
3.1. Evaluation of the analytic solution
The performance evaluation of the analytic solution for Gaussian to dark-hollow Gaussian (DHG) beam shaping is carried out analogously to Sec. 2.1. The parameters of the in- and output intensity distributions are w0 = 2.366mm, w1 = 2mm and n = 1. The two manually selected points of the Taylor series are R1 = 0.3mm and R2 = 4.2mm. By fixing the distance between Z1 = 0mm and Z2 = 16mm and the refractive index n2 = 1.5, the remaining parameter Z2 that adjusts the relative positioning in z-direction is optimized for an identical merit function as defined in Sec. 2.1. Figure 7(a) shows the variation in OPL and Fig. 7(b) shows the spatial mapping error for 3 different solutions of different order. All calculations have been repeated for 200 rays between 0 and 4 mm along r-direction (which corresponds to ≈ 99.7% encircled input energy). The evaluation is performed for 7th (k=3), 11th (k=5) and 15th (k=7) order Taylor polynomials for g, and f has been truncated to 21th order.
The root mean square values of the variation in OPL are 0.053μm (k=3), 0.064μm (k=5) and 0.025μm (k=7). The root mean square values of the mapping error are 0.792 × 10−2 (k=3), 0.776 × 10−2 (k=5) and 0.164 × 10−2 (k=7), respectively. As in the previous design example, the variation in OPL and the mapping error reduce very well with an increasing Taylor series order, confirming the convergence and high accuracy of the derived analytic solution. Furthermore, the maximum collimation error after refraction at the second plane surface is evaluated for the 3 different solutions. The maximum collimation error is 0.744mrad (k=3), 0.063mrad (k=5) and 0.058mrad (k=7). If required, an even higher accuracy can be achieved by increasing the order of the Taylor polynomials and/or by fine-tuning the radial positions R1 and R2.
Within the scope of this work, it has been shown how the optical design of refractive laser beam shaping systems based on two plano-aspheric lenses can be described by means of a new mathematical formulation using functional differential equations. The derived solution scheme has been successfully verified for two specific design tasks and allows an increasingly accurate calculation of the lens profiles described by high order Taylor polynomials. In case of the flattened Lorentzian beam shaping system, a maximum collimation error value of 0.265mrad has been achieved. If required, a higher accuracy seems feasible with this new design approach.
In case of the dark-hollow Gaussian output intensity distribution, a maximum collimation error value of 0.058mrad has been achieved. Furthermore, the design method established a direct link between the mathematical structure of the mapping function and the mathematical function describing the first optical surface (fractional exponents). This deeper insight provided by the optical design approach could prove beneficial when it comes to designing and manufacturing of such optical systems. The description of the lens surfaces with analytic functions allows imposing boundary conditions (e.g., for the local curvature of the surfaces) directly in the design process to meet manufacturing requirements.
Particularly noteworthy is the successful implementation of a two-point Taylor series expansion that could be also of interest for other optical designs, for example in imaging systems . Here, a multi-point Taylor series approach could be used to generalize the optical design for more than two surfaces that allows the simultaneous coupling of more than three ray-bundles.
A related functional differential equation based optical design approach has been already applied successfully to design a system with movable free-form optics for a new class of solar concentrators . This already demonstrated possibility to incorporate movement in the mathematical description could be useful for the design of dynamic refractive laser beam shaping systems. In such a system, the relative motion of the lenses will be used to dynamically change the shape or profile of the laser beam: a promising functionality for the laser materials processing industry. This dynamic optical design concept will be addressed in future work.
The work reported in this paper is supported in part by the Research Foundation - Flanders (FWO-Vlaanderen) that provides a post-doctoral grant to Fabian Duerr, and in part by the IAP-BELSPO grant IAP P7-35 photonics@be, the Industrial Research Funding (IOF), Methusalem, and the OZR of the Vrije Universiteit Brussel.
References and links
1. F. M. Dickey, S. C. Holswade, and D. L. Shealy, Laser Beam Shaping Applications (CRC Press, 2005). [CrossRef]
3. J. A. Hoffnagle and C. M. Jefferson, “Design and performance of a refractive optical system that converts a Gaussian to a flattop beam,” Appl. Opt. 39, 5488–5499 (2000). [CrossRef]
4. V. Oliker, “Optical design of freeform two-mirror beam-shaping systems,” J. Opt. Soc. Am. A 24, 3741–3752 (2007). [CrossRef]
5. L. Romero and F. Dickey, “Lossless laser beam shaping,” J. Opt. Soc. Am. A 13, 751–760 (1996). [CrossRef]
7. J. A. Hoffnagle and C. M. Jefferson, “Beam shaping with a plano-aspheric lens pair,” Opt. Eng. 42, 3090–3099 (2003). [CrossRef]
8. H. Ma, Z. Liu, P. Jiang, X. Xu, and S. Du, “Improvement of Galilean refractive beam shaping system for accurately generating near-diffraction-limited flattop beam with arbitrary beam size,” Opt. Express 19, 13105–13117 (2011). [CrossRef] [PubMed]
9. J. Rubinstein and G. Wolansky, “Intensity control with a free-form lens,” J. Opt. Soc. Am. A 24, 463–469 (2007). [CrossRef]
10. Z. Feng, L. Huang, G. Jin, and M. Gong, “Designing double freeform optical surfaces for controlling both irradiance and wavefront,” Opt. Express 21, 28693–28701 (2013). [CrossRef]
11. F. Duerr, P. Benítez, J. C. Miñano, Y. Meuret, and H. Thienpont, “Analytic design method for optimal imaging: coupling three ray sets using two free-form lens profiles,” Opt. Express 20, 5576–5585 (2012). [CrossRef] [PubMed]
12. F. Duerr, P. Benítez, J. C. Miñano, Y. Meuret, and H. Thienpont, “Analytic free-form lens design in 3D: coupling three ray sets using two lens surfaces,” Opt. Express 20, 10839–10846 (2012). [CrossRef] [PubMed]
13. D. L. Shealy and J. A. Hoffnagle, “Review: design and analysis of plano-aspheric laser beam shapers,” Proc. SPIE 8490, 849003 (2012). [CrossRef]
14. R. H. Estes and E. R. Lancaster, “Two-point Taylor series expansions,” NASA TMX-55759 (N67-23965) (1966).
16. F. Duerr, Y. Meuret, and H. Thienpont, “Potential benefits of free-form optics in on-axis imaging applications with high aspect ratio,” Opt. Express 21, 31072–31081 (2013). [CrossRef]