Stress analysis in solid-liquid parts of solidifying castings

In the paper, we present results of stress analysis in domains which are a mixture of solid and liquid phases. Such mixtures occur in solidifying castings and are a result of forming a structure with solid skeleton and filling of a liquid phase. In this structure, stress occurs due to the appearance of temperature gradients, different values of material properties for the solid and liquid phase, and the appearance of friction forces between the solidified part of the casting and the mold on a macroscopic scale. This can lead to casting defects, such as hot cracking. The results are obtained with the use of a authors computer program based on the Finite Element Method. The stress analysis takes into account the elastic-plastic state of considered computational area. The presented results are focused on the microscopic scale, for which a finite element mesh is created which imitates the growing grains of the metal alloy in the casting, on the basis of macroscopic parameters.


Introduction
The matter of stress analysis in an elastic-plastic state was favorite a scientific problem for many past years. Nowadays, research on that phenomenon come back and are subjected to further analysis.
Investigation of solidification problems is essential because of a variety of industrial applications and material processing technologies such as casting, welding, the growth of single crystals, purification of metals, and numerous other applications. In recent decades, vast research has been done to attain a correct understanding of this vital process. Determination of the stresses distribution during the solidification process helps to assess the quality of casting products and optimize the production circumstances. Undesirable stresses and deformations that occur in the solidification process are due to the high gradient of temperatures inherent in this process and the succeeding cooling [1]. High-temperature gradients and phase change lead to significant changes in thermomechanical properties. Moreover, in most studied cases the mechanical and thermal properties of a particular material at elevated temperatures are the major factors that make the solidification problem even more complicated and confusing. The mushy zone (the state between liquid and solid state) plays an essential role in determining the properties of a final product in the solidification of alloys [2]. The rate of solidifying as well as stress distribution in the mushy zone directly affects the morphology of crystallographic growth, hot cracking or gap, and shrinkage formation. In the present work, the stress analysis in solid-liquid parts of solidifying castings is presented. Stress analysis in a solid-liquid state is the preliminary studies related to a determination of the hot tearing criterion during the solidification process. To predict hot tearing must be taken into account many issues, including constitutive equations of alloys in the high-temperature range containing partially solidified state. As of now, tests have been conducted on several aluminum alloys in this state to obtain their mechanical properties, and the corresponding constitutive equations have been created. Unfortunately, rheological behaviors of the alloys have not been obtained what is mainly because the rheological properties are generally obtained using relations between steady-state stress and strain rate in tensile conditions, and the brittleness of the partially solidified alloys makes them difficult to obtain [3]. In our future work, we intend to exploit thermal stress analysis and a relation between temperature and a solid fraction.
In the second half of the last century, some analytical solutions of stress development during the solidifying process have been provided. Although many simplifying assumptions in these studies, they are handy for benchmarking the results of numerical methods. In the subsequent years, various numerical methods have been employed for the determination of stresses in solidifying bodies. Among others, there were used control-volume finite differences method (FDM), boundary elements method (BEM) or finite elements method (FEM) [4][5][6]. The using variety of numerical algorithms, there are still analyzed almost every issue, i.e., technical or engineering problems.

Methods
Stress analysis is a general term used to describe the quantities stress and strains. It also is known as structural analysis. Dependency between stress and strain characterizes the stress model, got from an experiment, and by the relationship of strain and displacement, got from geometric considerations. Stress is related to yield through the physical relationship: where σ is a stress tensor, C is a stiffness tensor, and ϵ is a tensor of elastic deformation. Physical properties that appear in elasticity tensor can depend on temperature. The relationship from (1) is called generalized Hook's law. In turn, the strain is related to displacement through the Cauchy relations [7]: where q is the displacement vector. This work focuses on stress distribution in three-dimensional case, so the stiffness tensor becomes: where: also, E is Young's modulus, is Poisson's ratio.
In this case, the stress tensor has the following structure: While the tensor of elastic deformation is given by: where = 2 . And the displacement vector is equal: where u, v, w are displacements in the direction of x, y, z axes, respectively. Strains can be divided into elastic and plastic range: where is vector of elastic strain and is vector of plastic strain. Substituting above relationship into generalized Hooke's Law, we get: Transition from elastic state into plastic is possible only, when stress achieve level described by function called plasticity condition. Here, the Huber-Mises-Hencky condition is used, which states that material goes into plastic state, when energy achieves level specific for material. H-M-H condition can be expressed as: where f is a plasticity surface, F is a plasticity function and is yield stress. In H-M-H condition, plasticity function contains only second invariant of stress deviator D. So, the above condition takes form: where 2 is second invariant of stress deviator D.
Because plastic strains forces nonlinear relation between displacements and stress, we have to use linearization. Moreover, appearance of plastic strains forces using incremental applying of external forces. This modifies the standard elastic formulation: into incremental form: where K is global stiffness matrix given by: R is the load vector of external forces, F is equivalent force of stress: and ∆ is displacement increment. The B matrix depends on derivatives of shape functions used in Finite Element Formulation and C ep is elasto-plastic stiffness tensor The procedure for solving is performed in iterative way, until convergence is achieved [ Stress value depends on state in which given finite element is. In order to calculate stresses following procedure is utilized: 1. Assuming elastic behavior, compute the strain increment, the trial stress increment and the trial stress by: 2. Check if is beyond plasticity surface. If not then is considered as true stress. If yes than correct stress with the formula: where = √ 2 and is the mean trial principle stress.

Numerical simulations
Simulation tools become indispensable for engineers who are interested in tackling increasingly larger problems or the ones who are interested in searching larger phase space of process and system variables to find the optimal design. Advances in hardware allow not only to solve the larger tasks (using more detailed grids), but also to describe the problem more accurately. Increasing capacity of computer memory makes it possible to consider growing problem sizes. At the same time, increased precision of simulations triggers even greater load. There are several ways to tackle these kinds of problems. For instance, one can use parallel computing [9], someone else may use accelerated architectures such as GPUs (Graphics Processing Units) [10], while another person can use the special organization of computations [11]. We have used parallel processing to assure effective time duration for simulations. During the implementation, we used TalyFem and PETSc libraries, which allowed us to split structures, such as matrices and vectors, into many computing nodes [12] [13]. The plasticity model described in previous section was used in analysis of stress distribution of solidifying casting. Castings made from alloys solidifies not in fixed temperature, but rather in temperature interval. This results in the presence of different zones: fully liquid, fully solid and mushy zone made from mixture of solid and liquid material. The presented simulations tried to imitate presence of different zones, by using domain presented in Fig. 1.

Fig. 1. Computational domain with tetrekaidecahedrons.
This domain is made from tetrakaidecahedrons that represent grains in casting and the area that fills space between them. While tetrakaidecahedrons can fully fill domain, it was chose to scale them down. This way it was possible to introduce the mushy zone. In simulations, the radii of sphere described the tetrakaidecahedron was equal to 0.08 mm. In order to prepare the geometry, a Python script was prepared that allows to generate domain using basic parameters (number of tetrakaidecahedrons in each direction, size of tetrakaidecahedron) that later is meshed by the Gmsh software.
The simulations were made for aluminum alloy in temperature above 800 K. Because of temperature, material properties were equal to: Young's modulus 1.4·10 7 Pa in grain and 1.25·10 7 Pa in mushy zone, the yield stress was equal to 2.7·10 4 Pa in grain and 2.5·10 4 Pa in mushy zone. The Poisson coefficient was equal to 0.35 in both materials. The domain was supported on three surfaces that went through origin of coordinate system. Each surface constrained one direction (i.e. surface parallel to x plane constrained x direction). Load was prescribed only in x direction on surface opposite to constrained x plane and was equal to 15 kPa.  The results are presented on Figs. 2 and 3. We can see that the domain deforms as a single region, but in places of connection between mushy and grain zones the accumulation of stress occurs.