82 | Numerical analysis of rapid solidification of NiTi alloy: Influence of boundary conditions Numerična analiza hitrega strjevanja NiTi zlitine: vpliv robnih pogojev Matej Zadravec 1,* , Primož Ternik 2 , Rebeka Rudolf 1 and Milan Svetec 3 1 Faculty of Mechanical Engineering, University of Maribor, Smetanova 17, 2000 Maribor, Slovenia 2 Ternik Primož - Zasebni raziskovalec, Bresterniška ulica 163, Bresternica 3 PAZU, Lendavska ulica 5a, Murska Sobota * Corresponding author e-mail: matej.zadravec@um.si Abstract: The objective of the present work is to analyze influence of temperature boundary conditions on the solidification process of NiTi alloy. Alloy is taken as an incompressible fluid where the heat is transferred by conduction and convection, including the thermal phase change phenomenon. The last one is modelled by the improvement procedure, so-called enthalpy- porosity formulation, where the liquid-solid mushy zone is treated as a porous zone with porosity equal to the liquid fraction. The numerical model is based on the finite volume method in body fitted coordinates with a SIMPLER scheme to couple the pressure and velocity fields. Simulation presents solidification for two cooling cases on the cylindrical part, where in the first case the adiabatic boundary condition is considered and in the second case the convective cooling over the wall is present. The results are presented for the velocity and temperature field as well as for the NiTi mass fraction during the solidification process. Results shows that velocity and temperature field is strongly affected by the different cooling condition on the cylindrical wall and therefore the solidification process of the alloy. Analysis shows that cooling on the cylindrical part is one of the major parameters for alloy solidification and therefore should not be neglected. Key words: Solidification, Finite volume method, Phase change. Povzetek: Delo obravnava vpliv temperaturnih robnih pogojev na strjevanje NiTi zlitine s pomočjo numerične simulacije. Zlitina je obravnavana kot nestisljiva tekočina, kjer se toplota prenaša s prevajanjem in prestopom, vključno s pojavom fazne spremembe. Numerično modeliranje fazne spremembe je izvedeno z izboljšanim pristopom, imenovanim formulacija entalpija-poroznost. Pri tem pristopu je območje vmesnega faznega stanja tekoče-trdno obravnavano kot porozno območje s poroznostjo, ki je enaka deležu tekoče faze. Numerični model temelji na metodi končnih volumnov, pri čemer sta tlačno in hitrostno polje povezana s SIMPLER algoritmom. Numerična simulacija podaja strjevanje NiTi zlitine v geometriji epruvete, ki je v prvem primeru po cilindričnem delu oboda idealno izolirana (adiabatni pogoj) in v drugem primeru se preko te površine hladi zaradi okoliškega zraka. V rezultatih so predstavljeni hitrostno in temperaturno polje, kakor tudi deleži posameznih faz NiTi zlitine med procesom strjevanja. Iz rezultatov je razviden precejšen vpliv različnega robnega pogoja na hitrostno in temperaturno polje, kar posledično vpliva tudi na strjevanje zlitine. Analiza je pokazala, da je ohlajanje preko stene epruvete eden izmed vplivnejših parametrov strjevanja zlitine in ga zato ne smemo zanemariti. Ključne besede: strjevanje, metoda končnih volumnov, fazna sprememba. Izvirni znanstveni članek TEHNIKA – materiali Datum prejema: 15. maj 2015 ANALI PAZU 4/ 2014/ 2: 82-88 www.anali-pazu.si | 83 NUMERICAL ANALYSIS OF RAPID SOLIDIFICATION OF NiTi ALLOY: INFLUENCE… 1. Introduction NiTi alloys belong to a group of shape memory materials which are known as functional materials. This definition is closely connected with two main properties: shape memory behaviour and superelastic effect [1]. It must be mentioned that the first observation of the shape memory behaviour of materials goes back to 1932, when Swedish researcher Arne Olander first noticed this phenomenon on a sample of gold and cadmium (Au-Cd). As late as in 1962, William J. Buehler and his colleagues in the "Naval Ordnance Laboratory" observed the shape memory effect in the alloy of nickel and titanium. This alloy was named NiTinol (Nickel-Titanium Naval Ordnance Laboratory) [2]. Shape memory describes the effect of restoring the original shape of plastically deformed material by heating it (Fig. 1). This phenomenon results from a crystalline phase change known as thermoelastic martensit transformation [3]. Superelastic effect is referred to as a reversible martensitic transformation, which is not due to change in temperature, but to the change in stress state. A sufficiently large load causes the transformation of austenite into martensite (Fig. 2). Unloading returns the transformation of the austenite and the original shape [2]. The reasons why the shape memory alloy NiTi was implemented in medicine or orthodontic treatment are the following: they have an excellent memory properties, excellent mechanical properties, good corrosion resistance and excellent biocompatibility. The most important characteristics for orthodontic application of NiTi alloys are superelasticity and biocompatibility. NiTi wires were used in the first (initial) phase of orthodontic treatment [4,5], which deals with alignment and de-rotation of teeth, as well as correction of vertical and horizontal discrepancies by levelling dental arches. In order to produce physiological response in the periodontal ligament (PDL) and bone, initial aligning wires should apply light continuous force on the teeth. The magnitude of the force caused by the orthodontic wire is largely dependent on the material properties (modulus of elasticity) [6,7]. This study investigated the possibility of production simplicity with the melting in further continuous vertical casting of Ni-Ti binary alloy with general composition 50,68051 at.% Ni and 49,31949 at.% Ti (or in wt.% 55,749516 Ni and 44,250483 Ti). For this purpose in the first step solidification and casting processes, modelling with the numerical approach was used. The geometry of casting device was cylinder like solidification device where the effect of thermal boundary conditions on solidification process was observed. 2. Numerical model 2.1. Governing equations The time-dependent NiTi alloy solidification using ANSYS Fluent was modelled by the enthalpy-porosity technique [6]. In this technique, the melt interface is not tracked explicitly. Instead, a liquid fraction, which indicates the fraction of the cell volume that is in liquid form, is associated with each cell in the computational domain. The liquid-solid mushy zone is a region in which the liquid fraction values lie between 0 and 1 and its temperature ranges between the liquidus (T l ) and solidus (T s ) temperatures. The mushy zone is modelled as a “pseudo” porous medium in which the porosity decreases from 1 to 0 as the material solidifies. When the material has fully solidified, the porosity becomes 0 and hence the velocities also drop to 0. The enthalpy of the material is computed as the sum of the sensible enthalpy (h) and the latent heat (∆H): Deformation Strain (%) Fig. 1. Schematic presentation of the shape memory behaviour [1-2]. Austenite→ Martensite Avstenite←Martensite Stress (10 6 N/mm 2 ) Fig. 2. Superelasticity curve in the tension deformation diagram [1-3]. (1) 84 | ANALI PAZU, 4/ 2014/ 2, str. 82-88 Matej ZADRAVEC et al. where is reference enthalpy, reference temperature, specific heat at constant pressure, liquid mass fraction and latent heat of fusion. The liquid fraction, , is defined as: Finally, for the solidification problem, the energy equation reads as: where is enthalpy (see Eq. 1), is density and is velocity. The conservation equation of mass and momentum are decoupled from the one of the thermal energy. These equations are solved using a segregated solver with the second order accurate upwind scheme. 2.2. Geometry and boundary conditions Geometry of experimental tube considered in the simulation is presented in Fig.3. Due to the symmetrical geometry, the NiTi alloy solidification was modelled as an axisymmetrical problem. The upper horizontal and lower curved walls were maintained at the constant temperatures, whereas the other boundaries are considered once to be adiabatic in nature and second case was as they are exposed (cooled) to an air with temperature of 25°C. The temperature of the curved wall was lower than the solidus temperature in order to allow the solidification and provide the phase change. The physical properties of the melt are given in Table 1. Table 1. NiTi melt thermal properties used in the analysis. Such a geometrical set-up will further allow us to perform the directional solidification experiments, in which the solidification processing parameters (e.g. temperature gradient and growth rate) can be independently controlled [8] so that one may study the dependence of microstructural parameter on either temperature gradient at constant growth rate or growth rate at constant temperature gradient for the constant initial NiTi solute composition. To perform a second simulation case with convective boundary conditions on walls it was necessary to define the heat transfer coefficient from interior part of the geometry (solid or liquid alloy) to the outer air. The schema of experimental tube with corresponding temperature profile behaviour through a tube presented on Fig. 4 shows convective and conductive heat transfer mechanism through a tube. (2) (3) Fig. 3. Geometry and boundary conditions. Density Variable (from experiments) Specific heat 0.7988 J/kg K Latent heat 24200 J/kg Thermal conductivity 21.5 W/m K Viscosity 0.00574 Pa s Solidus temperature 750 K Liquidus temperature 1310 K | 85 NUMERICAL ANALYSIS OF RAPID SOLIDIFICATION OF NiTi ALLOY: INFLUENCE… First, we have to define the outer properties of the media, which was air with dimensionless Prandtl (eq.4) and Grasshoff (eq.5) number to obtain the Rayleigh number (eq.6) and finally the Nusselt number (eq.7). In general, we can use the Nusselt number (eq.8) to define the convective heat transfer coefficient ( ) on the outer side of the tube and accounting for the heat conduction in the tube wall, we finally get the overall heat transfer coefficient (eq.9). Because the heat transfer coefficient is dependent from the inner temperature of the solid or liquid alloy there exist the temperature dependency as presented in Fig. 5. 2.3. Numerical accuracy assessment Final numerical results should be grid independent. To establish grid independent results we performed a detailed analysis using for different meshes. With each grid refinement the number of elements in a particular direction is doubled and the minimum element size is halved. Finally we apply Richardson’s extrapolation technique, which is a method for obtaining a higher- order estimate of the flow value (value at infinite grid) from a series of lower-order discrete values [9-11]. For a general variable ϕ the grid-converged value according to the Richardson extrapolation is given as: where is obtained on the finest grid and is the solution based on the next level coarse grid, is the ratio between coarse to fine grid spacing and is the expected order of accuracy. The variation of liquid fraction and temperature with grid refinement is given in tabulated form in Table 2. Fig. 4. Schema for heat transfer in simulation CASE 2. (4) (5) (6) (7) (8) (9) Fig. 5. Heat transfer coefficient on the cooling wall for CASE 2. (10) Table 2. Effect of mesh refinement upon the liquid fraction and temperature at (at x=0 and z=0.1) for a steady-state numerical analysis. MI MII MIII MIV Num. el. 4700 18800 75200 300800 / / 0.5410 0.5401 0.5396 0.5393 0.5391 0.06% 1052.95 1052.43 1052.17 1052.03 1051.90 0.02% 86 | ANALI PAZU, 4/ 2014/ 2, str. 82-88 Matej ZADRAVEC et al. The “percent” numerical error as given in Table 2 is a quantification of the relative difference between the numerical predictions and the extrapolated value obtained with Richardson’s extrapolation technique. Results of calculations of numerical accuracy indicate that as the mesh is refined, there is a consistent improvement in the accuracy of the predicted values, and the agreement between predictions obtained with mesh MIII and extrapolated value is extremely good. Based on this, the simulations in the remainder of the paper were conducted on mesh MIII which provided a reasonable compromise between high accuracy and computational effort. 3. Results From the temperature plots shown in Fig. 6 is obvious that the movement due to the natural convection is rather slow and reduces as the time progresses. Due to that, one can conclude that in the case of NiTi solidification in an experimental tube, the heat conduction is the sole transport mechanism. To compare effect of different boundary condition on walls we can see that in the case 2 where the convective boundary condition is used the velocities are expected to be higher and also present at lather times of solidification. On the temperature contours plots (Fig. 7) the conductive mechanism can be observed because contours are mainly lines. The exception is the region in the vicinity of the lower curved wall, where contours of the temperature are slightly curved. However, even in this region, the heat conduction is prevailing mechanism and curved contours results from the curved geometry of the lower wall. When we compare results for both boundary conditions we see in the case 2 that the contours of the temperature are slightly curved over the whole geometry height because of the convective boundary condition, which causes increasing cooling at the walls. During the solidification the NiTi liquid mass fraction changes as is shown with contours on Fig. 8. Contours are lines signifying the unidirectional solidification process in the experimental tube. Furthermore, the solidification front moves upward, which is consistent with the time evolution of the temperature field. From this contours plots we can also clearly observe the effect of convective boundary conditions on the solidification front, which is much higher than in the case with adiabatic boundary condition. Fig. 6. Contours of vertical velocity component in the middle of the geometry for different solidification time (A- adiabatic wall & C-convective wall boundary condition). | 87 NUMERICAL ANALYSIS OF RAPID SOLIDIFICATION OF NiTi ALLOY: INFLUENCE… 4. Conclusions In the present study, the effect of different boundary conditions circumstances in simulation of time- dependent phase-change of NiTi alloy was studied by numerical means. Solidification process was modelled by the enthalpy-porosity method. From numerical results we can observe slight movement of the alloy in the vertical direction, and is also reduced as the solidification time progresses. Temperature field in alloy clearly shows, that the main energy transport mechanism in case 1 with adiabatic walls is conductive heat transfer, but in the case 2 with non-adiabatic boundary conditions on the walls the convective mechanism is also present. From solidification point of view we can conclude that solidification of process in an experimental tube is almost unidirectional in both cases and solidification front moves upward, where there is a large impact of convective boundary condition on solidification. From the present study we see that we should not neglect the convective heat transfer, because this will cause totally different solid liquid distribution in calculated domain. Literature 1. M. Thier, D. Treppmann, D. Drescher, C. Boureaul. Fig. 7. Contours of temperature for different solidification time (A-adiabatic wall & C-convective wall boundary condition). Fig. 8. Contours of liquid volume fraction for different solidification time (A-adiabatic wall & C-convective wall boundary condition). 88 | ANALI PAZU, 4/ 2014/ 2, str. 82-88 Matej ZADRAVEC et al. Transformation characteristics and related deformation behaviour in orthodontic NiTi wire. Journal of Materials Science: Materials in Medicine. Volume 3, Number 3, 229-233, DOI: 10.1007/BF00713455 (2004). 2. D. Cimprič. Shape memory alloys: Seminar work. Ljubljana: Fakulteta za matematiko in fiziko, oddelek za fiziko, 2007. http://mafija.fmf.uni-lj.si/ seminar/files/2006_2007/SMA.pdf. 3. D. Stoeckel. The Shape Memory Effect – phenomenon, alloys and applications. Proceedings: EPRI, Freemont, Calofornia, 1995, p.1-13. 4. B. Coluzzi, A. Biscarini, L. Di Massoa, F.M. Mazzolai, N. Staffolanib, M. Guerrab, M. Santorob, S. Ceresara, A. Tuissi. Phase transition features of NiTi orthodontic wires subjected to constant bending strains. Journal of Alloys and Compounds, Volume 233, Issues 1-2, 15 January 1996, p. 197- 205. 5. J. Ferčec, M. Kos, M. Brunčko, I. Anžel, B. Glišić, E. Marković, R. Rudolf. Comparison of NiTi orthodontic archwires and a determination of the charasteristic properties. Materiali in tehnologije, ISSN 1580-2949, 2014, 48, 1, p. 99-104. 6. J. Ferčec, I. Anžel, R. Rudolf. Stress dependent electrical resistivity of orthodontic wire from the shape memory alloy NiTi. Materials & design, ISSN 0264-1275, Mar. 2014, vol. 55, p. 699-706. 7. J. Ferčec, B. Glišić, I. Šćepan, E. Marković, D. Stamenković, I. Anžel, J. Flašker, R. Rudolf, Determination of Stresses and Forces on the Orthodontic System by Using Numerical Simulation of the Finite Elements Method. Acta Phys Po A 2012;122:659-665. 8. V.R. Voller, C. Prakash, A Fixed-Grid Numerical Modelling Methodology for Convection-Diffusion Mushy Region Phase-Change Problems. Int. J. Heat Mass Transfer 30, 1709 (1987). 9. M. Gündüz, H. Kaya, E. Çadırlı, N. Maraşlı, K. Keșlioğlu , B. Saatçi, Effect of solidification processing parameters on the cellular spacings in the Al–0.1 wt% Ti and Al–0.5 wt% Ti alloys. J. Alloys Compd. 439, 114-127 (2007). 10. O. Turan, A. Sachdeva, N. Chakraborty, R.J. Poole, Laminar natural convection of power-law fluids in a square enclosure with differentially heated side walls subjected to constant temperatures, J. Non- Newt. Fluid Mech. 166, 1049-1063 (2011). 11. I. Biluš, M. Morgut, E. Nobile, Simulation of sheet and cloud cavitation with homogenous transport models. Int. J. Simul. Model. 12, 94-106 (2013).