Stress mirror polishing for future large lightweight mirrors: design using shape optimization

This study proposes a new way to manufacture large lightweight aspherics for space telescopes using Stress Mirror Polishing (SMP). This technique is well known to allow reaching high quality optical surfaces in a minimum time period, thanks to a spherical full-size polishing tool. To obtain the correct surface' aspheric shape, it is necessary to define precisely the thickness distribution of the mirror to be deformed, according to the manufacturing parameters. We first introduce active optics and Stress Mirror Polishing techniques, and then, we describe the process to obtain the appropriate thickness mirror distribution allowing to generate the required aspheric shape during polishing phase. Shape optimization procedure using PYTHON programing and NASTRAN optimization solver using Finite Element Model (FEM) is developed and discussed in order to assist this process. The main result of this paper is the ability of the shape optimization process to support SMP technique to generate a peculiar aspherical shape from a spherical optical surface thanks to a thickness distribution reshaping. This paper is primarily focused on a theoretical framework with numerical simulations as the first step before the manufacturing of a demonstrator. This two-steps approach was successfully used for previous projects.


Introduction
Stress Mirror Polishing (SMP) is a manufacturing method allowing to obtain the required shape of optical surfaces using deformation and elasticity properties of materials. In this approach, the mirror is deformed from a sphere, as we want to obtain a quasi parabola, to the desired inverse-shape, spherically polished using full-size tools, before to take the correct form when the applied constraints are released. This technique has been developed for forty years on several type of mirrors such as Keck segments [1] or the Spectro-Polarimetric High-contrast Exoplanet Reseach (SPHERE) instrument on the Very Large Telescope (VLT) [2]. Contrary to the classical aspherical polishing process where small size tools introduce some high-order spatial frequencies errors on the optical surface, the SMP process uses large tool which gives an extreme surface quality due to the continuous deformation and polishing behaviors.
Moreover, due to the large area of the polishing tools, SMP process [3] allows a considerable saving in time of manufacturing and is perfectly suited to produce a large number of primary mirror segments in a short time for future space and ground-based observatories such as Thirty Meter Telescope [4,5] and LUVOIR [6] projects. SMP techniques have also a significant benefit for the manufacturing of large lightweight aspheric primary mirrors for space telescopes, especially as they can be applied on a thin shell mirror and combined with an active support [7] or a rigid assembly using a sandwich structure with two facesheets and a Carbon Fiber Reinforced Composite (CFRC) honeycomb [8]. This could be perfectly suited for manufacturing of ultra high surface quality primary mirrors, as the ones foreseen in future high-contrast imaging missions for exoplanets detection and characterization, or simply to have a gain in term of manufacturing time for more classical mirrors.
We investigated the possibility to produce such a lightweight aspheric primary mirror for space telescopes using classical elasticity theory combined with new tools. They are provided by Finite Element Analysis (FEA), such as "shape optimisation", to go further the classical use described by Izazaga-Pérez [9]. In the first section, we introduce rotational symmetry aspherics, known as conicoids, and the method to produce them using Active Optics and SMP more specifically by giving an example for a 3 th -order spherical Zernike aberration. We will discuss the definition of shape optimization and the strategy adopted to implement it. Then, an algorithm has been developed allowing to optimize the deformation obtained, using data provided by Finite Element Model (FEM) outputs which has been implemented in Python language.

Definition
Primary mirrors are one of the most sensitive part of space telescopes. The increasing size of diameter allows collecting a large amount of photons and detecting fainter objects, while the surface quality is crucial to obtain the best images of stars, planets or exoplanets [10]. For these reasons, large diameter and high surface quality are requested on most space and ground-based future observatories. Nowadays most of the telescope optical designs use primary parabolic mirrors which are more convenient to focus rays from infinite objects without introducing any additional complexity. The sphere and the parabola are two peculiar solutions of what is commonly called conicoids.
Let us define Eq.(1) of a conicoids aspheric optical surface with rotational symmetry by : where R c is radius of curvature of aspherics with a conic constant C and r is the radial coordinate from the center to the semi-diameter a. Note that if C = 0 the optical surface is a sphere and a parabola if C = −1.
By developing expression in expansion series, we obtain for the first terms Eq.(2): Retrieving the asphere z a from the sphere z s with radii of curvature R a and R s respectively allows us to determine the deformation z d to be applied during the polishing process as: Replacing Eq.(2) in Eq.(3) we obtain (Eq.(4)): R s is defined as the first parameter of the substrate dimensions. R a is defined in order to minimize the removal volume defined by Eq.(5) as described by Unti [11]: where ∆z is the difference in abscissa between the sphere and aspheric (Eq.(6)): with is the shifted distance from the origin as shown in the Fig. 1.
Definitions of the best-fit sphere (BFS) minimizing the volume removal to obtain the targeted aspheric surface.
Unti [11] demonstrated for a perfect parabolic mirror with a conic constant equal to -1, the best-sphere radius of curvature minimizing the volume removal (Eq. (7)) is equal to: The best-fit sphere radius can also be determined by minimizing the Root-Mean Square (RMS) deviation between the targeted optical surface and the best-fit sphere as used in the optical software ZEMAX [12]. The deformation z d is determined by substituting the radius of the sphere R s (Eq.(8)) in the Eq.(3) by the best-fit sphere radius R BF S .

Active Optics and 3 r d order spherical aberration generation
Stress mirror polishing developed by Lemaître [13] is derived from Schmidt works to aspherize his refractive correctors using full-size spherical tool to avoid any local retouch [14]. This technique relies on Active Optic methods using polishing pressure as an uniform load distribution.
Active optics method based on elasticity's theory of materials is usually used to correct the Wave-Front Error (WFE) of an optical system by applying the adjusted deformation on one or several surfaces. Single or combined deformation modes can be generated in this purpose with one or several actuators. For example SPHERE, the Very Large Telescope (VLT) instrument dedicated to exoplanets detection, is assisted by this method as a means of compensating small deformations due to external constrains to the instrument, such as gravity load or temperature gradient, at long timescale [15].
For years, they have also been deployed for many applications such as Variable Curvature Mirrors (VCMs) [16] implemented in the delay line system of the the Very Large Telescope Interferometer (VLTI) in order to compensate optical path differences by changing focii of mirrors also called zoom mirrors [17]. Usually deformation modes are generated with a combination of pressure, forces and bending moments. However pure aspheric shapes are mostly difficult to execute on classical mirrors with quasi-constant thickness distribution. For this particular reason Variable Thickness Distribution (VTD) classes have been developed by Lemaître [16] to produce better first-order modes, by adapting the thickness distribution to the required deformation shape, the load distribution and the boundary conditions. In our particular case, VTD is the starting point to generate aspherics using stress polishing. Third-order spherical aberration (SA3) can be generated by applying a uniform polishing pressure P p on the optical surface with a pushing central force F c in reaction, as developed by Lemaître [18], while the mirror's edge is maintained with a device discussed in the following paragraph. In that specific configuration, not only pure 3 r d order spherical aberration is produced but also curvature and some other high-order spherical aberrations which have to be removed.
The dimensionless thickness distribution T D S A3 which produces 3 r d -order spherical aberration is described [18] by the following Eq.(9): with ρ = r/a ∈ [0; 1] ; a the mirror's semi-diameter and ν the Poisson's ratio. It is important to note that, in this solution, the thickness is zero at the mirror's edge and infinite at the center as shows in Fig. 2. In practice, primary mirrors have to be polished with an acceptable thickness at the edge according to the large mirror's diameter and finite thickness at the center.
To obtain the real variable thickness distribution T S A3 (Eq.(10)), it is necessary to multiply the normalized dimensionless thickness distribution T D S A3 by a scaling factor T 0 [18] and taking into account the sphere deviation and the edge thickness e. Scaling factor is depending on intrinsic parameters of the mirror's material such as Poisson's ratio ν and Young's modulus E but also external loads q and the sagitta of SA3 that expected to be generated A S A3 : with

Shape optimization and Pure 3 r d -order spherical aberration
Shape optimization is a very efficient tool as long as you know how to use it properly. Botkin [19] has defined the basis of such optimization to go beyond mass reduction problem using mainly size optimization of a model in which the geometry remains unchanged. Shape optimization, has since been developed [20] and is currently used in a large number of various applications such as fluid-structure interaction to find the optimal shape of a plane wing by minimizing drag [21]. Very few examples are referenced in literature as noticed Lund [22] for the fluid-structure coupling especially but the same could be said for field in optics. One should notice that it is not uncommon to encounter an amalgam between shape optimization and topology optimization. From an objective function, shape optimization refers to the relocation of points belonging to a given shape while topology optimization determines if a model element has to be removed or not based on the level of its involvement in the model rigidity.
The shape optimization is based on the variation of design variables depending on an objective function, as mentioned before, which is usually a scalar subject to design constraints. The objective function can be a combination of design responses to provide a single final response. The main issue of using shape optimization is to generate the correct shape basis vectors. Indeed, each selected model gridpoints, also called nodes, will move along these basis vectors directions. NASTRAN software integrates this functionality as described by Kodiyalam [23] but can also be generated by a self-programmed code. In this case, the model meshing should be adapted for both having enough space to relocate nodes and sufficient number N of nodes to converge. Because of the large amount of design variables, the optimizer is not able to perform all the possible options but it can use design sensitivity analysis in order to determine the optimal set of design variables to reach the goal.
Design sensitivity computes all the rates of change [24], also called partial derivatives coefficients λ i j , for each i-th design variables v i to anticipate the influence on the j-th model response r j given by Eq. (11):

Problem Formulation
In our case we want to find the optimal variable thickness distribution with the aim of producing a pure 3 r d -order spherical deformation defining by Zernike "Peak-to-Valley"(PtV) polynomials [25] as: with δ n,m = is the root-mean square value of the 3 r d -order spherical aberration linked to the PtV value using δ n,m coefficient.
In this context, the goal of the optimization is to reduce the deviation of the deformation from the theoretical SA3 surface, given in Eq. (12), by minimizing the threshold β(∆) with ∆ = δ 1 , δ 2 , . . . , δ N and subject to the design constraints ζ i (Eq.(13)): where δ i (Eq. (14)) is the deviation between model deformation responses z i and the SA3 theoretical shape for each gridpoints i of the optical surface z i S A3 , displayed in Fig. 3, defines as:  Fig. 4 is showing the Beta optimization principle and the design constraints ζ i which impose δ i values to be lower than the threshold β(∆) while it is reducing.

Shape optimization algorithm
Shape optimization method has to start with a model close to the optimal form such as the thickness distribution described in Eq.(9) to avoid large gridpoints displacements and local minima. Finite Element Analysis (FEA) is operated for the first calculation with the initial model parameters to get the input data: optical surface displacements as design responses and variable thickness distribution as design variables layered in Fig. 5.
Thereafter, the input data files are used to generate the bulk data entries with the aid of a self-programed PYTHON algorithm that will be embedded in the optimization input file. This program allows to calculate the maximum difference between the theoretical desired shape and the FEM displacements δ max = Max(δ i ) in order to define the minimum threshold to implement. Since the difference δ max is small, a multiplying factor c 3 has to be applied to assure an objective function β ob j large enough in order to prevent small values depending on the model unit (Eq.(15)): As mentioned, the threshold has to be higher than the maximum difference δ max and has to be used as design constraint for each optical surface gridpoints. Eq.(16) describes the design constraint ζ i implemented: where c 1 and c 2 are the two coefficients used to manage the shape optimization. c 2 is slightly higher than the difference δ max to avoid a constraint violation before running the optimization algorithm. The design constraint sensitivity is handled by the scaling factor c 1 depending on the design model.
As regards design variables selected on the thickness distribution boundary, they are all model gridpoints in this particular numerical simulation. They are compiled in a DVGRIDs entry defining the direction and magnitude for a given change in a design variable.
Once the optimization is terminated, particularly for the first one, due to the deformation of the elements induced by the gridpoints displacements, it is indispensable to refine the Finite Element Model (FEM) to get more accurate results. With a projection on the Zernike orthonormal basis, the optical surface is decomposed into a linear combination of weighted orthogonal polynomials.
Due to the ill-posed problem such as non-adapted optimization coefficients or an initial design very far from the target, it is rare to converge in one step of optimization. In most cases, the optimization process needs a certain number of iterations, but usually less than a dozen. The optimization range will be performed on the entire mirror's surface a but it can be narrowed to the interested area.

Initial design
This first study is limited to a two-dimensional axisymmetric problem using Finite Element Analysis (FEA) software. In this approach, CTRIAX6 elements are required to make the calculation possible. However, CTRIAX6 element defines isoparametric and axisymmetric triangular elements with mid-side gridpoints which do not enable a large deformation due to the elements distorsion. Therefore, remeshing the model at each step is critical to pursue the optimization process. The initial design of the "tulip-form" mirror based on the T S A3 variable thickness distribution is shown in Fig. 6  The preliminary parameters introduced for the initial design are listed in Table 1. We use a Zerodur vitro-ceramic substrate from Schott company as primary mirror material due to its low coefficient of thermal expansion perfectly adapted to space mirror stability. The mechanical properties of Zerodur are the following: Young's modulus= 91000 MPa, shear modulus=37295.08 MPa and density Poisson's ratio=0. 22.
The outer mechanical diameter for SMP a and the final outer mechanical diameter r f are 1080 mm and 1030 mm respectively as shown in the Figure 6. The clear aperture has been chosen for more convenience and set at 1000 mm for a quasi-parabolic mirror with a conic constant around -0.967. The planned central hole is 277 mm diameter which will be machined after the polishing process. A starting mirror's edge thickness has been established at 2.9 mm in order to be able to support the mirror during the deformation. The total mass of the mirror before the central hole drilling is 25.29 kg. A compromise between large gridpoint displacements and good accuracy in the first calculation was reached with around 500 CTRIAX6 elements and 1200 gridpoints in the finite element model. These numbers should be increased when the objective function decreases in order to converge more effectively. A combination of loads is applied in the FEA including a nominal polishing pressure of 31.5g/cm 2 on the optical surface, associated to a central pressure of 9.51M Pa at the "tulip-form" basis to generate pure 3 r d -order spherical aberration. Entire model is placed under 9.81m.s −2 gravity as shown in Fig. 7 and stays fixed during all the optimization process.
An articulated and movable border along the radial axis is disposed at the top edge of the outer mechanical diameter. The purpose of this boundary condition is to avoid introducing radial moments M r and radial tensile forces N r and to allow a higher sagitta of mirror's deflection and a pure 3 r d -order spherical aberration.

Validity of the FEM
Finite Element Model (FEM), and particularly Finite Element Analysis (FEA), is commonly used to study new opto-mechanical configuration before entering in the manufacturing phase of a mirror. In this approach the FEA outputs are cross-checked with some demonstrator experimental results in a second phase. Most of the time a very high correlation is noticed, however the last phase is mandatory to perfectly adjust the model. This procedure has been successfully used to design, simulate and polish toric mirrors of the SPHERE instrument for the VLT [26,27] or to modify the design of the warping harnesses in order to gain in performance on sky [15]. The same approach led to the study and development of active optics dedicated to space applications through the MADRAS deformable mirror [28,29] and is currently applied for the design and manufacturing of off-axis parabola for the WFIRST mission [30].
This specific background in combining FEA with opto-mechanical fabrication, acquired during previous projects, give us a good confidence in the numerical results presented in the following section with respect to the up-coming prototyping phase.

Results and Discussions
The goal of the optimization is to deform the mirror by producing 14.4 µm RMS of 3 r d -order spherical aberration during polishing process to obtain the targeted aspheric, which is an ellipsoid with a conic constant of -0.96682, for a full aperture diameter of 1080mm. Optimization parameters, c 1 , c 2 and c 3 coefficients have to be updated for each re-meshing step depending on design responses. DVGRIDs data entries are provided in function of design model but the  displacement of each gridpoints is decreasing as the final shape approaches. Because of the axisymmetrical model we are working on, we produce only axisymmetrical modes such as Focus and spherical aberrations. Fig. 8 shows the deformation resulting from the initial calculation before the first optimization. At first sight, the mirror has the shape of the planned 3 r d -order spherical aberration.
By subtracting the theoretical SA3 from the FEA deformation we can notice a large variation notably at the center of the mirror with 9.4 µm PtV surface error (red solid curve in Fig. 9). The global deviation from the theoretical shape along the mirror's radius appears to be 5 th -order spherical aberration (SA5) as shown in Fig. 10 depicting the three first theoretical Zernike spherical polynomials. The edge of the mirror is not moving along the optical axis, Z-axis, due to the adopted boundary condition. The pupil is taking into account the radial displacement of the mirror in order order to have the right Zernike's decomposition.
The initial calculation before the optimization process shows a deformation map without Piston (Fig. 11a) around 1.8 µm RMS, essentially SA5 with 1.4 µm RMS, as noticed previously. A projection of residuals on Zernike modes between SA7 (36 th Zernike polynomials) and 561 th Zernike polynomials indicates a high-order value of 84.1 nm RMS (Fig. 11b) with a maximum deflection at the center. Given the mirror's axisymmetrical geometry and the load condition, it is not necessary to use more than 561 Zernike polynomials to define the optical surface. Indeed, the spatial frequency between two peaks in the residual deformation is much lower in comparison with the radial component of the 561 th Zernike polynomial.    A first optimization has been computed on the entire mechanical radius [0mm ; 540mm] letting gridpoints on the variable thickness distribution free to move along the Z-axis. A significant reduction of the difference between desired theoretical shape and FEA deformation can be noticed leading to a residual of 6 µm PtV (Fig. 9) corresponding to a residuals sag. of 2.1 µm RMS on the entire pupil. Two extrema are present and contribute to the high-order spherical aberrations which can be visualised on the residuals error maps (Fig. 11).
The first error map (Fig. 11c) shows the residuals of the Zernike decomposition when the Piston aberration, which we do not need to consider, is substracted. This value decreases down to 52.9nm RMS by removing the 36 first Zernike polynomials by a numerical subtraction (Fig. 11d). These high-order polynomials need to be minimized as much as possible in the optimization process because of the complexity of applying the high-order removal procedure afterward.
The next optimizations are performed on the reduced pupil [138.5mm ; 540mm] to take into account the central hole and thus release constraints on the optimization process. The deviation from the theoretical SA3 is free to evolve within the central hole while the optimization is focused on reducing the difference in the pupil of interest (Fig. 12). As the threshold has to be decreased to reach a very low value of a few nanometers, the model unit has to be scaled in order to better control the optimization process. In this particular case, it was modified to use micrometers instead of millimeters to reach a better sensitivity on the optical surface coordinates.
By focusing the analysis on the reduced pupil it allows us to observe a substantial improvement in the aim of fitting the SA3 deformation given by the FEM model and the targeted theoretical SA3. However, some difficulties in the optimization process have to be noticed for a radius at distances larger than 400mm from the center. It is illustrated by a second regime (black and blue curves in Fig. 13) retaining any optimization in this radius range. The residual error caps at around 20nm RMS on the drilled optical surface. Indeed, this part of the mirror has a very small thickness, starting from 9.9mm to 2.9mm, which makes the model responses most sensitive to the design variables specially in this area. This leads to consider the optimization on the very limited pupil [138.5 ; 500mm] without the optical surface extension.  In this instance, a last optimization is performed on the final pupil [138.5;500mm] (green curve in Fig. 13) adding a supplementary degree of freedom in the optimization for gridpoints located in the range [500;540mm]. Naturally, it leads to a peak in the non-optimized area to ensure a optimal result in the optimized one. Finally, the deformation map of the optical surface (Fig. 14a) shows a residual error, without piston, about 16.5 nm RMS in the reduced pupil including 14 nm RMS of high-order spherical modes after the 36 first Zernike's polynomials (Fig. 14b). After verification by changing the model units twice and refining it with very small elements (150 µm between two gridpoints), it was concluded that these residuals are effectively attributable to the FEA deformation and not to the accuracy of the FEA sampling (Fig. 15).
Indeed, by taking a minimum mesh size (MMS) (distance between two gridpoints) between 0.15mm and 2.5mm, the residual error without the 36 first Zernike polynomials (R36) ranges from 14.03 nm RMS to 14.14 nm RMS. This numerical difference of 1 Angstrom is comparable to the dimension of a single atom. In our case, a MMS of 1.75mm has been chosen to make a trade-off between precision of the results and reduction of computational time and load. In this manner, a duration of 30 minutes for each iteration in the optimization process has been performed with a total of five in our particular case.
The optimized final variable thickness distribution in Fig. 16 (red curve) is different from the initial theoretical contour for many reasons. The first reason lies in the VTD definition which is valid for a zero thickness at the mirror's edge. In our case, the mirror needs to have a sufficient thickness at the edge to be maintained and a finite height at the center to respect manufacturing and loads specifications, leading to shift and cut the theoretical curve. Another reason explaining this VTD deviation is the useful surface used for the optimization. At the beginning, the entire optical surface has been exploited for the shape optimization but it has been reduced to the clear aperture to release some limitations and also accelerate the minimization process. In addition, polishing and central pressures have been slightly adjusted to cancel focus and residual SA3 completely. In fact, the variation of spherical aberrations according to the loads is more sensitive for lower than higher Zernike modes. In this case, a minor modification of loads impacts very low orders leaving spherical residuals from SA5 not affected.
This adjustment can also be used, during manufacturing by SMP as an iterative process, to help converging towards former specifications. After the first polishing run, the optical part is usually slightly out of spec on very low orders. Adjusting polishing pressure and loads allow to compensate for over or under deformation of the optical part. A second and a third run are generally sufficient to reach the correct shape letting high spherical modes not impacted.
The maximal Von Mises stress within the mirror during the deformation is equal to 8.09 MPa (Fig. 17) which is acceptable for standard use in this particular case (limited to 14 MPa). The maximal stress is located at the "tulip-form" basis owing to the important pressure applied in a relatively small area.

Conclusions
Thanks to this approach, combining SMP and shape optimization, we show a new manner to produce a conicoid, rotationaly symmetric aspheric, from a 2D finite element model by imposing pre-defined loads configuration and boundary conditions. The shape optimization method will assist the SMP process to reach high quality and very precise surface shapes in a minimum of iterations. In this study, DVGRIDs data entry have been defined in a basic way but they can be described with more degrees of freedom or even as a combination of orthogonal basis such as Zernike polynomials. It is also possible to define other theoretical symmetrical or non-symmetrical shapes to produce more complex optics with a combination of high-order modes.
The main following step consists of polishing the thin mirror using the shape optimization method developed in this paper. The last step of the development will associate this large (>1m class) and thin (around 20mm) mirror to a rigid lightweight support in order to obtain very stiff reflective optics for space applications. In fact, sandwiches structures offer higher stiffness for lower mass ratio compared to other lightweighted designs [8]. In this case, the parallelization of the SMP technique using shape optimization and the rigid support manufacturing will participate to a considerable gain of fabrication time. Each part can be produced simultaneously on two different manufacturing procedures to be assembled in a final stage.