- Email: [email protected]
Three-dimensional periodic Voronoi grain models and micromechanical FE-simulations of a two-phase steel
Computational Materials Science 24 (2002) 513–519 www.elsevier.com/locate/commatsci
Three-dimensional periodic Voronoi grain models and micromechanical FE-simulations of a two-phase steel Mikael Nyg ards, Peter Gudmundson
*
Department of Solid Mechanics, Royal Institute of Technology (KTH), 100 44 Stockholm, Sweden Received 5 July 2001; received in revised form 18 December 2001; accepted 18 December 2001
Abstract A three-dimensional model is proposed for modeling of microstructures. The model is based on the finite element method with periodic boundary conditions. The Voronoi algorithm is used to generate the geometrical model, which has a periodic grain structure that follows the original boundaries of the Voronoi cells. As an application, the model is used to model a two-phase ferrite/pearlite steel. It is shown that periodic cells with only five grains generate representative stress–strain curves. Ó 2002 Elsevier Science B.V. All rights reserved. Keywords: FEM; Microstructure; Voronoi cells; Plasticity; Pearlite; Ferrite
1. Introduction It is desirable to design microstructures to get special properties, such as high strength, ductility or toughness. To do this design properly, the role of microstructure must be thoroughly understood, since most macroscopic behaviors can be related to representative parts of the microstructure. Micromechanical modeling gives the possibility to investigate certain microstructural properties in more detail, and relate those to the macroscopic behavior of the material. The main advantage with modeling is that the influence of various parameters easily can be investigated. A corresponding
*
Corresponding author. Tel.: +46-8-790-7548; fax: +46-8411-2418. E-mail address: [email protected] (P. Gudmundson).
experimental study would be both expensive and time consuming. The morphology of a three-dimensional microstructure decides to a large extent the modeling strategy. If the microstructure is elongated in one direction, a three-dimensional micromechanical model that is similar to the real material is fairly straightforward to generate [1,2]. When the microstructure is random in three dimensions, the problem can be approached by a simplification that still captures the features of interest in the material, like distributing particles to investigate clustering [3]. In polycrystalline materials, the characteristic feature is the grain structure. Grains may be simplified geometrically to capture certain features, for example as cubes [4] or as elongated twodimensional structures [5]. Alternatively, the microstructure can be mapped, and a corresponding
0927-0256/02/$ - see front matter Ó 2002 Elsevier Science B.V. All rights reserved. PII: S 0 9 2 7 - 0 2 5 6 ( 0 2 ) 0 0 1 5 6 - 8
514
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
microstructure can be simulated [6,7]. However, to model more realistic microstructures, it is desirable to use an algorithm to create the grain structure, and to ensure its randomness. For this purpose the Voronoi algorithm has been shown to be successful. The topology of the Voronoi tessellation, see for example [8], is similar to the microstructure of real materials. Therefore, it is advantageous to use for modeling of materials at a microstructural scale. The three-dimensional Voronoi algorithm has been applied to materials science by other researchers in the past. Local stresses within the generated Voronoi structure have been reported by e.g. Kozaczek et al. [9], Kumar et al. [10], Wu and Guo [11]. While Roberts and Garboczi [12] used the Voronoi structure to study the elastic moduli of a cellular solid. Statistical properties of the Voronoi structure have been investigated by Kumar et al. [13]. The three-dimensional modeling by Voronoi cells offers the unique possibility to automatically generate unit cells and it can therefore capture statistical data. Hence, effects of morphology may be studied. In the present model, a concept for modeling microstructures is presented in Section 2. A periodic unit cell that follows the original boundaries of a generated Voronoi structure is loaded by periodic boundary conditions. Thereby, grains are not divided to fit within a square unit cell, which is advantageous when the grain topology is of importance, also edge effects are avoided. In Section 3, a validation with an orthotropic elastic material model is performed. As an application, a two-phase ferrite/pearlite steel is modeled with an elastic–plastic analysis. In the model, an initial ferrite phase is created by the Voronoi algorithm, thereafter, the properties of the adaptive mesh is used to distribute the pearlite. Hence, a two-phase irregular microstructure is created that is similar to the real microstructure.
2. Numerical model On a microstructural level, polycrystalline materials consist of a more or less random grain structure. In this study the so-called Voronoi al-
gorithm is used to model microstructures in three dimensions. 2.1. Periodic Voronoi structure The idea behind a micromechanical model based on the Voronoi algorithm is to model a representative part of the microstructure and to predict the macroscopic behavior of the material. The cell should be large enough to catch the macroscopic behavior, and it should be kept as small as possible to save computation time. To avoid undesired edge effects, a periodic cell loaded by periodic boundary conditions has been chosen as the representative unit cell. The Voronoi algorithm utilizes a set of sites to partition space into regions; one per site. The region for a site s consists of all points closer to s than to any other site. In practice, the distances between the site s and all its closest neighbors are divided in half by planes perpendicular to a line connecting the sites. In this manner, a three-dimensional structure resembling a real microstructure is created. In this work, a public domain software called Qhull [14] has been used to create the initial Voronoi structure. The program is given a random set of sites as input, and it outputs the edge points, surfaces and regions of the Voronoi cells. A periodic Voronoi structure results if the set of sites show a periodic symmetry. This property has been utilized in the present study to generate periodic Voronoi structures. Firstly, a preset number of sites (grains within a periodic cell) is generated within a unit cube. The location of each site is set by a random generator with uniform distribution in the cube. The site locations within the cube are then copied to horizontally, vertically and diagonally located cubes. Lastly, the Voronoi algorithm is applied to the sites within 3 3 3 ¼ 27 cubes. A unique periodic cell can be chosen from the center of the created Voronoi structure. In Fig. 1, a periodic cell containing five grains is shown. It is evident that each grain is represented by irregular convex polyhedrons. By considering the nature of the Voronoi algorithm, there are two regions on each side of every surface, every edge has three regions around it, and each edge corner has four regions sharing that corner. In other words, there
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
515
periodic cell and the average strain tensor eij can be expressed as ubi uai ¼ eij ðxbj xaj Þ;
Fig. 1. A typical periodic cell containing five grains.
ð1Þ
where xj denotes the position. This condition may with aid of the finite element method be implemented as periodic boundary conditions. However, a requirement for such a formulation is a periodic mesh. Hence, every outer boundary node must have at least one exact equivalent point on the opposite side of the cell. For complicated structures it often becomes difficult to create a periodic mesh. If an equivalent point is missing, the displacement of the exact position can be estimated by the corresponding finite element shape function f. Eq. (1) may then be reformulated as ubi fi ðuai k Þ ¼ eij ðxbj xaj Þ;
ð2Þ
In this study, an adaptive mesh generator has been selected. However, a drawback with this approach is that the mesh can be very coarse in some regions, while it in other regions becomes very fine. The disadvantage with a fine mesh is that is costs computational time, while a coarse mesh instead results in a less accurate calculation. There are several public domain and commercial mesh generators available [15]. Since the Voronoi structure is composed of mixed long and short boundaries, one important task is to generate elements with low element aspect ratios. For this purpose, a software called QMG [16] was chosen. The mesh algorithm is a quadtree/octatree-based algorithm developed by Mitchell and Vavasis [17], for meshing in two or three dimensions. The algorithm theoretically guarantees a bounded aspect ratio on the created mesh.
where fi ðuai k Þ denotes the displacement of the periodicity point. It is here expressed in terms of the finite element shape function ðfi Þ and the corresponding nodal point displacements uai k . The shape function fi assigns coefficients for a linear combination that describes how the three closest nodes are used to calculate the displacement of the exact position. Since the QMG mesh generator does not create a periodic node set, the latter formulation of the periodic boundary conditions must be used. The kinematic constraints are formulated by connecting two equivalent surfaces, i.e. outer boundaries of the periodic cell (as seen in Fig. 1). The node positions of the surface that has the most nodes (dense surface) are constrained by the corresponding nodes of the opposite surface (sparse surface), according to Eq. (2). Thus, all nodes on the dense surface are expressed in terms of node positions of the sparse surface. However, this may result in loose nodes on the sparse surface, those need to be constrained since loose nodes result in a weakening of the model. Thus, any loose nodes of the sparse surface are constrained by the corresponding node positions of the dense surface.
2.3. Periodic boundary conditions
2.4. Finite element simulation
The relation between the displacement vectors uai and ubi at two equivalent boundary points of the
The finite element simulations have been performed by the ABAQUS finite element program
will always be four grain boundaries which meet at each corner point. Due to periodicity, this is also valid for the outer faces. 2.2. Finite element mesh
516
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
3. Results
[18]. The average strain measures: exx , eyy , ezz , cxy , cxz and cyz are represented by degrees of freedom for dummy nodes located outside the model. The coupling between average strains and the nodal displacements is defined in Eq. (2). The reaction forces of the dummy nodes may by use of the principle of virtual work be identified as the corresponding effective stresses rxx , ryy , rzz , sxy , sxz , syz times the volume of the periodic cell. Thereby, a complete stress–strain response is evaluated from nodal displacements and nodal forces. In practice, loading is applied by constraining a combination of average stresses or average strains over the periodic cell. The interpretation of results is simplified since integration over the cell is not needed. Average stresses and strains are simply evaluated from the dummy nodes. Model generation is done in Matlab [19], where the aim is to generate an ABAQUS input file. The model generates 10-noded tetrahedronal elements with four Gauss points. From the finite element simulation, two different results may be evaluated. 0
C out
0:9942 0:4953 B 0:4953 0:9947 B B 0:4954 0:4956 ¼B B 0:0008 0:0008 B @ 0:0006 0:0004 0:0001 0:0000
0:4954 0:4956 0:9948 0:0005 0:0006 0:0001
0:0008 0:0008 0:0005 0:7455 0:0000 0:0003
3.1. Validation of the periodic boundary conditions For validation of the model, four different periodic cells have been generated and meshed. All cells contain five grains and several thousand of elements. For homogeneous material properties, periodic boundary conditions should result in homogeneous stress and strain fields, otherwise the formulation is invalid. In this work, the boundary conditions (Eq. (2)) are tested by an elastic analysis. A stiffness tensor with cubic symmetry is prescribed for each element in the cell, c11 ¼ 1 GPa, c12 ¼ 0:5 GPa and c44 ¼ 0:75 GPa. The cell is then loaded by the six average strain components (one at the time, while all others are zero). The output is used to reconstruct the stiffness tensor. If the formulation is exact, the stiffness tensor should be recovered. The output stiffness tensor C out as generated by the cell shown in Fig. 1 (cell 1) became
0:0006 0:0004 0:0006 0:0000 0:7457 0:0003
First, stresses and strains at a microscale within the periodic cell can be studied by contour plots. Secondly, average stresses and average strains can be evaluated from the dummy nodes, as nodal displacements and reaction forces.
1 0:0001 0:0000 C C 0:0001 C C GPa: 0:0003 C C 0:0003 A 0:7460
It is evident that the periodic boundary condition formulation is an approximation. In Table 1 the number of elements, nodes and the three stiffness out in out in out in component ratios C11 =C11 , C12 =C12 and C44 =C44 for each generated periodic cell are tabulated. The
Table 1 Number of elements, nodes and elastic validation for five unit cells Cell
Elements
Nodes
out in =C11 C11
out in C12 =C12
out in C44 =C44
1 2 3 4
71652 59649 73804 51411
105857 89059 110764 76831
0.9942 0.9941 0.9922 1.0042
0.9907 0.9907 0.9886 0.9910
0.9943 0.9942 0.9921 1.0012
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
517
out out diagonal terms C11 and C44 deviates less than 1% in comparison to the correct solution, while the out error in C12 is about 1%. Hence, it can be expected that the approximate treatment of the periodic boundary conditions in Eq. (2) results in an error of the order of 1%. In comparison to uncertainties in material properties, the accuracy was judged to be sufficient.
3.2. Modeling of a two-phase steel In this study a two-phase ferrite/pearlite steel is used as model material. The material contains 75% ferrite and 25% pearlite, a photograph can be seen in Fig. 3a. Each ferrite grain was assumed to have cubic anisotropy (c11 ¼ 233 GPa, c12 ¼ 124 GPa and c44 ¼ 117 GPa [20]). The direction of each ferrite grain is assigned from a uniform distribution, by prescribing the Euler angles. The pearlite structure, on the other hand, was assumed to be elastic isotropic (E ¼ 202 GPa and m ¼ 0:3). For modeling of plastic properties, it is assumed that the material properties of the single phase materials can be used to describe the constituents of the two-phase materials. Therefore, the yield stress and plastic hardening for the ferrite phase and the pearlite structure were determined from tensile tests, by manufacturing of ferrite and pearlite test specimens. The specimens were fabricated to have properties as similar as possible as in the two-phase structure, with respect to grain size, lamellar spacing etc. In the FE model, the measured hardening of ferrite and pearlite seen in Fig. 2 is used, and it is assumed to be isotropic. More details about the experiments are found in [21]. A major task in modeling two-phase materials is the configuration of the two phases [1]. There are several different possibilities to place the two phases; both or none of the phases can be continuous throughout the periodic cell. As the periodic cell is used to study the two-phase ferrite/ pearlite steels, the grain structure is divided into the two different phases. As a consequence of the adaptive mesh algorithm there are smaller elements close to boundaries. This fact was utilized to distribute the pearlite. With this localization cri-
Fig. 2. Measured plastic properties of ferrite and pearlite.
terion, the pearlite structure will be located along the boundaries and form a fairly irregular phase, which is more or less continuous. The resulting phase configuration is justified from a materials science viewpoint, which can be seen by comparing Fig. 3a and b. Fig. 3a shows a photograph of the material, while Fig. 3b shows the numerical model. In the latter figure, the displayed surface is a cut, through the middle of the cell shown in Fig. 1; the ferrite material is gray while the pearlite material is black. 3.2.1. Experiments Since the pearlite distribution should be physically motivated, the model is used in comparison to experimental data. Uniaxial tensile tests have been simulated by the four generated cells. Even though the three-dimensional models in this case only contain five ferrite grains, they well represent the two-phase materials at higher strains (Fig. 4). It is also apparent that the deviation between different cells is a fairly small, despite the random direction of the individual ferrite grains. The results from the model fall within the experimental scatter at large strains. However, around the yield point none of the models is able to reproduce the experimental data. In this region heterogeneous deformation mechanisms such as L€ uder band formation are of importance, which is not considered here.
518
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
Fig. 3. Pearlite configuration, containing 25% pearlite (black color) and 75% ferrite (gray color), shown in (a) a photograph of the material and (b) the numerical model.
Fig. 4. Comparison of stress–strain curves as generated by the three-dimensional model and experimental results.
4. Discussion and conclusions A three-dimensional model of microstructures at the grain level was presented. By use of a peri-
odic grain structure and periodic boundary conditions, representative results could be obtained with a limited number of grains. It is obvious that the combination of the finite element method and
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
the Voronoi algorithm constitutes an effective tool for microstructural investigations. Since the finite element mesh can be generated in an automatic fashion, the proposed model can easily be utilized to perform statistical analyses for several different materials. There is scope for some improvements of the finite element model. Firstly, the approximation of periodic boundary conditions in Eq. (2) can be avoided if a truly periodic finite element mesh can be generated. This is certainly possible to achieve, but it is not an easy task. Secondly, a more uniform mesh can be achieved if the shortest Voronoi edges are removed or elongated before meshing. This would reduce the total number of elements in the model, and thereby shorten the computational time. Thirdly, if an anisotropic grain structure is desirable, the periodic cell can simply be elongated in one direction before the mesh is created. In the material model, certain improvements may be of interest, especially if the local stress and strain fields within the unit cell are of interest. The crystal structure may be accounted for by consideration of crystal plasticity. Moreover, strain gradient plasticity theories could be used to account for characteristic internal length scales. The model was compared to experimental uniaxial tensile tests of two-phase steels. It was showed that the macroscopic stress–strain curves of ferrite and pearlite could be used to characterize the phases. Acknowledgements This work was supported by the Swedish Foundation for Strategic Research, SSF, through the Brinell Centre. References [1] H.J. B€ ohm, F.G. Rammerstorfer, F.D. Fischer, T. Siegmund, Microscale arrangement effects on the thermomechanical behavior of advanced two-phase materials, J. Eng. Mater. Technol. 116 (1994) 268–273. [2] Z. Xia, W.A. Curtin, P.W.M. Peters, Multiscale modeling of failure in metal matrix composites, Acta Mater. 49 (2001) 273–287. [3] C.I.A. Thomson, M.J. Worswick, A.K. Pilkey, D.J. Lloyd, G. Burger, Modeling void nucleation and growth within
[4]
[5]
[6]
[7]
[8]
[9]
[10]
[11]
[12]
[13]
[14]
[15] [16] [17] [18] [19] [20]
[21]
519
periodic clusters of particles, J. Mech. Phys. Solids 47 (1999) 1–26. E. Nakamachi, C.L. Xie, M. Harimoto, Drawability assessment of bcc steel sheet by using elastic/crystalline viscoplastic finite element analyses, Int. J. Mech. Sci. 43 (2001) 631–652. R. Parisot, S. Forest, A.-F. Gourgues, A. Pineau, D. Mareuse, Modeling the mechanical behavior of a multicrystalline zinc coating on a hot-dip galvanized steel sheet, Comput. Mater. Sci. 19 (2000) 189–204. N. Lippmann, T. Steinkopff, S. Schmauder, P. Gumbsch, 3d-finite-element-modelling of microstructures with the method of multiphase elements, Comput. Mater. Sci. 9 (1997) 28–35. M. Li, S. Ghosh, O. Richmond, H. Weiland, T.N. Rouns. Three dimensional characterization and modeling of particle reinforced metal matrix composites: Part I. Quantitative description of microstructural morphology, Mater. Sci. Eng. A 265 (1999) 153–173. S. Fortune, Voronoi diagrams and delaunay triangulations, in: D.Z. Du, F. Hwang (Eds.), Computing in Euclidean Geometry, World Scientific, Singapore, 1992, pp. 193– 233. K.J. Kozaczek, B.G. Petrovic, C.O. Ruud, S.K. Kurtz, A.R. McIlree, Microstructural modelling of grain-boundary stresses in alloy 600, J. Mater. Sci. 30 (1995) 2390–2400. S. Kumar, S.K. Kurtz, V.K. Agarwala, Micro-stress distribution within polycrystalline aggregate, Acta Mech. 114 (1996) 203–216. M.S. Wu, J. Guo, Analysis of a sector crack in a threedimensional voronoi polycrystal with mictrostructural stresses, J. Appl. Mech. 67 (2000) 50–58. A.P. Roberts, E.J. Garboczi, Elastic moduli of model random three-dimensional closed-cellular solids, Acta Mater. 49 (2001) 189–197. S. Kumar, S.K. Kurtz, J.R. Banavar, M.G. Sharma, Properties of a three-dimensional poisson-voronoi tesselation: A monte carlo study, J. Stat. Phys. 67 (1992) 523– 553. C.D. Barber, D.P. Dobkin, H.T. Huhdanpaa. The Quickhull algorithm for convex hulls, ACM Trans. on Mathematical Software, December 1996. http://www-users.informatik.rwth-aachen.de/roberts/software.html. http://www.cs.cornell.edu/home/vavasis/qmg-home.html. S. Mitchell, S. Vavasis, Quality mesh generation in higher dimensions, SIAM J. Comp. 29 (2000) 1334–1370. Hibbitt, Karlsson & Sorensen Inc., Pawtucket, RI, USA. ABAQUS 5.8, 1998 edition. The MathWorks Inc., Natick, Massachusetts 01760, USA. Matlab 5.1, 1997 edition. K.-H. Hellwege, A.M. Hellwege, (Eds.), Landolt-B€ ornstin, Group III: Elastic, piezoelastic, piezooptic and electrooptic constants of crystal, Crystal and Solid State Physics, Vol. I, 1966, Spring-Verlag, Berlin. M. Nyg ards, P. Gudmundson, Micromechanical modeling of ferritic/pearlitic steels, Mater. Sci. Eng. A 325 (2002) 435–443.
Three-dimensional periodic Voronoi grain models and micromechanical FE-simulations of a two-phase steel Mikael Nyg ards, Peter Gudmundson
*
Department of Solid Mechanics, Royal Institute of Technology (KTH), 100 44 Stockholm, Sweden Received 5 July 2001; received in revised form 18 December 2001; accepted 18 December 2001
Abstract A three-dimensional model is proposed for modeling of microstructures. The model is based on the finite element method with periodic boundary conditions. The Voronoi algorithm is used to generate the geometrical model, which has a periodic grain structure that follows the original boundaries of the Voronoi cells. As an application, the model is used to model a two-phase ferrite/pearlite steel. It is shown that periodic cells with only five grains generate representative stress–strain curves. Ó 2002 Elsevier Science B.V. All rights reserved. Keywords: FEM; Microstructure; Voronoi cells; Plasticity; Pearlite; Ferrite
1. Introduction It is desirable to design microstructures to get special properties, such as high strength, ductility or toughness. To do this design properly, the role of microstructure must be thoroughly understood, since most macroscopic behaviors can be related to representative parts of the microstructure. Micromechanical modeling gives the possibility to investigate certain microstructural properties in more detail, and relate those to the macroscopic behavior of the material. The main advantage with modeling is that the influence of various parameters easily can be investigated. A corresponding
*
Corresponding author. Tel.: +46-8-790-7548; fax: +46-8411-2418. E-mail address: [email protected] (P. Gudmundson).
experimental study would be both expensive and time consuming. The morphology of a three-dimensional microstructure decides to a large extent the modeling strategy. If the microstructure is elongated in one direction, a three-dimensional micromechanical model that is similar to the real material is fairly straightforward to generate [1,2]. When the microstructure is random in three dimensions, the problem can be approached by a simplification that still captures the features of interest in the material, like distributing particles to investigate clustering [3]. In polycrystalline materials, the characteristic feature is the grain structure. Grains may be simplified geometrically to capture certain features, for example as cubes [4] or as elongated twodimensional structures [5]. Alternatively, the microstructure can be mapped, and a corresponding
0927-0256/02/$ - see front matter Ó 2002 Elsevier Science B.V. All rights reserved. PII: S 0 9 2 7 - 0 2 5 6 ( 0 2 ) 0 0 1 5 6 - 8
514
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
microstructure can be simulated [6,7]. However, to model more realistic microstructures, it is desirable to use an algorithm to create the grain structure, and to ensure its randomness. For this purpose the Voronoi algorithm has been shown to be successful. The topology of the Voronoi tessellation, see for example [8], is similar to the microstructure of real materials. Therefore, it is advantageous to use for modeling of materials at a microstructural scale. The three-dimensional Voronoi algorithm has been applied to materials science by other researchers in the past. Local stresses within the generated Voronoi structure have been reported by e.g. Kozaczek et al. [9], Kumar et al. [10], Wu and Guo [11]. While Roberts and Garboczi [12] used the Voronoi structure to study the elastic moduli of a cellular solid. Statistical properties of the Voronoi structure have been investigated by Kumar et al. [13]. The three-dimensional modeling by Voronoi cells offers the unique possibility to automatically generate unit cells and it can therefore capture statistical data. Hence, effects of morphology may be studied. In the present model, a concept for modeling microstructures is presented in Section 2. A periodic unit cell that follows the original boundaries of a generated Voronoi structure is loaded by periodic boundary conditions. Thereby, grains are not divided to fit within a square unit cell, which is advantageous when the grain topology is of importance, also edge effects are avoided. In Section 3, a validation with an orthotropic elastic material model is performed. As an application, a two-phase ferrite/pearlite steel is modeled with an elastic–plastic analysis. In the model, an initial ferrite phase is created by the Voronoi algorithm, thereafter, the properties of the adaptive mesh is used to distribute the pearlite. Hence, a two-phase irregular microstructure is created that is similar to the real microstructure.
2. Numerical model On a microstructural level, polycrystalline materials consist of a more or less random grain structure. In this study the so-called Voronoi al-
gorithm is used to model microstructures in three dimensions. 2.1. Periodic Voronoi structure The idea behind a micromechanical model based on the Voronoi algorithm is to model a representative part of the microstructure and to predict the macroscopic behavior of the material. The cell should be large enough to catch the macroscopic behavior, and it should be kept as small as possible to save computation time. To avoid undesired edge effects, a periodic cell loaded by periodic boundary conditions has been chosen as the representative unit cell. The Voronoi algorithm utilizes a set of sites to partition space into regions; one per site. The region for a site s consists of all points closer to s than to any other site. In practice, the distances between the site s and all its closest neighbors are divided in half by planes perpendicular to a line connecting the sites. In this manner, a three-dimensional structure resembling a real microstructure is created. In this work, a public domain software called Qhull [14] has been used to create the initial Voronoi structure. The program is given a random set of sites as input, and it outputs the edge points, surfaces and regions of the Voronoi cells. A periodic Voronoi structure results if the set of sites show a periodic symmetry. This property has been utilized in the present study to generate periodic Voronoi structures. Firstly, a preset number of sites (grains within a periodic cell) is generated within a unit cube. The location of each site is set by a random generator with uniform distribution in the cube. The site locations within the cube are then copied to horizontally, vertically and diagonally located cubes. Lastly, the Voronoi algorithm is applied to the sites within 3 3 3 ¼ 27 cubes. A unique periodic cell can be chosen from the center of the created Voronoi structure. In Fig. 1, a periodic cell containing five grains is shown. It is evident that each grain is represented by irregular convex polyhedrons. By considering the nature of the Voronoi algorithm, there are two regions on each side of every surface, every edge has three regions around it, and each edge corner has four regions sharing that corner. In other words, there
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
515
periodic cell and the average strain tensor eij can be expressed as ubi uai ¼ eij ðxbj xaj Þ;
Fig. 1. A typical periodic cell containing five grains.
ð1Þ
where xj denotes the position. This condition may with aid of the finite element method be implemented as periodic boundary conditions. However, a requirement for such a formulation is a periodic mesh. Hence, every outer boundary node must have at least one exact equivalent point on the opposite side of the cell. For complicated structures it often becomes difficult to create a periodic mesh. If an equivalent point is missing, the displacement of the exact position can be estimated by the corresponding finite element shape function f. Eq. (1) may then be reformulated as ubi fi ðuai k Þ ¼ eij ðxbj xaj Þ;
ð2Þ
In this study, an adaptive mesh generator has been selected. However, a drawback with this approach is that the mesh can be very coarse in some regions, while it in other regions becomes very fine. The disadvantage with a fine mesh is that is costs computational time, while a coarse mesh instead results in a less accurate calculation. There are several public domain and commercial mesh generators available [15]. Since the Voronoi structure is composed of mixed long and short boundaries, one important task is to generate elements with low element aspect ratios. For this purpose, a software called QMG [16] was chosen. The mesh algorithm is a quadtree/octatree-based algorithm developed by Mitchell and Vavasis [17], for meshing in two or three dimensions. The algorithm theoretically guarantees a bounded aspect ratio on the created mesh.
where fi ðuai k Þ denotes the displacement of the periodicity point. It is here expressed in terms of the finite element shape function ðfi Þ and the corresponding nodal point displacements uai k . The shape function fi assigns coefficients for a linear combination that describes how the three closest nodes are used to calculate the displacement of the exact position. Since the QMG mesh generator does not create a periodic node set, the latter formulation of the periodic boundary conditions must be used. The kinematic constraints are formulated by connecting two equivalent surfaces, i.e. outer boundaries of the periodic cell (as seen in Fig. 1). The node positions of the surface that has the most nodes (dense surface) are constrained by the corresponding nodes of the opposite surface (sparse surface), according to Eq. (2). Thus, all nodes on the dense surface are expressed in terms of node positions of the sparse surface. However, this may result in loose nodes on the sparse surface, those need to be constrained since loose nodes result in a weakening of the model. Thus, any loose nodes of the sparse surface are constrained by the corresponding node positions of the dense surface.
2.3. Periodic boundary conditions
2.4. Finite element simulation
The relation between the displacement vectors uai and ubi at two equivalent boundary points of the
The finite element simulations have been performed by the ABAQUS finite element program
will always be four grain boundaries which meet at each corner point. Due to periodicity, this is also valid for the outer faces. 2.2. Finite element mesh
516
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
3. Results
[18]. The average strain measures: exx , eyy , ezz , cxy , cxz and cyz are represented by degrees of freedom for dummy nodes located outside the model. The coupling between average strains and the nodal displacements is defined in Eq. (2). The reaction forces of the dummy nodes may by use of the principle of virtual work be identified as the corresponding effective stresses rxx , ryy , rzz , sxy , sxz , syz times the volume of the periodic cell. Thereby, a complete stress–strain response is evaluated from nodal displacements and nodal forces. In practice, loading is applied by constraining a combination of average stresses or average strains over the periodic cell. The interpretation of results is simplified since integration over the cell is not needed. Average stresses and strains are simply evaluated from the dummy nodes. Model generation is done in Matlab [19], where the aim is to generate an ABAQUS input file. The model generates 10-noded tetrahedronal elements with four Gauss points. From the finite element simulation, two different results may be evaluated. 0
C out
0:9942 0:4953 B 0:4953 0:9947 B B 0:4954 0:4956 ¼B B 0:0008 0:0008 B @ 0:0006 0:0004 0:0001 0:0000
0:4954 0:4956 0:9948 0:0005 0:0006 0:0001
0:0008 0:0008 0:0005 0:7455 0:0000 0:0003
3.1. Validation of the periodic boundary conditions For validation of the model, four different periodic cells have been generated and meshed. All cells contain five grains and several thousand of elements. For homogeneous material properties, periodic boundary conditions should result in homogeneous stress and strain fields, otherwise the formulation is invalid. In this work, the boundary conditions (Eq. (2)) are tested by an elastic analysis. A stiffness tensor with cubic symmetry is prescribed for each element in the cell, c11 ¼ 1 GPa, c12 ¼ 0:5 GPa and c44 ¼ 0:75 GPa. The cell is then loaded by the six average strain components (one at the time, while all others are zero). The output is used to reconstruct the stiffness tensor. If the formulation is exact, the stiffness tensor should be recovered. The output stiffness tensor C out as generated by the cell shown in Fig. 1 (cell 1) became
0:0006 0:0004 0:0006 0:0000 0:7457 0:0003
First, stresses and strains at a microscale within the periodic cell can be studied by contour plots. Secondly, average stresses and average strains can be evaluated from the dummy nodes, as nodal displacements and reaction forces.
1 0:0001 0:0000 C C 0:0001 C C GPa: 0:0003 C C 0:0003 A 0:7460
It is evident that the periodic boundary condition formulation is an approximation. In Table 1 the number of elements, nodes and the three stiffness out in out in out in component ratios C11 =C11 , C12 =C12 and C44 =C44 for each generated periodic cell are tabulated. The
Table 1 Number of elements, nodes and elastic validation for five unit cells Cell
Elements
Nodes
out in =C11 C11
out in C12 =C12
out in C44 =C44
1 2 3 4
71652 59649 73804 51411
105857 89059 110764 76831
0.9942 0.9941 0.9922 1.0042
0.9907 0.9907 0.9886 0.9910
0.9943 0.9942 0.9921 1.0012
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
517
out out diagonal terms C11 and C44 deviates less than 1% in comparison to the correct solution, while the out error in C12 is about 1%. Hence, it can be expected that the approximate treatment of the periodic boundary conditions in Eq. (2) results in an error of the order of 1%. In comparison to uncertainties in material properties, the accuracy was judged to be sufficient.
3.2. Modeling of a two-phase steel In this study a two-phase ferrite/pearlite steel is used as model material. The material contains 75% ferrite and 25% pearlite, a photograph can be seen in Fig. 3a. Each ferrite grain was assumed to have cubic anisotropy (c11 ¼ 233 GPa, c12 ¼ 124 GPa and c44 ¼ 117 GPa [20]). The direction of each ferrite grain is assigned from a uniform distribution, by prescribing the Euler angles. The pearlite structure, on the other hand, was assumed to be elastic isotropic (E ¼ 202 GPa and m ¼ 0:3). For modeling of plastic properties, it is assumed that the material properties of the single phase materials can be used to describe the constituents of the two-phase materials. Therefore, the yield stress and plastic hardening for the ferrite phase and the pearlite structure were determined from tensile tests, by manufacturing of ferrite and pearlite test specimens. The specimens were fabricated to have properties as similar as possible as in the two-phase structure, with respect to grain size, lamellar spacing etc. In the FE model, the measured hardening of ferrite and pearlite seen in Fig. 2 is used, and it is assumed to be isotropic. More details about the experiments are found in [21]. A major task in modeling two-phase materials is the configuration of the two phases [1]. There are several different possibilities to place the two phases; both or none of the phases can be continuous throughout the periodic cell. As the periodic cell is used to study the two-phase ferrite/ pearlite steels, the grain structure is divided into the two different phases. As a consequence of the adaptive mesh algorithm there are smaller elements close to boundaries. This fact was utilized to distribute the pearlite. With this localization cri-
Fig. 2. Measured plastic properties of ferrite and pearlite.
terion, the pearlite structure will be located along the boundaries and form a fairly irregular phase, which is more or less continuous. The resulting phase configuration is justified from a materials science viewpoint, which can be seen by comparing Fig. 3a and b. Fig. 3a shows a photograph of the material, while Fig. 3b shows the numerical model. In the latter figure, the displayed surface is a cut, through the middle of the cell shown in Fig. 1; the ferrite material is gray while the pearlite material is black. 3.2.1. Experiments Since the pearlite distribution should be physically motivated, the model is used in comparison to experimental data. Uniaxial tensile tests have been simulated by the four generated cells. Even though the three-dimensional models in this case only contain five ferrite grains, they well represent the two-phase materials at higher strains (Fig. 4). It is also apparent that the deviation between different cells is a fairly small, despite the random direction of the individual ferrite grains. The results from the model fall within the experimental scatter at large strains. However, around the yield point none of the models is able to reproduce the experimental data. In this region heterogeneous deformation mechanisms such as L€ uder band formation are of importance, which is not considered here.
518
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
Fig. 3. Pearlite configuration, containing 25% pearlite (black color) and 75% ferrite (gray color), shown in (a) a photograph of the material and (b) the numerical model.
Fig. 4. Comparison of stress–strain curves as generated by the three-dimensional model and experimental results.
4. Discussion and conclusions A three-dimensional model of microstructures at the grain level was presented. By use of a peri-
odic grain structure and periodic boundary conditions, representative results could be obtained with a limited number of grains. It is obvious that the combination of the finite element method and
M. Nyg ards, P. Gudmundson / Computational Materials Science 24 (2002) 513–519
the Voronoi algorithm constitutes an effective tool for microstructural investigations. Since the finite element mesh can be generated in an automatic fashion, the proposed model can easily be utilized to perform statistical analyses for several different materials. There is scope for some improvements of the finite element model. Firstly, the approximation of periodic boundary conditions in Eq. (2) can be avoided if a truly periodic finite element mesh can be generated. This is certainly possible to achieve, but it is not an easy task. Secondly, a more uniform mesh can be achieved if the shortest Voronoi edges are removed or elongated before meshing. This would reduce the total number of elements in the model, and thereby shorten the computational time. Thirdly, if an anisotropic grain structure is desirable, the periodic cell can simply be elongated in one direction before the mesh is created. In the material model, certain improvements may be of interest, especially if the local stress and strain fields within the unit cell are of interest. The crystal structure may be accounted for by consideration of crystal plasticity. Moreover, strain gradient plasticity theories could be used to account for characteristic internal length scales. The model was compared to experimental uniaxial tensile tests of two-phase steels. It was showed that the macroscopic stress–strain curves of ferrite and pearlite could be used to characterize the phases. Acknowledgements This work was supported by the Swedish Foundation for Strategic Research, SSF, through the Brinell Centre. References [1] H.J. B€ ohm, F.G. Rammerstorfer, F.D. Fischer, T. Siegmund, Microscale arrangement effects on the thermomechanical behavior of advanced two-phase materials, J. Eng. Mater. Technol. 116 (1994) 268–273. [2] Z. Xia, W.A. Curtin, P.W.M. Peters, Multiscale modeling of failure in metal matrix composites, Acta Mater. 49 (2001) 273–287. [3] C.I.A. Thomson, M.J. Worswick, A.K. Pilkey, D.J. Lloyd, G. Burger, Modeling void nucleation and growth within
[4]
[5]
[6]
[7]
[8]
[9]
[10]
[11]
[12]
[13]
[14]
[15] [16] [17] [18] [19] [20]
[21]
519
periodic clusters of particles, J. Mech. Phys. Solids 47 (1999) 1–26. E. Nakamachi, C.L. Xie, M. Harimoto, Drawability assessment of bcc steel sheet by using elastic/crystalline viscoplastic finite element analyses, Int. J. Mech. Sci. 43 (2001) 631–652. R. Parisot, S. Forest, A.-F. Gourgues, A. Pineau, D. Mareuse, Modeling the mechanical behavior of a multicrystalline zinc coating on a hot-dip galvanized steel sheet, Comput. Mater. Sci. 19 (2000) 189–204. N. Lippmann, T. Steinkopff, S. Schmauder, P. Gumbsch, 3d-finite-element-modelling of microstructures with the method of multiphase elements, Comput. Mater. Sci. 9 (1997) 28–35. M. Li, S. Ghosh, O. Richmond, H. Weiland, T.N. Rouns. Three dimensional characterization and modeling of particle reinforced metal matrix composites: Part I. Quantitative description of microstructural morphology, Mater. Sci. Eng. A 265 (1999) 153–173. S. Fortune, Voronoi diagrams and delaunay triangulations, in: D.Z. Du, F. Hwang (Eds.), Computing in Euclidean Geometry, World Scientific, Singapore, 1992, pp. 193– 233. K.J. Kozaczek, B.G. Petrovic, C.O. Ruud, S.K. Kurtz, A.R. McIlree, Microstructural modelling of grain-boundary stresses in alloy 600, J. Mater. Sci. 30 (1995) 2390–2400. S. Kumar, S.K. Kurtz, V.K. Agarwala, Micro-stress distribution within polycrystalline aggregate, Acta Mech. 114 (1996) 203–216. M.S. Wu, J. Guo, Analysis of a sector crack in a threedimensional voronoi polycrystal with mictrostructural stresses, J. Appl. Mech. 67 (2000) 50–58. A.P. Roberts, E.J. Garboczi, Elastic moduli of model random three-dimensional closed-cellular solids, Acta Mater. 49 (2001) 189–197. S. Kumar, S.K. Kurtz, J.R. Banavar, M.G. Sharma, Properties of a three-dimensional poisson-voronoi tesselation: A monte carlo study, J. Stat. Phys. 67 (1992) 523– 553. C.D. Barber, D.P. Dobkin, H.T. Huhdanpaa. The Quickhull algorithm for convex hulls, ACM Trans. on Mathematical Software, December 1996. http://www-users.informatik.rwth-aachen.de/roberts/software.html. http://www.cs.cornell.edu/home/vavasis/qmg-home.html. S. Mitchell, S. Vavasis, Quality mesh generation in higher dimensions, SIAM J. Comp. 29 (2000) 1334–1370. Hibbitt, Karlsson & Sorensen Inc., Pawtucket, RI, USA. ABAQUS 5.8, 1998 edition. The MathWorks Inc., Natick, Massachusetts 01760, USA. Matlab 5.1, 1997 edition. K.-H. Hellwege, A.M. Hellwege, (Eds.), Landolt-B€ ornstin, Group III: Elastic, piezoelastic, piezooptic and electrooptic constants of crystal, Crystal and Solid State Physics, Vol. I, 1966, Spring-Verlag, Berlin. M. Nyg ards, P. Gudmundson, Micromechanical modeling of ferritic/pearlitic steels, Mater. Sci. Eng. A 325 (2002) 435–443.