Topological Optimization of Dynamic Characteristics for Orthotropic Material Structure Using Shape Derivative and Augmented Lagrangian Method

This paper presents a new level set method for topology optimization of dynamic structure of orthotropic materials using the shape derivative analysis and augmented Lagrangian method. The design boundary of the structure is embedded implicitly into the zero level set of a higher dimensional scalar function, which is mathematically represented as a Hamilton-Jacobi partial differential equation (PDE). The design sensitivity of the dynamic structure is obtained by the combination of the level set representation and the shape derivative method. In doing so, the evolution of the design boundary is advanced iteratively in terms of the solutions of the Hamilton-Jacobi PDE using explicit time marching schemes. Some typical numerical examples are applied to demonstrate the validity of the present method.


Introduction
Structural optimization, in particular the shape and topology optimization, has been identified as one of the most challenging tasks in structural design.By combining the finite element method with various optimization algorithms, investigation has made a considerable progress in the structural design engineering during the last few decades.The optimization method has been widely received by automotive, aerospace, nuclear and other high-tech industries, which is used to satisfy the design requirements for maximum performance, minimum weight, cost efficiency, and environmental consideration.Various techniques and approaches including ground structure method (Prager 1974, Rozvany 1979, Hagishita and Ohsaki 2009), homogenization method (Bendsøe and Kikuchi 1988), solid isotropic material with penalization (Zhou and Rozvany 1991, Mlejnek 1992, Bendsøe and Sigmund 1999) and level-set method (Wang etc. 2003, Allaire etc. 2004) have been developed.Other methods, such as the current flexible building block method (Kim etc. 2006) , the unit cell approach (Wang 2004), the radial-basis-function (RBF) level-set method (Wang 2006) can be considered nearly as evolutions of the approaches mentioned above, thus they will not be discussed here.The essential idea of all these methods is to transform the design problem into an optimal material distribution problem so that the configuration designed can be measured quantitatively by an objective function.Also, there are some differences among these approaches in their representations and modeling schemes.
It is not easy to design a mechanical structure with the desired natural frequency and mode shape.As for parametric structural optimization problems, many methods were presented for sensitivity analysis of the dynamic structures, such as these reported by Fox and Kapoor 1968, Nelson 1976, Seyranian etc.1994and Luo etc.2006 etc.For non-parametric structural optimization problems, a technique for topology optimization of linear elastic continua based on homogenization theory was investigated by Bendsøe & Kikuchi 1988, and it has been applied to topology optimization problems of the dynamic structures by Diaz andKikuchi 1992 andTenek andHagiwara 1994 etc.By employing this technique, the material distribution was formulated with parameters of periodic microstructure.Another few papers explored the nonparametric boundary optimization problems of linear elastic continuum.Based on the gradient method in a Hilbert space, Inzarulfaisham &Azegami 2004 andMeske etc. 2006 proposed a numerical method called the traction method and they applied it to maximum lower-order natural frequencies under a volume constraint.Later, Seonho etc. 2006 developed a topological shape optimization method which was applied to power flow problems at high frequencies and the necessary design gradients were computed by using an efficient adjoint sensitivity analysis method.Allaire and Jouve 2005 applied the level-set method to the shape and topology optimization, and they maximized the first eigenfrequency under the specified weight constraint in the three-dimensional space.
This paper explores the feasibility employing level set method to optimize an orthotropic material structure.The augmented Lagrangian multiplier with inequality constraint is extended into the topology optimization of the design structure.By the combination of the shape derivative, level set and augmented Lagrangian multiplier method, the sensitivity of the maximum lower order natural frequency is investigated theoretically.The optimization boundary of orthotropic material structure is embedded implicitly into the zero level set of a higher dimensional scalar function, which is mathematically represented as a Hamilton-Jacobi type partial differential equation (PDE).The movements of the design boundary are advanced iteratively in terms of the solutions of the Hamilton-Jacobi PDE using explicit time marching schemes.The 2D numerical examples are employed to demonstrate the validity of the current method.

Topological optimization theories
In order to obtain the dynamic topological structure of an orthotropic material, the key problem is to develop a formula of the maximum lower order natural frequency.The constructive relationship is the foundation of stress equation and strain equation.

Constitutive relationship
As we know, the constructive relationship of the orthotropic material under the plane stress state is given as  are the elements of stiffness matrix, stress matrix and strain matrix in the local coordinate system (or the material coordinate system).From the material coordinate system to the global coordinate system， the transformation equation of the elastic parameters is and  is the angle measured from the global coordinate system to the material coordinate system.If the material coordinate system is inconsistent with the global coordinate system, the stress-strain relationship in the global coordinate system can be given as follow (3)  and Q' ij are the elements of stress matrix, strain matrix and stiffness matrix in the global coordinate system.
Equation ( 3) is the constructive relationship of 2D orthotropic material.After obtaining the formula of the constructive relationship, the sensitivity of the maximum lower order natural frequency will be derived.

Augmented Lagrangian method
It is very difficult to obtain a reasonable optimal result with a fixed Lagrangian multiplier in the structural topology optimization.In this paper, the augmented Lagrangian method with slack variables is explored to solve the optimization problem with inequality constraints.
Generally speaking, the augmented Lagrangian method can be classified into equality constraint and inequality constraint problems.In the case of equality constraints, the quadratic penalty function is adopted to penalize constraint violations by squaring the infeasibility (Chen and Shapiro 2006).For the inequality constraint problem, it can be converted into an equality constraint by introducing the slack variables.Let us describe the issue firstly , and the augmented Lagrangian equation with equality constraints at the kth iteration can be obtained (Nocedal and Wright 1999).
The following equation can be Equation ( 6) is a convex quadratic function with respect to slack variables j r , so the optimal value of the slack variables is zero.By considering equations ( 6) and (8), j r value at kth iteration is expressed as follow Equation ( 9) solves the problem of slack variables j r .By substituting j r into equation ( 6), the next issue is an optimization problem of equality constraint as follow 10) is an augmented Lagrangian equation with equality constraints.We can apply the iteration formula of equality constraints to here, and then we have In this way, the inequality constraint problem can be optimized by the augmented Lagrangian method.For more details, readers can refer to Bertsekas 's book published in 1982.In the next section, this method will be applied to the topological optimization of maximum lower-order nature frequency.

Maximum lower-order natural frequency
The topological optimization structure of maximum lower-order natural frequency with inequality constraints is defined as follow . e(u) is strain tensor at any point ] is the stiffness matrix as shown in equation ( 3).Utilizing the standard and augmented Lagrangian multiplier methods for equation ( 13), we have Here v stands for both the adjoint variable and Lagrangian multiplier.By solving shape derivative of equation ( 14) with respect to t (Sokolowski etc. 1991), the following formula, according to equation ( 11), can be obtained Based on equation ( 12), the formula to update the Lagrange multiplier is As we know, the problem of maximum lower-order natural frequency is a self-adjoint problem (Meske etc. 2006) , so v u  .Let us define sensitivity of the objective functional as follows Λρuv e(u) Qe(u) The shape derivative of Lagrangian functional can be rewritten as The equation ( 18) will be changed into Equation ( 19) is an evolution speed of the frequency optimization problem.This velocity can make equation ( 20) not more than zero, so the objective function always decreases until 0 , the shape designed is the optimized result.In next section, level-set method will be investigated in detail.

Level-set method
After obtained the evolution velocity of the design boundary, the PDE of Hamilton-Jacobi equation (Osher and Sethian 1988 ) will be introduced in this section.
The level-set method was developed by Osher andFedkiw 2003 as well as Sethian 1999 to solve the problems such as tracking, computer vision, crack propagation in solid material, image processing, boundary evolution simulating in fluid mechanics etc.The topology optimization can be expressed as a dynamic evolution process of the level-set function varying with pseudo time t.The embedded function allows its surface to move up and down on a fixed coordinate system without changing its surface topology structure, and the optimization shape embedded in level-set function can automatically modify the topology structure by boundary merging and breaking.We can track the topological change of the design structure by checking the level-set surface.This process can be completed by solving PDE with the evolution velocity n V and a set of initial value.The Hamilton-Jacobi PDE is generally solved by using the upwind based on the fixed Eulerian grids.Therefore, in the following paragraph we can employ the level-set method to capture the optimization structure on a fixed mesh.
Let us define a closed subset In the level-set method, the design boundary is embedded implicitly into the zero level-set of a higher-dimensional level-set surface.The movement of the higher-dimensional embedded function is governed by the Hamilton-Jacobi PDE, the evolution velocity and a set of initial value.The Hamilton-Jacobi equation can be given by differentiating equation 0 ) ), ( ( Here  is the gradient operator.In the three-dimensional space, let us take Equation ( 23) demonstrates that only can the normal velocity n V provide the contribution for the boundary evolution.Furthermore, we simplify the second term of equation ( 23) The Hamilton-Jacobi equation is rewritten as

Now, we have velocity n
V to make the objective function decrease.The optimization boundary embedded into the level-set function can automatically be updated in the iteration process.In next section, the optimal numerical algorithm will be described in detail.

Optimum numerical algorithms
Once the velocity n V and Hamilton-Jacobi PDE (21) have been obtained, the finite difference method with respect to t can be obtained 0 Here n V is the velocity at time n t , and  denotes the gradient operator.For the 3D case, the formula can be written as . t  depends on the Courant-Friedreichs-Lewy (CFL) condition.For the 2D case, the CFL condition can be obtained Here  and  are the absolute values of the evolution velocities in x and y directions.x  and y  are lengths of sides over the entire Cartesian grid. 1 0   α , the common nearly-optimal value is 9 .0   .Equation (28) indicates that the numerical waves must be at least as fast as the physical waves.In order to make the level-set function regular, we reinitialize it periodically by solving the following equation The accurate and robust numerical method was developed by Osher and Fedkiw 2003.To simplify this paper, readers can refer to Osher and Fedkiw's book for more details.

Numerical examples
In this section, several topological optimization examples of the orthotropic material are explored by employing the method presented in this paper.The objective function of the optimization problem is the maximum first-order natural frequency.The numerical examples are the cantilever beams of 2D structure whose geometric parameters and boundary conditions are shown in Fig. 2. The material parameters are 1 E =135  10 9 , 2 E =28.8  10 9 , 12 G =4.47  10 9 and 12  =0.33, and mass fraction is no more than 0.4(M 0 =0.4).The initial design structure and initial value of level-set function are shown in Fig. 3 and Fig. 4, respectively.The length and width of the design domain are 0.1m and 0.05 m, and the thickness is 1 mm.When material directions are 0 0 , 30 0 , 45 0 , 60 0 and 90 0 , figures 5, 10, 12, 14, and 16 are some topological optimal results and figures 9, 11, 13, 15 and 17 are the final level-set surface.When material direction is 0 0 , Figures 6, 7, 8 are volume, augmented Lagrange multiplier and the first-order natural frequency variable with iteration number, respectively.If the material directions are 0 0 and 90 0 , the optimization structures are symmetric.When the material directions are 30 0 , 45 0 and 60 0 , the optimal results are not symmetric and the result of 30 0 material direction is the most non-symmetrical in the three cases.This result illustrates that the material performance as well as material direction can affect the topological optimal result.
The optimization structures shown in Figures 5, 10, 12, 14, and 16 are basically coinciding with the paper of Allaire and Jouve 2005, which proves our method to be very valid.

Conclusions
The feasibility of employing level set method to optimize an orthotropic material structure is explored theoretically.The augmented Lagrangian multiplier with inequality constraints is extended into the topology optimization of the design structure.By the combination of level-set, shape derivative and augmented Lagrangian method, the sensitivity formulae of the topological optimization of dynamic structure is derived.The optimization boundary of the orthotropic material structure is embedded into the zero level-set of a higher-dimensional scalar function.The movements of the boundary are iteratively promoted by the solutions of PDE.The evolution velocity of moving boundary can be controlled by CFL condition and sensitivity.Examples studied are employed to illustrate the validity of the present optimization method.The conclusions indicate this investigation will provide an important foundation for the advanced development topology optimization of the complex constitutive relationship material.
inequality constraints.By defining slack variables j r , equation (4) can be changed into the following form

A
linear elastic continuum is in an open domain n R Ω  , n=2, 3, with Dirichlet boundary D Γ as well as Neumann boundary N Γ , and by zero level-set is Lipschits continuous.We introduce an different parts of the design domain as shown in Fig.1.