Saeid Nezamabadi

Hamid Zahrouni

Julien Yvonnet

Michel Potier-Ferry

A Multiscale Finite Element Approach for Buckling Analysis of Elastoplastic Long Fiber Composites

come

INTRODUCTION

The development of numerical multiscale analysis has made possible new studies about effective properties of heterogeneous nonlinear materials containing dispersions of multiple phases in the microstructure. Thus, one alternative to meshing the whole structure, including heterogeneities, is to use a multiscale finite element procedure, also called in the literature computational homogenization, or nested finite element procedures (see, among many others, [START_REF] Smit | Prediction of the mechanical behaviour of nonlinear heterogeneous systems by multi-level finite element modelling[END_REF][START_REF] Feyel | 2 Multiscale approach for modelling the elastoviscoplastic behaviour of long fiber SiC/Ti composite materials[END_REF]. In this framework, a nonlinear finite element problem is introduced under each material point (for example, integration points) to obtain the macroscopic stresssstrain relationship, without restriction on the deforma-tion level or the local nonlinear constitutive laws. Each local nonlinear microscopic problem is associated to a representative volume element (RVE). The macroscopic stresses are then obtained by direct averaging of the microscopic stresses, found by solving the local problem [START_REF] Terada | A class of general algorithms for multi-scale analysis of heterogeneous media[END_REF][START_REF] Ghosh | A multilevel computational model for multi-scale damage analysis in composite and porous media[END_REF][START_REF] Feyel | A multilevel finite element method (F E 2 ) to describe the response of highly non-linear structures using generalized continua[END_REF][START_REF] Yvonnet | The reduced model multiscale method (R3M) for the non-linear homogenization of hyperelastic media at finite strains[END_REF].

Nevertheless, modeling multiscale buckling of heterogeneous materials remains a challenging task. The instability phenomena in the framework of multiscale homogenization analysis of heterogeneous materials are of structural and/or material types. These phenomena can occur on both macroscopic as well as microscopic scales and may influence each other. [START_REF] Karlsson | An investigation of localization in a porous elastic material using homogenization theory[END_REF] The present work is devoted to the microbuckling analysis of long fiber composites. A multiscale finite element method (FE2) is combined with the asymptotic numerical method (ANM) to study the elastoplastic instability which may occur in structures at both macroscopic and microscopic scales. The fiber is described by a linear material constitutive law, while the matrix phase is described by a nonlinear Ramberg-Osgood relationship. The stress field is then obtained via the total mechanical strain without any history dependence. Large strains are considered, which induce geometrical nonlinearities in both cases. The ANM framework allows obtaining complex response curves involving limit points in loading and displacement to be obtained. In the present path following procedure, adjustment of the step length is naturally automatic because the validity range of the asymptotic solution is a posteriori estimated depending on the local nonlinearity of the response branches. Numerical examples show the effectiveness of the proposed approach by investigating microscopic and macroscopic instabilities of long fiber composite structures in compression. KEY WORDS: asymptotic numerical method, nonlinear homogenization, multiscale finite element method, long fiber composite, plastic microbuckling phenomena in the heterogeneous materials. They studied the overall behavior of porous elastic media where they found that the homogenized material looses ellipticity, although the matrix material remains elliptic. Afterward, the instability of the heterogeneous materials originated numerous papers. (e.g., [START_REF] Triantafyllidis | On the comparison between microscopic and macroscopic instability mechanisms in a class of fiber-reinforced composites[END_REF][START_REF] Geymonat | Homogenization of nonlinearly elastic materials: microscopic bifurcation and macroscopic loss of rank-one convexity[END_REF][START_REF] Ohno | Microscopic symmetric bifurcation condition of cellular solids based on a homogenization theory of finite deformation[END_REF]Okumura et al., 2002;[START_REF] Miehe | Computational homogenization analysis in finite elasticity: material and structural instabilities on the micro-and macro-scales of periodic composites and their interactions[END_REF][START_REF] Gong | On the stability of kelvin cell foams under compressive loads[END_REF][START_REF] Lopez-Pamies | Second-order estimated for the macroscopic response and loss of ellipticity in porous rubbers at large deformations[END_REF][START_REF] Debotton | Neo-Hookean fiber reinforced composites in finite elasticity[END_REF][START_REF] Michel | Microscopic and macroscopic instabilities in finitely strained porous elastomers[END_REF]Yvonnet et al., 2007;Okumura et al., 2008). Fiber microbuckling is one of the most famous instability phenomena in heterogeneous materials. It is one of the main mechanisms of long fiber composite compressive failure. [START_REF] Rosen | Mechanics of composite strengthening, Fiber Composite Materials[END_REF] presented the first approach of this type and showed that in an elastic context, microbuckling stress is approximately equal to the composite shear modulus. As established by Budiansky (1983), matrix physical nonlinearity and fiber initial waviness have to be considered to correctly predict the failure level. One can find many finite element microbuckling studies in the literature, for instance, [START_REF] Kyriakides | On the compressive failure of fiberreinforced composites[END_REF] and [START_REF] Lee | Compressive response and failure of fiber reinforced unidirectional composites[END_REF]. In [START_REF] Drapier | A nonlinear numerical approach to the analysis of microbuckling[END_REF][START_REF] Drapier | Towards a numerical model of the compressive strength for long fibre[END_REF][START_REF] Drapier | A structural approach of plastic microbuckling in long fibre composites: Comparison with theoretical and experimental results[END_REF], a multiscale approach is considered, but the microscopic model is solved analytically and the assumption of periodicity along the fiber is a bit restrictive. However, to our best knowledge, this paper represents the first attempt to model microbuckling by a multiscale finite element method, which is of fundamental interest for structure analysis of fiber composites. In the present paper, we study the plastic microbuckling of long fiber composites using multiscale finite element analysis (FE 2 ) combined with an asymptotic numerical method (ANM). We consider an initial fiber curvature and a nonlinear behavior of the matrix by choosing a Ramberg-Osgood relation [START_REF] Chen | Plasticity for Structural Engineers[END_REF]. The ANM is based on the expansion of nonlinear problems in the form of power series that are truncated at rather large orders. This allows one to obtain a series of linear problems. Consequently, within the framework of homogenization, one can construct the localization tensors using the superposition principal. These tensors are the same for all the linear problems deduced from the perturbation procedure. Hence, an explicit constitutive relation at the macroscopic level is computed for each asymptotic step. Furthermore, the length of steps is naturally adaptive and is estimated a posteriori using the previously computed terms of the series. This is not the case for the classical iterative algorithms based on Newton-Raphson procedure. This advantage allows one to follow complex response curves.

The layout of this paper is as follows: In Section 2, we present the variational formulations and the coupling equations between the macroscopic and microscopic problems. We also show how to adapt a nonlinear constitutive law to the ANM framework. In Section 3 details of the perturbation procedure applied to the multiscale problem are given. Resolution and path-following strategies are also presented. In Section 4 we present two main numerical examples involving instability phenomena in the long fiber-reinforced composites.

MACROSCOPIC AND MICROSCOPIC VARIATIONAL FORMULATIONS

Equilibrium and Coupling Equations for Micro and Macro Problems

We present in this section the variational equations needed to model the macroscopic and the periodic microscopic problems. Similar formulations have been presented in a previous paper [START_REF] Nezamabadi | A multilevel computational strategy for handling microscopic and macroscopic instabilities[END_REF] in which we considered a heterogeneous material involving large displacement and an elastic and linear constitutive law for the different phases. In the present paper we deal with similar problems but taking into account the plasticity behavior of the matrix at the microscopic level.

Let Ω ba a domain in R d , d being the space dimension associated with a macroscopic structure and ∂Ω its external boundary, both in their reference configuration. In the framework of a total Lagrangian formulation, the weak form of the macroscopic equilibrium equation can be written as follows:

Find ū ∈ S(Ω) satisfying the essential boundary conditions, i.e. ū = ũ on ∂Ω u , with ũ being the prescribed displacements and essential boundary, and S(Ω) the space of sufficiently smooth functions, such that

Ω P : δ F dΩ = λT ext (δ ū) in Ω , ( 1
)
where the notation (.) denotes the macroscopic quantities, P is the first Piola-Kirchhoff tensor associated with a point X of the macroscopic structure in its initial configuration, ū denotes the macroscopic displacement field, F = ∇ X ū + I is the macroscopic deformation gradient tensor, and I is the second-order identity tensor. The term T ext (δū) represents the virtual work of external load, λ is a loading parameter, and δū ∈ S 0 (Ω), S 0 (Ω) being the space of all displacement fields vanishing on ∂Ω u and sufficiently smooth on Ω.

At this scale, the constitutive relationship between P and F is unknown. In the context of multiscale finite element analysis, the macroscopic stresses are obtained by solving a local nonlinear finite element problem related with the periodic microstructure.

For the microscopic scale, the material is assumed to be heterogeneous with a periodic microstructure, characterized by a representative volume element (RVE) that occupies a domain ω ∈ R d in its reference configuration. Let u(X) be a displacement at a microscopic point X in the undeformed configuration. The weak form associated with the microscopic problem is:

Find u ∈ S(ω) satisfying the microscopic boundary conditions, i.e., Eq. ( 5), such that

ω P : δF dω = 0 in ω , ( 2
)
with δu ∈ S 0 (ω), S(ω), and S 0 (ω) being defined as previously. In (2), P and F denote, respectively, the first Piola-Kirchhoff stress tensor and the deformation gradient tensor at a microscopic point X.

Problems ( 1) and ( 2) are coupled through two main relations. First, as the constitutive relation is not explicitly given at the macroscopic level, the effective stress tensor P is obtained by considering an average value of the microscopic stress field over the RVE for a given point X of the macroscopic structure. This relation is expressed as follows:

P = P(X) = 1 |ω| ω P(X)dω , ( 3
)
where |ω| represents the volume of the considered RVE.

The second relation concerns the mean value of the microscopic deformation gradient assumed in the form

F = F(X) = 1 |ω| ω F(X)dω . ( 4
)
This relation is deduced from the boundary conditions imposed on the RVE [START_REF] Miehe | Computational micro-to-macro transitions for discretized micro-structures of heterogeneous materials at finite strains based on the minimization of averaged incremental energy[END_REF]. Generally, three main boundary conditions can be considered for the RVE as linear deformations, uniform tractions, or periodic constraints. In the proposed work, we consider periodic conditions on the boundary of the RVE, which can be recast into the following form:

u + -u -= ( F -I)(X + -X -) on ∂ω , ( 5
)
where u is the microscopic displacement, and the exponents + andare associated with node indices on opposite sides of the RVE. Note that the boundary conditions depend on the macroscopic deformation tensor F.

Microscopic Constitutive Relations

At the microscopic scale, we assume that the constitutive relations are known in each phase of the RVE. For the applications that we target in the present work, a linear constitutive relation is considered for the fiber microconstituent and a nonlinear constitutive law is assumed for the matrix. In fact, we consider plastic behavior for the matrix without taking into account elastic unloading. Deformation theory of plasticity is used, which is convenient for problems where the physical nonlinearity is more important than the effect of irreversible process and the history of the loading [START_REF] Hencky | Zur theorie deformationen und der hierdurch im material hervorgerufnen nebenspannungen[END_REF]. Let us remark that Eq. ( 2) can be expressed in terms of the second Piola-Kirchhoff stress field S and the Green-Lagrange strain γ as follows:

ω S : δγ dω = 0. (6)
The stress tensors S and P are related through P = FS and the Green-Lagrange strain is given by

γ = 1 2 (F T F -1) = 1 2 ∇ X u + ∇ X u T + ∇ X u T ∇ X u . ( 7
)
The linear constitutive relation for the fiber can be expressed as

S = C : γ, (8)
where C refers to the classical fourth-order elastic tensor.

Concerning the matrix, we choose an elastoplastic constitutive law based on the Ramberg-Osgood relation [START_REF] Chen | Plasticity for Structural Engineers[END_REF]Abaqus Theory Manual, 2008;[START_REF] Zahrouni | Asymptotic numerical method for nonlinear constitutive laws[END_REF], which is written in the 3D case in the following form:

Eγ =(1+ν)S d -(1-2ν)P I+ 3 2 α S eq σ y n-1 S d , (9)
where E, ν, α, n, and σ y denote, respectively, the Young modulus, Poisson's ratio, yield offset, hardening component, and yield stress. P = -1/3S : I is the equivalent hydrostatic stress, S d is the stress deviator defined by S d = S + P I. S eq is the von-Mises equivalent stress defined as follows:

S eq = 3 2 S d : S d . ( 10
)
As Eq. ( 9) is highly nonlinear, it is convenient to set it into a regular and quadratic form which is well desired in the framework of ANM [START_REF] Potier-Ferry | Traitement des fortes non-linéarités par la méthode asymptotique numérique[END_REF][START_REF] Zahrouni | Asymptotic numerical method for nonlinear constitutive laws[END_REF][START_REF] Abichou | Asymptotic numerical method for problems coupling several nonlinearities[END_REF]Zahrouni et al., 2004). As the hardening exponent n is not an integer, this law is not analytic for null stress. This is why a regularization procedure and a differential relation are introduced in this model, allowing expansions into power series [START_REF] Zahrouni | Asymptotic numerical method for nonlinear constitutive laws[END_REF]. For this purpose, the von-Mises equivalent stress is slightly modified and recast in the following form:

S 2 eq = 3 2 S d : S d + η 2 σ y 2 , ( 11
)
where η is a regularization parameter. For η = 0, the initial constitutive law ( 9) is recovered. Furthermore, two additional variables κ and ζ are introduced to transform the power law into a differential equation:

κ(S eq ) = 3 2 α S eq σ y n-1 = 3 2 α 3 2σ y 2 S d : S d + η 2 n-1 2 , ( 12
)
ζ 2 = S 2 eq σ 2 y = 3 2 σ y 2 S d : S d + η 2 . ( 13
)
These two variables are linked as follows:

κ = 3 2 α ζ n-1 . ( 14
)
If we carry out a differentiation of Eq. ( 14), one obtains a relation which is more convenient for the asymptotic expansions:

ζ dκ = (n -1) κ dζ. ( 15
)
In this way we obtain a general problem with a quadratic nonlinearity with respect to the new variables. Moreover, to keep the same initial slope as for nonregularized law, the first member of Eq. ( 9) is multiplied by (1 + αη n ) [START_REF] Zahrouni | Asymptotic numerical method for nonlinear constitutive laws[END_REF]. The nonlinear constitutive law leads then to the five following equations:

E(1+αη n )γ =(1+ν)S d -(1-2ν)P I +κS d , ζ 2 = 3 2 σ y 2 S d : S d + η 2 , ζ dκ = (n -1)κdζ, (16)

ASYMPTOTIC NUMERICAL ALGORITHM

The micro and macro formulations are described by the nonlinear system of Eqs. (1,(5)(6)(7)(8)16). The regularization procedure and the differential form introduced in the nonlinear constitutive law have allowed a regular and quadratic framework to be obtained, which is preferred using ANM algorithms. The present section gives some details about the path following procedure within ANM.

We show that the perturbation procedure allows one to transform the nonlinear problem into a sequence of linear ones. Because of coupling equations, at each asymptotic order, the solution of the resulting linear problems is obtained by using the superposition principal. Details of this procedure are given in [START_REF] Nezamabadi | A multilevel computational strategy for handling microscopic and macroscopic instabilities[END_REF] and will be explained once more in what follows.

Perturbation Technique

Solution of the nonlinear multiscale problem [Eqs. (1,(5)(6)(7)(8)16)] is sought using the asymptotic numerical method. This consists of expanding the main macroscopic and microscopic variables of this problem into a power series with respect to a path parameter. In this way, the nonlinear problem is transformed into a sequence of linear ones having the same tangent operator. If these variables are represented by a mixed vector Λ = (ū, P, u, P, ...), the perturbation technique applied to Λ and to the load parameter λ leads to

Λ(a) = Λ 0 + N p=1 a p Λ p , ( 17
)
λ(a) = λ 0 + N p=1 a p λ p , ( 18
)
where a is a path parameter to be defined, N is the truncation order of the series, and Λ 0 and λ 0 denote a known initial solution. By substituting Eqs. ( 17), ( 18) into (1, 5-8, 16), one obtains a sequence of linear problems. For example, for order p, Eq. ( 1) is written as:

Ω Pp : δ F dΩ = λ p T ext (δū) . ( 19
)
Concerning the microscopic variables, if we set the Green-Lagrange strain γ into the following form:

γ = 1 2 ∇ X u + ∇ X u T + 1 2 ∇ X u T ∇ X u = γ l (u) + γ nl (u, u), (20)
and

δγ = γ l (δu) + γ nl (δu, u) + γ nl (u, δu) = γ l (δu) + γ s nl (u, δu) . ( 21
)
Equation ( 6) can then be written for order p as ω (S p : (γ l (δu) + γ s nl (u 0 , δu))

+S 0 : γ s nl (u p , δu))dω = 0. ( 22
)
Boundary conditions at the microscopic level are expressed at order p by

u + p -u - p = Fp (X + -X -) on ∂ω . ( 23
)
Constitutive law is defined for each phase of the microscopic structure. In the case of the linear law (8), order p leads to the following relation:

S p = C : γ p , (24)
and in the case of the nonlinear constitutive law ( 16), order p gives

S p = C t : γ p + S res p , (25)
where

γ p = γ l (u p )+γ s nl (u 0 , u p )+ p-1 i=1 γ nl (u i , u p-i ), (26)
and C t and S res p are defined in Appendix A.

Resolution Strategy

We have presented in Section 3.1 the perturbation technique which transforms the nonlinear problem to be solved into a sequence of linear ones. The aim of the present section is to give some details on how the resulting linear problems are solved in the context of a multiscale framework. Because of the coupling conditions between macroscopic and microscopic problems, we start by solving partly the linear problems at the microscopic level for each order. The obtained solution allows one to construct the localization tensor numerically and then to obtain the macroscopic tangent modulus. At this stage, one can solve the linear problem of the macroscopic equilibrium to obtain the displacement field ūp . This displacement is used to compute the macroscopic deformation gradient Fp and then to finish the computation of the variables at the microscopic level.

The linear problems at order p of the microscopic level (5, 6, 16) can be set in a generic form as follows:

L(u p , δu) = F nl p (δu) in ω, (27)
u + p -u - p = Fp (X + -X -) on ∂ω. ( 28
)
The operator L corresponds to the tangent stiffness of the microscopic structure. It is the same for all the linear problems resulting from the perturbation procedure. The right-hand side operator F nl p is given and depends only on the solutions of the previous linear problems i (1 ≤ i ≤ p). Details about L and F nl p are given in Appendix B. Note that the boundary conditions (28) depend on the local macroscopic deformation gradient, Fp , which is unknown. As the problems ( 27) are linear, we can use the superposition principle. Therefore, u p can be expressed as a linear combination of the solution modes obtained by imposing independently homogeneous deformation fields (see Appendix C). The solution u p can be then expressed in the form

u p (X) = A(X) : Fp + u nl p (X), (29)
where A is a third-order tensor deduced from solutions of the linear problems obtained by the superposition principle, and u nl p is a vector associated with the right-hand term of Eq. ( 27).

Deriving (29) with respect to X and using ∇ X u p = F p , we obtain

F p = A ,X : Fp + u nl p,X , (30)
where A ,X is a fourth-order tensor, which is identified as a localization tensor, and u nl p,X is a second-order tensor. To achieve the homogenized relationship, we consider the asymptotic expansion of P = FS. After some elementary operations, one obtains at order p the following relation:

P p = L : Fp + P nl p , ( 31
)
where expressions of operators L and P nl p are detailed in Appendix D. The effective stresses can be then obtained at each order p by averaging Eq. ( 31):

Pp = L : Fp + Pnl p . ( 32
)
By introducing (32) into Eq. ( 19), we obtain the final form of the linear problems at order p at the macroscopic level:

Ω Fp : L : δ FdΩ= λ p T ext (δū)- Ω Pnl p : δ FdΩ. ( 33
)
Solution ūp of Eqs. ( 33) and ( 29) are used to finish the computation of the variables at the microscopic level for each integration point of the macroscopic structure.

Path-Following Strategy

The path-following strategy is of critical importance, mainly in the presence of instability phenomena. This is the case of the present work where instability may occur in both macroscopic and microscopic levels. Nonlinear problems involving instability behavior are generally solved by using Newton-Raphson iterative algorithms associated with arc-length techniques to follow complex branches involving limit points in loading and displacement responses. To optimize the computation time, one has to use step-size adjustments. Indeed, for a fixed step size, a very small step length can lead to a high number of steps; however, a large step length allows reduction of the number of steps but can lead to divergence or to a very large number of iterations. Furthermore, the linear iterative algorithms converge sometimes to an unwanted path in the case of problems involving strong nonlinear solution response [START_REF] Eriksson | On step size adjustments in structural continuation problems[END_REF]. The best strategy must use large step lengths in the regions with very smooth response and small step lengths in the regions with strong nonlinearity, as for limit points or quasi-bifurcation situations. In the general case, the best strategy is not obvious [START_REF] Ramm | Strategies for tracing the nonlinear response near limit points[END_REF][START_REF] Riks | Some computational aspects of the stability analysis of nonlinear structures[END_REF][START_REF] Wriggers | A quadratically convergent procedure for the calculation of stability point in finite element analysis[END_REF][START_REF] Eriksson | Derivatives of tangential stiffness matrices for equilibrium path descriptions[END_REF][START_REF] Crisfield | Non-Linear Finite Element Analysis of Solids and Structures[END_REF][START_REF] Battini | A modified corotational framework for triangu-lar shell elements[END_REF][START_REF] Kouhia | Stabilized forms of orthogonal residual and constant incremental work control path following next term methods[END_REF].

In the framework of ANM, the power series allows one to obtain a large part of the solution branch by decomposing only one tangent stiffness matrix. This procedure can be considered as a high-order predictor which generally does not need any correction phase. As the series has convergence radius, a continuation procedure can be easily implemented, considering that the relative difference between the solutions of two consecutive orders N and N -1 must remain small enough as compared to a critical value δ. If we consider the discretized form of the displacement field, one can determine a validity range of the truncated series using the following approximation:

a max = δ ū1 ūN 1 N -1 , ( 34
)
where δ is the accuracy parameter and . refers to the norm of the vector. In this way, the whole solution is obtained in a stepby-step manner considering each end step as the starting point for the next one. Furthermore, the step length is naturally adaptive, depending upon the local nonlinearity of the response curve, since it is a posteriori computed using the terms of the series. Details of these procedures are given in many references [START_REF] Cochelin | A path following technique via an asymptotic numerical method[END_REF][START_REF] Zahrouni | Computing finite rotations of shells by an asymptotic-numerical method[END_REF][START_REF] Cochelin | Méthode Asymptotique Numérique[END_REF]. This technique to compute the validity range is not unique; [START_REF] Assidi | Regularization and perturbation technique to solve plasticity problems[END_REF] proposed a procedure which consists in computing the validity range for all the variables expanded into power series and keeps the smallest value of a max . This ensures that all the equations of the considered problem are solved with a high accuracy. To improve the validity range of the asymptotic solution, we propose that the power series representation be replaced by rational fractions named Padé approximants. Experience shows that this technique allows one to reduce the computation time by about 50%. Details of this technique are not given in this work. Readers can refer to the following references: [START_REF] Najah | A critical review of asymptotic numerical methods[END_REF][START_REF] Elhage-Hussein | A numerical continuation method based on Padé approximants[END_REF][START_REF] Lahmam | High-order predictor-corrector algorithms[END_REF][START_REF] Boutyour | Bifurcation points and bifurcated branches by an asymptotic numerical method and padé approximants[END_REF][START_REF] Aggoune | High-order predictioncorrection algorithms for unilateral contact problems[END_REF][START_REF] Aggoune | Asymptotic numerical methods for unilateral contact[END_REF] At last, to solve the resulting linear problems ( 33) and ( 27), we need to define a constraint equation as for the path-following technique proposed with the classical Newton-Raphson algorithms. The computation can be managed by imposing a force, a displacement, or a combination of both [i.e., Riks technique [START_REF] Riks | Some computational aspects of the stability analysis of nonlinear structures[END_REF]]. In the present work, a simple procedure is chosen by considering the projection of the macroscopic displacement field onto its tangent: a = (ū(a) -ū0 ) .ū 1 .

(35)

Note that Eq. ( 35), which uses the macroscopic variable ū, can be replaced by anyone using local variables as the displacement vector of the discretized RVE considered at an integration point of the structure. We can also combine many variables which do not exhibit the same behavior during the macroscopic loading.

NUMERICAL EXAMPLES

The accuracy and the efficiency of the proposed technique were studied in the elastic framework in [START_REF] Nezamabadi | A multilevel computational strategy for handling microscopic and macroscopic instabilities[END_REF]. The influence of microbuckling of the fiber with the initial imperfection on the macrostructure is assessed here. For this purpose, first, the plastic microbuckling mechanism is described briefly. Afterward, the compression of a fiber-reinforced composite material with the elastoplastic nonlinear constitutive law in the matrix is considered. The proposed problem has been discretized using a two-dimensional finite element in the plane stress framework. The macrostructure has been meshed with four node quadrangular elements, whereas the microstructure has been meshed with nine node quadrangular elements.

In the context of ANM, the truncation order N of the series and the accuracy parameter δ, which allows limiting the length of each asymptotic step in a continuation procedure, are important. In all the tests presented here, we have solved the nonlinear problems using ANM with Padé approximants and N = 15, but δ = 10 -6 for the first test and δ = 10 -8 for the second one.

Plastic Microbuckling Mechanism

In this part we focus on the plastic microbuckling mechanism of the long fiber-laminated composite. A bidimensional representation of a laminate is used (see Fig. 1), where e 2 is the 90 • direction corresponding to the loading direction. The imperfection magnitude (v 0 ) is constant through the thickness. This mechanism was studied in several papers [START_REF] Kyriakides | On the compressive failure of fiberreinforced composites[END_REF][START_REF] Drapier | A nonlinear numerical approach to the analysis of microbuckling[END_REF][START_REF] Drapier | Towards a numerical model of the compressive strength for long fibre[END_REF][START_REF] Drapier | A structural approach of plastic microbuckling in long fibre composites: Comparison with theoretical and experimental results[END_REF]. Here we wish to revisit this mechanism by considering a microstructure with one fiber using the periodic boundary conditions at microscopic scale under plane stress conditions. The imperfection wavelength of a fiber is equal to 200.π μm. The mechanical properties employed are shown in Table 1 which corresponds to a T300/914 composite. The parameters of the nonlinear behavior of the matrix are deduced from the work presented in [START_REF] Wisnom | The effect of specimen size on the bending strength of unidirectional carbon fibre-epoxy[END_REF]. The microstructure is de-

Fiber T300

Matrix 914 (isotropic) (isotropic)

E f = 240 Gpa E m = 4.5 Gpa ν f = 0.3 ν m = 0.4 Fiber volume fraction: f =0.6 n = 6 Fiber diameter: d f = 10μm σ y = 115 MPa
formation driven by a macroscopic compression mode;

F -I = F11 -1; F12 ; F21 ; F22 -1 = [0; 0; 0; -1] [see Eq. (5)].
Figure 2 shows the responses under compression mode for different imperfection magnitudes. These results are similar to ones obtained by [START_REF] Kyriakides | On the compressive failure of fiberreinforced composites[END_REF]. The stable and unstable responses toward prescribed strain can be distinguished. In the microstructure, the path following is well adopted and we can follow the response curve with a reasonable number of steps, despite the abrupt change of curvature at maximum stress.

Interaction between Micro-and Macro-Instabilities

To show the influence of plastic microbuckling on the macrostructure response, we consider the macrostructure shown in Fig. 3 which presents a rectangular plate clamped on three edges and submitted to a force distribution λP on the top edge. Because of the symmetry of the problem, only half of the structure is discretized. We use the same microstructure as for the previous example. The imperfection magnitude of fiber is 3 μm. We use our multiscale procedure to solve this micro-macro problem. In Fig. 5, the displacement of a point located on the top of a macroscopic element shown in Fig. 4 versus the loading parameter is presented. This figure shows the influence of microscopic instability on the macroscopic one: initially, in the portion of the curve preceding point A, the macroscopic behavior is mainly linear. However, we can observe an abrupt change in the portion of the curve between point A and point B. Note also the same phenomenon between the points B, C, D, and E. This is mainly caused by the occurrence of microscopic instabilities at different points. It is worth noting that thanks to the high order predictor of ANM, we can detect these instabilities. The deformed shapes of the microstructure at the integration point IP I shown in Fig. 4 at points A, B , C, D, and E (see Fig. 5) are presented in Fig. 6. There is a sig- nificant difference between these deformed microstructural shapes despite the small variation of macroscopic load, which illustrates the relationship between local and global instabilities. Furthermore, we observe that the microscopic buckling modes are the same for all points after the buckling of macrostructure (B, C, D, and E) which shows that we follow the same bifurcation branch.

In Fig. 7, the macroscopic stress-strain diagrams are presented for the different integration points shown in imum stress at the same time, which explains the shape of the load-displacement curve near the maximum load (see Fig. 7). This phenomenon consequentially induces a reduction of asymptotic step length. Nevertheless, ANM allows this type of response curves to be obtained, despite the appearance of successive instabilities.

CONCLUSION

The present paper has been devoted to study the interaction between microscopic instabilities and macroscopic ones in the framework of long fiber composites.

We have considered an elastic linear behavior for the fiber and a nonlinear constitutive relation for the matrix based on the total deformation theory. To this end, a Ramberg-Osgood relationship was implemented. Initial curvatures were taken into account to model geometric imperfections of the fiber. Nonlinear problems were solved using an efficient numerical technique, in the context of heterogeneous materials, combining multiscale finite element analysis (FE 2 ) and an asymptotic numerical method (ANM). Using ANM permits one to construct only once the localization and homogenization tensors for each asymptotic step. Consequently, the tangent stiffness matrix is also the same for each step. The length of each step is a posteriori estimated from the computed terms of the series. This leads to a naturally adaptive step length algorithm which allows investigations of instabilities phenomena with high accuracy and efficiency.

As numerical examples, we have shown first the plastic microbuckling mechanism by considering the compression mode of a microstructure in which we have consid- ered only one fiber with initial imperfection. We have then studied the compression of a fiber-reinforced composite with the microstructure defined above to assess the influence of the plastic microbuckling on the macrostructure behavior. The obtained results illustrate the capability of the proposed algorithm to deal with instabilities in the context of heterogeneous materials involving nonlinear constitutive relations. To our knowledge, a similar micro-macro analysis had never been proposed to account for the compressive failure of long fiber composites. Most of the numerical studies in the literature concern microscopic modeling that is more or less similar to that in Section 4.1 [see [START_REF] Kyriakides | On the compressive failure of fiberreinforced composites[END_REF], and [START_REF] Lee | Compressive response and failure of fiber reinforced unidirectional composites[END_REF]]. In [START_REF] Drapier | A structural approach of plastic microbuckling in long fibre composites: Comparison with theoretical and experimental results[END_REF], the material properties were homogenized only transversally and not in the fiber direction, which permitted to account for structural effects across the thickness of a plate, but this analysis could not be coupled to any macro-structural problem, contrary to the one proposed in this paper.

However, accumulation of asymptotic steps is often observed which may be explained by the multiple instabilities and the loss of ellipticity close to the maximum load. To overcome this issue, one possible way of extension of the method is the use of second-order computational homogenization which allows obtaining in this context a well posed problem at the macroscopic level [START_REF] Michel | Microscopic and macroscopic instabilities in finitely strained porous elastomers[END_REF]Kouznetsova et al., 2004).

where B is a fourth order tensor depending on the initial solution u 0 . Now let us express S p from equation ( 24) for the linear law and equation ( 25) for the nonlinear constitutive law. In these equations, we need to identifify γ p . Hence, we expand the Green-Lagrange strain (20) that gives at order p:

γ p = 1 2 ∇ X u p + ∇ X u T p + ∇ X u T 0 ∇ X u p + ∇ X u T p ∇ X u 0 + γ nl p , (B.4)
where γ nl 1 = 0, As each problem at order p is linear, its solution u p (X) can be decomposed into modes associated with different boundary condition responses. Here we give an illustration in the 2D case, though the methodology can be readily extended to three dimensions. We can write where ũ(ij) (X) are the solutions of the following problems:

γ nl p = 1 2 p-1 i=1 ∇ X u T i ∇ X u p-i ,
L(ũ (ij) , δu) = 0 in ω, ũ(ij)+ -ũ(ij)-= X (ij)+ -X (ij)-on ∂ω, (C.3) with X (11)+ -X (11)-= 1 0 0 0 (X + -X -), X (12)+ -X (12)-= 0 1 0 0 (X + -X -), X (21)+ -X (21)-= 0 0 1 0 (X + -X -), X (22)+ -X (22)-= 0 0 0 1 (X + -X -), (C.4) and u nl p (X) is the solution of the problem:

L(ũ (ij) , δu) = F nl p (δu) in ω.

(C.5)

The solution u p (X) can thus be expressed by: u p (X) = A(X) : Fp + u nl p (X), (C.6)

where A(X) is a third-order tensor defined by A ijk = ũ(jk) i .

APPENDIX D.

We expand P = FS which leads to:

  S d = S + P I,

FIG. 1 :

 1 FIG. 1:The microstructure of long fiber composite with the imperfection magnitude.

FIG. 2 :

 2 FIG. 2:Macroscopic stress-strain diagrams for the different imperfection magnitudes.

Fig. 4 .FIG. 5 :FIG. 6 :

 456 FIG. 4: Mesh topology of the macrostructure

FIG. 7 :

 7 FIG. 7:The macroscopic stress-strain diagram for the integration points (I, II, and III) shown in Fig.4.

S

  ∀p =2, . . . , N.(B.5) As γ p is linear with respect to ∇ X u p , equation (B.4) can be rewritten in the form:γ p = B : ∇ X u p + γ nl p = B : F p + γ nl p . (B.6)Then S p is expressed by:S p = C (r) : B : F p + S nl p , (B.7)where C (r) is equal to C for the linear constitutive law and C t for the nonlinear relation, and the residual stress vector, S nl p , is defined in the case of the linear law as:S nl p = C : γ nl p , (B.8)and in the case of the nonlinear constitutive law as:S nl p = C t : γ nl p + S res p . (B.9)By substituting Eq. (B.7) in (B.1), we obtain Eq. (27). The explicit forms of the operators involved in the resolution of the sequence of linear problems in equation (27) are given by: L(u p , δu)= Ω S 0 : δγ p +F p : B : C (r) : δγ 0 dΩ, (Bi : δγ p-i +S nl p : δγ 0 dΩ. (B.11) APPENDIX C.

P 1 =

 1 F 0 S 1 + F 1 S 0 , P p = F 0 S p +F p S 0 + S p-i = P * nl p , we can express P p as:P p = F 0 C (r) : B : F p +S nl p +F p S 0 +P * nl p , (D.3)which can be rewritten in the form:P p = H (r) : F p + F 0 S nl p + P * nl p . (D.4) Using (30), we obtain: P p =H (r) : A ,X : Fp +H (r) : u nl p,X +F 0 S nl p +P * nl p . (D.5) By setting P nl p = H (r) : u nl p,X + F 0 S nl p + P * nl p and L = H (r) (X) : A ,X (X), we finally obtain: P p = L : Fp + P nl p . (D.6)

  have investigated for the first time the instability Université Paul Verlaine de Metz, Laboratoire de Physique et Mécanique des Matériaux, FRE CNRS 3236, Ile du Saulcy 57045, Metz Cedex 01 France 2 Université Paris-Est, Laboratoire Modélisation et Simulation Multi Echelle, FRE CNRS 3160, 5 Bd Descartes, 77454 Marne-la-Vallée Cedex 2, France

	A MULTISCALE FINITE ELEMENT APPROACH FOR
	BUCKLING ANALYSIS OF ELASTOPLASTIC LONG
	FIBER COMPOSITES
	Saeid Nezamabadi 1 , Hamid Zahrouni 1 , Julien Yvonnet 2
	& Michel Potier-Ferry 1
	1

TABLE 1 :

 1 Mechanical characteristics for the T300/914 composite.

APPENDIX A.

In the preceding sections, we have developed the asymptotic numerical algorithm for problems involving nonlinear constitutive laws. By comparison with the elastic case, we have replaced the elastic stiffness tensor C by the tangent modulus one C t to compute the tangent stiffness matrix. Furthermore, because of this nonlinear law, stress field involves, at order p, residual vector S res p . This latter leads to additional terms in the right hand side F nl p . Hence, at order 1, we can write the constitutive relation in the following form:

One can show easily that this equation is provided from the expression detailed below:

3) where is the fourth-order identity tensor and C 1 , C 2 and C 3 are constants defined as follows:

where

At order p, the constitutive relation is expressed as

Matrix C t is the same at each order; however, residual stresses are expressed at order p as follows:

where

APPENDIX B.

Expansion of (6) gives The tensors δγ 0 and δγ p are given by expressing δγ [Eq. ( 21)]:


Determination of the crop coefficient (Kc) for optimal use of water in agriculture

In a global context of climate change and demographic pressure, the world is facing the scarcity of water resources. This includes the optimal use of this resource in agriculture.

The practice of irrigation and drainage is one of the strategies to be adopted to cope with and adapt to this threat. The determination of the crop coefficient is very important in sizing an irrigation system. Many crops in Africa, especially basic ones (cassava, millet)

and other medicinal plants have not yet seen their Kc determined.

Irrigation is the action of providing a plant with a quantity of water that can compensate for the water deficit in the growth of plant.

Requirement: specific water need of plant.

This water need is represented by the evapotranspiration operated by the plant (ETM)

Experimental determination of Crop ETM by using daily water balance during crop growth

For water monitoring : More than two lysimeters will be used for simultaneous measurement of:

• irrigation (water supply)

• deep percolation (drainage) • evaporation

Experimental determination of ET0 by installing a meteorological station

A complete weather station will be used to measure daily the parameters used in the Penman-Monteith calculation of ET0.

These parameters are :

• the minimum and maximum temperatures,

Experimental determination of ET0 by installing a meteorological station

The Penman-Monteith calculation of ET0.

(5)

  felmaks.one@gmail.com Presented at the RSPH Conference : « What Is the Future of Water and Public and Health? »,

  , R = Runoff D = drainage,The plant will be fed by irrigating without rain (P = 0), hence ΔS = 0The system will be installed to prevent runoff (R = 0

Fig

  Fig. 1: Lysimeter

•

  Fig. 3: Localized and sprinkling irrigation


Burkina Faso, a poor country located in the Sudano Sahelian zone, is progressively making progress in terms of access to water, energy and agricultural development.

On the other hand, the DRC, a naturally rich country located in the most watered part of Africa, still faces the challenges of access to drinking water, weak development of agriculture, ... After a quick analysis, we realize that Burkina Faso is investing a lot of resources in water engineering. What are water engineering applications in which DRC can invest?

• https://www.populationdata.net/pays/republique-democratique-du-congo/, on January 11th, 2018

• https://www.populationdata.net/pays/burkina-faso/, on January 12th, 2018

• http://www.banquemondiale.org/fr/results/2013/04/12/water-sanitationresults-profile, on January 13th, 2018

• http://oneabf.com/notre-projet-dentreprise/, on January 13th, 2018 • Mortality-rate: 13,27 ‰ (2015)

• Child mortality rate: 88,62 ‰ (

• 53% as rate of access to drinking water (83% at Kinshasa)

•

Access to electricity: 13,5%

• Agriculture=less than 10% of GDP • ...

BURKINA FASO STRATEGIES

Many studies focused mainly on water and energy:

SECTORS AND APPLICATIONS OF WATER AND ENVIRONMENTAL ENGINEERING CHALLENGE FOR DRC

Wastewater treatment plant


Bangalore Ravi Kiran

Bogdan Stanciulescu
email: bogdan.stanciulescu@mines-paristech.fr

Jesús Angulo
email: jesus.angulo@mines-paristech.fr

Unsupervised clustering of hyperspectral images of brain tissues by hierarchical non-negative matrix factorization

Keywords: Hyperspectral Image, Hierarchical clustering, Brain tissue, Non-negative matrix factorization

Hyperspectral images of high spatial and spectral resolutions are employed to perform the challenging task of brain tissue characterization and subsequent segmentation for visualization of in-vivo images. Each pixel is a high-dimensional spectrum. Working on the hypothesis of pure-pixels on account of high spectral resolution, we perform unsupervised clustering by hierarchical non-negative matrix factorization to identify the purepixel spectral signatures of blood, brain tissues, tumor and other materials. This subspace clustering was further used to train a random forest for subsequent classification of test set images constituent of in-vivo and ex-vivo images. Unsupervised hierarchical clustering helps visualize tissue structure in in-vivo test images and provides a inter-operative tool for surgeons. Furthermore the study also provides a preliminary study of the classification and sources of errors in the classification process.

Introduction

Medical Hyperspectral imaging (MHSI) is a noninvasive and non-ionizing modality used in robust tissue identification, characterization achieved reasonable classification and localization accuracy in of various kinds of cancerous tissues on account of vast improvements in spectral and spatial imaging resolutions (LF14). Discriminating cancerous and normal tissues with Near Infra-Red (NIR) spectroscopy has been an active topic in the past two decades. (PYF + 07) study the use of MHSI for identifying residual tumor in a resection bed and to indicate regions requiring more resection. (LHW + 14) use HSI in-vivo to achieve tumor identification in tumor bearing mice. (LWL11) study localization of tumor tissue on the tongue of human patients using a tunable HSI camera. In this work we study the characterization and identification of brain tissues by performing joint clustering and un-mixing. We also demonstrate the utility of the hierarchical clustering as a visual aid for inter-operative(IoP) procedures (SC03) in the resection of cancerous tissue. Though one should note that the standard in IoP brain tumor surgery modality has been MRI, it is precise but a costly technology. Hyperspectral images are cheaper and furthermore faster to acquire.

The primary goal of this paper is to address the common problem of lack of an exact ground truth labels in the hyperspectral images in-vivo. This produces a challenging problem of hyperspectral image segmentation and classification. Our solution consists in generating unsupervised clustering that can serve as a stable identification and labeling, for a subsequent semi-supervised classification. Further we note that, it is unrealistic to obtain a well localized marker of tumors, in-vivo, during a surgical procedure, in view of different constraints of space and time. Finally we also do not know the exact spectrum of cancerous and normal tissue at this stage of our study. We do note that from previous studies on different tumor types, normal and cancerous tissue spectra are close (see subsection 4), and are difficult to distinguish. We thus characterize the tissue by clustering their spectrum in similar subspaces. This helps visualize different structures in the tissue as well as locate tumors and related structures. Unsupervised clustering is rarely followed, for such problems, but given that we operate at high spatial and spectral resolutions, we stand to gain some useful information. We use the clustering of ex-vivo tumors with surrounding tissue, blood and other structures as a target model we predict the location of tumors in in-vivo images. Finally, the clustering will also be used as tool by surgeons and pathology to label the tissues to generate ground truth data.

Figure 1: We classify the input images into three categories: In-vivo images (during surgery), Ex-vivo images (extracted tumor after surgery) and material cubes (carefully windowed regions containing only surgical material). Black plastic rings markers highlight the locations of healthy tissue and a tumor tissue in a scene. The ex-vivo images serve as localized target that ensure the presence of tumor tissue in a localized window. We partition the operations into a training and test set, after which we perform hierarchical clustering on the test images. The H2NMF clustering algorithm also calculates the endmembers for each cluster. We use clusters at a particular level to train a random forest classifier. We show the classifier prediction on other ex-vivo and in-vivo images in the test set. The flow-diagram at the bottom illustrates all the steps in our study.

The primary contributions of the paper include the identification of pure-spectrum or pure-pixels of materials using VNIR camera images of brain tumors by applying the hierarchical clustering by 2-rank NMF approximation (GKP15). Pure-pixels are found at the leaves of this hierarchical structure and correspond to rank-1 matrices. This is discussed in section 2. We use the unsupervised clustering as supervised labels to identify these pure-pixels in the in-vivo images by training a decision forest. See figure 1 for a summary of the work-flow. The organization of the paper is as follows. We describe briefly the imaging setup in this section. A brief introduction to endmember extraction is provided in section 2, the separable non-negative matrix factorization(NMF) and the hierarchical NMF decomposition algorithm (GKP15) employed in this study.

Imaging Setup: The medical imaging demonstrator consists of two camera pairs covering different parts of the spectrum at different spatial resolutions, a VNIR (visible-to-near-infrared) and NIR (near-infrared cameras) cameras. The VNIR camera captures 826 spectral bands, in the spectral range of 400-1000 nm with a spectral resolution of 2-3 nm, and each pixel has a dimension of 0.1287 mm x 0.1287 mm. While the NIR camera captures 172 spectral bands, in the spectral range of 900-1700 nm, each pixel has a dimension of 0.48 mm x 0.48 mm. We shall work with the VNIR camera images mainly since they consist of precise spectral information. The increment in wavelength is 1.9285 nm for the current setup. Images are calibrated using a standard white and dark reference.

In this study we work on hyperspectral images from a VNIR camera with a very high spatial and spectral resolution. Images were obtained of the brain tissue and brain tumors, in-vivo (brain tissue open during surgery) and ex-vivo (extracted brain tumor on surgical table). Images were obtained over 20 operations. The images from both cameras were obtained during the surgical procedure. The information of the type of tumor and other associated patient information was missing during the current paper's study. This will be discussed more in detail. The images were taken during the surgical procedure (in-vivo) and right after the procedure outside the body (ex-vivo) with the extracted tumor placed under the camera. In this paper we restrict our study to the VNIR camera images due to their high spectral and spatial resolution.

Endmember Extraction

Hyperspectral unmixing is the problem of inferring which are the materials in the scene (visible directly or visible as a mixture) by the general process of dimensionality reduction. This basically consists in decomposing a hyperspectral image X = [x 1 , x 2 , ...x n ] with n pixels, where each pixel x i is p-dimensional spectrum. The two elements of the decomposition are a set of r reference spectra or endmembers and W = [w 1 , w 2 , ...w r ] and their relative proportions at The 2-rank approximation divisional step is used is to arrive with a unit rank matrices and the leaves of the hierarchical decomposition. The labels in the segmentations are consistently colored to match with the tree nodes. They help show the persistence of a subspace/clustering visually. This tree structure gives a visual idea of the structure of the tissue, labels that are close in terms subspaces branch. By viewing the endmembers (spectra) at different levels, we can infer visually that blood is a dominant material in the scene (blue peak around 400nm). The shallow clusters correspond to submatrices of low approximation error which close to rank-1 and thus correspond to a single tissue/material. each pixel(abundances

) H = [h 1 , h 2 , ...h n ],
where, w i is a p-dimensional spectrum and h i is r-dimensional coefficient weights. The spectrum at a pixel is a convex composition of endmembers. This imposes two conditions: i) the abundances need to be positive, ii) the sum of abundances by each endmember sum to unity at each pixel (MBDC + 14). The goal of the hyperspectral unmixing(HU) or endmember extraction process is to obtain the spectrum that corresponds to "pure" materials in the scene. This can be a complex problem based on the physical interactions between the materials, rays and the reflections. We assume the Linear mixing model(LMM) for its inherent simplicity and popularity in the remote sensing community. The model assumes that at a given pixel location the incident ray interacts only with one material in the true site of the pixel. Given the spectral and spatial resolution of the VNIR camera, this assumption seems to be a reasonable.

Separable NMF: Endmember extraction is well studied as a problem of non-negative matrix factorization problem since it respects the endmember constraints. It is written as a matrix factorization problem X ≈ W H, where this approximation is often evaluated by minimizing the Frobenius norm of the residual X -W H 2 F . Independent of how the decomposition is calculated, geometrically the NMF generates the vectors h i (endmembers) whose simplicial cone contains all the data points {x i } ∈ X. There can be many cones possible and the decomposition is nonunique. In case of high spatial and spectral resolution of HSI cubes, one can make the hypothesis that there exists pixels in the image which themselves form the bounding vertices of the a convex cone, and thus directly the endmembers. Physically this means that given a pixel spectrum we can infer the material it contains and be sure that its contribution is from a single material. This is called the separability condition (GL14) under which we obtain a unique and tractable decomposition of the NMF problem. Geometrically, the separability condition reformulates the NMF problem as that of finding the extreme rays of the conical hull of a finite set of vectors. Given a HSI cube X ∈ R n×p , we can write its separable NMF decomposition as

X = W H = W [I r , H ]Π
(1) where W ∈ R m×r , H > 0, Π is a permutation matrix, which rotates the r different endmembers. Columns of W are the endmembers, while the entry H(i, j) are the abundances/scalar coefficient of the ith endmember at the jth pixel in X. Geometrically this corresponds to the convex hull of the columns of X. Hierarchical Two-rank NMF: The hierarchical two rank NMF (H2NMF) first introduced in the paper (GKP15) aims at obtaining a rank-1 sub-matrix by hierarchical application of a NMF decomposition. The decomposition step consists of: i) projecting the data points on to the first two eigen vectors, and thus a 2-d cone, by calculating the singular value decomposition(SVD). ii) the direction matrix is approximated by a separable NMF decompositions, namely the SPA (GL14) of the data points under decomposition.

This yields pure-pixels in the HSI cube approximating the first principle component. The goal is to finally arrive with a set of leaf sub-matrices of rank one. This corresponds to pure-pixel or a single material. Though this hierarchical method is stated, the stopping condition for the growth of the hierarchy is not provided. Instead this is a parameter of the algorithm that determines the number of endmembers being extracted. At each step of this divisional hierarchical clustering step, the approximation error minimized is (GKP15):

E k = σ 2 1 (X(:, K 1 i )) + σ 2 1 (X(:, K 2 i )) -σ 2 1 (X(:, K i )) (2)
where σ 1 refers to the first singular value of the matrices, and K i corresponds to the indices of the parent and K 1 i , K 2 i are indices of the child clusters in the binary split. This approximation error corresponds to the rank approximation of input cluster or matrix. The split is decided by choosing line after the projection on to a 2-d cone. This is chosen by parameter δ * (see algorithm 1) so as to have a balanced set of points across the split (GKP15) as well a clusters with elements that are homogeneous. Cluster Endmember: The H2NMF algorithm associates an endmember with each cluster. The endmember of a cluster is the spectrum best approximating the first eigen value of submatrix, we follow (GKP15)'s measure. Mean removed spectral angle (MRSA) is defined for any two pairs of spectra x, y as:

φ(x, y) = 1 π arccos (x -x) T (y -ȳ) x -x 2 y -ȳ 2 ∈ [0, 1] (3)
where for a vector z ∈ R m , z = ∑ m i=1 z i . The MRSA provides the a measure to compare the pure-pixel in the input sub-matrix and its endmember spectrum produced by a rank-2 NMF decomposition. We use MRSA to evaluate the homogeneity of the hierarchical clustering at different levels.

We use the H2NMF algorithm (see Algo. 1) to provide a hierarchical subspace decomposition of the ex-vivo cubes to recover the pure-pixels by varying the level of decomposition. The decomposition is robust to noise and calculates quickly, given the huge sizes of the datasets (X's size is several 100000 pixels × 800 wavelengths). Secondly we visualize the hierarchical clustering structure which provides us a way to localize different tissue based on their paths in the hierarchy. In figure 2 we can observe from the clustering and the endmembers that, blood (endmember spectrum in blue, peaking around 400 nm for r = 4, 8, 12) at various levels of the hierarchy. Continuing further r > 8 we see that the tissues themselves are seen to have different layers and structuring of these layers.

Experiments and analysis

We partition the surgeries (numbered 1-20) into two disjoint sets: first the training set over which we perform unsupervised clustering and subsequent classifier training, second the test set over which we predict the cluster labels. We predict both on in-vivo and ex-vivo images from different operations in the test set to evaluate the validity of the different cluster labels. The training set consists of ex-vivo images, material cubes and certain windows from In-vivo images that contain synthetic surgical materials such as cotton, plastic cables and straps, as demonstrated in figure 1. The training set also contains supplementary tissue information such skull and surrounding tissues. Ideally this process captures the dictionary of materials to be clustered in the same subspace.

Training set : The training set is a union of pixels from images of the tumor ex-vivo and dictionary of materials such as plastic and cotton, from the training subset of operations. We denote the training images by X train ∈ R N×p , where N are the total number of pixels in the training set images together, and p = 826 is the dimensionality of the spectra at each pixel. The ex-vivo and material cubes in the test set will be clustered together so as to separate the different pure-pixels corresponding to different materials The training set does not contain any in-vivo images since images containing the target tumor samples shall remain unlabelled and not part of a class. Test set : The test set consists of both types of images: Images ex-vivo of the tumor, and images in-vivo from operations not in the training set. We refer to the different test cubes/images as X i and their corresponding classification map by the random forest as Y i . We use algorithm (1) described in section 2 to perform hierarchical clustering of the training-images. The principal goal here is to visualize the tissue structure as decomposed by the hierarchical clustering and aid a surgeon or expert in annotating these images. Another important goal is to determine the depth of the hierarchy at which we are able to separate the tumor, blood, tissues and materials and extract their endmembers or spectral signatures.

Cluster Evaluation

In this section we compare the performance of H2NMF, the clustering method used in our study, with respect to other clustering algorithms: hierarchical K-means(HKM), hierarchical spherical Kmeans(HSKM). For the cases of NMF when used as a clustering method, one does not explicitly use a distance function or dissimilarity (Kua14). The pre- Algorithm 1 H2NMF (GKP15)

1: procedure H2NMF(X ∈ R p×n + ) 2: K 1 ← {1, 2, ..., n} and K i ← / 0 for 2 ≤ i ≤ r 3: (K 1 1 ), K 2 1 ) = 2-rank-Split(X, K 1 ) 4:
while k < r do Iterate until r-clusters are reached 5:

j = arg max i=[1,2...r] E k from eq(2)
pick cluster j with largest error 6:

K j = K 1 j , K k = K 2 j
update cluster j into j and k 7:

(K 1 j , K 2 j ) = splitting(X, K j ) 8: (K 1 k , K 2 k ) = splitting(X, K k ) 9:
return K i cluster indices forming hierarchy 10: procedure SPLITTING(X ∈ R p×n

+ , K ⊂ {1, 2, ...n}) 11: [W, H] = rank2NMF(X(:, K ))
Projection onto 2-d cone 12:

x(i) = H(1,i) H(1,i)+H(2,i) 13:
Compute split parameter δ * parameter trades-off balanced split and cluster homogeneity 14:

K 1 = {K (i)|x(i) ≥ δ * } 15: K 2 = {K (i)|x(i) < δ * } 16: return K 1 , K 1
return 2-sets of indices corresponding to two sub-matrices of the input-matrix X. diction strength (TW05) is a well suited measure for cluster validation for NMF-based clustering because it does not rely on a distance function. Here we introduce a method similar to the one described in (Kua14) to evaluate the NMF. Given matrix X we split them into X train , X test over which we obtain NMF decompo-sitions, giving us

arg min X train -W train H train 2 F arg min X test -W test H test 2 F (4)
Now the prediction error can be obtained by solving the non-negative least-squares problem for

H predict = arg min X test -W train H 2 F (5)
We then compare the maximum abundance values in H test , H predict to check which endmember index contributed the most, and use this as a measure of prediction strength. This is because the H2NMF algorithm provides decomposition with one dominant endmember. This is a variant on the classical prediction strength measure (TW05). We demonstrate the results for the three hierarchical clustering algorithms. The prediction strength provides us a way to calculate the optimal number of clusters. Figure 4 demonstrates the prediction strength calculated over random samples of the training set X train . This measure could also provide a reasonable guide in determining which clusters to split, aside the approximation error. We plot the maximum value of MRSA from equation (3) for different number of clusters in in figure 4. This demonstrates how far the pure-pixel in cluster in a given level of hierarchy from its rank-2 NMF approximation. The figure also shows the plot of the average approximation error in eq. 2 over different leaf sub-matrices for the three hierarchical clustering algorithms. Both of these measures provide us a tool to asses if the spectrum and sub-matrix at a leaf are well approximated by a pure-pixel. Finally, we evaluate the hierarchical clusters intrinsically by using the gap statistic (TW05) which using the average pairwise dissimilarity per cluster to evaluate the number of clusters graphically, by looking for an elbow.

We now train the resulting clustering on the training set X train . We use the Matlab's implementation of the random forest to train the Random forest(RF) with the weak learner set to axis-aligned-hyperplane. We first cut the subspace hierarchy by choosing the number of clusters r, and then use the clusters as labels to train the RF-classifier.

The predictions of the RF-classifier are shown in figure 5 when trained on two different levels of clustering. We remark here that the actual subspace hierarchy obtained by a sequence rank-2 NMF decompositions, is not completely approximated by the RFclassifier. Nevertheless we do see a reasonable segmentation. The rings shown in the in-vivo images correspond to physical markers that were placed and localize healthy and tumor tissues for evaluation purposes. We do not obtain a single cluster label that might correspond to tumor. To improve the classification it would be ideal to integrate the structure of the subspace hierarchy. This can be done by learning a random forest for each level in the hierarchy. This would of course be costly. A more principled solution would be to minimize an appropriate loss function that approximates the H2NMF subspace hierarchy.

The resulting segmentation from the classifier provides firstly the stable pure signatures of different materials in the ex-vivo tumors. The results are verified basically by checking the presence of unique label (or sets of labels) in a ring (1)) in the in-vivo images. The tissue samples within these rings have been verified posteriori by pathology. One ring is assured to surround a cancerous tissue sample, while the other a healthy one. In this setup there was estimate of the depth cancerous and its uniform exposure to the VNIR-camera. This is tough to ensure during critical surgeries. The color-map of the classification result in figure 5 correspond to the color-map of the unsupervised clustering in figure 3. The semantic meaning of the labels is temporary and the color-map further requires verification by a surgeon or pathology expert, to which finally a more sensible name shall be attached. But it is evident that we are able to identify plastic surgical materials, cotton.

Issues

Spectral proximity between tissues : (SMP + 13) study differentiation between benign and malignant tumors in canine mammary tumors, and detect tu-mors using a imaging spectrometer. Ovarian tissue characterization and tumor differentiation (UBS + 01), were performed using NIR spectroscopy. The results showed that both slope and intensity of reflectance spectra may have the ability to discriminate normal from abnormal conditions of the ovary. Slopes between 510 nm and 530 nm were most discriminatory for ovaries, while 630 nm to 900 nm in the study on kidneys (Pes07). The differences between malignant and benign tissues in human breast tissues, have been attributed to the metabolic difference due to the presence of more oxy-hemoglobin, lipids and water. These elements were used as the spectral markers considered for the additive construction of the spectrum. (KCT + 10). This proximity in the spectral response of normal, benign and malignant tumor plays an crucial role in classification of tumors pixels. In our study, due the lack of information on the actual spectrum of tumor tissues we are unable to provide a good division of the subspace and the subsequent refinement of decision boundary between normal and possible cancerous tissue.

Specular reflections : Specular reflections are a prime problem in all our clustering and classification steps and contributes to mis-classification. Specular reflection occurs when the light source can be seen as a direct reflection on the surface of the tissue under study. The pixels affected by specular reflection contain mainly the light source and may or may not contain sufficient information on the reflectance of the tissue under study. These reflected rays may also be suspect of undergoing multiple reflections. An example is shown in figure 6.

Conclusion

In this paper we applied and studied the hierarchical rank-2 NMF clustering and endmember extraction to segment brain tissue structure in an unsupervised setting. Subsequently we have tried to train a random forest to learn the structure of the subspace clusters.

Subspace learning and Spatial features : Though the classifier can't discriminate between normal tissues and tumor spectra, it provides a good segmentation of the tissues and surgical materials in the scene. The H2NMF hierarchy provides a hierarchy of low rank sub-matrices that contain pure-pixels, and in this case correspond to different brain tissues. This cluster hierarchy can further be used by a surgeon in loop to search more pertinent subspaces and feed the results back to the classifier. We have seen the effect Figure 4: Left: the variation of the maximum MRSA(mean-removed spectral angle) to evaluate the homogeneity of clusters in terms of endmember signature. This is evaluated across the three hierarchical clustering algorithms: H2NMF, HSKM, HKM. We observe that all the three algorithms provide satisfactory results. Center: figure also shows the variation approximation error in equation ( 2) with increasing leaves for three hierarchical clustering algorithms. Right: Plot demonstrating the mean prediction strength calculated using the NMF decomposition. The red bars show the variance of the prediction strength. Optimal clusters are chosen based on the maximum of the prediction strength. Here for example it is 14.

Figure 5: Results of the classification on the test set which consists of Ex-vivo images and In-vivo images test-set. Two sets of classification maps are shown, the top image tile corresponding to learning the labels by cutting the subspace tree at r = 18 clusters and the bottom tile at r = 24. We observe that there are different tissue structures and materials found in the original clustering in figure 3. It would be ideal to have a cluster/endmember corresponding to the tumor, though we remark that it is not obvious to extract them or check for their existence. We can also notice that there are sets of labels co-occurring in the marker rings though no trivial or consistent co-occurrence that can ensure a tumor detection. We can also notice in the figure that certain materials have been approximated by the same subspace cluster: the cotton and the plastic. This can also be visually confirmed by look at the subspace tree. Discriminating materials in this subspace tree structure is a crucial analysis tool. As already mentioned, multiple classifiers are needed to approximate well the subspace clustering hierarchy. This would merge the two levels of classification in this image into a more precise classification map. of learning classifiers on different levels of the hierarchical clustering. Though this produces reasonable segmentations, it does not approximate the subspace hierarchy exactly. We will study how to define loss functions to encode this structure. Given that the difference in spectrum between tumor tissues and normal tissues are small, and prone to noise and variations across patients, building a robust spatial structure descriptor is important. One of the results of our study points to the fact that the tissue structures surrounding a tumor is a key feature. There is already evidence that in the micro-scale these normal and cancerous tissues have a different toplogical arrangement (Dvo03). This structural information of tumor is useful to obtain a better detector for cancerous tisues in hyperspectral images, since spectrum alone is not sufficient to classify them robustly. Spatio-spectral features have been used in (LHW + 14), though our aim in the future is to use the spatial features extracted by scattering-transform (BM13) on the different abundance maps extracted by the hierarchical clustering algorithm. This enables us to perform a principled search for features across various scales.

Figure 2 :

 2 Figure2: Segmentation Hierarchy showing a tumor ex-vivo decomposed at different levels r = 1 to 32 of the hierarchical rank-2 decomposition (H2NMF) algorithm. The 2-rank approximation divisional step is used is to arrive with a unit rank matrices and the leaves of the hierarchical decomposition. The labels in the segmentations are consistently colored to match with the tree nodes. They help show the persistence of a subspace/clustering visually. This tree structure gives a visual idea of the structure of the tissue, labels that are close in terms subspaces branch. By viewing the endmembers (spectra) at different levels, we can infer visually that blood is a dominant material in the scene (blue peak around 400nm). The shallow clusters correspond to submatrices of low approximation error which close to rank-1 and thus correspond to a single tissue/material.

Figure 3 :

 3 Figure 3: H2NMF clustering on the training set X train which consists of input ex-vivo tumor tissue HSI cube and on materials. We have in this image the ring on a cotton, and a cable in scene with open cranium and healthy tissues. This labeling with the corresponding HSI cubes serves as the training set for constructing the random forest classifier.

Figure 6 :

 6 Figure6: An in-vivo image and its classification map with the specular reflection marked in dotted red ellipse. This phenomena is not restricted to the optical range.

Acknowledgment : This work has been supported in part by the European Commission through the FP7 FET Open programme ICT-2011.9.2, by the European Project HELICoiD "HypErspectral Imaging Cancer Detection" under Grant Agreement 618080.


A L Zimpeck

C Meinhardt

L Artola

G Hubert

F L Kastensmidt

R A L Reis

Impact of Different Transistor Arrangements on Gate Variability

This paper evaluates a set of complex cells with different transistor arrangements that implement the same logic function. These cells were evaluated under nominal conditions and with gate variability at layout level. The purpose is to verify what topology is more appropriate to increase the robustness of cells regarding the process variability issues. Results emphasize the importance of investigating the effects caused by process variability in FinFET technologies, as the electrical characteristics of circuits suffer significant changes. In general, the best choice is to use the network that the transistor in series is as far as possible to the output node. However, a trade-off needs to be done due to performance and power consumption penalties.

Process Variability on Nanotechnologies

For continuing with the technology scaling on sub-22nm, novel materials and new devices architectures had to be adopted in integrated circuit designs [START_REF] Pradhan | Investigation on Asymmetric dual-k spacer (ads) Trigate wavy FinFET: A novel Device[END_REF]. Multigate devices were introduced with the objective of overcome obstacles encountered in standard Complementary Metal Oxide Semiconductor (CMOS) technologies and keep scaling [START_REF] King | FinFETs for Nanoscale CMOS Digital Integrated Circuits[END_REF]. FinFETs gained prominence for presenting an excellent short-channel effects (SCE) control, reduced leakage currents, high driving capability, better yield, and the fabrication process compatible to the conventional CMOS [START_REF] Mitra | Robust System Design to Overcome CMOS Reliability Challenges[END_REF].

The emerged technologies raised relevant topics related to reliability and robustness that to need to be investigated. The process variability still as one of the most significant challenges in nanometer regime becomes even more meaningful in designs based on FinFET technology. The sub-wavelength lithography introduced difficulties in transferring the small geometric patterns required by the layout of nanometer technologies to the substrate surface [START_REF] Dadgour | Grain-Orientation Induced Work Function Variation in Nanoscale Metal-Gate Transistors -Part I -Modeling, Analysis and Experimental Validation[END_REF]. The side effect is that this manufacturing step impacts directly on the transistor threshold voltage (V TH ).

FinFET technologies use high-k/metal gate stack to improve the gate control on the channel region [START_REF] Ortiz | High-k Organic, Inorganic, and Hybrid Dielectrics for Low-Voltage Organic Field-Effect Transistors[END_REF]. Thereby, the Metal Gate Granularity (MGG) has been identified as a source of statistical variability. The orientation of the grains used in the gate is not considered ideal as a unique metal uniformly aligned. The metal gate having different work-functions (ϕm) randomly aligned that implies in higher work-functions fluctuations (WFF) [START_REF] Dadgour | Grain-Orientation Induced Work Function Variation in Nanoscale Metal-Gate Transistors -Part I -Modeling, Analysis and Experimental Validation[END_REF][START_REF] Dadgour | Statistical Modeling of Metal-Gate Work-Function Variability in Emerging Device Technologies and Implications for Circuit Design[END_REF]. Moreover, the imperfections caused by process variations can influence the fin height, fin width, fin-to-fin similarity, and the resistance MOL (Middle of Line).

All these factors can compromise entire blocks of cells because they modify the transistor structure. Therefore, the electrical properties of the circuits also suffer deviations [START_REF] Pedram | Thermal Modeling, Analysis, and Management in VLSI Circuits: Principles and Methods[END_REF]. The impact on performance and power metrics can accelerate the circuit degradation besides introducing errors in the circuit functionality. The consequence is the parameter yield loss leading to a higher cost of fabrication. This emphasizes the importance of creating new design guidelines able to deal with the challenges imposed by the process variability in nanometer technologies.

Many related works studied the challenges and the impact of variability in the circuits and devices at nanometer technologies [START_REF] Dadgour | Statistical Modeling of Metal-Gate Work-Function Variability in Emerging Device Technologies and Implications for Circuit Design[END_REF][START_REF] Pedram | Thermal Modeling, Analysis, and Management in VLSI Circuits: Principles and Methods[END_REF][START_REF] Wang | Statistical Variability and Reliability in Nanoscale FinFETs[END_REF][START_REF] Mutlu | Statistical Methods for the Estimation of Process Variation Effects on Circuit Operation[END_REF][START_REF] Cao | Mapping Statistical Process Variations Toward Circuit Performance Variability: An Analytical Modeling Approach[END_REF][START_REF] Wu | Impact of Surface Orientation on Vth of FinFET[END_REF]. The influence of the variation of the main geometric parameters on the currents I ON and I OFF of FinFET transistors considering a set of predictive FinFET technologies sub-22nm was analyzed in [START_REF] Meinhardt | Predictive Evaluation of Electrical Characteristics of sub-22nm FinFET Technologies under Device Geometry Variations[END_REF]. In [START_REF] Zimpeck | Impact of PVT Variability on 20nm FinFET Standard Cells[END_REF], the impact of Process, Voltage and Temperature (PVT) variations in timing and power on a subset of circuits using a commercial standard cell was studied. Different transistor sizing techniques were applied in FinFET circuits with PVT variations in [START_REF] Zimpeck | FinFET Cells with Different Transistor Sizing Techniques against PVT Variations[END_REF]. However, there are only few approaches in the literature to mitigate the variability effects. The replacement of traditional inverters by Schmitt Triggers (ST) into near-threshold full-adder architecture was proposed by [START_REF] Dokania | Circuit-Level Design Technique to Mitigate Impact of Process, Voltage and Temperature Variations in Complementary Metal-Oxide Semiconductor Full Adder Cells[END_REF] as a process variability mitigation technique. Their results showed a significant reduction in timing and power deviation at 16-nm bulk CMOS technology. The similar analysis was done in [START_REF] Toledo | Impact of Schmitt Trigger Inverters on Process Variability Robustness of 1-Bit Full Adders[END_REF].

Transistor arrangement optimization is a technique used to design faster circuits [START_REF] Butzen | Transistor Network Restructuring against NBTI Degradation[END_REF], but also to deal with Bias Temperature Instability (BTI) effects [START_REF] Da Silva | CMOS Logic Gate Performance Variability related to Transistor Network Arrangements[END_REF] or to improve design robustness against permanent and transient faults. The adoption of complex gates reduces the transistor number that is correlated to the area and also, reduces the delay and power consumption. Some logic cells can be designed using different kinds of transistor networks [START_REF] Da Rosa | Switch Level Optimization of Digital CMOS Gate Networks[END_REF]. It is well known that different transistor combinations, that implement the same logic function, present different electrical and physical characteristics as well as distinct behavior under process variations [START_REF] Da Silva | CMOS Logic Gate Performance Variability related to Transistor Network Arrangements[END_REF].

In this context, this work evaluates the electrical characteristics of a set of complex gates with different transistor arrangements considering the ideal behavior and a realistic behavior, i.e., with WFF variations. The more suitable topology to process variability mitigation is analyzed. Moreover, gains and penalties regarding performance and power consumption were investigated alongside variability robustness. The set of information provided here helps the designers to choose the most appropriate physical layout depending on the target application.

Section II shows the methodology adopted in this work. A review about the 7nm FinFET technology from ASAP and the layout techniques are presented in Section III. Section IV and V propose a description of simulation and a discussion of the results.

Methodology

This work evaluates the impact of process variability on different transistor arrangements. Fig. 1 shows the design flow adopted in this research. First, a set of six complex gates of 3 to 5-inputs was chosen to perform this analysis. This choice was made due to the cells present serial transistors allowing more than one topology type. All complex cells are designed with symmetric sizing such that NFET and PFET transistors have the number of fins equal to three. The purpose of the next stage was dedicated to the development of the schematic and the physical layout for all investigated complex gates. The layouts were submitted to the physical verification flow composed by DRC (Design Rule Check) and LVS (Layout Versus Schematic) steps. Both steps are based on the technology rules of the 7-nm FinFET Predictive Process Design Kit (ASAP7) developed at Arizona State University in partnership with ARM Ltd [START_REF] Clark | ASAP7: A 7-nm FinFET Predictive Process Design Kit[END_REF]. The main parameters of this technology are summarized in the Table I. The behavior of each complex gate is verified to certificate that it implements your function correctly. Nominal conditions were used as a reference to the variability evaluation. After setting the parasitic extraction configurations with the options desired, the parasitic wire capacitances and resistances from the layout were extracted. A new netlist was generated where each net contains one subckt with the RC tree structure modeling the net and the connections between the parasitic networks.

Process variability was taken through two thousand Monte Carlo simulations carried out in HSPICE with the WFF parameter modelled as a Gaussian function with a 3σ deviation of 5% from nominal values [START_REF] Meinhardt | Predictive Evaluation of Electrical Characteristics of sub-22nm FinFET Technologies under Device Geometry Variations[END_REF] [START_REF] Alioto | Variations in Nanometer CMOS Flip-Flops: Part II-Energy Variability and Impact of Other Sources of Variations[END_REF]. The mean (μ), the standard deviation (σ) and the normalized standard deviation (σ/μ) were analyzed for all complex gates.

The choice of the transistor topology with less impact of variability effects is given observing the deviation, i.e., the σ/μ relation. An arrangement is pointed as the best choice to mitigate the variability impact if it has the lowest value for the σ/μ relation. All complex gates drive a FO4 inverter and this has two inverters used as the load connected to each input.

Layout Description based on ASAP7

All the complex gates were designed using the ASAP7 Process Design Kit (PDK) that considers a not yet available technology node for academic use. This PDK was chosen because admits realistic design conjectures regarding the lithography step and the current technology competencies [START_REF] Clark | ASAP7: A 7-nm FinFET Predictive Process Design Kit[END_REF]. Table II shows some essential layers and the design rules considered in this PDK.

Table II -Some key layers, widths and pitches [START_REF] Clark | ASAP7: A 7-nm FinFET Predictive Process Design Kit[END_REF] FinFET technologies have a width quantization characteristic [START_REF] Alioto | Variations in Nanometer CMOS Flip-Flops: Part II-Energy Variability and Impact of Other Sources of Variations[END_REF]. With a 27 nm fin pitch, highdensity layout design is achieved with three fins for each PFET and NFET devices [START_REF] Chava | Standard Cell Design in N7: EUV vs. Immersion[END_REF]. Fig. 3 shows the layout view of the AOI21 complex gate with the serial transistor close to the output. The difference of the far topology is the serial transistor connected to supply voltage. The cell height is set to 7.5 tracks of metal 2 (M2) that correspond to 0.27 μm. The supply/ground rails have a tall approximately equal to 1.22 tracks of M2. The total area of complex cells with three, four and five inputs is 84.78 nm 2 , 101.74 nm 2 and 118.67 nm 2 , respectively.

This PDK has the manufacturing process composed by front end of line (FEOL), middle of line (MOL) and back end of line (BEOL) as shown the Fig. 3. The layouts were developed with a continuous diffusion layer and every gate must have at least one other gate surrounding the horizontal axis. The fin layer polygons should have an equal length along the X-axis. The Source-Drain Trench (SDT) connects the active area to the LISD layer. The Local-Interconnect Gate (LIG) is used to connect the gate terminal and the Local-Interconnect Source-Drain (LISD) connects the source and drain of the transistors. The function of V0 is to join the LIG and LISD to the BEOL layers. The metal 1 (M1) is used for intra-cell routing and short connections. The power rails are made using M1 and LIG that are connecting at V0 in regular intervals. For obtain success in the LVS verification step is necessary to add a TAP cell in the design of each cell.

Investigation of Arrangement Effects

In this work, the analysis is divided into three main parts. First, a comparison of electrical characteristics of each complex gate with different networks considering the typical behavior and under WFF variations is made. The second part verifies the more suitable topology to gate variability mitigation as well as the improvements and penalties on delays and power consumption.

Fig. 4 and Fig. 5 show the difference between nominal values and the mean of the results obtained by Monte Carlo simulations for power consumption and considering the worst-case of the propagation delay, respectively. The worst-case delay is the higher propagation time from all timing arches of each function. The red bars indicate the increase caused in the electrical characteristics due to process variability. The power consumption suffers around 4-5.5% of deviation while the effects on the worst case of the propagation delay were of up to 7.5% from the nominal values. series far to the output node result in more power consumption and longer delays in the majority of cases. This outcome is very similar to what happens in the analysis of timing and aging as related in the literature [START_REF] Da Silva | CMOS Logic Gate Performance Variability related to Transistor Network Arrangements[END_REF]. III summarizes the mean (μ) and the normalized standard deviation (σ/μ) of the worst-case delay and the power consumption for all complex gates analyzed. Also, to show the process variability impact dimension on energy and the energy consumed per switching event, results are evaluated considering the Power-Delay-Product (PDP), where the PDP is calculated for each Monte Carlo simulation. The metrics of intermediate topology are between the close and far values, and then, they were omitted. III enforce the high sensibility to process variability of these cells on 7-nm technology, inserting wide deviations on power and delay, and consequently on PDP. When the power deviation is analysed, all complex gates found as the best cases the options where transistors in series are placed as far as possible from the output node ((b) alternatives presented in Fig. 1). Moreover, the OAIs complex gates are around 7-9.5% more sensible to WFF variations than the AOIs. On the other hand, the propagation time is most impacted by more than 30% of deviation independently of the function or transistor topology. The OAI21, AOI211 and AOI221 gates have attenuation in the effects of variability when the transistors in series are closer to the output node.

All results in Table

If considering the only the on delay, the choice of the far topology is disadvantageous for the OAI21, AOI211 and AOI221. This choice increases in the delay mean value of about 53.2% compared to the expected nominal behavior. In this casethe far topology brings few advantages with a better trade-off. Adopting the close topology in these cells, it is possible to attenuate the effects of process variability until 2.3% besides improving the performance of OAI21 in 14.5%.

Otherwise, for the OAI211 cell, the far topology is promising, presenting a variability robustness gain of 4.9% on average delay compared to close topology with a decrease of 14.7% in delay mean value with this topology. The far topology also becomes the OAI221 6.2% more robust to variability effects. However, in comparison to close topology, this cell presents an increase of 8.7% in the delay mean values. Although with some penalties, the far topology is a good choice to improve the variability tolerance of this cell.

However, the adoption of the far topology increases in the power robustness variability for all complex gates and this effect is reflected on the PDP, as shown Fig. 6. For the majority of investigated complex gates, the far topology remains the best choice. All the AOI complex gates achieve a decrease in the process variability effects using the far topology. For the AOI221, the improvement in relation to the close topology reaches to 6.9%.

Fig. 6 PDP Normalized Standard Deviation considering WFF variations

This benefit of far topology is due to the high improvement on power robustness compared to the close topology. Fig. 7 illustrates the power evaluation and highlights the main benefits and penalties in term of power of using the far topology instead of the close topologye. The complex gates with three and four inputs present improvements that can reach until 7.6% and 8.6%, respectively. However, a penalty in the power consumption mean is observed. The AOI211 and the OAI211 gates consume around 10% and 8% more power compared to nominal values, respectively. In general, the AOIs present a little more gain in power robustness variability, but they have a higher increase in the power consumption compared to the OAIs complex gates. Fig. 7 Gains and penalties related to power in adopting the far topology

Conclusions

This work presented an evaluation of the influence of process variability on a set of complex cells designed in 7nm FinFET ASAP technology. The objective was verified what is the best arrangement to mitigate the effects of process variability. The improvements and penalties also are analysed. For ensure the improvement in the power robustness variability, the best alternative is adopted the far topology in the designs. On the other hand, the choice of the far topology to increase the delay robustness variability is only advantageous for the OAI211 and OAI221 complex gates. For the others complex cells analysed, the close arrangement is the better alternative to mitigate the process variability. The improvements in variability caused by the power are reflected on the PDP. Then, far topology seems to be the best choice to mitigate the process variability effects.

Fig. 1

 1 Fig. 1 Design flow of the experiments

Fig. 2

 2 Fig. 2 Complex gates logically equivalent with series transistor (a) close, (b) far and an (c) intermediate place to the output node

Fig. 3

 3 Fig. 3 Layout of AOI21 Complex Gate in 7nm FinFET Technology (ASAP7) with serial transistor close to the output

Fig. 4 Fig. 5

 45 Fig. 4 Difference between power consumption nominal values and under WFF variations

Table I -

 I 7nm FinFET device parameters[START_REF] Clark | ASAP7: A 7-nm FinFET Predictive Process Design Kit[END_REF]

	Parameter		7nm
	Nominal Supply Voltage 0.7 V
	Gate Length (L G )	21nm
	Fin Width (W FIN )	6.5nm
	Fin Height (H FIN )	32nm
	Oxide Thickness (Tox)	2.1nm
	Channel Doping	1x10 22 m -3
	Source/Drain Doping	2x10 26 m -3
	Work	NFET	4.3720
	Function	PFET	4.8108

  The arrangements with transistors in

	Layer	Lithography	Width/ drawn (nm)	Pitch (nm)
	Fin	SAQP	6.5/7	27
	Active	EUV	54/16	108
	Gate	SADP	21/20	54
	SDT/LISD	EUV	25/24	54
	LIG	EUV	16/16	54
	VIA0-VIA3	EUV	18/18	25
	M1-M3	EUV	18/18	36
	VIA4-VIA5	LELE	24/24	34
	M4-M5	SADP	24/24

Table III -

 III Mean and standard deviation of the worst-case delay, power consumption and PDP

	Metrics	AOI21 Close Far	OAI21 Close Far	AOI211 Close Far	OAI211 Close Far	AOI221 Close Far	OAI221 Close Far
	Delay	μ (ps) σ/μ (%)	4.2 34.3	6.4 33.8	4.4 32.7	6.4 33.4	9.3 33.6	10.0 34.3	9.4 31.6	8.0 35.0	11.6 35.0	12.0 35.6	9.9 34.4	10.7 32.2
	Power	μ (nW) 274.2 297.8 255.6 270.8 278.5 306.4 302.8 308.0 308.0 328.5 393.4 307.5 σ/μ (%) 24.0 22.1 26.0 24.2 27.8 25.4 29.7 27.3 28.9 27.3 31.2 29.8
	PDP	μ (aJ) σ/μ (%)	2.3 27.0	2.6 26.0	2.1 26.7	2.2 27.8	3.4 29.7	4.0 28.1	3.0 30.5	3.2 31.4	4.7 30.3	5.2 28.2	3.9 31.8	4.1 31.4

Acknowledgements

This research is sponsored by the Brazilian National Council for Scientific and Technological Development (CNPq).