Imitative model of the formation of shrinkage cavities and macroporosity

INTRODUCTION

Computational modeling of casting processes is an effective tool for reducing time and resources when starting the production of cast parts. The main task in developing the technological process (achieving dense casting) is solved by predicting the porosity obtained from the modeling.

The physical schematization and the mathematical model of the formation of macroporosity and shrinkage cavities are presented in works [1, 2]. According to this model, the formation of shrinkage cavities is considered a process of liquid metal flow under the influence of pressure differences that arise in the biphasic zone of the casting during solidification. The formation of closed shrinkage cavities and macroporosity occurs as a process of creating new separation surfaces, i.e., due to the rupture of the liquid metal's continuity under tensile stresses that exceed its strength [3].

The main difficulties in modeling the solidification process of casting are associated with determining the continuously changing boundary between the metal and the environment (open shrinkage cavity) and the appearance of new separation surfaces within the metal: macroporosity and closed shrinkage cavities. Existing methods for constructing an arbitrary free surface of a moving liquid, such as front-tracking, adaptive meshes, and the volume-of-fluid (VOF) method, provide high accuracy in determining the curved surface of the liquid metal, but are complicated to implement and require substantial computational resources, hence they are more frequently used for scientific purposes [4].

To solve practical problems in foundry production, an imitative simulation model based on the step-by-step method and the casting's mass balance is widely used, where a single physical process of liquid metal flow under the influence of solidification forces is decomposed into several processes occurring independently of each other. The imitative model is significantly simpler to implement and sufficiently adequate. Many years of experience using this model as part of the "PoligonSoft" computational casting process modeling system [5] show that with the correct choice of model parameters, the prediction of porosity in casting is almost always experimentally confirmed.

Previously, a step-by-step method was proposed [6] to determine the shape of the shrinkage cavity, taking into account the capillary feeding of interdendritic spaces above the free surface of the liquid metal, applicable to the finite difference model of casting. This article is dedicated to the development of a numerical model that simulates the processes of shrinkage porosity formation by considering the capillary forces that arise in the biphasic zone of the solidifying casting, using a finite element mesh, and its implementation as part of the "PoligonSoft" computational casting process modeling system.

FOUNDATIONS OF THE POROSITY MODEL

The formation of shrinkage porosity occurs in the presence of a forming solid structure. The sizes of the interdendritic spaces determine the capillary forces acting in the biphasic zone. Since the magnitude of these forces can exceed the metallostatic pressure and the surrounding pressure, it can be assumed that in the formation of internal shrinkage cavities and macroporosity, capillary forces play an important role. Considering the capillary forces allows for the formulation of a model in which the formation of macroporosity depends on the dispersion of the dendritic framework, making the model more adequate and opening up possibilities to predict the pore size.

Consider the crystallization process of the liquid metal poured into the mold, which cools due to heat dissipation into the environment, both through the mold walls and directly from the free surface of the liquid metal. The casting, represented by its mesh model, is divided into elemental volumes: tetrahedrons with nodes located at their vertices. This mesh model is used to solve the thermal problem: calculating the temperature distribution in the casting block using the finite element method.

Suppose the temperature in the elemental volumes associated with node iii of the mesh is known from solving the corresponding heat conduction equation. We also assume that the free surface of the liquid metal (mirror) is flat and initially located at the upper edge of the casting. The phase composition of the metal in the elemental volume is characterized by the equilibrium volumetric proportions of the liquid phase f_L and solid phase f_S, whose calculation, based on the chemical composition of the alloy and the known temperature, can be performed with the help of one of the known thermodynamic databases.

Following the accepted representations of the model, we consider that if the liquid phase fraction satisfies the condition f_L>f_L^∗∗​, the solid phase does not form an immobile dendritic structure and moves along with the liquid metal mirror [1, 5, 6]. From the moment a rigid dendritic structure forms f_L= f_L^∗∗​​, the movement of the solid phase becomes impossible.

Fig. 1 - Schematic representation of the stages of solidification of the casting, where M is the surface of the liquid metal; L is the liquid metal; S is the solid phase; F is the mold; a) initial state; b) formation of shrinkage cavity; c) crystallization of the closed volume of liquid metal in a thermal node with the formation of an internal shrinkage cavity or dispersed porosity.

In general terms, the crystallization process of the casting goes through the following stages: formation of shrinkage cavities during the crystallization of open thermal nodes, formation of closed thermal nodes, and formation of internal shrinkage cavities or dispersed porosity (Fig. 1).

A thermal node is considered open if the liquid metal is in contact with the environment (Fig. 1b). If there is no direct contact of the liquid metal with the environment, the thermal node is considered closed (Fig. 1c).

Based on the crystallization sequence of the cast piece presented in Fig. 1, we will consider the mathematical models of shrinkage defect formation at each stage.

Cristalización de nodo térmico abierto

The crystallization of the liquid metal is accompanied by the contraction of the metal. In an open thermal node, crystallization does not cause a pressure drop if the contraction is compensated by the descent of the free surface (liquid metal mirror). The liquid metal surface is free and can move if it does not have a fixed solid phase framework, i.e., f_L>f_L^∗∗​​ (Fig. 2).

Fig. 2 - Scheme of crystallization of an open thermal node in the presence of a liquid metal mirror, where: a) standard model of "PoligonSoft"; b) new porosity model; 1 - zone of dry dendrites; 2 - zone fed by the liquid metal due to the capillary effect; S/L - zone of continuous dendritic structure, f_L ≤ f_L^∗∗

The descent of the liquid metal mirror must compensate for the metal contraction at each time step. The movement of the mirror is determined by the expression:

where V_Ω is the contraction volume, S_M​ is the area of the liquid metal mirror.

The contraction volume at this time step is equal to:

where v_Ωj​ is the contraction volume at the node of the finite element mesh; N is the number of mesh nodes within the thermal node (hereafter, within the zone).

According to the model conditions, in the mesh nodes located above the mirror, where there is no fixed dendritic framework, i.e., f_L>f_L^∗∗​, there should be no metal, and therefore, the actual fractions of the liquid phase g_L​ and solid phase g_S should be 0, and the porosity fraction g_P should be 1.

The volume of the casting assigned to each node of the finite element mesh is equal to one-fourth of the sum of the volumes of the elements to which that node belongs. Thus, the descent of the mirror below that node leads to the exclusion of the metal from the calculation, whose volume V_m​ may be greater than the contraction volume V_sh . The error introduced by this operation in the calculation of the casting mass increases with the size of the finite element mesh elements. To eliminate this error and maintain the casting mass constant in the model, a proportional reduction of the metal volumes in the nodes above the liquid metal mirror is adopted according to the expression: ϕV_m=V_sh , where ϕ is the proportionality coefficient. The fraction of liquid and solid phases in the nodes above the liquid metal mirror decreases according to the expression: g'=(1-ϕ)g where g and g' are the phase fractions at the beginning and end of the time step, respectively.

Fig. 2 schematically shows the principle of interaction of the biphasic zone of the casting with the liquid metal mirror, implemented in the macroporosity model of the "PoligonSoft" software and the new model proposed in this work. In the "PoligonSoft" software model, the movement of the liquid metal mirror in a fixed dendritic framework leads to the drying of the interdendritic spaces and the formation of macroporosity (see Fig. 2a, zone 1), whose volumetric fraction is numerically equal to f_L​ [7].

In the new model, it is assumed that the dendritic framework above the liquid metal level is completely filled with liquid metal due to the capillary effect. In the capillary feeding zone g_L=f_L  and f_g=f_S  (see Fig. 2b, zone 2). In the zone where f_L ≤ f_L^∗∗​ there is no free surface of the liquid metal.

Fig. 3 - Scheme of macroporosity formation in the biphasic zone at the boundary of an open thermal node, where L - liquid metal; S - solid phase; F - mold; 1 - zone of dry dendrites (macroporosity); r - radius of curvature of the liquid metal surface; S/L - zone of continuous dendritic structure.

The formation of a continuous dendritic framework (zones where f_L ≤ f_L^∗∗​​) around the liquid metal hinders the liquid metal's contact with the environment (see Fig. 3a). The free surface of the liquid metal, when situated within the dendritic framework, loses the ability to move freely. Due to the influence of capillary forces acting in the dendritic framework, the contraction of the metal during crystallization is only partially compensated by the change in the liquid metal level, leading to a decrease in pressure in the thermal node. The pressure distribution is determined by the expression:

(1)

where P_atm​ is the external pressure at the time of the thermal node formation, h is the height of the liquid metal column from the highest point in the thermal node, where f_L ≥ f_L^∗∗; E  is the compressibility modulus of the melt; V_Ω  is the contraction volume that has occurred in the thermal node at the current time step; V_P is the porosity volume produced in the current time step due to the change in the melt level; V_L is the volume of liquid metal in the thermal node.

Due to the resulting depression, the liquid metal is drawn towards the center of the thermal node, drying the interdendritic spaces on its periphery, leading to the formation of macroporosity (zone 1 in Fig. 3b). The pressure drop in the thermal node and the volume of porosity formed depend on the capillary forces acting in the dendritic framework. Conditionally, the expression defining the equilibrium between the forces drawing the liquid metal from the periphery to the center of the thermal node and the capillary forces opposing this process can be written as follows:

(2)

r is the radius of curvature of the menisci in the interdendritic spaces. To estimate r, the following expression is used:

where λ_II​ is the spacing between the secondary dendrite arms. Equation (2) allows determining the value of V_P and assigning porosity to the boundary nodes of the finite element mesh.

FORMATION OF CAVITIES IN CLOSED THERMAL NODES

From a certain point, due to the decrease in the liquid phase fraction, the boundaries of the thermal node become impermeable, and the metal contraction during crystallization is no longer compensated by the change in the liquid metal level in the dendritic framework, and the thermal node closes (see Fig. 4). This leads to an intense pressure drop in the thermal node, which is determined by the following expression:

(3)

The last term in the expression determines the pressure drop due to contraction. After some time, when the pressure at some point in the thermal node falls to a critical value P_crit, it becomes "energetically" favorable to form a new free surface of the liquid metal in the free liquid metal zone (f_L>f_L^∗∗). It should be noted that to form a new separation surface, some work is required, so P_crit<0​. [3, 8]

Fig. 4. Scheme of shrinkage cavity formation (b) in a closed thermal node (a), where L - liquid metal; S - solid phase; F - mold; M - free surface (mirror) of the liquid metal; V_P - volume of the shrinkage cavity

The appearance of a new flat separation surface fully compensates for the accumulated contraction in the thermal node from the moment of its isolation. Therefore, the location of the free surface of the liquid metal can be determined from the condition of equality between the contraction volume and the volume of the cavity formed: V_P=V_Ω .

With the appearance of the free surface of the liquid metal, the thermal node reopens in the sense that crystallization contraction will henceforth be compensated by the descent of the liquid metal's free surface. In the node, a shrinkage cavity will form according to the algorithm described earlier (see Fig. 2).

FORMATION OF MACROPOROSITY IN CLOSED THERMAL NODES

If a closed thermal node has a fixed solid phase framework everywhere, i.e., f_L ≤ f_L^∗∗​​, the formation of a flat free surface of the liquid metal is impossible.

The rupture of the liquid metal's continuity due to the pressure drop occurs through the formation of macroporosity. The first pore forms when the pressure falls below the critical value P_crit​. The most likely place for the pore to initiate is the point where the pressure is minimal, i.e., where the driving force for the formation of the separation surface Φ_nuc is maximum.

Once the pore formation process has begun in the zone, it is likely that the contraction occurring in the next time step will be compensated by the growth of the existing pores. The condition for this process is  P^k+1<-Pσ^k​ in those nodes of the finite element mesh where there is porosity. Here, k is the number of the time step. The driving force of this process is characterized by the function:

As the liquid metal crystallizes, the radius of the channels in the solid phase framework decreases, the capillary pressure at the pore edge increases, and consequently, the pressure in the liquid metal decreases, which can ensure the growth of the pore (see Fig. 5). If the pressure required for the growth of the existing pore is less than P_crit, it is more likely that a new porosity zone will form. The relationship between Φ_poro​ and Φ_nuc determines which of the processes, creation of a new porosity zone or development of the existing one, will predominate.

When Φ_poro​ ≥ Φ_crit , the existing porosity zone will develop; when Φ_poro​< Φ_nuc , a new zone will form.

Fig. 5 illustrates the algorithm for modeling pore formation in the dendritic structure.

Fig. 5. Formation of macroporosity in a closed thermal node; a) scheme of pore formation in an immobile dendritic structure; b) forces acting at the pore-liquid metal interface; c) porosity distribution in the finite element nodes. L - liquid metal; S - solid phase; V_P - pore; r - radius of curvature of the meniscus; λ_I y λ_II - distances between the primary and secondary dendrite arms.

The pore formed in node iii initially occupies a volume V₁, assigned to that node. If this volume is less than the contraction volume that needs to be compensated, the pore "grows" at the expense of neighboring nodes in the finite element mesh with a total volume V₂ and so on (see Fig. 5c). The pore growth stops when the negative pressure in the thermal node balances with the capillary forces.

(4)

where :

P_atm​ is the atmospheric pressure, h is the height of the liquid metal column; P_σ are the surface tension forces; P_ε is the pressure drop due to contraction.

Based on equation (4), the pressure balance for node i at the boundary between the liquid metal and the pore can be written as:

(5)

where g^k_L and g^k+1_L​ are the actual fractions of the liquid phase in the node before and after the formation of the pore (at time steps k and k+1); v_ω(j) is the volume associated with node number j in the list ω​. The sum is taken over the set of nodes ω(j) that participated in the formation of the pore.

The growth of the porosity zone continues as long as the left side of equation (5) is greater than the right side. To meet this equality, in addition to the sum of the neighboring nodes of node i where the pore originated, the sums of the neighbors of these neighbors must be included –ω₁,ω₂,...ω_n:

(6)

Here ⟨f_L^k+1⟩_ω is the average value of the liquid phase fraction in the set of nodes ω(j).

V_1​=g_l,iv_i is the contraction volume compensated by the drying of the volume in the node where the pore originates; V₂ is the contraction volume that can be compensated by the removal of all the liquid metal from the neighboring nodes of the node where the pore originates (from its "star"), i.e.:

V₃ is calculated for the nodes that form a star for the nodes considered in V₂, and so on. The proportionality coefficient ϕ characterizes the degree of drying of the nodes in the set ω_n that are located on the periphery of the formed porosity zone. The solution of the balance equation (6) allows determining the nodes in which, for pressure balance reasons (4), porosity should be assigned.

Below is one of the possible solutions to equation (6). Formalizing the sequence of pore propagation through the nodes of the finite element mesh, we write:

Equation (6) takes the form:

(7)

Based on this solution, porosity g_P=g_L is assigned to the nodes in the set ω₁,...ω_n-1 . In the nodes of the set ω_n​ porosity g_P=ϕg_L is assigned.

RESULTS AND DISCUSSION

The model presented in this article has been implemented in the "PoligonSoft" casting simulation system for trial tests.

Figure 6 presents the results of the simulation of the solidification process of a large working blade of a gas turbine (GTU) made with the ChS70 alloy using the technology described in the work [9].

Fig. 6. Result of the porosity modeling in the working blades of the ST-20 gas turbine; a) exterior view of the casting block; b) P_crit= -1 MPa, λ_II​ = 30 μm, E = 200 MPa;
c) P_crit= -0.1 MPa, λ_II​ = 300 μm, E = 2000 MPa; d) = -0,1 MPa, = 30 μm, E = 2000 MPa.
1 - cast part; 2 - ceramic shell; 3 - thermal insulation.

The calculations were performed considering f_L^∗∗=0.7 and different values of the parameters P_crit , λ_II​ and E, whose range of variation was chosen based on the following considerations:

The formation of shrinkage pores occurs when the pressure in the molten metal drops to a critical value P_crit ; the pressure in the melt at the time of shrinkage pore formation is negative, which can be estimated from the Laplace equation P=-2σ/r

Since the surface tension σ of molten nickel is approximately 1.7 N/m and the diameter of the pores observed in this casting ranges from 3.5 to 60 µm, it can be assumed that the pressure in the melt at the time of pore formation can vary from -0.1 to -1 MPa.

The meniscus radius r depends on the size of the interdendritic spaces, which are determined by the cooling rate of the melt in the biphasic zone. It was assumed that during the solidification of castings of this size, possible interdendritic spacings range from 20 to 300 µm.

The compressibility modulus of the melt can be estimated using the formula E=a²ρ_l , where a is the speed of sound in the melt; ρ_l is the density of the melt. Based on known data about the speed of sound in pure molten metals, taken from [10], it can be assumed that E=3*10⁴÷10⁵ MPa. It should be noted that in alloys, the speed of sound is significantly lower than in pure metals.

In the proposed model, the compressibility modulus characterizes the pressure drop process in the thermal node. Under ideal conditions, the rate of pressure drop in a closed thermal node is proportional to E. During the crystallization of a real casting, the solid metal crust surrounding the thermal node may not be airtight, which reduces the rate of pressure drop in the node. There is also the possibility of deformation of the crust under the action of the pressure difference between the environment and the inside of the thermal node, which also reduces the rate of pressure drop, as part of the crystallization contraction is compensated by deformation.

This indicates that the effective compressibility modulus of the melt must be significantly lower than the theoretical estimate. Within the framework of this model, the listed phenomena are not considered, and therefore, the compressibility modulus, like the critical pressure, is a fitting parameter that must be determined based on experimental data.

As expected, the volumetric fraction of porosity depends on the choice of model parameters P_crit and λ_II, which determine the moment when the pressure reaches the point at which porosity begins to form.

In the variants shown in Figures 6a and 6b, the porosity in the turbine blade does not exceed 0.2%. In Figure 6c, the porosity in the blade reaches 1.5%. The location of the porosity zone practically does not depend on the model settings and is determined by the casting geometry and the pouring technology.

In general, the simulation results are quite consistent with the results of the metallographic study of the castings and with the simulation in the standard model of the "PoligonSoft" system [9].

To assess the sensitivity of the model to the parameters  P_crit , λ_II​ and E, the solidification process of St.3 steel castings was also simulated. The solidification took place in a permanent mold made of St.3 on a sand base (Fig. 7).

Fig. 7. Result of porosity simulation in a casting made of St.3. a) exterior view of the casting block; 1 - casting; 2 - St.3 metal mold; 3 - sand; b) Standard PoligonSoft model;
c) P_crit= -1 MPa, λ_II​ = 300 μm, E = 2000 MPa;
d) P_crit= -20 MPa, λ_II​ = 300 μm, E = 2000 MPa;
e) P_crit= -20 MPa,  λ_II​= 300 μm, E = 200 MPa;
f) P_crit = -20 MPa,  λ_II​= 10μm, E = 200 MPa;
g) P_crit = -20 MPa,  λ_II​= 1 μm, E = 200 MPa

The thermophysical properties of the casting were calculated using the Computherm thermodynamic database [11] based on the alloy's chemical composition. The initial temperature of the metal was 1550°C, and the mold and sand were at 500°C. Cooling was conducted in ambient air at 20°C through convection with a heat transfer coefficient of 10 W/m²/K and radiation, with an emissivity factor for both the mold and the metal of 0.8.

For comparison, calculations were performed in the "PoligonSoft" system with standard configurations used for simulation in the casting production of the FGUP NPCG "Salut" company.

The simulation results are presented in Figure 7. The solidification of the casting in all presented variants corresponds to the scheme shown in Figure 1. In the first stage, the crystallization of the open thermal node occurs, compensating for contraction by the displacement of the casting mirror. In the second stage, a closed thermal node is formed where, depending on the set model parameters, a closed cavity (Figures 7c and 7d) or dispersed macroporosity (Figures 7e-7g) is formed.

The realization of one of these variants depends on the rate of two processes: the pressure drop in the thermal node due to crystallization contraction and the growth of the solid phase fraction leading to the formation of a continuous dendritic structure.

A rapid pressure drop in the thermal node to a critical value P_crit before the formation of a continuous dendritic structure leads to the formation of a free surface and an internal shrinkage cavity. This mode of porosity formation is favored by a large compressibility modulus E, a low critical pressure P_crit and a high value of f_L^∗∗. As can be seen in Figures 7c and 7d, a compressibility modulus value of 2000 MPa is quite large, ensuring a rapid pressure drop in the thermal node and the formation of a casting mirror at any critical pressure (down to -20 MPa). Therefore, assigning a higher value of E does not lead to a change in the porosity prediction.

With a low compressibility modulus, the pressure in the thermal node drops slowly and reaches a critical value when the continuous dendritic structure has already formed. In this case, the accumulated crystallization contraction manifests as dispersed macroporosity (Figures 7e, 7g). The volumetric fraction and size of the formed pores depend on the amount of liquid phase and the dendritic structure parameters, i.e., the spacing between the secondary dendrite arms.

REFERENCES

[1] Zhuravlev V.A., Sukhikh S.M. Computer modeling of the formation of distributed porosity and shrinkage cavities during the crystallization of alloys in ingots // Izvestiya AN SSSR. Metals. No. 1, 1981, pp. 80-84.
[2] Zhuravlev V.A., Bakumenko S.P., Sukhikh S.M. et al. On the theory of the formation of closed shrinkage cavities during the crystallization of alloys in large volumes // Izvestiya AN SSSR. Metals. No. 1, 1983, pp. 43-48.
[3] Zhuravlev V.A. On the macroscopic theory of alloy crystallization // Izvestiya AN SSSR. Metals. No. 5, 1975, pp. 93-99.
[4] Dawei Sun and Suresh V. Garimella. Numerical and experimental investigation of solidification shrinkage // Numerical Heat Transfer, Part A, 52, 2007, pp. 145-162.
[5] Casting simulation system "PoligonSoft", trademark of "CSoft Development", http://csdev.ru.
[6] Monastyrsky V.P. Modeling of macroporosity and shrinkage cavity formation in castings // Russian Foundryman. No. 10, 2011, pp. 16-21.
[7] Tikhomirov M.D. Fundamentals of casting process modeling. Shrinkage problem. - Moscow: Appendix to the journal "Foundry Production". - 2001. No. 12, pp. 8-14.
[8] Kent D. Carlson, Zhiping Lin, C. Beckermann, G. Mazurkevich, and Marc C. Schneider. Modeling porosity formation in aluminum alloys. Proceedings MCWASP-XI, ed. Ch.-A. Gandin, M. Bellet (Warrendale, PA: The Minerals, Metals & Materials Society, 2006), pp. 627-634.
[9] Monastyrsky V.P., Monastyrsky A.V., Levitan E.M. Development of casting technology for large gas turbine blades for power plants using "Poligon" and ProCAST systems // Foundry Production. No. 9, 2007, pp. 29-34.
[10] Physical quantities: Reference manual / Babichev A.P., Babushkina N.A., Bratkovsky A.M. et al.; Ed. Grigoriev I.S., Meilikhov E.Z. - Moscow: Energoatomizdat, 1991. - 1232 pp.
[11] Thermodynamic database for iron-based superalloys: PanIron 5.0, CompuTherm, LLC, USA.

Translated by  A.J. Camejo Severinov
‍Original text in Russian
‍