Shape and topology sensitivity analysis and its application to structural design

The incremental sensitivity analysis associated with variation of structure material parameters, shape or topology variation is generally discussed by analyzing the evolution of potential and complementary energies, or arbitrary functionals of state fields. The concept of configuration and sensitivity generalized forces is used in presenting the sensitivity derivatives. The general reciprocity relations are derived for the case of potential or complementary energy variations. The topology variations in bar structures related to introduction of elements and introduction of inclusions and voids in plates are discussed, and the sensitivity forces are derived.

on the boundary portion S T and by the displacements u 0 = λ k u 0 (x) on the boundary portion S u , where x denotes a certain point. The loading process is now induced by varying λ s and λ k with the structure parameter vector fixed. The structure transformation process occurs for varying vector s and fixed loading factors λ s and λ k . For the elastic structure the stationary conditions of the potential energy = (u, λ s , s) or of the complementary energy = σ , λ k , s specify solution of the boundary value problem, thus δ = ,u u, λ s , s · δu = 0 ( 1 ) provides the equilibrium equations in a class of kinematically admissible fields u k and variations δu. Here ,u is the Frechet derivative of the potential energy functional with respect to the displacement vector at the actual state u k = u and δ denotes its variation. The dot between symbols denotes the scalar product. For a class of statically admissible stress fields σ and variations δσ , there is where ,σ is the Frechet derivative with respect to σ . Now (2) provides the compatibility conditions satisfied at the actual stress state σ . When the loading or displacement factors vary, the structure follows the loading path C l through the constitutive equilibrium states. The incremental equations along the loading path follow from (1) and (2), namely ,uuu and ,σ σσ whereu l ,σ l andλ s ,λ k are the increments or rates with respect to a selected evolution parameter. The product with respect to all indices of lower rank tensor applies for symbols put in juxtaposition. On the other hand, when the loading parameters are fixed but the structural vector s varies, the structure states evolve along the transformation path C s . In view of (1) and (2), we have ,uuu s + ,usṡ · δu = 0 ( 5 ) and ,σ σσ whereu s ,σ s are the displacement and stress increments induced by the structural vector incrementṡ. The variation of potential and complementary energies along the transformation path arė ,s = ,u · u s +˙ ,s ·ṡ = ,s ·ṡ = −S ·ṡ, ,s = ,σ · σ s + ,s ·ṡ = ,s ·ṡ = S ·ṡ.
In view of (5) the following equalities can be written for incremental paths of two structural parameters s 1 and s 2 ,uuu or −Ṡ 1 ·ṡ 2 = ,suu We have thereforeṠ 1 ·ṡ 2 =Ṡ 2 ·ṡ 1 (13) and the reciprocity relations occur for the increments or rates of configurational forces and increments or rates of conjugated structural parameters. These reciprocity relations are valid for elastic linear or nonlinear materials. Such relations were already extensively discussed by Herrmann and Kienzler [3,4] and Kienzler and Herrmann [5], but the proof presented here is general and provides the expressions for the second-order sensitivities of energy functionals. Here the structure evolution vector s can present material parameter, shape and topology variations. Let us note that s can now be a scalar or vector parameter of structure modification. The reciprocity relations for rates of loading and generalized conjugate displacements are derived from the equilibrium Eq. (3) along the loading path. For two loading paths associated with load factorsλ s 1 andλ s 2 we have ,uuu Along the loading path the variation of potential energy is expressed as follows: where q = − ,λ is the generalized displacement conjugate to the load factor λ s . The second variation of in view of (3) equals ,λ = ,uuu λ ·u λ + 2 ,uλ ·u λλs + ,λλ λ s 2 = − ,uuu λ ·u λ + ,λλ λ s 2 = −qλ s .
The reciprocity relations are now expressed in the form The relations (13) and (18) apply to both linear and nonlinear elastic structures. In particular, when the load factor λ s occurs linearly in the potential energy, the second terms in (16) and (17) vanish. Consider now an arbitrary functional and its variation induced by the structure evolution along the transformation path whereu s satisfies the incremental equilibrium Eq. (3). The variationu s can be specified by solving (3) for each specific case of structure variation in a direct sensitivity approach. However, a more effective adjoint structure approach can be proposed by considering an augmented functional where the multiplier u a is a kinematically admissible displacement field. The variation of G along the transformation path now is expressed as follows: Let us specify an adjoint problem providing the adjoint field u a . The first-order sensitivity of the functional (19) is now presented in a simple formĠ whereŜ = G ,s − ,us u a is the sensitivity force conjugate to s and depending on both primary and adjoint fields. The reciprocity relations (13) now cannot be proved for a general case. In the following, the sensitivity expressions will be derived by applying the adjoint structure approach. The details of derivation of the second-order sensitivity for arbitrary state functionals and of the potential energy have been discussed by Mróz [6].

Topological sensitivity analysis for bar structures
In this section we shall derive the sensitivity forces associated with variation of topological parameters for bar structures. In the structural context this parameter corresponds to introduction or removal of a load carrying element or support. When the sensitivity derivative does not vanish for the vanishing cross-section of the element, it provides an assessment of the potential force versus the cost gradient and induces structure modification. The infinitesimal and finite topology modifications can then be considered, cf. Bojczuk and Mróz [7][8][9].

Topological derivative for bar structures in regular states
Consider now a plane bar structure made of n − 1 rectilinear elements and loaded in-plane. We assume that this structure can be modified by introduction of additional n-th element of length l. The sensitivity of functional G with respect to introduction of this element can be expressed as follows: where S A = E n A n and S I = E n I n are the extensional and flexural stiffnesses of the virtual element; E n , A n and I n denote Young's modulus, cross-section area and moment of inertia of virtual bar and are the topological derivatives of functional G with respect to extensional and flexural stiffnesses of virtual element.

Now, consider the functional of generalized displacements of the form
where analyzed structure consists of n elements including one virtual element, and f denotes function of stiffnesses moduli S A , S I of virtual element and of vector U of generalized displacements, whose components u and w are the longitudinal and transverse displacements and θ is the rotation angle. Introduce an adjoint structure of the same configuration and boundary conditions as the considered structure and respectively loaded by distributed longitudinal and transverse forces and distributed moments, namely In view of (28) the first variation of functional G with respect to introduction of virtual element can be presented as follows: Taking into account the virtual work equation next the complementary virtual work equation and the relations the sensitivity of the functional G with respect to introduction of virtual bar now takes the form so it is expressed in terms of axial strain fields ε, ε a and curvature fields κ, κ a in the virtual element of primary and adjoint structures. Moreover, N , N a , and M, M a are the conjugate normal forces and bending moments, and superscript a denotes quantities in the adjoint structure. Now, we express fields of strains and curvatures of the virtual element by kinematic quantities on its ends. Assume that displacements in longitudinal and transverse directions and rotation angle at the end 1 of the virtual element n in local coordinate system (x n , y n ) are u 1 , w 1 , θ 1 , while at the end 2 respectively u 2 , w 2 , θ 2  Fig. 1). Using transformation formulae, the relations between generalized displacements in local coordinate systems (x i , y i ), (x j , y j ) related to bar elements i, j and in the coordinate system (x n , y n ) are of the form where u 1i , w 1i , θ 1i are generalized displacements of the point 1 in coordinate system (x i , y i ), u 2 j , w 2 j , θ 2 j are generalized displacements of the point 2 in coordinate system (x j , y j ), while α 1 , α 2 are the angles between the axes x i , x j and the axis x n . Taking into account that δS A = const, δS I = const along virtual element and that virtual element is unloaded, the equilibrium equations for this element are of the form Integrating the first Eq. (35) in view of the boundary conditions where l denotes length of the virtual bar, we obtain the constant field of virtual strain, namely Next, integrating the second Eq. (35), we obtain curvature, transverse displacement and rotation angle of the virtual element in the form In order to determine constants A, B, C and D the boundary conditions on the ends of virtual element can be used, namely Then, we get Now, the field of the virtual curvature takes the form Analogously, in the virtual element of the adjoint structure, we have Moreover, curvatures at the ends 1 and 2 of the virtual element of the primary and adjoint structures are Now, taking into account (37), (42) 1 and (43), the sensitivity (33) of the functional G with respect to introduction of virtual bar can be finally presented in the form When the functional G expresses the potential energy , namely the problem becomes self-adjoint and sensitivity (44) takes the form

Topological derivative of buckling load for bar structures
Analogously as in Sect. 3.1, consider now plane bar structure made of n − 1 rectilinear elements and loaded in-plane. We assume that this structure can be modified by introduction of additional n-th element of length l, extensional stiffness S A = E n A n , located transversely to one of the existing elements and playing role of the elastic support simply connected at the ends. The critical state condition for buckling now is where κ = w and θ = w is the curvature and rotation angle of the eigenmode w, N denotes normal force, E, A, I are the Young's modulus, cross-sectional area and moment of inertia. Moreover, ε n is the virtual strain and taking into account that where w n denotes the eigenmode deflection of the structure at the supported point, and in view of (37), we have Now, calculating variation of (47) with respect to introduction of virtual support along the path of critical states in view of (49) and assuming that S A = const, we obtain Using the virtual work equation along the path of critical states and taking into account that the third term under the integral in (50) related to redistribution of normal forces vanishes in this case (cf. Mróz and Bojczuk [10]), the topological derivative of critical value of load parameter with respect to introduction of virtual element finally can be presented as follows:

Topological sensitivity analysis for plates
In this section we shall derive the sensitivity forces associated with variation of topological parameters for plates. Here, the composite material design is considered and the topological modification is related to introduction of inclusions of new materials or introduction of voids.

Topological derivative with respect to introduction of inclusions
Consider now an elastic plate, whose middle surface occupies the domain A⊂R 2 , with the boundary Γ = Γ u ∪ Γ T , Fig. 2. The plate is loaded by tractions T = T 0 (T 0 i = σ i j n Γ j , i, j = 1, 2) on the boundary portion Γ T and by body forces p 0 1 , p 0 2 applied respectively in the plate domains A ξ and is the unit vector normal to the boundary Γ T . The stress and strain states occurring in the plate domain can be presented in the vector form as σ = [σ 11 , σ 22 , σ 12 ] T , ε = [ε 11 , ε 22 , 2ε 12 ] T . Moreover, displacements u = u 0 are specified on the boundary portion Γ u .
The topological derivative of the functional G with respect to the hole or inclusion area is defined as follows (cf. Sokołowski andŻ ochowski [11,12], Bojczuk and Mróz [13]) where x(x ∈ A) is an arbitrary position in the plate domain, in which the derivative is specified, denote respectively the domain of inclusion (void) and the reference domain, ξ is the expansion parameter, which specifies the size of the modification and A ξ = A − A 0 . Let us consider the response functional of the form where u 1 and u 2 are the displacement fields respectively in A ξ and A 0 domains. In order to specify the sensitivity derivative with respect to introduction of inclusion the variational approach is applied here. Now, the first variation of functional (54) in the expansion process of inclusion located at the material point x is expressed as follows where n = [n 1 , n 2 ] T is the unit vector normal to the interface Γ ξ , δϕ = [δϕ 1 , δϕ 2 ] T is the interface transformation vector, [ f ] = f 1 − f 2 denotes "jump" of the quantity f on the interface Γ ξ , and comma preceding an index denotes partial derivative. Following the previous derivation for plates cf. Bojczuk and Mróz [9] and the general methodology of sensitivity analysis cf. Dems and Mróz [14], the variations of state fields can be eliminated by introducing an adjoint plate structure of the same form, as the primary plate, but with induced body force fields, namely and satisfying the following boundary conditions where σ a , u a , ε a are the state fields in the adjoint structure. Taking into account continuity conditions on the interface Γ ξ where n, t denote respectively normal and tangential direction to this interface and using virtual work equation next, the complementary virtual work equation we finally obtain Let us introduce a general curvilinear coordinate system, where the first of coordinates ξ attains constant value on the interface Γ ξ and the second θ varies in the interval 0; 2π . In the particular case let us choose elliptical coordinates, which are related with the local coordinate system x 10 , x 20 coinciding with the ellipse semi-axes ξ a, ξb by relationships Then, taking ξ as the design parameter, (61) takes the form To express the sensitivity derivative with respect to the inclusion area A = πξ 2 ab, the following incremental formula is used d A = 2πξabdξ , thus When the last term in (54) disappears, [p 0 ] = 0, while the two first terms correspond to the total strain energy U , we get Now, let us consider the particular case, namely topological sensitivity derivative with respect to introduction of void. Taking into account that on the free boundary of the void σ nn = 0, σ nt = 0 the general formula (64) now takes the form When the functional G represents the strain energy U and p 0 = 0, the formula (65) can be rewritten as follows In order to use derived formulae for topological derivative, the stress or strain distributions on the interface should be known. In many cases, they can be analytically determined using methods of the elasticity theory (cf. Mura [15]; Muskhelishvili [16]), or they can be specified numerically. Let us consider an infinite plate under plane state of stress and containing a circular inclusion made of different material, with its center introduced at a point x. We assume that Young's modulus, Kirchhoff's shear modulus and Poisson's ratio are denoted respectively for matrix by E 1 , G 1 , ν 1 and for inclusion by E 2 , G 2 , ν 2 . The displacement and stress components appearing in the matrix, in the coordinate system n, t, are (cf. Muskhelishvili [16]) and where σ 1 , σ 2 denote the principal stresses in the neighbourhood of the point x for plate without inclusion (hole), r, θ are the polar coordinates, R is the radius of inclusion (Fig. 3) and β 1 , γ 1 , δ 1 are the unknown constants. The displacement and stress components expressed by (68), (69), can be presented on the boundary between matrix and inclusion i.e. for r = R, in the form The displacement and stress components occurring in the inclusion are (cf. Muskhelishvili [16])

Application of shape and topology sensitivity analysis to structural design
Consider now general optimization problem defined as follows where μ (μ ≥ 0) is the Lagrange multiplier and s i , i = 1, 2, ..., n −1 are the design parameters. The sensitivity derivatives G ,s i are specified by the formulae derived in the paper. The cost derivatives C ,s i can easily be derived for the assumed cost function expressed in terms of design parameters. The optimal values of the parameters s i , i = 1, 2, ..., n − 1 and the multiplier μ are determined in the incremental process of gradient optimization. Next, using the sensitivity derivative with respect to topological parameter s n the condition of topology modification can be formulated. When the condition is satisfied, standard optimization with respect to s i , i = 1, 2, ..., n should be performed again.

Formulation of topology modification conditions for bar structures
Consider first topology modification by introduction of virtual element into a bar structure in the regular state. Assume, that the functional G is expressed by (27) and the cost of structure is where C (n−1) denotes cost of the structure without additional bar, while c, S A , l are respectively unit cost, extensional stiffness and length of virtual element. Now, the condition of topology modification by introduction of a new bar, in view of (44), takes the form and it corresponds to non-positive value of topological derivative of the Lagrangian (78), where μ ≥ 0 is the Lagrange multiplier calculated for the unmodified structure.
Consider now the problem of the form (77), but with buckling constraint cost λ − λ cr ≤ 0, i.e. G ≡ λ, G 0 ≡ λ cr , and λ cr denotes critical value of load parameter λ. We assume that the bar structure is modified by introduction of additional n-th element of unit cost c, length l, extensional stiffness S A , which plays the role of the elastic, transverse support simply connected on the ends. Now, the topological derivative of the Lagrangian L = C + μ (λ − λ cr ), in view of (52), provides the condition of new support introduction in the form where w n is the eigenmode deflection of the beam at the supported point, θ = w denotes the slope and μ ≥ 0, as previously, is the Lagrange multiplier calculated for the unmodified structure.

Formulation of topology modification conditions for plates
In the case of the problem (77) for plates, the condition of introduction of an infinitesimally small inclusion (hole) at the point x, takes the general form (Bojczuk and Mróz [8]) where T L ,A 0 (x), T C ,A 0 (x), T G ,A 0 (x) are the topological derivatives at the point x respectively of Lagrangian, cost functional C and functional G, and this last derivative is expressed by (64) -(67) and (76). So, an inclusion (hole) should be introduced at a point, where the condition (83) is satisfied and T L ,A 0 (x) attains minimal value (cf. bubble method, Eschenauer et al. [17]). This condition enables localization of introduced new structural element.
However, in order to accelerate the optimization process, finite topology modifications can be applied. Now, the problem consists in introduction of finite holes of unknown size and shape together with introduction of finite changes of other boundaries. It is assumed that domains of relatively small values of the topological derivative of Lagrangian, which is expressed by (83), should be eliminated. It can be done using level-set method by removal of all domains, where this topological derivative is smaller than an adequately chosen negative iso-value (cf. Burger et al. [18], Wang et al. [19]). However, also alternative approach can be applied, where the specially constructed design quality function is used, which contains the scaling factor controlling amount of removed domain related to the respective iso-value (cf. Bojczuk and Mróz [8]).

Concluding remarks
The present paper provides the general framework of specification of sensitivity derivatives of potential or complementary energies and of arbitrary state functionals, especially displacement functionals. Here, also topological derivatives with respect to introduction of additional elements into bar structures and inclusions or voids into plates are taken into account. These sensitivity derivatives can be regarded as configurational or sensitivity forces generating effective redesign or optimal design procedures. The reciprocity relations naturally follow from the second-order sensitivity expressions of energies and increments of configurational forces.
Open Access This article is distributed under the terms of the Creative Commons Attribution License which permits any use, distribution, and reproduction in any medium, provided the original author(s) and the source are credited.