The numerical simulation of temperature fields during metal casting is a critical tool for modern foundry engineering. Its primary objective is to solve the solidification heat transfer model using numerical methods, thereby reconstructing the temperature distribution throughout the casting process on a computer. This capability allows for the prediction of defect formation—such as shrinkage porosity, cavities, slag inclusions, and hot tears—by identifying their probable locations and root causes. Consequently, it provides a scientific basis for optimizing casting parameters and gating system design, ultimately leading to significant improvements in the quality and yield of cast components. In recent years, alongside advancements in computational power and numerical algorithms, numerous commercial software packages have emerged for casting simulation. While powerful, analyses conducted solely within these software environments can sometimes obscure the underlying physical phenomena governing heat transfer at critical boundaries. This work, therefore, employs a fundamental numerical approach to dissect the influence of interfacial thermal resistance, a key physical parameter, on the solidification behavior of magnesium alloy sand castings.
The thermal process in sand castings is governed by the classical heat conduction equation. In a two-dimensional Cartesian coordinate system, the governing differential equation is expressed as:
$$
\frac{\partial T}{\partial t} = \alpha \nabla^2 T = \alpha \left( \frac{\partial^2 T}{\partial x^2} + \frac{\partial^2 T}{\partial y^2} \right)
$$
Here, \( T \) represents temperature, \( t \) is time, \( x \) and \( y \) are spatial coordinates, and \( \alpha \) is the thermal diffusivity. The thermal diffusivity is defined by the material’s intrinsic properties: \( \alpha = \lambda / (\rho c_p) \), where \( \lambda \) is the thermal conductivity, \( \rho \) is the density, and \( c_p \) is the specific heat capacity.
The geometry under investigation is a classic T-shaped casting, commonly studied for its sensitivity to shrinkage defects due to varying section thicknesses. The computational domain encompasses not only the casting but also a significant portion of the sand mold to accurately capture the heat sink effect. The key dimensions for the model are summarized in Table 1.
| Parameter | Value (cm) |
|---|---|
| Base Width | 70 |
| Base Height | 12 |
| Riser Width | 30 |
| Riser Height | 60 |
| Mold Thickness (All sides) | 10 |
The total computational domain is a square region of 82 cm × 82 cm. This domain is discretized using a uniform structured grid with a spatial step of \( \Delta x = \Delta y = 1 \) cm, resulting in 82 × 82 = 6724 computational cells. The material assignment—distinguishing between AZ91 magnesium alloy and the silica sand mold—is made within this grid based on the geometrical boundaries defined in Table 1.
The thermophysical properties of the materials are crucial inputs for an accurate simulation. The values selected for AZ91 magnesium alloy and the green sand mold are listed in Table 2. Note that the alloy’s solidification range is defined by its liquidus and solidus temperatures.
| Material | Density, \( \rho \) (g/cm³) | Specific Heat, \( c_p \) (J/g·K) | Thermal Conductivity, \( \lambda \) (J/cm·s·K) | Initial Temp. (K) | Liquidus/Solidus Temp. (K) |
|---|---|---|---|---|---|
| AZ91 Magnesium Alloy | 1.81 | 1.02 | 0.51 | 871 | 708 / 693 |
| Silica Sand Mold | 1.60 | 0.27 | 0.0025 | 293 | N/A |
The boundary conditions are defined in terms of thermal resistance, \( R \), at the interfaces. Thermal resistance (in s·cm²·K/J) is the inverse of the interfacial heat transfer coefficient (IHTC). A higher resistance indicates a greater impediment to heat flow. The baseline conditions for the simulation are defined in Table 3. Three distinct interfacial resistances are considered: between the casting and the mold (\( R_{c/m} \)), between the casting and the ambient air at the riser top (\( R_{c/a} \)), and between the outer mold surface and the ambient air (\( R_{m/a} \)).
| Interface | Thermal Resistance, \( R \) (s·cm²·K/J) |
|---|---|
| Casting / Mold (\( R_{c/m} \)) | 1500 |
| Mold / Air (\( R_{m/a} \)) | 5000 |
| Casting / Air (\( R_{c/a} \)) | 120 |
To solve the transient heat conduction problem, a direct finite difference method based on an explicit forward-time scheme is implemented. This method applies the principle of energy conservation directly to each discrete control volume (grid cell). The change in internal energy of a cell over a time step \( \Delta t \) is equal to the net sum of heat fluxes entering through its faces.
Considering a cell indexed (\( i, j \)) with volume \( \nu = \Delta x \Delta y \cdot 1 \) (for unit depth), its temperature at the new time level \( n+1 \) is calculated from its temperature at the previous time level \( n \) by:
$$
T_{i,j}^{n+1} = T_{i,j}^{n} + A_x \cdot TNS_{i,j}^{n}
$$
Where:
- \( T_{i,j}^{n} \): Temperature of cell (i,j) at time step \( n \).
- \( A_x = \frac{\Delta t}{\rho c_p \nu} \): A material-dependent coefficient.
- \( TNS_{i,j}^{n} \): The “Total Net Sum” of heat flow rates into the cell (i,j) during the time step from its four neighboring cells (North, South, East, West) and from boundary interfaces, if applicable.
The heat flow between two adjacent internal cells, for example from cell (i,j) to (i+1,j), is calculated using Fourier’s law in discrete form:
$$
Q_{i \rightarrow i+1} = \frac{\lambda_{eff} \cdot A_{face} \cdot (T_{i,j}^{n} – T_{i+1,j}^{n})}{\Delta x}
$$
Here, \( A_{face} \) is the area of the shared face, and \( \lambda_{eff} \) is an effective conductivity, often taken as the harmonic mean for dissimilar materials. For a boundary cell, the heat flux is governed by the interfacial thermal resistance. For instance, the heat flow from a casting boundary cell to the mold interface is:
$$
Q_{boundary} = \frac{A_{face} \cdot (T_{cast}^{n} – T_{mold}^{n})}{R_{c/m}}
$$
The explicit time-marching scheme requires adherence to stability criteria. For a 2-D explicit formulation, the stability condition is given by the Fourier number constraint:
$$
Fo = \alpha \frac{\Delta t}{(\Delta x)^2} \leq \frac{1}{4}
$$
Given the material properties and grid size (\( \Delta x = 1 \) cm), a time step of \( \Delta t = 0.02 \) s was selected, which comfortably satisfies this condition for both materials, ensuring a stable and accurate numerical solution. The simulation code was implemented from the ground up, and the resulting temperature field data was processed for visualization and quantitative analysis.
Prior to investigating parametric effects, the fundamental validity of the model was checked against general experimental trends. The simulation captures the characteristic temperature history at the casting/mold interface: a rapid initial drop followed by a progressively slower cooling rate as the temperature gradient decreases. This trend aligns perfectly with experimental measurements reported in literature for different alloy systems, confirming that the direct finite difference approach is fundamentally sound for simulating the thermal history of sand castings.
The interfacial thermal resistance between the casting and the mold (\( R_{c/m} \)) is arguably the most critical parameter in simulating sand castings, as it dictates how efficiently heat is extracted from the metal into the mold. A parametric study was conducted by varying \( R_{c/m} \) while keeping other resistances at their baseline values. The temperature distributions at a fixed solidification time of 4000 seconds are profoundly different, as conceptually summarized in Table 4.
| \( R_{c/m} \) Value (s·cm²·K/J) | Dominant Heat Path | Solidification Sequence | Temperature Gradient Description |
|---|---|---|---|
| 3000 (High) | Casting/Air (Riser Top) | Strong directional solidification from riser downwards. | Steep gradient. Five distinct, strong temperature bands from cold riser to hot base. |
| 1500 (Baseline) | Mixed: Riser Top & Mold Sides | Directional from riser, but base cooling enhanced. | Moderate gradient. Five bands visible, but the hottest region (base) is smaller. |
| 750 (Low) | Mold Sides become very effective. | Nearly simultaneous cooling from riser top and mold sides/base. | Shallow gradient. The hottest central region shrinks further; distinct banding weakens. |
| 150 (Very Low) | Mold Sides are highly effective. | Strong simultaneous cooling from all surfaces. | Very shallow gradient. The casting approaches a near-uniform temperature. No distinct hot spot. |
Quantitatively, the temperature evolution at a strategic “corner” node within the thick base section of the T-casting reveals critical insights. The cooling curves for different \( R_{c/m} \) values show significant divergence during early solidification. With high resistance (\( R_{c/m} = 3000 \)), the corner cools very slowly initially because heat extraction through the adjacent mold wall is poor. With very low resistance (\( R_{c/m} = 150 \)), the corner cools rapidly from the start. However, as solidification progresses and the internal temperature gradients diminish, all curves tend to converge toward similar cooling rates, indicating that internal conductive resistance eventually dominates over the interfacial resistance. This behavior can be related to the governing equation where the driving force (\( \Delta T \)) decreases. The initial condition is heavily influenced by the boundary condition \( \frac{\partial T}{\partial n} \propto \frac{1}{R} \Delta T \), but at later times, the internal Laplacian term \( \nabla^2 T \) becomes the rate-limiting factor.
In many sand castings with open risers, the top surface of the riser is exposed to air. The thermal resistance at this casting/air interface (\( R_{c/a} \)) controls the heat loss by radiation and convection from this crucial “hot-top.” The study varied \( R_{c/a} \) while keeping \( R_{c/m} \) constant at 1500 s·cm²·K/J. The impact on the temperature field is summarized in Table 5.
| \( R_{c/a} \) Value (s·cm²·K/J) | Riser Top Heat Loss | Solidification Sequence | Temperature Gradient Description |
|---|---|---|---|
| 120 (Low / Baseline) | Very efficient. | Strong directional solidification from riser downwards. | Steep gradient. Five strong temperature bands, coldest at riser. |
| 480 (Moderate) | Less efficient. | Weakened directional solidification. Base cooling becomes relatively more important. | Moderate gradient. Fewer cold bands, hotter riser and larger hot zone in base. |
| 840 (High) | Inefficient. | Directional solidification from riser is lost. | Shallow gradient. Most of casting is at high temperature. No clear cold riser. |
| 1200 (Very High) | Highly insulated. | Simultaneous solidification from all mold surfaces. Riser acts as a hot spot. | Very shallow gradient. Entire casting remains hot for longer. Defect risk in riser neck/base. |
The temperature history at the center of the riser highlights this effect. During the very initial moments, all curves overlap because the heat loss is initially dominated by the massive temperature difference. However, as the surface temperature drops, the imposed resistance begins to govern the flux. The cooling rate at the riser node decreases dramatically as \( R_{c/a} \) increases. A high \( R_{c/a} \) essentially insulates the riser top, nullifying its intended function as an efficient heat sink and a feeder. This leads to a reversal of the desired thermal gradient, potentially causing shrinkage defects in the casting body rather than in the riser. This is a critical finding for the design of riser covers or insulating toppings in sand castings.
The resistance at the outer mold-air interface (\( R_{m/a} \)) was also varied across a wide range. Surprisingly, the simulated temperature distributions within the casting at 4000 seconds showed negligible visual difference. This indicates that for a typical sand casting with a reasonable mold thickness, the primary thermal barrier is the casting/mold interface itself, combined with the low thermal conductivity of the sand. The additional resistance of the mold-air boundary has a minimal secondary effect on the casting’s internal thermal history because the mold itself acts as a significant thermal capacitor and resistor. The heat flux bottleneck occurs at the first interface (\( R_{c/m} \)) and within the mold body. Therefore, while insulating the outer mold wall might slightly change the mold’s own temperature profile, it does not significantly alter the heat extraction dynamics from the casting in this configuration. This conclusion simplifies modeling efforts for sand castings, as precise knowledge of external mold convection coefficients is often not critical.
Through a fundamental numerical analysis of a T-shaped AZ91 magnesium alloy sand casting, this work elucidates the distinct roles of various interfacial thermal resistances. The casting/mold resistance (\( R_{c/m} \)) is the paramount parameter, directly controlling the early-stage cooling rates and the establishment of the initial thermal gradient. Lowering this resistance promotes more uniform cooling and can eliminate isolated hot spots. The casting/air resistance (\( R_{c/a} \)) at exposed surfaces like riser tops is equally critical for directional solidification control. Increasing it (e.g., with an insulating sleeve) can severely impair riser efficiency and reverse the thermal gradient, leading to defective sand castings. Conversely, the mold/air resistance (\( R_{m/a} \)) has a negligible impact on the casting’s solidification pattern for molds with adequate thickness. These insights underscore the importance of accurately characterizing or selecting the casting/mold and casting/air interface conditions in numerical simulations to reliably predict solidification patterns and optimize the quality of magnesium alloy sand castings. The direct finite difference method proves to be a transparent and effective tool for developing this fundamental understanding.