## 1.

## Introduction

Construction of a freeform surface to reallocate the radiant flux is a substantial challenge in terms of nonimaging optics. The basic concept is to convert the radiant flux emitted from a point source onto a prescribed illumination, usually in Cartesian coordinate. Parkyn^{1} firstly pioneered a tessellation method to map the cell area between a point source in polar grid with a particular rectangular illumination in Cartesian grid. The cells of both grid topologies were sized so that they share common energy. The methodology was based on the assumption that when the polar cell on the light source was exactly mapped to the corresponding rectangular target cell, the flux of the angular cell would be fully transferred to the target cell accordingly. Each corresponding pair of cells on these grids has an associated surface normal vector, and then governs the surface shape using extrinsic differential geometry. The major concern about equal flux grid method lies in a fact that the constructed surface, either reflective or refractive, is unlikely to guarantee the integrability condition in the overall surface, and thus leads to a discontinuous boundary. The cumulative deviations from the ideal surface normal become more pronounced at the far corners of the output grid. The discontinuous boundary could lead to a problem. The problem is that a small amount of radiant flux will be distorted by the sharp cliff. Hence, sufficient partition cells shall be numbered to confirm the estimated accuracy for surface construction.

To alleviate the discontinuous issue of freeform, Wang et al.^{2} introduced the variable separation mapping strategy to take into account the normal deviation of the curve below a present value. The surface curve was regenerated with tangent vectors along the curve perpendicular to the normal vector, thus eliminating the accumulated normal deviation. The proposed scheme constructs a smooth surface by interpolating the discontinuity surface with non-uniform rational basis spline (NURBS).^{3} On the other hand, the integrability condition $\eta =N\xb7(\nabla \times N)=0$ is necessary in order to make a smooth freeform surface, where $N$ is the surface normal.^{4} To fulfill the integrability condition, Fournier^{3} interpolated the ellipsoids, obtained from Oliker’s algorithm with appropriate grid partition, by using the NURBS where parameters of the NURBS were determined by minimizing the residual curl $\eta $. Although the surface can be made smooth and the resulted radiant flux can be less distorted via the above mentioned approaches, the amount of energy deviation is still difficult to be estimated before the freeform surface is constructed. Furthermore, the computational costs of the constrained optimization for the NURBS’ parameters may rise significantly and the convergence is not guaranteed as the number of grid partition increases.

In additional to the tessellation method, the freeform surface design problem is also related to the solution of the well-known Monge-Ampere (MA) equation. In 1972, Schruben first derived the MA equation obtained from reflector design. Beginning in 1980s until the present, Oliker et al. formulated the freeform problem into a highly nonlinear MA type of partial differential equation (PDE).^{5}6.7.8.^{–}^{9} The MA PDE also appears in other research fields, such as differential geometry, optimal control and mass transportation, meteorology, and geostrophic fluid.^{10}11.^{–}^{12} Although many works about existence and regularity of MA PDE have been conducted, numerical computation of freeform surface reconstruction by solving the MA PDE remains a challenge in practice. We are curious if there exists a systematic way to construct the freeform in which not only the surface is guaranteed to be smooth, but also the correspondence between light source and target plane is ensured under a controllable energy deviation.

In this article, we formulate the ray forward and inverse propagation with negligible energy deviation within the MA PDE framework. By following the works of Caffarelli and Oliker^{13}, the freeform considered here is assumed to be convex. The model is based on local energy conservation in which the flux correspondence associated with each local domain is prescribed. We take advantage of mathematical analysis of the MA PDE and successfully formulate the freeform design problem into the popular finite element paradigm. The formulation could maintain the ability to control the light beam and guarantee a smooth freeform surface. Moreover, this article describes this in detail by going through the different modeling steps of the model. All the main restrictions influencing the target distribution, collective solid angle, or source distribution are discussed, alongside the viability of realizing an optical surface for different specific requests.

## 2.

## Freeform Model in Localized PDE

To develop a method for freeform construction subject to the source-target information, we must deduce an adequate, physically realistic model of ray propagation in the free space. The model can then be employed to understand the reflective (refractive) forward process of rendering illuminance according to a known optical configuration. Also, the model is applicable to the inverse algorithm of estimating the freeform surface from a given source-target corresponding condition. First of all, the ray transportation must satisfy the energy conservation among coordinate transformation.^{14}^{,}^{15} The radiant intensity $I$ ($\text{Watts}/\mathrm{sr}$) of a point source (in polar domain ${\mathrm{\Omega}}_{\theta}$) and prescribed illuminance distribution $E$ ($\text{Watts}/{\mathrm{m}}^{2}$) of the target region (in Cartesian domain ${\mathrm{\Omega}}_{x}$), in one dimension, can be related as

## (1)

$${\int}_{{\mathrm{\Omega}}_{\theta}}I(\theta )\mathrm{d}\theta ={\int}_{{\mathrm{\Omega}}_{x}}E(x)\mathrm{d}x.$$The goal is to design a freeform surface by analyzing the potential and limitation of using the PDE. In this article, the freeform surface is assumed to be continuous and convex. By partitioning the global domain ${\mathrm{\Omega}}_{\theta}$ and ${\mathrm{\Omega}}_{x}$ into many local subdomains ${\mathrm{\Omega}}_{\theta}^{i}$ and ${\mathrm{\Omega}}_{x}^{i}$, so that each corresponding pair satisfies the equal flux $i=1,2,\dots ,n$, we are able to ensure the convergence of the integration. As long as the number of partitions is sufficiently large, the freeform can be considered to be a flat surface. As a result, a reduced but simplified representation of a freeform can be derived on each local domain. Hereafter, a flat surface criterion to the freeform in each local domain is called a small local planar approximation.

The representation of the local freeform is indeed a nonlinear differential equation and shall be derived in the following. First, let us introduce some notations. Global coordinate frame is denoted by ${[]}^{0}$ with basis vectors $\{{\mathbf{e}}_{1}^{0},{\mathbf{e}}_{2}^{0}\}$ and origin point $\mathbf{O}$. The target plane $\mathbf{T}$ is expressed by a plane equation ${\mathbf{C}}_{\mathbf{0}}\xb7{[\mathbf{X}]}^{0}-{d}_{0}=0$, where ${\mathbf{C}}_{\mathbf{0}}={[\begin{array}{cc}{c}_{1}& {c}_{2}\end{array}]}^{T}$, ${d}_{0}$, and ${[\mathbf{X}]}^{0}$ are the normal vector, shift and spatial variables, respectively. Although the laws of reflection (or refraction) are invariant with coordinate system, solution of freeform problem requires that the relations derived from these laws can be expressed in a coordinate system appropriate to the geometry of the given problem. Hence, we define a local coordinate system, by which the $i$’th local coordinate of any vector $\mathbf{V}$ is denoted by ${[\mathbf{V}]}^{i}$ with corresponding orthogonal basis $\{{\mathbf{e}}_{1}^{i},{\mathbf{e}}_{2}^{i}\}$ and local origin ${\mathbf{O}}^{i}$. Coordinate transformation from $i$’th to $j$’th coordinate system can be conducted via a sequence of operators, ${\mathscr{H}}_{i}^{j}={\mathcal{T}}^{j}{\mathcal{R}}_{i}^{j}$, where ${\mathcal{R}}_{i}^{j}$ and ${\mathcal{T}}^{j}$ represent the rotation of basis vectors and translation: ${\mathcal{T}}^{j}(\mathbf{V})={[\mathbf{V}]}^{j}+{[{\mathbf{O}}^{i}-{\mathbf{O}}^{j}]}^{j}$, respectively. By means of subsequent intermediate local coordinates, the light source $\mathbf{O}$ and the target plane $\mathbf{T}$ can be linked as ${[\mathbf{O}]}^{i}={\mathscr{H}}_{0}^{i}{[\mathbf{O}]}^{0}$ and ${\mathbf{C}}_{0}\xb7{\mathscr{H}}_{i}^{0}{[\mathbf{X}]}^{i}={d}_{0}$, respectively.

The scenario of our freeform construction is schematically in Fig. 1. The segmental freeform surface to be constructed on $i$’th local coordinate system can be analytically presented as a function ${[\mathbf{A}]}^{i}=[u,f(u)]$ with surface normal as

As shown in Fig. 1(a), while a ray from point source $\mathbf{O}$ is hitting on the point $\mathbf{A}$ at the freeform reflector, the ray is guaranteed to be transported toward the target plane $\mathbf{T}$. The refraction law in vector form is

## (3)

$$\frac{{n}_{2}}{{n}_{1}}\frac{{[\overrightarrow{\mathbf{XA}}]}^{i}}{\Vert {[\overrightarrow{\mathbf{XA}}]}^{i}\Vert}=\frac{{[\overrightarrow{\mathbf{OA}}]}^{i}}{\Vert {[\overrightarrow{\mathbf{OA}}]}^{i}\Vert}-\left[\right(\frac{{[\overrightarrow{\mathbf{OA}}]}^{i}}{\Vert {[\overrightarrow{\mathbf{OA}}]}^{i}\Vert}\xb7{[\mathbf{N}]}^{i})\phantom{\rule{0ex}{0ex}}-\sqrt{{n}_{2}^{2}-{n}_{1}^{2}+{n}_{1}^{2}{(\frac{{[\overrightarrow{\mathbf{OA}}]}^{i}}{\Vert {[\overrightarrow{\mathbf{OA}}]}^{i}\Vert}\xb7{[\mathbf{N}]}^{i})}^{2}}]{[\mathbf{N}]}^{i},$$^{7}the refraction index can be defined by ${n}_{1}=1.0$ and ${n}_{2}=-1.0$. Here, we size each cell of subdomain pair (${\mathrm{\Omega}}_{\theta}^{i}$ and ${\mathrm{\Omega}}_{x}^{i}$) to allow equal flux transfer. In order to derive the freeform surface via the preceding equations, we must introduce an intermediate variable $u$ in the computational domain (${\mathrm{\Omega}}_{u}^{i}$) to link the computation in both source and target domains (${\mathrm{\Omega}}_{\theta}^{i}$ and ${\mathrm{\Omega}}_{x}^{i}$). By changing the variable, Eq. (3) can be rewritten as

## (4)

$${\int}_{{\mathrm{\Omega}}_{u}^{i}}I(u)\frac{\mathrm{d}{\theta}^{i}}{\mathrm{d}u}\mathrm{d}u={\int}_{{\mathrm{\Omega}}_{u}^{i}}E(u)\frac{\mathrm{d}{x}^{i}}{\mathrm{d}u}\mathrm{d}u.$$Obviously, the integral form of energy conservation, Eq. (1) can be represented by a differential equation in summation:^{7}

## (5)

$${E}_{{\mathrm{\Omega}}_{u}^{i}}\xb7\frac{\mathrm{d}{x}^{i}}{\mathrm{d}u}-{I}_{{\mathrm{\Omega}}_{u}^{i}}\xb7\frac{\mathrm{d}{\theta}^{i}}{\mathrm{d}u}=0\phantom{\rule[-0.0ex]{1em}{0.0ex}}\text{for}\text{\hspace{0.17em}}\text{\hspace{0.17em}}u\in {\mathrm{\Omega}}_{u}^{i}.$$Let $({u}_{h}^{i},0)$ be the coordinate of the hitting point where the incident ray $\overrightarrow{\mathbf{OA}}$ with light source ${[\mathbf{O}]}^{i}=({L}_{x}^{i},{L}_{z}^{i})$ hits on $i$’th local $x$-axis. The $x$-coordinate of the hitting point ${u}_{h}^{i}$ can be expressed as two separate functions of $u$ and ${\theta}^{i}$ by trigonometric formula:

## (6)

$${u}_{h}^{i}({\theta}^{i})={\rho}_{0}\frac{\mathrm{sin}({\theta}^{i})}{\mathrm{cos}({\theta}^{i}+{\varphi}_{0}^{i})},$$## (8)

$$\frac{\mathrm{d}{\theta}^{i}}{\mathrm{d}u}=\frac{\mathrm{d}{\theta}^{i}}{\mathrm{d}{u}_{h}^{i}}\frac{\mathrm{d}{u}_{h}^{i}}{\mathrm{d}u}.$$Similarly, suppose the incident ray $\overrightarrow{\mathbf{OA}}$ is reflected onto the point ${[\mathbf{X}]}^{i}=({x}^{i},{z}^{i})$ at the target plane $\mathbf{T}$. The variables ${x}^{i}$ and $u$ related by the law of reflection can be expressed as

## (9)

$$\frac{\Vert {[\overrightarrow{\mathbf{OA}}]}^{i}\Vert}{\Vert {[\overrightarrow{\mathbf{XA}}]}^{i}\Vert}(u-{x}^{i})=2({[\overrightarrow{\mathbf{OA}}]}^{i}\xb7{[\mathbf{N}]}^{i})[-{f}^{\prime}(u)]-(u-{L}_{x}^{i}),$$## (10)

$$\frac{\Vert {[\overrightarrow{\mathbf{OA}}]}^{i}\Vert}{\Vert {[\overrightarrow{\mathbf{XA}}]}^{i}\Vert}[f(u)-{z}^{i}]=-2({[\overrightarrow{\mathbf{OA}}]}^{i}\xb7{[\mathbf{N}]}^{i})[-{f}^{\prime}(u)]-[f(u)-{L}_{z}^{i}].$$Dividing Eqs. (9) by (10), the equation that defines a local freeform surface can be obtained:

## (11)

$$\frac{u-{x}^{i}}{f(u)-{z}^{i}}=\frac{2({[\overrightarrow{\mathbf{OA}}]}^{i}\xb7{[\mathbf{N}]}^{i})-(u-{L}_{x}^{i})}{-2({[\overrightarrow{\mathbf{OA}}]}^{i}\xb7{[\mathbf{N}]}^{i})[-{f}^{\prime}(u)]-{[f(u)-{L}_{z}^{i}]}^{i}}.$$If we fix the source and target plane position, the freeform surface defined in Eq. (11) can also be considered as a parameterized surface where the parameter is the $x$-coordinate of the $i$’th coordinate system. The freeform parameters depend on the external geometry at which each optical route O-A-X constitutes an epipolar plane. In order to transform the underdetermined Eq. (11) to a determined problem, we replace ${z}^{i}$ by ${x}^{i}$ through the equation ${\mathbf{C}}_{0}\xb7{\mathscr{H}}_{i}^{0}{[\mathbf{X}]}^{i}={d}_{0}$:

## (12)

$${x}^{i}(u)=\frac{{C}_{0}^{i}{f}^{\prime}(u)+{C}_{1}^{i}f(u)+{C}_{2}^{i}u+{C}_{3}^{i}}{{C}_{4}^{i}{f}^{\prime}(u)+{C}_{5}^{i}f(u)+{C}_{6}^{i}u+{C}_{7}^{i}}+R(u,f),$$## (13)

$${\alpha}_{1}^{i}{f}^{\prime \prime}+{\alpha}_{2}^{i}{f}^{\prime}+{\alpha}_{3}^{i}f+{\alpha}_{4}^{i}=0,$$The construction of freeform surface by numerically differential formulation has an advantage in that we can regularize the smoothness tolerance by setting the partition number (${\mathrm{\Omega}}_{\theta}^{i}$ and ${\mathrm{\Omega}}_{x}^{i}$) in a tradeoff between computational efficiency and surface smoothness. In order to ensure the numerically estimated accuracy of freeform surface, sufficient partition cells shall be numbered in both domains. Here, we impose the edge rays ${R}^{i-1}$ and ${R}^{i}$ in the polar cell toward the corresponding boundary $[{b}^{i-1},{b}^{i}]$ in the Cartesian cell, so the radiant flux is fully transferred without stray loss. The edge ray,${R}^{i}$, reflected by the freeform surface is propagating toward the corresponding point ${b}^{i}$. We shall set a reasonable right end point of the intermediate cell ${\mathrm{\Omega}}_{u}^{i}$, by which the freeform possesses convexity in the local domain. To solve the nonlinear PDE in Eq. (13), we employ the vanished moment method proposed by Feng and Neilan and use Newton’s method to find the adequate root function. The linearized equation at each iteration step can be solved by the finite element method. Here, we skip some straightforward but tedious expansion and rearrangement. Figure 2 demonstrates the key steps of the freeform construction process. The algorithm of freeform construction in this session can be summarized as following:

## Algorithm

Freeform surface via local reconstruction (FS-LRA).

Input:Ωθ: Domain of light source: Ωx: Domain of target; O: Light source; T; Target plane equation; Global coordinate frame []0; O1: initial point of freeform surface: |

Output: f(x): freeform surface function; |

1. Partition the source and target domain into Ωθi and Ωxi for i=1,2,…,n; |

2. Setup position and orientation of the freeform (i.e., determine global and the first local coordinate frame []1): |

3. fori=1 to ndo |

4. Determine Ωui; |

5. Solve equation (13) on Ωui; |

6. Setup the i+1-th coordinate frame []i+1; |

7. Apply transformation operator Hii+1 to the light source and target plane; |

8. end for |

## 3.

## Freeform Surface via Local Reconstruction Algorithm: Error Evaluation

In order to validate the proposed freeform scheme, we first constructed a surface by the previous algorithm and compared it with a known parabolic reflector; both are aimed to deliver a collimation beam as shown in Fig. 3. A point source was placed at the focal point $(x,y)=(0,10)$ of an ideal parabolic reflector $y=-{x}^{2}/12+13$. Here, we limit the emitted angle from the point source within the range of ${\mathrm{\Omega}}_{\theta}=[-119\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{deg},119\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{deg}]$, and the radiant flux reflected by the ideal parabolic reflector was directed toward the given target plane ${\mathrm{\Omega}}_{x}=[-10,10]$ with a uniform and parallel propagation. Figure 3(a) exhibits the surface comparison between the curve constructed via proposed algorithm and ideal parabolic curve. The surface error ${e}_{s}=\phantom{\rule{0ex}{0ex}}{\mathrm{max}}_{x\in {\mathrm{\Omega}}_{x}}\Vert {S}_{\mathrm{FF}}(x)-{S}_{P}(x)\Vert $, maximum deviation between the constructed freeform surface (${S}_{FF}$) and the ideal parabolic curve (${S}_{P}$), is on the order of ${10}^{-4}$ to ${10}^{-6}$, which is extremely small and quadratically inverse to the number of partition cells. The results are in agreement with the theoretical prediction that the error comes from the finite element method with Hermite element.^{16}17.^{–}^{18}

## 4.

## Construction Algorithm: Two Design Types

## 4.1.

### Uniform Illuminance

After validating our freeform algorithm with a comparative parabolic surface by forward testing, we now aim to achieve a uniform illuminance with a geometric correspondence by the freeform surface via local reconstruction algorithm (FS-LRA). As shown in Fig. 4, a 20-m-wide target plane was placed beneath a point source with 6-m distance. The emitted angle from the point source was confined for reflection and refraction within ${\mathrm{\Omega}}_{\theta}=[-30\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{deg},30\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{deg}]$. With a mere 80 partition cells and 960 partition cells in both source-target subdomains, a freeform reflector and refractor can be generated in about 1.1 and 10 s, with 2.3 GHz Intel® dual core computer. The study employed a ray tracing simulation on the constructed freeform surface with a commercial tool LightTools®, in which light rays from the source were incidental to the target plane. The root-mean-square (RMS) nonuniformity is merely 1.05% (reflector) and 2.10% (refractor), respectively. It is noted that the accuracy or partition number of FS-LRA is highly dependent on the tolerance of the freeform surface. For uniform illuminance request in both reflect/refractive cases, a monotonic curve is accessible with the prescribed simple target surface. Because the mathematical form in the refraction case is more complex than that in reflection case, 960 partitions are required for convergence to achieve the design uniform illumination pattern, but only 80 partitions are required for the constructed reflector to achieve the same uniformity.

## 4.2.

### Linear Varying Illuminance

In this example, we extended our technique to create a prescribed color illumination by a pair of independent sources with different colors. Figure 5 shows the schematic diagrams of the system layout for both reflect/refractive cases. With a 20-m target plane for both cases, we aimed to create a linearly varying international commission on illumination (CIE) chromatic $(x,y)$ along the $x$-direction associated with two color sources. The distance between source and target plane is 10 and 70 m, respectively. The light sources leave the center at 10 m and 15 m, respectively. The different configuration conditions of reflection and refraction design are to achieve desired illuminance performance. As shown in Figs. 5(c) and 5(d), the illuminances feature a high linearity along the $x$-direction. There are only about 5.58% and 1.52% RMS nonuniformity with the perfect chromatic line for reflection and refraction cases, respectively.

## 5.

## Conclusion

In this article, we proposed a freeform model based on local energy conservation associated with a series of coordinate transformation. The major advantage of this method is to exploit the well-proven existence and continuity of an elliptic nonlinear differential equation of the Monge-Ampere type. We successfully formulated the freeform design into a finite element paradigm with Newton’s iteration. The key maneuver is to unify the variables ($\theta $ and $x$) in both domains (${\mathrm{\Omega}}_{\theta}^{i}$ and ${\mathrm{\Omega}}_{x}^{i}$) into a common intermediate variable $u$ in the computational domain (${\mathrm{\Omega}}_{u}^{i}$). After imposing the small local planar approximation, which is absolutely valid with a sufficient number of partition cells, a complete freeform reconstruction algorithm can be developed. In addition, the segmental differential formulation is capable of taking numerical error into account in each step. The proposed surface parameterization guarantees the surface smoothness with no restriction on the desired target distribution, collective solid angle, or source distribution, respectively. This technique still leaves many opportunities open and clearly more research must be carried out to explore its potential in full. First, at this moment, the issue was tackled in a one-dimensional case; complete treatment in consideration of twist deformation in freeform surface is underway. Second, point source approximation is another subject to be addressed when the freeform structure was placed in proximity of the light source such as light emitting diodes applications. The preliminary results presented in this article, however, indicate that parameterizing the freeform surface via local freeform PDE may create another route in the field of freeform optics, as it has the potential to take freeform design into many nonimaging applications.

## Acknowledgments

This work was financial supported by the National Science Council of Taiwan government under contract NSC101-2622-E-009-001-CC2 and 100-2115-M-009-001. The first author and second author contributed to this work equally.

## References

## Biography

**Yu-Lin Tsai** is a PhD candidate in the Department of Photonics at National Chiao Tung University. He received his MA degree in the Institute of Mathematical Modeling and Scientific Computing from National Chiao Tung University, Hsinchu, Taiwan, in 2010. His research interests include nonimaging optical design and freeform design problems.

**Ming-Chen Chiang** is an engineer in the Department of Wireless Product of Universal Scientific Industrial. He received his MA degrees both in IMMSC (Institute of Mathematical Modeling and Scientific Computing) and Institute of Communication Engineering from National Chiao Tung University, Hsinchu, Taiwan, in 2011. His interests include the freeform design problems and communication channel analysis.

**Ray Chang** is a master student in the Institute of Display at National Chiao Tung University. He received his MA degree from the Institute of Display from National Chiao Tung University, Hsinchu, Taiwan, in 2013. His research includes nonimaging freeform design problems.

**Chung-Hao Tien** is an associate professor in the Department of Photonics at the National Chiao Tung University (NCTU) where he has been a faculty member since 2004. He received his PhD degree in the Institute of Electro-Optical Engineering from NCTU in 2004. His research interests include computational imaging, free form nonimaging, display and\ lighting optics.

**Chin-Tien Wu** is an associate professor in the Department of Applied Mathematics at National Chiao Tung University. He received his PhD degree in Scientific Computation and Mathematical Modeling from the University of Maryland in 2003. His research interests include scientific computing and freeform design problems.