Laurent Chougrani
email: laurentchougranil.chougrani@poly-shape.com

Jean-Philippe Pernot
email: jean-philippe.pernot@ensam.eu

Philippe Véron
email: philippe.veron@ensam.eu

Stéphane Abed
email: s.abed@poly-shape.com

Parts internal structure definition using non-uniform patterned lattice optimization for mass reduction in additive manufacturing

Keywords: Additive manufacturing, 3D modelling, Shape optimization, Lattice structures, Finite element analysis, Mass reduction, Bio-mimicry, Topology optimization, Variable neighborhood search

Today, being able to generate and produce shapes that fit mechanical and functional requirements and having as low as possible mass is crucial for aerospace and automotive applications. Besides, the rise of new additive manufacturing technologies has widened the possibilities for designing and producing complex shapes and internal structures. However, current models, methods and tools still represent a limitation to that new horizon of printable shapes. This paper addresses the way internal lattice structures can be generated and optimized to reduce the mass of a product. A new framework is introduced that allows the modeling and optimization of non-uniform patterned lattice structures. Using non-uniform structures, additional degrees of freedom are introduced and allow the definition of a wide variety of shapes which can better fit the requirements. First, a non-uniform patterned lattice structure is generated using the results of an initial finite element analysis. This initial structure is then optimized while iteratively removing the beams considered as useless with respect to a user-specified mechanical criteria. At each iteration, the lattice structure is sent to a finite element solver that returns the von Mises stress map used to drive the simplification process. Here, the simulations are performed on the wireframe lattice structures to speed up the optimization loops. Once this process is completed, the final structure is no longer fully patterned, but it is re-organized to reduce the mass while satisfying the mechanical criteria. This approach is illustrated with examples coming from our prototype software.

Introduction and motivations

Additive manufacturing (AM) is a relatively new but yet wildly developing industry which changes drastically the way of designing parts. Classical manufacturing processes like milling or machining tend to use physical tools to remove material, extracting the desired structure from it. Differently, AM processes, especially Laser Beam Melting (LBM or DMLS) and Electron Beam Melting (EBM), tend to start from barely nothing, e.g., metal powders, and build the desired structure out of it, using energy beams (nonmaterial tools) to locally fuse and solidify the powder [START_REF] Gebhardt | Understanding additive manufacturing[END_REF]. Consequently, those new technologies introduce a totally different way of designing, but also new production constraints and limitations as well as new terminologies. Thus, the notion of lattice structures has become popular, and refer numeric methods to produce lightweight products. To reduce the mass of an object, two main ways have been investigated: lattice structure and topology optimization. Both technics lead to complex shapes that can be hard or impossible to manufacture without additive manufacturing. Today, there already exist examples of industrial success stories involving mass reduction through AM. For example, the Thales/ ESA's deployment mechanism structure for which the mass has been reduced by 5, the costs reduced by 4, the number of components reduced by 10 thanks to a combination of topology optimization and lattice insertion (Fig. 1). Costs for such parts is high, it is still a niche market, essentially for aerospace, automotive and energy industrial fields. Moreover, as the final geometry might be convoluted, the computation time for optimization and the number of iterations can be very high.

In this paper, a new way to generate and optimize lattice structures is proposed. The contribution is threefold:

(1) a non-uniform lattice structure generation technique is introduced, based on an analogy to crystallography; (2) a new optimization framework is developed to remove lattice beams according to mechanical criteria, thus decreasing the mass of the structure; (3) a new stress-based criterion is used to compare different results from the proposed algorithm. The resulting structure satisfies to user-specified requirements such as mass reduction and rigidity.

Even though researchers and industrials have investigated this research field, there still exist some limitations on the level of freedom the designers have when treating complex parts. The limitations of the existing approaches are discussed in Sect. 2. The proposed framework is introduced in Sect. 3 and the different steps are detailed in Sect. [START_REF] Zhou | The simp-srv method for stiffness topology optimization of continuum structures[END_REF]. The implementation and the results are presented in Sect. 5. Section 6 concludes this paper and highlights future perspectives of this work.

Related work

In recent years, different approaches have been introduced to define rigid, light and printable structures. They can be classified according to three main categories presented in the following sections.

Topology optimization

Continuous topology is a way of modifying structures to get lightweight parts while optimizing material layout within a given design space with the goal of maximizing the performance of the system [START_REF] Duysinx | Topology optimization of continuum structures with local stress constraints[END_REF]. Unfortunately, topology optimization can be time consuming when the required results needs to be of a high resolution (i.e., high number of voxels). Some other limitations can be discussed depending on the method. For instance, SIMP methods [START_REF] Zhou | The simp-srv method for stiffness topology optimization of continuum structures[END_REF] may lead to noisy structures (Fig. 14a) and is the most used method in the industry due to its good computation time performances. Level set and/or BESO/ESO methods [START_REF] Allaire | Topology and geometry optimization of elastic structures by exact deformation of simplicial mesh[END_REF][START_REF] Huang | A further review of ESO type methods for topology optimization[END_REF] can significantly reduce noisy results but need a re-meshing step that can be time consuming in the long run. Topology optimization is also quite a sensitive analysis [START_REF] Deng | Topology optimization of total femur structure: application of parametrized level set method under geometric constraints[END_REF]. However, those limitations can be dealt with through seeing that topology results could be interpreted as inter-connected structures, plates or beams. This observation combined with AM capacities has also led to lattice structure studies [START_REF] Wang | A unit cell approach for lightweight structure and compliant mechanism[END_REF]. In this paper, advantages of both technics are combined in a unified framework which optimizes both the topology of the object and a non-uniform patterned lattice structure inserted in it. Such structures (inter-connected beams) can be found in nature in both organic and inorganic structures. Since AM allows to generate such structures with much more resolution than any other manufacturing process, it opens the field of biomimicry, i.e bio-inspired designs.

Bio-mimicry

Nature has always been a source of esthetic experience for the mankind and of inspiration for artists and designers [9]. Actually, natural shapes created by natural forces resemble topology optimization results. In 2005, Witzel [START_REF] Witzel | Finite-element model construction for the virtual synthesis of the skulls in vertebrates: case study of diplodocus. The Anatomical Reord, special issue: Finite Element Analysis in Vertebrate[END_REF][START_REF] Witzel | International design conference[END_REF] managed to artificially compute a reptile skull through topology optimization. Through this study, they have shown that seeking for shape ideas in nature can be valuable (Fig. 2).

This seek for nature comparison and inspiration is called bio-mimicry and it is getting more and more investigated [START_REF] Avsec | Biomimicry and 3d printing[END_REF]. Some approaches tend to look for patterns that could be used to define lattice structures. Among them, Voronoi diagrams, Delaunay triangulation or fractals are examples of patterns that can be found in tortoise shells, soap bubble agglomerations, giraffe's skin patterns, and many others. Those schemes can be of interest in the AM context because they can be easily generated and propagated through a volume. However, they Fig. 2 Topology optimization for bio-mimetic simulation of a skeleton [START_REF] Witzel | Finite-element model construction for the virtual synthesis of the skulls in vertebrates: case study of diplodocus. The Anatomical Reord, special issue: Finite Element Analysis in Vertebrate[END_REF][START_REF] Witzel | International design conference[END_REF] Fig. 3 a Examples of atomic lattice structures [START_REF] Rousseau | Cristallographie gTomTtrique et radiocristallographie[END_REF]; b insertion of a crystallographic quadratic type lattice within a sphere represent a difficulty in terms of mechanical models as they require too many elements to be described properly both in terms of geometry and behavior. Non-living structures using lattice can be found in nature, the most famous being crystals [START_REF] Rousseau | Cristallographie gTomTtrique et radiocristallographie[END_REF]. Those lattices have been historically the first to be used in 3D printing for the simple reason that they consist of the duplication of a unitary cell (Fig. 3a) Patterned lattice often refers to crystallographic structures, replacing atoms with nodes and inter-atomic links with beams (Fig. 3b). Such structure is commonly used due to its symmetry and its defined topology, which allows to quickly compute a mechanical response through some finite element models. Crystallographic lattice are defined in terms of inter-atomic planes. Following this analogy, it is then possible to stretch or compress the distances between those planes to increase the set of accessible shapes (Fig. 4). This can be implemented in different ways [START_REF] Hadi | Design configurations and creation of lattice structures for metallic additive manufacturing. 14eme Colloque National AIP PRIMECA[END_REF]. In this paper, as the structure needs to be optimized with respect to mechanical stress, the distances between the crystallographic planes will be defined so as to follow as much as possible the von Mises stress distribution. Thus, high stress density will contract the nearby lattice elements and low stress density will stretch it.

Nevertheless, once an elementary cell has been selected (Fig. 4a), and even though the distances between the crystallographic planes can be adjusted, the topology of the structure is fixed. This reduces the degree of freedom in the definition of the lattice. In some case, it would be interesting to let more freedom also in the definition of the topology of the lattice. The framework proposed in this paper allows such local modifications to optimize the topology as well as the geometry of the lattice structures.

Mechanical model

Another approach to define and optimize lattice structures is to make use of so-called homogenized materials [15]. This consists in looking for a fully filled material with similar properties to the lattice structure that is studied, and use this new material for finite element analysis as a solid material. However, this type of approach owns several shortcomings. For instance, as this equivalent material is set for a given lattice structure, the lattice topology (i.e., its connectivity matrix) cannot be changed throughout the optimization process (Fig. 5). This would mean that the equivalent material needs to be known for the new topology, thus leading to the need of an exhaustive equivalent material library. Anyhow, such an equivalent material would smooth the mechanical properties and the influence of the geometric structure.

Furthermore, this method would also not be directly reliable for un-patterned lattice structures, as the equivalent material would become a set of equivalent materials connected to each other. Such an approach has not yet been tested in the literature, but would likely lead to some difficulties.

Fig. 5 Example of an optimized lattice structure using homogenized materials and preserving the topology (source Autodesk)

Non-uniform patterned lattice

As a conclusion, patterned lattice structures are widespread in current AM processes and are generally preferred to other lattice structures. Unfortunately, there are still gaps between the needs for generating non-uniform patterned lattice structures (Fig. 4c), the actual generation capabilities of 3D modelling software (Fig. 4b) and the capacities of modelling the mechanical response of those structures. This paper introduces a new framework to generate an initial non-uniform patterned lattice structure and optimize its topology while removing beams according to mechanical criteria. The generation and optimization modules are decoupled in the sense that the optimization could also be performed on a user-specified lattice structure that would not have been generated with our approach. The generation module defines the topology and geometry of the lattice, whereas the optimization module acts on its topology to try to reduce the mass while smoothing the stress distribution. During the optimization loops, the mechanical simulations are performed using a technic that derives from the methods of Alzahrani and Wang [START_REF] Wang | A unit cell approach for lightweight structure and compliant mechanism[END_REF]16]. A specific criterion is used to compare the optimality of the resulting structures in terms of stress distribution. The new framework is introduced in Sect. [START_REF] Duysinx | Topology optimization of continuum structures with local stress constraints[END_REF] and the details of the different modules are given in Sect. 4. Section 5 discusses the implementation and results.

Bio-inspired optimization framework

Going further into our bio-mimetic comparison, nature's ultimate issue is minimizing energy consumption. Any organ or body part needs to be fed with nutrients, oxygen and other elements. So, considering that the amount of energy generated per amount of time is fixed for a creature, it needs to distribute this energy among all of its body parts. The percentage of energy distributed to a particular body part depends on the contribution of the body part to the whole body capacity to sustain itself. In this paper, the beams of the lattice can be compared to bones of a skeleton, as their main purpose is to give rigidity to the creature's body (all other purposes of bones, such as fluid circulation, will be neglected to make it simpler). Inspired by this analysis, an optimization framework is defined, that looks for every beam (comparatively to bones in a skeleton), and sorts the contribution of each beam to the structure's stress contribution to decide if "it should be fed", meaning that it should remain, or "starve", meaning that the beam should be suppressed. This is close to what topology optimization does but applied to lattice structures. This principle has been experimented on human bones by manipulating the amount of stress received by a bone and looking at the evolution of the bone's shape and mass over time [START_REF] Howell | an experimental study of the effect of stress and strain on bone development[END_REF]. This is based both on Darwin's theories of evolution [START_REF] Darwin | The origin of species by mean of natural selection[END_REF] and bone structure observations [19]. Those studies infer that in terms of density (that can be related to beams diameter), bones tend to reduce their thickness and re-organize their inner structure as their stress load varies. Of course, such and evolution is done over a relatively long period of time. In this paper, this bio-inspired framework does not modify the beam's radii, but it will be considered that if the beam contribution is under a certain amount, this beam can directly be removed. In the future, further studies could act on radii evolution, thus, leading to a more accurate bio-inspired method. The proposed optimization framework consists of different steps presented in Fig. 6 and detailed in the next section. As previously introduced, to explore a wider range of available lattice structures, the inter-reticular planes can be stretched. To do so, a spatial scalar field is needed. von Mises field will be used here to position the crystallographic planes according to the initial stress distribution (Fig. 6b,c). Then, a lattice pattern is selected based on its mechanical behavior and its printability (Fig. 6d), it is then propagated within the volume (Fig. 6e). A first calculation is performed to ensure that the limit criterion (von Mises maximum stress) is not already reached (Fig. 6f). Finally, the optimization loop is run until the limit criterion is reached (Fig. 6g). As explained, this loop consists of suppressing beams that contribute the less to the mechanical "stress absorption" (this will be discussed later on) of the structure, according to our bio-mimicry statements as well as to a statistical evaluation of the considerate beam contribution among the whole beam population. This Fig. 6 Framework for patterned lattice structure optimization process stops when a certain limit is reached (here, a von Mises maximum value) and the structure is saved.

Algorithms and methods

Crystallographic plane definition

The planes that define the cell's boundaries can be manipulated in such ways that their relative positions can vary non-uniformly in three dimensions. Thus, new designs can be tackled with a wider range of lattice fibers and possible curvature evolutions (Figs. 4,7). To drive the position of the crystallographic planes, the von Mises scalar field is used as an attractive or repulsive field. Therefore, the more the von Mises stress is important, the more the distance between the planes is small. Of course, many other distribution laws could be defined and used. Figure 7b, c shows two distributions which can be obtained through that proposed method.

Cell topology definition

Lattice can be created out of crystallographic structures by replacing atoms with nodes and atom links with beams. Thus, several patterns can be imagined and compared in terms of mechanical behavior. To quickly illustrate the way an initial pattern can be identified, results for two extreme cases are presented: the so-called grid pattern (Fig. 8b) and face-centered cubic (FCC) pattern (Fig. 8c).

The topology of the grid pattern (resp. FCC pattern) can be encoded using a connectivity matrix C Grid (resp. C FCC ) coming from the Graph theory [START_REF] Bondy | Graph theory with applications[END_REF]. Considering a set of (n + 1) nodes P 0 -P n , one can create a set of C 2 potential beams, denoted , ranging from B 00 to B nn and containing all the combinations of node connections. The effective lattice structure  does not necessarily contain all the elements of  (it generally does not). An element of  contained in  is said to be "activated". Thus, the elements c ij of the connectivity matrix  of dimension ((n + 1) × C 2 ) are defined as stated in the following equation:

( 1 if B ij ∈ 
been applied on a unitary cell of size 5 mm and have been loaded with a shear stress of 10 N. The radius of the beams is 1mm. Results are displayed in Table 1.

One can notice that the FCC pattern offers a good compromise between a great reduction of the max displacement (about ten times less) and a relatively low increase of the mass (about two times more). So, even if our objective is to minimize as much as possible the mass, starting with a FCC-based lattice structure is a good idea. Such a choice

C ij = 0 if B ij ∉  (1)
can allow a great variety of orientations and a high possibility of combinations while securing the rigidity aspects.

Even if the Graph theory also allows for -1 values, it is not needed for this purpose. Using this equation, the following two elementary matrices C Grid and C FCC can be obtained and used to define the whole lattice structure by looping through the volume. For the sake of readability, rows containing only 0 are removed. 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 Anyhow, since our optimization process removes beams, it can converge to a solution close to the one that would have been obtained while using the Grid pattern. Finally, one can mention other researches that follow mechanical criteria to insert multiple cell types at different places within the structure [START_REF] Watts | Exploiting the design freedom of RM, Rapid Manufacturing Research Group[END_REF]. In this case, the difficulty lies on the interconnectivity between cells. Even if this can be an interesting ⎡ ⎤ way of improving computation times, it will not be studied ⎢ 1 0 0 1 0 0 0 0 0 0 0 0 0 0 0 ⎢ ⎢ 1 0 0 0 1 0 0 0 0 0 0 0 0 0 0 ⎢ ⎢ 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0 ⎢ in this paper.

⎢ ⎢

Lattice deployment

⎢ 0 1 0 0 0 1 0 0 0 0 0 0 0 0 0 ⎢ ⎢ 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0 ⎢ Given the bounding box of the volume to fill in with cells,

C Grid = ⎢ ⎢ (2
) 0 0 1 0 0 0 1 0 0 0 0 0 0 0 0 ⎢ ⎢ 0 0 0 1 0 0 0 1 0 0 0 0 0 0 0 ⎢ ⎢ 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0 ⎢ the lattice is fully known through both the number of needed cells per dimension and the deformation fields. Thus, two ⎢ ⎢ cases can be distinguished: ⎢ 0 0 0 0 0 1 1 0 0 0 0 0 0 0 0 ⎢ ⎢ 0 0 0 0 0 1 1 0 0 0 0 0 0 0 0 ⎢ ⎢ ⎢ ⎣ 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0 Using our approach, the topology of the initial lattice structure is defined by the connectivity matrix whereas the geometry of the lattice structure is given by the position of the crystallographic planes. Of course, the choice of a pattern impacts the mechanical behavior of the entire structure. To go further in this understanding, those two patterns have Finally, the simplified structure is cleaned to remove self-hanging beams (Sect. 4.4.4).

FE analysis

P (t) = � B (t) × s , ∀t ∈ [0;1] (4)
Once the lattice structure is generated and a load case is i j j=0

With (n i + 1) control points positioned with respect to the local maximums of the von Mises stress that have been computed before. Those control points act like attractors on the curve, contracting planes around positions of high von Mises stress. The number of control points can be adjusted depending on the stress distribution. Adding control points would result in a more contracted lattice structure, even though their number generally remains below four; otherwise, other control curves could be useful such as NURBS. The initial and final control points are positioned on the boundary of the bounding box to cover the entire space. Finally, the user specifies the number N i of cells per direction which is used to position the crystallographic planes according to points P i (t k ) of the Bezier curve located using a uniform law: 1 defined, the simulation model is sent to a finite element solver that returns the average von Mises stress in each beam. This stress distribution provides the contribution of each beam to the stress absorption. Here, the beams are modeled as a wireframe structure (1D elements) to reduce computation times. This simplification brings up some limitations that are not discussed in this paper.

Beam deletion

To optimize the topology of the lattice structure, it is needed to understand how the structure behaves when beams are removed. For that matter, one can look at the structure in terms of inter-connected unit cells. This concept can adapt to any type of pattern, so this discussion will address a random beam neighborhood to remain general. It will only depend on the previously introduced connectivity matrix C. A lattice node is defined as a spatial point where beams interconnect with each other. It can be represented in terms of its neighborhood

t k = k × Δ i , ∀k ∈ {0, … , N i } and Δ i = i (5)
such as in Fig. 9 where the different loads are represented per beam [START_REF] Wang | A unit cell approach for lightweight structure and compliant mechanism[END_REF]. Since the volume to fill in is defined by a triangle mesh whose shape may be complex, a trimming step is required to trim the cells which go outside. Thus, once the planes have been positioned using one of the two technics, it is possible to sort out the intersections between the planes and the triangles. Cells are then created at the intersection between six planes, they are Using the energy conservation principle, assuming the hypotheses of linear static mechanics and that beams have the same radii, one can go back to a simplified approximation of the virtual work principle. In 2D, one can write that the isolated system is such as (Fig. 9): ⎧ ∑ 5 i i=0 f 1 = 0 filled in with a given pattern (Sect. 4.2) and potentially ⎢∑ 5 f i = 0 (6) trimmed if they are intersecting with triangles. Follow-

⎨ i=0 2 ⎢ ∑ 5 i
ing this strategy, it is ensured that every beam is inside the triangle mesh.

i=0 f 3 = 0 ⎩ j circulating through another node. Then, through conservation of energy, the total amount of redirected stress needs to be of 100% and the local redirected amount of stress cannot overcome the deleted amount of stress. This statement leads to the following Eq. 9:

n � e(i) = 1 with e(i) ∈ [-1;1]. i=0 (9)
Fig. 9 Example of a node neighborhood and associated loads [START_REF] Wang | A unit cell approach for lightweight structure and compliant mechanism[END_REF] For each beam removal, Eq. 6 leads to a new stress repartition within the cells of the whole structure. The way the stress distribution evolves strongly relies on the lattice pattern and the load case. In the proposed approach, only lowstressed beams are removed so that the amount of stress that was circulating within the suppressed beam is redirected into the other surrounding beams. Consequently, the stress levels of the surrounding beams rise. Let 8(i) be the variation of stress on node connected to a removed beam. Then, each node in the remaining structure will see its stress level varying and become f �(i) such that Eq. 6 still remains and becomes :

Thus, c(i) represents the percentage of stress that each remaining beam will have to sustain. This phenomenon leads to stress homogenization within the structure. It appears clearly that removing one low-stressed beam will not highly disturb the structure. However, either removing one highstressed beam or too many low-stressed beams could lead to a high amount of stress redistribution, which could lead to the structure disruption.

It is then chosen to remove beams according to a criterion that will select a certain amount of elements based on their stress level and remove them. This criterion will not be discussed here, and so is the overall algorithm.

Optimization criterion

After the iterative removal process, a lattice structure is considered as optimal if all the remaining beams are at the same level of the von Mises stress (this level being defined by the user as the von Mises limit criterion). On that particular matter, as 1D beams are used, one could argue that von Mises criterion does not apply. The finite element solver used:

⎧ ∑ n �(i)
Cast3m [START_REF]Cast3m general notice[END_REF], proposes the same operator called "VMIS" to ∀k ∈ {0, … , 5} ⎢ i=0 f 0 = 0 ... [START_REF] Deng | Topology optimization of total femur structure: application of parametrized level set method under geometric constraints[END_REF] compute an equivalent constraint criterion, no matter if ele-⎨ ∑ n �(i) ments are 3D, 2D or 1D. In the 1D case (beams) the criterion

⎢ i=0 f k = 0 ⎩ is computed such as: f ( EFFX 2 2 2 2 (10) VMIS = SECT + MOMX * DX TORS + MOMY * DY INRY + MOMZ * DZ INRZ
With n the number of beams on a neighborhood. This stands for a 3D structure (i.e., 3D neighborhood). Then, for each node in the remaining structure:

with respect to the Cast3m notations which can be found in [START_REF]Cast3m general notice[END_REF][START_REF]Cast3m, up to date von mises notice and website[END_REF]. In the proposed approach, this stress level directly corresponds to a user-specified limit. From that point, one can define a scalar value A that will be used as a performance

f �(i) i (8) j = f j + c(i) × 8(i),
where ë(i) can be discussed regarding the selected pattern (the different incoming angles, radius of each beams, and the applied load case). Thus, its value will not be defined here. c(i) represents the percentage of stress that was circulating indicator for the algorithm. This scalar value is the area under the curve representing the distribution of the stress level in the lattice structure (Fig. 10). For an optimal structure, i.e., a structure for which the remaining beams are at the same level of the von Mises stress, this area is simply equal to: through the node connected to the deleted beam that is now could use a criterion that mimics those used in topology optimization, i.e., based on either a level set method or a homogenization of the compliance matrix [START_REF] Allaire | Conception optimale de structures[END_REF]. However, an adaptation has to be performed to use this criterion on a lattice structure. A result of such a method on 3D elements is discussed in Sect. 5. In this paper, a statistical data criterion are used to sort low-stressed beams. To be effective, the adopted criterion must be acceptable in terms of computation times and lead to a non-degenerated structure, i.e., a structure which can still sustain the load case. Here, the so-called cut criterion C cut corresponds to the difference between the von Mises mean stress M VM and the von Mises stress standard deviation SD VM :

A optimal = N B ⋅ L VM (11)
C cut = M VM -SD VM ( 12
)
where N B is the number of beams and L VM is the limit von Mises criterion. To compare lattice structures composed of a different number of beams, it is preferable to rescale the von Mises stress and the number of beams over the interval [0; 1] as shown in Fig. 10. This criterion helps evaluating how much the stress is equally distributed in the lattice structure during the optimization loop. Consequently, the closest A is from 1, the closest the structure is to the optimal one (Fig. 10). Thus, at each step of the optimization process, knowing the stress repartition in the structure, one can remove the lattice elements that are unnecessary for the structure rigidity. However, removing the elements one by one would lead to high computation times making it non-suitable for industrial applications. It is thus needed to use a criterion to remove sets of beams and not the beams one by one. Ideally, one Using a statistical approach leads to a non-absolute criterion, each beam gets a notation relatively to other beams and to the maximum stress required. It means that even if the absolute simulated stress values deviate from the experiment, it is still accurate in terms of relative stress distribution, and one can still sort out which beams need to be removed. Moreover, since wireframe simulation models are far faster to compute than 3D models, this method can compute dense lattice structures quickly and as the criterion is relative, one only needs to calibrate a safety coefficient based on expert knowledge to stop the optimization process. Then, if mandatory, a 3D finite element analysis can be used to get more accurate results. The proposed criterion gives both the position of the beam within the whole population in term of von Mises stress (through standard deviation) and its distance to mean stress. For an homogeneous population of beams, this leads to a removal of less than 10% of the total population at each iteration [START_REF] Teufelhart | optimization of strut diameters in lattice structures[END_REF]. Furthermore, the number of beams that will be removed is proportional to the total number of beams remaining within the structure. This allows to auto-regulate the beam population. This process directly uses the von Mises (approximated von Mises due to wireframe models) distribution. It affects the connectivity matrix C and the final lattice structure is not anymore a pattern. Figure 11 illustrates two stress distributions for the two crystallographic plane arrangements of Fig. 7. A pure shear load case of 10 N is applied at the extremity and no beams are removed. Even for the same FCC elementary cells, and for the same number of beam elements, the stress distribution can vary depending on how the crystallographic are distributed. Figure 11b shows that the nonuniform adaptive lattice has a more homogeneous stress repartition (see decile of von Mises stress). Applying our performance criterion would lead to a value A = 0.027 for the uniform lattice and A = 0.030 for the non-uniform lattice (11% gain) (Fig. 7).

Cleaning

Once the beams have been iteratively removed and the stop criterion has been reached, self-hanging beams are removed. They correspond to beams that sustain no stress. Thus, beams that have a stress level of 0 MPa (±c) are removed. As a remark, one can understand that as 1D beams are used, the stress level of a beam B ij is computed as: 1

removing of unnecessary beams can then not be based on their stress level but on the fact that one of their nodes is not stressed. This issue could be address more effectively using the C matrix to sort out overhanging beams (nodes with only one beam connected) and delete them.

Implementation and results

The proposed framework has been implanted on our prototype platform developed in C# and based on MVVM WPF technology. It creates non-uniform patterned lattice structures that are then optimized using Cast3m as a finite element solver, integrated in our optimization loop. To validate the proposed approach, the case study of Fig. 7 will be used.

A load case of 10 N perpendicular to the structure main axis is applied at the extremity, and the user-specified maximum von Mises criterion L VM is set to 40 MPa. After few iterations, the structure has evolved towards an optimized structure that suits the technical requirements and allows an effective mass reduction. The stress has been homogenized within the structure that is not degenerated (Figs. [START_REF] Rousseau | Cristallographie gTomTtrique et radiocristallographie[END_REF][START_REF] Hadi | Design configurations and creation of lattice structures for metallic additive manufacturing. 14eme Colloque National AIP PRIMECA[END_REF].

VMIS(B ij ) = 2 ⋅ (VMIS(N i ) + VMIS(N j )) (13)
Then, for an overhanging beam, the stress level on one node would be 0, and on the other node would it have a certain value, leading to a beam stress level superior to 0. The Here, our performance indicator is A = 0.11 for the uniform lattice and A = 0.14 for the non-uniform one with a relative gain of 27% (Fig. 12). This illustrates a better von Mises stress distribution throughout the beam population.

As one can see comparing Figs. 11, 12, the optimization process tends to smooth the curve of stress repartition, illustrating the homogenization process. Figure 12 also contains the decile repartition of von Mises stress for the uniform and non-uniform lattice structures. Comparing Figs. [START_REF] Witzel | International design conference[END_REF][START_REF] Avsec | Biomimicry and 3d printing[END_REF] shows that through the overall diminution of beam elements, the stress level in each of those elements has raised up toward the user-specified von Mises limit criteria. It corresponds to the homogenization phenomenon.

Figure 13 shows that the denser parts of the resulting structures are located on the upper and lower planes near the fixations, and that the center area's density has been significantly lowered. Such a behavior is consistent with what usual optimization leads to, like IPN structures for example.

Table 2 gives output values, the ratio column provides the percentage of relative values from non-uniform to uniform structure. The gain from initial to final volume is of 85.6% for the non-uniform structure and 81.1% for the uniform, showing that with the same initial number of beam, their spatial repartition is of importance for the optimization process. Finally, Table 2 compares all the results and sorts the obtained gain for each value (gain from non-uniform to uniform lattices). Mass was the main objective of the optimization process and can be improved, for this test case, by 33% only by stretching the lattice structure. One can also see that mean stress is higher (and so closer to the limit criterion) when using a non-uniform lattice.

Figure 12 shows that the stress distribution is better than for the initial lattice structure (Fig. 11) as most of the beams have a stress level under [mean stress + standard deviation]. The performance indicator A is still smaller than 1, but its value has been multiplied by a factor 5 compared to the initial structure (from 0.03 to 0.14). In this work, only the removal of beams was authorized, but allowing other transformations could lead to an even better distribution and could help raising the A value further more. Among those foreseen transformations, one could cite beam radius evolution, beam switching and beam to plate transformation.

To further analyze the results, a classic topology optimization based on SIMP method has been run on the initial fulfilled volume. The results are shown in Fig. 14. One can see that the general shape of our optimized lattice structure tends to keep material in the same place as the classic topology optimization does.

Finally, one can see that the way the crystallographic planes are distributed is important for the optimization process. Thus, alternative plane distribution rules still have to be studied so as to improve even more the performance indicator and reduce the mass.

Conclusions and future works

In this paper, a non-uniform patterned lattice structure optimization framework has been proposed and validated. The optimization process starts with an initially non-uniform structure and iteratively removes beams according to the von Mises stress distribution within the whole structure. At the end, the structure becomes un-patterned, the stress is better distributed and the mass is reduced. Compared to Fig. 2, results of Fig. 14a, b illustrate that the proposed approach also tends to generate nature-like structures.

A first remark one would have would be to use directly the constraints tensor field instead of the von Mises scalar field. This would allow to differentiate tension from compression, and to enable a better space redistribution of the crystallographic planes. For the sake of simplicity, and as a first attempt, von Mises criterion has been chosen precisely for being a simple scalar instead of nine scalars embedded in a tensor, but future works will extend this particular point. Such method has already been tested and has shown interesting results [START_REF] Tam | Stress line additive manufacturing (slam) for 2.5-d shells[END_REF].

Similarly, using 1D finite elements does not allow to go up to a complete optimization, but to get a lattice shape containing as less beams as possible. Then, it is necessary to run a 3D finite element calculation to ensure that the technical requirements are met. This is an alternative to common approaches that still needs to be improved but allows to change lattice topology during optimization contrary to homogenized material methods (Fig. 5). This widens the design possibilities.

Finally, through Figs. 12, 14 one can see that there is still potential for pushing the von Mises front up to the limit criterion and get the A criterion as close as possible to 1. As already said, this improvement could be done by allowing more transformations. Instead of only removing beams, beams could be switched, new beams could be added, beam's radii could be changed and so on [START_REF] Teufelhart | optimization of strut diameters in lattice structures[END_REF]. Moreover, Fig. 14 and Table 2 illustrate that working on non-uniformly distributed crystallographic planes improves significantly the performance indicator A and leads to a significant mass reduction when compared to a uniform lattice structure. Thus, the way those planes are located throughout space still needs to be thought. It can be seen as a pre-optimization process. This method could give access to a wider range of reachable geometries and topologies for complex structures, improving material space repartition in a continuous topology-related fashion, thus reducing significantly the mass of systems through lattice structures.

Finally, the comparison between the adopted simulation models and the corresponding experiments was not tackled and still remains an issue for this method to be efficient.

Fig. 1

 1 Fig. 1 Satellite solar panel deployment mechanism (courtesy of Thales Alenia Space, ESA and Poly-Shape): a non-additive manufacturing design, b additive manufacturing redesign

Fig. 4

 4 Fig. 4 Example of lattice structures: a elementary cell, b uniform lattice deployment, c adaptive lattice deployment

Fig. 7

 7 Fig. 7 Insertion of lattice structure: a von Mises field. b Uniform lattice insertion. c Non-uniform lattice insertion

Fig. 8

 8 Fig. 8 Cell topology definition: a FCC crystallographic atoms and links. The grid pattern only uses corner nodes P0-P7 (b), and the FCC pattern uses both corner and centered nodes P0-P14 (c)

Fig. 10

 10 Fig. 10 Rescaled distribution of the stress level for the current and optimal structures

Fig. 11

 11 Fig. 11 Initial A criterion for both classic and stretched lattice structures (Fig. 7), von Mises stress value per element sorted form highest to lowest. a Uniform lattice, b non-uniform lattice

Fig. 13 Fig. 12

 1312 Fig. 13 Isometric views of the final structures. a Uniform lattice. b non-uniform lattice

Fig. 14

 14 Fig. 14 Comparison between a SIMP continuous topology, b nonuniform lattice structure and c uniform lattice structure

Table 1

 1 Estimated result values for both structures (for a same given nxnxn lattice cell unitary size)

	i	j,n

-

  Uniform lattice In this case, a simple loop based on the (dx, dy, dz) unit cell's dimensions, and the (X, Y, Z) bounding box's dimensions is used. The user gives the number N i of cells per axis i ∈ {x, y, z} through which the unit cell's dimensions can be computed:In this paper, the mass of the structure is reduced step by step while removing beams of an initially patterned lattice as defined in Sect. 4.2. The removal is performed according to a mechanical criterion, i.e., the distribution of von Mises stress.

			4.4 Optimization loop
	dx = X ⎧ N x ⎢ dy = Y ⎨ N y ⎢ dz = Z ⎩ N z	(3)	Manufacturing constraints are not fully addressed here, as the FCC pattern is printable regarding its stretched dimensions in x, y and z. Further work would be needed to consider this printability constraints. The method is discussed in the next
	-Non-uniform lattice Here, a deformation field is needed to drive the crystallographic plane stretching. In the proposed approach, Bezier curves are used to define this deformation field along the three directions of the	sections. The overall idea is to sequentially run a FE analysis (Sect. 4.4.1) to define the beams with a low stress, and then remove those beams (Sect. 4.4.2). The removal step is per-formed while following an optimization criterion (Sect. 4.4.3).
	bounding box. For a given direction of the bounding
	box i ∈ {x, y, z}, the Bézier curve is defined as follows
	[22]:
	ni

Table 2

 2 Comparison of the maximum stress and associated volume for the initial/final and uniform/non-uniform structures

		Uniform	Non-uniform	Ratio (%)
	Nber of beam
	Initial	117,376	117,376	0.00
	Final	15,680	13,514	16.03
	Max VM (MPa)
	Initial	9.3	12.8	+27.3
	Final	40	39.6	1.01
	Volume (mm 3	)
	Initial	5092	4998	1.88
	Final	960.6	721.1	33.2
	Mean VM (MPa)
	Initial	1.11	1.21	+8.26
	Final	4.4	5.6	+21.43
	Std. deviation
	Initial	1.8	1.88	+4.26
	Final	15.2	24.2	+37.23


Clément Dupont

PhD Maximilien Vermandel
email: maximilien.vermandel@chru-lille.fr

PhD Henri-Arthur Leroy

MD Mathilde Quidet

MD Fabienne Lecomte

MSc Nadira Delhem

PhD Serge Mordon

PhD Nicolas Reyns

INtraoperative photoDYnamic Therapy for GliOblastomas (INDYGO): study protocol for a phase I clinical trial (Short title: Intraoperative Photodynamic Therapy for GBM)

Keywords: NCT number: NCT03048240 EudraCT number: 2016-002706-39 5-ALA, photodynamic therapy, fluorescence-guided resection, glioblastoma, neurosurgery, medical device, laser therapy Protocol title, protocol identifying number

for providing the extraperitoneal distension balloons (

Pr. Nicolas Reyns, MD, PhD, Principal Investigator, Department of Neurosurgery, CHU Lille, INSERM OncoThAI U1189, +33320446721 nicolas.reyns@chru-lille.fr Dr. Maximilien Vermandel, Co-Principal investigator, PhD, Department of Medical Physics, CHU Lille, INSERM OncoThAI U1189, +33320446721 maximilien.vermandel@chru-lille.fr Name(s) and address(es) of the clinical laboratory(ies) and other medical department(s) or institutions involved in the research. Institut National de la Santé Et de la Recherche Médicale, 1, Avenue Oscar Lambret, 59000 Lille, France

I.

Rationale and Background Information: Glioblastoma (GBM) is a brain tumor with a particularly poor prognosis. The median overall survival is approximately 15 months with the current standard of care (SOC) [START_REF] Iacob | Current data and strategy in glioblastoma multiforme[END_REF][START_REF] Ostrom | The epidemiology of glioma in adults: a "state of the science" review[END_REF][START_REF] Stupp | Radiotherapy plus concomitant and adjuvant temozolomide for glioblastoma[END_REF] . Although GBM is a rare neoplastic disease with low prevalence (3-5 /100,000 persons), it remains the most frequent primary malignant brain tumor in adults 4,[START_REF] Fisher | Epidemiology of brain tumors[END_REF] . Currently, no existing therapeutic agent is able to stop GBM progression. Complete resection is rarely feasible since isolated tumor cells usually infiltrate the surrounding brain. Adjuvant therapies to improve local control are thus highly expected.

Currently, the European Society of Medical Oncology (ESMO) guidelines provide recommended treatment options according to the grade of the glioma [START_REF] Stupp | Highgrade glioma: ESMO Clinical Practice Guidelines for diagnosis, treatment and followup[END_REF] . Given the invasive nature of GBM, local relapse is unavoidable. Radiotherapy (RT) applied to the tumor bed with concomitant and adjuvant chemotherapy with temozolomide (TMZ) complete the surgery and prolong the progression-free survival (PFS) and the overall survival (OS). Although these adjuvant treatments statistically improve outcomes, the benefits remain modest. An increase of less than 3 months in the median OS with RT and TMZ versus RT alone has been observed [START_REF] Stupp | Radiotherapy plus concomitant and adjuvant temozolomide for glioblastoma[END_REF] . Therefore, the management of newly diagnosed patients should include a maximal resection when achievable; the extent of resection (EOR) is of prognostic value. To optimize the EOR, surgery commonly goes along with an intraoperative assistance to diagnosis residual tumor tissues.

Fluorescence-guided resection (FGR) with 5-aminolevulinic acid (5-ALA) as a prodrug uses the fluorescence phenomenon to distinguish healthy tissue from tumor tissue in order to support the surgical procedure. Due to a defective regulatory process in tumor cells, exogenous 5-ALA induces the accumulation of a photosensitive molecule, protoporphyrin IX (PpIX), in tumor cells only. When the surgery cavity is illuminated under blue light (375-440 nm), the tumor tissues shine red (635 nm) (see figure 1). An increase in the EOR and an improvement in PFS has been observed when 5-ALA and is used for FGR [START_REF] Stummer | Fluorescence-guided surgery with 5-aminolevulinic acid for resection of malignant glioma: a randomised controlled multicentre phase III trial[END_REF] . However, it has been shown that extended resections performed using FGR carry the risk of temporary impairment of neurological function [START_REF] Stummer | Counterbalancing risks and gains from extended resections in malignant glioma surgery: a supplemental analysis from the randomized 5-aminolevulinic acid glioma resection study. Clinical article[END_REF] . Furthermore, relapse occurs in the cavity margin in 85% of cases [START_REF] Petrecca | Failure pattern following complete resection plus radiotherapy and temozolomide is at the resection margin in patients with glioblastoma[END_REF] . The improvement of local control by EOR and additional therapy delivered intraoperatively is a key point in increasing progression-free survival and overall survival [START_REF] Pavlov | Combining intraoperative carmustine wafers and Stupp regimen in multimodal first-line treatment of primary glioblastomas[END_REF] . In this context, intraoperative photodynamic therapy (PDT) appears relevant in treating areas bordering the resection cavity while maintaining a maximal but safe resection. PDT is based on the photochemical reaction activated by light between a photosensitive molecule and oxygen. This reaction creates reactive oxygen species, including singlet oxygen, which leads to cell death. The selective accumulation of PpIX in tumor cells only makes this treatment local [START_REF] Valdes | delta-aminolevulinic acid-induced protoporphyrin IX concentration correlates with histopathologic markers of malignancy in human gliomas: the need for quantitative fluorescence-guided resection to identify regions of increasing malignancy[END_REF][START_REF] Johansson | Protoporphyrin IX fluorescence and photobleaching during interstitial photodynamic therapy of malignant gliomas for early treatment prognosis[END_REF][START_REF] Johansson | 5-Aminolevulinic acidinduced protoporphyrin IX levels in tissue of human malignant brain tumors[END_REF] . In addition to 5-ALA being used for intraoperative diagnosis, previous human experiences using 5-ALA as a prodrug for PDT have been reported and highlight the positive role of PDT for the treatment of GBM [START_REF] Johansson | Protoporphyrin IX fluorescence and photobleaching during interstitial photodynamic therapy of malignant gliomas for early treatment prognosis[END_REF][START_REF] Stummer | Long-sustaining response in a patient with nonresectable, distant recurrence of glioblastoma multiforme treated by interstitial photodynamic therapy using 5-ALA: case report[END_REF][START_REF] Beck | Interstitial photodynamic therapy of nonresectable malignant glioma recurrences using 5-aminolevulinic acid induced protoporphyrin IX[END_REF] . Additionally, an immune response seems to be triggered during PDT treatment, which would allow a delayed therapeutic action on GBM [START_REF] Maeding | Boosting Tumor-Specific Immunity Using PDT[END_REF][START_REF] Wachowska | Immunological aspects of antitumor photodynamic therapy outcome[END_REF] .

Illumination of the resection cavity remains challenging since the device has to suit the shape of the cavity. In this context, a deformable balloon is expected to meet this requirement. Several devices have been designed for this purpose in the past, such as an inflatable balloon with different light sources [START_REF] Wilson | Instrumentation and light dosimetry for intraoperative photodynamic therapy (PDT) of malignant brain tumours[END_REF][START_REF] Dwyer | Optical integrating balloon device for photodynamic therapy[END_REF][START_REF] Moseley | In vitro light distributions from intracranial PDT balloons[END_REF][START_REF] Eljamel | ALA and Photofrin fluorescence-guided resection and repetitive PDT in glioblastoma multiforme: a single centre Phase III randomised controlled trial[END_REF][START_REF] Schmidt | Light-emitting diodes as a light source for intraoperative photodynamic therapy[END_REF] , an optical fiber coupled to a reflector, or a laser light source coupled to the surgical microscope, enabling a step by step illumination [START_REF] Akimoto | Preliminary clinical report on safety and efficacy of photodynamic therapy using talaporfin sodium for malignant gliomas[END_REF] . Moseley et al. described the only clinical study using 5-ALA as a precursor for the treatment of GBM with PDT [START_REF] Moseley | In vitro light distributions from intracranial PDT balloons[END_REF][START_REF] Eljamel | ALA and Photofrin fluorescence-guided resection and repetitive PDT in glioblastoma multiforme: a single centre Phase III randomised controlled trial[END_REF] . For this study, two photosensitizers, 5-ALA and Photofrin® (Pinnacle Biologics, Chicago, IL), were considered, and repetitive PDT was evaluated. 5-ALA was used for the FGR procedure, and PpIX was activated during the first illumination. Following the procedure, further illuminations (post-operative) were applied to activate Photofrin®, which has a longer half-life than 5-ALA. Thus, none of these previous studies used 5-ALA as a primary photosensitizer precursor to provoke a PDT effect on areas bordering the resection cavity, as proposed in the present study.

To be more acceptable for both the surgeons and the patients, intraoperative PDT has to be fully embedded within the usual surgical workflow. Consequently, the concept of this new therapeutic approach was to deliver 5-ALA PDT intraoperatively early after FGR.

This study adopts a repurposing strategy, using 5-ALA in a therapeutic manner for glioma, and hence develops a model for future therapeutic strategies for treating high-grade gliomas combined with the current SOC in addition to intraoperative 5-ALA PDT. Intraoperative PDT takes advantage of the presence of the photosensitizer, which is already approved and administered to the patient for the FGR procedure [START_REF] Leroy | Fluorescence guided resection and glioblastoma in 2015: A review[END_REF] . Indeed, even if no fluorescence is observed after FGR, remaining tumor cells around and deeper in the cavity are still present and take up PpIX at a concentration, limiting their visualization through the microscope. These remaining cells that are not microscopically visible are the main target of the PDT treatment that is delivered by means of a new, recently designed lighting device [START_REF] Dupont | A novel device for intraoperative photodynamic therapy dedicated to glioblastoma treatment[END_REF] .

The proposed treatment strategy can fit into the current SOC and can be applied at the initiation of care (see figure 2) [START_REF] Stupp | Highgrade glioma: ESMO Clinical Practice Guidelines for diagnosis, treatment and followup[END_REF] . Thus, intraoperative 5-ALA PDT is more ethically acceptable compared with the use of interstitial PDT on the resectable tumor.

II. Study Goals and Objectives

The current INDYGO study aims to assess the feasibility and safety of the intraoperative 5-ALA-mediated PDT procedure in patients with a newly diagnosed GBM that is accessible for complete surgical removal with contrast enhancement. The primary endpoint is to assess the feasibility of intraoperative PDT early after the surgical resection of GBM without unacceptable and unexpected toxicities. The secondary endpoints are the evaluation of the following:

-progression-free survival, -overall survival, -treatment response evaluated by MRI every three months until tumor relapse, -toxicity evaluated every three months until tumor relapse and follow-up of serious unexpected events potentially linked to patient death, -quality of life measured every three months until tumor relapse.

To validate the study, at least seven out of the ten patients included in the study must be free of unacceptable and unexpected toxicities.

An immunological study has been added to the workflow to explore biological markers of the PDT-induced immune response.

III.

Study Design This study is an interventional, single group, pilot study and single center (university hospital of Lille, France). Enrollment of 10 patients is expected. All patients with newly-diagnosed probable GBM based on clinical and radiological criteria and in whom surgery is achievable could be recruited. The main exclusion criteria are contraindications to 5-ALA or MRI procedures.

IV.

Methodology A simplified workflow of the study is presented in figure 3. Once a patient is included in the INDYGO trial, MRI (T1, T1-weighted, T2, FLAIR, diffusion/perfusion sequences) and blood samples for immunology measurements are acquired. Next, 5-ALA is orally administered 4 hours before the surgery (20 mg/kg), and a standard surgical resection is performed, followed by the FGR procedure. Once there is no remaining fluorescence or potentially infiltrating functional areas, an intraoperative MRI is performed to assess the quality of the resection. After the MRI examination, the PDT procedure starts. During the setup of the device, a quality control evaluation of the laser is performed. The lighting device is placed inside the surgical cavity and is inflated with diffusing liquid until it reaches the cavity boundaries. The treatment time is deduced from the quantity of diffusing liquid injected in the device. The treatment sequence is composed of five fractions of illumination separated by two-minute breaks to re-oxygenate brain tissues. Soon after the PDT treatment is completed, a second intraoperative MRI and blood samples are acquired to assess potential early and acute effects of PDT [START_REF] Leroy | MRI assessment of treatment delivery for interstitial photodynamic therapy of high-grade glioma in a preclinical model[END_REF] .

Figure3: Flowchart of the INDYGO trial

After surgery, patients undergo an SOC, including [START_REF] Stupp | Highgrade glioma: ESMO Clinical Practice Guidelines for diagnosis, treatment and followup[END_REF] RT and concomitant TMZ, followed by adjuvant TMZ. MRI and blood samples are then acquired every three months alongside patient follow-up until recurrence.

Alongside this study, biological samples are collected and enable an exploratory study. This collection will be used to design future immune monitoring for each patient. Several parameters will be explored, such as measurements of cytokines, growth factors, immune cell profiles and peripheral blood mono-nucleated cells.

V. Discussion

The main difficulty of this clinical trial relies on the use of a new medical device simultaneously with the evaluation of a drug repurposed for this condition. This issue was raised due to a comprehensive assessment of this new medical device [START_REF] Dupont | A novel device for intraoperative photodynamic therapy dedicated to glioblastoma treatment[END_REF] . Thus, once the safety of the medical laser system combined with the medical device was evaluated in a preclinical stage, a risk management procedure was delivered to an institutional review board for ethical consideration and to the French National Medicine Agency for legal consideration. Additionally, a comprehensive analysis of the previous use of 5-ALA for the purpose of PDT associated with an efficacy and safety preclinical trial in a rodent model [START_REF] Leroy | MRI assessment of treatment delivery for interstitial photodynamic therapy of high-grade glioma in a preclinical model[END_REF] validated the repurposing of the drug for therapy. Managing both the medical device and the drug issues remain specific for PDT, and both have to be considered before being approved for human use.

VI. Trial Status

The INDYGO clinical trial is still in the recruitment phase.

VII. Safety Considerations

For each patient recruited, all the adverse events are recorded from the beginning of PDT treatment to relapse, up to 24 months. This collection of adverse events is reviewed by an independent safety monitoring board twice: in the middle and at the end of the study. The independent safety monitoring board is composed of three neurosurgeons from different university hospitals in France. Additionally, as the primary endpoint of this trial is the safety of PDT, this safety board composes also the data monitoring committee.

VIII. Follow-up

Several monitoring visits are planned after the PDT treatment, as recommended with the SOC6. RT and TMZ sequences are still administered to all patients. The INDYGO trial utilizes blood sample collections and MRI imaging every 3 months until recurrence. A quality of life questionnaire has been established to investigate the impact of PDT treatment on the patient's quality of life (EORTC QLQ-C301). This collection is filled in every 3 months until relapse, by the signature of a consent form, and is assessed up to 24 months. This questionnaire includes a global health status and a quality of life scale.

IX.

Data Management and Statistical Analysis Since this is a pilot feasibility study, in agreement with the statistical board, the number of subjects is not calculated on the basis of the statistical hypothesis. A total of 10 patients will be recruited over a planned period of one year. The feasibility objective is to have a minimum of seven patients included who have been fully treated with PDT, without unacceptable and unexpected toxicity. Furthermore, this number of 10 patients seems compatible with the average number of patients recruited annually in the investigative center for the surgical treatment of GBM.

Figure 1 :

 1 Figure 1: Fluorescence-guided resection under visible light (a) and blue light (b). Residual tumor cells shine red due to the fluorescence property of PpIX.

Figure 2 :

 2 Figure 2: Workflow of the intraoperative PDT procedure

http://groups.eortc.be/qol/eortc-qlq-c30

The authors would like to thank Estelle GALLET and Laurent JEAUNEAU from B. Braun (

All the data are described individually using case reports. If possible, a Kaplan-Meier curve will be established for the PFS. The analysis of the data will be carried out with the methodological support platform of Lille University Hospital.

X.

Quality Assurance All the medical observations are kept in the patient's file; the data regarding the study are reported on case report forms (CRF) according to good clinical practice. Any deviation from the protocol must be noted, and the reason for deviation must be documented. The data collection will be exhaustive and checked regularly by a clinical research nurse, according to the protocol procedures. The monitoring of the trial is carried out according to the monitoring plan. A planning meeting with the principal investigator has been held before the start of the trial. During the trial, several checkpoints will be defined, including the presence of signed informed consent forms obtained by the neurosurgeon, respect of the inclusion and exclusion criteria, the reporting of any adverse event and the monitoring of all steps of the patient's follow-up. At the end of the trial and once the final analysis is completed and validated, all the files and data will be sealed and archived, according to specific procedures, in a secure location.

XI.

Expected Outcomes of the Study The main expected result of this INDYGO study is the feasibility and safety assessment of intraoperative PDT. This is a crucial requirement to consolidate data prior to a multicenter randomized clinical trial. As secondary objectives, the study aims to demonstrate the potential effects of PDT in terms of progression-free survival and to provide new data on the effectiveness of PDT by radiological biomarkers. From the patient side, the translational impact is highly relevant. Local control is highly expected, and PDT brings new solutions to its optimization. Adding PDT to the GBM SOC is expected to significantly increase progression-free survival.

XII. Duration of the Project

Enrollments are scheduled for a period of 12 months, and follow-up is scheduled for 24 months.

XIII. Project Management

Pr. Reyns is the principal investigator of the study and will submit the final report. Dr. Vermandel is the scientific manager of the project. Pr. Mordon is the scientific expert. They will coordinate the publication plan to communicate trial results.

XIV. Ethics

The authors state that they have obtained appropriate institutional review board approval and have followed the principles outlined in the Declaration of Helsinki for all human or animal experimental investigations. In addition, for investigations involving human subjects, informed consent are collected from the participants involved. Access to the data are in agreement with the French National Agency regulating Data Protection (Commission Nationale Informatique et Libertés, CNIL).


RESS OR SAS: ELEMENTS OF CHOICE FOR FINE PARTICLE PRODUCTION

INTRODUCTION

Supercritical fluid technology presents a new and interesting route for particle generation. Most of the drawbacks of the traditional routes, milling or liquid crystallisation, are avoided. Supercritical processes, RESS and SAS, give micro-or even nano-particles with narrow size distribution and can also be used to achieve micro-encapsulation or surface coating of an active substance particle with a polymer.

The bio-availability of pharmaceutical molecules depends on their absorption by the gastrointestinal tract which is governed by their dissolution rate. Micronisation techniques and especially those using SC-CO 2 which lead to an increase in the specific surface area (typically several tens of m 2 /g for a non-porous solid) can significantly contribute to this bio-availability improvement.

It is known that RESS is usable with CO 2 -soluble molecules while SAS can process nonsoluble molecules. However, the choice is not so simple in reality. Very often, the solubility of pharmaceuticals is very low in pure CO 2 at moderate pressures (10 to 20 MPa) and temperatures (313 to 333 K). These operating conditions can be changed by the use of a cosolvent or by tuning the operating P or T. To know whether these possibilities exist or not, a good knowledge of the phase equilibria is required.

In this paper we will discuss how to determine the right strategy to be applied when facing such a problem and how to choose between a RESS or a SAS process.

OPTIMISATION OF THE RESS PROCESS

RESS is a two-step process : after having solubilised a substance in a SCF, the mixture is suddenly depressurised in a nozzle causing fast nucleation and fine particle generation. In designing a RESS process, it is necessary to have the best possible understanding of what happens upstream from the nozzle that is to say during the extraction step. Therefore, it is important either to collect data from the literature when they exist, or to perform experiments or modelling, about the solute solubility in the SCF [START_REF] Sauceau | A novel apparatus for accurate measurements of solid solubilities in supercritical phases[END_REF].

The key parameters of this extraction step are obviously the operating T and P. The flow rate of the fluid may also play an important role since a thermodynamic equilibrium may or may not be reached in the extraction autoclave. In fact, the kinetics of the dissolution must not be neglected and diffusional limitations can occur. Additional problems may be encountered when the solid to be extracted is not a pure component. Fractionation of the load can lead to a variation in the composition of the particles upon depressurisation.

To understand what will happen during fluid expansion, a Mollier diagram (pressure P vs. molar enthalpy h) is very useful. An acceptable assumption is to consider that the depressurisation is an iso-enthalpic process (in fact, there is a variation of kinetic energy of the fluid, which should be taken into account). Therefore, the expansion step will be a vertical line (dashed arrow) on this diagram. Such a drop in pressure implies an important decrease in temperature as shown on figure 1 which is the Mollier diagram for pure CO 2 . In the example shown on this figure, the starting point as well as the arrival point, are chosen to be in the single-phase zone of the diagram. However, to ensure such conditions, the temperature before expansion must be high (380 K in the example given), which may be incompatible with the solute to be extracted. A lower pre-expansion temperature may lead to a condensation or freezing of the SCF (the triple point of CO 2 is at 0.52 MPa).

To compensate this drop in temperature and to avoid potential problem of clogging, a heated nozzle or a heating device just upstream from the nozzle can be used. The postexpansion pressure is also of the utmost importance.

Nucleation of the solute occurs during the sudden expansion of the solution and is due to the mechanical perturbation, which propagates at the speed of sound leading to very uniform conditions (2) and therefore to a very narrow particle size distribution.

Indeed, the expansion time depends not only on the speed of sound (which varies with the operating conditions and the nature of the fluid) but also on the length of the capillary nozzle.

← S L V → S S L L L V V S V
However, an approximation between 10 -4 to 10 -6 second has been proposed. It is during this very short period of time that solubility falls by several orders of magnitude leading to large supersaturation ratios. The driving force of the nucleation process is the difference between the chemical potential of the solid in both phases (fluid and solid) which is related to the activity of the solute in the solution and at equilibrium. It is therefore the supersaturation which drives the nucleation. This supersaturation comes from the drop in density, which is very large especially near the critical point. The value of the enhancement factor (ratio of the solubility of a given solid in SCF over the theoretical solubility in an ideal gas) can be as high as 10 5 or 10 6 . These values give an idea of the large supersaturation ratios obtained, provided the fluid after expansion can be considered to be an ideal gas. This is the case when the expansion pressure is close to atmospheric pressure. However, the expansion, may not be complete and can be controlled. It may be advantageous to limit the depressurisation and to let a pressure of a few MPa. This will diminish the supersaturation ratio but the recycling of the fluid might be facilitated. In fact relatively small changes in pressure may cause a dramatic fall in density and therefore in solubility. Again, the h,P diagram shown in fig. 1 is of great help to determine which conditions can be used.

The most obvious drawback of RESS is the lack of solubility of several families of molecules in the CO 2 . To overcome this problem one way is to change the supercritical fluid. However, this is seldom possible since the few other candidate molecules (N 2 O, light hydrocarbons,...) are much more hazardous than is CO 2 [START_REF] Raynie | Warning Concerning the Use of Nitrous Oxide in Supercritical Fluid Extractions[END_REF].

To overcome this disadvantage, another possibility is to use the RESS process with a cosolvent being previously dissolved in the CO 2 . In this case the advantage of not using any chemical solvent is lost but it might be interesting to use a co-solvent using such low-toxicity solvents like acetone or ethanol.

The diagram h,P of the binary mixture is different from the pure component diagram. As an example, the diagram for the mixture (95 % CO 2 -5% Ethanol) shows that the two-phase envelope is considerably enlarged by the addition of 5 % co-solvent (figure 2). This has to be taken into account when determining the operating conditions. The maximum mole fraction has to be determined to maintain a single phase fluid after depressurisation.

The choice of the proper co-solvent is not trivial and the affinity two-by-two supercritical fluid, co-solvent and solute are to be considered. A mixture of co-solvents may also be used . Moreover, it is not possible to predict the co-solvent effect from the solubility of the solute in the liquid co-solvent. Therefore, solubility measurements of the solid in the mixture SCF-cosolvent are highly recommended (4 ).

OPTIMISATION OF THE SAS PROCESS

In the SAS process the substance to be micronised is dissolved in a solvent. The SCF acts as an anti-solvent that causes the crystallisation of the solute.

We will limit the discussion to semi-continuous process where solvent and anti-solvent are co-injected at the nozzle level. This process is sometimes referred as ASES or SEDS. Several parameters have an influence on the size and the morphology of the obtained particles. As with the RESS process the operating T and P upstream and downstream from the nozzle are obviously key-parameters. The role of the nozzle is to favour droplet formation which is rather different from the RESS case in which the main role of the nozzle is to provoke a huge drop in pressure. Nozzle geometry is also an important parameter in controlling the process. In addition, the concentration of the solute in the solution and the relative flow rates of the solution and the anti-solvent are also to be taken into account.

The volumetric expansion of the organic solvent due to the dissolution of the anti-solvent has been recognised as an important parameter [START_REF] Chang | Solvent Expansion and Solute Solubility Predictions in Gas-Expanded Liquids[END_REF]. The volume expansion can be predicted by using an equation of state with adequate mixing rules. Several authors have measured this expansion depending on the operating conditions and on the nature of the solvent. These data are of the utmost importance for the batch mode processes where the particle are produced from a liquid phase rich in solvent.

Fig. 2. P,h diagram of a CO 2 -Ethanol mixture. Mole fractions: CO 2 =95% and Ethanol=5%

(Calculation : PR EOS) For the semi-continuous mode, the expansion curves are less pertinent since the phenomena involved are rather different. The particles are formed from the drops of the solution which collide into the fluid phase. As the antisolvent flow rate is generally much higher than those of the solution, one can use the P,h diagram to see how many phases are present at the equilibrium.

On figure 2, the diagram of a 95 % CO 2 -5% Ethanol mixture, may correspond to a practical ratio of the anti-solvent over the solution of about 20. This diagram shows the twophase L+V envelope which is to be avoided in order not to have a liquid phase in the expansion chamber. Such a liquid phase might re-solubilise the particles formed. In these semi-continuous processes which are those most studied today because of their potential industrial applications, the mass-transfer and nucleation occur at the surface of solution droplets in a continuous flow of antisolvent. It is therefore necessary to understand the fate of these droplets. There is a competition between two phenomena which occur simultaneously : • the anti-solvent effect i.e. the dissolution of the anti-solvent into the solution. The consequence is the swelling of the solution droplets. • the evaporation of the solvent into the anti-solvent. The consequence is the shrinking of the solution droplets.

A key-question must therefore be addressed : do the droplets of solvent in a flow of antisolvent swell or shrink ? Shekunov et al. [START_REF] Shekunov | Optical characterisation and mechanism of antisolvent precipitation in turbulent flow[END_REF] showed results obtained with optical methods, image analysis, laser interferometry and PIV. They showed that droplets of ethanol in a turbulent flow of supercritical CO 2 do not swell as expected but on the contrary seem to retract as long as the ethanol dissolves ("evaporates") in the CO 2 . Depending on the pressure (below or above the critical mixture point) the behaviour of the droplet changes : at low pressures (below the mixture critical pressure) the droplets of ethanol decrease regularly. Nucleation and growth of crystals occur at the droplet interface leading to agglomerated particles. At higher pressures, a liquid ethanol jet appears from each droplet. Nucleation occurs within the shrinking ethanolrich droplets and the mixing regime has a great importance. No swelling have been observed whatever the operating conditions.

Lora et al. [START_REF] Lora | Simulation of the semi-continuous supercritical antisolvent recrystallization process[END_REF] have recently given a simulation of a semi-continuous SAS with CO 2 as anti-solvent, toluene as the solvent and naphtalene and phenanthrene as solutes. They showed that three issues must be addressed : phase equilibria, hydrodynamics and mass transfer. Their calculations showed that the solvent droplets undergo a rapid swelling followed by a slower shrinking phase. CO 2 dissolution in the solvent and solvent evaporation into CO 2 increase with the CO 2 solubility in the solvent : toluene was compared with acetone and DMSO. Crystallisation of the two tested solutes varies : the anti-solvent effect is sufficient to provoke phenanthrene crystallisation while solvent evaporation is involved in naphtalene crystallisation. The variation of the flow rate ratio solvent/anti-solvent is therefore important for the latter solute while it is ineffective on the first solute.

Another study with the CO 2 -toluene system (without a third component to be crystallised) showed the dependence of droplet fate on the respective densities of the solvent and the antisolvent [START_REF] Werling | Numerical Modeling of Mass Transfer in the Supercritical Antisolvent Process: miscible conditions[END_REF]. When the anti-solvent density is lower than that of the solvent -the most common case-the droplets can swell. However, they shrink in the opposite case. The extent of swelling (or shrinking) as well as the lifetime of the droplets are very sensitive to the operating conditions.

To optimise the production of particles from the droplets, and particularly when very fine particles are expected it is advantageous to work in the single-phase zone at high pressures where higher mass transfer and higher supersaturation ratios can be obtained. In this zone, the average lifetime of a droplet is short, implying rapid mass transfer, high supersaturation ratios and finally very small particles. The knowledge of the phase equilibria of the binary system solvent-antisolvent is compulsory before designing a semi-continuous SAS process.

CONCLUSION

Supercritical-assisted particle formation has made a lot of progress in the recent years. Both RESS and SAS processes continue to undergo fundamental and applied research and have benefited from the recent years' advances. Several controversial issues still remain about : the influence of operating parameters on the characteristics of the particle produced (size, morphology, polymorphism), the comprehension of the fluid dynamics, the nucleation phenomenon, the crystal growth under these specific conditions, the particle agglomeration in the jet, the method of scaling up, etc.

Because of its simplicity, RESS will remain the first process to be tested while the universality of SAS will ensure future developments for very different type of materials. Very often, semi-continuous SAS process and RESS process can compete for pre-industrial particle generation. A thorough study based on the phase equilibria and the solubility in pure CO 2 and in CO 2 + cosolvent within a broad range of T and P gives most of the required data for choosing the "right" process.

Figure 1 .

 1 Figure 1. RESS expansion in the pressure-enthalpy diagram for pure CO 2


Klervi Rannou

Gilles Thuret

Philippe Gain

Jean-Charles Pinoli

Yann Gavet
email: gavet@emse.fr

A novel quantitative methodology for age evaluation of the human corneal endothelium

Keywords: quality control, human corneal endothelium, reduced minkowski functionals

The human corneal endothelium regulates the cornea transparency. Its cells, that can not regenerate after birth, form a tesselated mosaic with almost perfect hexagonal cells during childhood, becoming progressively bigger and less ordered during aging. This study included 50 patients (in 10 decades groups) and 10 specular microscopy observations per patient. Five different criteria were measured on the manually segmented cells: area and perimeter of the cells as well as reduced Minkowski functionals. All these criteria were used to assess the probability of age group membership. We demonstrated that the age evaluation is near the reality, although a high variability was observed for patients between 30 and 70 years old.

INTRODUCTION 1.1 Medical context

The cornea is the transparent tissue in the front of the eye. Its two main roles are to refract the light on the retina, and to protect the eye from the outside. It is constituted of three cell layers; among these layers, the innermost called endothelium is dedicated to a constant deswelling of the cornea, keeping it clear during lifetime. The number of endothelial cells (ECs) per surface unit, called endothelial cell density (ECD in cells/mm 2 ) is thus crucial to maintain corneal transparency. Below 300 to 500 cells/mm 2 , the endothelial function is insufficient and the cornea becomes edematous, causing an irreversible vision loss. The ECD of an healthy adult with a normal cornea varies between 2000 and 3000 cells/mm 2 . The corneal endothelium has two particularities: it is a monolayer that can be easily observed with specular microscopes (a routine imaging device first developed in the 60's), and its cells do not regenerate. [START_REF] Stocker | The endothelium of the cornea and its clinical implications[END_REF] The latter is of importance because EC loss caused by several ocular severe diseases or traumatism (accidental or iatrogenic) can result in corneal blindness. This explains why ophthalmologists routinely check the endothelial quality using specular microscopes.

Two parameters are universally used: ECD (the most important) and EC morphology. Cell morphology is assessed only by the polymegethism (variations in cell size, given by the coefficient of variation of cell area = standard deviation/mean) and the pleomorphism 2 (variations in cell shape, namely percentage of hexagonal cells). The ECD is an objective and easy to interpret criterion, but the two others are not intuitive and therefore virtually not used by ophthalmologists. The physiological age-related decrease in ECD has been extensively studied in the past and a mean annual cell loss of 0.6% occurs during adulthood. [START_REF] Bourne | Central corneal endothelial cell changes over a ten-year period[END_REF] Nevertheless inter-individual variations of ECD are important. In the present article we studied the link between the age and a morphological analysis of ECs. We propose a method for predicting the age from the endothelial parameters with a view of detecting outliers liable to indicate a specific health status.

Corneal endothelium observation and segmentation

The images of the corneal endothelium observations have been taken in vivo using a small field non-contact specular microscope (SP3000, Topcon), see Figs.1. For each decade (0 to 10, 11 to 20, ..., 91 to 100) five patients have been observed. They had no history of corneal disease and benefited from a routine corneal examination. For each patient, five images per eye were taken: one in the central part, and the others at 3 to 4 millimeters from the center in the superior, nasal, inferior and temporal parts of the cornea. A total of 500 images have been manually segmented by an expert ophthalmologist in order to avoid the bias induced by an automatic segmentation methods. The segmentation forms closed cells with thin contours.

CHARACTERIZATION OF THE MOSAICS

With some experience, the expert can easily distinguish between normal and pathological endothelial mosaics. In this article, two families of characterizations are presented: the first one will deal with all cells taken separately (called individual cell characterization), and the second one considers the mosaic as a whole (called global mosaic characterization). The measures performed in the first case are the mean ECD and the mean cell perimeter. In the second case, the reduced Minkowski functionals are introduced.

Individual cell characterization

The statistical distribution of the cells areas and perimeters was previously reported. [START_REF] Rannou | Comparison of corneal endothelial mosaic according to the age: The corimmo 3d project[END_REF] It logically showed that the mean ECD decreased with age. In addition, the standard deviation was high after 50 years. The perimeter had the same kind of distributions.

For each age group, the mean and standard deviation is represented in Fig. 2. These graphs show that the mean monotonically increases with regard to the age. Roughly speaking, this defines a bijection between the age and the measure, and consequently a specific measurement allows to evaluate an age.

Global mosaic characterization

The reduced Minkowski functionals 5 are a way to characterize point processes or sets. To perform these measurements in the corneal mosaic, the ECs are simplified to their geometric center, forming a point field. This is mathematically defined as follows. Let Φ be a finite point field within W ⊂ R 2 with N points, Φ r the dilation of Φ by a closed ball b r of radius r > 0. Let A(r), U (r) and χ(r) be the Minkowski functionals of Φ r (area, perimeter and Euler number). The reduced Minkowski functionals a(r), p(r) and e(r), are the normalized Minkowski functionals of Φ r :

a(r) = A(r) πN r 2 (1) u(r) = U (r) 2πN r (2) e(r) = χ(r) N . (3)
An illustration of the method is proposed in Fig. 3.

The Minkowski functionals are measured for all patients, and the mean value for each age group is presented in Fig. 4. These figures clearly show an evolution of each criterion with the age. To quantify this visual evidence, a specific value is introduced. It is a particular radius value defined for each functional (r a , r u and r e for area, perimeter and Euler number, respectively), such that:

∀r ≥ r a , a(r) ≤ 0.3, (4)
∀r ≥ r u , u(r) ≤ 0.3, (5)
∀r ≥ r e , e(r) = 0.

Notice that the value 0.3 is empirically chosen. These values are represented in Fig. 5. These graphs show that the mean values are strictly increasing curves. Again, a bijection can be made between the age group and the mean r value, for the three criteria area, perimeter and Euler number. Three other critera are thus defined to evaluate the age of a cornea.

Evaluation of the corneal age group

So far, five criteria have been presented: mean area and perimeter, and the values r a , r u and r e deduced from the reduced Minkowski functionals. From each criterion, the age group g can be deduced by searching for the closest mean value. Let µ i and σ i be the mean and standard deviation of one criterion, for instance the cells density, of the i-th group. Let m be the same criterion measured on one endothelium (the mean value is computed on the 5 images of the same eye). The estimated age group of the endothelium for this criterion is: Furthermore, an interval of the age groups to which the endothelium could belong is illustrated in Fig. 6 and given by:

g = arg min i=1,...,10 |m -µ i |. (7) • • • • • • • • • • • • (a)
I = [min(E), max(E)] , ( 8
)
where E is:

E = {i = 1, . . . , 10 | m ∈ [µ i -σ i , µ i + σ i ]} . ( 9
)
Notice that in Fig. 6, the values µ represents a continuous extrapolation of the discrete values µ i . The representation of these results (the evaluated age g and the age group interval I) is shown as an example (Fig. 7) for the reduced Minkowski functionals, and only for the 5th age group. This kind of graph is a bit hard to understand, especially if lots of measures are taken into account. Euler number e 7: Age group g and interval I of each endothelium the 5th age group according to the criteria r a , r u and r e . There is not specific ordering in the endothelium indices (in abscissa). The age group 5 (in ordinates) should be included in the interval if the measure would have been perfect.

Probability of age group membership

Let N = 5 be the number of evaluation criteria. For each endothelium e and for each age group i, the numbers n g,i and n I,i are defined as follows:

n g,i = N j=1 1{g j = i} (10) n I,i = N j=1 1{i ∈ I j } (11)
where g j and I j are the estimated age group and interval for the j-th criterion (j ∈ [1; N ]). For the N criteria (j denotes the j-th criterion), n g,i is the number of times the endothelium e is associated to the age group i, and n I,i is the number of times that the interval I j contains the age group i.

From these values, a probability of age group membership is constructed. If G i is the i-th age group (G i = [(i -1) × 10; i × 10[ in years), then the probability of the endothelium e to belong to G i is given by:

P [e ∈ G i ] = 1 (c g + c I ) c a × n g,i N + c I × n I,i N j=1 n I,j , (12)
where c g = 2 and c I = 1 are two positive numbers arbitrarily defined such that more weight is attributed to the age evaluation, and less to the interval.

RESULTS

The probability of group membership is reproduced in the tables 2, 3, 4 and 5. A color scale is applied to enhance the values (see table 1 for the scale). Only a few groups are represented to simplify the presentation of the results.

Except for the patients between 0-30 years, the variability in ECD and the other measured parameters is high (see standard deviations on Figs. 2 and5). This simply reflects the reality observed on patients. The characterization method is a good way to reflect the evolution of the corneal endothelium, compared to an average evolution. The prediction of the age was more efficient for extreme ages (younger, 0-20 years old, and older people, 80-100 years old.) than for the rest of the population. Nevertheless, some old patients may appear much younger from a "corneal endothelium point of view". For example, in the 5-th age group, some patients appear to be classified either in the younger age groups or in the older age groups, denoting a high disparity of evolution of the corneal endothelium.

CONCLUSION AND PERSPECTIVES

As it has been previously stated, this article presents a novel methodology for age evaluation of the corneal endothelium. The results are quite good for these five criteria. Some other criteria are under development and evaluation, and it is planned to published these new results, with an analysis between all these criteria. Furthuer analysis will determine whether the other criteria will efficiently help determining the corneal age.

Using the methodology presented in this article, a second series of hundred patients with one central observation in each cornea will be analyzed. It is also planned to perform a long term study (at least for ten years) with a high number of patients (to be determined) that could link the corneal age with other clinically determined ages.

  (a) Specular image of the corneal endothelium of a 5-year old subject. The cell density is high, and the shapes are regular hexagons. (b) Manual segmentation of the cells of Fig.1a. (c) Specular image of the corneal endothelium of a 92-year old subject. The cell density is low and the shapes are irregular. (d) Manual segmentation of the cells of Fig.1c.

Figure 1 :

 1 Figure 1: Examples of corneal endotheliums observed by specular microscopy.

Figure 2 :

 2 Figure 2: Mean (µ) and standard deviation (σ) of the areas and perimeters of the cells per age group.

  Φ point field. (b) Φr, dilation of Φ with r = 0.5. (c) Φr, dilation of Φ with r = 1. (d) Φr, dilation of Φ with r = 1.5.

Figure 3 :

 3 Figure 3: Point field Φ that schematically represents the geometric centers of the cells, and several dilations with different sizes. The reduced Minkowski functionals are performed on the different binary sets.

Figure 4 :

 4 Figure 4: Mean reduced Minkowski functionals for each age group. The distance in µm is the dilation size.

Figure 5 :

 5 Figure 5: Mean (µ) and standard deviation (σ) of r a , r u and r e per age group.

Figure 6 :

 6 Figure 6: Interval age group definition.

Table 1 :

 1 Color scale employed in the following tables, in percentage.

	0-9% 10-19% 20-29% 30-39% 40-49% 50-59%
	60-69% 70-79% 80-89% 90-99%	100%

Table 2 :

 2 Probability of age group membership for the 1st age group.

ACKNOWLEDGMENTS

The authors thank the French National Research Agency for financial support (ANR-12-TECS-004, CorImMo 3D). Table 3: Probability of age group membership for the 5th age group.