**1. Introduction**

The multiscale homogenization procedure for microstructural material is an important topic in the science and engineering of materials. According to Nguyen et al. [1], computer modeling is the most efficient technique for this purpose, and it has been a subject of interest in recent years. One of the main challenges for multiscale modeling is in developing numerical models that predict the behavior of microstructures with high accuracy and efficiency. Many works on this subject have been developed over the years [2–5]. According to Pundale et al. [6], microstructural materials are considered to be composites for the numerical modeling approach, where the challenge relies on coupling the individual behavior of each constituent to obtain an overall response.

Achenbach and Zhu [7,8] carried out numerical studies with multiphase media using a unit cell model to study the variation of the interface parameters and to estimate the effect of distribution of stresses in both the matrix and in the fibers. Moorthy and Ghosh [9] developed a finite element (FE) model in conjunction with Voronoi cells to examine the small deformation of arbitrary heterogeneous

material in a 2D model. The influences of the shape, size, orientation, and distribution of inclusions on the micro- and macroscopic responses were also investigated. The results indicated that the orientation of every grain did not play as important a role as the size of the inclusions in the response of the overall properties.

Effective properties of macrostructural materials using synthetic geometries (i.e., regular geometries such as cylinders, spheres, and ellipsoids) of different sizes and different mechanical properties were studied by Yao et al. [10], who used a 2D formulation in the boundary element method (BEM) to demonstrate that the BEM formulation may be more suitable for interface analysis when compared to the FE method. Using the following two approaches, Zheng et al. [11] studied a solid whose macrostructure comprised fluid-filled pores: the superposition method and a multi-subdomain method 2D BEM model using subregions. It was shown that the multi-subdomain method (with subregions) was more efficient and accurate for determining the effective properties. Gitman et al. [12] discussed the periodicity of the material, and they established the no-wall-effect concept, which refers to the inability of inclusions to penetrate through the sample borders within a representative volume element (RVE). According to Gitman [12], the periodicity of a material considers that a material experiencing no-wall effects due to the RVE is a representative volume and should, therefore, represent any part of the material. Despite the many experimental and computational methods devoted to the study of the effective elastic modules of nodular cast iron (NCI) [6,13–15], it is important to highlight that computational methods are much more attractive from the point of view of saving time and reducing the operational costs.

Numerical homogenization was used based on the concept of the RVE for determining the effective properties. Through this method, the statistical data of many RVEs can be obtained, and by averaging the values of each distribution of spherical inclusions, it is possible to find the global response of the material [16]. According to Zohdi [16], this technique is more trustworthy than performing individual direct simulations. Buroni and Marczak [13] developed a homogenization procedure using a 2D BEM as a numerical approach to determine the effective modulus of NCI. In this work, the inclusions were considered to be voids, which were discretized with a single special hole element. This new numerical approach demonstrated good accuracy and low computational cost. Carazo et al. [17] performed an RVE study to determine the effective properties of NCI. The authors imposed rectangular and hexagonal geometric shapes of the RVE; the numerical results were compared to analytical expressions, and the authors concluded that the graphite fraction has the largest influence on the effective modulus. They realized that the effective properties are independent of the RVE shape. Fernandino et al. [15] recently used the FE method and a multiscale analysis to determine the effective elastic proprieties of NCI where a complete experimental characterization was introduced and used to evaluate the numerical results.

In this paper, the main objective is to develop a homogenization process for NCI. Nodular cast iron, also known as spheroidal graphite cast iron (SGI), was chosen to develop a numerical model of homogenization. The main reason is the morphology of the microstructure, which contains nearly spherical inclusions of graphite embedded in a homogeneous ferritic matrix. The nodular graphite has a high carbon content, while the matrix medium has a ferritic, pearlitic, ferritic-pearlitic, or austenitic structure [18]. For the numerical homogenization process, two approaches were employed: the first one models the nodules as a synthetic geometry, and the second one models them as their real shapes obtained via microtomography. A numerical routine using the BEM for the homogenization analysis was supported by experimental studies, such as microtomography by X-rays (micro-CT), tensile tests, hardness tests, microhardness tests, micrographic analyses, and scanning electron microscopy equipped with energy-dispersive X-ray (SEM-EDX) analyses. The main goal of this investigation is to determine the influence on the effective properties obtained when considering the nodules as synthetic and real geometries in the homogenization process.

#### **2. Numerical Model**

A multiscale numerical routine was written based on the BEM as the numerical approach for determining the effective Young's modulus of heterogeneous materials with linear elastic behavior. The graphite was modeled as a second constituent by using the subregion methodology provided by BEM. The subregion formulation takes into account the stiffness of nodular inclusions within the ferritic homogeneous matrix of NCI. To accelerate the convergence of the homogenization process, periodic boundary conditions were also implemented in the numerical routine, particularized to plane stress (PT).

#### *2.1. The BEM and Subregions*

The BEM is essentially based on the domain's boundary, which is partitioned in a coarse mesh where the differential governing equation is numerically integrated [19,20]. This numerical approach has all the unknown variables only on the boundary, meaning that no information from inside the domain is required. Figure 1 illustrates a boundary value problem considering both an RVE by sub-region and a no-wall-effect concept [12]. The ferritic matrix is denoted by the subscript 1, and the graphitic nodules are indicated by the subscript 2.

**Figure 1.** BEM using subregions (or multi-domains): (**a**) RVE for NCI and (**b**) RVE for a unitary cell.

The governing equation for a linear elastic material is presented in Equation (1):

$$u(\mathbf{x}) = \int\_{\Gamma\_i} \mathcal{U}\_{ij}(\mathbf{x}, \mathbf{y}) t(\mathbf{y}) t(\mathbf{y}) ds\_y - \int\_{\Gamma\_i} T\_{ij}(\mathbf{x}, \mathbf{y}) u(\mathbf{y}) ds\_y \tag{1}$$

where *U* and *T* are the components of displacement and tractions containing the fundamental solutions for the nodal displacement (*u*) and tractions (*t*), respectively. The fundamental solutions for a 2D elastic model are as follows:

$$\mathcal{U}l\_{\vec{l}\vec{j}}(x,y) = -\frac{1}{8\pi\mu} \left[ (3-4\nu)\log\frac{1}{r}\delta\_{\vec{l}\vec{j}} + r\_{\vec{l}}r\_{\vec{j}} \right] \tag{2}$$

$$T\_{\vec{\eta}}(\mathbf{x}, \mathbf{y}) = -\frac{1}{4\pi(1-\nu)r} \left\{ \frac{\partial r}{\partial n} \left[ (1-2\nu)\delta\_{\vec{\eta}} + 2r\_{\vec{r}}r\_{\vec{\eta}} \right] + (1-2\nu) \left( n\_{\vec{r}}r\_{\vec{\eta}} - n\_{\vec{r}}r\_{\vec{\eta}} \right) \right\} \tag{3}$$

where *μ* and *ν* are the shear stress and Poisson's ratio, respectively; *r* is the distance between the source point and a field point; and *δ* stands for the Kronecker delta. Once the external boundary conditions are imposed, the system is integrated over the total number of nodes (*Tn*), resulting in a linear system of equations:

$$\sum\_{i=1}^{T\_n} \left[ H\_{\vec{i}\vec{j}} u\_i \right] \;=\sum\_{i=1}^{T\_n} \left[ \mathcal{G}\_{\vec{i}\vec{j}} t\_{\vec{i}} \right] \to \left[ H \right](\mu) \;=\,\left[ \mathcal{G} \right](t) \tag{4}$$

where *H* and *G* are the arrays of displacement and tractions, respectively. Every subregion is discretized by discontinuous quadratic elements. Therefore, the BEM equations for each subregion were implemented according to Brebbia and Dominguez [19] and Aliabadi [20].

For the homogenization process to succeed, the influence of the boundary conditions imposed on the problem must be taken into account in order to increase the accuracy and convergence. In this sense, three classical methods may be used to impose the boundary conditions [21]: (1) a linear displacement boundary condition (Dirichlet condition); (2) a constant traction boundary condition (Neumann condition); and (3) a periodic boundary condition (PBC). Some works [21–23] emphasize the importance of imposing the PBC due to its efficiency and accuracy for the final results of analyses in microstructural material. In this work, a PBC that considers the Dirichlet condition was implemented to ensure good accuracy in the homogenization procedure.

#### *2.2. Constitutive Law*

The constitutive linear equation of an elastic material is written in index notation as:

$$
\sigma\_{ij} = \mathbb{C}\_{ijkm} \varepsilon\_{km} \tag{5}
$$

where *Cijkm* represents the components of the fourth-order constitutive tensor, which relates the stress and strain tensors, and it is related to the mechanical properties according to Equation (6):

$$\mathbf{C}\_{\text{ijkm}} = \frac{E}{1 - \nu^2} \left[ (\mathbf{1} - \nu) \mathbf{I}\_{i\bar{j}\bar{k}m} + \nu \left( \mathbf{I}\_{i\bar{j}} \times \mathbf{I}\_{i\bar{j}} \right) \right] \tag{6}$$

where *E* and *ν* denote the Young's modulus and Poisson's ratio, respectively. In our problem using subregions, these parameters are given by:

$$E = \begin{cases} \begin{array}{c} E\_{matrix} \quad \text{and} \quad \nu\_{matrix} \quad \text{if} \\\ E\_{module} \quad \text{and} \quad \nu\_{module} \quad \text{if} \end{array} \quad \begin{array}{c} (\mathbf{x}, \mathbf{y}) \in \Omega\_{\text{matrix}} \\\ (\mathbf{x}, \mathbf{y}) \in \Omega\_{\text{module}} \end{array} \tag{7}$$

where *Ematrix* and *νmatri*<sup>x</sup> are the Young's modulus and Poisson's ratio, respectively, for the ferritic matrix, while *Enodule* and *νnodule* are the Young's modulus and Poisson's ratio for the graphite nodule.

#### *2.3. The Plane Stress (PT) and Plane Strain (PS) Hypotheses and the Mean Field Theory*

The PT and PS hypotheses are two considerations that should be taken into account to decrease the computational costs as well as to minimize the number of variables for a 2D analysis. The macroscopic properties of a heterogeneous material are typically determined using Equation (8), where the microand macro-couplings are implicitly defined. Equations (9) and (10) establish the relationship between the infinitesimal macroscopic stress and strain tensors and the averages of the stress and strain in the RVE [24]:

$$
\langle \langle \sigma\_{\rm ij} \rangle\_{\Gamma} = \langle \mathcal{L}\_{\rm ijkm}^\* \langle \varepsilon\_{km} \rangle\_{\Gamma} \tag{8}
$$

$$
\varepsilon\_{ij}^{\mathcal{M}} = \langle \varepsilon\_{ij} \rangle\_{\Gamma} \tag{9}
$$

$$
\sigma\_{i\dot{j}}^{\mathcal{M}} = \langle \sigma\_{i\dot{j}} \rangle\_{\Gamma} \tag{10}
$$

where ‹› stands for the average stress taken over the boundary Γ. It is the spatial average operator, which is defined as follows:

$$
\langle \bullet \rangle \, \, = \, \frac{1}{|\Gamma|} \int\_{\Gamma} \bullet d\Gamma \, \tag{11}
$$

when • is the mechanical property to be evaluated. According to Zohdi et al. [25], a specimen can be computationally tested to permit the solution to the elastic problem for determining the stress and strain micro fields. Solving the system of equations shown in Equation (8) determines the constants for the matrix of elasticity. The microstructural sample material is subjected to the same boundary conditions, resulting in a field of tensions and deformations that are distributed uniformly as an entire homogeneous body [13,24].

#### **3. Methodology and Materials**

In this work, the graphite nodules were modeled numerically using two different procedures. The first one employs a statistical analysis, assuming that nodular graphite has a synthetic geometry (circular shape), and the second procedure assumes that the nodular graphite has the real geometry obtained by X-ray microtomography. The morphological characterization and mechanical properties were determined through experimental tests, which were used to support and validate the proposed numerical homogenization.

#### *3.1. X-Ray Microtomography*

To validate the proposed methodology, digital images were acquired using micro-CT. The samples were scanned using a ZEISS X radial 510 Versa, high-resolution, micro-CT scanner system (Versa 510, Bregnerødvej, Birkerød, Denmark). This system generates X-rays with a cone-beam geometry, and a 3D digital object can be visualized from the acquired radiograms. The small samples of NCI have dimensions of 6 mm × 7 mm × 24 mm. Each sample was then mounted on a rotating cylindrical support inside the X-ray chamber and positioned so that the X-ray lines were perpendicular to their axis of rotation. Then, the scanning parameters were adjusted and scanned for 15 s. The geometry obtained by micro-CT was straightforwardly exported to generate a BEM mesh, making it possible to use the real geometry of the microstructure.

Through the micro-CT images, it was also possible to determine the shape coefficient *Kf* (*Kf* = *D*1/*D*2) and the characteristic parameter *r*/*d* (*r* = (*D*<sup>1</sup> + *D*2)/4). Here, *D*<sup>1</sup> and *D*<sup>2</sup> are the smallest and the largest diameter, respectively, measured inside the graphite nodule; *r* is the average radio; and *d* is the average distance measured between two graphite nodules [26,27].

To obtain an approximate representation of the characteristic parameters of the NCI microstructure, a statistical analysis was implemented by random selection of graphite nodules over the analysis zone. All data were tabulated in histograms.

#### *3.2. Microstructure of NCI*

A micrographic analysis was used to characterize the spatial distribution of graphite nodules and the average nodule size. Micrograms were obtained from NCI blocks with dimensions of 5 mm × 5 mm × 5 mm. The surface was polished and etched with nital at 3%, and the samples were then mounted on a sample holder and analyzed in a vacuum according to ASTM E112-96 standard. Three micrograms were obtained with scales of 400, 100, and 20 μm, and the resulting images were used to measure the amount, area, and overall size of the graphite nodules.

#### *3.3. Hardness Tests*

Vickers and Brinell hardness testing was performed to measure the hardness of both the ferritic matrix and the graphite nodules in line with ASTM E92-16 standard. By using the micro-indentation methodology [28], it was possible to determine the Young's modulus for both the matrix and the graphite nodules. In this work, the measurements were carried out with an FM 700 Vickers microhardness apparatus (FM 700 Vickers, Karola Miarki 12, Gliwice, Polska), equipped with a glass camera with a 100× magnification capability and which had a load of 10 g. The average hardness of the ferritic bulk was determined using a Brinell hardness (HB) test, in line with ASTM E92-16, with a semiautomatic HB tester (ZHU250 HB, Kennesaw, GA, USA). Afterwards, the tensile strength of NCI was determined through the standard correlations between the HB hardness and tensile strength. Figure 2a,b depict the shallow indentations made in the microhardness tests, revealing that the graphite nodule is softer than the matrix due to the difference in highlights between the areas of the hardness

impression. The HB test is presented in Figure 2c. In this case, the spherical indentations are larger than the Vickers hardness HV indention, meaning that the hardness value only indicates the toughness of the bulk material.

**Figure 2.** Hardness tests: (**a**) HV matrix indentation; (**b**) HV graphite nodular indentation; and (**c**) HB indentation.

The mathematical model discussed in [28] consists of identifying the geometric parameters of the indented surface. Then, from the geometrical parameters of the hardness impressions and the mechanical and geometrical properties of the indenter, it is possible to determine the Young's modulus according to Equation (12):

$$E = \frac{1 - \upsilon^2}{\frac{1}{E\_r} - \frac{1 - \upsilon\_i^2}{E\_i}} \tag{12}$$

where *ν<sup>i</sup>* is the Poisson's ratio of the indenter, *ν* is the Poisson's ratio of the sample, *Ei* is the Young' modulus of the indenter, *Er* is the reduced Young's modulus of the sample, and *E* is the Young's modulus of the sample. The parameters used for estimating the effective Young's modulus are listed in Table 1.

**Table 1.** Parameters used for estimating the experimental Young's modulus.


The values of *Er* were extracted by using Equation (13):

$$E\_r = \frac{\varpi P\_{\text{max}} \sqrt{\pi}}{2 \text{h}\_\text{s} \sqrt{A}} \tag{13}$$

where  is a constant that depends on the geometry of the indenter; *Pmax* is the maximum load; *hs* is the amount of sink-in, which is assumed to be approximately equal to the final depth *hf*; and *A* is the projected area of the elastic contact. Due to the uncertainty in the Vickers system and the structurally inhomogeneous sample, a statistical test was developed using micro-indentation patterns containing 50 indentations over the analysis zone.

#### *3.4. Tensile Tests*

Tensile testing was performed on sub-size test samples, which were obtained by machining the NCI blocks according to the specifications in ASTM E8/E8M-09 standard, using an Instron 8801 test machine (Instron, Glenview, IL, USA) [±100 kN (22,500 lbf)] with a 2620-601 dynamic gauge extensometer (Instron, Glenview, IL, USA). The Young's modulus, tensile strength, and ultimate tensile strength were determined from the experimental stress-strain curve.

#### *3.5. Scanning Electron Microscope Equipped with EDX*

The quantitative compositional information for NCI was determined by using a scanning electron microscope system equipped with energy dispersive X-ray analyzer system (JEOL JSM 6610 SEM-EDX, Musashino, Tokyo, Japan). The samples were polished, treated, and placed in a vacuum chamber for visualization of the surface. The samples were analyzed at 100×, and the microstructural examinations were carried out for a representative sample with a large number of graphite inclusions.

#### *3.6. Computational Homogenization*

The homogenization concept is based on a multiscale analysis that establishes the relationships between macro- and micro-scale entities [29]. The micromechanical analysis typically focuses on the concepts of RVEs; the analysis involves an average representative volume extracted from the micro-heterogeneity in an infinitesimal-volume sample [1,12,24]. This sample contains a reasonable representation of the microstructure's material; therefore, the material properties are obtained through average volumes of the respective micro fields obtained from an RVE [12,24]. The numerical homogenization was developed following two procedures.

#### 3.6.1. Case I: Graphite Nodules Considered as Synthetic Geometry

The synthetic homogenization procedure was developed using the characteristic parameters obtained from statistical analyses based on morphological characterizations. There is no definitive way in which to establish an RVE [16].; the most frequently used procedure currently involves either fixing the nominal RVE size and then increasing the number of graphite nodules (Figure 3) or increasing the size of the RVE so that the number of heterogeneities changes by a fixed percentage.

**Figure 3.** The RVE with different sizes of graphite nodules.

In Figure 3, RVE1 ⊂ RVE2 ⊂ RVE3 ⊂ RVE4. The discretization strategy of the microstructure was assumed to be a square sheet for the 2D RVE case and a cube for the 3D RVE case. In both cases, a constants dimension, containing a specific number of nodules, was defined. The characteristic parameter obtained with the procedure described in Section 3.1 was used to model NCI. The nodules/matrix ratio of the sample was obtained by morphological characterization using image analysis during the binarization process, which establishes the percentage of each phase. The homogenization model employed samples with the following numbers of nodules: 4, 6, 8, 10, 15, 20, 25, 30, 40, 45, 50, 55, and 60.

#### 3.6.2. Case II: Graphite Nodules Considered as Real Geometry

Figure 4 presents the homogenization procedure used to estimate the effective Young's modulus considering a real morphology. This procedure involves the generation of several RVEs with different sizes and in random coordinates over the sample (Figure 4a) [12]. The RVE sizes were obtained from numerical tests, where constant displacements were imposed as a boundary condition. Furthermore, the micrograms obtained by micro-CT were used to obtain the forms and distribution of the nodules (Figure 4b). In the sequence, a 2D mesh that represents the inclusions of the irregular geometry inside the ferritic matrix was generated (Figure 4c).

**Figure 4.** The RVE for NCI: (**a**) sample, (**b**) random RVE, and (**c**) 2D mesh.

The boundary was discretized with a discontinuous quadratic element. Moreover, the dimensions of the RVE were progressively increased until the effective Young's modulus was reached.

#### **4. Results and Discussion**

The results obtained from the experimental and numerical tests described in the previous section are presented below.

#### *4.1. X-Ray Microtomography*

According to the micro-CT images, the graphite nodules exhibited a shape that is approximately a circle. In this sense, the shape parameter (*Kf*) and the characteristic parameter (*r/d*) were estimated through statistical testing. Figure 5 illustrates the nodule distributions in an NCI sample. The results for *Kf* and the *r/d* ratio are presented in a histogram graph.

**Figure 5.** Geometric parameters for NCI: (**a**) histograms of *Kf* and (**b**) *r*/*d* ratio.

According to Figure 5a, it is possible to verify that the graphite nodules present a dominant circular shape, as most of the samples present *Kf* values between 0.75 and 1.0. However, some inclusions also present a vermicular shape. Due to this accentuated characteristic, in this work, the homogenization procedure using synthetic geometry will be developed considering the graphite as a circular shape. The *r*/*d* parameter was 0.3 ± 0.01, as indicated in Figure 5b, and the error was sufficiently low for the model representation of NCI. These parameters were coded into the numerical model to represent NCI.

#### *4.2. Microstructure of NCI*

The micrograms obtained from the laser confocal system were used to determine the NCI grain size. Figure 6 presents the microstructure of a sample of NCI after etching with a 3% nital solution. In Figure 6a, the ferritic-pearlitic structure of the sample is observed. This structure is surrounding the nodular graphite, forming a so-called "bullseye," and it contains certain alloying elements, such as magnesium and silicon, and small quantities of other elements that do not act during the nucleation process. The presence of pearlite lamellae is due to the high cooling rate of the iron from the eutectoid temperature. Figure 6b depicts the matrix structure with the "bullseyes" and the mixture of fine lamellar pearlite and a ferrite matrix structure. Finally, Figure 6c presents the pearlite as a thin lamellar structure, which is the result of a high cooling rate after casting. Graphite is surrounded by perlite, with the presence of some ferrite [30].

**Figure 6.** Metallographic structure of NCI: (**a**) ferritic-pearlitic structure, (**b**) "bullseye" structure, and (**c**) lamellar pearlite structure.

Table 2 summarizes the result for NCI that corresponds to the nodule diameter for a grain size number (G) of 4.62 (See ASTM E112-96). The average grain size of NCI obtained from the statistical count is illustrated in Figure 7.


**Table 2.** Experimental results from the laser micrographs.

**Figure 7.** Grain size number (G) for GGG-40 obtained by the 3D laser microscopic system.

From Figure 7, the nodular graphite size indicates a value of G of around 5. In addition, ASTM E112-96 standard presents a relationship between G and the average diameter. For this study, the average diameter is approximately 64 ± 10 μm, which is in good agreement with related works [18,31,32].

#### *4.3. Hardness Tests*

The effective Young's modulus was determined from the Vickers micro-hardness tests, with the associated results presented via the histograms in Figure 8a,b for the Young's modulus of the matrix and the nodules, respectively.

**Figure 8.** Histogram depicting the results for the Young's modulus of (**a**) the matrix and (**b**) the nodules.

The mean and standard deviations for the matrix and nodules are 241 ± 19 GPa and 29.17 ± 1.1 GPa, respectively. There is a significant difference between the Young's modulus of the matrix and that of the graphite nodules; the main reason for this is the low strength of the graphite. To ensure that the present results are valid, some studies regarding the determination of the Young's modulus for NCI are introduced and presented in Table 3.

#### *4.4. Tensile Tests*

Tensile tests were performed, and the stress versus strain curve is presented in Figure 9. The mechanical properties, including the yield stress (*σ*0.2), ultimate tensile strength (*σu*), Young's modulus (*E*), and Poisson's ratio (*ν*), were obtained from the experimental curve, and they are reported in Table 3. The Young's modulus was 171 GPa, averaged over the four tests. The values of the modulus of elasticity that others obtained for NCI are 187 GPa [33], 172 GPa [15], and 179 GPa [34], demonstrating good agreement with the present values.

**Figure 9.** Stress-strain behavior of NCI determined from a tensile test.


**Table 3.** Experimental results from the tensile and micro-hardness testing.

#### *4.5. Scanning Electron Microscope Equipped with EDX*

During the formation of the nodular graphite in the casting process, additional elements, such as silicon and magnesium, are required to alter the solidification mechanism and chromium to improve the tensile strength. The chemical composition is presented in Table 4.


**Table 4.** The chemical composition of NCI (wt %).

Nodular cast iron is formed by approximately 3.6% carbon, and it makes up 92.8% per mass of iron. The remaining chemical compounds were measured by SEM-EDX equipment (JEOL JSM 6610 SEM-EDX, Musashino, Tokyo, Japan). Despite obtaining information regarding all percentages of constituents (Table 4), it is important to note that this type of equipment is not able to detect small percentages of chemical constituents with good accuracy. For those constituents with high percentages (Fe, Si, and C), the values measured through the use of SEM-EDX equipment demonstrated good agreement with those provided by the foundry industry from which these samples were acquired. The high amounts of carbon are responsible for the formation of nodular graphite and for the lamellar perlite produced by the carbon remainder that did not form graphite during nucleation. According to Table 4, the high percentage of Si joined to the Mg explains the nodular graphite's round shape in the ferritic matrix. The small amounts of Cr and Cu account for the formations of the pearlite structure and the high tensile strength.

#### *4.6. Computational Homogenization*

#### 4.6.1. Case I: Graphite Nodules Considered as Synthetic Geometry

#### Two-Dimensional Model

Figure 10 illustrates the convergence curve for *E*\*, where a wide dispersion for each computational test is evidenced. According to this graph, when a few nodules are considered, *E*\* varies considerably, indicating a positional dependence, which affects the effective properties significantly. It is possible to observe that as the number of nodules increases, *E*\* stabilizes for each set of samples. The elasticity modules, as well as the statistical data, are presented in Table 5.

In Table 5, *E\** is the average value for each sample, (*Emax* − *Emin*) is the difference between the maximum and minimum values, (SD) is the standard deviation, (DVM) is the percentage deviation of the difference between maximum and minimum values with respect to the average, and (DM) is the percentage deviation of the standard deviation from the average value. Based on Table 5, it is possible to realize that *E\** is reached with 35 nodules, and after this increment, the Young's modulus no longer varies significantly (less than 5%) except for particular points that affect the general result generated, as can be seen in Figure 11a. Figure 11a depicts the results for the *E*\* microstructures containing 35 nodules in each analysis; the maximum, minimum, average, and standard deviations of the group of samples tested are shown. Figure 11b presents the corresponding histograms of the results, where most of the analyses yielded values for *E\** of 181 GPa. The value of *E*\* for the NCI obtained from the experimental data is 171 GPa.

**Figure 10.** Convergence curve of the RVE. Values of the effective Young's modulus for samples with different random distributions of nodules.

**Table 5.** *E*\* for computational tests of the NCI microstructure with various random distributions of nodules.


**Figure 11.** Statistical studies for average value: (**a**) *E*\* for 20 microstructures with 35 randomly distributed graphite nodules, and (**b**) the histogram of *E*\*.

The slight difference between the experimental and numerical results is approximately 8%. This difference is perhaps due to the 2D microstructural model employed and its simplification. According to Rodrígues et al. [14], a difference exists in the effective properties when using 2D and 3D RVEs. The authors claim that hydrostatic effects on the deformation of a 2D nodule are different from those effects that occur in a real 3D solid.

#### Three-Dimensional Model

In this work, the same strategy from [13] and [16] was used to study the effective property in an NCI RVE. Figure 12 presents the convergence curve of *E*\* for a 3D model. Furthermore, Table 6 lists the (*Emax* − *Emin*), the standard deviation, the percentage change in the difference between the maximum and minimum values with respect to sock (DVM) and the mean values of each of the samples for *E*\*, and the difference between the maximum and minimum values and the standard deviation (DM).

**Figure 12.** Convergence curve of the RVE. Values of the effective Young's modulus for samples with different random distributions of nodules.

**Table 6.** The *E*\* for computational test of the NCI microstructure with various random distribution of nodules.


It can be observed that an *E*\* is reached for 30 inclusions, where there are increases of less than 5%. Figure 13a presents the result for *E*\* for an RVE with 30 random distributions of inclusions, while Figure 13b presents the corresponding histograms of the results, where it can be appreciated that the data for *E*\* were 186 GPa.

**Figure 13.** Statistical studies for average values: (**a**) *E*\* of material for 50 microstructures with 30 randomly distributed graphite nodules, and (**b**) the histogram of *E*\*.

Now, for the 3D analysis, it can be seen that the Young's modulus was 186 GPa due to the hydrostatic effects that generate a higher stiffness in the material.

#### 4.6.2. Case II: Graphite Nodules Considered as Real Geometry

The homogenization was performed through the use of the real geometry of the NCI microstructure obtained by micro-CT acquisition. The results of the homogenization process are presented in Figure 14. It can be observed that as the number of inclusions increases, the effective Young's modulus is stabilized, reaching values of 175 GPa for an RVE with 30 inclusions. The RVE size in this case was 1.5 mm2. The present methodology demonstrated good accuracy for determining the effective Young's modulus in relation to the values obtained with the experimental results, where the difference is approximately 8%. The error in the Young's modulus value can be justified by the 2D model approach. Despite considering the real geometry of the graphite, a difference of 8% remains between the values of the effective Young's modulus obtained experimentally and numerically. In this sense, perhaps a 3D model that takes into account the real geometry should be more realistic, thereby minimizing this difference. Table 7 summarizes the results for the effective Young's modulus obtained through homogenization (real and synthetic geometries) and the experimental procedure and those available in the literature.

**Figure 14.** *E*\* for NCI.

In this work, three numerical models for determining the effective Young's modulus were considered. Table 7 summarizes the results for the effective Young's modulus obtained for each numerical model, the Young's modulus obtained by tensile test. Considering the approach of the graphite nodule as a synthetic geometry, both 2D and 3D models resulted in values greater than those obtained experimentally and those available in the literature. However, the 2D model that considers the real geometry of a graphite nodule resulted in good agreement with the values obtained experimentally and with [33]. In this sense, it may be possible to claim that those models that consider the real geometry of graphite nodules are more accurate for determining the effective properties.



#### **5. Conclusions**

In this paper, a homogenization process for determining the effective Young's modulus of NCI was performed. A numerical strategy for the homogenization process was written for a 2D and a 3D elastic model using all features provided for the BEM. The homogenization process was performed considering the nodular graphite modeled with synthetic and real geometries. Periodic boundary conditions were imposed to increase the convergence to reach an RVE. The real geometry of the microstructure was obtained from micro-CT images, and a special subroutine was developed to generate the BEM mesh. The graphite nodules were always considered as inclusions, and their geometries were modeled as subregions. In this sense, experimental studies were also performed to identify NCI's mechanical properties and morphological structure. From the results of the homogenization process, it is possible to verify that the modeling of the graphite as inclusions instead of voids significantly increases the accuracy. Regarding the graphite nodules, two different approaches were considered. The first one used the real geometry in which the BEM mesh was built based on micro-CT images of NCI. In the other approach, the graphite nodules were considered with a synthetic geometry; that is, the graphite shapes were considered to be perfectly round. The numerical results computed by the two approaches had nearly the same accuracy, whereas the modeling considering the real geometry of the graphite had higher computational efficiency than the model considering the synthetic geometry. The authors are, thus, devoting efforts to 3D BEM modeling considering the real geometry of NCI.

**Author Contributions:** A.B. and C.A. contributed with the numerical modelling and the experimental tests for the material characterization. A.M. and R.L. contributed with micro CT test.

**Acknowledgments:** The authors thank Capes and FAPDF for financial support.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


© 2018 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article*
