Residual Stresses around Individual Graphite Nodules in Ductile Iron: Impact on the Tensile Mechanical Properties

Residual stresses in ferritic ductile iron castings have been studied for decades. However, little attention has traditionally been given to the local residual stresses which may arise in the microstructure as a result of the thermal contraction mismatch between the matrix and the graphite nodules during solid-state cooling. Recent synchrotron X-ray measurements performed by the authors have demonstrated that in the ferritic phase these local stresses can be in the order of 100-150 MPa, hence of the same order of magnitude as the material macroscopic yield stress. This suggests that they might have a significant influence on the mechanical properties of ductile iron components. However, no systematic research appears to have been conducted so far to investigate this aspect. The present work takes a first step in this direction by presenting an integrated theoretical analysis which addresses both the formation of these local residual stresses at the microscopic level and their role during mechanical loading at the macro-scale.


Introduction
Ductile iron (DI) is nowadays widely used in key industrial sectors like transport, marine, and energy production, accounting for as much as 25 % of the total casting production worldwide [1]. One of the main reasons for the great success of this material is its microstructure, which is formed by spheroidal graphite particles (nodules) embedded in a continuous metallic matrix. This composite microstructure provides a combination of strength, ductility and toughness which is quite unique within the cast iron family, and, together with the low price and excellent castability of the material, makes DI ideal for several advanced engineering applications [2].
As with many other cast alloys, DI parts may contain macroscopic residual stresses which arise due to the presence of mechanical constraints that hinder the thermal contraction of the material during solidification and solid-state cooling. This problem is well-known, and today reliable tools exist that allow predicting and coping with it in almost all practical cases. On the other hand, the intrinsic composite nature of DI suggests that another type of residual stress may form, at much smaller length scales, as a consequence of the thermal contraction mismatch between the nodules and the matrix. As thoroughly discussed in [3] [4], this issue has not received much consideration in the past, probably because of two main reasons. First, direct measurements of the local stress state around particles which normally do not exceed approximately 150 μm in size were not possible until recently. Second, any reliable theoretical estimate was also prevented by the very large uncertainty existing on the mechanical properties of the nodules, which in some reports were considered as voids while in others were even regarded as rigid spheres.
To shed light on this controversial topic, a combined set of experimental and theoretical analyses have been performed by the present authors in recent times. In detail, a novel synchrotron X-ray diffraction method has been used to characterize a ferrite grain neighboring a nodule located below the surface of a thin sample, see Fig. 1 (a), which was extracted from a DI wind turbine main shaft [5]. The measurements have revealed the existence of both plastic deformation and a complex residual elastic strain field close to the nodule, as Fig. 1 (b) shows. In parallel, a theoretical model has been developed to provide a quantitative description of the thermo-elastic behavior of the nodules [6]. This model has been used to explain the X-ray findings and uncover the important effect of the conical sectors forming the internal structure of the nodules, which drive residual stress concentrations in the surrounding ferrite regions [7].
As can be seen in Fig. 1 (c), taken from [7], the maximum magnitude of the local residual stresses, associated with the elastic strain measured via X-ray diffraction, is likely in the order of the macroscopic proof stress of DI. This suggests that these stresses might have a significant influence on the mechanical properties of the material. The aim of the present paper is to provide a preliminary investigation of this aspect, by presenting an integrated theoretical analysis which addresses the role played by the local residual stresses on the standard tensile properties of DI.
(a) (b) (c) Figure 1. (a) Setup of the synchrotron X-ray diffraction measurements performed in [5] and (b) corresponding residual elastic strain in one of the ferrite grains neighboring the nodule, reported as a function of the z-distance from the nodule-matrix interface. Note that the strain refers to the Fdirection, which is perpendicular to the diffraction planes sketched in (a). (c) Von Mises stress in the matrix as a function of the radial distance from the nodule, as predicted by the model developed in [7], which was validated against the X-ray diffraction data.

Numerical Model
A two-scale hierarchical approach is adopted to study the role of the local residual stresses upon tensile loading. As Fig. 2 shows, it consists of 3 main steps, which are carried out using two separate models -one at the macro-scale and one at the micro-scale -connected in a sequential manner.
Macro-scale: Solidification and solid-state cooling. First, solidification and solid-state cooling of the entire DI main shaft considered in [5] are simulated at the macro-scale using MAGMA5 [8], a commercial code widely adopted in the foundry industry. The aim is to predict the thermal history at the specific location from which the sample used in the X-ray diffraction analysis was extracted, indicated with a red dashed arrow in Fig. 3 (a). The specific parameters used to setup the simulation can be found in [9], whereas the simulation outcome, in terms of cooling curve predicted at the abovementioned location, is reported in Fig. 3   Micro-scale: Solid-state cooling and tensile loading. The thermo-mechanical interaction between the nodules and the matrix is simulated at microstructural level by means of the 3D periodic unit cell model depicted in Fig. 4 (a). The model is composed of a single spherical nodule embedded in a cube of matrix material, whose side-length is adjusted to provide the 11.5 % graphite volume fraction of the DI under examination. The highly-inhomogeneous internal structure of the nodule is discretized explicitly following the approach described in detail in [6]. Accordingly, the nodule is subdivided into 80 conical sectors radiating from the center of the particle as shown in Fig. 4 (b). Each sector is subdivided further into a bulk region, composed of graphite platelets oriented perpendicularly to the radial direction, and a superficial layer characterized by smaller and highly misoriented graphite crystals.
The model is implemented in the finite element software ABAQUS. Symmetry boundary conditions are applied to the three faces of the cell intersecting the nodule, whereas standard planeremain-plane constraints combined with zero average surface traction requirements are imposed on the other faces. Following [10], a frictionless contact condition is applied between the nodule and the matrix, due to the experimental evidence of the negligible interface strength reported in [11].
The formation of local residual stresses is simulated by prescribing a uniform temperature variation over the entire unit cell according to the curve in Fig. 3 (b). For simplicity, it is assumed that solid-state cooling starts at the upper inter-critical temperature, which defines the stress-free reference state. The eutectoid reaction is modelled as a ferrite halo that nucleates at the nodulematrix interface and then grows in radial direction. In this respect, both the nodule growth due to carbon diffusion and the matrix volume expansion associated with the austenite-to-ferrite transition are accounted for. During the subsequent cooling in the ferritic field, it is assumed that the carbon mass in the nodule remains constant, and the nodule-matrix interaction is governed by the mismatch in thermal contraction only. Further details and motivations for all these assumptions can be found in [7]. Once the temperature of the model has reached 20 °C, mechanical tensile loading begins. This is simulated by prescribing a uniform x-displacement of the cell face whose normal points in the negative direction of the x-axis, while keeping the previously defined boundary conditions active on the other faces. During deformation, the "uniaxial" stress � and strain ̅ are calculated as volume averages of the normal components in the x-direction of the microscopic stress and strain , i.e.
Note that, due to the particular boundary conditions applied, the volume averages of all the other stress components vanish during the period of loading.

Constitutive Models
Matrix. It is assumed that the constitutive behavior of the matrix can be described by the isotropic elastic-visco-plastic model proposed by Peric [12]. The basic equations of the model are: -Additive strain decomposition where , , , ℎ denote the total, elastic, visco-plastic and thermal strain, respectively. -Linear elastic law where is the Cauchy stress and ℂ is the 4 th order isotropic elastic stiffness tensor. -Von Mises-type yield function where ℎ is the so-called stress threshold, is the equivalent visco-plastic strain which obeys the evolution law ̇=, is the 2 nd order identity tensor and tr(⋅) denotes the trace operator. -Associative flow rule: where and are the strain-rate sensitivity and the viscosity-related parameter respectively. To complete the model, the exponential strain hardening law is selected, where 0 , ∞ and * are material parameters. One of the key features of Peric's model is its behavior in the limit ℎ → 0 and → 0. Indeed, in the first case the model reduces to the well-known Norton-Hoff model, which is used to describe power-law creep at high temperature. In the second case instead, the model reduces to the standard

468
Science and Processing of Cast Iron XI rate-independent J2-flow theory of plasticity, which is normally adopted to describe the behavior of metals at room temperature. In such a condition, ℎ defines the yield stress evolution, with 0 and ∞ representing the initial yield stress and the ultimate tensile strength respectively. As a consequence, Peric's model is well suited for the present analysis, where applicability during both solid-state cooling and room temperature loading is required.
Concerning the values of the model parameters, the Young's modulus is assumed to be given by where = −0.81, ∆ is the difference between the actual temperature and 20 °C, is the Young's modulus at 20 °C and is the matrix melting temperature [13]. In the present work, and are taken equal to 205 GPa and 1773 K, respectively. Poisson's ratio and the hardening parameter * are assumed constant instead and equal to 0.29 [13] and 0.1 [7] respectively. The values of the remaining parameters are prescribed according to the temperature-dependent data reported in Fig. 5. Figure 5. Temperature dependence of selected material parameters. (a) Parameters 0 and ∞ of Peric's model [7]. (b) Parameters and of Peric's model [7]. (c) Engineering thermal strain experienced upon cooling from a temperature T to 20 °C of ferrite with 2.0 wt% Si [7], of graphene [14] and of fine grain isotropic graphite IGI110 [14].

(a) (b) (c)
Nodule. The mechanical response of the nodule is assumed to be linear elastic. Following the analysis of [6], the bulk region of each conical sector, formed by graphite platelets stacked on each other, is assumed to have the elastic symmetry characteristic of rhombohedral graphite. The associated elastic moduli, which are valid for a perfect stack of graphene sheets, are reported in Table 1, where the 3-direction denotes the axis perpendicular to the stack. In order to account for the decrease in stiffness due to the presence of defects such as twins and bending planes within the graphite platelets, the elastic moduli of Table 1 are multiplied by a factor 0.15. As demonstrated in [7] and supported by the findings of [15], this method provides a precise prediction of the local residual stresses occurring around the nodule. The superficial layer of the sectors is assumed to be isotropic instead, with Young's modulus and Poisson's ratio equal to 10 GPa and 0.2 respectively, for the reasons discussed in [6]. No temperature dependence is taken into account, as this quantity seems to affect the elastic properties of graphite only to a limited extent within the temperature range considered [16].
Concerning the thermal contraction, the proper values of graphene are employed for the platelets forming the bulk region of the sectors, whereas the values of the isotropic fine-grain graphite IGI-110 are used for the superficial layer. Both sets of data, taken from [17], are plotted in Fig. 5 (c).

Results and Discussion
The tensile stress-strain curve for DI as predicted by the two-scale numerical analysis which takes the local residual stresses into account is shown in Fig. 6. In the same figure, an analogous curve obtained by simulating only the loading stage, i.e. by assuming that the DI microstructure is initially stress-free at room temperature is plotted for comparison. It can be seen that, at low stress levels, the two curves exhibit seemingly linear elastic regions which are almost overlapping. This qualitative impression is confirmed by the very close values of the secant modulus E 0÷100 listed in Fig. 6, which represent the average slope of the curves between 0 and 100 MPa. On the other hand, substantial differences between the two curves are visible in the transition zone from the macroscopically elastic to the macroscopically plastic regime. In particular, when the local residual stresses are neglected, an abrupt yield point is observed, whereas when they are accounted for, a more gradual transition is recorded. In relation to this, the 0.2% proof stress R p0.2 turns out to be approximately 40 MPa higher in the latter case as compared to the former, which might sound surprising, as the residual stresses are often believed to weaken the material by promoting earlier plastic flow in the matrix. Finally, once the transition to the macroscopically plastic regime is completed, both curves display the same hardening behavior for the remaining part of the deformation range analyzed in the present work.
To shed light on the mechanisms controlling the shape of the curves shown in Fig. 6, particularly during the early deformation stages, the formation and subsequent evolution of the plastic region in the unit cell is recorded stepwise up to ̅ = 0.30%. The corresponding graphical representation, reported in Fig. 7, allows two crucial observations. First, plastic flow in the matrix begins earlier when the local residual stresses are included in the simulation. Nevertheless, up to ̅ = 0.12% ( � ≈ 200 MPa) the plastic zone remains limited and localized close to the nodule, meaning that it hardly has any impact on the overall stiffness of the unit cell. Second, a comparison with Fig. 6 reveals that the pronounced stiffness decrease at the macroscopic level, which marks the transition from the elastic to the plastic regime in the stress-strain curves, occurs due to the formation of a continuous plastic band crossing the entire unit cell. The spread of this plastic band is delayed when the residual stresses are accounted for, which explains why the associated stressstrain curve is characterized by a larger value of R p0.2 .
It should be emphasized that the present findings do not allow stating that the local residual stresses are always beneficial to the tensile properties of DI. Indeed, the way plasticity propagates Figure 6. Comparison of the tensile stress-strain curves of DI as predicted by the unit cell model considering local residual stresses (red line with square markers) and as calculated assuming an initial stress-free microstructure at room temperature (blue line with triangular markers).

470
Science and Processing of Cast Iron XI Figure 7. Comparison of the propagation of plasticity in the unit cell upon tensile loading as predicted either assuming an initial stress-free microstructure at room temperature (top sequence) or considering local residual stresses (bottom sequence). The symbol denotes the volume fraction of the matrix which is in the plastic regime.
upon loading depends on the particular geometry of the model. This means that if more realistic representations of the DI microstructure had been used instead of the simple periodic unit cell of Fig. 4, the outcome of the analysis might have been different. However, the present results are of use since 1) they indicate that the local residual stresses may have an impact on the mechanical properties of DI and 2) they justify further investigations into this scarcely investigated topic.

Conclusions
In the present work an integrated modelling approach based on two sequentially connected models -one at the macro-scale and one at the micro-scale -was developed and used to investigate the effects of the local residual stresses during tensile loading of ductile iron. A comparison of simulations performed with and without considering the local residual stresses indicate that, upon loading, plastic flow in the ferrite matrix begins earlier when residual stresses are accounted for. Nevertheless, the subsequent growth of the plastic region seems to be slowed down by the presence of residual stresses, and this leads to an increment of the macroscopic proof stress in the order of 10%. Despite biased by the simple geometry used to describe the material microstructure, the present findings demonstrate that the residual stresses developing at microstructural level may have a significant impact on the mechanical properties of ductile iron.

Acknowledgements
Part of this work has been supported by the Strategic Research Center "REWIND -Knowledge based engineering for improved reliability of critical wind turbine components", Danish Research Council for Strategic Research, grant no. 10-093966. The authors gratefully acknowledge Vestas Wind Systems A/S for providing the geometry and the casting details of the wind turbine main shaft used in the simulations.