## Abstract

The SQM is a versatile methodology for designing freeform optics for irradiance redistribution. Recently it has often been applied in design of freeform lenses, mirrors and diffractive optical elements. Still, many questions regarding theory and performance of optics designed with the SQM are open. Here we investigate theoretically plano-freeform refractive lenses designed with the SQM when an incident collimated beam must be transformed into a beam illuminating with prescribed irradiances a large number of pixels on a flat screen. It is shown that a lens designed for such task with the SQM operates as a multifocal lens segmented into subapertures with focal lengths providing accurate control of the irradiance distribution between pixels. These subapertures are patches of hyperboloids of revolution. Two different designs are possible, one of which defines a concave lens. Eikonal function for such lenses is also derived. As a proof of concept, we numerically analyze performance of a plano-freeform lens designed with the SQM for transforming a uniform circular parallel light into an image of A. Einstein represented by gray values at ≈ 38K pixels.

© 2017 Optical Society of America

## 1. Introduction

The irradiance control problem(ICP) shown in Fig. 1 is a problem of design of an optical element for redirecting and reshaping the irradiance profile of an incident light beam so that a prescribed irradiance distribution is generated on a given target. Such design problems arise in many illumination and nonillumination tasks, for example, in automotive and architectural lighting, indoor wireless communications, beam splitting, focusing, metal cutting, etc. In nonimaging optics, problems of this type form one of the two main groups [1]. When the shapes of the input or output beams and irradiance profiles lack particular symmetries, it is natural to turn to freeform optics whose design is not constrained by a priori assumed symmetries. However, a freeform solution to the ICP is hard to obtain because, in general, even in the geometrical optics(GO) approximation such solution must be a differentiable function satisfying a special nonlinear second order partial differential equation(PDE) of Monge-Ampère type with an unusual boundary condition [2]; see section 2. A mathematically rigorous theory for solving such PDE problems is not yet available. Heuristic (usually iterative) approaches to problems of ICP type have been proposed [3, 4], often, with a posteriori justification by computed examples.

The supporting quadric method(SQM) developed in [5, 6] and its more recent extensions [7,8] is a highly versatile, physically transparent and mathematically rigorous methodology for solving optical design problems requiring freeform surfaces. In the SQM framework problems such as the ICP are formulated in geometrical optics(GO) and solved as variational problems with guaranteed convergence to the mathematically rigorous solution. At the same time, dealing directly with the complicated PDEs is avoided. The obtained solution generates a continuous freeform optical surface. In fact, two kinds of surfaces are obtained. The surfaces of the first kind are concave, while the surfaces of the second kind may be neither convex nor concave. When the surface is concave, known results in mathematics guarantee that such a surface is two times differentiable everywhere except, possibly, for a set of zero surface area (sets of “measure zero” in mathematical terms). Additional efforts are usually needed to implement the solution of the variational problem into a numerical procedure [6, 9–11]. Though the SQM has often been used to design freeform optics [9, 12–17], many general questions concerning geometry and diffractive properties of optics designed with the SQM still remain open.

In this paper we provide answers to several of such questions for a refractive plano-freeform lens, designed with the SQM, in the practically important “semi-discrete” case when a collimated beam is transformed into an output beam radiating at a discrete set consisting of an arbitrary number of pixels in a plane with prescribed irradiance concentrations; see Fig. 3(a). Such lens operates as a multifocal lens segmented into subapertures. Each subaperture is a piece of a hyperboloid of revolution with its own focal length determined by the SQM so that accurate control of the irradiance distribution between foci is achieved. It is shown here that the support (=“footprint”) of each subaperture on the planar side of the lens is a convex set whose boundary has a special geometric stucture which is explicitly described. We also obtain a simple expression for the eikonal for each subaperture and of the entire lens. In addition, an explicit relation between the discontinuities in the field of normals across the boundaries of adjacent subapertures and the distance between the foci of these subapertures are established.

As a proof of concept we present the results of analysis of a plano-freeform lens transforming a uniform circular parallel light beam into an image of high complexity, namely, the photographic monochromatic image of A. Einstein represented by gray values on a grid with ≈ 38K pixels. For the designed lens each pixel is a focus for the corresponding subaperture creating the required irradiance. Explicit determination of each of the subapertures allows us to estimate the diffractive spot sizes and assess diffraction effects on the quality of the image generated by the lens more realistically without calculating Fresnel-Kirchhoff integrals for each pixel. These results extend substantially the preliminary results in [18]. In order to make the presentation reasonably self-contained, we included in section 2 a brief review of the GO formulation of the ICP. Because the results in this paper are very tightly connected with the SQM we felt that a quick introduction into the basic concepts behind this method would make the reading of the paper easier; see section 3. Mathematical details are presented in Appendices A and B.

## 2. Problem statement and analytic formulation

Consider the ICP of designing a plano-freeform lens with capability to redirect and reshape a collimated light beam so that a pre-specified irradiance distribution is produced on a given target set in a plane; see Fig. 1. In analytic form problems dealing with redirection and redistribution of energy require determination of a lens performing an optical transformation with controlled Jacobian. In the GO approximation, the required transformation (often called *ray mapping*) is derived using the Snell and the energy conservation laws along tubes of light rays [19] propagating between the source and target regions,
$\overline{\mathrm{\Omega}}$ and
${\overline{T}}_{d}$, respectively, as in Fig. 1. (Everywhere in the paper ovebar over sets means that the sets are closed in the sense of set theory.) Denote this transformation by *P _{d}*. We recall the expressions for

*P*and its Jacobian [2].

_{d}It is assumed that the planar side of the lens is inactive and the light rays change direction only due to the right (active) side. That is, the thin lens approximation is in effect. For brevity, we refer to the active surface of the lens as the lens *R*. In Cartesian coordinates *x* = (*x*_{1}, *x*_{2}), *z* as in Fig. 1, the lens *R* is the graph of a smooth (sag) function *z*(*x*),
$x\in \overline{\mathrm{\Omega}}$. In vector form, *P _{d}* =

*P*+

*d*

**k**, where

**k**is a unit vector in direction of the

*z*axis,

*d*=

*const*> 0, $P:\overline{\mathrm{\Omega}}\to \overline{T}$ is the

*refractor*map with $\overline{T}$ being the projection of ${\overline{T}}_{d}$ on the plane

*α*; thus,

*P*is the projection of

*P*on

_{d}*α*. Often

*P*is more convenient to use than

*P*.

_{d}Put *M*(*z*) := [1 + (1 − *n*^{2})|∇*z*|^{2}]^{1/2}, where ∇*z* := grad*z* and *n* is the refractive index of *R*. It was shown in [2] that

*I*defined in Ω and output irradiance

*L*defined in

*T*(Ω) is obtained using the energy conservation law [2], given by where

_{d}*J*(

*P*) is the matrix of the Jacobian of

_{d}*P*and det is the determinant. Noting that constancy of the vector

_{d}*d*

**k**in

*P*implies

_{d}*J*(

*P*) ≡

_{d}*J*(

*P*), the condition (2) is transformed to

*Thus, analytically, the ICP in the GO formulation consists in finding a function*$z\in {C}^{2}(\mathrm{\Omega})\cap {C}^{1}(\overline{\mathrm{\Omega}})$

*satisfying both equations {Eq. (3), Eq. (4)} with the map P defined by Eq. (1).*[Here ${C}^{2}(\mathrm{\Omega})\cap {C}^{1}(\overline{\mathrm{\Omega}})$ denotes the set of functions one time continuously differentiable in $\overline{\mathrm{\Omega}}$ and twice continuously differentiable in Ω.]

An explicit expression for *J*(*P*) in terms of the sag *z* is given in [2]. It is a strongly nonlinear second order PDE. Since it will not be used it in this paper, we do not reproduce it here.

The problem {Eq. (3), Eq. (4)} is not a standard boundary value problem in PDE’s theory. A mathematically rigorous theory for solving this problem for
$z\in {C}^{2}(\mathrm{\Omega})\cap {C}^{1}(\overline{\mathrm{\Omega}})$ does not exist. Because of strong nonlinearities in Eq. (3) and unusual boundary condition Eq. (4) this problem is hard to analyze directly and, even more so, to solve numerically. The Eq. (3) includes a term with the Hess(*z*) under determinant, and is referred to as an equation of Monge-Ampère type. The known results on solvability of Monge-Ampère equations are not applicable to {Eq. (3), Eq. (4)} because of the particular nonlinear terms in Eq. (3) and unusual constraint Eq. (4).

If the requirement that $z\in {C}^{2}(\mathrm{\Omega})\cap {C}^{1}(\overline{\mathrm{\Omega}})$ in Eqs. (3) and (4) is relaxed then a mathematically rigorous theory for solving an appropriately modified variant of the problem {Eq. (3), Eq. (4)} in the class of “weak” solutions (explained in section 3) can be developed and this was done in [16] using the SQM.

## 3. SQM, semi-discrete ICP and weak solutions

In this section we recall the main steps of the SQM.

#### Step 1

The use of SQM begins with the following basic fact crucial to our construction. Consider a parallel light beam with cross section
$\overline{\mathrm{\Omega}}$ propagating in direction **k** of the *z*-axis. Let the refractive index of the desired lens be *n* > 1. Let, as in Fig. 1, *α* and *α*′ be two parallel planes at a distance *d* > 0 from each other. Fix a number
$\overline{f}<0$. Assume that the target
${\overline{T}}_{d}$ on the plane *α*′ consists of only one arbitrarily picked point
$(\overline{p},d)$ with
$\overline{p}=({x}_{1}(\overline{p})$,
${x}_{2}(\overline{p}))\in \alpha $. Consider a refracting conic surface
${F}_{\overline{p},\overline{f}}$ which is a graph of the sag function

*α*and passes through $\overline{p}$ and $(\overline{p},d)$. The center of ${F}_{\overline{p},\overline{f}}$ is $(\overline{p},h)$, where $h:=\frac{\overline{f}n}{{n}^{2}-1}+d$, and its eccentricity is

*n*[20]. The left and right foci of ${F}_{\overline{p},\overline{f}}$ are, respectively,

**k**hitting ${F}_{\overline{p},\overline{f}}$ from the left refracts and after refraction passes through $(\overline{p},d)$ as in Fig. 2. Thus, $P(x)=\overline{p}$ and ${P}_{d}(x)=(\overline{p},d)$ for any

*x*∈

*α*, that is, $(\overline{p},d)$ is a caustic point of the refracted light beam. This means (and also can be checked by substituting ${z}_{\overline{p},\overline{f}}$ into Eq. (1)) that $J({P}_{d}(x))=J(P(x))=J(P({z}_{\overline{p},\overline{f}}(x)))\equiv [0]$ for any

*x*∈

*α*; in particular, for $x\in \overline{\mathrm{\Omega}}$. Here, [

^{0}] is the zero matrix. Thus, the inverses ${P}_{d}^{-1}$ and

*P*

^{−1}of the maps

*P*and

_{d}*P*are not defined in this case as point-to-point maps.

However, we can define
${P}_{d}^{-1}$ and *P*^{−1} as set valued maps by setting

Clearly,
${P}_{d}^{-1}(\overline{p},d)={P}^{-1}(\overline{p})=\overline{\mathrm{\Omega}}$ and below we discuss only the map *P*^{−1}. For any *p* ∈ *α* put

*I*is the irradiance of the input beam. Thus, $G({z}_{\overline{p},\overline{f}},\overline{p})=\mu $. The conservation law is satisfied here since the total radiation distributed on $\overline{\mathrm{\Omega}}$ is equal to the total irradiance at $(\overline{p},d)$. Therefore, the lens ${F}_{\overline{p},\overline{f}}$ defined by Eq. (5) with

*x*restricted to $\overline{\mathrm{\Omega}}$ now solves the ICP in this case. But instead of solving the problem {Eq. (3), Eq. (4)} we solved for the lens ${F}_{\overline{p},\overline{f}}$ defined by ${z}_{\overline{p},\overline{f}}$ the problem

*δ*is the Dirac delta function. Note that the formula for

*F*in Eqs. (6) implies that the lens ${F}_{\overline{p},f}$ also solves this problem with any

^{r}*f*< 0, that is, the solution of the ICP is not unique.

#### Step 2

Suppose now that
$\overline{\mathrm{\Omega}}$ is as before,
$\overline{T}=\{{p}_{1},\dots ,{p}_{K}\}$, where *p _{i}*,

*i*= 1, …

*K*,

*K*∈ ℕ,

*K*≥ 1, are distinct points on

*α*, and ${\overline{T}}_{d}=\{({p}_{1},d),\dots ,({p}_{K},d)\}$ as in Fig. 3(a). For

*K*hyperboloids with sag functions ${z}_{{p}_{i},{f}_{i}}$,

*f*< 0,

_{i}*i*= 1, …,

*K*, defined as in Eq. (5) with $\overline{p},\overline{f}$ replaced by

*p*,

_{i}*f*, define

_{i}*p*

_{1}, …,

*p*

_{5}is in Fig. 3(b). According to Eq. (7), the lens

*R*defined as the graph of

_{z}*z*(

*x*), $x\in \overline{\mathrm{\Omega}}$, consists of pieces of hyperboloids ${F}_{{p}_{i}.{f}_{i}}$,

*i*= 1, …,

*K*, such that a ray in direction

**k**incident on

*R*is intercepted by the “nearest” ${F}_{{p}_{i}.{f}_{i}}$ and refracted so that it passes through the focus (

_{z}*p*,

_{i}*d*). Thus, in this case the lens is an “envelope” defined by Eq. (7) of a family of suitable hyperboloids.

The lens *R _{z}* defines the map
${P}_{z}:\overline{\mathrm{\Omega}}\to \overline{T}$ and its inverse
${P}_{z}^{-1}:\overline{T}\to \overline{\mathrm{\Omega}}$. Namely,

*z*as in Eq. (7) the total irradiance “delivered” by the lens to (

*p*,

_{i}*d*) is defined as

Let *μ*_{1} ≥ 0, …, *μ _{K}* ≥ 0 be such that

*f*

_{1}< 0, …,

*f*< 0 such that for hyperboloids ${F}_{{p}_{i},{f}_{i}}$,

_{K}*i*= 1, …,

*K*, the lens

*z*defined by Eq. (7) satisfies the conditions

*f*

_{1}, …,

*f*are determined with a provably convergent iterative algorithm [16].

_{K}#### Step 3

The last step of the SQM is the proof that a solution to the ICP for a given distributed irradiance pattern *μ* on
$\overline{T}$ is a limit of a sequence of solutions to semi-discrete problems when discrete distributions converge to *μ*; see [16] for details. In practice, it is the semi-discrete version of the ICP defined by Eqs. (9) that has to be actually solved.

In *summary*, at the end of step 2 the described methodology provides a framework for designing a freeform multifocal lens solving the semi-discrete ICP. This lens is defined by Eq. (7) and consists of subapertures
${F}_{i}=(x,z(x))\equiv (x,{z}_{{p}_{i},{f}_{i}}(x))$,
$x\in {\overline{\mathrm{\Omega}}}_{i}$, *i* = 1, …, *K*, each of which is a part of *R _{z}* over
${\overline{\mathrm{\Omega}}}_{i}$. The subapertures “split” the entire beam into parts

*B*,

_{i}*i*= 1, …,

*K*, with cross sections ${\overline{\mathrm{\Omega}}}_{i}$ and ${\cup}_{i=1}^{K}{\overline{\mathrm{\Omega}}}_{i}=\overline{\mathrm{\Omega}}$. For each ${p}_{i}\in \overline{T}$ the rays of

*B*emitted by

_{i}*x*∈ Ω

*, where Ω*

_{i}*is the interior of ${\overline{\mathrm{\Omega}}}_{i}$, are refracted by*

_{i}*F*so that they are focused at (

_{i}*p*

_{i},

*d*). The behavior of the subaperture “boundary” rays emitted by $x\in \partial {\overline{\mathrm{\Omega}}}_{i}$ is discussed in section 4, subsection 4. It follows from the discussion above that the SQM produces a continuous in $\overline{\mathrm{\Omega}}$ solution

*z*of the semi-discrete ICP. Because

*z*is not necessarily differentiable, it is called a “weak” solution.

## 4. Properties of weak solutions to semi-discrete ICP

Let *z* be a weak solution defined by Eq. (7) of the ICP Eqs. (9). Assume that
$\overline{\mathrm{\Omega}}$ is convex.

#### 1. Concavity and differentiability

By construction, *z* is continuous in
$\overline{\mathrm{\Omega}}$. Because the minimum in Eq. (7) is taken over a set of concave functions, the function *z* is also concave [21]. This is the Keplerian configuration when refracted rays crossover [22]. Replacing min by max in Eq. (7) we obtain new solutions of ICP which are continuous and piece-wise concave but globally may be neither convex nor concave (the refracted rays in this case don’t crossover). Most of the results below are valid for both types of solutions. For definiteness, for the rest of the paper and without further reminding we discuss only the concave solutions defined by Eq. (7).

Because for each *i* the function
$z(x)\equiv {z}_{{p}_{i},{f}_{i}}(x)$,
$x\in {\overline{\mathrm{\Omega}}}_{i}$, and
${z}_{{p}_{i},{f}_{i}}$ is infinitely differentiable on the entire plane *α*, the function *z* is infinitely differentiable in Ω* _{i}*. Denote by ∂Ω

*the boundary of ${\overline{\mathrm{\Omega}}}_{i}$. For*

_{i}*x*∈ ∂Ω

*the “one-sided” partial derivatives are defined but the usual derivatives of*

_{i}*z*of the first and higher orders may not be defined on ${\mathrm{\Omega}}_{ij}:={\overline{\mathrm{\Omega}}}_{i}\cap {\overline{\mathrm{\Omega}}}_{j}$ for

*i*≠

*j*and Ω

*≠ ∅.*

_{ij}#### 2. Focal function

It is convenient to introduce the “focal” function of the lens defined by Eq. (7). Let
$f:\overline{T}\to (-\infty ,0)$ and *f _{i}*:=

*f*(

*p*). The pair (

_{i}*z*,

*f*) is called a

*refractor pair*if (7) and the relations

*i*= 1, …,

*K*. The function

*f*is called the

*focal*function of the lens

*R*. Note that the last Eq. defines

_{z}*f*implicitly. It follows from Eq. (7) that

_{i}*d*−

*z*(

*p*) = −

_{i}*f*/(

_{i}*n*− 1) = the distance from the apex of ${F}_{{p}_{i},{f}_{i}}$ to the focal plane

*α*′.

Fix an *i* ∈ {1, …, *K*}. A direct check shows that for any
$x\in {\overline{\mathrm{\Omega}}}_{i}$ the expression in square brackets in Eq. (10) is equal to *nf _{i}*/(

*n*

^{2}− 1), that is, the supremum in Eq. (10) is attained at each $x\in {\overline{\mathrm{\Omega}}}_{i}$. By Eq. (7) for any

*j*∈ {1, …,

*K*},

*j*≠

*i*and

*x*∈ Ω

*we have $z(x)={z}_{{p}_{j},{f}_{j}}(x)<{z}_{{p}_{i},{f}_{i}}(x)$. This implies that the expression in square brackets evaluated at such*

_{j}*x*is strictly less than

*nf*/(

_{i}*n*

^{2}− 1), that is, for a given

*i*the supremum in Eq. (10) is attained only in ${\overline{\mathrm{\Omega}}}_{i}$. Since

*z*is differentiable in Ω

*, Eq. (10) implies that*

_{i}*P*(

*x*) =

*p*at any

_{i}*x*∈ Ω

*, and*

_{i}*x*∈ ∂Ω

*provided ∇*

_{i}*z*(

*x*) is defined as a one-sided gradient at

*x*. When $x\in {\overline{\mathrm{\Omega}}}_{ij}$ such value of ∇

*z*(

*x*) is not unique and the map

*P*is not single valued. Property 4(c) in subsection 4 implies that the sets Ω

*have zero areas.*

_{ij}#### 3. Normals on boundaries of subapertures

In Appendix A, it is shown that if the points *p _{i}* and

*p*in the target $\overline{T}$ are sufficiently close to each other and

_{j}*F*∩

_{i}*F*≠ ∅ then for $x\in {\overline{\mathrm{\Omega}}}_{ij}$ the normals

_{j}**n**

*(*

_{i}*x*) to

*F*and

_{i}**n**

*(*

_{j}*x*) to

*F*at any

_{j}*x*∈ Ω

*satisfy the inequality*

_{ij}*C*> 0 which can be chosen the same for all points in $\overline{\mathrm{\Omega}}$ and $\overline{T}$. This property shows that a higher resolution of the target reduces the discontinuities of the normals. It also suggests that for an irradiance distribution on $\overline{T}$ with continuous and positive density the solution

*z*of the ICP has continuous first derivatives in Ω.

#### 4. Boundaries of subapertures

For *p _{i}* ≠

*p*the set Ω

_{j}*is the projection of Γ*

_{ij}*:=*

_{ij}*F*∩

_{i}*F*. Assume Γ

_{j}*≠ ∅ and not a single point. We claim:*

_{ij}*4(a)*Γ

*is a hyperbola in a plane parallel to*

_{ij}**k**;

*4(b)*each ${\overline{\mathrm{\Omega}}}_{i}$ is a convex set;

*4(c)*the boundary $\partial {\overline{\mathrm{\Omega}}}_{i}$ is a union of straight line segments and possibly pieces of $\partial {\overline{\mathrm{\Omega}}}_{i}\cap \partial \overline{\mathrm{\Omega}}$. Proofs of these properties are in Appendix B. A simple example of a segmented beam is shown in Fig. 4. In this example, the freeform lens splits a uniform circular beam shown in Fig. 4(a) so that the

*p*

_{1}, …,

*p*

_{7}in the plane

*α*′:= {

*z*= 150

*mm*} are illuminated with irradiances indicated in Fig. 4(b). This lens consists of seven subapertures.

## 5. Eikonal and paraxial approximation

#### 1. Derivation of the eikonal

The eikonal (=phase function) does not explicitly enter the strong or the weak GO formulations of the ICP. However, the sag function *z* of the lens in the GO solution can be converted to the phase shift due to the lens by multiplying it by (*n* − 1)*k*, where *k* = 2*π*/*λ* is the wave number [23]. For the rotationally symmetric case and when the input beam is Gaussian and output beam is uniform, Hoffnagle [23] has shown that in paraxial limit the eikonal of the GO solution is identical to the one derived using Fourier optics. Here we derive the eikonal of a GO solution to the ICP without any symmetry assumptions.

The lens *R _{z}* designed with the SQM is defined by Eq. (7) and the wavefront is generated by subapertures
${F}_{{p}_{i},{f}_{i}}$,

*i*= 1, 2, …,

*K*. From the point of view of GO each ${F}_{{p}_{i},{f}_{i}}$ focuses the light into the focus (

*p*,

_{i}*d*). Because the eikonal is determined up to a constant, we may assume that the region $\overline{\mathrm{\Omega}}$ in the plane

*α*: {

*z*= 0} is the planar side of the lens

*R*. The material of the lens is assumed to have refracting index

_{z}*n*which is also the eccentricity of each ${z}_{{p}_{i},{f}_{i}}$,

*i*= 1, 2, …,

*K*. Thus the geometry of the lens

*R*is consistent with the refractive properties of the material.

_{z}Fix an *i* ∈ {1, …, *K*}. We construct the eikonal
${\mathrm{\Psi}}_{{p}_{i},{f}_{i}}$ for the subaperture *F _{i}*. The reference phase is taken as
$n{z}_{{p}_{i},{f}_{i}}({p}_{i})=n\left(\frac{{f}_{i}}{n-1}+d\right)$. The phase shift at

*x*due to the interval between

*z*= 0 and $z=\frac{{f}_{i}}{n-1}+d$ which “includes” the lens, is calculated at the plane $z=\frac{{f}_{i}}{n-1}+d$ as

*R*as

_{z}The sag *z* is continuous in
$\overline{\mathrm{\Omega}}$. Hence the discontinuities in
${\mathrm{\Psi}}_{{R}_{z}}$ are due to discontinuities in the focal function *f* (*p*),
$p\in \overline{T}$. On the other hand, by Eq. (5)*z*(*p _{k}*)−

*d*= −

*f*/(

_{k}*n*−1),

*k*=

*i*,

*j*, and then |

*f*−

_{i}*f*| = (

_{j}*n*− 1)|

*z*(

*p*) −

_{i}*z*(

*p*)|. By Lemma 6.1 in [16], for any

_{j}*p*, ${p}_{i}\in \overline{T}$ we have $|z({p}_{i})-z({p}_{j})|\le |{p}_{i}-{p}_{j}|/\sqrt{{n}^{2}-1}$. Thus,

_{i}*p*−

_{i}*p*| decreases. By property 4, section 4, the discontinuities may occur only at hyperbolas or points defined by

_{j}*F*∩

_{i}*F*. Such sets have zero areas. In some applications [17] the hyperboloids can be selected so that

_{j}*f*

_{1}= … =

*f*. In these cases ${\mathrm{\Psi}}_{{R}_{z}}$ is continuous over the entire aperture.

_{K}For future reference we record here an expression for the optical path length(OPL) of a lens designed with SQM. Clearly, for each *i* = 1, …, *K* and all
$x\in {\overline{\mathrm{\Omega}}}_{i}$ we have

*we obtain*

_{i}#### 2. Paraxial approximation of the eikonal

Applying the standard approximation for the square root in *z _{p,f}* and taking into account that

*f*< 0 we obtain for each

*i*= 1, 2, …,

*K*

*x*−

*p*| << |

_{i}*f*|, $x\in {\overline{\mathrm{\Omega}}}_{i}$ (as is the case in paraxial approximation), then the term $O\left(\frac{{|x-{p}_{i}|}^{4}}{{|{f}_{i}|}^{4}}\right)$ is negligible and the approximate expression for the eikonal (denoted below by ${\widehat{\mathrm{\Psi}}}_{{R}_{z}}$) is

_{i}*f*.

_{i}## 6. Specific computational design - image of A. Einstein

The photographic image of A. Einstein in Fig. 5(c) from the internet is given as a set
$\overline{T}=\{{p}_{1},\dots ,{p}_{K}\}$ of *K* = 37, 904 pixels with gray values *μ*_{1}, …, *μ _{K}* representing brightness of each pixel in the range from 0 (dark) to 255 (maximal brightness). The pixels form a rectangular 184 × 206 array. In our example the given gray values were treated as prescribed irradiances. The goal was to determine a plano-freeform lens which intercepts, redirects and reshapes a circular, uniform, collimated light beam into illumination pattern in the plane

*α*′= {

*z*= 150mm} producing the rectangular array with gray values

*μ*

_{1}, …,

*μ*. The schematic in Fig. 3(a) applies in this example. The results are shown graphically in Fig. 5 (the pictures are not to scale).

_{K}#### Overview of parameters and results

This ICP was formulated as a semi-discrete problem in the form given by Eqs. (9) and solved numerically in the SQM framework on a standard desktop computer. The numerical results were obtained using a system of algorithms and computational codes developed jointly by the present author and Dr. B. Cerkasskiy.

The dimensions of the image
$\overline{T}:61\text{mm}\times 68.3333\text{mm}$ with distance of 0.3333mm between adjacent pixels in *x*_{1} and *x*_{2} directions; refracting index *n* = 1.4607 (fused silica at wavelength *λ* = 532nm); the input beam with cross section
$\overline{\mathrm{\Omega}}$ was modeled with 515, 431 point sources of equal strength distributed uniformly (in radial direction and on each azimuthal circle) on a disk of radius 20.05mm centered at (0, 0) on the reference plane *α* := {*z* = 0}; the target set is in the plane *α*′; the curved side of the designed lens, denoted, again, by *R _{z}* is concave. Its peak-to-valley width (in

*z*− direction) is 10.397666mm; the distance from the plane

*α*′ to the apex of the lens is 105.055008mm. No despecling of the synthetic image and no smoothing by interpolation or any other method was performed on the synthetic image nor on the lens.

Denote by
$\widehat{G}\left(z,{p}_{i}\right)$ the irradiance at *p _{i}* obtained by numerically simulated ray tracing through the designed lens. Put
$\widehat{G}:=(\widehat{G}\left(z,{p}_{1}\right),\dots ,\widehat{G}\left(z,{p}_{K}\right))$ and

*μ*:= (

*μ*

_{1}, …,

*μ*). Tables 1 and 2 below shows values of several validation parameters for our example. The notations are: the symbol ║ · ║ indicates the RMS value and $\overline{\mu -\widehat{G}}$ is the mean value. The small value of Δ/║

_{K}*μ*║ shows the high accuracy of performance of the designed lens.

In the next table we present results of calculated optical path values OPL* _{i}*,

*i*= 1, …,

*K*of the subapertures. Note that our computational algorithm finds the solution satisfying the irradiance requirements on the target and the OPL

*are defined by the computed sag but not constrained a priori; the values in the table are obtained from the numerical solution. In this table: OPL := (*

_{i}*OPL*

_{1}, …,

*OPL*),

_{K}*δ*:= ║

*OPL*║,

*β*:= ║(

*OPL*

_{1}−

*δ*, …,

*OPL*−

_{K}*δ*)║. Table 2 shows that the differences in the OPL

*for different subapertures is relatively small.*

_{i}The values
$\widehat{G}\left(z,{p}_{i}\right)$, *i* = 1, …, *K*, were also used to generate the Fig. 5(a) to allow visual comparison with the original photograph. All light rays incident on the lens were refracted and redistributed on the target rectangle
$\overline{T}$. Thus, the theoretical energy efficiency of the designed lens is ≈ 100%. The fine details in Fig. 5 (a) demonstrate that our approach to ICPs can be used in problems in which high resolution of the reconstructed image is required.

#### Diffraction effects

Since a lens designed with the SQM functions as a multifocal lens, it is natural to attempt to evaluate the mutual diffraction effects due to subapertures. The SQM allows to determine each subaperture of a designed lens and, in principle, one can calculate the Fresnel-Kirchhoff integral for each subaperture. However, for a freeform lens with a large number of subapertures this approach may not be practical. The supports of the subapertures while convex, don’t have special symmetries and efficient and accurate analysis of diffraction effects deserves a separate investigation. Below, we present some observations in the case of the image of A. Einstein.

We are aware of only two papers discussing diffraction effects for freeform optics designed with the SQM. In [17] the sizes of diffraction spots were estimated for diffractive optical elements in the Fresnel regime under a priori assumption that diffraction spots are identical squares and in [18] a refracting lens with five foci was considered in the Fraunghofer regime under assumption that supports of subapertures and diffraction spots are circular. Their radii were estimated by the prescribed irradiances using energy conservation. Fig. 4(a) shows that the supports of subapertures may be neither square not circular. Diffraction effects, mainly, for freeform mirrors producing one-dimensional patterns have also been discussed in [24].

Here, we use the diameter and widths in the *x*_{1} and *x*_{2} coordinate directions of a subaperture support to generate “one-number” estimates allowing us to identify subapertures with diffraction spots overlapping with adjacent spots in *all or in coordinate directions*. Such estimates are certainly approximate but may be useful for initial analysis of designs. For pixels in the image of A. Einstein our estimates of diffraction spot sizes are constructed as follows. Recall that
${\overline{\mathrm{\Omega}}}_{i}$ denotes the planar support of the subaperture *F _{i}* with focus (

*p*,

_{i}*d*),

*i*= 1, …,

*K*. Denote by

*τ*the diffraction spot formed by

_{i}*F*at (

_{i}*p*,

_{i}*d*) on the focal plane

*α*′. In the SQM for a given solution we determine for each

*i*= 1, …,

*K*the supports ${\overline{\mathrm{\Omega}}}_{i}$, their diameters

*D*, and widths ${D}_{{x}_{1}i}$ and ${D}_{{x}_{2}i}$ in coordinate directions. By analogy with the one-dimensional Fresnel diffraction [25] the diffraction spot diameter for each

_{i}*τ*is calculated as DSD

_{i}*:=*

_{i}*f*/

_{i}λ*D*, where

_{i}*f*is the focal length of

_{i}*F*and

_{i}*λ*is the wavelength immediately after the lens. The latter is calculated from the given wavelength of the input beam and refractive index of the lens. Similarly, the diffraction spot widths in

*x*

_{1}and

*x*

_{2}directions are calculated as ${f}_{i}\lambda /D{S}_{{x}_{1}i}$ and ${f}_{i}\lambda /D{S}_{{x}_{2}i}$.

The chart in Fig. 6 shows the cumulative distribution of pixels vs. DSD* _{i}* for the image in Fig. 5(a) produced by the lens in Fig. 5(b). It allows to estimate the number of pixels with DSD

*≥*

_{i}*γ*, where

*γ*∈ (0, 0.76). For example, for

*γ*= 0.3333 there are ≈ 6, 000 pixels (≈ 16%) with DSD

*> 0.3333. As noted earlier, the adjacent pixels in the image Fig. 5(a) are ≈ 0.3333mm from each other in the*

_{i}*x*

_{1}and

*x*

_{2}directions. Thus, diffraction spots of pixels with DSD

*larger than 0.3333 may overlap with diffraction spots of neighboring pixels. Similarly, overlaps may occur in*

_{i}*x*

_{1}and

*x*

_{2}directions for pixels with $D{S}_{{x}_{1}}$, $D{S}_{{x}_{2}}\ge 0.3333$. The inter-relations between diffraction spots of pixels in the image can be seen in the color maps in Fig. 7 showing the distributions of

*DSD*, $D{S}_{{x}_{1}}$, $D{S}_{{x}_{2}}$ for all pixels in the image Fig. 5(a) (left, middle and right images in Fig. 7, respectively).

Because the irradiance of the input beam is nearly uniform, the subapertures illuminating pixels with small prescribed irradiances have small areas and consequently large diffraction spots. For example, for the image in Fig. 5(c) the gray values for pixels with DSD* _{i}* > 0.3333 are in the range between 1 and 34. These pixels are shown in black in Fig. 5(a) (and in Fig. 5(c)) and in magenta and black in Fig. 7. On the average, the gray values in that range are on the order of 0.08% of the overall maximum. The majority of such pixels are near the boundary of the image; cf. with Fig. 5(a)). Thus, diffraction spots of large size have limited effect on the image quality.

## 7. Summary

The irradiance control problem arising in many applications is considered in the case when the freeform refracting lens transforms a collimated light beam into a beam generating a prescribed irradiance distribution at a given set of pixels in the focal plane. The analytical GO formulation of such design problem leads to a second order nonlinear PDE and boundary condition which can not be solved directly and reliably by presently available methods. By contrast, the supporting quadric method(SQM) provides a mathematically rigorous and physically transparent framework for determination of weak solutions defining freeform lenses producing the prescribed irradiance redistribution, while avoiding the difficulties connected with the approach involving PDEs. New geometric and diffractive properties of freeform refracting lenses designed with the SQM are established and investigated in this paper. These properties, combined with the fact that such lenses are segmented into subapertures which are pieces of hyperboloids of revolution, show that their geometry is very simple and they operate as multifocal lenses. This facilitates their performance evaluation. The particularly simple structure of subapertures of such lenses should be also useful at the fabrication stage. Explicit expressions for the eikonal of such a lens and of its subapertures are also derived. In an illustrative example, we analyze numerically the properties of a freeform lens designed with the SQM for transforming the input irradiance of a uniform collimated beam into an irradiance distribution representing the photographic image of A. Einstein. In this example the image consists of 37, 904 pixels with given gray values which were used as prescribed irradiances. Performance of the lens was validated by simulated ray tracing, comparison of prescribed irradiances with approximate irradiances determined by ray tracing, and by approximate estimates of diffraction effects due to the subapertures of the lens at each pixel; see Fig. 5, Tables 1, 2 and Fig. 6, 7. Most of results extend to one- and two-component optics designed with the SQM for applications with targets in the near- and far-field.

## 8. Appendices

## Appendix A

Let
$\overline{x}\in {\mathrm{\Omega}}_{ij}$. For *k* = *i*, *j* put

*z*the unit normal to its graph is the vector (−∇

*z*, 1)[1 + |∇

*z*|

^{2}]

^{−1/2}. Applying Eq. (11) we obtain ${\tilde{S}}_{k}{\mathbf{n}}_{k}=(\overline{x}-{p}_{k},{S}_{k})$,

*k*=

*i*,

*j*, and then

*z*(

*p*) −

_{k}*d*= −

*f*/(

_{k}*n*− 1),

*k*=

*i*,

*j*, imply that By Lemma 6 in [16], Combining this with Eq. (17) we obtain Thus, if |

*p*−

_{i}*p*| ≤

_{j}*σ*for some small

*σ*> 0 then Since ${S}_{i}^{2}-{S}_{j}^{2}={f}_{i}^{2}-{f}_{j}^{2}+({n}^{2}-1)[{(x-{p}_{i})}^{2}-{(x-{p}_{j})}^{2}]$ we conclude that |

*S*−

_{i}*S*|≤

_{j}*C*

_{0}|

*p*−

_{i}*p*| for some positive constant

_{j}*C*

_{0}that can be chosen the same for all

*p*

_{1}, …,

*p*. A similar estimate holds also for $|{\tilde{S}}_{i}-{\tilde{S}}_{j}|$. Eq. (16) implies now the property claimed in section 4, subsection 3.

_{K}## Appendix B

For *p _{i}* ≠

*p*the set Ω

_{j}*is the projection of Γ*

_{ij}*:=*

_{ij}*F*∩

_{i}*F*. Assume Γ

_{j}*≠ ∅ and not a single point. We claim:*

_{ij}*4(a)*Γ

*is a hyperbola in a plane parallel to*

_{ij}**k**;

*4(b)*each ${\overline{\mathrm{\Omega}}}_{i}$ is a convex set;

*4(c)*the boundary $\partial {\overline{\mathrm{\Omega}}}_{i}$ is a union of straight line segments and possibly pieces of $\partial {\overline{\mathrm{\Omega}}}_{i}\cap \partial \overline{\mathrm{\Omega}}$. We begin with

*4(a)*.

For any
$p\in \overline{T}$ and *f* < 0 the surface *F _{p}*,

*is a hyperboloid of revolution with axis passing through*

_{f}*p*and parallel to

**k**. The intersection of

*F*,

_{p}*with the plane {*

_{f}*z*=

*m*} for any real

*m*is either empty, or a point or a circle. It follows from Eq. (5) that for a fixed

*m*such that $0\le m\le {\mathrm{max}}_{{\mathrm{\Omega}}_{ij}}{z}_{{p}_{k},{f}_{k}}(x)$,

*k*=

*i*,

*j*, the Eqs.

*α*. These circles may degenerates into a point if

*m*=

*d*+

*f*/(

_{k}*n*− 1). The set of all {

*x*} satisfying both Eqs. in (18) is the set Ω

*. Subtracting the equalities in (18) from each other, we obtain*

_{ij}*x*is a function of

*m*only, provided

*p*,

_{i}*f*,

_{i}*p*,

_{j}*f*,

_{j}*d*,

*n*are fixed. On the other hand, for each $m\in [0,ma{x}_{{\mathrm{\Omega}}_{ij}}{z}_{{p}_{i},{f}_{i}}(x))$ the two circles defined by Eqs. (18) intersect at two distinct points

*x*(

*m*). At $m=ma{x}_{{\mathrm{\Omega}}_{ij}}{z}_{{p}_{i},{f}_{i}}(x)$ such circles intersect at one point. The point $\overline{x}$ corresponding to $m=ma{x}_{{\mathrm{\Omega}}_{ij}}{z}_{{p}_{i},{f}_{i}}(x)$ is defined uniquely and it divides Ω

*into two connected subsegments. Each of these subsegments is parametrized by $m\in [0,\phantom{\rule{0.2em}{0ex}}ma{x}_{{\mathrm{\Omega}}_{ij}}{z}_{{p}_{i},{f}_{i}}(x)]$ and, in particular, we can differentiate (20) with respect to*

_{ij}*m*on each of this subsegments. Then, that is, any

*x*∈ Ω

*lies on a straight line in plane*

_{ij}*α*forming a constant angle with the vector

*p*−

_{i}*p*. Consequently, the curve Γ

_{j}*(whose projection on*

_{ij}*α*is Ω

*) lies in a plane parallel to the*

_{ij}*z*–axis which intersects ${F}_{{p}_{i},{f}_{i}}$ and ${F}_{{p}_{j},{f}_{j}}$. That is, Γ

*is a hyperbola in that plane.*

_{ij}It remains to prove *4b* and *4c*. Denote by
${\widehat{F}}_{{p}_{k},{f}_{k}}$ the convex body bounded by
${F}_{{p}_{k},{f}_{k}}$, *k* = *i*, *j*. Then the sets
${\widehat{F}}_{{p}_{k},{f}_{k}}\cap \overline{\mathrm{\Omega}}$, *k* = *i*, *j*, and
${\widehat{F}}_{{p}_{i},{f}_{i}}\cap {\widehat{F}}_{{p}_{j},{f}_{j}}$ are convex. Therefore, each
${\overline{\mathrm{\Omega}}}_{i}$ is the intersection of convex sets and therefore convex. The boundary
$\partial {\overline{\mathrm{\Omega}}}_{i}$ is a union of linear segments and possibly pieces in
$\partial {\overline{\mathrm{\Omega}}}_{i}\cap \partial \overline{\mathrm{\Omega}}$. Thus, the proofs of the claimed properties are now complete.

## Funding

This work was partially supported by the US-Israel Binational Science Foundation (BSF) under Grant 2010217.

## Acknowledgments

This paper is a substantially expanded version of an invited talk presented at the SPIE meeting, San Diego, CA, August 28, 2016. It is a pleasure to thank Professor Roland Winston for very useful discussions of the subject of this paper. The paper was written during Winter-Spring of 2016 while the author was a Visiting Professor at the Israel Institute of Technology (Technion), Haifa, Israel.

## References and links

**1. **R. Winston, J. C. Miñano, P. Benítez, N. Shatz, and J. Bortz, *Nonimaging Optics* (Elsevier, 2005).

**2. **V. I. Oliker, “Differential equations for design of a freeform single lens with prescribed irradiance properties,” Opt. Eng. **53**(3), 031302 (2014). [CrossRef]

**3. **H. Ries and J. Muschaweck, “Tailored freeform optical surfaces,” J. Opt. Soc. Am. A **19**(3), 590–595 (2002). [CrossRef]

**4. **R. Wu, L. Xu, P. Liu, Y. Zhang, Z. Zheng, H. Li, and X. Liu, “Freeform illumination design: a nonlinear boundary problem for the elliptic Monge-Ampére equation,” Opt. Lett. **38**(2), 229–231 (2013). [CrossRef] [PubMed]

**5. **L.A. Caffarelli and V.I. Oliker, “Weak solutions of one inverse problem in geometric optics,” J. Math. Sci. **154**(1), 37–46 (2008).

**6. **V. I. Oliker, “Mathematical aspects of design of beam shaping surfaces in geometrical optics,” in *Trends in Nonlinear Analysis*, M. Kirkilionis, S. Krömker, R. Rannacher, and F. Tomi, eds. (Springer, 2003), pp. 191–222.

**7. **V. I. Oliker, “Geometric and variational methods in optical design of reflecting surfaces with prescribed irradiance properties,” Proc. SPIE **5942**, 594207 (2005). [CrossRef]

**8. **V.I. Oliker, “Designing freeform lenses for intensity and phase control of coherent light with help from geometry and mass transport,” Arch. Rational Mech. Anal. **201**, 1013–1045 (2011). [CrossRef]

**9. **S. Kochengin and V. I. Oliker, “Determination of reflector surfaces from near-field scattering data II. Numerical solution,” Numer. Math. **79**(4), 553–568 (1998). [CrossRef]

**10. **V. I. Oliker and B. V. Cherkasskiy, “Controlling light with freeform optics: recent progress in computational methods for optical design of freeform lenses with prescribed irradiance properties,” Proc. SPIE **9191**, 919105(2014). [CrossRef]

**11. **P. M. M. de Castro, Q. Mérigot, and B. Thibert, “Far-field reflector problem and intersection of paraboloids,” Numer. Math. **134**(2), 389–411 (2016). [CrossRef]

**12. **F. R. Fournier, Freeform reflector design with extended sources (Ph.D. Diss., University of Central Florida, 2010).

**13. **D. Michaelis, P. Schreiber, and A. Bräuer, “Cartesian oval representation of freeform optics in illumination systems,” Opt. Lett. **36**(6), 918–920 (2011). [CrossRef] [PubMed]

**14. **W. Song, D. Cheng, Y. Liu, and Y. Wang, “Free-form illumination of a refractive surface using multiple-faceted refractors,” Appl. Opt. **54**(28), E1–E7 (2015). [CrossRef] [PubMed]

**15. **C. Canavesi, *Subaperture conics and geometric concepts applied to freeform reflector design for illumination* (Ph.D. Diss., University of Rochester, 2014).

**16. **V. I. Oliker, J. Rubinstein, and G. Wolansky, “Supporting quadric method in optical design of freeform lenses for illumination control of a collimated light,” Adv. in Appl. Math. **62**, 160–183 (2015). [CrossRef]

**17. **L. L. Doskolovich, M. A. Moiseev, E. A. Bezus, and V. I. Oliker, “On the use of the supporting quadric method in the problem of the light field eikonal calculation,” Opt. Express **23**(15), 19605–19617 (2015). [CrossRef] [PubMed]

**18. **M. N. Ricketts, R. Winston, and V. Oliker, “Diffraction effects in freeform optics,” Proc. SPIE **9572**, 957200 (2015).

**19. **M. Born and E. Wolf, *Principle of Optics*, 7th (expanded) edition (Cambridge University, 1999). [CrossRef]

**20. **R. Luneburg, *Mathematical Theory of Optics* (University of California, 1964).

**21. **R. T. Rockafellar and R. J-B Wets, *Variational Analysis*, Grundlehren der Math. Wiss. 317 (Springer, 2009).

**22. **D. L. Shealy and J. A. Hoffnagle, “Geometrical methods,”), *Laser Beam Shaping: Theory and Techniques*, 2-nd edition, F.M. Dickey, ed. (CRC, 2014), Chap. 6.

**23. **J.A. Hoffnagle, “A new derivation of the Dickey-Romero-Holswade phase function,” Proc. SPIE **5876**, 587606 (2005). [CrossRef]

**24. **S. Zwick, R. Feßler, J. Jegorov, and G. Notni, “Resolution limitations for tailored picture-generating freeform surfaces,” Opt. Express **20**(4), 3642–3653 (2012). [CrossRef] [PubMed]

**25. **V. A. Soifer, V.V Kotlar, and L.L. Doskolovich, *Iterative Methods for Diffractive Optical Elements Computation* (CRC, 1997).