A computational homogenization approach for ... - Wiley Online Library

4 downloads 6352 Views 2MB Size Report
a Department of Civil Engineering, International University - VNU HCM, Viet Nam ... and Civil Engineering, University of Technical Education - HCMC, Viet Nam.
Accepted Article

INTERNATIONAL JOURNAL FOR NUMERICAL METHODS IN ENGINEERING Int. J. Numer. Meth. Engng 2017; 00:1–23 Published online in Wiley Online Library (www.onlinelibrary.wiley.com). DOI: 10.1002/nme.5561

A computational homogenization approach for limit analysis of heterogeneous materials Canh V. Lea∗ , Phuong H. Nguyenb , H. Askesc , D.C. Phamd

a b

c

Department of Civil Engineering, International University - VNU HCM, Viet Nam Faculty of Applied Mechanics and Civil Engineering, University of Technical Education - HCMC, Viet Nam University of Sheffield, Department of Civil and Structural Engineering, Sheffield S1 3JD, United Kingdom d Institute of Mechanics - VAST, 264 Doi Can, Hanoi, Vietnam

SUMMARY

The macroscopic strength domain and failure mode of heterogeneous or composite materials can be quantitatively determined by solving an auxiliary limit analysis problem formulated on a periodic representative volume element (RVE). In this paper, a novel numerical procedure based on kinematic theorem and homogenization theory for limit analysis of periodic composites is developed. The total velocity fields, instead of fluctuating (periodic) velocity, at microscopic level are approximated by finite elements, ensuring that the resulting optimization problem is similar to the limit analysis problem formulated for structures, except for additional constraints which are periodic boundary conditions and averaging relations. The optimization problem is then formulated in the form of a standard second-order cone programming problem to be solved by highly efficient solvers. Effects of loading condition, RVE architecture and fiber/air void volume fraction on the macroscopic strength of perforated and fiber reinforced composites are studied c 2017 John Wiley & Sons, Ltd. in numerical examples. Copyright ⃝ Received . . .

KEY WORDS: Limit analysis, computational homogenization, periodic composite, second-order cone programming

1. INTRODUCTION

Composites and heterogeneous materials are widely used in many structural elements in the fields of civil and mechanical engineering, and the estimation of their effective properties plays an important role in structural design and safety assessment of structures under quasi-static loadings. ∗ Correspondence

to: Canh V. Le, Department of Civil Engineering, International University - VNU HCM, Viet Nam, e-mail: [email protected]

This article has been accepted for publication and undergone full peer review but has not been through the copyediting, typesetting, pagination and proofreading process, which may lead to differences between this version and the Version of Record. Please cite this article as doi: 10.1002/nme.5561 This article is protected by copyright. All rights reserved.

Accepted Article 2

CANH V. LE ET AL.

The computation of the collapse load of a structure can be performed using the elastic-plastic incremental analysis method. However, for the determination of the ultimate strength the plastic limit analysis approach has been proved to be more effective [1, 2], i.e. without performing a cumbersome series of incremental elastic-plastic analysis and without knowing the history of the whole loading path. In order to correctly capture the effects of inhomogeneities, limit analysis of micro-structures, based on upper/lower bound theorems and homogenization theory, has been developed and gained an increasing attention in recent years. Homogenization theory in limit analysis has been firstly introduced to determine a macroscopic strength criterion for fiber-reinforced composite materials in [5–7]. A computational homogenization limit analysis based on finite elements and linear programming has been proposed in [8] to calculate the macroscopic strength of composites governed by Tresca’s criterion. The behavior of cylindrical porous material was studied in [9, 10] using Gurson’s model together with coupled kinematic and static theorems as well as homogenization theory. Based on finite elements and primal-dual interior point optimization algorithm, the assessment of the macroscopic strength criterion and stability analysis of soils reinforced by stone columns were performed in [11–15]. For periodic plates, a theoretical study on the overall homogenized Love-Kirchhoff strength domain of a rigid perfectly plastic multi-layered plate was reported by [16, 17], and numerical determination of the macroscopic bending strength criterion was presented in [18, 19]. By means of homogenization techniques and kinematic limit analysis in conjunction with non-linear programming, the plastic limit loads and failure modes of periodic composites governed by ellipsoid yield criteria can be determined [20–25]. Using the elastic stresses of the periodic microstructure, a static direct method in combination with homogenization was proposed in [26–30] for 2 and 3-dimensional limit analysis of periodic metal-matrix composites. In the method, the strong form of the equilibrium equations was transformed into the so-called weak form, and to be satisfied locally in an average sense using approximated virtual displacement fields. Based on linear matching method, limit analysis of reinforced concrete structures has also been studied in [31–35]. Practical multi-scale modelling and homogenization can be found in [3]. Direct approaches for limit analysis of structures and materials lead to a problem of constrained optimization involving either linear or nonlinear programming, and hence the development of efficient optimization algorithms to enable solutions to be obtained is one of the major research directions in the field. From a mathematical point of view, linear programming is very attractive, and has been widely applied to limit analysis problems [38–40]. Although powerful software based on interior point algorithms is available for the resolution of the resulting optimization problem, a large number of constraints generated in the linearization process would be needed in order to provide accurate solutions, thereby increasing the computational cost. Attempts have been made to solve problems involving exact convex yield functions using non-linear programming algorithms [41–46]. In recent years, thanks to the development of such a efficient primal-dual interior point algorithm, implemented in a commercial optimization software Mosek [54], limit analysis problems have gained increasing attention. In fact, most commonly used yield criteria can be cast in the form of conic or semi-definite constraints, and optimization problems involving such constraints can be solved efficiently by such an algorithm, evidenced by recent works [13–15, 18, 19, 47–52]. The objective of the present paper is to develop a novel computational homogenization procedure for kinematic limit analysis of periodic materials. The total velocity fields at microscopic level This article is protected by copyright. All rights reserved.

Accepted Article

COMPUTATIONAL HOMOGENIZATION APPROACH FOR LIMIT ANALYSIS

3

are approximated by finite elements, instead of the fluctuating (periodic) velocity as in [21–25]. Consequently, the auxiliary limit analysis problem of microstructure is similar to a limit analysis problem formulated for structures. The plastic dissipation or the objective function can be easily transformed into a form involving a sum of Euclidean norms. Additional constraints, periodic boundary conditions and averaging relations can be treated naturally as equalities. The corresponding optimization problem is then formulated in the form of a standard second-order cone programming problem, for which highly efficient solvers are available. Note that the study of the overall macroscopic behavior of composite materials and RVE existence using an incremental computational homogenization approach based on total displacement approximation was reported in [36, 37]. In [10], the total velocity fields were also used in combination with linear displacement boundary condition for limit analysis of cylindrical porous material. In this paper the upper bound limit analysis of periodic materials is formulated in terms of the full velocity fields and periodic boundary conditions. The main contribution of the present work is to handle the corresponding finite element based discrete optimization problem using conic programming, enabling the efficient computation of a large number of points to describe a yield surface.

2. FUNDAMENTALS OF LIMIT ANALYSIS

In this section, fundamentals of limit analysis theory are recalled; more details can be found in [40]. Consider a rigid-perfectly plastic body of area Ω ∈ R2 with boundary Γ, and subjected to body forces f and to surface tractions g on the free portion Γt of Γ. The constrained boundary Γu is fixed and Γu ∪ Γt = Γ, Γu ∩ Γt = ⊘. Let X denote an appropriate space of statically admissible stress states, whereas Y is an appropriate space of kinematically admissible velocity states. For smooth fields σ and u, the classical form of the equilibrium equation can always be transformed to a variational form as a(σ, u) = F (u), ∀ u ∈ Y (1) where the external work rate and internal work rate can be expressed respectively in linear and bilinear forms as ∫ ∫ T F (u) = f u dΩ + gT u dΓ (2) Ω



Γt

σ T ϵ(u) dΩ

a(σ, u) =

(3)



where ϵ(u) = [ ϵxx ϵyy γxy ]T are strain rates. Note that in the framework of limit analysis, the strain rates are always of concern, not strains, and hence for simplicity we denote ϵ (without a dot) as strain rates and this rule is also applied for other kinematic quantities. This article is protected by copyright. All rights reserved.

Accepted Article 4

CANH V. LE ET AL.

The exact collapse load multiplier can be obtained by solving one of the following optimization problems λ =

max {λ− | ∃ σ ∈ B : a(σ, u) = λ− F (u), ∀ u ∈ Y }

max min a(σ, u) σ∈B u∈C = min max a(σ, u) u∈C σ∈B = min D(u)

(4) (5)

=

(6) (7)

u∈C

where the plastic dissipation rate D(u) is expressed in terms of σ and u as D(u) = max a(σ, u), σ∈B

(8)

C = {u ∈ Y | F (u) = 1}

(9)

B = {σ ∈ X | ψ(σ(x)) ≤ 0}

(10)

with ψ(σ) is the yield criterion. Problems (4) and (7) are, respectively, known as static and kinematic principles of limit analysis, for which the stress or displacement field must be discretized, respectively note that in order to obtain a strict lower bound on the actual collapse load multiplier the strong form of Equation (1) must be used and satisfied everywhere in the problem domain). On the other hand, the mixed formulations (5) and (6) require the approximation of both stress and displacement fields, and therefore mixed finite elements can be used. However, in the present work, the kinematic form (7) is of particular interest. Most commonly used yield criteria can be expressed as ψ(σ) =



σT Pσ − 1

(11)

where P is a coefficient matrix consisting of the strength properties of materials. For von Mises criterion, applied to isotropic materials and plane stress condition, P is given by 

 1 −1/2 0 1   P = 2  −1/2 1 0  σ0 0 0 3

(12)

where σ0 is uniaxial yield stress. On the other hand, for anisotropic materials Hill’s yield criterion is often used, and P is expressed by   P=

αzx + αxy

−αxy

−αxy 0

αxy + αyz 0

0



 0  3ηxy

(13)

where αzx , αxy , αyz and ηxy are material constants characteristic of anisotropy, and determined by

This article is protected by copyright. All rights reserved.

Accepted Article

COMPUTATIONAL HOMOGENIZATION APPROACH FOR LIMIT ANALYSIS

2αzx

=

2αxy

=

2αyz

=

ηxy

=

1 1 1 + − σ2pz σ2px σ2py 1 1 1 + − σ2px σ2py σ2pz 1 1 1 + − σ2py σ2pz σ2px 1 3τ2pxy

5

(14) (15) (16) (17)

in which σpx , σpy and σpz are the uniaxial yield stresses in three orthotropic directions, and τpxy is the shear yield stress. In the framework of a limit analysis problem, only plastic strain rates are considered and are assumed to obey the normality rule. Therefore, the power of dissipation can be formulated as a function of strain rates as ∫ √ D(ϵ) = ϵT Θ ϵ dΩ (18) −1

where Θ = P



.

3. HOMOGENIZATION THEORY

Consider a heterogeneous macro-continuum V ∈ R2 with associated micro-structures Ω ∈ R2 , as shown in Fig. 1. As for homogenization theory [4], one of the key assumptions is that the mechanical problem contains two separate scales: macroscopic (or overall) and microscopic (or local) scales. This means that the microscopic scale is small enough, compared with the one of the structure, for the heterogeneities to be separately identified, and the macroscopic scale is large enough for the effects of the heterogeneities to be smeared-out. The connection between the two scales involves two stages: down-scaling (or localization) and up-scaling (or globalization). Downscaling is the macro-to-micro translation, in which the macroscopic quantities are transformed to the RVE as boundary conditions. The inverse procedure, the micro-to-macro translation by which the microscopic properties originating in the RVE are averaged at the macroscopic level, is called up-scaling. The macroscopic quantities are linked to microscopic ones by the average relations 1 E ≡ ⟨ϵ⟩ = |Ω| 1 Σ ≡ ⟨σ⟩ = |Ω|

∫ ϵ dΩ

(19)

σ dΩ

(20)



∫ Ω

where ϵ and σ denote the local strain rate and stress terms; E and Σ are the average strain rate and average stress terms over the RVE; the operator ⟨·⟩ stands for the volume average of a field over the RVE, and |Ω| is the area of the RVE. For a RVE located at a sufficiently large distance from the boundary of the heterogeneous body and subjected to boundary conditions that would produce uniform fields of strain rate E and stress This article is protected by copyright. All rights reserved.

Accepted Article 6

CANH V. LE ET AL.

‫ݔ‬ଶ

‫ݔ‬ଵ

ܺଶ

ܺଵ

Heterogeneous material

Down-scaling

RVE

ܺଶ

Up-scaling

ܺଵ

Homogenized material

Figure 1. Homogenization of periodic composites

Σ in the associated homogeneous body, the microscopic strain rate and stress fields conform at the

microscopic level to the periodicity of the geometry, or in other words ϵ and σ are periodic fields. This leads to the fact that at the microscopic scale in the RVE, the local mechanical fields can be split into its average and a fluctuating term which accounts for the presence of heterogeneities as ˜ (x) u(x) = E · x + u

(21)

˜ (x) ϵ(x) = E + ϵ

(22)

˜ (x) σ(x) = Σ + σ

(23)

˜ and σ ˜ are the fluctuation terms of the local velocity, strain rate and stress, ˜ (x), ϵ where u ˜ is periodic, respectively. Periodicity requirements imply that i) the local fluctuation velocity u and ii) the traction vectors σ · n are opposite on opposite sides of the RVE boundary (where the external normal vectors n are also opposite) [4], or anti-periodic. There are three classical boundary conditions commonly used in homogenization theory: (i) linear displacements, (ii) constant tractions and (iii) periodicity conditions. However, only periodic boundary conditions will be considered in this paper. Direct consequences of Eq. (20) are ˜ ⟩ = 0; ⟨σ ˜⟩ = 0 ⟨ϵ

(24)

For any admissible velocity and stress field respectively satisfying periodic and anti-periodic conditions mentioned previously, the following relation holds Σ : E = ⟨σ : ϵ⟩

(25)

This virtual work principle is termed Hill-Mandel’s macro-homogeneity equality, and plays an important role in homogenization theory. This article is protected by copyright. All rights reserved.

Accepted Article

COMPUTATIONAL HOMOGENIZATION APPROACH FOR LIMIT ANALYSIS

7

4. COMPUTATIONAL HOMOGENIZATION FOR KINEMATIC LIMIT ANALYSIS

4.1. Limit analysis based on homogenization theory

Consider a ductile composite, all the constituents of which are rigid-perfectly plastic and obey the associative normality rule. The effective properties of such a composite can be determined by solving a specific auxiliary problem formulated on the RVE of periodic micro-structure, taking account the connection between macro and micro-quantities. In [26–30] the elastic stresses of the periodic microstructure were calculated according to different types of prescribed loading conditions, either prescribed stress Σ or prescribed strain E, and then macroscopic strength domains were approximated using the static direct method. Recently, Bleyer et al. [18] proposed a novel numerical procedure for homogenization of periodic plates, in which kinematic theorem was combined with the prescribed curvatures , and prescribed moments were employed in the static formulation. However, in this paper we will present a novel computational homogenization limit analysis based on the kinematic approach described in [20–25]. On the boundary of the RVE, the external prescribed loading condition considered are macroscopic strains, which are treated as variables. The set of macroscopic stresses which can be carried by the heterogeneous material characterizes the strength of the homogenized material. In the present study, the kinematic formulation for limit analysis of composite materials will be reformulated in terms of total strain rates (or total velocities), and hence the resulting formulation is different from those presented in [20–25], in which the fluctuation strain rates/velocities are the problem variable fields. Let the total strain rates/velocities be the problem variable fields, the computational homogenization formulation for kinematic limit analysis of a periodic micro-structure can be expressed as λ

+

= min

∫ √

ϵT Θ ϵ dΩ

(26a)



s.t

F (u) = ΣT E = 1

(26b)

u(x) − E · x periodic ∫ 1 ϵ dΩ E= |Ω| Ω

(26c) (26d)

The problem statement of Eqs. (26) is closely similar to a limit analysis problem formulated for structures, the main difference is that periodic boundary conditions and averaging relations are taken into account, i.e. Eqs. (26c) and (26d) (To avoid redundancy, Eq. (26c) should not be enforced for corner nodes, see details in Section 4.2). Note that the stress should be divergencefree at microscopic scale, and hence no body forces appear in the microstructural limit analysis formulation. Furthermore, in the present formulation only continuous velocity fields are considered although it may be possible to consider discontinuous fields and take into account the discontinuities contribution in the dissipated power.

4.2. Discrete formulation by finite elements In the kinematic formulation, plastic strain rates must be approximated using discretization methods, e.g. finite element method (FEM). The problem domain is discretized into elements according to This article is protected by copyright. All rights reserved.

Accepted Article 8

CANH V. LE ET AL.

Ω ≈ Ω1 ∪ Ω2 ∪ · · · ∪ Ωnel and Ωi ∩ Ωj = ⊘, i ̸= j , where nel is the number of elements. Then, the finite element approximation of velocity and plastic strains at the microscopic scale in the RVE can be expressed as uh = Nd; ϵh = Bd

(27)

where N are the shape functions, B is the strain-displacement matrix and d are corresponding nodal velocity values. By means of finite element discretization and Gaussian integration technique, the objective (or the plastic dissipation) function of the problem (26) can be rewritten as DF EM =

NG √ ∑ ξi (Bi d)T ΘBi d

(28)

i=1

where ξi is the integral weight at the ith Gaussian integral point and N G is the total number of integration points over Ω. Similarly, constraint (26d) can be formulated as 1 ∑ ξi Bi d |Ω| NG

E=

(29)

i=1

Now, attention will be focussed on the treatment of periodic boundary condition of the RVE. To carry out this task, the boundary Γ of the RVE will be decomposed into two parts such that ∪ Γ = Γ+ Γ− with outward normals n+ = −n− at two associated points x+ ∈ Γ+ and x− ∈ Γ− , as shown in Fig. 2. The key kinematical constraint for the limit analysis of the micro-structure is that the velocity fluctuation must be periodic on the boundary of the RVE. This means that for each pair {x+ , x− } of boundary material points, we have

4

݊ଵି

‫ݔ‬ଶା Gା ଶ

3

݊ଵା

Gଵା

Gଵି

‫ݔ‬ଵି

1

‫ݔ‬ଶି Gି ଶ

2

‫ݔ‬ଵା

Interior nodes ା Nodes on G

ି Nodes on G

Corner nodes

Figure 2. Periodic boundary condition

˜ (x+ ) = u ˜ (x− ) u

(30)

This article is protected by copyright. All rights reserved.

Accepted Article

COMPUTATIONAL HOMOGENIZATION APPROACH FOR LIMIT ANALYSIS

9

By combining Eqs. (21) and (30), one obtains uΓ+ − uΓ− = u2 − u1 , uΓ+ − uΓ− = u4 − u1 1

1

2

2

(31)

where uΓ+ and uΓ− (i = 1, 2) are respectively the nodal velocities of nodes on the boundary Γ+ i and i i − Γi , and uc (c = 1, 2, 4) is the nodal velocities of the corner nodes. Note that Equation (31) ensures the periodicity of the fluctuation field at all points on the boundary of the RVE except the corner nodes, and kinematic conditions for these corner nodes are implicitly enforced in Equation (29). By assembling the global matrix form, constraints (31) may be cast as (32)

Cd = 0

in which C is a link-topology matrix, consisting of {0, 1, −1} only. Finally, problem (26) becomes λ+

= min

NG √ ∑ ξi (Bi d)T ΘBi d

(33a)

i=1

s.t

ΣT E = 1

(33b)

Cd = 0

(33c)

E=

1 |Ω|

NG ∑

ξi Bi d

(33d)

i=1

4.3. The solution of the discrete problem The problem for limit analysis of periodic micro-structure (33) is a non-linear optimization problem with equality constraints. Its objective function is not differentiable everywhere while powerful optimization algorithms require their gradients to be available everywhere. In the following, the problem will be formulated in the form of a second-order cone programming (SOCP), making sure that it can be solved using highly efficient primal-dual interior-point solvers. For plane stress problems, Θ is a positive definite matrix (shown in Section 2 for both von Mises’ and Hill’s criteria). The objective or the plastic dissipation function in (33) can be rewritten straightforwardly in a form involving a sum of norms as DF EM =

NG ∑

ξi QT Bi d

(34)

i=1

where ||.|| denotes the Euclidean norm appearing in the plastic dissipation function, i.e, ||v|| = (vT v)1/2 , and Q is the so-called Cholesky factor of Θ. By introducing a vector of additional variables ρi ρi = QT Bi d

(35)

and auxiliary variables t1 , t2 , . . . , tN G , the problem (33) can be cast in the form of second-order cone programming as This article is protected by copyright. All rights reserved.

Accepted Article 10

CANH V. LE ET AL.

λ+

= min

NG ∑

ξi ti

(36a)

i=1

s.t

ΣT E = 1

(36b)

Cd = 0

(36c)

NG 1 ∑ ξi Bi d E= |Ω|

(36d)

||ρi || ≤ ti , i = 1, 2, . . . , N G

(36e)

i=1

which involving equalities and conic constraints.

5. NUMERICAL APPLICATIONS

To illustrate the performance of the described computational homogenization strategy for limit analysis of periodic micro-structures, a number of benchmark plane stress problems, in some of which other numerical solutions are available for comparison, will be considered. The problem (36) is formulated within the Matlab environment and solved using the highly efficient primal-dual interior point solver, Mosek [54]. In all examples, square RVEs of a × a (a = 1 mm) are used. For a given set of macroscopic stresses Σ, load multipliers, λ+ , can be obtained by solving the discrete optimization problem (36). Therefore, points on the macroscopic strength domains are determined by Σy = λ+ Σ (37)

5.1. Perforated materials governed by von Mises’ criterion

Perforated materials can be treated as a special composite, and its load-bearing capacity plays an important role in the design of many engineering structures. The RVE of the perforated materials with a rectangular or circular hole at its center, as shown in Fig. 3, is considered. The RVE is subjected to a pair of macroscopic stress components (Σ11 , Σ22 ) in the plane (x1 , x2 ) as shown in Fig. 3. The matrix material for the rectangular perforated RVE is aluminium Al with yield stress σ0 = 137 MPa, while it is mild steel St3S with yield stress σ0 = 273 MPa for the circular one. These problems have been examined by Li in [21, 22] using a kinematic finite element formulation in combination with an iterative solution algorithm,and by Zhang et al. [55] using a quasi-lower bound approach. Finite element meshes of three-node triangular elements (T3) for the computational models are presented in Fig. 4. Numerical macroscopic strength domains of the perforated material with circular or rectangular hole for different values of α, α = 00 and α = 450 , are shown in Figs. 5 and 6. It is evident that the present numerical solutions are in good agreement with those obtained by Li et al. [22] and by Zhang et al. [55]. Note that the quasi-static method presented in [55] may provide a higher solution than the actual macroscopic strength. This article is protected by copyright. All rights reserved.

Accepted Article

11

COMPUTATIONAL HOMOGENIZATION APPROACH FOR LIMIT ANALYSIS

‫ݔ‬ଶ ‫ܮ‬ଶ

‫ܮ‬ଵ

‫ݔ‬ଶ

Sଵଵ ‫ݔ‬ଵ

Sଵଵ ‫ݔ‬ଵ

‫ݎ‬

ܽ

Sଶଶ

a

Sଶଶ

Figure 3. RVEs of the perforated material: rectangular (left, L1 × L2 = 0.1 × 0.5 mm) and circular (right, r/a = 0.25) holes

0.8

(b) 1752 T3 elements

Figure 4. RVEs of the perforated material: finite element meshes.

Result by Li et al. Result by Zhang et al. Present method

0.8 0.6

0.4

0.4

0.2

0.2 Σ22/σ0

Σ22/σ0

0.6

(a) 2038 T3 elements, L1 × L2 =

0.1 × 0.5 mm

0

−0.2

−0.2

−0.4

−0.4

−0.6

−0.6

−0.8

−0.5

0 Σ11/σ0

(a) α = 00

0.5

Result by Li et al. Result by Zhang et al. Present method

0

−0.8

−0.5

0 Σ11/σ0

0.5

(b) α = 450

Figure 5. Macroscopic strength domain of the perforated material: circular hole

The macroscopic uniaxial strength for various values of α and different sizes of the hole is plotted in Figs. 7 and 8. As expected the effective strength of the perforated materials decreases significantly as the size of microscopic hole increases. Associated failure mechanisms and plastic dissipation distribution of the perforated materials are also presented in Figs. 9 and 10, where the color contours show the plastic dissipation at limit state. This article is protected by copyright. All rights reserved.

Accepted Article 12

CANH V. LE ET AL.

1

0.4 0.3 0.2 0.1 Σ22/σ0

Σ22/σ0

0.5

0

−0.5

−1 −0.4

−0.2

0 −0.1 −0.2 −0.3

0 Σ11/σ0

0.2

−0.4 −0.4

0.4

−0.2

(a) α = 00

0 Σ11/σ0

0.2

0.4

(b) α = 450

Figure 6. Macroscopic strength domain of the perforated material: rectangular hole

1

0.8

0.8

0.6

0.6

0.4

Experimentals by Litewka et al. Present method Result by Li et al.

0.2

0 0

Σ11/σ0

Σ11/σ0

1

20

40

α

60

(a) L1 × L2 = 0.1 × 0.5

80

Experimentals by Litewka et al. Present method Result by Li et al.

0.4

0.2

0 0

20

40

α

60

80

(b) L1 × L2 = 0.1 × 0.7

Figure 7. The macroscopic uniaxial strength for various values of α and different sizes of the hole: rectangular hole.

It can be observed that for all cases the yield zone extends from the hole’s boundary towards the boundary of the RVE. Next, the efficacy of the proposed computational homogenization for kinematic limit analysis of periodic materials will be studied. Three meshes of elements were considered for the computation of the effective properties of the RVE with rectangular hole (L1 × L2 = 0.1 × 0.5). From Table I, we can see that the present numerical homogenization limit analysis can produce solutions very quickly by the cone-based optimization algorithm, just taken less than 6 seconds to solve the problem of up to 40000 variables. The macroscopic strength domains of the perforated samples under three independent macroscopic loads (Σ11 , Σ22 , Σ12 ) are plotted in Figs. 11 and 12. As expected, for the RVE with rectangular hole the macroscopic strength in Σ22 direction is greater than that in Σ11 direction, while they are equivalent for the RVE with circular hole. This article is protected by copyright. All rights reserved.

Accepted Article

13

COMPUTATIONAL HOMOGENIZATION APPROACH FOR LIMIT ANALYSIS

Present method, α = 00 Result by Li et al., α = 00

1

Present method, α = 450 0

Result by Li et al., α = 45

11

Σ /σ

0

0.8

0.6

0.4

0.2

0 0

0.1

0.2

0.3

0.4

0.5

r/a

Figure 8. The macroscopic uniaxial strength for various values of α and different sizes of the hole (r/a): circular hole. 30 18

25 16 14

20

12

15

10 8

10 6 4

5

2

(a) α = 00

(b) α = 450

Figure 9. The failure mode of the perforated material: rectangular hole (L1 × L2 = 0.1 × 0.5). 16 12 14 10

12 10

8

8 6 6 4

4 2

(a) α = 00

2

(b) α = 450

Figure 10. The failure mode of the perforated material: circular hole (r/a = 0.25).

This article is protected by copyright. All rights reserved.

Accepted Article 14

CANH V. LE ET AL.

Table I. On the performance of the proposed method

Present periodic solutions sdof CPU time (s∗ ) Σ11 /σ0 642