1. Introduction

In many applications an optical system is required to redirect the light from a source and shape its output irradiance distribution so that a prescribed field is generated on a given target. Systems with such capabilities are needed in numerous areas of modern science and technology. Laser technology is one such area in which the beam generated by a laser has to be transformed so that some given region of space is illuminated with an irradiance pattern different from that of the source [1,2]. Similar tasks must be solved in nanotechnology, photonics, photolithography, directed-energy, material processing, astronomical instruments, sensing systems, solar collection, optical engineering metrology, buildings and roads illumination, and many others.

A systematic approach to beam shaping (and to many other optical design problems) is to assume a priori that the solution and the data are rotationally symmetric (RS) [1–6]. In fact, most of the currently used conventional optical devices rely on lenses and mirrors which are rotationally or rectangularly symmetric or can be described by explicit formulas. The principal drawback of this approach is that the inherent limitation of assumed special symmetry excludes or makes highly inefficient utilization of so designed optics in applications requiring control of irradiance produced by nonsymmetric sources and distributed over areas of irregular shape. By contrast, a freeform optical surface is a surface designed without any a priori assumptions of available symmetry or explicit formula. Consequently, freeform surfaces provide significantly more degrees of freedom to create optical systems for a wide variety of applications. In addition, such systems can be more compact, energy efficient and light. The goal of this paper is to present a systematic and mathematically rigorous approach to design of a pair of plano-freeform optical surfaces for beam shaping without any a priori symmetry assumptions. At the same time it will be shown that our approach is engineeringly viable and practical.

A typical configuration of a beam shaping optical system is shown schematically on Fig. 1. The spatial coordinates are (x, z) with x = (x₁, x₂) ∈ α, where α is the plane z = 0. We work in the geometrical optics (GO) framework and use the associated notions of wavefronts, rays, radiances and irradiances. In analytic formulation of the beam shaping problem (BSP) the lens R₁ is represented as a graph of some function z(x), x ∈ Ω̄, where Ω̄ is the cross section of the input beam. In practice, a lens has the second side, say, to the left of R₁. We assume that that side is inactive, that is, it does not change the directions of incoming light rays; this comment applies also to the lens R₂ with obvious modifications. Everywhere in the paper Ω denotes a bounded open domain on α and Ω̃ its closure. The same convention is used to denote the target set T̄_d in the plane z = d. Points of the set T̄_d to be illuminated by the output beam produced by the lens R₂ are denoted by (p, d), where p = (p₁, p₂) ∈ T̄ ⊆ α and T̄ is the projection of T̄_d on α. This two-lens system defines a ray mapping P_d : Ω̄ → T̄_d; P_d(x) = (p = P(x), d). The graph of the function w(p), p ∈ T̄, is the lens R₂. The radiance I₁ of the input beam and the irradiance I₂ of the output beam are related by Eq. (1) below. The design task is to determine the functions z and w realizing the required optical transformation of rays and radiances.

 

Fig. 1 A system of two plano-freeform lenses for beam shaping; only the active (non-planar) sides of the lenses are shown. The two lenses R₁ and R₂ are required to transform a collimated input beam with an arbitrary radiance I₁(x) into a collimated output beam with a prescribed irradiance I₂(p) at a given region on the plane α′ : z = d; n = const is the normalized refractive index of R₁ and R₂ (as explained at the beginning of Subsection 2.1). From Ref. [7].

Download Full Size | PDF

Using the GO setting and applying the Snell and energy conservation laws, the relation between the radiance I₁ and irradiance I₂ the plane z = d can be expressed as

The PDE for design of beam shaping optics with two refractive lenses is not among such special cases. In such circumstances the usual approach to a given design task is to use an ad hoc heuristic technique in each specific case [11,12].

In this paper we avoid direct use of PDE methods. Instead, our approach is based on the mathematically rigorous geometric supporting quadric method (SQM) [9] combined with the variational (optimal transport) theory [7].

It is worthwhile noting that freeform optics can be useful even in the RS case. For example, in [1] the authors designed a two-lens beam shaper under assumption that the input and output radiance patterns are RS. Their approach requires an initial circularization and filtering of the input beam from the laser before it can be accepted by the beam shaper; see Fig. 2. Clearly, a solution unconstrained by RS requirement will provide a useful improvement of the process by eliminating the step (a)→(b) and possibly (b)→(c); see Fig. 2.

 

Fig. 2 Design of a two-lens system for transforming a non-uniform laser beam into a uniform intensity RS beam requires first circularization (a)→(b) and filtering (b)→(c) of the beam from the laser. Then a RS beam shaper can convert it into a “flat-top” (c)→(d); from Ref. [13] with permission.

Download Full Size | PDF

2. Variational formulation and solution of two-lens beam shaping problem (BSP)

2.1. Statement of BSP

The geometry of the BSP is shown in Fig. 1 and below we refer to this figure. It is assumed that (after a normalization, if necessary) the media to the left of R₁ and to the right of R₂ have constant refractive index n, while the media between R₁ and R₂ has refractive index 1. In other words, n denotes the relative refractive index. Because collimated beams propagate longer distances without dispersion and distortion [1] it is physically natural to assume (and we do so) that the input beam is collimated and propagates in direction k. The output beam is also required to be collimated, propagating in the same direction, and of constant phase. It is assumed that the lenses R₁ and R₂ as well as the media are lossless, that is, no energy is lost during the propagation of light. An individual light ray of the input beam propagating in the positive z-direction k is tagged by a point x ∈ Ω̄. The light rays of the input beam are refracted by the lens R₁ and produce a beam (not collimated!) which is intercepted and refracted by R₂, producing the collimated output beam of direction k. The output beam, after lens R₂, in addition to being collimated, must produce a prescribed irradiance distribution μ on a given target domain T̄_d, which is in the plane z = d parallel to α and located at the distance d > 0 from α. Thus, T̄_d is the cross section of the output beam by the plane z = d.

Similar to [1] we consider here the case when the refractive index n > 1. The case 0 < n < 1 is similar but requires additional considerations. We plan to treat this case in a separate paper. In [1] the authors, in addition to the refractive index n, use the distance between the lenses as the given design parameter. In the radially symmetric case such distance is taken as the distance along the common axes of revolution for the two lenses. When radial symmetry is not assumed this parameter is not well defined. Instead, when the lenses are freeform surfaces, we use the parameter β = l − nd, where l = nz(x) + t(x) + n(d − w(P(x)) is the optical path length (OPL) [14]. Thus, the refractive index n, the distance d in Fig. 1 and OPL l are constant design parameters. It is known that β = 0 is physically inconsistent with n ≠ 1 [7]. Physical considerations show that only the case β < 0 is of interest and this is assumed in this paper.

It is convenient to recall here the expression for the ray mapping in case when the function z is differentiable. Decompose the ray mapping P_d : Ω̄ → T̄_d as P_d(x) = P(x) + dk, where P : Ω̄ → T̄ is the orthogonal projection of P_d on the plane α. We refer to P as the refractor map. Put p = P(x) and $M(x)=\sqrt{1+(1-n^2){\left|\nabla z(x)\right|}^2}$, (∇ ≡ grad). It was shown in [7,15] that

It is important to note that the factor M(x) in the denominator of Eq. (2) precludes application of the Monge principle as in [17]. Replacing M(x) ≈ 1 corresponds to paraxial approximation.

2.2. The supporting quadric method (SQM)

In recent years the SQM has often been applied to design freeform mirrors, lenses and diffractive optical elements [18–23]. It provides a viable and rigorous alternative to the approach requiring solution of the complex nonlinear PDE. With the SQM combined with optimal transport theory the solution of the BSP is found in the class of functions called weak solutions (to be defined below) which are continuos, convex but not necessarily two times differentiable everywhere. Consequently, the SQM requires re-formulation of the BSP in geometric terms not involving explicit use of Jacobians. The basic principles of the SQM have been described before [9,22,24]. However, its particular implementation for the BSP requiring two refracting surfaces is presented here for the first time. In order to make our presentation reasonably self-contained we need to recall the geometric/optical properties of two-sheeted hyperboloids. Eventually, the required lenses are obtained as envelopes of families of sheets of such hyperboloids.

The notations here are as in Fig. 1. Put also

 

Fig. 3 Hyperbola H^l maps the segment [a, b] into F^r.

Download Full Size | PDF

Similarly, for an arbitrary fixed (x, ξ) ∈ α × ℝ the graph of the function

The freeform lenses for beam shaping designed here are envelopes of families of suitable sheets of hyperboloids of the type (5) and (7). Such envelopes are defined as follows.

Definition 2.2 Let Ω and T be two bounded domains on the plane α and Ω̄, T̄ their closures. Denote by C(Ω̄) and C(T̄) the set of continuous functions on Ω̄ and T̄, respectively. Let n > 1 and β < 0. A pair (z, w) ∈ C(Ω̄) × C(T̄) is called a two-lens (TL) system of type A (or Keplerian) if

An illustration of Eqs. (9)–(10) in 2D and T̄ consisting of five points p₁, ..., p₅ is in Fig. 4, where for brevity we put F_pi = H(x, p_i) + w(p_i), i = 1, ..., 5. The points Q_i = (p_i, w(p_i)) are their second foci lying on the lens R₂. Because β < 0, by Lemma 2.1 for a fixed p ∈ T̄ the hyperboloid H(x, p) + w(p) transforms the plane front into spherical front with center at (p, w(p)). According to Eq. (9), the lens R₁ defined as the graph of z(x), x ∈ Ω̄, consists of pieces of hyperboloids F_pi,fi, i = 1, ..., K, such that a ray in direction k incident on R₁ is intercepted by the “nearest” F_pi,fi and refracted so that it passes through the focus Q_i. Thus, in this case the lens R₁ is an “envelope” defined by Eq. (9) of a family of suitable hyperboloids and the lens R₂ is an envelope defined by Eq. (10).

 

Fig. 4 The function z in Eq. (9) defines the first lens R₁. The second lens R₂ is defined by Eq. (10).

Download Full Size | PDF

We call the branches H of hyperboloids in the definition (2.2) “supporting” hyperboloids. Obviously, at each point of a TL system of type A and B there exists at least one supporting hyperboloid. In Fig. 4 the supporting hyperboloids of R₂ are not shown.

Definition 2.3 For a TL system in Λ(Ω̄, T̄) the refractor map P_z,w : Ω̄ → T̄, possibly multivalued, and its inverse $P_{z,w}^{-1}$ : T̄ → Ω̄ are defined by equations

We now formulate precisely the weak form of the beam shaping problem. First we define, following [7], the light transfer from Ω̄ to T̄ for TL systems. Let I₁ be a nonnegative integrable function on α with spt(I₁) ⊆ Ω̄, where $\text{spt}(I_1):=\overline{\left\{x\in \overline{\mathrm{\Omega }}\left|I_1(x)\ne 0\right|\right\}}$. Physically, I₁ is the intensity of the radiance distribution of the input beam. Denote by ℬ(T̄) the set of subsets of T̄ measurable with respect to the usual area measure on α. Let (z, w) ∈ Λ(Ω̄, T̄) and P_z,w the refractor map defined by (z, w). Then for any set τ ∈ ℬ(T̄) the set $P_{z,w}^{-1}(\tau )$ is measurable and the function

2.3. Optimal transport and discrete minimization problem

Our next task is to describe an algorithm for solving the discretization of Eq. (17). We will concentrate on the case when (z, w) ∈ Λ_A(Ω̄, T̄) since both lenses in this case are convex. First we note that solution of Eq. (17) is equivalent to a solution of an optimal mass transportation problem [7]. In order to outline the connection, recall the formulas (9)–(10) and define the following class of admissible pairs of functions

Theorem 2.5 Suppose n > 1 and β < 0 and I₁ and μ as above and satisfy Eq. (18). Put

Furthermore, this pair (z_min, w_min) is a weak solution of Eq. (17) and conversely any weak solution of Eq. (17) is a minimizing pair of Eq. (21).

Remark 2.6 It is important to note that in Theorem 2.5 it is assumed that Adm_A(Ω̄, T̄) ≠ ∅. This imposes a constraint on separation of Ω̄ and T̄ in the direction perpendicular to the z axis. For example, if Ω̄ = T̄ then this condition is satisfied. Alternatively, geometric arguments show that condition Adm_A(Ω̄, T̄) ≠ ∅ can be achieved by adjusting the OPL and the distance d in a way similar to the way it was done in [27]. In practice this property can be verified by a trial run of the code on the given data.

In view of 2.5, to solve Eq. (17) we need to solve the minimization problem (20)–(21). For that we use a computational algorithm based on linear programming theory [28]. The discrete version of Eqs. (20)–(21) is stated and solved as follows.

Step 1. Discretization. Let M, N ∈ ℕ, M ≥ 1, N ≥ 2. Subdivide Ω̄ into M measurable subsets ${\omega }_i^M\subset \overline{\mathrm{\Omega }}$, i = 1, ..., M, and T̄ into N measurable subsets ${\tau }_j^N\subset \overline{T}$, j = 1, ..., N, such that for each i = 1, ..., M diameters ${\omega }_i^M<1/M$ and for each j = 1, ..., N the diameters ${\omega }_j^N<1/N$, and for 1 ≤ t, k ≤ M, 1 ≤ h, l ≤ N,

Step 3. Refining the subdivisions of Ω̄ and T̄ one obtains a sequence of solutions of discrete problems guaranteed to converge to a weak solution of Eqs. (20), (21).

Remark 2.7 It was proved in [7], Lemma 4, that the slopes of designed surfaces do not exceed $1/\sqrt{n^2-1}$ . We also note that the above theory guarantees collimation of the shaped beam for lenses designed exactly. In practice, possible collimation errors depend on the numerical accuracy of the design, fabrication and propagation distance. Some related results are contained in [15]. We hope to further address these issues in a separate publication.

3. Design examples

As a first example, we consider the calculation of two plano-freeform lenses transforming the circular incident beam with a uniform radiance I₁(x) = 1/(πR²), $x\in \mathrm{\Omega }=\left\{(x_1,x_2)|x_1^2+x_2^2\le R^2\right\}$ into a rectangular beam with a uniform radiance I₂ (p) = 1/(w₁w₂), p ∈ T_d = {(p₁, p₂) | |p₁| ≤ w₁/2, |p₂| ≤ w₂/2}. The propagation directions of the incident and of the outgoing beams coincide with the direction of the z axis. The lenses were designed for the following parameters: radius of the incident beam R = 10 mm, the side lengths of the outgoing rectangular beam w₁ = 50 mm and w₂ = 25 mm, the distance between the tips of the lenses along the z axis t(0) = w(0) − z (0) = 30 mm, and the refractive index of the lens material n = 1.5. Without loss of generality, we assume z(0) = 0.

It follows from the results in Subsections 2.2 and 2.3 that the surfaces of the plano-freeform lenses z(x) and w(p) are the solution of a linear programming problem (22)–(23) for the Keplerian configuration. For the Galilean configuration one needs to solve the linear programming problem obtained from Eqs. (22)–(23) by reversing the inequality in Eq. (22) and maximizing the functional in Eq. (23). In both cases, the function z(x), x ∈ Ω was defined on a regular grid with M = 10000 points ${\overline{z}}_i^M$, i = 1, . . ., M. The function w(p), p ∈ T_d was defined on a 100 × 100 rectangular grid with N = 10000 points $w_j^N$, j = 1, . . ., N. The corresponding linear programming problems were solved using an in-house implementation of the interior point algorithm taking into account the special sparse form of the constraints matrix given by Eq. (22). The calculation time on a modern PC was about two hours. The calculated two-lens systems are shown in Fig. 5. In the Keplerian configuration, both lenses are convex [Figs. 5(a), 5(b)], whereas for the Galilean configuration [Figs. 5(c), 5(d)], only the second lens is convex [Figs. 5(c), 5(d)].

 

Fig. 5 The designed systems of two plano-freeform lenses in a Keplerian (a, b) and a Galilean (c, d) configurations to transform an incident circular beam into a rectangular outgoing beam. The dashed lines show the z axis.

Download Full Size | PDF

In order to estimate the performance of the designed two-lens systems, we simulated them in the commercial ray-tracing software TracePro in the following way [31]. First, the calculated optical surfaces represented by the sets of points ${\overline{z}}_i^M$, i = 1, . . ., M and $w_j^N$, j = 1, . . ., N were exported to a computer-aided design software Rhinoceros [32], in which 3D-models of the two plano-freeform lens systems were created (Fig. 5). Then, these 3D-models were used in TracePro for the simulations. Figures 6(a) and 6(b) show the calculated normalized radiance distributions generated by the Keplerian lens system in the planes z = 10 cm and z = 100 cm. The presented results were obtained by tracing 1 000 000 rays and demonstrate the formation of a rectangle with the required size.

The fact that the size of the rectangle and the radiance distribution remain nearly unchanged at different distances from the element confirms that the wavefront of the beam generated by the lenses is nearly plane. Indeed, the calculated root-mean-square deviation of the optical path length from a constant value at z = 10 cm does not exceed 1 nm. The normalized root-mean-square deviations (NRMSD) of the numerically simulated radiance distributions from a constant value amount to 12.1% and 13.8% at z = 10 cm and z = 100 cm, respectively.

Note that the radiance smoothly decreases toward the boundaries of the rectangle [see the cross-sections in Figs. 6(a) and 6(b)]. This decrease is due to Fresnel losses, which were neglected when designing the lens surfaces. Indeed, in the Keplerian configuration the deflection angles of the rays upon refraction by the freeform surfaces of the lenses vary from zero for the central ray going to the center of the generated rectangle to 35.8° for the rays going to the corners of the rectangle. Since the Fresnel losses increase with an increase in the deflection angle [33], the generated radiance smoothly decreases toward the boundaries of the rectangle in Figs. 6(a) and 6(b) Because of such large deflection angles, the total Fresnel losses on the freeform surfaces of the lenses amount to 24.6% of the incident beam flux.

 

Fig. 6 Normalized radiance distributions generated in the planes z = 10 cm (a, c) and z = 100 cm (b, d) by the designed Keplerian (a, b) and Galilean (c, d) lens systems calculated in TracePro. The radiance cross-sections along the coordinate axes are shown at the top and at the right of the figures.

Download Full Size | PDF

Figures 6(c) and 6(d) show the calculated normalized radiance distributions generated by the Galilean lens system at the same distances z = 10 cm and z = 100 cm. The radiance distributions generated in this case are near-uniform and have better quality comparing to the previous case. The NRMSD of the numerically simulated radiance distributions from a constant value amount to 3.9% and 4.1% at z = 10 cm and z = 100 cm, respectively. The better uniformity of the generated rectangular radiance distribution can be explained by decreased deflection angles leading to lower Fresnel losses. Indeed, in the Galilean configuration, the deflection angles of the rays do not exceed 26°. In this case, the total Fresnel losses on the freeform surfaces of the lenses amount to 10.1% of the incident beam flux.

As a second, more sophisticated example, we considered the design of a two-lens system transforming an incident circular beam with uniform radiance into a star-shaped beam with uniform radiance. The vertices of the star (of the region T_d) lie on a circle with radius R_st = 24.2 mm. We assume the distance between the tips of the lenses along the z axis t(0) = w(0) − z (0) = 50 mm, with the rest geometric and simulation parameters (R, n, M, N) remaining the same as in the previous case. In this example, we consider only the Galilean configuration since it provides a higher quality of the generated distribution. The plano-freeform lenses calculated by solving a linear programming problem Eqs. (22)–(23) at M = N = 10000 are shown in Fig. 7. Figure 8 shows the calculated normalized radiance distributions generated by the two-lens system at z = 10 cm and z = 100 cm. The simulation results confirm the formation of a star-shaped beam with the required size. The NRMSD of the calculated distributions from a constant value amount to 5.1% and 8.4% at z = 10 cm and z = 100 cm, respectively. As in the previous example, the estimated root-mean-square deviation of the optical path length from a constant value at z = 10 cm does not exceed 1 nm.

 

Fig. 7 (a) The calculated system of two plano-freeform lenses transforming an incident circular beam into a star-shaped beam. The dashed line shows the z axis. The first plano-freeform lens with the surface z(x) (b) reshapes the input radiance, whereas the second lens with the surface w(p) (c) reshapes the phase.

Download Full Size | PDF

 

Fig. 8 Normalized radiance distributions generated by the two-lens system shown in Fig. 7 in the planes z = 10 cm (a) and z = 100 cm (b) calculated in TracePro. The radiance cross-sections along the coordinate axes are shown at the top and at the right of the figures.

Download Full Size | PDF

4. Summary

It is shown theoretically and by computed simulation examples that beam shaping with two plano-freeform refractive lenses is feasible even when the lenses are not required to be rotationally symmetric. Utilization of freeform lenses improves dramatically the efficiency of a laser output applied for illumination of a given area with a prescribed irradiation distribution. Since our approach is not constrained by the requirement of rotational symmetry of any component of the designed system, the number of possible applications is also substantially increased in comparison with designs relying on systems with all components rotationally symmetric [1]. In this paper such capabilities are demonstrated for nonisotropic parallel input and output beams. Spherical input with parallel output beams are also possible and will be considered in separate publications.