Void growth and morphology evolution in fcc bi-crystals are investigated using crystal plasticity nite element method. For that purpose, representative volume element of bi-crystals with a void at the grain boundary are considered in the analysis. Grain boundary is assumed initially perpendicular/coaxial with the straight sides of the cell. Fully periodic boundary conditions are prescribed in the representative volume element and macroscopic stress triaxiality and Lode parameter are kept constant during the whole deformation process. Three di erent pairs of crystal orientations characterized as hard-hard, soft-soft and soft-hard has been employed for modelling the mechanical response of the bi-crystal. Simulations are performed to study the implications of triaxiality, Lode parameter and crystallographic orientation on slip mechanism, hardening and hence void evolution. The impact of void presence and its growth on the heterogeneity of lattice rotation and resulting grain fragmentation in neighbouring areas is also analysed and discussed. Keywords: Crystal plasticity; Bi-crystals; Void growth; Stress Triaxiality; Lode parameter; Unit cell calculations. 1. Introduction One of the common challenges to overcome when using advanced metals and alloys in en- gineering applications is their insucient ductility (Fourmeau et al., 2013; Basu et al., 2017). Strengthening of materials is usually achieved by introduction of precipitates, second phase par- ticles or new grain boundaries while at the same time these additions are possible sites of dam- age/fracture initiation, cf. (Yerra et al., 2010). Predominant mechanisms of ductile failure in polycrystaline metals are nucleation, growth and coalescence of micro-voids (small scale voids), cf. (Benzerga and Besson, 2001). Accounting for these mechanisms in modelling, a family of macroscopic approximations were proposed, namely prominent models like Gurson-types models (Gurson, 1977; Tvergaard and Needleman, 1984), their extension in later works by taking into Corresponding author. Email addresses: mdakshin@ing.uc3m.es (Manjunath Dakshinamurthy), kkowalcz@ippt.pan.pl (Katarzyna Kowalczyk-Gajewska), gvadillo@ing.uc3m.es (Guadalupe Vadillo) Preprint submitted to Elsevier October 26, 2020 account many aspects such as void shape, size, volume fraction and distribution of voids, as well as distinct features of the constitutive model of a virgin material such as strain and kinematic hardening, viscoplasticity and plastic anisotropy. For respective contributions see extensive review by Besson (2010). Such type of approach was also used to formulate the micro-scale Gurson-type yield condition for single crystal (Han et al., 2013; Paux et al., 2015). In the case of Han et al. (2013) such condition is the result of micromechanical analysis based on variational estimates due to de Botton and Ponte Casta~ neda (1995). Note that, there are number of papers in which such variational approach is directly applied to estimate the yield surface of single crystal with voids (Mbiakop et al., 2015; Song and Casta~ neda, 2017). On the micro-scale, two failure modes have been found to be operational in parallel: cleavage and dimple fracture (Papaefthymiou et al., 2006). Studies of Kadkhodapour et al. (2011a,b); Lani et al. (2007); Furn emont et al. (2007) have reported di erent types of damage mechanisms based on their experimental work with the help of SEM and light optical microscopy. By experimental observations, nucleation of voids is related to the fracture and decohesion of second phase parti- cles or other precipitates. The next steps in understanding the failure process are intergranular fracture along high angle boundary (HAB) and microvoid-induced transgranular fracture. For example, such observations have been made for aluminium alloys AA7XXX subjected to the so- lution treatment along the processing route, which resulted in partial recrystallization (Dorward and Beerntsen, 1995; Deshpande et al., 1998). For such alloys it was observed that precipitates are often grouped along HABs which promotes intergranular fracture along them. Importance of HAB as a possible location of void growth has been also con rmed by molecular dynamics simulations (Bringa et al., 2010). Additionally, as seen in Morere et al. (2000), concentration of precipitates is also observed along subgrain boundaries, characterized by a low misorientation angle, together with associated intersubgranular failure. It is also worth to note that the e ect of crystal orientation and the grain boundary was observed under shock loading conditions for copper bi-crystal by Perez-Bergquist et al. (2011). They observed formation of voids both at the grain boundary and within each crystal, with the number and shape of voids highly dependent on the crystallographic orientation. Although this type of loading is beyond the scope of this paper, these experiments con rm the role of crystal anisotropy and grain boundaries in the failure of crystalline material. Although there is quite a lot numerical studies of single crystals with voids, (O'Regan et al., 1997; Potirniche et al., 2006; Yerra et al., 2010; Srivastava and Needleman, 2015a; Selvarajou et al., 2019) very few numerical studies have been dedicated to the studies on void in bi-crystals. One can mention earlier works by Liu et al. (2009) in which the fcc bi-crystal under uniaxial tension was investigated to assess the e ect of boundary inclination and the crystallographic orientation and misorientation on the void growth and the plastic deformation distribution around the voids. Liu et al. (2007, 2010) studied mechanisms of coalescence of voids located at adjacent grains in bi-crystal unit cell using nite element calculations. Recently, Jeong et al. (2018) considered the unit cells with a void inside a grain, at a grain boundary and at a triple junction. The e ect of crystal orientation on the ow strength and growth rate of the void was discussed under prescribed boundary conditions for constant stress triaxialities. Even though studies mentioned above shed light on the behaviour of void at boundary on bi-crystals, there are still open questions, like the e ect of neighbouring grain orientation on void shape, the e ect of stress state on void growth or void collapse behaviour among others. 2 Appreciating the important role of a grain boundary on void behaviour, as noted by Deshpande et al. (1998); Bringa et al. (2010); Morere et al. (2000), the main purpose of the current work is to study the factors a ecting void growth and morphology evolution in bi-crystals, such as stress triaxiality and Lode parameter and the relative orientation between grains, focusing speci cally on voids originating in high angle grain boundaries. The response of bi-crystals containing a void at the grain boundary are also compared with the response of a single crystal with a void. The paper is constructed as follows, after this introductory section, the problem formulation including crystal plasticity framework, representative volume element and prescribed boundary conditions are presented in Section 2. In Section 3 results of FE simulations are compared from single crystals with bi-crystals and from di erent bi-crystals RVEs. The paper is closed by con- clusions. Details of implementation of user de ned multi-point constraints (MPC) subroutine for controlling stress state (stress triaxiality and Lode parameter) and periodic boundary conditions (PBC) are included in the Appendix. Table 1: Slip systemsfem ;en gof the FCC crystal structure en em en em 1(111) [101] 7(111) [011] 2(111) [011] 8(111) [101] 3(111) (110] 9(111) [110] 4(111) [101] 10(111)[110] 5(111) [011] 11(111)[101] 6(111) [110] 12(111)[011] 2. Problem formulation Three dimensional nite element calculations are carried out to model the response of voids, located at the grain boundary of a bi-crystal, under triaxial loading conditions using a unit cell model. The unit cell is modelled as an fcc single crystal with the 12 potentially active slip systems taken to be f111gh110i(Table 1). Finite element analyses are carried out using a rate dependent crystal plasticity constitutive relation in the large deformation framework. Description of the crystal plasticity model and material parameters used for calculations are given below. The plastic parameters of the material analyzed in this work correspond to annealed OFHC copper (Kalidindi et al., 1992). 2.1. Constitutive model for a single crystal Description of kinematics follows the standard large strain deformation framework. The mul- tiplicative decomposition is applied, according to which the deformation gradient Fis the product of an elastic part Feand a plastic part Fp, namely: F=FeFp: (1) 3 Table 2: Material parameters for crystal plasticity model (Han et al., 2013; Kalidindi et al., 1992) Parameter Notation Value Reference shear rate _ 0(s1) 0.001 Inverse of strain rate sensitivity p1 20 Initial critical resolved shear stress ^ 0(MPa) 20 Initial strain hardening modulus h0(MPa) 180 Saturated critical resolved shear stress ^ s(MPa) 117 Strain hardening exponent p2 2.25 Self-hardening coecient q ( = ) 1 Latent hardening coecient q ( 6= ) 1.4 Elastic constants c11(GPa) 199 c12(GPa) 136 c44(GPa) 105 It is assumed that the plastic deformation occurs by dislocation glide on slip systems, therefore the evolution of the plastic part follows the equation _Fp=0 @NSX =1_ fM 1 AFp;fM =em en (2) fM being a Schmid tensor where the slip directions em and planesen are speci ed in the interme- diate con guration. It should be noted that the plastic deformation described by Eq.2 is volume preserving, det Fp= 1, sinceem ·en = 0. The current orientation of crystallographic directions and planes can be found using the elastic part of the deformation gradient as follows n(t) = (FeT)1en ;m(t) =Feem : (3) Because the face centred cubic (fcc) crystal structure is considered in the study, the dislocations glide happens on the f111gplanes inh110idirections. For the NS= 12 slip systems, the slip planes and directions are listed in Table 1. The shear rate _ on a given slip system is related the stress state by the phenomenological visco-plastic type power law (Asaro and Needleman, 1985) _ = _ 0 ^ p1 sign( ): (4) with the resolved shear stress =eT ·fM (5) where the parameter _ 0is the reference shear rate, p1is the inverse of the strain rate sensitivity parameter and eTis the Mandel stress. When the elastic strains are small as compared to the 4 inelastic ones, the Mandel stress can be approximated by the IInd Piola-Kirchho stress pushed forward to the intermediate con guration eS, namely: eT= (FeTFe)|{z} IeSFpSFpT(6) The stress tensor eSis calculated with help of the hyperelastic law assuming the Kirchho -type function of free energy density per unit volume in the reference con guration eS=@ @Ee=C ·Eewhere Ee=1 2 FeTFeI (7) withCthe elastic IVth order sti ness tensor. The critical resolved shear stress ^ , obeys the following evolution rule (Kalidindi et al., 1992) _^ =h0NSX =1q  1^ ^ sp2 j_ j (8) whereh0is the initial hardening rate, ^ sthe saturation resistance, p2the hardening exponent and q the cross hardening coecient matrix. The model parameters for elasticity and plasticity are given in Table 2. The constitutive equations mentioned above are implemented into a user subroutine UMAT in ABAQUS/Standard (2019). The implementation of the crystal plasticity model follows the fully implicit time integration procedure with the consistent material tangent operator presented by Meissonnier et al. (2001). 2.2. Representative Volume Element and boundary conditions In this section, representative unit cells used for studying the behaviour of a void in both single crystal and bi-crystal cases are described brie y. For studying the void response in a single crystal, the representative volume element consisting of a cubic cell with a initial spherical void at the center as shown in gure 1 is considered and initial void volume fraction f0= 0:0044 is used. The initial void volume fraction is selected to be small, so that we can study the void shape changes and its evolutions at large deformations. Void volume fraction is calculated as f0=Vvoid=Vcell, whereVcellis the initial volume of the cell and Vvoidis the initial volume of the void. The complete matrix-void zone is discretized into 90000 C3D8R elements. For the bi-crystal case, representative unit cell consists in a cube with a void inside divided into two equal parts, each half having half spherical void with f0= 0:0022, as shown in gure 4. Finite element mesh and number of elements in the mesh are the same as for the single crystal case. Depending on the cases to be studied, di erent crystallographic orientations are assigned to each half of the bi-crystal matrix and will be discussed in detail in the next section. Initially, unit cells are aligned such that the edges of the unit cell are parallel to the global XYZ axes. In the case of the bi-crystal cell, interface (grain boundary) between the two halves of the unit cell is perpendicular to the global X axis. The notation adopted here denotes each orientation by the crystallographic directions, aligned with the main loading direction. The main loading direction is parallel to the Xaxis and the secondary loading direction is parallel to Z 5 orYaxis depending on the selected Lode parameter. Six crystallographic orientations are used in this study. Orientations are given in terms of the Bunge Euler angles and crystallographic planes parallel to each loading direction. Pole gures for the orientations considered are shown in gure 1. For simplicity, each 6 orientations will be henceforth called O1 - O6 as shown in Table 3. Macroscopic Cauchy stress ijcomponents are related to the microscopic Cauchy stresses  ij through the relations: ij=1 VcellR VcellijdVcell (9) and the de nition of macroscopic equivalent stress eqv, triaxiality Tand Lode parameter Lare given by: eqv=q 3 20ij0ij;T=h eqv;L=2213 13(10) with 0ij=ijhij;h=ii 3(11) and1,2and3the principal macroscopic true stress components. Following Srivastava and Needleman (2015a), an equivalent strain is calculated based on the following numerical volume average of local logarithmic strain  "ij, namely: h"iji=1 VmatrixZ Vmatrix"ijdVmatrix (12) whereVmatrix denotes the volume of bi-crystal excluding the void. It must be stressed that, in general, this volume averaged value h"ijiis not equal to the macroscopic logarithmic strain as long as deformation is non-uniform. The de nition of the equivalent strain is then proposed as follows: "eqv=r 2 3h"0ijih"0ijiwhereh"0 iji=h"iji"hij;"h=h"iii 3(13) This value is used only for the purpose of graphical presentation of results. Proper boundary conditions are important for studying symmetric (O1, O2, O3) and non symmetric (O4, O5, O6) crystallographic orientations, as crystals tend to come shear strain com- ponents. If shear strain components are constrained, cells deformation behaviour is a ected and hence void growth. In order to avoid this e ect in our RVEs, fully periodic boundary conditions are applied to the cell in all three directions, thereby allowing crystal to deform freely without constraints. The main idea of periodic boundary conditions is that the relative displacement be- tween opposing nodes is the same, i.e., the nodal displacement of opposite nodes is coupled, but still allowing global strain in this direction. Note that by prescribing periodic boundary condi- tions we are analysing a unit cell which corresponds to the con guration representing a bi-crystal laminate with voids periodically distributed along every second grain boundary. In addition, the values ofTandLare prescribed throughout the numerical simulations as constant values. This is achieved by implementing a user de ned multi point constraints subroutine in the framework of ABAQUS/Standard (2019). Detailed description of periodic boundary conditions and multi point constraints implementations in ABAQUS/Standard (2019) are given in Appendix A. 6 In this work four values of stress triaxiality ( T= 0, 1=3, 2=3, 1) and three values of Lode parameter ( L= 0, 1,1) with fully periodic boundary conditions are considered. Some of the combinations of triaxiality and Lode parameters analyzed corresponds to general loading scenarios used in experiments and simple loading conditions, like T= 1=3,L=1 (uniaxial tension), T= 2=3,L= 1 (biaxial tension) and T= 0,L= 0 (pure shear stress). Since anisotropic material model is assumed in each case orientation of principal stress axes is explicitly speci ed with respect to material anisotropy axes. It should be underlined that in the considered examples dependence of the void growth on the Lode parameter cannot be studied separately from the e ect of plastic anisotropy. Finally, as our material model is rate dependent, a value of _ "eqv=_ 0in the range (100150) is also assured in all the simulations performed in this work. Table 3: Orientation of single crystal matrix for di erent RVEs in terms of Euler angles. The Euler angles are de ned with respect to global axis (X,Y,Z) RVE Bunge Euler angles ()Crystallographic orientation on global ( 1;; 2) coordinates X-Z O1 (0,0,0) [100] -[001] O2 (0,45,0) [100] -[011] O3 (45,0,0) [110] -[001] O4 (45,54.73,0) [110] -[111] O5 (0,45,54.73) [111] -[1 10] O6 (45,35.26,90) [111] -[ 112] 2.3. Single crystal without voids Before studying the void response in single crystal and bi-crystal, we will rst study the response of single crystal without void for di erent loading conditions (Lode parameters L= 0;1;1) and the orientations given in Table 3 . Equivalent stress-strain curves are plotted for all the Lode parameters and orientations considered (see gures 2). De nitions of equivalent stress and strains are as given by equations 10 and 13. Note that yielding of a single crystal without void is not a ected by triaxiality parameter T since resolved shear stress in Eq. 5 depends only on the deviatoric part of stress. This leads to the conclusion that, in particular, an equivalent strain-stress curve for L=1 (resp.L= 1) and a selected orientation is the same for uniaxial tension (resp. uniaxial compression) and biaxial compression (resp. biaxial tension). Additionally, due to the form of Eq. 4 and 5, the same response in terms of equivalent stress and strain is obtained when we only change the sign of stress. In particular, the same response is found for uniaxial tension ( L=1 andT= 1=3) and uniaxial compression ( L= 1 andT=1=3) along the same crystallographic direction. However, in view of the convention assumed in our calculations for direction of primary loading, for L=1 (andT= 1=3) we have to do with uniaxial tension along crystallographic direction parallel to X axis, while for L= 1 (andT=1=3) with uniaxial compression along crystallographic direction parallel to Z axis. Therefore, if the orientation is not symmetric with respect to this change, like e.g. O4, O5 or O6, we may observe di erent responses for L= 1 andL=1 in Fig. 2. It can be said that in this case the e ect of Lode parameter is only apparent, and in fact the 7 Figure 1: Cell containing a initially spherical void at its center and pole gures representing the initial orientations studied in this work. observed di erence is due to material anisotropy. Similar e ects are observed when the quadratic Hill criterion is used. Such e ects were thoroughly discussed by Benzerga and Besson (2001). Nevertheless in the present case also a direct e ect of the Lode parameter is observed due to non- linearity of the model. The best example to demonstrate such e ect is orientation O1 for which three principal directions of stress are equivalent from the point of view of crystal symmetry. In this case as expected the same stress-strain curve is found for L= 1 andL=1, however, the response is di erent for L= 0. It is stressed that this di erence is found only when p1>1 in Eq. (4) and increases with increased non-linearity. For Lode parameter L=1, the value of equivalent stress is higher for O5 and O6 followed by O4 and least for O1, O2, O3, as shown in gure 2(a). Similarly for Lode parameter L= 1, the highest value of equivalent stress is for O4 followed by O5, O6 and least for O1, O2, O3 as shown in gure 2(b). Finally, as shown in gure 2(c) for L= 0, the highest value of equivalent stress is for O4, O5, O6 orientations and least for O1, O2, O3. 8 From this data, we categorized the crystal orientations into hard and soft orientations. For a given Lode parameter, the orientation with higher values of equivalent stress are considered as hard orientations and orientations with lower values of equivalent stress are considered as soft orientations. In reference to classical crystal plasticity studies, it can be said that soft (resp. hard) orientations are such which have high (resp. low) Schmid-like factor under given loading conditions. (a) L=-1 (b) L=1 (c) L=0 Figure 2: Equivalent stress - equivalent strain for single crystals without voids for di erent Lode parameters and di erent crystal orientations.The presented curves are the same for all triaxiality values considered. 3. Results and Discussion Analytical studies on the void growth in an anisotropic rigid-plastic matrix described by the quadratic criterion performed by (Benzerga and Besson, 2001) for spherical voids, and continued by (Keralavarma and Benzerga, 2010; Keralavarma et al., 2020) for ellipsoidal voids, revealed that 9 this growth is governed by the formula _f f(1f)3 h1 ~2 eq(N)sinh3 hm 1 (14) in which ~eq(N) depends on the matrix plastic anisotropy, changing with the loading directions Nwith respect to the main anisotropy axes, mis the mean stress, 1the yield stress in the selected direction, while his the so-called net anisotropy invariant calculated as a single value for a given material. Accordingly an exponential void growth is a ected by the triaxiality Tandh, while stress direction Nmodi es mainly the proportionality coecient through ~ 2 eq. In analogy to this study Paux et al. (2015) proposed an approximate condition for the voided single crystal described by the non-quadratic regularized Schmid law, bearing mathematical similarities with the power-law rate-dependent model studied here. Authors have assessed the coecient 1=has equal to 0:506 for the fcc crystals, so it can be assumed as constant here. Therefore, qualitatively, mainly the dependencies on triaxiality and the plastic anisotropy are expected to be observed in the present analysis. Nevertheless, as it will be shown, they are additionally modi ed by the fact that the void is located at the bi-crystal boundary and not in the homogeneous matrix. At this moment such additional dependencies can be analyzed only numerically. In this section results from di erent numerical simulations are presented. We will focus on the joint e ect of stress state ( T;L) and matrix crystallographic orientations (i.e. material anisotropy) on void volume fraction evolution and void morphology. We will compare results from single crystals with bi-crystals and results from di erent bi-crystal RVEs. In this work results of single crystal are used as reference values to compare with bi-crystals. Calculations of single crystals with voids are performed for four di erent crystallographic orientations O1, O2, O4, O6, 3 Lode parameters ( L=1;0;1) and 4 values of triaxialities ( T= 0;1=3;2=3;1). Based on the previous studies on the void growth in fcc single crystals Srivastava and Needleman (2015b), at this low and moderate triaxiality level, porosity growth is expected to be rather slow, void shape change signi cant and Lode parameter to have an important in uence on void behavior in this respect. It should be added that for the advanced regime of the considered processes the changing shape of a void and a unit cell may additionally a ect the dependencies observed in the calculations (Keralavarma and Benzerga, 2010). 3.1. Void growth in a bi-crystal For studying the e ect of grain orientation on the void growth in bi-crystals, three di erent representative volume elements (RVE) representing three di erent microstructures are created. The di erence between each RVE are the crystallographic orientation assigned to the matrix of each half cell. Di erent orientations are assigned for half cell 1 and half cell 2, such that three di erent bi-crystal RVEs are formed. The orientations sets for the simulations of bi-crystals are selected from the cell calculations of single crystals without voids for di erent Lode parameters. By using the hard and soft orientations data, three di erent RVEs of bi-crystals representing; 1.soft-soft; 2. hard-hard; 3. soft-hard are created. Orientations used for soft-soft RVE are O1 and O2 for half cell 1 and half cell 2, respectively. Similarly for hard-hard RVE, orientations used are O6 and O4 for half cell 1 and half cell 2, respec- tively. For soft-hard case, the orientations considered are O1 and O4 for half cell 1 and half cell 2. 10 It can be veri ed that for all three cases the grain boundary between two halves of bi-crystal is a high angle boundary. The misorientation angle and misorientation axis speci ed by the common crystallographic direction of two orientations is for soft-hard 56 :6andh0:590;0:769;0:245i, for hard-hard is 35 :6andh0:598;0:800;0:046i, while for soft-soft is 45andh100i. Note that for soft-soft bi-crystal the grain boundary is a twist boundary for which its unit normal is coaxial with misorientation axis. In the remaining two cases these two vectors are inclined with respect to each other. For bi-crystals without voids, gs. 3(a) and 3(b) represent equivalent stress- equivalent strain curves for soft-hard orientation when di erent Lode parameter values ( L=1 andL= 1) are prescribed in the whole cell. Green lines represent stress-strain relations for (soft) half cell 1 and red lines represent stress strain relation for (hard) half cell 2. Blue lines represent stress -strain curves when the whole bi-crystal is considered. The presented curves are the same for any triaxiality value. For comparison purposes, equivalent strain is calculated in all the cases considering the volume of the whole cell. Note that due to imposed periodic boundary condition analyzed bicrystal without a void represent the behaviour of the laminate with alternating layers with O1 (soft) and O4 (hard) crystallographic orientation. Thus the predicted stress and strain elds within each half cell are uniform, however, considerably di erent from each other and from the overall behaviour. For L= 1, g. 3(c) presents stress-strain curves for soft-soft orientation. ForL=1, g. 3(d) shows stress-strain relations for the hard-hard case. For L=1 andL= 0, for the soft-soft bi-crystal case, stress-strain curves overlaps for both half cells and the full cell. ForL= 1 andL= 0, stress-strain curves for hard-hard bi-crystal are similar as for L=1 . These last cases are not shown in the document for the sake of brevity. For bi-crystals with a void, gure 4 shows FE mesh and the 3 RVEs used in the present study with their respective initial pole gures. The projection of misorientation axis is also marked in the pole gure. The initial void volume fraction of each half cell is f0= 0:0022. Void volume fraction is calculated as f=Vvoid=Vcell, whereVcellis the volume of the whole cell and Vvoidis the volume of the void in the half cell. Calculating the void volume fraction for symmetric orientations O1, O2, O3 is straight forward, as the cell faces remain straight and the current volume of the cell is simply given by (Ux+Dx) ·(Uy+Dy) ·(Uz+Dz), whereDx;Dy;Dzare the initial lengths of the cell in X,Y,Z direction respectively and Ux;Uy;Uzare the displacement applied at each direction. For non symmetric orientations, the initially straight sides of the unit cell along which fully periodic conditions are applied becomes curved. Since the cell faces are not straight it is not direct to calculate the current volumes of the half cells and the void. To face this problem, the following approach has been adopted in this work: The 3D (X,Y,Z) coordinates of all the nodes on the surface of the void (initially spherical) is extracted from the nite element calculations for all the time steps. By using these coordinates, void volume fraction of each half cell is calculated at each time increment with the help of convex hull algorithms (convex hull of a given set of points is the smallest convex polygon that contains all the points). 3.1.1. Void growth in a soft-soft bi-crystal. First we will present some results related to soft-soft bi-crystal RVE. This RVE has O1 ori- entation for half cell 1 and O2 orientation for half cell 2. Evolution of void volume fraction f with respect to the equivalent strain "eqvfor Lode parameter L=1 is shown in gure 5(a) for 11 (a) Soft-hard bi-crystal, L=-1 (b) Soft-hard bi-crystal, L=1 (c) Soft-soft bi-crystal, L=1 (d) Hard-hard bi-crystal, L=-1 Figure 3: Equivalent stress - equivalent strain in half cell 1, half cell 2 and the whole bi-crystal without void for soft-hard, soft-soft and hard-hard orientation for di erent Lode parameters. The presented curves are the same when any value of triaxiality is prescribed in the whole cell. di erent triaxialities ( T). Continuous lines represent void volume fraction evolution in half cell 1 and dotted lines in half cell 2 and di erent colours represents di erent triaxiality values. For triaxilityT= 1,T= 2=3 andT= 1=3, we can see that void volume fraction fincreases with "eqv in both half cell 1 and half cell 2, whereas for triaxiality T= 0, void volume fraction decreases belowf0in both half cell 1 and half cell 2. This behaviour of void volume decreasing below f0was also observed for porous single crystal by Srivastava and Needleman (2015a). In the case of low triaxialities, simulations are terminated when the outer surface of the voids comes in contact with each other. In general void grows faster at moderate triaxialities ( T= 1;2=3) compared to low triaxialities ( T= 1=3;0), this trend has been previously observed by Srivastava and Needleman (2015a); Yerra et al. (2010); Srivastava and Needleman (2013) for porous single crystals. Com- paring void volume fraction evolution between the two half cells, we observe that di erence in void volume evolution between two half cells is hardly noticeable. This results from the fact that for this case direction X of major stress 16=2=3coincides with the same crystallographic 12 Figure 4: Di erent bi-crystal RVEs and respective pole gures representing the orientations used for studying the e ect of matrix orientation on the void behaviour. Blue colour represents orientations in half cell 1, green colour represents orientation in half cell 2 and red colour represents the projection of misorientation axis . direction [001] in both half cells. As seen from gure 5(b), for Lode parameter L= 1, void grows slower than for L=1 (note that the vertical axis of both gures has di erent scale). We can see di erences in void growth behaviour of half cell 1 and half cell 2 for all triaxialities considered. In this case direction Z of a minor stress 36=1=2varies between both half cells (see Table 1). As seen in Fig. 3(c) the stress response in two half cells is then di erent, so the level of equivalent stress in half cell 1 is lower than in half cell 2. Following relation (14), void grows faster in half cell 1 (continuous lines) compared with void growth in half cell 2 (dotted lines) at moderate triaxialities ( T= 1;2=3). However, at low triaxialities ( T= 0;1=3) void volume in half cell 1 decreases faster as compared to void in half cell 2. To nd a source of such behviour the local triaxialities have been calculated for the Cauchy stress averaged over each half cell separately. It has been found that the local triaxility in half cell 1 for case T=0 is negative and equal to 0:2, while for half cell 2 it is higher and equal to 0:23. Similarly, respective values for the case T=1/3 are 0 :22 and 0:44, correspondingly. For the remaining cases local triaxialities per half cell are approximately equal to the macroscopic ones. All the values are collected in the table B.6 For Lode parameter L= 0, the behaviour of the void in both half cell 1 and half cell 2 are similar to the case of Lode parameter L= 1, hence plot of f"eqvis not presented. Among the Lode parameters considered, void growth rate is in general higher for Lode parameter L=1, followed by L= 0 andL= 1. 13 (a) L=-1 (b) L=1 Figure 5: Void volume evolution at di erent stress triaxialities and Lode parameters for the soft-soft bi-crystal. Matrix orientation of half cell 1 is O1 and half cell 2 is O2. (a) L=-1 (b) L=1 Figure 6: Void volume evolution at di erent stress triaxialities and Lode parameters for the hard-hard bi-crystal. Matrix orientation of half cell 1 is O6 and half cell 2 is O4. 3.1.2. Void growth in a hard-hard bi-crystal. For hard-hard bi-crystal RVE, half cell 1 has orientation O6 and half cell 2 has orientation O4. Similarly to previous comparison, void volume fraction evolution for di erent Lode parameters and triaxialities are plotted in gures 6. For Lode parameter L=1, void growth in half cell 1 is quite similar to void growth in half cell 2 for triaxiality ( T= 0;1=3;1). For triaxiality T= 2=3, we observe a di erent void growth behaviour in both half cells. Void in half cell 1 grows faster compared with void in half cell 2. 14 For Lode parameter L= 1, we can see distinction in void growth between half cell 1 and half cell 2 at moderate triaxialities ( T= 1;2=3), as void in half cell 2 grows faster than void in half cell 1. For the remaining triaxialities ( T= 0;1=3), there is no signi cant di erence in the void growth in two half cells. Similarly to the soft-soft bicrystal, when calculating local triaxilities per each half cell it has been found that the local triaxilities in half cell 2 are higher than in half cell 1. ForL= 0, the di erence in void growth in half cell 1 and half cell 2 is small. Void in half cell 1 grows faster than void in half cell 2. Plots are not presented for L= 0 for brevity. As in the case of soft-soft bi-crystal RVE, void growth rate is higher for Lode parameter L=1, followed byL= 0 and nally by L= 1. 3.1.3. Void growth in a soft-hard bi-crystal. For soft-hard bi-crystals, half cell 1 has orientation O1 (soft) and half cell 2 has orientation O4 (hard). Evolution of fversus"eqvfor half cell 1 and half cell 2 for di erent values of Tand Lare given in gures 7. The way results are presented in gures 7 are di erent as in previous sections to clearly show the changes occurring at low triaxialities ( T= 0;1=3), as the amount of void growth is small compared to moderate triaxialities ( T= 1;2=3). ForT= 0, as we can see from gure 7(a), void volume fraction in both half cell 1 and half cell 2 decreases below the initial void volume fraction f0. For half cell 1, void volume fraction increases slightly in the early stages of deformation, but starts dropping as deformation progresses, whereas void volume fraction fin half cell 2 drops faster than fin half cell 1 for all Lode parameters. Among the 3 Lode parameters considered, void volume fraction drops faster for L= 1, followed byL= 0 and least for L=1 in both half cells. Considering T= 1=3, as shown in gure 7(b), at L=1, we can see that the volume of the void in half cell 1 increases, whereas in half cell 2, void volume fraction grows slightly up until "eqv= 0:2, but later on the void volume fraction decreases. For L= 0, void volume fraction in half cell 1 increases up to "eqv= 0:3, but for"eqvgreater than 0.3 it starts decreasing, and during the nal stages of deformation, void volume decreases below the initial value f0= 0:0022. For void in half cell 2 even though we see a slight increase in funtil"eqv= 0:2, at higher "eqvvoid volume starts decreasing below f0. ForL= 1, void in half cell 1 and void in half cell 2 behaves similarly to the case of L= 0. The behaviour of the void in half cell 1 is similar for both triaxialities T= 2=3;1, as shown in gures 7(c) and 7(d) as we observe void growth. For void in half cell 2, we observe void growth forT= 1 andL=1;0;1, but forT= 2=3 void growth is very small and gradually during the course of deformation, we observe void growth saturation (i.e. no change in f) forL= 0 and L=1. ForT= 2=3;1, void in half cell 1 grows faster than void in half cell 2 irrespective of Lode parameters. In general void grows faster for L=1, followed by L= 0 andL= 1. The di erence in the behaviour of the void in both half cells are signi cant. For moderate triaxialities, the e ect of Lode parameter is more prominent for void in half cell 1 (soft) compared to void in half cell 2 (hard). For all triaxialities and Lode parameters values considered, falways grows faster (or shrinks slower) in half cell 1 (soft orientation) compared with half cell 2 (hard orientation).This general trend is following the formula (14) since the level of equivalent stress is much higher in the half cell 2 than in half cell 1, as seen in Figs. 3(a),3(b). Additionally, when calculating the average triaxialities per each half cell we have found almost two times higher values 15 in half cell 1 (soft) than in half cell 2 (hard). The values are included in Table B.5. (a) T=0 (b) T=1/3 (c) T=2/3 (d) T=1 Figure 7: Void volume evolution at di erent stress triaxialities and Lode parameters for the soft-hard bi-crystal. Matrix orientation of half cell 1 is O1 and in half cell 2 is O4. 3.2. Deformed Void shape in a bi-crystal In this section, we will present results on the in uence of matrix orientation and stress state on void shapes and accumulated shear distribution =P Rt 0j_ jdton void shape evolution in di erent bi-crystal RVEs. 16 (a) T=0 (b) T=1 Figure 8: Deformed void shape and distribution of accumulated shear for soft-hard bi-crystal at stress triaxiality T= 0;1 and equivalent strain "eqv= 0:45 forL= 0;1;1: 3.2.1. Deformed void shape in a soft-hard bi-crystal For soft-hard bi-crystal, at triaxiality T= 0, gure 8(a) shows the cross section of the deformed cell in two di erent planes (XY and XZ) at "eqv= 0:45 and the contours showing distribution of the accumulated shear for Lode parameters L= 0;1;1. As can be clearly seen from the plots, the distribution of and void shape in half cell 1 and half cell 2 are di erent for all the Lode parameters considered. For L= 0, the deformed cross-section of the cell along XY plane shows the void evolving into an ellipsoid on both half cells with the major axis of ellipsoid along 17 X axis (direction of the maximum applied stress). Looking at deformed cross section along XZ plane, we see the void in half cell 2 attening into an elliptical crack like shape, with its major axis along X axis and the crack like surface arising from collapsing of void located away from the grain boundary, whereas void shape in half cell 1 is an oblate ellipsoid. The accumulated shear values are higher in half cell 2 compared with half cell 1 and the highest values of the accumulated shear is observed near the collapsed void surface in half cell 2, with a maximum value of 4.71. From the deformed cell cross-section along XY axis, for L= 1, we see the void shape in both half cell 1 and half cell 2 to be circular at "eqv= 0:45. As in the case of L= 0, the void in half cell 2 collapses into a crack in XZ plane, with its tip along X axis and away from the grain boundary, whereas for void in half cell 1, the void evolves into an oblate ellipsoid. The value of the accumulated shear is higher in half cell 2 as compared to half cell 1 with a maximum value of 3.24. ForL=1, the void shape in both half cell 1 and half cell 2 are very similar, and void evolves into an ellipsoid with the major axis along X. At "eqv= 0:45, we do not see void collapsing, but at higher strains void tend to collapse in both half cell 1 and half cell 2. Unlike L= 0;1, forL=1, void tends to collapse along Y axis and void radius evolves in the order Rx> Rz> Ry. Similar to the cases of L= 0;1, the values of accumulated shear is higher in half cell 2. For the 3 Lode parameters analysed, the values of accumulated shear is higher for L=1, followed by L= 0 and least by L= 1. Similar analysis has been performed for T= 1 at Lode parameters L= 0;1;1 as shown in gure 8(b). As we discuss previously from gure 7(d), at T= 1 void growth is observed for all Lode parameters and the shape of the voids in both half cell 1 and half cell 2 are ellipsoid in all the cases, but there is a clear di erence in the amount of growth within each half cell, as well as the void growth along each principal loading direction. Void radius is di erent along each direction and can be characterized as follows: for "eqv= 0:45 and for L= 0, void radius is of the order Rx> Ry> Rzand forL= 1, order is Rx=Ry> Rz, nallyL=1 void radius is of order Rx>Rz>Ry. Accumulated shear values are higher in half cell 2 compared to half cell 1 and the value of the accumulated shear is higher for L=1, followed by L= 0 and least in L= 1. 3.2.2. Deformed void shape in a hard-hard bi-crystal In this section we will study void shape and accumulated shear distribution in the hard-hard bi-crystal, where the matrix orientation for half cell 1 is O6 and for half cell 2 is O4. Figure 9 shows the cross section of the deformed cell in XY and XZ planes, for L= 0;1;1 andT= 0 when"eqv= 0:45. For a hard-hard bi-crystal, we observe that the whole cell is deforming in such a way that it bends around the grain boundary. Such behaviour is understood when observing deformation of a respective bicrystal without a void for which in addition to axial deformation also shearing in opposite directions in two half cells is found. This shearing in XZ plane is also seen in Fig. 9. This is because of the anisotropic nature of the cell and the incompatibility between both half cells when hard-hard crystals comes in contact. This deformation behaviour of the cell has a huge in uence on the way the void shape evolves. For L= 0, void shows elliptical shape with the major axis along X axis of the cell on both XY and XZ plane. Void shape is not symmetric along X axis of the cell and void radius along Z axis is higher in the upper part of the cell compared to the lower part along the grain boundary. The trend for void collapsing is clearly evident from the shape of the void. Void tends to collapse in half cell 2 faster than in half cell 1. 18 Figure 9: Deformed void shape and distribution of accumulated shear for a hard-hard bi-crystal at stress triaxiality T= 0 and equivalent strain "eqv= 0:45 forL= 0;1;1 Similarly, for L= 1 and along XY plane, we observe void shape being circular and along XZ plane void shape is irregular and it is unsymmetrical along X axis. Similar to the case of L= 0, void radius on the upper part of the cell in Z axis is higher than void radius in the lower part of the cell and void growth is observed along X and Y axis, whereas void collapse along Z axis. Finally forL=1, void shape is elliptical in both half cells and in both XY and XZ planes. Void radius evolution is of the order Rx>Rz>Ryand void tend to collapse along Y axis. The location of the maximum value of accumulated shear is close to the void in half cell 2 for all 3 Lode parameters considered, and away from void we observe a homogeneous distribution of accumulated shear and similar values in both half cell 1 and half cell 2. At "eqv= 0:45 andT= 0, the maximum is found in L= 1 and has the value 4 :42. ForT= 1 (results not shown here) and all Lode parameters considered, we observe elliptical void shapes in both XY and XZ planes. Void shape is similar in both half cells and the accumulated shear distribution is similar to the case T= 0. 19 3.3. E ect of neighbouring grain orientation (a) soft-hard, soft-soft and single crystal with orientation O1 (b) soft-hard, hard-hard and single crystal with orientation 04 Figure 10: soft-hard, soft-soft, hard-hard bi-crystal and single crystal RVEs used for studying the e ect of neigh- bouring grain orientation on the void behaviour. Blue colour represents orientations in half cell 1 and green colour represents orientation in half cell 2. In this section, three di erent comparisons are made to study the e ect of neighbouring grain orientation on void growth, 1. a matrix with hard orientation neighbour, 2. a matrix with soft orientation neighbour, and 3. neighbour grain with the same orientation (i.e., single crystal). For that end, rst we will compare void evolution in half cell 1 of the soft-hard bi-crystal and soft-soft bi-crystal. As shown in gure 10(a), soft-hard bi-crystal and soft-soft bi-crystal has the same orientation for half cell 1 i.e., O1, but orientation of half cell 2 in soft-soft bi-crystal and soft-hard bi-crystal is di erent. Second, we will show results of void evolution in half cell 2 of hard-hard bi-crystal and soft-hard bi-crystal. For hard-hard bi-crystal and soft-hard bi-crystal orientation of half cell 2 is the same (i.e., O4 orientation), but orientation of half cell 1 in hard-hard bi-crystal and in soft-hard bi-crystal is di erent (see gure.10(b)). Similarly, void behaviour is compared between half cell in a bi-crystal and half cell in a single crystal. The orientation of the single crystal is the same as half cell in the bi-crystal. For the case of the soft-hard bi-crystal, void in half cell 1 is compared with void in single crystal with the same orientation, i.e., O1, whereas for void in half cell 2 comparison is made with void in single crystal with the same orientation i.e., O4. By comparing the void behaviour in these cases, we can better understand the e ect of neighbouring grain orientation on void behaviour and its morphology. 20 3.3.1. Void growth in soft-hard, soft-soft bi-crystals and single crystal (a) L=-1 (b) L=1 (c) L=-1 (d) L=1 Figure 11: Comparison of void volume fraction-equivalent strain evolution for half cell 1 in soft-hard, soft-soft bi- crystals and in a single crystal with the same orientation as half cell 1 at di erent stress triaxialities T= 0;1=3;2=3;1 and Lode parameters L= 1;1. In the case of soft-hard bi-crystal, soft-soft bi-crystal and single crystal, we will compare the e ect of hard, soft and identical neighbouring matrix, respectively, on the void evolution in a soft matrix (half cell 1 is always O1 orientation). Figures 11 show void volume fraction-equivalent strain at di erent triaxialities and Lode parameters L=1;1. Figures 11(c) and 11(b) shows void volume evolution in half cell 1 for soft-soft and soft-hard bi-crystals. Figures 11(a) and 11(d) show void evolution in half cell 1 of soft-hard bi-crystal and of a single crystal. Single crystal calculations are performed for the full cell, but for comparison with bi-crystals, only half of the 21 cell is taken into account. For Lode parameter L=1, void in half cell 1 of soft-hard bi-crystal (continuous lines) grows faster than in soft-soft bi-crystal (dotted lines) for triaxialities T= 1;2=3 as shown in gure 11(c). For triaxiality T= 1=3, void volume fraction in soft-hard bi-crystal increases slowly, whereas for soft-soft bi-crystal, void stops growing (i.e. void volume saturates after initial growth). For triaxiality T= 0, void volume fraction in both bi-crystals drop below f0from the beginning, and void in soft-soft bi-crystal drops faster than in soft-hard bi-crystal. ForL= 1 andT= 1;2=3, as we can see from gure 11(b), void in half cell 1 of soft-hard bi-crystal grows faster than void of soft-soft bi-crystal and the di erence in void volume evolution between both bi-crystals is smaller compared to L=1. ForT= 1=3, the amount of void growth in both bi-crystals is very less and void volume fraction in soft-hard bi-crystal is higher than in soft-soft bi-crystal. For T= 0, void volume in both bi-crystals drops below the initial volume f0, dropping faster in soft-soft bi-crystal than in soft-hard bi-crystal. Void evolution in half cell 1 for soft-hard bi-crystal and single crystal is presented in g- ures 11(a), 11(d). For L=1 and for all triaxialities considered, void growth in half cell 1 in a soft-hard bi-crystal is higher than void growth in a single crystal. Interestingly, for L=1;if we compare the void volume evolution of single crystal with half cell 1 of soft-soft bi-crystal shown previously in gure 11(c), we observe that the response of the void is the same for both cases at all triaxialities considered. This is explained by the same reason as almost the same response in two half cells seen in Fig. 5(a). ForL= 1 andT= 1 we observe that void growth in half cell 1 of soft-hard bi-crystal is much higher than void growth in a single crystal. For T= 2=3, 1=3, void growth in half cell 1 of soft-hard bi-crystal is higher than the single crystal case, but the di erence in value of fis not as big as forT= 1. ForT= 0, void volume fraction decreases below f0for both cases and void in half cell 1 of soft-hard bi-crystal drops slightly faster than void in single crystal. For L= 1 void behaviour in half cell 1 of soft-soft bi-crystal is closer to the void behaviour of the single crystal case. The observed variation in the void growth for O1 orientation placed in these 3 di erent con g- urations for a given Lode parameter can be correlated with the hierarchy of an average triaxiality in the corresponding half cell 1 demonstrated in Table B5 and B6. Note that for a single crystal case its value is equal to the imposed triaxility. 22 3.3.2. Void growth in soft-hard, hard-hard bi-crystals and single crystal (a) L=-1 (b) L=1 (c) L=-1 (d) L=1 Figure 12: Comparison of void volume evolution-equivalent strain for half cell 2 in soft-hard, hard-hard bi-crystals and and single crystal with the same orientation as half cell 2 at di erent stress triaxialities T= 0;1;1=3;2=3 and Lode parameters L= 1;1. In this section, the behaviour of the void in half cell 2 of a soft-hard bi-crystal is compared with the behaviour of the void in half cell 2 of hard-hard bi-crystal and single crystal with O4 orientation. The e ect of a softer, harder and identical neighbouring matrix on the evolution of a void located in a hard matrix (half cell 2 is always O4 orientation) is compared and discussed. Figures 12(a) and 12(b) show void volume fraction evolution with equivalent strain for voids of half cell 2 in a soft-hard (continuous lines) and in a hard-hard bi-crystal (dotted lines). For Lode parameter L=1, void volume fraction in half cell 2 of a hard-hard bi-crystal grows faster 23 than half cell 2 in a soft-hard bi-crystal for triaxialities T= 1;2=3;1=3. The di erence in void volume fraction evolution of both bi-crystals is less remarkable for T= 1=3 than forT= 2=3;1. For T=0, void volume fraction for both bi-crystals drop below initial void volume fraction f0and void volume fraction of half cell 2 in soft-hard bi-crystal drops faster than half cell 2 in hard-hard bi-crystal. As we can see from the gure, void grows faster when the half cell has hard neighbour compared to soft neighbour. Figure 13: Cross section showing deformed void shape and accumulated shear for a (soft) single crystal, a (hard) single crystal, soft-hard bi-crystal, soft-soft bi-crystal and hard-hard bi-crystal for T= 0,L= 1 and"eqv= 0:45. The soft orientation in half cell 1 is always O1 and the hard orientation in half cell 2 is always O4 for con gurations containing soft or hard component, respectively (see gure 10(a) and 10(b)). Figure 12(b) shows void volume fraction evolution for Lode parameter L= 1. ForT= 1, void volume fraction in both bi-crystals initially increases at the same rate, but at "eqv>0:15 void volume in a hard-hard bi-crystal increases faster than in a soft-hard bi-crystal. For T= 2=3, we observe contrasting void behaviour in both bi-crystals: as fof a hard-hard bi-crystal gradually increases at constant rate, fof a soft-hard bi-crystal start decreasing below initial void volume fractionf0. ForT= 1=3 andT= 0, void volume fraction decreases in both bi-crystals below initial void volume fraction f0, dropping faster in half cell 2 of a soft-hard bi-crystal than in a hard-hard bi-crystal. Similarly, gures 12(c) and 12(d) show comparisons of half cell 2 void evolution for soft-hard bi-crystal (continuous lines) and for a single crystal with the same orientation as half cell 2 (dotted lines). Single crystal calculations are performed considering the full cell, but for comparison with 24 bi-crystals, only the behaviour of half cell is considered. For L=1 and comparing soft-hard bi-crystal with a (hard) single crystal, we observe that, at T= 1;2=3, void growth in the single crystal is higher than in the soft-hard bi-crystal. At T= 1=3;0, void evolution fin the soft-hard bi-crystal is very similar to the single crystal case. For L= 1 and at T= 1, void growth in the single crystal is higher than in the soft-hard bi-crystal. For T= 2=3, we observe void growing in di erent way compared to other cases, as void volume fraction in both soft-hard bi-crystal and single crystal increases initially, but starts decreasing during nal stages of deformation and void in single crystal drops faster than in soft-hard bi-crystal. For T= 1=3;0 void volume decreases below initial void volume fraction for both cases and void volume fraction in single crystal drops faster than void in soft-hard bi-crystal. Overall, for L=1;1 and atT= 1, void growth is higher in a hard-hard bi-crystal, followed by a (hard) single crystal and nally by a soft-hard bi-crystal. For T= 1 andL=1, void evolution growing faster for orientation O4 than orientation O1 were also found by Ling et al. (2016) for a voided single crystal. At T= 2=3,fis higher in a hard-hard bi-crystal and the void behaviour of a soft-hard and a (hard) single crystal has similar trend. At T= 1=3;0, void volume fraction decreases faster in a single crystal, followed by a soft-hard bi-crystal and nally by a hard-hard bi-crystal. The tendency for void collapse is higher in the single crystal case. For all triaxialities tested, the behaviour of the hard-hard bi-crystal is in general closer to the response of the single crystal case. The calculated average triaxialities per each half cell included in Table B4 and B5 conform with the exposed tendencies. An additional role in shaping the void evolution is played by the speci c deformation of a unit cell in the case of the hard-hard bi-crystal as demonstrated in the next section. 25 3.4. E ect of neighbouring grain orientation on the deformed void shape. Figure 14: Distribution of misorientation angle with respect to the initial orientation for soft-hard bi-crystal for T= 0,L= 1 and equivalent strain 0 :45, Pole gures 111 at the bottom present the spread of the current crystal orientations for each Gauss point in the cell (pole gure on the left is a density plot, pole gure on the right { discrete point plot. It is clearly evident from the results presented in the previous section that void growth is a ected by neighbouring grain orientation and, depending on the Lode parameter and triaxiality values, the behaviour of the void is either closer to a single crystal or its respective bi-crystal counterpart. In this section, the e ect of neighbouring grain orientation on the void shape is presented and discussed. As in the previous section, the cases that are going to be compared are: 1. Void in half cell 1 in a soft-hard bi-crystal, void in half cell 1 in a soft-soft bi-crystal and in a single crystal with the same orientation as half cell 1 (O1). 2. Void in half cell 2 in a soft-hard bi-crystal, void in half cell 2 in a hard-hard bi-crystal and in a single crystal with the same orientation as half cell 2 (O4). As noted by Srivastava and Needleman (2013, 2015a), changes in void shape in a initially spherical void are more signi cant at low triaxialities than at moderate or high triaxiality values. We now focus on studying the in uence of neighbouring grain orientation on the deformed void 26 shape. We choose to present only one case here, i.e., triaxiality value T= 0 andL= 1 as a representative example. Deformed void shapes are compared in di erent RVEs at "eq= 0:45. Fig. 13 shows the deformed void shape in single crystals with O1 and O4 orientation respec- tively, soft-soft, soft-hard and hard-hard bi-crystals. By comparing deformed void shape in half cell 1(left-half) in a soft-hard bi-crystal and in a single crystal with O1 orientation, we observe that void in half cell 1 is ellipsoidal in both RVEs at both XY and XZ planes, with void radius in orderRxRy>Rz. In the soft-soft bi-crystal case, we observe that void within half cell 1 is also ellipsoidal but the order of void radius is Rx>Ry>Rz. In a soft-hard bi-crystal void shape in half cell 2 (right-half) is quite similar to the void shape in a single crystal with O4 orientation in both XY and XZ planes. As we see from gure 13 for the aforementioned RVEs, void tends to collapse into a "crack like" shape along X direction, and void radius for both RVEs are in order RxRy> Rz. For a hard-hard bi-crystal, void evolves into irregular shape before collapsing along XZ plane, but void collapse occurs at rate slower than for single crystal. By comparing accumulated shear in half cell 1 among soft-hard, soft-soft and single crystal RVE, the maximum value of accumulated shear is observed in single crystal with O1 orientation (soft orientation) in the vicinity of the void in XZ plane. Similarly comparing accumulated shear in half cell 2 among soft-hard, hard-hard and single crystal RVE, the maximum value of is obtained in single crystal with O4 orientation (hard orientation) in XZ plane. Overall, for T= 0,L= 1 and at "eq= 0:45, values are higher in (hard) single crystal, followed by hard-hard bi-crystal, soft-hard bi-crystal, (soft) single crystal and nally soft-soft bi-crystal. Figure 15: Distribution of misorientation angle with respect to the initial orientation for soft-hard bi-crystal for T= 0,L=1 and equivalent strain 0 :45, Pole gures (111) at the bottom present the spread of the current crystal orientations for each Gauss point in the cell (pole gure on the left is a density plot, pole gure on the right-discrete point plot.) 27 Figure 16: Distribution of misorientation angle with respect to the initial orientation for hard-hard bi-crystal for T= 0,L= 1 and equivalent strain 0 :45, Pole gures 111 at the bottom present the spread of the current crystal orientations for each Gauss point in the cell (pole gure on the left is a density plot, pole gure on the right-discrete point plot.) Figure 14, presents the in uence of void on the heterogeneity of lattice rotation in the soft-hard bi-crystal, at T= 0 andL= 1. Let us stress that initially all elements within each half cell have the same orientation. While for a bi-crystal without a void, lattice rotates uniformly within each half cell (results are con rmed by calculations, but not presented here for brevity), strong variation is observed when the voids are present at the grain boundary. What is not obvious at rst is that heterogeneity is more pronounced within hard half cell as initially homogeneous orientation in half cell is fragmented into two distinct orientations after deformation, as visible in the left pole gure. This is probably caused by the fact that di erent set of slip systems are operating in the top and bottom part of the half cell and the lattice is rotating in opposite directions in these two parts. Note that the largest misorientation angle ( >45o) is observed close to the void boundary, where the void is distorted most. Figure 15 presents analogical results for the same bi-crystal but forL=1 andT= 0. Here, also the larger misorientation angles are found in the hard 28 half cell. Note that in both cases observed, di erences in misorientation angles between elements belonging to the hard half cell are often higher than 15o, which suggests, appearance of new high angle boundaries and possibility of grain re nement. Smaller di erences in misorientation angles in two half cells are observed for hard-hard bi-crystal in gure 16. Nevertheless, also in this case the presence of void leads to the non-uniform lattice rotation leading to an increase of a lattice curvature. 4. Summary and Conclusions Three-dimensional nite element cell calculations are carried out with the aim of analysing the response of voids at the grain boundary of a bi-crystal RVE. The e ect of crystal orientation on porosity evolution are studied for three di erent cases: soft-hard, soft-soft and hard-hard RVEs subjected to di erent loading conditions. Unit cells consisting of cubic cells divided in two halves with a spherical void at the grain boundary with total void volume fraction in the cell f0= 0:0044 and volume of void in each half cell f0= 0:0022 are analysed. The responses of bi-crystal RVEs are studied for four stress triaxiality values: T= 0;1=3;2=3;1, three Lode parameters L=1;0;1 with a strictly speci ed orientation of loading directions with respect to the initial bi-crystal orientation, and a dimensionless equivalent strain rate _ "eqv=_ 0in the range (100 { 150). In order to select soft and hard orientations rst calculations have been performed for the case of single crystal without voids. For a given Lode parameter, the orientation with higher values of equivalent stress are considered as hard orientations and orientations with lower values of equivalent stress are considered as soft orientations. The evolution of void volume fraction strongly depends on the overall stress triaxiality and more speci cally on the triaxiality for the stress averaged over each half cell for all bi-crystal RVEs considered. When the level of strain increases the Lode value of macroscopic loading importantly modi es this dependencies. Moreover the Lode value has a prominent e ect on the evolving shape of void, which especially seen at low triaxialities. Even though all bi-crystal RVEs considered have high angle grain boundaries, signi cant di erence in void growth between the two half cells occurs specially for soft-hard bi-crystal case. The di erent behaviour of the void evolution in each half cell can be attributed to the di erence in the strength (soft-hard) between the two initial orientations (O1 and O4) and the incompatibility between the two half cells. In the soft-hard bi-crystal case, at high triaxialities, void grows faster in softer crystal compared to harder crystal counterpart, whereas at low triaxialities, void tends to collapse faster in harder half cell compared to neighbouring softer crystal. Void growth is higher for higher triaxialities and void growth is faster forL=1, followed by L= 0 andL= 1. Void shape evolution in bi-crystals depends on triaxiality, Lode parameter and initial orien- tation of each grain. At high triaxialities, void tends to evolve into spherical/ ellipsoidal shapes irrespective of Lode parameter and grain orientation. At low triaxialities, void shape is greatly a ected by the initial bi-crystal orientation with respect to the assumed loading conditions spec- i ed by the Lode parameter as the void tend to evolve into ellipsoidal shapes for soft half cells, and into crack like shapes for hard half cells. Void shape in half cell 2 of soft-hard bi-crystal looks more similar to hard (O4) single crystal, rather than half cell 2 of hard-hard bi-crystal. Void shape in half cell 1 of soft-hard bi-crystal, half cell 1 of soft-soft bi-crystal and soft single crystal looks similar. 29 If the volume and morphology of voids with the same matrix orientation, but di erent neigh- bouring cell orientation are compared we observe that void tend to behave like the single crystal case (same orientation as bi-crystal half cell) if the di erence in strength between two half's of the bi-crystal is small, i.e., in the cases soft-soft and hard-hard. When the di erence in strength between the two half cells is higher (i.e., soft-hard) void behaviour tend to move away from the respective single crystal behaviour. At high triaxialities, the profound e ect the orientation of the neighbouring grain has on the evolution of the void volume becomes even more important than the orientation of the grain itself. Irrespective of Lode value, for void in half cell 1 of soft-soft, soft-hard bi-crystals and soft (O1) single crystal, void growth is faster in the soft grain when the neighbouring grain has hard orientation. For void in half cell 2 of soft-hard, hard-hard bi-crystal and hard (O4) single crystal, void grows faster in the hard grain if the counterpart has hard orientation. Slip system activity also depends on triaxiality, Lode parameter and initial crystal orientation and far away from the void becomes similar to a fully dense crystal under the same imposed loading conditions. Depending on crystal orientation we have di erent Schmid factors, which in turn leads to di erent slip activity. Higher the value of stress triaxiality, higher will be the amount of slip activity, and hence higher plastic deformation in the matrix around the void, leading to higher void growth. For a soft-hard bi-crystal, at high triaxiality, the maximum value of accumulated shear occurs around the void in the softer part, and away from the void, slip activity in the harder part is higher. At low triaxiality, slip system activity is signi cantly higher in the harder part and the maximum value of accumulated shear is observed around the void also in the harder part, irrespective of Lode value. Moreover, in soft-hard bi-crystal the heterogeneity of slip activity, especially within the harder crystal increases leading to a non-uniform lattice rotation and consequently grain fragmentation. Slip activities dependence on Lode parameter is clearly evident from our studies. Slip activity is higher for L=1, followed by L= 0 and least forL= 1, irrespective of triaxiality value. For hard-hard and soft-soft bi-crystals, slip activity away from the void is similar in both half cells and around the void maximum slip activity occurs. The crystal plasticity model used in the present study is standard. Observed dependencies can be modi ed when enhanced strain gradient formulations Wul ngho et al. (2013); Stupkiewicz and Petryk (2016); Ling et al. (2018), which introduce size e ects, are applied. For example, as recently demonstrated by Ling et al. (2018) using the reduced micromorphic nite strain crystal plasticity model, the void growth is slow down and the void coalescence delayed when the void size is small as compared to the intrinsic length scale of the material. It happens due to the modi cation of accumulated shear eld close to the void. Similarly, the strain gradient formulation may help to account for the grain boundary e ect on the resistance to plastic ow, as demonstrated by Wul ngho et al. (2013). This may a ect the results obtained in the context of bicrystal. Therefore, in the future more advanced crystal plasticity formulations can be considered for the analysis of the grain boundary e ect on the void growth. Furthermore, to complement numerical results presented in this work and in order to extract de nite conclusions, more numerical simulations should be developed and experimental campaigns focus on studying the evolution of intergranular voids in fcc bi-crystals should be designed and conducted in the future. 30 Acknowledgements The research leading to these results has received funding from the European Union's Horizon 2020 Program (Excellent Science, Marie-Sklodowska-Curie Actions) under REA grant agreement 675602 (Project OUTCOME). We are thankful to Dr. Jose A. Rodr guez-Mart nez of University Carlos III, Madrid, for numerous fruitful discussions. 31 Appendix A. Prescribed Boundary Conditions Figure A.17: Periodic boundary conditions applied to the 3D cell. Boundary conditions applied for a representative volume element is described in detail in this section. In the undeformed cell, the edges of the cell are aligned along the coordinate axes (X,Y,Z), and the origin of the reference coordinate system corresponds to the node a6 (see Figure A16). The faces of the cell are initially straight and the nodes on the opposite faces are connected in order to obtain a periodic response from the cell. Apart from the periodic response, in the present work, constant values of stress triaxiality and Lode parameter are applied to the cell. The procedure to prescribe loading to keep constant values of macroscopic stress triaxiality and Lode parameter throughout the deformation history follows the work of Vadillo et al. (2016). If the ratio of principal macroscopic true stress RandQare de ned in the form: R=2 1;Q=3 1; (A.1) triaxialityTand Lode parameter L(Eq. (10)) can be written as: T=p 2(R+Q+ 1) 3p (1R2) + (1Q2) + (R2Q2);L=2RQ1 1Q(A.2) Since the macroscopic principal true stresses ( 1,2,3) and the normal components of macro- scopic strain rate in the basis of stress principal directions ( _ "11, _"22, _"33) are equal to the volume 32 average values in a cell (Hill (1967)), the total rate of deformation work in the whole cell _Wcan be written as: _W=V1_"11+V2_"22+V3_"33 (A.3) _Wcan be also expressed in terms of the transformed rates of deformation and forces, as: _W=V(I)_"(I)+V(II)_"(II)+V(III)_"(III) (A.4) considering the transformation (A.5) 0 @_"(I) _"(II) _"(III)1 A=N0 @_"11 _"22 _"331 A;0 @(I) (II) (III)1 A=N0 @1 2 31 A; with N=0 @A11A12A13 A21A22A23 A31A32A331 A; and: A11=1p 1 +R2+Q2;A12=Rp 1 +R2+Q2;A13=Qp 1 +R2+Q2 A21=Rp 1 +R2;A22=1p 1 +R2;A23= 0: (A.6) A31=Qp (1 +R2) (1 +R2+Q2);A32=RQp (1 +R2) (1 +R2+Q2) A33=(1 +R2)p (1 +R2) (1 +R2+Q2) If in the transformed coordinate system, the three imposed incremental boundary conditions are prescribed as stress uniaxial: (II)= 0;(III)= 0; _"(I)= _"dummy (A.7) the three prescribed boundary conditions in 1,2,3 directions, using N1relations are therefore: R12= 0;Q13= 0; _"11+R_"22+Q_"33= _"dummyp 1 +R2+Q2 (A.8) Values of _"11;_"22;_"33and hence stress states ( TandL) can be controlled in the cell by xing properly values R and Q and prescribing incremental strain on an added dummy node. The multi- point constraints given in equations A.6, A.8 are implemented into ABAQUS/Standard (2019) via a user de ned subroutine MPC. Due to the anisotropy of the problem, and in order to determine a proper deformation behavior of the cell, periodic boundary conditions should be adopted in the external surfaces of the cell in all three directions. In this work, the general ideas of the implementation of periodic boundary conditions into a Finite Element solver given in (Segurado et al., 2002) is adopted. For the implementation of periodic boundary conditions, it is necessary to couple the displacements of opposite external nodes, on which an average macroscopic strain is allowed. 33 For simplicity and easy implementation, nodes on the surfaces of the cell are categorized into 3 di erent groups: corner nodes, edge nodes and surface nodes (see gure A16 for nomenclature details) The nodea2 is connected to the dummy node by the MPC subroutine, through which the average strains are applied to the system. This means that any deformation that is applied to the cell through the MPC has to be applied only connecting a2 to the added dummy node. In order to prescribe displacements in the outer surfaces which ensure periodicity: uk+ iuk i=Eijxj;i;j= 1;2;3; xj=xk+ jxk j;t+t= 0; t=:n (A.9) beingEijthe macroscopic strain tensor and xjthe position vector, the following set of relations should be applied to the external nodes of the cell in order to avoid constraint linkages. These conditions are implemented in ABAQUS/Standard (2019) by the use of the command *Equation . 1.For nodes on the faces of the cell: uTop(x;y;z )uBottom(x;y;z ) =ua5(x;y;z )ua8(x;y;z ) uRear(x;y;z )uFront(x;y;z ) =ua8(x;y;z )ua4(x;y;z ) uLeft(x;y;z )uRight(x;y;z ) =ua8(x;y;z )ua7(x;y;z ) 2.For nodes on the edges of the cell: ue7(x;y;z )ue3(x;y;z ) =ua8(x;y;z )ua4(x;y;z ) ue3(x;y;z )ue1(x;y;z ) =ua5(x;y;z )ua8(x;y;z ) ue1(x;y;z )ue5(x;y;z ) =ua4(x;y;z )ua8(x;y;z ) ue12(x;y;z )ue11(x;y;z ) =ua8(x;y;z )ua7(x;y;z ) ue11(x;y;z )ue10(x;y;z ) =ua5(x;y;z )ua8(x;y;z ) ue10(x;y;z )ue9(x;y;z ) =ua7(x;y;z )ua8(x;y;z ) ue8(x;y;z )ue4(x;y;z ) =ua8(x;y;z )ua4(x;y;z ) ue4(x;y;z )ue2(x;y;z ) =ua8(x;y;z )ua7(x;y;z ) ue2(x;y;z )ue6(x;y;z ) =ua4(x;y;z )ua8(x;y;z ) 3.For nodes on the corner of the cell: ua6 x=ua7 x;ua3 x=ua2 x;ua7 x=ua3 x;ua1 x=ua4 x; ua3 y=ua4 y;ua1 y=ua2 y;ua5 y=ua1 y;ua6 y=ua5 y; ua1 z=ua4 z;ua4 z=ua2 z;ua7 z=ua6 z;ua3 z=ua2 z; ua8 x;y;z= 0;ua4 x;y= 0; 34 Appendix B. Average triaxialities and Lode values in each half cell of the voided bi-crystal Table B.4: Average local Triaxiality and Lode parameter for each half cell of hard-hard voided bi-crystal Applied global TandLLocalTandL(half cell 1) LocalTandL(half cell 2) T= 0,L=1 T= 0,L=1 T= 0,L=0:8 T= 0,L= 0 T= 0,L= 0 T= 0,L= 0 T= 0,L= 1 T=0:2,L= 1 T= 0:23,L= 1 T= 0:33,L=1T= 0:33,L=1 T= 0:33,L=1 T= 0:33,L= 0 T= 0:33,L= 0 T= 0:3,L= 0 T= 0:33,L= 1 T= 0:22,L= 1 T= 0:4,L= 1 T= 0:66,L=1T= 0:66,L=1 T= 0:66,L=1 T= 0:66,L= 0 T= 0:66,L= 0 T= 0:66,L= 0 T= 0:66,L= 1 T= 0:53,L= 1 T= 0:72,L= 1 T= 1,L=1 T= 1,L=1 T= 1,L=0:86 T= 1,L= 0 T= 1,L= 0:1 T= 1,L=0:1 T= 1,L= 1 T= 1,L= 1 T= 1,L= 1 Table B.5: Average local Triaxiality and Lode parameter for each half cell of soft-hard voided bi-crystal Applied global TandLLocalTandL(half cell 1) LocalTandL(half cell 2) T= 0,L=1 T= 0,L=1, T= 0,L=0:8 T= 0,L= 0 T= 0,L= 0 T= 0,L= 0 T= 0,L= 1 T= 0,L= 1 T= 0,L= 1 T= 0:33,L=1T= 0:64,L=1 T= 0:14,L=1 T= 0:33,L= 0 T= 0:6,L= 0 T= 0:3,L= 0 T= 0:33,L= 1 T= 0:48,L= 1 T= 0:22,L= 1 T= 0:66,L=1T= 0:85,L=1, T= 0:5,L=1 T= 0:66,L= 0 T= 1,L= 0 T= 0:4,L= 0 T= 0:66,L= 1 T= 0:86,L= 1 T= 0:51,L= 1 T= 1,L=1; T= 1:43,L=1, T= 0:65,L=0:86 T= 1,L= 0 T= 1:61,L= 0:1, T= 0:62,L=0:1 T= 1,L= 1 T= 1:3,L= 1 T= 0:78,L= 1 35 Table B.6: Average local Triaxiality and Lode parameter for each half cell of soft-soft voided bi-crystal Applied global TandLLocalTandL(half cell 1) LocalTandL(half cell 2) T= 0,L=1 T= 0,L=1 T= 0,L=0:8 T= 0,L= 0 T= 0,L= 0 T= 0,L= 0 T= 0,L= 1 T=0:2,L= 1 T= 0:23,L= 1 T= 0:33,L=1T= 0:33,L=1 T= 0:33,L=1 T= 0:33,L= 0 T= 0:33,L= 0 T= 0:3,L= 0 T= 0:33,L= 1 T= 0:22,L= 1 T= 0:4,L= 1 T= 0:66,L=1T= 0:66,L=1 T= 0:66,L=1 T= 0:66,L= 0 T= 0:66,L= 0 T= 0:66,L= 0 T= 0:66,L= 1 T= 0:66,L= 1 T= 0:66,L= 1 T= 1,L=1 T= 1,L=1 T= 1,L=0:86 T= 1,L= 0 T= 1,L=0:11 T= 1,L= 0:11 T= 1,L= 1 T= 1,L= 1 T= 1,L= 1