M A T H E M A T I C A L M O D E L I N G OF H E A T TRANSFER DURING W A T E R SPRAY COOLING A N D C O N T R O L L E D S L O W COOLING OF S T E E L T U B E S . by Steve Mwenifumbo B . A S c , University of British Columbia, 1999 A THESIS SUBMITTED IN PARTIAL F U L F I L L M E N T OF T H E REQUIREMENTS FOR T H E D E G R E E OF M A S T E R OF APPLIED SCIENCE in T H E F A C U L T Y OF G R A D U A T E STUDIES D E P A R T M E N T OF M E T A L S A N D M A T E R I A L S ENGINEERING We accept this thesis as conforming to the required standard T H E UNIVERSITY OF BRITISH C O L U M B I A February 2003 © Steve Mwenifumbo In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of HeA-aXs aw^ M ^ W I Q A S The University of British Columbia Vancouver, B C , Canada Date tAoLv. G | ol Bto^cnrei/iV Abstract Mathematical models were developed which incorporate heat flow, phase transformation kinetics, and temperature and composition dependent thermal-physical properties to predict through-thickness temperature evolution during the manufacturing of steel tubular products under industrial cooling conditions. The two industrial cooling conditions investigated were an 'in-line' water quench, used to control austenite grain growth, and controlled slow cooling, applied at the end of the process to produce the desired rmcrosrructure and mechanical properties. This investigation combined computer modeling with laboratory measurements on industrial samples (i.e., AISI 5130 as-rolled steel tubes). A 2-D finite-difference axisymmetric accelerated cooling model (ACM) and a 2-D finitedifference circumferential controlled slow cooling model (CSC) were coded to run on a personal computer. The A C M model can handle multi-stage cooling with varying boundary conditions including radiation, natural convection, and forced convection. The model has the capacity to handle temperature and microstructural dependent thermal-physical properties associated with different phases that can form in the material. The C S C is capable of simulating the effects of tube spacing, forced air-cooling, and various ambient conditions. The basic heat diffusion parts of both the A C M and CSC models, including the routines associated with the liberation of latent heat, were validated against a benchmark commercial F E M software package, ABAQUS™. In addition, the results of the C S C model were verified experimentally by acquiring thermal data from a tube instrumented with thermocouples during the slow cooling of 5130 steel tubes. The CSC model predictions agreed reasonably well with the measured data, capturing the various factors affecting heat-transfer between tubes and their surroundings. A sensitivity analysis conducted with the C S C model indicated that the radiation exchange between different components (e.g., adjacent tubes, tubes and the fume hood, tubes and the floor) plays a significant role in heat-transfer. In addition, several factors influenced the effective emissivity within the model, including loose scale formation on the tubes, geometric assumptions in the model to approximate radiation view factors, and assumptions related to the temperature of the surrounding furnace environment. A B A Q U S is a trademark of Hibbitt, Karlsson and Sorensen Inc. ii Table of Contents Abstract ii Table of Contents iii List of Figures vi List of Tables x List of Symbols xi Acknowledgements 1.0 2.0 Introduction 1 References 4 Literature Review 6 2.1 Introduction 6 2.2 Numerical Methods 7 2.2.1 Finite-Difference Methods (FDM) 9 2.2.2 Alternate Direction Implicit Method (ADI) 12 2.2.3 Finite-Element Methods (FEM) 13 2.3 3.0 xv Boundary Conditions 14 2.3.1 Convection 14 2.3.2 Radiation 21 2.4 Thermal-Physical Properties 22 2.5 Phase Transformations 24 2.5.1 Ferrite/Pearlite 26 2.5.2 Bainite 27 2.5.3 Martensite 28 2.6 Process Modeling 29 2.7 Summary 31 References 34 Scope and Objectives 41 3.1 Objectives 41 3.2 Methodology 41 References 44 iii 4.0 Mathematical Modeling 46 4.1 General P D E and Solution Approach 47 4.1.1 Accelerated Cooling Model 48 4.1.2 Controlled Slow Cool Model 58 4.2 5.0 6.0 66 4.2.1 A C M Model 73 4.2.2 C S C Model 73 4.3 Coupling to Microstructure Models 74 4.4 Material Properties 75 4.5 Latent Heat of Evolution 76 4.6 Summary 77 References 80 Experimental Tests 81 5.1 Material 81 5.2 Controlled Slow Cooling Bed Simulator 81 5.3 Sample Instrumentation 84 5.4 Controlled Slow Cool Test Runs 85 5.4.1 87 Measured Cooling Curve Data References 93 Model Validation 94 6.1 Numerical Validation 94 6.1.1 Accelerated Cooling Model 95 6.1.2 Controlled Slow Cool Model 99 6.2 7.0 Implementation of Boundary Conditions Experimental Analysis 103 References 109 Discussion 110 7.1 Sensitivity Analysis 110 7.1.1 Stability and Convergence 112 7.1.2 Temperature Profiles 117 7.1.3 Fitting Parameters 127 7.2 Industrial Applications: 128 References 130 iv 8.0 Summary and Recommendations 131 8.1 Summary and Conclusions 131 8.1.1 Sensitivity Analysis - Solution Parameters 133 8.1.2 Sensitivity Analysis - Process Parameters 134 8.2 Recommendations for Future Work 135 Appendix A : Heat-Transfer Finite-Difference Equations-ACM Model (Interior Node) 137 Appendix B: Heat-Transfer Finite-Difference Equations - C S C Model (Interior Node) 140 Appendix C: Experimental Validation Cooling Curves - CSC Model 144 v List of Figures Figure 1.1 Schematic of the seamless tube making process. Figure 2.1 Schematic of thermo-mechanical processing model components and coupling. Figure 2.2 7 Schematic depicting heat flow of an axisymmetrical differential element in cylindrical coordinates. Figure 2.3 10 Schematic of ADI directions for solving the parabolic equations in an axisymmetric cylindrical coordinate case. Figure 2.4 12 Boiling curve representing the relationship between the heat flux and the surface temperature during boiling heat-transfer. Figure 2.5 3 20 Shape factor associated with radiation exchange between elemental surface of area dA[ and dAj. 22 Figure 2.6 Temperature dependent thermal-physical properties 23 Figure 3.1 Schematic of the methodology used in this study. 43 Figure 4.1 Schematics of different cooling systems: (a) water-quench stands and (b) slow cooling bed. Figure 4.2 46 Schematic depicting the relationship between temperature and microstructure in cooling process. 47 Figure 4.3 Schematic of axisymmetrical geometry - A C M model. 48 Figure 4.4 Schematic of ADI directions for Accelerated Cooling FD model. 50 Figure 4.5 Schematic of spatial discretization. 50 Figure 4.6 Schematic node numbering. 51 Figure 4.7 3-D and 2-D schematic of control volume. 51 Figure 4.8 Schematic of A C M model interior node and associated control volume. 52 Figure 4.9 Schematic of A C M model outer surface node and associated control volume. Figure 4.10 53 Schematic of A C M model inner surface node and associated control volume. Figure 4.11 54 Schematic of A C M model outer surface corner node and associated control volume. Figure 4.12 55 Schematic of A C M model inner surface corner node and associated control volume. 56 vi Figure 4.13 Schematic of A C M model surface end node and associated control volume. 57 Figure 4.14 Schematic of CSC geometry. 59 Figure 4.15 ADI implicit directions for Controlled Slow Cool F D model. 61 Figure 4.16 Schematic of spatial discretization and node numbering. 62 Figure 4.17 3-D & 2-D schematic of control volume with respect to spatial nodes. 62 Figure 4.18 Schematic of CSC model interior node and associated control volume. 63 Figure 4.19 Schematic of C S C model inner surface node and associated control volume. Figure 4.20 64 Schematic of CSC model outer surface node and associated control volume. 65 Figure 4.21 Spray boiling curves. 70 Figure 4.22 Schematic representing boiling curve regimes. 71 Figure 4.23 Schematic illustration depicting the view factor approximation. 72 Figure 4.24 Depiction of tube motion. 74 Figure 5.1 Schematic of controlled slow cooling bed system used to cool tubular specimens. 82 Figure 5.2 Photograph of the experimental test rig. 82 Figure 5.3 Photograph of the furnace hood and schematic of thermocouple positions. 83 Figure 5.4 Photograph of furnace hood pre-heat. 84 Figure 5.5 (a) Photograph of instrumented 5130 steel specimen depicting thermocouple positions, (b) schematic of corresponding CSC model nodal/thermocouple positions. 85 Figure 5.6 Photograph of experimental test run # 6 (without furnace hood). 86 Figure 5.7 Furnace hood thermocouple readings - experimental Test 1. 87 Figure 5.8 Measured cooling curves - Test 1. 88 Figure 5.9 Furnace hood thermocouple readings - experimental Test 2. 88 Figure 5.10 Measured cooling curves - Test 2. 89 Figure 5.11 Furnace hood thermocouple readings - experimental Test 3. 89 Figure 5.12 Measured cooling curves - Test 3. 90 Figure 5.13 Measured cooling curves - Test 4. 90 Figure 5.14 Furnace hood thermocouple readings - experimental Test 5. 91 Figure 5.15 Measured cooling curves - Test 5. 91 Figure 5.16 Measured cooling curves - Test 6. 92 vii Figure 6.1 Schematic of accelerated cooling processing conditions (cases 1 & 2). Figure 6.2 Case 1: Comparison of A C M model to A B A Q U S F E M model - no latent heat evolution. Figure 6.3 98 Case 2: Comparison of A C M model to A B A Q U S F E M model - latent heat evolution. 99 Figure 6.4 Schematic of controlled slow cool processing conditions (cases 1 & 2). Figure 6.5 Case 1: Comparison of CSC model to A B A Q U S F E M model - no latent heat evolution. Figure 6.6 101 101 Case 2: Comparison of CSC model to A B A Q U S F E M model -latent heat evolution. Figure 6.7 97 102 Predicted cooling curve and experimental data - 5130 Test 4 (without furnace hood). 105 Figure 6.8 Microstructure evolution - 5130 Test 4 (without furnace hood). 106 Figure 6.9 Predicted cooling curve and experimental data - 5130 Test 5 (with furnace hood). 107 Figure 6.10 Microstructure evolution - 5130 Test 5 (with furnace hood). 108 Figure 7.1 Comparison of predicted cooling curves for an outer surface node on a tube experiencing a complex cooling regime. Figure 7.2 Comparison of ATmax values for an outer surface node experiencing an external quench. Figure 7.3 114 Simulated cooling curves of an outer surface node experiencing natural convection for several different mesh configurations. Figure 7.4 119 Simulated through-thickness cooling curves for an externally quenched tube traveling at 0.05 m/s ( A C M model). Figure 7.8 118 Simulated through-thickness cooling curves for an externally quenched tube traveling at 0.02 m/s ( A C M model). Figure 7.7 116 Predicted through-thickness cooling curves for tubes of varying wall thickness experiencing an external quench. Figure 7.6 115 Enlarged view of simulated cooling curves of an outer surface node experiencing natural convection for several different mesh configurations. Figure 7.5 113 120 Simulated through-thickness cooling curves for an externally quenched tube traveling at 0.1 m/s ( A C M model). viii 120 Figure 7.9 Simulated through-thickness cooling curves for an externally quenched tube traveling at 0.15 m/s ( A C M model). 121 Figure 7.10 Predicted cooling curves for internal/external quench. 122 Figure 7.11 Comparison of predicted cooling curves of an external node submitted to Figure 7.12 an external quench with and without phase transformation kinetics. 123 Predicted cooling curves of an outer surface node undergoing natural 124 convection at various inter tube spacing values. Figure 7.13 Effects of furnace temperature on an outer surface node under slow cooling conditions. 125 Figure 7.14 Comparison of tube end effects using the A C M model. 126 Figure 7.15 Comparison of emissivity values for predicted cooling curves. 127 Figure A.1 Schematic of ADI directions. 137 Figure A.2 Schematic of A C M model interior node and associated control volume. 138 Figure B.l Schematic of ADI directions. 140 Figure B.2 Schematic of C S C model interior node and associated control volume. 141 Figure C l Predicted cooling curve and experimental data - 5130 Test 1 (furnace hood). Figure C.2 144 Predicted cooling curve and experimental data - 5130 Test 2 (furnace hood). Figure C.3 145 Predicted cooling curve and experimental data - 5130 Test 3 (furnace hood). Figure C.4 146 Predicted cooling curve and experimental data - 5130 Test 4 (without furnace hood). Figure C.5 147 Predicted cooling curve and experimental data - 5130 Test 5 (furnace hood). Figure C.6 148 Predicted cooling curve and experimental data - 5130 Test 6 (without furnace hood). 149 ix List of Tables Table 2.1 Approximate values of convection heat-transfer coefficients. 16 Table 2.2 Constants C and n for natural convection on a horizontal cylinder. 16 Table 2.3 Values of C and n for flow over a horizontal cylinder. 18 Table 2.4 Values of C and m for flow normal to a single circular cylinder. 18 Table 4.1 Thermodynamic properties of dry air at atmospheric pressure. 68 Table 4.2 Heat-transfer correlations for spray boiling data. 69 Table 4.3 Thermal-physical property data for an AISI 5130 low alloy steel. 75 Table 4.4 Enthalpies of transformation in steel. 77 Table 5.1 Nominal composition of AISI 5130 alloy steel 81 Table 5.2 As-rolled 5130 alloy steel test matrix. 86 Table 6.1 Steel properties used in the numerical validation of the A C M and C S C models. Table 6.2 95 Case 1: Processing parameters used in A C M validation (thermal diffusion model - no latent heat). Table 6.3 Case 2: Processing parameters used in A C M validation (thermal diffusion model - latent heat). Table 6.4 96 Case 1: Processing parameters used in CSC validation (thermal diffusion model - no latent heat). Table 6.5 96 100 Case 2: Processing parameters used in CSC validation (thermal diffusion model - with latent heat). 100 Table 6.6 As-rolled 5130 alloy steel test matrix. 103 Table 6.7 CSC model parameters. 103 Table 6.8 Material properties for 5130 alloy steel. 104 Table 6.9 Thermodynamic properties of dry air at atmospheric pressure. 104 Table 7.1 Sensitivity analysis: solution parameters. 111 Table 7.2 Sensitivity analysis: processing parameters. 111 Table 7.3 Sensitivity analysis: boundary conditions for htmutai comparisons. 112 X List of Symbols Chapter 2 A„ Koistinen-Marburger coefficient b Kinetic parameter representing nucleation and growth rates C, m, n Empirical constants for convection over a horizontal cylinder c Specific heat J/kg/K D Diameter of the rod or tube m d Prior austenite grain size p r Internal energy J/kg/m g Gravitational acceleration m/s Gr Grashof number h Heat-transfer coefficient W/m /K hov Overall heat-transfer coefficient W/m /K h T Heat-transfer coefficient at radiation boundary W/m /K h c Heat-transfer coefficient at convection boundary W/m /K Hi Enthalpy of transformation of phase i J/kg k Thermal conductivity W/m/K M Martensitic transformation start temperature °C Mf, Martensite finish temperature °C s rn, n 2 2 2 2 Spatial positions Nu Nusselt's number, hDolk* Pr Prandlt's number, p c^k p Temporal dependence Q Latent heat of transformations J/kg/m q Heat flux W/m q Heat flux - radial direction W/m q Heat flux - axial direction W/m Re Reynolds number, [p po(l + Ra Rayleigh number r Radius m T Temperature K T Surface temperature K r z s a x) A/lvRTWT'a) xi 2 Tomb Ambient temperature K T Mean film temperature K Too External/environmental temperature K t Time s m X Volume fraction of phase i Xf Volume fraction of ferrite t Xf Temperature dependent equilibrium volume fraction of ferrite X Volume fraction of pearlite z Axial length >eq p V m Differential operator. p Density kg/m ^ Angle in the circumferential direction, radian V Average velocity m/s H Absolute viscosity N-s/m B Spatial coefficient of propagation yi, Influence parameter cr Stefan Boltzmann's constant J/K /m /s v Kinematic viscosity of air N-s/m e Emissivity 4 xii 2 Chapter 4 Bi Biot number, hR/ks c Specific heat J/kg/K D Outer diameter m E Stored energy J/kg/m p 2 F Fourier number g Gravitational acceleration m/s h Heat-transfer coefficient W/m /K h Heat-transfer coefficient for convective heat flow W/m /K hov Overall heat-transfer coefficient W/m /K Heat-transfer coefficient for radiation W/m /K Enthalpy of formation J/kg Hi Enthalpy of transformation of phase i J/kg k Thermal conductivity W/m/K Al Characteristic length m Pr Prandtl's number, \i c^k Q Latent heat generated per unit volume W/m 2 qoD Energy entering or leaving the node outside diameter W/m 2 qm Energy entering or leaving the node inside diameter W/m 2 r Radius m T Temperature K Tamb Ambient temperature K Ti Initial temperature K T Mean film temperature K Tnua Maximum temperature K T Surface temperature K T External/environmental temperature K Vfmax Maximum volume fraction At Time increment s AV Volume increment L 0 c m s x 2 2 2 2 xiii 2 Ax Mid-plane distance m Xi Volume fraction of phase i V Differential operator p Density P Spatial coefficient of propagation v Kinematic viscosity of air N-s/m p Absolute viscosity N-s/m e Emissivity a Stefan Boltzmann's constant A Austenite B Bainite F Ferrite M Martens it e P Pearlite kg/m J/K /m /s 4 xiv 2 Acknowledgements I would like to express my sincere thanks to my supervisors Dr. Mary A . Wells, Dr. Steve Cockcroft, and Dr. Daan Maijer for their guidance, support and patience during my time at U B C . In addition, I would like to acknowledge the Timken Bearing and Steel Company, Ohio, and the U.S. Department of Energy for funding this project. I am indebted to a host of people at U B C , in particular Gary Lockhart and Dianfeng L i for their technical assistance and help during the experimental tests. In addition, I would like to thank the staff at Timken, in particular Jeff Ives, Daqing Jin, and Robert Kolarik, for informative discussions on model development, facilitation of industrial trials, and project management. Much appreciation is given to all the people in the Department of Metals and Materials Engineering at U B C who have supported me during this project. Specifically, thanks are given to all the other graduate students for the many stimulating discussions. Finally, many thanks to my good friend Kira Procter and my family for the encouragement and support they provided when it was most needed. xv Chapter 1: Introduction 1 CHAPTER 1 INTRODUCTION Increasingly, manufacturing industries are turning to mathematical modeling to enhance product quality, to improve processing efficiency, and to increase the cost effectiveness of the manufacturing process. Steel companies are seeking to achieve these goals through the use of microstructural engineering models and the application of controlled thermo-mechanical processing (CTMP) to produce hot rolled products. [17] Mathematical models, when applied to the manufacturing process, link the basic principles of microstructural phenomena with heat and mass transfer. In order to calibrate these models to operating parameters, empirical and semiempirical relationships must be developed through laboratory experiments and industrial trials. [8] The C T M P of tubes and pipes involves the conversion of materials into value-added products through precise control of processing variables. The properties of metals and their alloys depend strongly on their processing history. For example, the distribution of phases, grain structure, alloy compositional segregation and defects in final commercial products depend strongly on the conditions under which the materials were fabricated. Precise control is required of thermal history, mechanical deformation and microstructure evolution during processing, if desired structures and properties are to be achieved. Recent advances in computer technology, sensor design, and numerical model simulations of materials and processes have allowed the application of C T M P to the production of steel tubes and pipes. Some anticipated benefits of CTMP include: • Reduced energy consumption. • Reduced emission ofgreenhouse gases and toxic wastes. • Reduced machining costs. • Reduced rework and scrap accumulation. • Increased overall cost savings. • Increased product performance. The Timken Company of Canton, Ohio, has recently carried out a development program to fully exploit the latest developments in both production and material technologies for seamless tube products. This five year program of research is directed towards the development of a Chapter 1: Introduction 2 system for controlled thermo-mechanical processing manufacturing and performance.' ' 9 of tubes and pipes for enhanced The program has been sponsored by the Department of Energy (DOE) and includes a number of industrial and institutional participants: The Timken Company, The University of British Columbia (UBC), Oak Ridge National Laboratories (ORNL), Idaho National Environment & Engineering Laboratories (INEEL), The National Research Council of Canada (NRC-IMI), US Steel, The Ford Motor Company, Daimler Chrysler Corporation, and The Friedrich Kocks Company. The methodology being employed for this program includes the development of a mathematical model to simulate an entire tube mill. The simulation will include components of heat-transfer and deformation for the cooling, heating, and rolling processes, as well as provide the expected metallurgical responses, based on the evolution of microstructure. One of the components of this program will use numerical modeling techniques to develop the Hot Tube Mill Model (HTMM) to optimize a virtual tube mill where the evolution of temperature, stress/strain, and microstructure can be simulated as they occur in the tube making process. In addition, the H T M M will be used to develop recipes for various steel grades and products. Overall, the objective is to yield a final product with enhanced properties while minimizing processing and energy costs. Two of the largest sectors in the domestic tube and pipe industry are seamless mechanical tubing (e.g., drive shafts, bearings) and seamless oil country tubular goods. [10] These products can be produced from a variety of steel grades, including 52100, 5130, 4320, and 8620, with seamless tube outer diameters (OD) ranging from 25.4mm to 154mm, and wall thickness ranging from 2.54mm to 50.8mm. Figure 1.1 depicts the proposed seamless tube production process at Timken where a billet enters the process, goes through seven stages, and exits as a finished product. There are two novel stages in the proposed Timken C T M P process. These include: 1) an in-line' heat treatment, and 2) controlled slow cooling. The 'in-line' heat treatment consists of an inductive heater or water quench stand (or a combination of both), and occurs after the product leaves the elongation mill. This treatment has a dual purpose. It removes any temperature profile variations developed from previous processing steps (i.e., conduction from the mandrel when first entering the piercing mill) by heating or cooling the tube. It also controls austenite grain growth through the use of a water quench stand. The controlled slow cooling stage in the C T M P process consists of a controlled slow cool bed located at the end of the process. The purpose of Chapter I: Introduction 3 the slow cool bed is to ensure the final product has the desired microstructural and mechanical properties. Rotary Hearth Furnace Piercing Mill Elongation Mill r^lM) Induction Heater Water Quench Stand Reducing Mill Rotary Sizing Mill 3 Slow Bed Cool Figure 1.1 - Schematic of the seamless tube making process. In order to model both the 'in-line' water quench stand and the controlled slow cooling bed, it is necessary to accurately describe the appropriate boundary conditions. These include boiling water heat-transfer associated with the 'in-line' water quench stands and the radiation view factors between adjacent tubes and their surroundings, associated with the controlled slow cooling bed. Another factor to consider is the latent heat that evolves when the austenite decomposes into various transformation products and its effect on the tube thermo-profile. In this research project, a set of mathematical models for simulating water quenching and controlled slow cooling of tubular products for the Timken C T M P process have been developed. The models will be integrated into the H T M M and be utilized by INEEL to populate a neural net that will be used for control purposes in the C T M P process. These models provide important insights into the tube manufacturing process leading to cost benefits and quality enhancements in tubular products. Chapter 1: Introduction 4 References 1. I. V. Samarasekera, D. Q. Jin and J. K. Brimacombe, "The Application of Microstructure Engineering to the Hot Rolling of Steel", Proceedings of 38 Mechanical Working and Steel th Processing Conference, Cleveland, O H , Vol. 34, 1996, pp 313-327. 2. E. Anelli, P. Farsetti, and G. Cumino, "New Process for On-line Normalizing of Seamless Pipes", Proceedings of 38 th Mechanical Working and Steel Processing Conference, Cleveland, O H , Vol. 34, 1996, pp 285-295. 3. C. A. Muojekwu, D. Q. Jin, V . H. Hernandez, I. V. Samarasekera and J. K. Brimacombe, "Hot-direct Rolling, Runout Table Cooling and Mechanical Properties of Steel Strips Produced from Thin Slabs", Proceedings of 38 Mechanical Working and Steel Processing th Conference, Cleveland, O H , Vol. 34, 1996, pp 351-366. 4. C. Dilg, J. Kirsch, W. Schutz, E . Amoris, and A. Streibelberger, "TMCP Process Including Two Accelerated Cooling Stages", La Revue De Metallurgie-CIT, 1995, pp 884-892. 5. E . Anelli, P. Farsetti, and G . Cumino, "Prediction of the Response to Quenching and Tempering of Low Allow Steels for Seamless Tubes", Proceedings of 36 Mechanical Working and Steel Processing Conference, Baltimore, M D , Vol. 32, 1994, pp 149-159. 6. E . Anelli, "Application of Mathematical Modeling to Hot Rolling and Controlled Cooling of Wire Rods and Bars", ISIJ International (Japan), Vol. 32, No.3, 1992, pp 440-449. 7. A. Kumar, C. McCulloch, E . B. Hawbolt, and I. V . Samarasekera, "Modeling Thermal and Microstructural Evolution on Runout Table of Hot Strip Mill", Materials Science and Technology, April 1991, V o l . 7, pp 360-368. 8. P. C. Campbell, E.B. Hawbolt, and J. K. Brimacombe, "Microstructural Engineering Applied to the Controlled Cooling of Steel Wire Rod: Part I. Experimental Design and Heat Transfer", Metallurgical Transactions A , Vol. 22A, No. 11, 1991, pp 2769-2778. 9. R.V. Kolarik, "Controlled Thermo-mechanical Processing of Tubes and Pipes for Enhanced Manufacturing and Performance", 1999, Timken internal report. Chapter 1: Introduction 5 10. "Controlled Thermo-Mechanical Processing (CTMP) of Tubes and Pipes for Enhanced Manufacturing and Performance", Office of Industrial Technologies, Energy Efficiency and Renewable Energy, U.S. Department of Energy, April 2000, Steel Project Fact Sheet. Chapter 2: Literature Review 6 CHAPTER 2 LITERATURE REVIEW 2.1 Introduction With increasing emphasis on product quality and process optimization, industry and academia have strived to relate the structure and mechanical properties of quenched and tempered steel products to heat treatment variables and steel composition. A common problem associated with thermo-mechanical processing is the development of residual stresses and distortion in the final product. This distortion and residual stress often make it difficult, and even impossible, to manufacture components to final product dimensions prior to heat treatment. Instead, an extensive amount of costly machining is required to remove excess material after hardening. In addition, removal of hardened surface material can often lead to a decrease in the mechanical properties of the products. 111 In order to minimize the impact of distortions and post- heat treatment machining, knowledge of the thermal history and microstructural evolution becomes crucial. In the past, practical experience was used to predict the final product properties. Unfortunately, these traditional methods have turned out to be inefficient and insufficient, because small changes in the thermo-mechanical processing can lead to major changes in final product properties. As a result, researchers are working to develop mathematical models capable of quantitatively predicting the evolution of temperature, deformation and microstructure during thermo-mechanical processing. These models may then be used to optimize the processing from the standpoint of component mechanical properties (including ease of machining) and dimensional tolerances. Ideally a comprehensive thermo-mechanical process model would include treatment of: 1) the heat-transfer processes, 2) the phase transformation processes, and 3) the deformation processes, and would also account for the interactions between these various phenomena (see Figure 2.1). This literature review will concentrate on the work performed to model the thermal fields and microstructure evolution, as it applies to hot deformed steel rod and tubular products, and exclude components associated with mechanical responses, as the latter are not the subject of this research program. The analysis of heat flow must be coupled to microstructural evolution Chapter 2: Literature Review 7 since there is a release of latent heat associated with the transformation event, and the final microstructure of the product is dependent on the thermal history experienced by the component. Induction heating Hoat transfer coefficient Temperature analysis Latent heat TTT diagram Composition Thermal strain Mechanical response Stress dependent transformation Phase transformation Transformation strain Transformation induced plasticity Figure 2.1 - Schematic of thermo-mechanical processing model components and coupling. 111 In the following sections, this review will examine: numerical methods used in mathematical modeling; the modeling of boundary conditions as they relate to 'in-line' quenching and controlled slow cooling of round products; material property dependencies; microstructure evolution as it relates to diffusion-controlled and martensitic type phase transformations; and, finally some of the research investigating thermal process modeling. 2.2 Numerical Methods Numerical simulations can provide detailed thermal and microstructural response data for a particular heat treatment, as well as being useful tools in improving the understanding of heattransfer mechanisms active in a process. The formulation of mathematical models describing controlled thermo-mechanical processing of tubes requires a set of partial differential equations (PDEs) to describe the heat-transfer encompassing heat conduction, heat accumulation and generation and/or consumption. For transient heat conduction analysis, the governing PDE is: t2J Chapter 2: Literature Review 8 V(k-VT) + Q = p c ^ dt (2.1) p where Vis the differential operator, k is thermal conductivity (W/m/K), T is temperature (K), Q is latent heat generated per unit volume associated with phase transformations (W/kg/m ), p is 3 density (kg/m ), c is specific heat (J/kg/K), and r is time (s). In the case of rods or tubular 3 p products the general PDE in cylindrical coordinates is: \_d_f yrr dr ury dr J d T^ 2 r d(/> 2 2 dz 2 [3] dT + Q = pc — p (2.2) J where r is the radius (m), <j> is the angle in the circumferential direction (radians), z is the axial length (m), and the thermal conductivity, k, is temperature independent. In analytical methods, the conducting body is treated as a continuum, where the temperature may be ascertained at any point or time through the successful solution of the PDE. ^ This enables the deduction of a large amount of practical information (i.e., the effect of [4 various parameters on temperature, etc..) owing to the closed and analytical nature of the results. Only a limited number of relatively simple geometrical shapes (e.g. cylinders, spheres, etc..) can be handled analytically. In addition, analytical solutions are limited to the application of boundary conditions that are typically very simple. Unfortunately, most heat conduction problems associated with thermo-mechanical processing possess no real analytical solution due to the presence of non-linearities. There are many sources for non-linearity in heat-transfer problems (temperature, time, and geometry dependence, etc..) that are generally associated with: • Material Properties - density (p), thermal conductivity (k), emissivity (e), specific heat (c ), etc... p Internal Heat Generation - latent heat, plastic deformation, etc... Boundary Conditions: a. Convection (temperature or geometry dependent) b. Radiation (function of view factors, emissivity, texture, T to the fourth power, etc..) Consequently, this absence of an analytic solution requires another approach to solve the PDE. Chapter 2: Literature Review 9 A variety of numerical methods can be utilized to solve the partial differential equations that describe the evolution of the temperature and microstructure. O f these, the numerical methods most suitable for computer programming are finite-difference approximation techniques and finite-element techniques. This section will introduce these techniques as well as some of [5] their limitations. 2.2.1 Finite-Difference Methods (FDM) The finite-difference method has been one of the most widely used numerical methods for decades. This method provides an approximate solution by dividing the continuous domain into a network or mesh of discrete points. The field unknowns are sought only at these discrete points rather than everywhere in the domain. [6] The basic principle of this numerical approach is to replace the derivatives appearing in the PDE with a series of piece-wise approximations. For example, the temporal derivative in Equation (2.2) may be approximated by Equation (2.3): j>p+\ _ rj,p 8T_ dt (2.3) At where the subscripts m and n represent the spatial positions and p the temporal dependence, and the finite-difference equations of the spatial second derivative terms are expressed as (explicit scheme will be discussed in a forthcoming section): 1 dT 2 8T j,p _ yrp , jp m+l,n ' m,n ' m-\,n 1 z dr 8T dz (Az) 2 1 1 (Ar) rpp _ 2fP m,n+\ m,n 2 1 I fP m,n-\ 2 1 K -2T^ dT 2 2 + Tf n+l (A*) < d<f> 2 2 Kumar et al. [7] . 1 r 1 m jp m+l,o _ 1 jip m-\,n 2(Ar) (2.4) (2.5) (2.6) and others used a conservation of energy approach to develop an energy [8] balance for differential elements. In this method, also known as the control-volume or heat balance approach, the finite-difference equation for any node may be obtained by formulating an energy balance at each node (see Figure 2.2). Chapter 2: Literature Review 10 Differential Element Primary Node Aq (m+l,n) Xm,n) z (m.rrt-1) ®—-Aq z + A z (m-lin) Aq r Figure 2.2 - Schematic depicting heat flow of an axisymmetric differential element in cylindrical coordinates. At time t, the rate at which energy is transferred to the environment or a neighboring node minus the rate at which energy leaves this node by conduction must be equal to the rate of change of internal energy for the node. For the case shown in Figure 2.2, this may be expressed as: + E.AV (2.7) where q is the energy entering or leaving the node, r and z are the subscripts referring to the radial and axial directions respectively, Q is the internal energy generation term, AV is the volume, and is the internal energy storage term. There are a variety of different schemes that can be used to solve the finite-difference approximations. In the explicit schemes the spatial derivative is expressed in terms of the temperatures at the current time step. Hence, the determination of a nodal temperature at some point in time is independent of temperatures at other nodes for the same time. Thus, if at the end of a certain time period, all the nodal temperatures are known, then each of the nodal temperatures at the end of the next moment may be explicitly found, node by node. Davis has [9] shown that the explicit method is both convergent and stable under the following criterion (twodimensional case - Cartesian coordinates): At <8 k (2.8) Although adhering to this criterion alleviates any instability, the selection of spatial and temporal increments places strong limitations on the explicit methods. For a given spatial Chapter 2: Literature Review 11 increment, the time increment must be compatible with the stability requirements. For example, decreasing the mesh size to improve the approximation of the spatial derivative will require a power law increase in the number of time increments and thus an increase in computational time. By employing an implicit scheme, rather than an explicit one, a reduction in the amount of computational time may be realized. In implicit schemes the spatial derivative is approximated in terms of the temperatures at the future time increment.^ 101 When this is substituted into the original PDE, it results in a difference equation that contains several unknowns. Therefore, since the future nodal temperature is expressed in terms of its current value and the future values of its neighbors' temperatures, the progression from one time step to the next requires a system of equations to be solved simultaneously. Despite the fact that this scheme utilizes somewhat more complicated algorithms, it eliminates some of the stability and convergence problems that occur with explicit schemes. Consequently, the foregoing restrictions on the spatial and temporal increments are eliminated and larger values of At may generally be used, resulting in reduced computational times. In addition, implicit schemes can often be set up to form tri-diagonal matrices. A tri- diagonal matrix is a square matrix that has all elements equal to zero, with the exception of a band centered on the main diagonal (bandwidth equal to 3). These tri-diagonal systems can be formulated to eliminate the storage of large numbers of useless zeros in the matrix. This spacesaving modification is advantageous because the resulting algorithm requires less computer memory, and may be solved in an efficient manner. Unfortunately, the simple implicit schemes are limited in that they are not second-order. An alternative scheme that is second order accurate in both space and time is the Crank-Nicolson method. [I1] Within this scheme the finite-difference approximations are developed at the midpoint of the time increment. In this scheme, the temporal derivative is approximated at t+At/2 and the spatial derivative is determined at the midpoint by averaging the difference approximations at the beginning (t) and at the end (t+At) of the time increment. [n] Although implicit techniques offer alternatives that guarantee stability, the direct application of implicit methods, such as the Crank-Nicolson technique, leads to the solution of m x n simultaneous equations which result in increased computation times. In addition, in two or three dimensions, matrix storage and computation time can become extremely large owing to the Chapter 2: Literature Review 12 fact that the matrix, in general, ceases to be tri-diagonal. [10] The alternating direction implicit method (ADI) is one way around this dilemma. 2.2.2 Alternating Direction Implicit Method (ADI) The ADI method, first proposed by Peaceman and Rachford represents an additional [12] method for solving the set of implicit finite-difference equations derived from transient parabolic equations. It has been shown to be approximately seven times faster than the Crank-Nicolson method and provides a means for solving parabolic equations in two spatial dimensions using tridiagonal matrices' . 121 In this method, the time increment is separated into two steps. illustration of this scheme, as proposed in Canale and Chapra, [10] An for an axisymmetrical heat conduction problem involving cylindrical coordinates is presented in Figure 2.3. r ! -•Z (b) 2 direction (a) Indirection n d Figure 2.3 - Schematic of ADI directions for solving the parabolic equations in an axisymmetric cylindrical coordinate case.' 101 In the 1 half-time step, the nodal temperatures are solved implicitly in the axial direction and explicitly in the radial direction using Equation (2.9): (JP+\I2 _2pp+-y2 , jp+yi jp N In the 2 nd _yrp 4-T p 1T 9 —T p 2M 2 ^ TP+V2 _JP (2.9) Ar half-time step, the nodal temperatures are implicitly solved in the radial direction and explicitly in the axial direction using Equation 2.10: m*\n m,n (Az) 2 m-\n nyn-1 m/i m,n-\ _,_ *• m,/H-1m,n-l m,n- 2(*) nnp+1 rpp*AJ2 (2.10) Chapter 2: Literature Review 13 Owing to its advantages of unconditional stability, second order accuracy, and a cost efficient solution algorithm employing tri-diagonal matrices for two dimensional transient heat conduction problems/ ' ' the ADI method is rather attractive in comparison to many of the 11 13 other finite-difference techniques. 2.2.3 Finite-Element Methods (FEM) In contrast to finite-difference techniques, the finite-element method divides the solution domain into a series of sub-domains called elements where the partial differential equations are reduced to a finite system of ordinary differential equations, which are then solved by matrix solution techniques. Originally developed to calculate stress in irregularly shaped objects and analyze structural problems in aircraft, 161 F E M has recently become a dominant numerical method for the simulation of thermal and stress problems. The reduction of the governing equation is typically performed using either the RayleighRitz method or the method of weighted residuals (MWR). The M W R incorporating a Galerkin procedure is the more common of the two, where element matrix equations are formulated by evaluating terms from general integral equations using elemental interpolation functions. 1141 This scheme is popular because of its general approach that permits a functional form of the dependent variables to be obtained for any equation regardless of its complexity/ 151 Although, the F E M is a powerful technique that can be applied to a wide variety of problems, memory storage issues can quickly become an issue for complex problems. In comparison, FDMs are simple to formulate, can easily be extended to two or three-dimensions, and require considerably less computational work and memory storage than FEMs for an equivalent number of nodes. with the user. Ultimately, the choice of which particular technique to use lies In general, simpler techniques are often used to solve geometrically simple problems, with the sophistication of the technique increasing with problem complexity. study performed by Thomas et a/. [131 In a for heat conduction problems with simple geometry, the ADI finite-difference method coupled with a specific latent-heat technique was found to be more cost effective and stabile, with accuracy similar to the finite-element methods. Chapter 2: Literature Review 14 2.3 Boundary Conditions In order to develop accurate thermal models, it is critical to have boundary conditions that accurately describe a process. The boundary conditions at the surface (t > 0) may be characterized in one of several different ways: Dirichlet type (specified temperature): q = -k^- = T(n ,n ,n ,t) on l 2 (2.11) J Neumann type (specified heat flux): q = -k—- = q(n ,n ,n ,t) on [ 2 (2.12) 3 Cauchy type (specified heat-transfer coefficient): q = -k^ on = h( ,n ,n ,tlT -Tj n] 2 3 (2.13) s Where q is the heat flux, — is the outward pointing derivative normal to the surface, and h is dn the heat-transfer coefficient. During the processing of rod and tubular products, the boundary conditions of convection (natural, forced-air, and spray-water) and radiation play a significant role. 2.3.1 Convection In the convective cooling of steel rod and tubular products, cooling occurs when fluid (gases or liquids) contacts the surface. The flow is typically treated as being normal to a single cylinder. To mathematically describe the process of heat convection, the basic laws of heat conduction must be coupled with those of fluid motion/ ' 5 The aforementioned process is sufficiently complex that an analytical description is virtually impossible. In addition, it is usually the average heat-transfer coefficient, rather than the local, that is desired for engineering applications. As a result, numerous studies have been performed to develop correlations for convective heat-transfer from cylindrical bodies in cross-flow conditions for a range of different Chapter 2: Literature Review fluids. 15 Most are empirical in nature and relate the Nusselt's number, Nu, to the Reynolds 116-201 number, Re, and Prandlt's number Pr (see Equations 2.14-2.16). Nusselt number: Nu =^D (2.14) k where h is the heat-transfer coefficient, D is the diameter of the rod or tube, and k the thermal conductivity. Reynolds number: R e „ = ^ (2.1S) where V is the average velocity, p the density, and p the absolute viscosity. Prandd's number: Pr = - ^ - (2.16) where c is the specific heat. One of these equations, as given by Krieth and Black/ is: 31 p Nu = CRe Pr x (2.17) y where C, x, and y are constants that depend on the magnitude of Re and on the cooling medium. In addition to these dimensionless parameters, the Rayleigh number, Ra, is often used for natural convection: Ra D = Gr Pr = ^ P ^ T ' ~ J ^ ' Pr D (2.18) where B is the spatial coefficient of propagation, g the gravitational acceleration, and Gr is the Grashof number. Chapter 2: Literature Review 16 The heat-transfer coefficient is a complex function of the composition of the fluid, the geometry of the solid surface, and the dynamics of the fluid motion past the surface. The typical range of values one might expect to encounter for h are quoted in Table 2.1. Table 2.1 - Approximate values of convection heat-transfer coefficients.' ' 3 Convection Mode and Fluid h (W/m'/K) Natural convection, air Forced convection, air Convection (natural, forced, boiling), water 5-25 10-200 20-100000 c Natural Convection In natural convection, flow around the body results from the development of temperature and density gradients. These gradients are the result of the heat-transfer between the body and the surrounding air.' To describe this behavior around cylindrical bodies, McAdams' ' as well 5] 20 as Morgan' ' suggested a correlation for the heat-transfer coefficient: 19 h = ^(CRa" ) (2.19) D where the heat-transfer coefficient and the Rayleigh number are based on the cylinder diameter and the values of C and n are presented in Table 2.2. Table 2.2 - Constants C and n for natural convection on a horizontal cylinder. Range of Rao io- -io10~ - 10 10 -10 10 -10 10 - 10 u j 2 2 2 4 4 7 7 12 c n 0.675 1.02 0.850 0.480 0.125 0.058 0.148 0.188 0.250 0.33 • T211 Golja 1 adopted a modified version of McAdams and Morgan's correlation for natural convection from a horizontal cylinder. It was postulated that since the Prandlt number does not change significantly with the air temperature, a constant value of 0.7 < Pr < 1 can be used and consequently the horizontal cylinder heats the surrounding air with a heat-transfer coefficient: h= _ 0.38 3/%_( ^ 2\ s \ amb) r D D T T l V l V 4 (2.20) Chapter 2: Literature Review 17 where v is the kinematic viscosity of air. Churchill and C h u , in another extensive study, derived analytically the following [18] relation for natural convection from a horizontal cylinder that is valid for a wide range of Rao'- r0.559V ' h= 0.6 + 0.387Ral 1 + D Pr J /16 -8/27 A 6 (2.21) where 0 <Pr < ao, JO' < Rao < 10 , and the properties are taken at the mean film temperature 5 12 T The choice of which equation to use for cross-flow over cylinders is subject to some conjecture, however, equations applicable to a wide range of flows are more practical from a computer standpoint. Forced-Air Convection Forced-air cooling has been used by numerous steel companies to achieve higher cooling rates than natural convection during the processing of steel rod and tubular products (e.g. Stelmor line). Because of the complexities associated with flow over cylinders, emphasis has been placed on developing empirical correlations to characterize a mean heat-transfer coefficient. Kreith and Black [3] developed an equation to characterize a mean heat-transfer coefficient in terms of the Reynolds number and the Prandtl number for flow over a horizontal cylinder: h = — CRe" Pr D ° 1/3 (2.22) where all the properties of air are calculated at the mean film temperature and the values of C and n are listed in Table 2.3 for various ranges of Reo- Chapter 2: Literature Review 18 Table 2.3 - Values of C and n for flow over a horizontal cylinder. Range of Re 0.4-4 4-40 40 - 4000 4000 - 40000 40000 - 400000 c 0.989 0.911 0.683 0.193 0.0266 D n 0.330 0.385 0.466 0.618 0.805 A study to measure the thermal response of an instrumented steel rod under controlled cooling conditions was performed by Campbell et alP ^ in the laboratory and on an operating 2 Stelmor line. Their results indicated that the correlation, published by Kreith and Black , [3] predicted reasonably well measured experimental values. Zhukauskas [16] performed a comprehensive study of flow normal to single tubes and banks of tubes. For single tubes the following correlation, used in conjunction with the values listed in Table 2.4 along with the mean film temperature, was developed: fp A/4 r h= CRe Pr" (2.23) m D D Pr. where 0.7 <Pr < 500, 1 <Re < 10 , n = 0.37 for Pr < 10, n = 0.36 for Pr > 10, the properties 6 D are calculated at the mean film temperature T , except for Pr (7/= 7^) and Pr (T=T ), and C and m x s S m are from Table 2.4. Table 2.4 - Values of C and m for flow normal to a single circular cylinder. Range of Reo 1-40 40- 1000 10 -2xl0 2xl0 - 10 3 5 5 6 Churchill and Bernstein 1 J C m 0.75 0.51 0.26 0.076 0.4 0.5 0.6 0.7 proposed the following correlation for forced-air convection over a horizontal tube which is based on the analysis of many data covering wide ranges of Pr and Reo'. Chapter 2: Literature Review 19 0.62Re£ Pr 2 0.3 + 1/3 1/4 1 + 2.82 x l O ) (2.24) 5 J where the properties are calculated at T and for all ReoPr > 0.2. m As with the correlations developed for natural convection, relations that cover a large range of Re are often better suited for modeling a variety of cooling conditions. Spray Water Cooling - Boiling Water Heat-transfer Spraying water on hot metal surfaces is a common cooling procedure and has been used extensively in a variety of manufacturing processes ranging from hot deformation processes, where the water is used to control microstructure, to heat-treating, where the water is used to cool parts down to room temperature to obtain final product properties. In spray water-cooling, water is conveyed onto the surface of the product being cooled via spray systems. The associated heat-transfer process is rather complex, often involves boiling water, and does not readily submit to fundamentally based mathematical approaches.' ' Various 23 studies have been performed to quantify the heat flux changes during the cooling process, " [24 29] including the influence of temperature, droplet distribution, thermal radiation, and thermalphysical properties. The distribution of droplet size depends on air pressure, water-flow rate, spray velocity and nozzle characteristics. Additionally, the velocity field within the spray has spatial dependencies (in the axial and radial directions) and is usually accompanied by plume oscillations. The thermal-physical properties of the water spray are dependent on the temperature and mole fraction of the constituents. As well as the aforementioned factors, the high temperature nature of the heat treatment process may introduce additional factors that can influence the prediction of temperature. Brimacombe et al. , ' summarized the relationships between the rate of heat extraction [26 by water sprays and the spray variables for a number of different experimental studies. another study, Mudawar and Valentine' 241 In experimentally determined the local quench curves for the single phase, nucleate boiling and transition boiling regimes of spray-cooled metallic surfaces at 400°C. In the range of 400 to 800°C, some of the most reliable published heat- transfer coefficient data for water cooling is that of Mitsutsuka,' ' which is based on seven 29 Chapter 2: Literature Review 20 independent data sources and is also in excellent agreement with more recent work by Ohnishi et a/. > [30 As discussed in Nagasaka et a/., [31] water spray cooling presents additional advantages, at high pressures, over traditional water quenching processes in that the film boiling regime is prevented because the temperature at the surface of the product decreases to a value near the water temperature in a very short period of time. Thus, at high pressures, the water spray may bounce back from the impinging surfaces causing a break down of the film layer and a subsequent decrease in the heat-transfer. To model spray water-cooling, a boiling curve is usually used to represent the relationship between the heat flux or heat-transfer coefficient and surface temperature. A schematic of a typical boiling water heat-transfer curve is presented in Figure 2.4. Forced Convection Regime Nucleate Transition Boiling Boiling Regime Regime Film Soiling Regime Figure 2.4 - Boiling curve representing the relationship between the heat flux and the surface temperature during boiling heat-transfer. [32] The heat-transfer process during spray cooling can be divided into four distinct regions: 1) a natural convection region where there is no boiling at the heated surface, 2) a nucleate boiling region where isolated bubbles form at nucleation sites and separate from the surface, 3) a Chapter 2: Literature Review 21 transition boiling region where partial film boiling occurs and the heat-transfer coefficient decreases as the surface temperature increases, and 4) a film-boiling region where heat-transfer from the surface to the liquid occurs by radiation through a stable film or vapor phase and the heat flux increases with increasing change in temperature. The critical heat flux (CHF) point is important in most industrial applications because it is the point of maximum heat flux. The minimum heat flux at elevated temperature is the Liedenfrost point. 2.3.2 Radiation Radiation cooling also plays a significant role in the processing of metals and, as a boundary condition, is highly non-linear and very difficult to accurately model. This is due to the fact that emissivity is dependent on temperature, surface finish, and view factors. [5] In addition, the temperature term is to the fourth power in the flux equation. Finally, very complex geometrical effects can influence the radiative heat-transfer. The exact character of radiant heat emission is not completely agreed upon, but it is known that the Stefan-Boltzmann law gives the rate at which energy is emitted per unit area of emitting surface:' 331 q = soT (2.25) 4 where q is the rate of energy emission, T the absolute temperature of the body, cr Stefan Boltzmann's constant, and sthe emissivity. In radiation, the heat exchange is proportional to the difference of the fourth power of the absolute temperature of the radiating bodies. Thus, for a given temperature difference the heat transferred is much greater at high temperatures than at low temperatures. The phenomenon of radiation is complicated by the fact that not all the radiation leaving one object is incident on another. Thus, a technique needs to be developed to handle exchanging radiant energy between multiple bodies. This direct flux of energy between one surface and another may be expressed in terms of the shape factor, which is a geometrical quantity that is found by determining the exchange between differential area elements in each surface and then performing simultaneous integrations over both surfaces (Figure 2.5). Chapter 2: Literature Review 22 dA/cos dt Figure 2.5 - Shape factor associated with radiation exchange between elemental surface of area </A, and dAj. In many instances in numerical modeling, ' 1 J it is customary to account for radiation occurring in the cooling process by linearizing the radiation equation so that the heat-transfer rate is proportional to a temperature difference where the effective heat-transfer coefficient, h , is r given by: (2.26) where T is the surface temperature and T s amo is the temperature of the environment. This in turn allows a Cauchy type boundary condition to be applied and an overall heat-transfer coefficient, hov, to be expressed that accounts for both radiation and convection at the boundaries: (2.27) 2.4 Thermal-Physical Properties Accurate estimation of the thermal-physical properties is of paramount importance when applying mathematical models to heat-transfer problems. By and large, all of these material parameters are both temperature and phase dependent.' 331 In solving the general heat conduction equation, the properties of importance are thermal conductivity, density, specific heat, and enthalpy of transformation. Figure 2.6 shows some of the temperature dependant behavior of selected thermal-physical properties for three different steel grades. Chapter 2: Literature Review 8200, 23 1 • , , r 8000 V Temperature (°C) Figure 2.6 - Temperature dependent thermal-physical properties. [31] Although, density is temperature and phase dependent, a constant averaged density value is often employed by researchers ' ' 17 8 31,341 for mathematical modeling, owing to the conservation of mass principle during thermal changes for problems involving static meshes. In the application of thermal-physical properties to mathematical models where there is a mixture of phases present, a large number of researchers' 8,351 apply a simple rule of mixtures for any property P: (2.28) Chapter 2: Literature Review 24 where the volume fraction and property of phase / are denoted by X and P, respectively. t For problems involving convection where fluid motion is involved, the additional property of fluid viscosity becomes extremely important. The temperature dependence of fluid properties means that there may be a significant variation in property quantities. Consequently, the accuracy of the results depends on the temperature chosen for the evaluation of the fluid properties. In most cases a mean film temperature is used. The mean film temperature is an arithmetic mean of the surface temperature, T , and the external fluid/environment temperature s S a (2.29) 2 2.5 Phase Transformations Phase transformations can be of significant importance to heat-transfer when considering cooling during thermo-mechanical processing. Most commercial steel products are produced using heat treatments where austenite is cooled continuously from a high temperature to a low temperature. As it passes through the transformation temperature range, the austenite decomposes into a mixture of different transformation products, each of which can form by different mechanisms. The interactions of these transformations may involve either 'hard- impingement' (adjacent particles come into contact) or 'soft-impingement' (overlapping of diffusion or thermal fields) and are significant in determining the final microstructure.' 371 Depending on the steel chemistry and thermal history, the final microstructure may be comprised of ferrite, pearlite, bainite, martensite or some mixture. These transformation products can be grouped into diffusion-controlled equilibrium transformation products (i.e. ferrite and pearlite) or martensitic type non-equilibrium products (i.e., bainite and martensite). The bainite transformation is not well understood and is often treated as either a diffusion-controlled or martensitic transformation. Hence, the treatment of bainite will be examined in a forthcoming section. The rate of heat generation/consumption during phase transformations is proportional to the transformation rate: m phases dt (2.30) Chapter 2: Literature Review 25 In Equation (2.30), p represents the density, H the enthalpy of transformation and X the volume t t fraction of phase /. The individual contributions, /, from the different transformations (N h e in p as total) are summed to arrive at the rate of volumetric heat generation or consumption at a given time. A large portion of the literature has been concerned with the quantification of phase transformation kinetics in steels, based on thermal history and chemistry. In the past, a significant number of these studies have used isothermal measurements, due to the relative simplicity of conducting and interpreting the results. 138-421 However, the resulting models tend to have limited flexibility because they are typically only valid for the specific material of interest and do not account for the evolution of non-equilibrium products. In addition, most industrial processes are usually non-isothermal and the prediction of microstructural evolution by nonisothermal cooling of steels is complicated by variations in both the driving force for the transformation and the diffusivity of the rate controlling species/ 431 The evolution of volume fraction during solid-state phase transformations, under isothermal conditions, is usually described using the classical Johnson-Mehl-Avrami- Kolmogorov (JMAK) theory/ " 38 121 (2.31) X(t) = l-exp(-bt ) n where X is a fraction transformed, t is the transformation time, and b and rt are empirically determined constants (in literature, b is often referred to as a kinetic parameter that represents the combination of nucleation and growth rates, whereas n is related to the geometry of the growing phase and the conditions of nucleation). Umemoto et a / / 4 4 1 proposed the following modification to the J M A K equation to incorporate the effect of austenite grain size on transformation kinetics: X(t) = l-exp (2.32) where d represents the prior austenite grain size and the value of m indicates the operational r nucleation mode, where m = 1 for nucleation at grain faces, m = 2 for nucleation at edges, and m - 3 for nucleation at corners. Due to limitations in applying isothermal kinetics to describe non-isothermal conditions, several authors have proposed an additivity principle to apply the J M A K equation over a range Chapter 2: Literature Review 26 of temperatures under continuous cooling conditions/ ' "* ' As proposed by Scheil 37 45 8 [45] and Cahn/ ' the additivity rule states that if r(X, T) is the time required to reach the transformed 46 fraction X, at a constant temperature, T, starting from X = 0, and X is attained in an isothermal process after time, /, then such a transformation is an additive reaction and can be written in general form as: Cahn suggested a less restrictive condition for additivity based on early site saturation. The reaction is additive whenever the rate is a function only of the amount of transformation and the temperature. (2.34) *%- = f(Xj) at Although, the Cahn site-saturation condition is sufficient for the additivity condition to be applied, it has experimentally been found to be insufficient to predict austenite decomposition/ ' 47 To address this insufficiency, Kuban et a// ' proposed a criterion, "effective-site-saturation", 47 where the later growth of the early nuclei dominate the transformation kinetics. Effective site saturation assumes that the nuclei formed in the initial stages of the transformation dominate the subsequent transformation kinetics. 2.5.1 Ferrite/Pearlite Often, JMAK theories have been proposed to describe the diffusional phase transformation kinetics for the decomposition of austenite into the equilibrium products, ferrite and pearlite/ ' ' " ' 38 4347 49 51 By employing the principle of additivity, which considers the equilibrium portion of transformed fraction at each temperature, this approach has been successful in predicting the transformation kinetics during continuous cooling/ ' 44 501 Numerous authors have developed equations, based on steel composition and cooling rate, for the empirically determined constants, b and n, as well as the start temperature for the transformations, for both ferrite and pearlite/ ' ' ' " 43 44 47 50 551 Unfortunately, it is difficult to make comparisons between the many studies that have been conducted on the austenite to pearlite Chapter 2: Literature Review 27 transformation, because differences in chemical composition and austenite grain size result in variations in the kinetics. 2.5.2 Bainite The transformation of austenite to bainite occurs when steel is cooled rapidly to below the bainite start temperature, B . s Bainite, in terms of its microstructure, can be defined as the product of a non-lamellar eutectoid reaction. In eutectoid steels it is a mixture of ferrite and cementiteJ 561 The austenite to bainite reaction has a dual nature; it is comprised of nucleation and growth transformation characteristics similar to the austenite to pearlite transformation, yet in other ways has characteristics of the austenite to martensite transformation. Numerous transformations. theories have been Rees and Bhadeshia [57] developed to explain the kinetics of bainite proposed a diffusionless displacive mechanism for bainite formation where bainite evolves through random nucleation of ferrite subunits to form bainite sheets. In this nucleation-controlled model, the effective nucleation depends on the current size and volume fraction of the bainite phase. In contrast, Lusk and L e e [581 developed a diffusional transformation model using non-additive kinetics to simulate the kinetics: (2.35) where Xf eq is the temperature dependent equilibrium volume fraction of ferrite and the expression X ,b =Yb(Xf + X ) reflects the influence of the existing ferrite and pearlite on the rate P En of formation of bainite. p, is the influence parameter and kb the temperature dependent term: k = - >2 ea b 029g e (A -T) Q/R{T+273) (2.36) el where g is the A S T M grain size number of austenite, Q the effective thermal activation energy, and nb the critical exponent that relates the driving force to the degree of undercooling. Umemoto et a / . [59] and Suehiro et a / . [531 bainite describing the rate constant, b. also developed diffusional transformational models for Minote et a / . [60] stated that a combination approach is possible, using a diffusion mechanism to make predictions where the process occurs at temperatures above 350°C and a displacive mechanism provides the most accurate predictions at temperatures below 350°C. The Avrami model has also been applied to predict bainite transformation kinetics despite the controversy regarding the different growth mechanisms. This Chapter 2: Literature Review 28 method is limited, as different rate constants have been reported and it is not possible to predict these constants. 2.5.3 Martensite Martensite is formed when austenite rapidly cools to a relatively low temperature. Martensite is a non-equilibrium single-phase structure that results from a diffusionless transformation of austenite.' ' It is athermal (i.e., independent of time), and is only a function of 61 the temperature to which the alloy is quenched or rapidly cooled. Very little literature exists to describe the kinetics of transformation for martensite. In general, the kinetics of transformation have been characterized by the empirical equation proposed by Koistinen-Marburger:' ' 62 X = l- xp[-A (M -T)] Q m (2.37) s where A is the Koistinen-Marburger coefficient and M is the martensitic transformation start m s temperature. Various authors A. m 18,59,63 ' have employed this equation assigning a value of 0.011 to The martensite start temperature, M , is dependent on the concentration of the austenites stabilizing alloying elements. Andrews has described the M with the quantitative formula: ' 164 s M , ( ° C ) = 512-453C-16.9M +150 -9.5M> + 217C - l \ . 5 C M n - 6 7 . 6 C O 2 (2.38) Although, the trend among researchers is to use the Andrew's equation, others' ' have 65 used the formula of Steven and Haynes:' ' 66 M (°C) = 561 - 474C - 17M - 21M> - 33Mw - 17Cr s Another less common approach, used by Cheng et al, (2.39) assumes the fraction of martensite between the M and the martensite finish temperature, M/, is a function of temperature s where: (2.40) resulting in: max m 1- T-M, M -M s Y fj (2.41) Chapter 2: Literature Review 29 where p is a constant (value of 2.5) and X™ represents the maximum amount of martensite (i.e., K already transformed martensite and austenite). 2.6 Process Modeling The mathematical modeling of cooling in thermal processing is generally comprised of two components: heat-transfer and microstructure evolution. Owing to the fact that the final evolved microstructure in rod or tubular product is dependent on its thermal history and viceversa, a large number of researchers couple these two components/ ' ' ' ' ' ' ' " 8 22 23 31 34 35 43 70 The 721 linking of the thermal response and microstructural evolution is typically performed through transformation kinetics and associated latent heats of transformation. In addition, these mathematical models require the availability of thermal-physical property and transformation kinetic data. Experimental tests performed in the laboratory and/or industrial plants are generally used for the validation and verification of the different mathematical models. A detailed discussion of a selection of published mathematical models, developed to predict thermal distribution and microstructural evolution during thermo-mechanical processing, is presented below. Early models, like those proposed by Inoue and Tanaka, [68] assumed constant thermal- physical properties and neglected the heat of transformation. More recent work by Rammerstorfer et al. \ [69 attempted to account for the latent heat evolved from phase transformations by adopting, in the temperature range where phase transformations take place, a modified value of specific heat. This approach simplified the calculations for computing the rate of heat evolved because it eliminated the need for calculating the microstructural evolution. Agarwal and Brimacombe 1701 modeled the heat flow and the microstructure evolution (austenite to pearlite transformation) during water quenching of eutectoid steel rods. In the model, finite-difference methods were adopted to numerically solve the heat balance equation. The rate of volume fraction of pearlite transformed was determined by invoking the additivity principle and adopting a two-parameter empirical equation. Isothermal transformation diagrams were used to estimate the empirical parameters in the kinetic equation. In an improvement to this model, Iyer et o/. [71] separated the incubation period from the transformation event. In this study, the experimentally determined kinetic parameters were utilized, rather than ones computed from published T T T diagrams. Campbell et a/./ ' ' 22 34 431 in an in-depth investigation, extended the phase transformation kinetics (including austenite to ferrite transformation) Chapter 2: Literature Review 30 database to hypo-eutectoid and eutectoid steels. A mathematical model incorporating heat flow, microstructural phenomena, and mechanical property relationships for steel rods experiencing controlled cooling for the Stelmor line was developed. These results compared reasonably well to measured laboratory tests and plant trials. Kumar et alS ^ investigated three different methods for modeling the thermal fields 7 associated with the austenite-to-pearlite phase transformation during heating and cooling o f hypo-eutectoid and eutectoid steels. Their methods included: 1) computing the rate o f the fraction transformed from IT diagrams, 2) computing the rate o f fraction transformed from an iron-carbon equilibrium diagram, and 3) adopting a modified heat capacity in which the latent heat value was assigned a value o f 0. A control-volume based finite-difference approach was adopted for the numerical solution o f the heat balance equation and the boundary conditions associated with cooling consisting o f air-cooling (combined radiative and convective heat flux) and water quenching (constant heat-transfer coefficient). A constant density and specific heat as well as a temperature dependent thermal conductivity were adopted for the calculations. The predicted results were validated against analytical solutions as well as experimental thermal response data. Under the experimental conditions examined, it was determined that there was very little difference between the methods owing to the relatively small contribution o f the heat source term. The second and third method, however, were determined to be easier to implement algorithmically. Nagasaka et al. [3l] developed a mathematical model, based on finite-element techniques, that incorporates phase transformations and the elastic-plastic stress behavior during water spraying o f steel bars. This involved the coupling o f diffusion-dependent and martensitic phase transformations with transient heat flow. The code, initially verified with published literature for pure iron, was subsequently validated against experimental test on bars o f 1035 carbon and nickel-chromium alloyed steels. The predicted results compared reasonably w e l l with the measured data o f temperature, microstructure, hardness and residual stresses. Thomas et al [23] developed three mathematical models to determine the temperature and thermal-stress distribution in steel cylinders subject to spray-water cooling. These models were based on finite-element methods and incorporated temperature-dependent material properties, latent heat generation, elastic-plastic behavior, and volumetric expansion associated with the formation o f martensite. The simulated results were then verified experimentally, using an AISI Chapter 2: Literature Review 31 4140 steel grade at an initial temperature of 1273 K, and shown to predict temperature and microstructural behavior reasonably well. Anelli [72] developed a comprehensive set of integrated mathematical models for simulating hot rolling and controlled cooling of wire rods and bars. Experimental tests were carried out on C - M n and eutectoid steels. These models were developed primarily for investigating relationships between process variables and the thermal and microstructure evolution in rods and bars, with the end goal being able to quantify the effects of rod diameter and cooling practices on microstructure and temperature along the Stelmor line. In a later investigation, Anelli et a/. ' extended this model to predict the through-thickness temperature [8 evolution and formation of microstructure constituents in tubular products for a range of steel grades. This model was again validated using laboratory and industrial trials. Hodgson et a/., [35] modeled the post deformation cooling behavior of hot worked steel for a range of spray cooling systems. The heat balance equation was solved using finite-difference techniques and the latent heat associated with phase transformations was incorporated into the thermal analysis as a function of the rate of transformation during cooling. The model was restricted to the transformation of ferrite, pearlite, and martensite. The rate of transformation for the diffusional products was implemented through an Avrami type equation, assuming additivity, while the fraction martensite formed was modeled as a function of the martensite start temperature. In this work a generalized equation was developed that relates the heat-transfer coefficient during spray cooling to the water flux and the surface temperature. The model was validated using both laboratory and mill data for cooling systems ranging from soft cooling to intense cooling. The application of the model demonstrated that the heat-transfer coefficient does require scaling when it is applied to the cooling of large components or when air-mist rather than hydraulic cooling is used. 2.7 Summary A review of the literature concerning the thermal-microstructural modeling of steel rod and tubular products has been performed and a summary of the key concepts is presented below. A variety of different numerical methods can be used to solve the partial differential equations that describe the evolution of the temperature and microstructure within thermomechanical processing. The choice of which particular technique to use ultimately lies with the Chapter 2: Literature Review 32 user, however, simpler techniques are normally used to solve simple problems, with the sophistication of the technique increasing with problem complexity. The most popular of these methods for computer programming are the finite-difference and finite-element methods. Within the finite-difference methods the PDEs may be solved using explicit and/or implicit techniques. For problems with simple geometry, ADIfinite-differencemethods incorporating latent heats of transformation are more cost-effective than thefinite-elementmethods. The development of thermal models requires boundary conditions that accurately describe the process. During post hot-deformation cooling, convective (natural, forced-air and spray-water cooling) and radiative boundary conditions are of particular significance. The complexity of convective heat-transfer generally necessitates the use of empirical correlations. For natural and forced-air convection these correlations relate the Nusselt's number (Nu), to the Reynolds number(/?e), Rayleigh's number (Ra), and Prandlt's number (Pr). In general, these correlations are specific to the particular convective phase, however, relations covering a large range of Re or Ra are often better suited for modeling a variety of cooling conditions. Spray-water cooling is generally describe through the use of spray-water boiling curves that represent the relationship between the heat flux or heat-transfer coefficient and surface temperature. Within these curves the heat-transfer process may be divided into four distinct regions: natural convection, nucleate boiling, transition boiling, and film boiling. In most industrial applications the critical heat flux is the most important point on the boiling curve, as it represents the maximum heat flux. The non-linearity and complexity of radiative heat-transfer makes it difficult to model. Shape functions (geometrical quantities), are often used to account for adjacent surfaces owing to the fact that not all the radiation leaving one object is incident to another. One approach to simplify problems involving radiation is to linearize the radiation equation so that the heattransfer rate is proportional to a temperature difference and may be written in terms of an effective heat-transfer coefficient. This effective heat-transfer coefficient allows a Cauchy-type boundary condition to be applied, and an overall heat-transfer coefficient to be expressed that accounts for radiation and convection. When applying mathematical models to heat-transfer problems it is important to have an accurate estimation of the thermal-physical properties. Temperature and phase dependencies of Chapter 2: Literature Review 33 the material parameters may add complexity to the heat flow problem. In situations where there is a mixture of phases present, a simple rule of mixtures is used. Transformation products can be grouped into diffusion-controlled equilibrium transformation products (i.e. ferrite and pearlite) or martensitic type non-equilibrium products (i.e., bainite and martensite). The transformation kinetics for ferrite, pearlite, bainite, and martensite reactions affect the thermal model through changes to the thermal-physical properties and the heat of transformation (proportional to the transformation rate). The two mutually coupled components associated with post hot-deformation cooling are temperature field calculations and phase transformations. In order to effectively model these, several key elements are required: 1. Numerical techniques for solving the PDEs, that accurately and efficiently solve the heatflow problem. 2. Boundary conditions that accurately describe the heat-transfer occurring at the surface of the product for the particular process. 3. Material property data that is both temperature and phase dependent. 4. Microstructure evolution routines describing austenite decomposition as it relates to diffusion-controlled and martensitic type phase transformations. 5. Appropriate coupling and integration of the individual components. 6. An important aspect when developing computer simulation models is to compare the results with experimental data. The literature has revealed that it is both necessary and possible to integrate thermal and microstructural models incorporating temperature and phase dependent material properties for cooling of steel rod and tubular products. There are, though, still a number of areas that require further development before these methods can be more widely applied. These areas include improvements in transformation kinetic routines and improvements in the characterization of boundary conditions (i.e., radiation and boiling water heat-transfer). Chapter 2: Literature Review 34 References 1. J. Rohde and A. Jeppsson, "Literature Review of Heat Treatment Simulations with Respect to Phase Transformation, Residual Stresses and Distortion", Scandinavian Journal of Metallurgy, Vol. 29, 2000, pp 47-62. 2. R.W. Lewis, K. Morgan, and O.C. Zienkiewicz, Numerical Methods in Heat-transfer, John Wiley & Sons, New York, 1981, pp 1. 3. F. Kreith and W.Z. Black, Basic Heat-transfer. Harper and Row, New York, 1980. 4. T. Shih, Numerical Heat-transfer. Hemisphere, New York, 1984. 5. A.J. Chapman, Heat-transfer. 4 ed., Macmillan, New York, 1984. th 6. W.J. Minkowycz et ai, Handbook of Numerical Heat-transfer, John Wiley and Sons, New York, 1980. 7. A . Kumar Singh and D. Mazumdar, "Comparison of Several Numerical Methods for Thermal Fields During Phase Transformation of Plain Carbon Steels", ISIJ International, Vol. 31, No. 12, 1991, pp 1441-1444. 8. E . Anelli, P. Farsetti, and G . Cumino, "Prediction of the Response to Quenching and Tempering of Low Alloy Steels for Seamless Tubes", Proceedings of 36 th Mechanical Working and Steel Processing Conference, Baltimore, M D , Vol. 32, 1994, ppl49-159. 9. P.J. Davis and P. Rabinowitz, Methods of Numerical Integratioa Academic Press, New York, 1975. 10. R.P. Canale and S. Chapra, Numerical Methods for Engineers: With Programming and Software Applications, 3 ed., McGraw-Hill, New York, 1998. rd 11. B. Carnahan, H . Luther, and J.Wilkes, Applied Numerical Methods, John Wiley and Sons, New York, 1969, pp 452. 12. D. Peaceman and H. Rachford, Journal of Society of Industrial Applied Mathematics, Vol. 3, 1955, pp 28-41. Chapter 2: Literature Review 35 13. B.G. Thomas, I.V. Samarasekera, and J.K. Brimacombe, "Comparison of Numerical Modeling Techniques For Complex, Two-Dimensional, Transient Heat-Conduction Problems", Metallurgical Transactions . B, Vol. 15B, 1984, pp 307-318. 14. K. Huebner, The Finite Element Method for Engineers. John Wiley and Sons, New York, 1975. 15. T. Hughes, The Finite Element Method, Dover Publications, New York, 2000. 16. Zhukauskas, "Heat-transfer from Tubes in Crossflow", Advances in Heat-transfer, Vol. 8, 1972, pp 93. 17. S. W. Churchill, and M . Bernstein, "A Correlating Equation for Forced Convection from Gases and Liquids to a Circular Cylinder in Crossflow", Journal of Heat-transfer, Vol. 94, 1977, pp 300. 18. S.W. Churchill and H.S. Chu, "Correlating Equations for Laminar and Turbulent Free Convection from a Vertical Plate", International Journal of Heat Mass Transfer, Vol. 18, 1975, pp 1323. 19. V . T Morgan, "Overall Convective Heat-transfer for Smooth Circular Cylinders", in T.F. Irvine and J.P. Hartnett, Eds., Advances in Heat-transfer. Vol. 11, Academic Press, New York, 1975, pp 199-264. 20. W.H. McAdams, Heat Transmission. 3 ed., New York, McGraw-Hill, 1954. rd 21. M . Golja, "Mathematical Model of Steel Hollow Shell Cooling after Rolling", Metalurgija, Vol. 36, No. 4, 1997, pp 229-234. 22. P.C. Campbell, E.B. Hawbolt, and J. K. Brimacombe, "Microstructural Engineering Applied to the Controlled Cooling of Steel Wire Rod: Part I. Experimental Design and Heat-transfer", Metallurgical Transactions A , Vol. 22A, No. 11, 1991, pp 2769-2778. 23. R. Thomas, M . Ganesa-Pillai, P.B. Aswath, K . L . Lawrence, and A . Haji-Sheikh, "Analytical/finite-element Modeling and Experimental Verification of Spray-cooling Process in Steel", Metallurgical and Materials Transactions A, Vol. 29A, 1998, pp 1485-1498. Chapter 2: Literature Review 36 24.1. Mudawar and W.S. Valentine, "Determination of the Local Quench Curve for SprayCooled Metallic Surfaces", Journal of Heat Treating, Vol. 7, 1989, pp 116-117. 25. D. L i , "Boiling Water Heat-transfer Study during D C Casting of Aluminum Alloys", M . A . Sc. Thesis, University of British Columbia, 1999. 26. J.K. Brimacombe, P.K. Agarwal, S. Hibbins, B. Prabhaker, and L . A . Baptista, "Spray Cooling in the Continuous Casting of Steel", Iron and Steel Engineer, pp 109-132. 27. M . Bamberger and B. Prinz, "Detenriination of Heat-transfer Coefficients during Water Cooling of Metals", Materials Science and Technology, Vol. 2, 1986, pp 410-415. 28. S. Segerberg, "Spray Quenching a Method with Great Flexibility", in Heat Treatment and Surface Engineering: New Technology and Practical Applications, Chicago, IL, 1988, pp 2830. 29. M.Mitsutsuka, Tetsu-to-Hagane, Vol. 69, 1983, pp 268. 30. A.Ohnishi, H. Takashima and M . Hariki, Trans. ISIJ, Vol. 27, 1987, pp B299. 31. Y . Nagasaka, J.K. Brimacombe, E.B. Hawbolt, I.V. Samarasekera, B. Hernandez-Morales, and S. E . Chidiac, "Mathematical Model of Phase Transformations and Elasto-Plastic Stress in the Water Spray Quenching of Steel Bars", Metallurgical Transactions A , Vol. 24A, 1993, pp 795-808. 32. M . Bamberger and B. Prinz, "Methods for Controlling Heat Flow during Quenching", Metallwissenschaft & Technik, Vol. 44, No. 2, 1990, pp 154-157. 33. F.P. Incopera and D.P. DeWitt, Fundamentals of Heat and Mass Transfer, John Wiley & Sons Inc., New York, N Y , 1990. 34. P.C. Campbell, E.B. Hawbolt, and J. K. Brimacombe, "Microstructural Engineering Applied to the Controlled Cooling of Steel Wire Rod: Part III. Mathematical Model - Formulation and Predictions", Metallurgical Transactions A , Vol. 22A, No. 11, 1991, pp 2791-2805. 35. P.D. Hodgson, K . M . Browne, R.K. Gibbs, T.T. Pham, and D . C . Collinson, "The Mathematical Modeling of Temperature and Microstructure During Spray Cooling", Chapter 2: Literature Review 37 Proceedings of the First International Conference on Quenching & Control of Distortion, Chicago, IL, 1992, pp 41-49. 36. J.P. Holman, Heat-transfer. McGraw-Hill, New York, 1981. 37. J. Christian, The Theory of Transformation in metals and alloys. Pergamon, Oxford, U K , 1975, pp 492. 38. M . Avrami, J. Chem. Phys., Vol. 7, 1939, pp 1103-1112. 39. M . Avrami, J. Chem. Phys., Vol. 8, 1940, pp 212-224. 40. M . Avrami, J. Chem. Phys., Vol. 9, 1941, pp 177-183. 41. W.A. Johnson and R.F. Mehl, "Reactions Kinetics in Process of Nucleation and Growth", Trans. AIME, Vol. 135, 1939, pp 416. 42. A . N . Kolmorgorov, "Izv.Akad. Nauk. SSSR, Ser. MatemNo. 3, 1937, pp 335. 43. P.C. Campbell, E.B. Hawbolt, and J.K. Brimacombe, "Microstructural Engineering Applied to the Controlled Cooling of Steel Wire Rod: Part II. Microstructural Evolution and Mechanical Properties Correlations", Metallurgical Transactions A , Vol. 22A, 1991, pp 2779-2790. 44. M . Umemoto, N. Komatsubara, and I. Tamura, "Prediction of Hardenability Effects from Isothermal Transformation Kinetics", Journal of Heat Treating, Vol. 1, No. 3, 1980, pp 5764. 45. E . Scheil, "Anlaufzeit der Austenitumwandlung, Archiv f. Eisenhuttenw, 8 (565), 1935. 46. J.W. Cahn, "Transformation kinetics during continuous cooling", Acta Met., Vol. 4, 1956, pp 572. 47. M.B. Kuban, R. Jayaraman, E.B. Hawbolt, and J.K. Brimacombe, "An Assessment of the Additivity Princeple in predicting Continuous-Cooling Austenite-to-Pearlite Transformation Kinetics Using Isothermal Transformation Data", Metallurgical Transactions A., Vol. 17A, 1986, pp 1493-1503. 48. T. Reti, T. Bell, Y . Sun, A. Bloyce, Material Science Forum 163-165, 1994, pp 673-680. Chapter 2: Literature Review 38 49. M . Militzer, R. Pandi, and E.B. Hawbolt, "Ferrite Nucleation and Growth during Continuous Cooling", Metallurgical and Materials Transactions A, Vol. 27A, 1996, pp 1547-1556. 50. E. B. Hawbolt, B. Chau, and J . K. Brimacombe, "Kinetics of Austenite-Ferrite and AustenitePearlite Transformations in a 1025 Carbon Steel", Metallurgical Transactions A , Vol. 16A, 1985, pp 565-578. 51. R. Vandermeer, "Modeling Diffusional Growth during Austenite Decomposition to Ferrite in Polycrystalline Fe-C Alloys", Acta. Metall. & Mater., Vol 38, N.12, 1990, pp 2461-2470. 52. E . B. Hawbolt, B. Chau, and J . K. Brimacombe, "Kinetics of Austenite-Pearlite Transformations in Eutectoid Carbon Steel", Metallurgical Transactions A, Vol. 14A, 1983, pp 1803-1815. 53. M . Suehiro, K. Sato, Y . Tsukano, H . Yada, T. Senuma, Y . Matsumura, "Computer Modeling of Microstructural Change and Strength of Low Carbon Steel in Hot Strip Rolling", Trans. I S I J , 27, 1987, pp 439-445. 54. A. Kumar, C. McCulloch, E . B. Hawbolt, and I. V. Samarasekera, "Modelling Thermal and Microstructural Evolution on Runout Table of Hot Strip Mill", Materials Science and Technology, April 1991, Vol.7, pp 360-368. 55. M . Militzer et al., "Modelling the Microstructure of Microalloyed Low-Carbon Steels",Thermec 1997, TMS, Warrendale, PA. 56. W.F. Smith, Structure and Properties of Engineering Alloys, McGraw-Hill Inc, New York, N Y , 1993. 57. G . I. Rees and H. K. D. H. Bhadeshia, "Bainite Transformation Kinetics, Part 1 Modified Model", Materials Science and Technology, Vol. 8, 1992, pp 985-993. 58. M . T. Lusk and Y . Lee, "A Global Material Model for Simulating the Quench Kinetics of Low Alloy Steels", Hungarian Scientific Society of Mechanical Engineering, Heat Treatment and Surface Engineering of Light Alloys: Proceedings of the 7th International Seminar of IFHT (Hungary), Sept. 1999, pp 273-282. Chapter 2: Literature Review 39 59. M . Umemoto, A. Hiramatsu, A. Moriya, T. Watanabe, S. Nanaba, N . Nakajima, G. Anan, Y . Higo, "Computer Modeling of Phase Transformation from Work-Hardened Austenite", ISIJ International, 32, 1992, pp 306-315. 60. T. Minote, S. Torizuka, A. Ogawa, M . Niikura, "Modeling of Transformation Behavior and Compositional Partitioning in Trip Steel", ISIJ International, 36, 1996, pp 201-207. 61. W.D. Callister, Materials Science and Engineering - An Introduction. John Wiley & Sons, New York, N Y , 2001. 62. D. P. Koistinen and R. E . Marburger, Acta Metall., Vol 7,1957, pp 59-60. 63. L. Zhou, X . Wei, H . Zeng, "A Cellular Automaton Approach to Modeling of Phase Transformation in Quenching of Steels", Proceedings of Heat Treating 18 Conference, Liu th Dai Memorial Symposium, A S M International, 1998, pp 700-704. 64. K.W. Andrews, JISI, Vol. 203, 1965, pp 721. 65. J.S. Kirkaldy, B.A. Thomson, and E . A . Baganis, in "Phase Transformations in Ferrous Alloys", A.R. Marder and J.I Goldstein eds., TMS-AIME, Warrendale. PA, 1984, pp 82-125. 66. W. Steven and A . G . Hayes, "The Temperature of Formation of Martensite and Bainite in Low Alloy Steels", JISI, Vol. 183, 1956. 67. H . Cheng, F. Jinag, and H.G. Wang, Journal of Material Processing Technology, Vol. 63, 1997, pp 568-572. 68. T. Inoue and K. Tanaka, "An Elastic-Plastic Stress Analysis of Quenching when Considering a Transformation", International Journal of Mech. Sci., Vol. 17, 1975, pp 361-367. 69. F.G. Rammerstorfer et al., "On Thermo-Elastic-Plastic Analysis of Heat Treatment Including Creep and Phase Changes", Computers and Structures, Vol.13, 1981, pp 771-779. 70. P.K. Agarwal and J.K. Brimacombe, "Mathematical Model of Heat Flow and AustenitePearlite Transformation in Eutectiod Carbon Steel Rods Transactions B, Vol. 12B, 1981, pp 121-133. for Wire", Metallurgical Chapter 2: Literature Review 40 71. J.K. Iyer, J.K. Brimacombe, and E.B. Hawbolt, "Prediction of the Structure and Mechanical Properties of Control-cooled Eutectoid Steel Rods", Mechanical Working and Steel Processing XXII (ISS), 1985, pp 47-58. 72. E . Anelli, "Application of Mathematical Modeling to Hot Rolling and Controlled Cooling of Wire Rods and Bars", ISIJ International (Japan), Vol. 32, No.3, 1992, pp 440-449. Chapter 3: Scope and Objectives 41 CHAPTER 3 SCOPE AND OBJECTIVES 3.1 Objectives Although much research has been done to model temperature and microstructure evolution in the thermo-mechanical processing of steels such as sheet, plate, and rod material/ " 1 1 5 1 limited work has been done with respect to modeling seamless steel tubing. Specifically, the process modeling of cooling conditions on tubes, relating heat treatment variables and steel parameters for process control purposes, has had little attention. The objectives for this research are: • To develop 2-D transient heat diffusion models to predict the variation in temperature under various cooling conditions for seamless tube during the following manufacturing processes: 1. Tubes subject to water spray quenching (i.e., as might occur between deformation stands to control austenite grain growth). 2. Tubes subject to controlled slow cooling, where radiation between adjacent tubes would be significant (i.e., as might occur at the end of the process to produce the final desired microstructure and mechanical properties). 3.2 Methodology The methodology applied to this research involved both the development of 2-D mathematical transient heat diffusion models and experimental tests at the Timken Bearing & Steel Company, Canton, Ohio. An evolutionary approach was adopted for model development as well as model validation. It started by building a simple model and then verifying it, with subsequent versions being built on the previously validated versions. The validation of the thermal diffusion models, incorporating latent heat additions, were first performed using a commercial finite-element package A B A Q U S and second through temperature measurements obtained from experimental tests. The industrial steel tube samples instrumented at U B C , and an experimental set-up designed and built at Timken, were used to measure the thermal history of 42 Chapter 3: Scope and Objectives each sample during cooling operations. Figure 3.1 depicts the overall methodology for this project. As illustrated in the diagram, the project has two major components: 1. The development of two 2-D transient heat diffusion models to predict the temperature distribution during cooling conditions of steel tubes. 2. Validation through a combination of comparison to benchmark codes (ABAQUS) and experimentally derived measurements (experimental tests of controlled slow cooling conditions). Although, the modeling activities involved the development of two models to describe the cooling conditions of rapid quenching and controlled slow cooling, the experimental portion of this study focused on the controlled slow cooling. Work is currently underway to validate the quench model and is beyond the scope of this thesis. The Timken Company provided the microstructural element (i.e., incorporating recrystallization kinetics, grain growth kinetics, and austenite decomposition kinetics for the formation of ferrite, pearlite, bainite and martensite) of the coupled 2-D transient heat diffusion models. Consequently, the proprietary nature of this research places the validation of the phase transformation routines, as well as the microstructural experimental analysis, outside the scope of this current study. The coupling, however, of the transformation kinetics routines and the associated evolved latent heat, as well as the incorporation of temperature and chemistry dependent performed by U B C . thermal-physical properties, was Chapter 3: Scope and Objectives 43 Chapter 3: Scope and Objectives 44 References 1. J. Rohde and A. Jeppsson, "Literature Review of Heat Treatment Simulations with Respect to Phase Transformation, Residual Stresses and Distortion", Scandinavian Journal of Metallurgy, Vol. 29, 2000, pp 47-62. 2. R. Thomas, M . Ganesa-Pillai, P.B. Aswath, K . L . Lawrence, and A. Haji-Sheikh, "Analytical/finite-element Modeling and Experimental Verification of Spray-cooling Process in Steel", Metallurgical and Materials Transactions A , Vol. 29A, 1998, pp 1485-1498. 3. I. V. Samarasekera, D. Q. Jin and J. K. Brimacombe, "The Application of Microstructure Engineering to the Hot Rolling of Steel", 38th Mechanical Working and Steel Processing Conference Proceedings, Vol. 34, 1996, pp 313-327. 4. C. A . Muojekwu, D. Q. Jin, V . H. Hernandez, I. V . Samarasekera and J. K. Brimacombe, "Hot-direct Rolling, Runout Table Cooling and Mechanical Properties of Steel Strips Produced from Thin Slabs", Proceedings of the 38th Mechanical Working and Steel Processing Conference, Cleveland, Ohio, Vol. 34, 1996, pp 351-366. 5. C. Dilg, J. Kirsch, W. Schutz, E. Amoris, and A. Streibelberger, "TMCP Process Including Two Accelerated Cooling Stages", La Revue De Metallurgie-CIT, 1995, pp 884-892. 6. Y . Nagasaka, J.K. Brimacombe, E.B. Hawbolt, I.V. Samarabekera, B. Hernandez-Morales, and S.E. Chidiac, "Mathematical Model of Phase Transformations and Elasto-Plastic Stress in the Water Spray Quenching of Steel Bars", Metallurgical Transactions A , Vol. 24A, 1993, pp 795-808. 7. P.D. Hodgson, K . M . Browne, R.K. Gibbs, T.T. Pham, and D . C . Collinson, "The Mathematical Modeling of Temperature and Microstructure During Spray Cooling", Proceedings of the First International Conference on Quenching & Control of Distortion, Chicago, IL, 1992, pp 41-49. 8. E . Anelli, "Application of Mathematical Modeling to Hot Rolling and Controlled Cooling of Wire Rods and Bars", ISIJ International (Japan), Vol. 32, No.3, 1992, pp 440-449. 9. A. Kumar, C. McCulloch, E . B. Hawbolt, and I. V . Samarasekera, "Modelling Thermal and Microstructural Evolution on Runout Table of Hot Strip Mill", Materials Science and Technology, April 1991, V o l . 7, pp 360-368. Chapter 3: Scope and Objectives 45 10. A. Kumar Singh and D. Mazumdar, "Comparison of Several Numerical Methods for Thermal Fields During Phase Transformation of Plain Carbon Steels", ISIJ International, Vol. 31, No. 12, 1991, pp 1441-1444. 11. P.C. Campbell, E.B. Hawbolt, and J. K. Brimacombe, "Microstructural Engineering Applied to the Controlled Cooling of Steel Wire Rod: Part I. Experimental Design and Heat Transfer", Metallurgical Transactions A, Vol. 22A, No. 11, 1991, pp 2769-2778. 12. P.C. Campbell, E.B. Hawbolt, and J.K. Brimacombe, "Microstructural Engineering Applied to the Controlled Cooling of Steel Wire Rod: Part II. Microstructural Evolution and Mechanical Properties Correlations", Metallurgical Transactions A , Vol. 22A, No. 11, 1991, pp 2779-2790. 13. P.C. Campbell, E.B. Hawbolt, and J.K. Brimacombe, "Microstructural Engineering Applied to the Controlled Cooling of Steel Wire Rod: Part III. Mathematical Model - Formulation and Predictions", Metallurgical Transactions A , Vol. 22A, No. 11, 1991, pp 2791-2805. 14. J.K. Iyer, J.K. Brimacombe, and E.B. Hawbolt, "Prediction of the Structure and Mechanical Properties of Control-cooled Eutectoid Steel Rods", Mechanical Working and Steel Processing XXII (ISS), 1985, pp 47-58. 15. P.K. Agarwal and J.K. Brimacombe, "Mathematical Model of Heat Flow and AustenitePearlite Transformation in Eutectoid Carbon Steel Rods Transactions B, Vol. 12B, 1981, pp 121-133. for Wire", Metallurgical Chapter 4: Mathematical Modeling 46 CHAPTER 4 MATHEMATICAL MODELING During the controlled thermo-mechanical processing of seamless tubes, a variety of different cooling systems may be used to cool the tubes. These may occur both 'in-line' (in between deformation processes) and/or at the end of the process (post-deformation). In-line cooling often involves a series of water-quench or forced-air stands (Figure 4.1a) acting on individual tubes moving through the process in the axial direction. Post-deformation cooling usually involves groups of tubes (stacked side-by-side and/or on top of each other) rolling on outer-diameters, moving perpendicular to the axial direction along cooling-beds (frequently incorporating hoods) - Figure 4.1b. (a) (b) Figure 4.1 - Schematics of different cooling systems: (a) water-quench stands and (b) slow cooling bed. The heat flow component, of the overall process mathematical model, is necessary to provide accurate predictions of the thermal response of the tubes as they are processed. Because latent heat can be evolved from the solid-state transformations that may occur in the process, the heat flow analysis has to be coupled to a microstructural model capable of quantitatively predicting the evolution in microstructure. In the present study, two heat-transfer models have been developed that are capable of predicting the thermal history within the seamless tubes during two of the processing stages: 1) during rapid water quenching, such as might occur between deformation stands to control austenite grain growth; and 2) during slow controlled cooling, such as might occur at the end of the process to produce the final desired microstructure. These models are referred to as the accelerated cooling model (ACM) and the controlled slow cool model (CSC), respectively. These models have been coupled with a microstructure evolution model (incorporating recrystallization kinetics, grain growth kinetics, and austenite decomposition kinetics) provided Chapter 4: Mathematical 47 Modeling by the Timken Company. The coupling occurring between the temperature and microstructure models is illustrated in Figure 4.2. Temperature and cooling rate dependent phase transformations Latent heat Figure 4.2 - Schematic depicting the relationship between temperature and microstructure in cooling process. 4.1 General PDE and Solution Approach The general governing equation for the transient heat conduction for both models can been expressed as : [1] V-kVT dT + Q = pc„ — (4.1) In Equation (4.1), the differential operator, V, is modified depending on whether the A C M or C S C model is being considered (Sections 4.1.1 and 4.1.2 respectively). The governing P.D.E. (Equation 4.1) is solved subject to initial and boundary conditions. The material properties and rate of latent heat generation vary significantly with temperature under typical process conditions. Hence, analytical solutions cannot be used and a numerical scheme must be adopted. A finite-difference approach has been used to solve the heat-transfer problems. The development of the A C M and CSC mathematical models is presented in the following sections. Chapter 4: Mathematical Modeling 48 4.1.1 Accelerated Cooling Model The A C M model was developed with the premise that the tubes enter the cooling stands separately, with their axial lengths parallel to the direction of motion. The decision to use axial symmetry in the A C M analysis originates from the assumption that the tubes experience uniform boundary conditions in the circumferential direction. Although axial symmetry has been postulated for the A C M model, inherent limits exits within this assumption. In particular, these limits stem from implant processing conditions (i.e. rotational speeds of the tubes proceeding down the processing line, the level conduction between tubes and processing-line supports, quench stand configurations, etc..) and will affect the level of axial symmetry within the tube during heat conduction and overall heat-transfer. The model allows for moving boundary conditions, multiple cooling zones, and on/off switches that can be applied to the outer diameter (OD) and/or the inner diameter (ID) of the tube. It incorporates a dynamic time-step approach to improve the speed and accuracy of results, and allows for both user defined thermal physical property data and heat-transfer coefficient data (e.g., boiling water curves). The model is specifically designed to predict thermal history and microstructure evolution. Geometry A rectangular section has been used to represent the axisymmetrical geometry of the seamless tube within the A C M , in keeping with the aforementioned assumptions (Figure 4.3). Figure 43 - Schematic of axisymmetrical geometry - A C M model. Chapter 4: Mathematical 49 Modeling This rectangular section consists of four sides that represent the outer diameter, inner diameter, and two tube ends. To increase versatility, the model permits user definition of the inner and outer diameter and the tube length. Heat Flow For the A C M model, the gradient and Laplacian operator within the transient heat conduction, Equation (4.1) are expressed in cylindrical coordinates as follows: V.Wr-I|.(/b-fl) |.(*fl) (4.2) + r or or dz dz The applicable boundary conditions for the solution of Equations (4.1) and (4.2) are defined as follows: 1) At the outside surface of the tube (r = R ), at times greater than 0 (P*0), a -k = KXT -T ) dT\ dr s =q amb @r = R , OD a (4.3) f>0 2) At the inner surface of the tube (r = Ri), at times greater than 0 (/>0), -k -L\ d = KXT -T ) s @r = R, t>0 (4.4) E0 @z = 0 t>0 (4.5) =q @z = z t>0 (4.6) =q amb ID dr 3) At the ends of the tubes {z-0 & zj): •A dz = K( -T ) =q T s amb z=0 -k*L = KST -T ) s amb Ef fi dz where h ov represents an overall heat-transfer coefficient that incorporates both radiative and convective heat-transfer. The initial condition is given in Equation (4.7) and is based on the assumption that the sample starts at an initial temperature T , corresponding to the temperature of the tube from the t previous processing stage: ( A^ i T r =T @Ri<r<R , 0 0<z<z fi t=0 (4.7) Numerical Formulation The finite-difference method has been employed to solve the axisymmetric heat-transfer problem as described by the governing partial differential equation, boundary conditions and Chapter 4: Mathematical Modeling initial conditions. 50 The finite-difference equations for the seamless tube have been developed based on the heat balance on each node in both the radial and axial directions. An alternating direction implicit method (ADI) is utilized to solve them. This method permits the use of a tridiagonal matrix solver, which is applied twice in each time increment. In the first half of the time increment the implicit direction is in the radial direction. In the second half of the increment the implicit direction is in the axial direction. Figure 4.4 depicts the two implicit directions that are used to determine the tri-diagonal equations. r r z (a) 1* direction (b) 2 direction n d Figure 4.4 - Schematic of ADI directions for accelerated cooling FDM model. The seamless tube has spatially discretized, with a subsequent heat balance being applied to a control volume. Figures 4.5 and 4.6 illustrate the spatial discretization and node numbering employed in the model. The model allows for user definition of the number of nodes in the radial and axial directions. However, for thermal validation purposes a mesh of 11 nodes in the radial direction and 101 nodes in the axial direction (1111 nodes in total) has been used. r n o d e Outer diameter Ar: Az -> Z Inner diameter Figure 4.5 - Schematic of spatial discretization. Chapter 4: Mathematical Modeling 51 r .A LA [ . A .A -•z Figure 4.6 - Schematic node numbering. Generalized Heat Balance Approach to Development of Finite-Difference Equations In order to increase the robustness of the model a generalized heat balance approach has been used. The heat balance approach allows for a variety of different boundary conditions to be applied to a control volume (see section 4.2). The heat balance, given below, also incorporates the latent heat evolved from any phase transformations that may occur in the process under study. QA 1B 1C + + D Q^V +A + (4.8) = E^AV In Equation (4.8) the q's represent the energy entering or leaving the control volume associated with conduction, convection or radiation. The subscripts A, B, C, and D refer to the particular face of the control volume on which the particular energy is entering or leaving (refer to Figure 4.7). The term Q represents the latent heat evolved, E is the stored energy, and AV is the volume. The control volume consists of two surface areas in the radial direction (A&C) and two in the axial direction (fi&D), with the lengths of Ar and Az maintaining constant values. control volume control volume A r D As TO H_I Figure 4.7 - 3-D and 2-D schematic of control volume. B 52 Chapter 4: Mathematical Modeling Within the A C M model six different types of control volumes are present (Figures 4.8 4.13). These control volumes correspond to the various nodal positions present in Figure 4.5. A detailed example derivation (interior node) of the heat-transfer finite-difference equations that have been established, based on a heat balance on each node by an energy analysis, is presented in Appendix A. The Fourier number, Fo, has been used to simplify the heat-transfer calculations: Fo,= kAt (4.9) pc Al where k is the thermal conductivity, At is the time increment, p is the density, c is the specific p heat, and Al is the characteristic length. In addition, for nodes lying on the tube surfaces, the Biot number, Bi, which is proportional to (thermal internal resistance)/(surface film resistance) has been used: Bi. = hAx (4.10) where h is the heat-transfer coefficient and Ax the mid-plane distance. The corresponding derived equations (Equations 4.11 - 4.22) for both steps in the alternating direction implicit method are presented below: Interior node: (r=m, z=n) control volume Outer Diameter Ar ..Inner z Diameter Ar Figure 4.8 - Schematic of A C M model interior node and associated control volume. 53 Chapter 4: Mathematical Modeling 1 ADI Step: (f st P+,/2 = t + At/2) - radial direction is implicit and axial direction is explicit: p (i-^C+^(C+,+C>e, 2 2pc Ar r 2r_ _ 2 nd ADI Step: (t p+l Fo. =r P+,/2 At Ar P l/2 m-\,n + + P \/2 r (4.11) + + At/2) - axial direction is implicit and radial direction is explicit: Ar Ar 2r +(i-^)C +e ,/2 ~2~) At (4.12) Outer surface node: (r=M, z=n) Figure 4.9 - Schematic of A C M model outer surface node and associated control volume. 1 ADI Step: ( r st p+7/2 = f + At/2) - radial direction is implicit and axial direction is explicit: (l-tt>,)C + l + Fo. + T ^ ) 2 + Q ^ 2pc + Fo Bi T„ r r p Ar) „. r -—\ + Bi Um pP+l/2 Fo. M,n r ~ 1 M-\,n (4.13) 54 Chapter 4: Mathematical Modeling 2 nd ADI Step: (f" = f * +/ + At/2) - axial direction is implicit and radial direction is explicit: m _ Ar J J l-Fo. \ M , „ r + BL 2 J V Fo. 1 Ar M,n ' ' M,, At 2pc f (4.14) Inner surface node: (r=M°, z=n) Figure 4.10 - Schematic of ACM model inner surface node and associated control volume. 1 ADI Step: (f st = r + At/2) - radial direction is implicit and axial direction is explicit: +m P (l-F )T„ +^{T* 0 z v *' M ,n 2 ' + T> M Ar 1 + Fo. 2 nd ADI Step: (f l-Fo, P+/ .1+1 =f P+//2 \ + M ,n-\) Q - ^ + Fo Bi T„ r R 2pp Fo. +JBL (4.15) + zl//2) - axial direction is implicit and radial direction is explicit: V* Ar + 2 P+l/2 M',n T Fo. Ar^ At f>e 2 f Chapter 4: Mathematical Modeling = (i + Fo.)rM°,n 2 55 M°,n+\ 2 (4.16) M°,n-\ Outer surface corner node: (r=M, z=N) Figure 4.11 - Schematic of ACM model outer surface corner node and associated control volume. 1 ADI Step: {f* st = f + At/2) - radial direction is implicit and axial direction is explicit: m r ^ + ^ K i + i V l H n d ADI Step: (t l-Fo. P+1 ' =f 1 |f p+l/2 M,N \ r B Ar J-P+l/2 _ MJi r M,N r M,N 1 m rpP+X/2 (4.17) U-\,N 1 + At/2) - axial direction is implicit and radial direction is explicit: Fo. Ar^ P+i/2 T 1 \ + (Fo,Bi,+Fo Bi,)T Fo. Ar l + Fo. Bi. +• 2 +Q^- J J U,N \T V P+ * M.N = [l + Fo,(l + * r , ) K j - F O J Z % M,N r 2 I 2 M-1.N +Q At + {Fo Bi +Fo Bi )T r r z 2 a> 2pc (4.18) Note: A similar set of equations has been developed for the outer surface corner node on the other end of the tube. Chapter 4: Mathematical Modeling 56 Inner surface corner node: (r=M°, z=N) AqiDs Aqros Figure 4.12 - Schematic of A C M model inner surface corner node and associated control volume. 1 ADI Step: (f st = f° + At/2) - radial direction is implicit and axial direction is explicit: +m [l-Fo,(l + 1 + Fo. 2 nd ».)]T^ + ^ ADI Step: (f l-Fo. +1 =f + +1/2 ~ 1 1 M',N r Ar + V + (Fo Bi r Fo. f r M',N M',N +Fo,BiX r LV Af,JV i Ar\ T r p + V 2 M°+\,N (4.19) + At/2) - axial direction is implicit and radial direction is explicit: Fo. M°,N M° ,N rpP+\J2 Y Ar Bi. + Q-^- + Z T r + Q - At + {Fo Bi +Fo,Bi,)r r r (4.20) Note: A similar set of equations has been developed for the inner surface corner node on the other end of the tube. m Chapter 4: Mathematical Modeling 57 Surface end node: (r=m, z=N) Figure 4.13 - Schematic of A C M model surface end node and associated control volume. 1 ADI Step: {f st +m =f + At/2) - radial direction is implicit and axial direction is explicit: P [1 - Fo (1 + Bi f * , z + z 2r.m,N 2 n d ADI Step: ( r p+/ = f* Fo. 2r.m,N m,N r = [l + Fo,{l + m 2j FOJZ^ L ' m,N + Fo Bi T„ z ' 0 + z * awt,iV : | Q At ' m,N 2 (4.21) J * "HJf + At/2) - axial direction is implicit and radial direction is explicit: T y p+ 2 +1 r +^L)T P Bi ))T™-Fo T™_ z z x : + (l - F o ^ j ^ + Fo Bi T 2 r z z a0 +Q At (4.22) Note: A similar set of equations has been developed for the inner surface corner node on the other end of the tube. Dynamic Time-Stepping Routine In order to increase execution efficiency within the program, a dynamic time-stepping or adaptive iteration subroutine has been added to the model to determine the appropriate time increment (At) needed for solution convergence. The subroutine begins by calculating the maximum temperature, Tmax, and volume fraction, Vfmax, difference for all the nodes and subsequently determines which flags to turn on (i.e., if time step needs to be increased, Chapter 4: Mathematical Modeling 58 decreased, or remain constant). If need be, the subroutine adjusts the time increment according to the following specified criteria: • If the maximum values are less than 1/3 the user defined T max or Vf values then it or Vf values then it max increases the time increment. • If the maximum values are greater than the user defined T max max decreases the time increment by a half. If the time increment is modified by either of the aforementioned criterion, then the subroutine restarts the solution techniques at the end of the last valid time increment (i.e., resets the time increment) using the newly modified time increment. Otherwise the current time increment is maintained and the solution has met convergence stability criterion. 4.1.2 Controlled Slow Cool Model The Controlled Slow Cool (CSC) Model has been developed with the premise that the tubes are entering the cooling regime side by side (no vertical stacking) with the axial length of the tube perpendicular to the direction of motion. Like the A C M , the C S C couples a heat diffusion model with a microstructure evolution model. However, in the C S C model the geometry of the analysis is different in that radial and circumferential conduction have been incorporated into the model and axial conduction has been ignored. Conduction in the axial direction was neglected on the postulate that the tube is experiencing relatively uniform boundary conditions in the axial length, except at the tube ends. Since the axial length is relatively long in comparison with the tube diameter and wall thickness, end effects are assumed to be minimal and; therefore, can be neglected. The model allows for multiple cooling zones that permit variable tube speeds, variable zone temperatures, and boundary conditions. The influence of adjacent tubes has also been incorporated into the model. As in the A C M , a dynamic time step approach has been included in the CSC model to improve computational speed and accuracy of results. Geometry A circular cross-section has been used to represent the geometry of the seamless tube within the CSC model, owing to the aforementioned assumptions (Figure 4.14). The circular cross-section consists of two sides (representing the inner and outer diameter of the tubes) to Chapter 4: Mathematical Modeling 59 which boundary conditions may be applied. To increase versatility, the model permits user definition of the inner and outer diameter. Figure 4 . 1 4 - Schematic of C S C geometry. Heat Flow For the C S C model, the gradient and Laplacian operator within the transient heat conduction Equation (4.1) are expressed in cylindrical coordinates as follows: rdr dr X J r 80 2 (4.23) 80 The applicable boundary conditions for the solution of Equations (4.1) and (4.23) are defined as follows: 1) At the outside surface of the tube (r = R„), at times greater than 0 (/>0), -k*4 dr r=R„ = ^(.T-T ) a =q OD @r = R , 0 r>0 (4.24) 2) At the inner surface of the tube (r = R,), at times greater than 0 (/>0), -k 8T\ 8r =0=q ID @r = R h f>0 (4.25) r=R, where h ff represents the effective heat-transfer coefficient due to radiation and convection. e Chapter 4: Mathematical Modeling 60 The initial condition is given in Equation (4.26) and is based on the assumption that the sample starts at an initial temperature, Ti, corresponding to the temperature of the tube from the previous processing stage: T(r,0\^=T, @R <r<R , t 0 0<0<2K, t=0 (4.26) Numerical Formulation The fmite-difference heat-transfer equations for the C S C model, similar to the A C M model, have been established by completing an energy balance on each nodal control volume. The finite-difference approximations have been used to calculate the heat balance in both the radial and circumferential directions. A n ADI method has been used to solve these equations. The time step is divided into two increments. As shown in Figure 4.15, in the first Vz time-step a tri-diagonal solver is used to solve the implicit equations for conduction in the radial direction. For the second A time-step, a full matrix solver is used to solve for conduction in the l ckcumferential direction. The decision to use an ADI method and divide the time step into two increments, rather than utilize one pass through a full matrix solver (i.e. solve conduction in both directions in one time step), was made to mmimize the computation time. Using this method, it was found that the computation time was reduced by 3 to 5 times, depending on the number of nodes, compared to one-pass through the full matrix solver. It is important to note that the solution difference between the two methods was found to be negligible. 61 Chapter 4: Mathematical Modeling (a) 1 direction - radiial st (b) 2 direction - circumferential nd Figure 4.15 - A D I implicit directions for Controlled Slow Cool F D model. The nodal volumes on which the energy balances have been completed are defined by spatially discretizing the seamless tube about a circumferential cross section. Nodal locations are defined based on an origin at the centerline of the tube. Figure 4.16 depicts the spatial discretization and node numbering. The model allows for user definition of the number of nodes in the radial and circumferential directions. For thermal validation purposes, a mesh of 5 nodes in the radial direction and 36 nodes in the circumferential direction (generating 180 nodes and 180 elemental control volumes, with node 1 beginning at the inner diameter and 0°) has been used. Chapter 4: Mathematical Modeling 62 y t AO Figure 4.16 - Schematic of spatial discretization and node numbering. Generalized Heat Balance Approach to Development of Finite-Difference Equations The heat balance approach used in the development of the finite-difference equations, presented below, incorporates the latent heat evolved from phase transformations and allows a variety of different boundary conditions to be applied to a control volume (see section 4.2). 4A+1 +<1C+<1D+ b = E.AV (4.27) The control volume is treated as having four surface areas through which heat can flow. Figure 4.17 shows the control volume consisting of two surface areas in the radial direction (A and C) and two in the circumferential direction (B and D), with the length of Ar and the angle A0maintaining constant values. Control Figure 4.17 - 3-D & 2-D schematic of control volume with respect to spatial nodes. Chapter 4: Mathematical Modeling 63 Within the CSC model, three different types of control volumes are present (Figures 4.18 - 4.20). These control volumes correspond to the various nodal positions present in Figure 4.16. A detailed example derivation (interior node) of the heat-transfer finite-difference equations that have been established based on a heat balance on each node by an energy analysis is presented in Appendix B. As with the A C M model, the CSC model uses both the Fourier and Biot numbers (Equations 4.9 and 4.10) to simplify the calculations within the model. The corresponding derived equations (Equations 4.28 - 4.33) for both steps in the alternating direction implicit method are presented below: Interior node: (r=m, 0=n) r control volume Inner Diameter Figure 4.18 - Schematic of C S C model interior node and associated control volume. 1 ADI Step: (f st +1/2 = f + At/2) - radial direction is implicit and circumferential direction is explicit: 2r. -^V r rn n t A 2 2 nd ADI Step: (f explicit: P+I = f* m + V 2 +(r +— m—\,n J m,n V (4.28) » 2, + At/2) - circumferential direction is implicit and radial direction is Chapter 4: Mathematical Modeling Fo. 2r_ 64 Ar m,n r 1 2 Ar^ j>P+\/2 m-\,n m+\,n ~2j At +0-^OC +e P ,/2 2 c f (4.29) 7n«er surface node: (r=Af, &=n) T control volume Outer Diameter Inner Diameter Figure 4.19 - Schematic of C S C model inner surface node and associated control volume. 1 ADI Step: (t st p+I/2 = f + At/2) - radial direction is implicit and circumferential direction is explicit: Ar) 1 + Fo. yrP+1/2 _ ~2J 2 nd ADI Step: (f* = r 1 P+m Fo. M'.n M°,n 2 ) M °- * 1 (4.30) + At/2) - circumferential direction is implicit and radial direction is explicit: l-Fo. Ar + Bi. P+V2 T M°,n Fo. Ar • B + 2 At Chapter 4: Mathematical Modeling - { 1 -I- Fn \ \T ° r F P+X 65 - F AC- ° r A * « ° , » rpP+\ (4.31) Ow/er Surface node: (r=M, 0=n) ]ODS control volume A* i i i i i i i 1 \ . i j- i i i ; i AO { AO Inner Diameter Figure 4.20 - Schematic of CSC model outer surface node and associated control volume. 1 ADI Step: (f explicit: st = f + At/2) - radial direction is implicit and circumferential direction is +1/2 (l-Fo^Vi, + T ^ ) + 2 f i +1 =f p+//2 ir ^ r ,\Or M,n l liJ r p+1 r _ ArV 2 i , EV, rU Fo Bi T„ P p + 1 / 2 (4.32) + At/2) - circumferential direction is implicit and radial direction is 1-Fo. L + C ' lyrP+l/2 M,n r 2 ADI Step: (f explicit: nd 2 f l + Fo. Q - * P «,» ' r A< ) ^ 77p+i _Fo. % ^ A * _,p+i 2" 7 A/,n-l J 2 pc (4.33) 66 Chapter 4: Mathematical Modeling Dynamic Time-Stepping Routine Owing to the fact that the structural configurations of the A C M and C S C models are similar, an analogous treatment of the dynamic time-stepping routine has been preformed. 4.2 Implementation of Boundary Conditions Within the models, different boundary conditions can be selected for the various cooling zones, each characterized by its length and the tube traveling speed. For the numerical analysis, the heat flow at the surface has been characterized with a Cauchy-type boundary condition: 9 = (4.34) K*.-TJ) where T (K) is the surface temperature, Tomb (K) is the ambient temperature (air or water), and h s the heat-transfer coefficient in W/m /K. 2 The boundary conditions considered include: no heat-transfer (default), radiation, natural air-cooling, forced air-cooling, and water spray cooling. Natural air-cooling has been used for all periods where the tubes were not in contact with the water sprays or the forced air jets. A n overall heat-transfer coefficient, hov, has been used to account for heat-transfer in those cases where there is a combination of convection (natural, forced, and water) and radiation: (4.35) where h is heat-transfer coefficient for convective heat flow and h is the heat-transfer c coefficient for radiation. r The methodology for evaluation of h and h are outlined in the c r following sections for each of the heat-transfer regions. Convection The values of the convective heat-transfer coefficient, h , in the Timken seamless tube c process have been characterized by a variety of different equations depending on the method of cooling. Chapter 4: Mathematical Modeling 67 Natural & Forced Air Convection Two equations are available within the models to describe the heat-transfer coefficient for convective losses due to natural convection. The first developed by Golja , models a horizontal [2] tube that heats the surrounding air: D 03^@L(T -T ) S (4.36) amb where D represents the outer diameter, k the thermal conductivity, /? the spatial coefficient of propagation, g the gravitational acceleration, v the kinematic viscosity of air, and the h values c for this correlation are around 11-19 W/m /K. 2 Kreith & Black [3] also developed an equation incorporating Prandtl's number Pr to characterize convective losses due to natural convection over horizontal cylinders: D c 0.6 + 0.387 (4.37) 1+ 0.56 Pr where the h values are around 15-22 W/m /K. 2 c The correlation developed by Golja has been used for the experimental validation of the controlled slow cooling bed (section 6.2). However, further experimental tests are required before making a choice on the final correlation. For forced air convection, the flow has been postulated to be similar to that of flow across a horizontal cylinder. The resulting heat-transfer coefficient characterizing the convective losses has been calculated from the Equation (4.38) : [3] 0.466 / x - h = 0.683|— (4.38) where h values are around 14-96 W/m /K. 2 c 68 Chapter 4: Mathematical Modeling The properties of air used for both natural convection and forced-air convection have been determined using a mean film temperature: (4.39) and the thermodynamic property data for dry air has been obtained from literature and is presented in Table 4.1. Table 4.1 - Thermodynamic properties of dry air at atmospheric pressure. Temperature T (K) Density 273 293 313 333 353 373 473 573 673 773 1273 1.252 1.164 1.092 1.025 0.968 0.916 0.723 0.596 0.508 0.442 0.268 Specific Heat P , C (kg/in) 3 P (J/kg/K) 1011 1012 1014 1017 1019 1022 1035 1047 1059 1076 1139 Thermal Conductivity Prandtl Number Pr gpVxio-* (W/m/K) Absolute Viscosity fixlO (N-s/m ) 0.0237 0.0251 0.0265 0.0279 0.0293 0.0307 0.037 0.0429 0.0485 0.054 0.0762 17.456 18.24 19.123 19.907 20.790 21.673 25.693 39.322 32.754 35.794 48.445 0.71 0.71 0.71 0.71 0.71 0.71 0.71 0.71 0.72 0.72 0.74 1.85 1.36 1.01 0.782 0.6 0.472 0.164 0.0709 0.035 0.0193 0.00236 k 6 (1/K/m ) 3 2 Spray Water Cooling To model the water spray cooling (accelerated water cooling), spray-boiling data has been used to characterize heat-transfer coefficients as a function of temperature. The data included a series of boiling curves that represented different spray configurations. Practically speaking, the heat-transfer coefficient is affected by the quench medium, nozzle diameter, spray velocity and volumetric flow rate. The data for these boiling curves has been acquired from two sources; literature' ' 4 61 and experimental tests performed at Timken . [5] The heat-transfer correlations representing the heat-transfer coefficient as a function of surface temperature for the boiling curve regions are presented in Table 4.2. 69 Chapter 4: Mathematical Modeling Table 4.2 - Heat-transfer correlations for spray boiling data. Boiling Curve now Rate (L/m /s) 2 TKRtest-R613R 9.7 T K R test - R64R1 13.45 T K R test-R610Rl 15.4 Curve A 2.02 Curve B Curve C Curve D Curve E 1.0 0.6 4.99 5.0 Heat-transfer Coefficient (W/m /K)* Surface Temperature 84983(l-(80-T)/55) 5926+4348147/(T-25) 114290-7571750/(T-25) 23077+1203 8463/(T-25) -36283+36376000/(T-25) 89091(l-(80-T)/55) 11250+4281250/(T-25) 105500-15982500/(T-25) -21226+23936350/(T-25) -40769+40450000/(T-25) 156556(l-(80-T)/55) 13889+7847225/(T-25) 69444-763800/(T-25) -33594+33754150/(T-25) 8.4287T+5398.5 1236.5T-466381 -86.48T+61703 0.0018T+1498.8 1.8295T+4901.7 674.39T-246063 -62.694T+43720 -0.0027T+1456.8 3.445T+2521.4 583.87T-212905 -54.565T+36182 0.0052T+1467.3 -0.0989T+20447 1489T-566490 -134.12T+108982 0.0013T+2199 38.46T+2535.3 1468.9T-556973 -136.84T+103220 0.0069T+1792.6 25 < T < 80 8 0 < T < 135 135<T<240 240 < T < 435 435 < T < 1000 25 < T < 80 80 < T < 240 240 < T < 340 340 < T < 870 8 7 0 < T < 1000 25 < T < 80 8 0 < T < 180 180<T<360 3 6 0 < T < 1000 288.15<T< 384.15 384.15 <T<399.15 399.15 <T<696.15 T > 696.15 288.15 <T<373.15 373.15 < T < 393.15 393.15 < T < 674.15 T > 696.15 288.15 < T < 371.15 371.15<T<390.15 390.15 <T<636.15 T > 636.15 288.15 < T < 394.15 394.15<T<416.15 416.15<T<796.15 T > 796.15 288.15 <T<391.15 391.15 < T < 411.15 411.15 < T < 741.15 T>741.15 2 Ref. [5] [5] [5] [4] [6] [4] [6] [4] [6] [4] [6] [4] [6] *T is in CC) for the Timken data and (K) for the literature data (curves A-E). In an attempt to account for the increase in the heat-transfer coefficient in the film-boiling regime, a higher spray flux has been assumed to modify the literature data curves . These data [6] curves have been generated using a series of polynomial regressions. Five boiling curves have been developed using literature data and three more from the experimental data obtained from Timken. The graph in Figure 4.21 depicts all eight boiling curves used to characterize the heattransfer coefficient associated with spray cooling. 70 Chapter 4: Mathematical Modeling J 200000 I I Boiling Curve Flow Rate Ref. (L/mVs) CM E R613R R64R1 R610R1 A B C D E 150000 H c o § I 1_ 100000 o 0 L. X 50000 H (0 0) •A • • 100 9.70 13.45 15.40 2.02 1.00 0.60 4.99 5.00 [5] [5] [5] [4,6] [4,6] [4,6] [4,6] [4,6] R613R R64R1 R610R1 —A B C D E 1000 Surface Temperature (K) Figure 4.21 - Spray boiling curves. These boiling curves consist of four different regions (see Figure 4.22): 1) a natural convection region where there is no boiling at the heated surface, 2 ) a nucleate boiling region where isolated bubbles form at nucleation sites and separate from surface, 3 ) a transition boiling region where partial film boiling occurs and the heat-transfer coefficient decreases as the surface temperature increases, and 4) a film-boiling region where heat-transfer from the surface to the liquid occurs by radiation through the vapor phase and the heat flux increases with increasing change in temperature. The critical heat flux (CHF) point is important in most industrial applications because it is the point of maximum heat flux. The heat flux is at a minimum at the Liedenfrost point. 71 Chapter 4: Mathematical Modeling a | X 3 55 CD log Temperature (K) Figure 4.22 - Schematic representing boiling curve regimes. w Radiation Radiation occurring during the cooling process has been accounted for by linearizing the radiation equation so that the heat-transfer rate is proportional to a temperature difference where the effective heat-transfer coefficient is given by: K=™(T,+T iT, T] ) amb + mb (4.40) where s is the emissivity of the steel and a Stefan Boltzmann's constant. This allows a Cauchytype boundary condition to be used. Owing to its dependence on the local temperature to the third power, the radiative heattransfer coefficient is highly non-linear. It must therefore be updated for each time step and a value of the local temperature at the time step being evaluated is required. One method of evaluating the radiative heat-transfer coefficient is to iterate within the time step. Evaluate h on r the basis of an initial value (i.e., the previous time-step temperature), solve the current time step, then re-evaluate h based on the new solution and compare successive solutions until a r convergence criterion is attained. However, this method is computationally intensive so an alternative approach has been used in the models. The method used has been to simply base the Chapter 4: Mathematical 72 Modeling evaluation of the radiative heat-transfer coefficient on the previous time step. By restricting the maximum temperature change at any node within one time step to sufficiently small value (as implemented by the dynamic time-stepping routine) the error introduced by this approximation is minimized and the overall computation time is reduced over that of an iterative approach. For the CSC model, the heat-transfer on the outer diameter surface of the tube due to radiation is more complex. To account for the influence of adjacent tubes, the heat-transfer coefficient due to radiation has been multiplied by a radiation shape factor (Fi)- This radiation shape factor has been approximated using a simplified view factor approach. In this approach, three different views have been considered at a given point/node on the tube surface: 1. The view of the furnace top. 2. The view of adjacent tubes. 3. The view of the furnace bottom. The view factor approximation is determined based on the distance between tube centers, radii of the tubes, and angular position of the nodes. A schematic representing the underlining methodology used ih determining the radiation shape factor is depicted in Figure 4.23. Top z Bottom Calculations: V= r * cos(a) u=L - v h = r * sin(a) Figure 4.23 zeta= SQRT(u +h) 8/2 = sin" (r/zeta) P = 90 - sin" (v/r) 2 1 1 2 fi y = cos^OVzeta) = y-8/2-0 appa = 180-fi-8 View factors: Top: F, = appa/180 Other Tube: F = 9/180 Bottom: F = fi/180 2 3 Schematic illustration depicting the view factor approximation. Chapter 4: Mathematical 73 Modeling The effective heat-transfer coefficient due to radiation is subsequently calculated using the following relation: K= +T XT amb 2 s +T ) 2 (4.41) a mb i=\ where F , is the radiation shape factor for each view / and T amb corresponding to each view i (e.g. T ambi temperature, and T is the ambient temperature - furnace top temperature, T ^ - adjacent tube surface amb ^ - furnace bottom temperature). amb 4.2.1 ACM Model Using a time on/off method applied to either the outer diameter (OD) and/or the inner diameter (ID) nodal control volumes, boundary conditions (adiabatic, radiation, convection, or radiation and convection) are implemented. The determination of the appropriate boundary condition is based on the location of the various nodes at a given time with respect to the specific cooling regimes. A default boundary condition of natural-convection/radiation has been applied to the tube ends and the nodes lying outside of the cooling regimes. 4.2.2 CSC Model Due to uncertainty in the exact configuration that the controlled slow cooling process would take, the model has been coded to accommodate multiple furnace zones. The zones have been configured as having separate tube velocities where the tube is assumed to be moving horizontally at a given velocity. It is also assumed that the tube rolls without slipping. Thus, in Figure 4.24, the center of the tube is assumed to move horizontally from C to C and the nodes move in a sigmoidal path. The appropriate zone and boundary conditions are determined by tracking the horizontal and angular position of the various nodes at any given time. The number of revolutions the node makes, along with its angular position, determines the location of these nodes within the cooling zones. Again, a time on/off method has been applied to the boundary conditions (radiation, convection, or radiation and convection) to the OD nodal control volumes. boundary condition has been assumed for the ID surfaces. A n adiabatic 74 Chapter 4: Mathematical Modeling co s= i6 Figure 4.24 - Depiction of tube motion. In addition, the model facilitates zones with varying top and bottom temperatures. This has several implications with regard to conductive and radiative heat-transfer, particularly the influences of neighboring tubes. In order to reduce the complexity of the process, no conductive heat-transfer between the tube and zone bottom (i.e., the support structure below the tubes) has been assumed, owing to the relatively small contact times and contact surface areas between the tube and the bottom surface. A default boundary condition of radiation has been applied to the nodes lying outside of the cooling regimes. 4.3 Coupling to Microstructure Models The microstructural transformation has been modeled using transformation kinetics subroutines supplied by the Timken Company. Two different transformation kinetics subroutines have been integrated into the models, with an option for the user to select the appropriate subroutine depending on steel grade. The first, being steel grade specific, has been developed through collaborative efforts between the Timken Company and the Colorado School of Mines. It utilizes a state variable approach in which the transformation kinetics are tailored for individual steel grades. The second, a chemistry dependent approach, has been developed to add flexibility to the simulation tool by utilizing steel chemistry to detenriine the transformation kinetics. This subroutine has been developed through a joint venture between the Timken Company, Questek, and Minitech Limited. By substituting the appropriate heat of transformation and fraction transformed for each transformation reaction, the total heat released during a given time step can be calculated. 75 Chapter 4: Mathematical Modeling 4.4 Material Properties To conduct the analysis of the seamless tube cooling, the thermal-physical properties of the steel (i.e., density, specific heat, conductivity, and emissivity) are required as inputs into the models. This data was obtained from literature , as well as from experimental data supplied by [7] the Timken Company . Table 4.3 shows typical model data for an AISI 5130 low alloy steel. [8] Table 4.3 - Thermal-physical property data for an AISI 5130 low alloy steel. Property Density (kg/m ) Symbol" PA,F,P,B,M 3 Thermal conductivity (W/m/K) k k kp,B,M CpA A F Value Source 7800 Timken 11.363+1.5586xlO" T 65.422-5.2176x 10" T+9.7673x 1 O^T 50.742-3.0567x 10" T+l. 1539x 10" T 2 2 2 7 AJSI 2 [7] 2 674.73-0.24371T+2.0265xl0^T (T > 787°C): -10034.5+5.9668T+5.2xl0 T (787°C < T< 769°C): 34754.5-31.9196T (769°C < T< 727°C): -11462.6+12.4346T CpF (727°C < T< 527°C): -4704.5+4.568T+1.1057xl0 T" (527°C < T): 503.13-0.13068T-51.702xl0 r +4.47xl 0"*T 449.5+0.4501T C P,B,M 0.75 £A,F,P,B.M 2 9 Specific heat (J/kg/K) [8] AISI l7J 2 9 2 AISI . [7] 5 2 2 AJSI Emissivity Timken The subscripts A, F, B, M refer to austenite, ferrite, bainite and martensite respectively. i;] P l8J Two different approaches have been utilized to incorporate this data into the models. The first uses a simple look-up table approach in which tabulated data, relating the thermal-physical properties to grade, structure and temperature, is entered into a model input data file. The second utilizes polynomial equations (supplied by Timken) to characterize the thermal-physical property data. In both of these cases, linear interpolation techniques have been used to retrieve the appropriate data for the given temperature or phase. In the application of thermal-physical properties to mathematical models where there is a mixture of phases present, a simple rule of mixtures for any property P has been applied: M phase i=i where the volume fraction and property of phase /' are denoted by X and P, respectively. t Chapter 4: Mathematical 76 Modeling For the numerical analyses, the densities of the steel grades have been held constant (i.e., independent of temperature), owing to the postulate of a static mesh and consequently allowing for the conservation of mass during thermal changes. 4.5 Latent Heat of Evolution The latent heat per unit volume, Q, represents the evolved/consumed heat energy when there is an exothermic or endothermic phase transformation (i.e., austenite decomposition). In the austenite temperature field, prior to any transformation, Q is set to zero during the numerical solution. The evolution of heat is considered as soon as the transformation reactions begin and the latent heat is evolved. The amount of heat released during each time step in the numerical calculations depends on the rate of austenite decomposition into the various transformation products during a given section of the cooling history. Thus, the heat generated during phase evolution has been calculated using a discretized form: (4.43) where p is the density, At the time increment, H the enthalpy of transformation of phase i, and t AXj the volume faction of phase i. There is a latent heat associated with the formation of the various transformation products: ferrite, pearlite, bainite and martensite. The latent heats vary with temperature due to the temperature dependence of the specific heat for each of the phases. The enthalpies of formation have been characterized using polynomial equations and the latent heats have been integrated into the model by way of transformation kinetics subroutines. enthalpies of transformation in steel are presented in Table 4.4. The values for the Chapter 4: Mathematical 77 Modeling Table 4.4 - Enthalpies of transformation in steel. Phase Transformation Enthalpy of Transformation (J/kg)* 221656.4-864.4T+1.9795T'0.001478T 3 2.917xl0 +114590T-148.8T + 0.06399T 3277373-10575T+11.545T 0.00424T 70651+225.23T0.3469T +6.775xl0" T 7 austenite - ferrite Source T<720 2 3 720 < T < 780 Timken [8J Timken [8] 2 3 austenite - pearlite 2 5 austenite - bainite 92000 austenite - martensite 83600 T>780 3 T is in degrees Celsius. 4.6 Summary Two 2-D transient heat diffusion models have been developed that are capable of predicting the thermal history within seamless tubes during two processing stages: 1) an accelerated cooling model (ACM) to model rapid water quenching, such as might occur between deformation stands to control austenite grain growth; and 2) a controlled slow cool model (CSC) to model slow controlled cooling, such as might occur at the end of the process to produce the final desired microstructure. Both models have been coupled with a microstructure evolution model (incorporating recrystallization kinetics, grain growth kinetics, and austenite decomposition kinetics) that has been provided by the Timken Company. These models have been coded in F O R T R A N to run on a personal computer and use finite-difference methods to solve the general transient heat conduction equation, with an alternating-difference-implicit scheme to solve the system of equations. A summary of the key aspects in the A C M and C S C models is presented below: The A C M model employs an axisymmetric configuration, owing to the premise that the tubes enter the cooling stands separately, with their axial lengths parallel to the direction of motion. The model has been designed as a flexible platform that allows for moving boundary conditions, multiple cooling zones, and on/off switches that can be applied to the outer diameter (OD) and/or the inner diameter (ID) of the tube. Chapter 4: Mathematical 78 Modeling Unlike the A C M model, the configuration of the CSC model is different. The (CSC) model employs a 2-D circumferential geometry where radial and circumferential conduction have been incorporated into the model and axial conduction has been ignored. This stems from the premise that the tubes enter the cooling regime side-by-side (i.e., no vertical stacking) with the axial length of the tube perpendicular to the direction of motion. Conduction in the axial direction has been neglected on the postulate that the tube is experiencing relatively uniform boundary conditions in the axial length, except at the tube ends. Since the axial length is relatively long in comparison with the tube diameter and wall thickness, end effects are assumed to be minimal and can be neglected. The model allows for multiple cooling zones that permit variable tube speeds, variable zone temperatures, and boundary conditions. For the purpose of increasing robustness in the A C M and CSC models, a generalized heat balance approach (incorporating the latent heat evolved from any phase transformations that may occur in the process) has been used. This approach allows for a variety of different boundary conditions to be applied to the control volumes. The boundary conditions describing the heat flow at the surface have been characterized by a Cauchy-type boundary condition. A n overall heat-transfer coefficient, h , accounts for ov heat-transfer in cases where there is a combination of convection (natural, forced, and water) and radiation. Natural and forced-air convection have been characterized by a variety of empirical correlations, that are functions of Re, Pr, and Ra numbers, to describe the heat-transfer of air flowing over horizontal tubes. Water spray cooling (accelerated water cooling) has been modeled using boiling curve data (for a variety of different spray configurations) to characterize the heat-transfer coefficients as a function of temperature. Radiation occurring during the process has been accounted for by linearizing the radiation equation so that the heat-transfer rate is proportional to a temperature difference. In the CSC model, the heat-transfer on the outer diameter surface due to radiation is more complex. Consequently, a geometric view factor approximation has been incorporated into the model to account for the influence of adjacent tubes and the surrounding environment. Chapter 4: Mathematical 79 Modeling A dynamic time-stepping subroutine has been added to the models to increase execution efficiency. The solution algorithm of this subroutine takes into account a maximum temperature and volume fraction difference in order to determine the appropriate time increment needed for solution convergence. The incorporation of the material property data has been conducted using two different approaches: 1) a simple look-up table approach in which tabulated data, relating the thermalphysical properties to grade, structure, and temperature, is entered into a model input data file, and 2) polynomial equations to characterize the thermal-physical property data. Linear interpolation techniques have been used to retrieve the appropriate data for the given temperature or phase. In those cases where there is a mixture of phases present, a simple rule of mixtures for any property P has been applied. The latent heat associated with the formation of the various transformation products (ferrite, pearlite, bainite and martensite) is temperature dependent. Hence the enthalpies of formation have been characterized using polynomial equations and the latent heats have been integrated into the model by way of the transformation kinetics subroutines. The models have been designed to predict the thermal history and microstructure evolution of tubes and are capable of simulating a wide range of cooling conditions (from the soft cooling associated with natural and forced-air convection to a hard quench utilizing water quench stands). 80 Chapter 4: Mathematical Modeling References 1. R.W. Lewis, K. Morgan, and O.C. Zienkiewicz, Numerical Methods in Heat Transfer, John Wiley & Sons, New York, 1981, pp 1. 2. M . Golja, "Mathematical Model of Steel Hollow Shell Cooling after Rolling", Metalurgija, Vol. 36, No. 4, 1997, pp 229-234. 3. F. Kreith and W.Z. Black, Basic Heat Transfer, Harper and Row, New York, 1980. 4. I Mudawar and W.S. Valentine, "Determination of the Local Quench Curve for SprayCooled Metallic Surfaces", Journal of Heat Treating, Vol. 7, 1989, pp 116-117. 5. D. Li, U B C - Timken Internal Report, 2000. 6. J.K. Brimacombe, P.K. Agarwal, S. Hibbins, B. Prabhaker, and L . A . Baptista, "Spray Cooling in the Continuous Casting of Steel", Iron and Steel Engineer, pp 109-132. 7. AISI Project. 8. D.Q. Jin, UBC-Timken Internal correspondence, 2000. 9. Y . Nagasaka, J.K. Brimacombe, E.B. Hawbolt, I.V. Samarasekera, B. Hernandez-Morales, and S. E . Chidiac, "Mathematical Model of Phase Transformations and Elasto-Plastic Stress in the Water Spray Quenching of Steel Bars", Metallurgical Transactions A , Vol. 24A, 1993, pp 795-808. Chapter 5: Experimental Tests 81 CHAPTER 5 EXPERIMENTAL TESTS Experimental trials were conducted in collaboration with the Timken Company at their facilities in Canton, Ohio. A controlled slow cooling bed simulator was constructed, in which steel tubes and a furnace hood were instrumented with thermocouples, and time-temperature histories were recorded during various slow cooling conditions. The simulator was designed to investigate the influences of inter-tube spacing and the surrounding cooling environment on cooling history and microstructure evolution. The data collected were used to verify the C S C model formulation. 5.1 Material The experimental work was carried out on a commercially significant low alloy, chromium steel grade, AISI 5130 that was supplied by the Timken Company. A typical chemical composition for AISI 5130 steel is shown in Table 5.1. Table 5.1 - Nominal composition of AISI 5130 alloy steel. ' 1 1 %c %Cr 0.3 0.95 %Mn 0.8 %Si 0.21 %S 0.04 %P 0.04 The experimental runs were conducted on as-rolled 5130 steel tubes with dimensions of 0.08575 m outer diameter (OD), 0.0138 m wall thickness, and 0.914 m tube length. 5.2 Controlled Slow Cooling Bed Simulator The experimental program was performed using a custom system designed and built at the Timken facilities. This system included: a gas furnace to preheat the steel tube samples to the desired initial temperature; a test rig used to hold the tubes; an insulated furnace hood with a removable front door to simulate the cooling bed; natural gas burners to preheat the cooling bed prior to sample insertion; and a data acquisition (DAQ) system that was used to record the thermal history within the samples and the surrounding cooling bed environment. A schematic of the controlled slow cooling bed apparatus is shown in Figure 5.1. Chapter 5: Experimental Tests 82 Instrumented tube sample Furnace Data Acquisition (DAQ) System Furnace Hood without door Test Rig / A Furnace Hood with door Data Acquisition (DAQ) System Data Acquisition (DAQ) System Figure 5.1 - Schematic of controlled slow cooling bed system used to cool tubular specimens. Different cooling conditions were examined by: 1) varying the distance between adjacent tubes, 2) including or excluding use of the furnace hood, and 3) using the furnace hood with or without the door. The test rig was designed to accommodate three tubes, less than a meter in length, and allowed for variation in tube spacing and rotation (shown in Figure 5.2). Figure 5.2 - Photograph of the experimental test rig. Chapter 5: Experimental Tests 83 The furnace hood consisted of a steel shell insulated with 0.102 m of FiberFrax* insulation (see Figure 5.3). The removable door was added to facilitate the introduction of the test rig into the furnace hood and to allow for boundary condition manipulation. The top of the furnace hood was instrumented with five type-K thermocouples as shown in Figure 5.3. Thermocouples 1,3, and 5 were used to determine the furnace hood temperature (located at the surface of hood insulation) and thermocouples 2 and 4 measured the furnace air temperature. Tube rn rn i Tube entry J • Thermocouple fS\ 4 i tf 5 ffi 3 H Distance above tube © 0.190 m ® 0.082 m 2 1®1 LJ LJ Insulation Tube entry J Figure 53 - Photograph of the furnace hood and schematic of thermocouple positions. Prior to the insertion of the test rig containing the experimental tubes, natural gas burners were used to preheat the furnace hood as shown in Figure 5.4. " F I B E R F R A X i s a trademark o f the Standard O i l Engineered Materials C o m p a n y , N i a g a r a Falls, N Y . 84 Chapter 5: Experimental Tests Figure 5.4 - Photograph of furnace hood preheat. 53 Sample Instrumentation The temperature history of the as-rolled tubes was obtained from thermocouples located in the center tube of the three-tube configuration. This center tube was instrumented with six type-K intrinsic thermocouples having a wire diameter of 0.25 mm and a stainless steel cladding of 1/16 inch in diameter. Clad thermocouple wires were fed along the center of the tube and brought to the surface through 5/16-inch holes. Center punching around the metal cladding secured the wires to the tube. Mullite insulation tubes were placed over exposed C H R O M E L A L U M E L * wires, and the wires were inserted into each of six 1/16-inch diameter flat bottom holes that were drilled at different depths and locations along its length (see Figure 5.5). Thermocouples 1,3, and 5 were positioned at a depth of approximately % of the wall thickness (ranging from 3.2 mm to 4.35 mm) and thermocouples 2, 4, and 6 at a depth of approximately A X of the wall thickness (ranging from 6.55 mm to 7.35 mm). The variation in thermocouple position along the length of the tube allowed an investigation into tube end effects. Care was taken to provide sufficient space between the thermocouples so they did not influence the thermal profile of the sample. The radial variation in thermocouple placement was done to measure any through-thickness thermal profile for the samples. Thermocouples were placed far enough away from the surface, so that large surface variations/irregularities (i.e., scale ' C H R O M E L - A L U M E L is a trademark of the Hoskins Manufacturing Company, Hamberg, MI. Chapter 5: Experimental Tests formation) did not impact on the temperature measurements. 85 The CSC model node positions that correspond to the thermocouple positions for a 1 7 X 9 node mesh model are depicted in Figure 5.5. Figure 5.5 - a) Photograph of instrumented 5130 steel specimen depicting thermocouple positions, b) schematic of corresponding C S C model nodal/thermocouple positions. 5.4 Controlled Slow Cool Test Runs In order to validate the C S C model, a series of experimental tests were run for the asrolled 5130 steel tube specimens at the Timken facilities in Canton, Ohio. Six tests were designed (involving three tubes) to simulate various controlled slow cooling conditions. The test matrix for the 5130 specimens is presented in Table 5.2. For each test, the sample tubes were first heated to the desired temperature (approximately 1125°C) in an electric furnace and held for about A hour at temperature to ensure L uniform heating. The samples were then quickly transferred onto the test rig with tongs. During cooling, the D A Q system recorded the sample temperature as a function of time at a frequency of 1 Hz. To ensure precise thermal histories for each sample, it was important that the thermocouples made good contact with the sample and had fast response time. To make certain Chapter 5: Experimental Tests 86 that the thermocouple wire made good contact with the sample, the wires were separately spot welded to the bottom of each hole using a capacitance discharge unit, thus forming a doubleintrinsic thermo-element junction. Such junctions provide a fast thermocouple response time and insure accurate positioning and contact at the bottom of the holes. [2] Following the controlled slow cooling runs, the specimens were cut for microstructural studies. Table 5.2 - As-rolled 5130 alloy steel test matrix. Test Tube Preheat (°Q Hood Doors Spacing (m) OD (m) Wall Thickness (m) Tube Length (m) 1 2 3 4 5 6 5130-1 5130-1 5130-2 5130-2 5130-3 5130-3 1145 1145 1145 1145 1145 1145 Yes Yes Yes No Yes No No No Yes No Yes No 0.419 0.203 0.203 0.203 0.419 0.419 0.08575 0.08575 0.08575 0.08575 0.08575 0.08575 0.0138 0.0138 0.01385 0.01385 0.01365 0.01365 0.914 0.914 0.914 0.914 0.914 0.914 Figure 5.6 depicts an actual cooling test run (excluding the use of the furnace hood). The presence of dark spots on the outer surface indicates the presence of scale formation. The darker regions at the ends of the tube indicate cooler sections of the tube associated with tube end effects. Figure 5.6 - Photograph of experimental test run # 6 (without furnace hood). In the series o f experimental runs, owing to the time consuming nature of tube instrumentation, each of the instrumented tubes was used twice in order to minimize costs and Chapter 5: Experimental Tests 87 variability within the geometry. Unfortunately, the reuse o f the tubes for the experimental runs gave rise to several unexpected results. First, since the thermocouple's wire diameter was quite small, not all o f the thermocouples survived subsequent tests. Second, each time the tubes were reused, the microstructure within the tubes changed. Finally, scale formation on the outer and inner surfaces o f the tubes from the previous tests may have modified some o f the heat-transfer properties. A larger design matrix was to be tested to overcome some o f these problems, however, time constraints at the Timken facilities prevented this. 5.4.1 Measured Cooling Curve Data Plots o f the resulting thermal histories for the various experimental runs (Tests 1 through 6 described in Table 5.2) are shown in Figures 5.7-5.16 and include the furnace hood thermocouple readings for the applicable tests. A discussion o f the basic trends and comparison with the model data w i l l be performed in the experimental analysis section (Section 6.2). i i i i I i i i i i i i i i 1 Ji iI iI iI Experimental Cooling Curves for an AISI 5130 Steel Alloy (Test 1 - Instrumented Furnace Hood) —i 400 i I I I I L Thermocouple 300 H -*-2 i_ 03 13 CO t_ <D 200 H CL 100 0 0 i I i I i i i i i i i i i i i i i i i i i i i—i—i—i—i—i—r 10 20 30 Time (min) Figure 5.7 - Furnace hood thermocouple readings - experimental Test 1. Chapter 5: Experimental Tests 1200 —i i i i iii i i i iiii i II Experimental Cooling Curves for an AISI 5130 Steel Alloy (Test 1 - Instrumented Tube) 30 Time (min) Figure 5.8 - Measured cooling curves - Test 1. 400 _i i i i i i i i_ J I I I L Experimental Cooling Curves for an AISI 5130 Steel Alloy (Test 2 - Instrumented Furnace Hood) ~i i i i i i i i | i i i i i—i—i—i—i—|—i—i—i—r 10 20 30 Time (min) Figure 5.9 - Furnace hood thermocouple readings - experimental Test 2. Chapter 5: Experimental Tests 1200 89 —I l I I l I I L_ Experimental Cooling Curves for an AISI 5130 Steel Alloy (Test 2 - Instrumented Tube) 1000 2 I 800 - Q_ E .CD 600 400 Time (min) Figure 5.10 - Measured cooling curves - Test 2. 600 j_i i i — i _1 I I I L_ Experimental Cooling Curves for an AISI 5130 Steel Alloy (Test 3 - Instrumented Furnace Hood) Thermocouple 200 n—i—i—i—r 10 20 30 Time (min) Figure 5.11 - Furnace hood thermocouple readings - experimental Test 3. Chapter 5: Experimental Tests 1000 J L_l I I L_J_ -J—I—I—III I I I I I I I I I Experimental Cooling Curves for an AISI 5130 Steel Alloy (Test 3 - Instrumented Tube) O o 800 i_ CD CD CL £ 600 400 n—i—i—i—r 0 10 n—i—i—i—i—r 20 Time (min) Figure 5.12 - Measured cooling curves - Test 3. Time (min) Figure 5.13 - Measured cooling curves - Test 4. 30 Chapter 5: Experimental Tests 91 500 _l I I I L_ _J I I I I I I I I I I I L Experimental Cooling Curves for an AISI 5130 Steel Alloy (Test 5 - Instrumented Furnace Hood) 0 10 20 30 Time (min) Figure 5.14 - Furnace hood thermocouple readings - experimental Test 5. 1200 I I I I L I I I I I I I I L_L_ Experimental Cooling Curves for an AISI 5130 Steel Alloy (Test 5 - Instrumented Tube) 1000 O o 1 800 ® CL E ,a> 600 400 Time (min) Figure 5.15 - Measured cooling curves - Test 5. Chapter 5: Experimental Tests Time (min) Figure 5.16 - Measured cooling curves - Test 6. Chapter 5: Experimental Tests 93 References 1. A S M Metals Handbook: Properties and Selection: Irons, Steels, and High-Performance Alloys. 10 ed., A S M International, Materials Park, OH, Vol. 1, 1990, pp. 179. th 2. G. Lockhart, Internal U B C Memo, July 2002. 94 Chapter 6: Model Validation CHAPTER 6 MODEL VALIDATION The A C M and CSC mathematical models were coded in F O R T R A N and developed to run on a personal computer (PC). An evolutionary approach was adopted for both model development and model validation. It started by building a simple model and then verifying it, with each subsequent version being built on the previously validated version. The thermal diffusion models, incorporating latent heat evolution, were first validated using a commercial finite-element package ABAQUS™ and second through experimentally derived temperature measurements (experimental tests of controlled slow cooling conditions) performed at the Timken Company, Canton, Ohio. Although, the modeling activities involved the development of two mathematical models: one to describe accelerated cooling and a second to describe controlled slow cooling conditions, the experimental tests of this study focused on controlled slow cooling conditions. Work is currently underway to validate (experimentally) the accelerated cooling model. The tests were conducted on industrial steel tube samples (instrumented at UBC) using an experimental facility designed and built at Timken. The Timken Company provided the microstructural routines (i.e., incorporating recrystallization kinetics, grain growth kinetics, and austenite decomposition kinetics for the formation of ferrite, pearlite, bainite and martensite) for the coupled 2-D transient heat diffusion models. The proprietary nature of this research places the validation of the phase transformation routines, as well as the microstructural experimental analysis, outside the scope of the current study. However, the coupling of the transformation kinetics routines and the associated evolved latent heat, as well as the incorporation of temperature and chemistry dependent thermal-physical properties, were performed at U B C . 6.1 Numerical Validation The numerical validation involved validating the thermal diffusion model as well as the algorithm associated with latent heat evolution. Owing to the lack of an analytical solution due to the complexity of the boundary conditions and the cylindrical nature of the geometry, the finite-difference models were compared with the results obtained using a benchmark commercial finite-element software package, ABAQUS™. The PC based F D M mathematical models were A B A Q U S is a trademark of Hibbitt, Karlsson and Sorensen Inc. Chapter 6: Model Validation 95 run on a Dell™ Dimension XPS T650r PHI machine with 256 Mb of R A M and a 650 M H z processor. The F E M mathematical models were run on a Silicon Graphics Inc. (SGI™) Origin 200 machine with 270MHz 64-bit MIPS RISC R12000™ microprocessors. The thermal-physical properties used in the analyses are presented in Table 6.1. For the purposes of thermal model validation the material properties were assumed to be homogenous throughout the calculation domain and equal to those for austenite. The characteristic enthalpies of transformation for the austenite-ferrite, austenite-pearlite, austenite-bainite, and austenitemartensite reactions are reported to be 75.2, 92.0, 92.0, and 83.6 kJ/kg, respectively (Nagasaka, et al.y* \ For validation purposes the larger enthalpy value of 92 kJ/kg, representing the pearlite l or bainite reactions, was adopted in the models to validate the algorithm for the evolution of latent heat. During this portion of the validation, the microstructural routines were turned off and a lumped sum approach to latent heat additions was used. In this approach, a linear incremental evolution of latent heat (92 kJ/kg x mass (kg) of simulated tube) in the 750°C to 700°C temperature range was applied. Table 6.1 - Steel properties used in the numerical validation of the A C M and CSC models. Property Symbol Values" Thermal conductivity Specific heat Density Emissivity Enthalpy of transformation k, (W/m/K) c , (J/kg/K) 11.363+1.5586"^ 467.1+0.1706T 7800 0.75 92 (700-750°C) p A (kg/m ) 3 £ H , (kJ/kg) f T is in degrees Celsius. 6.1.1 Accelerated Cooling Model The numerical validation of the A C M model involved the examination of several processing parameters; the geometry of the tube, the steel grade, the initial tube temperature, the boundary conditions, the spray boiling curve data (Section 4.2 - Figure 4.21) and the velocity at which the tube is traveling. Some cases are presented below as examples. Tables 6.2 and 6.3 summarize the parameters used for the comparisons between the A C M model predictions and the A B A Q U S F E M model. A sensitivity analysis was performed on the A C M model geometry and processing conditions in order to determine an appropriate mesh size to use in the numerical validation. For the analyses, a mesh of 11 nodes in the radial direction Dell is a trademark of the Dell Computer Corporation. SGI is a trademark of Silicon Graphics Inc. MIPS RISC R12000 is a trademark of MIPS Technologies Inc. Chapter 6: Model Validation 96 and 101 nodes in the axial direction (producing 1111 nodes and 1000 elemental control volumes) was used for both the F D M and F E M models. The temperature response of nodes located at the inner surface, mid-wall thickness and outer surface at the mid-axial position in the tubes were used as a basis for comparison between the two models. Figure 6.1 depicts the simulated boundary conditions considered in cases 1 & 2. Table 6.2 - Case 1: Processing parameters used in A C M validation (thermal diffusion model - no latent heat). Model Parameters Tube length (m) Value 0.2766 Outer diameter (m) 0.05222 Inner diameter (m) Type Convection Radiation Convection Type 0.03236 1417.59 298.15 293.15 0.02 400 Distance 0.025 m 0.2266 m 0.025 m Distance Adiabatic 0.1258 m Convection 0.025 m Adiabatic 0.1258 m Tinitial (K.) Tambient (K-) Twater Tube velocity (m/s) Time (s) OD boundary conditions ID boundary conditions Boiling Curve C C Boiling Curve C Table 6.3 - Case 2: Processing parameters used in A C M validation (thermal diffusion model - latent heat). Model Parameters Value Tube length (m) 0.2766 Outer diameter (m) 0.05222 Inner diameter (m) Tube velocity (m/s) Time (s) 0.03236 1417.59 298.15 293.15 0.02 400 Latent heat 92 kJ/kg evolved between (700-750°C) Tinitial (K) Tambient (K) Twater OD boundary conditions ID boundary conditions Type Convection Distance 0.025 m Radiation 0.2266 m Convection 0.025 m C Type Distance Boiling Curve Adiabatic 0.1258 m Convection 0.025 m Adiabatic 0.1258 m Boiling Curve C C Chapter 6: Model Validation 97 d = 0m i r I L = 0.2766 m J • < d = 0.1258 m , d = 0.1258 m d = 0.025 m Figure 6.1 - Schematic of accelerated cooling processing conditions (cases 1 & 2). The first case examines the thermal diffusion of heat in the absence of latent heat evolution under a variety of thermal processing conditions. As can be seen in Figure 6.2, the agreement between the results of the A C M model and the F E M model are relatively good. Both models predict the same behavior and there is only a slight variation of several degrees at the end of the cooling cycle. It has been postulated that the basis for the nodal temperature discrepancies between the F E M and F D M models may be the different ways in which the A B A Q U S F E M and the F D M mathematical models integrate in time. A B A Q U S uses a modified Newton-Raphson technique for solution of the problem in which so-called equilibrium iterations are conducted within a time-step. The F D M use a simple direct time-stepping approach. Chapter 6: Model Validation 98 1450 1400 ID-ACM MID - A C M 1350 g 1300 1250 3 l. ". mp 1150 H 1100 1050 ID-ABAQUS MID-ABAQUS I- • * 1 ".. !; OD - ABAQUS t * x * ti- .. Radiation Adiabatic I=I Convection 1 = 1 • i .• i i i i i i i * 1 = 1 i i r * • !; !: * * :! i " " * ' . i i i i i 1000 — 950 • • 1 * 1 1200 OD-ACM l' • 20 i • , i i i , i 40 Time (s) ID 60 80 10 Figure 6.2 - Case 1: Comparison of A C M model to A B A Q U S F E M model - no latent heat evolution. The second case, presented in Figure 6.3, compares the thermal predictions of the two models for identical processing conditions as described in the first case, but with the addition of a linear incremental evolution of latent heat in the 750°C to 700°C temperature range. Once again, the comparison was made on the basis of predictions for nodes located at the inner surface, mid-wall thickness and outer surface of the tube. The results for the second case also show good agreement between the A C M model and the A B A Q U S F E M model. Similar temperature differences are observed between the F E M and F D M models in this case, as in the first case, and are speculated to have arisen for the same reasons as discussed above. 99 Chapter 6: Model Validation 1450 1400 ID-ACM 1350 MID-ACM OD-ACM 1300 ID-ABAQUS 1250 im s E «-* H OD-ABAQUS 1200 •mm im 1150 a MID-ABAQUS Hi Start of latent heat evolution 1100 1050 1000 950 End of latent heat evolution 900 850 0 50 100 150 200 250 300 350 Time (s) Figure 63 - Case 2: Comparison of A C M model to A B A Q U S F E M model - latent heat evolution. 6.1.2 Controlled Slow Cool Model The CSC model was also validated against the commercial F E M model A B A Q U S using a number of problems selected to verify the various components of the model, including the basic heat-transfer model, the ability to account for the release of latent heat and the ability to handle a number of different processing configurations. The processing conditions examined include the geometry of the tube, steel grade, initial tube temperature, boundary conditions, adjacent tube spacing and the velocity at which the tube is moving. A selection of cases are presented below. Tables 6.4 and 6.5 summarize the processing parameters for the comparisons between the CSC model predictions and the A B A Q U S F E M model. As with the A C M model, a sensitivity analysis was performed on the C S C model geometry and processing conditions in order to determine an appropriate mesh size to use in the numerical validation. For the analyses, a mesh of 5 nodes in the radial direction and 36 nodes in the circumferential direction (producing 180 nodes and 180 elemental control volumes) was used for both the F D M and F E M models. The nodal comparisons of the models occurred at the inner surface, mid-radius and outer surface Chapter 6: Model Validation 100 nodes at the 0° position (see Section 4.1.2 - Figure 4.16). Please refer to Table 6.1 for material property data. The simulated cooling zones considered in cases 1 and 2 are depicted in Figure 6.4. Table 6.4 - Case 1: Processing parameters used in CSC validation (thermal diffusion model - no latent heat). Model Parameters Value Initial tube spacing(m) 0.1 Outer diameter (m) 0.05222 Inner diameter (m) BC Zone 1 Radiation 0.03236 1417.59 301.15 3600 Zone 2 Radiation Length (m) 8 8 T, p (K) 773.15 373.15 373.15 323.15 Tinitial (K.) Tambient (K.) Time (s) 0 Tbottom (K) T.(K) H T C (W/m /K) Velocity (m/s) 2 — — — — 0.005 0.01 Zone 3 Convection 1 — — 301.15 30 0.01 Table 6.5 - Case 2: Processing parameters used in CSC validation (thermal diffusion mod - with latent heat). Model Parameters Value Initial tube spacing(m) 0.1 Outer diameter (m) 0.05222 Inner diameter (m) 0.03236 1417.59 301.15 3600 92 kJ/kg evolved between (700-750°C) Zone 1 Zone 2 Zone 3 Radiation Radiation Convection Tinitial (K) Tambient (K) Time (s) Latent heat BC Length (m) 8 8 Ttop(K) 773.15 373.15 373.15 323.15 Tbottom (K.) Tcc(K) H T C (W/m /K) Velocity (m/s) 2 — — — — 0.005 0.01 1 — — 301.15 30 0.01 Chapter 6: Model Validation 0* 8m 1m (© Zone 1 Zone 2 Zone 2 BC: Radiation BC: Radiation BC: Convection Velocity = 0.005m/s Velocity = G.Olm/s Velocity = 0.01 m/s T T HTC = 30 W/nf/K top = 773.15 K Tbottom 373.15 K top = 373.15 K To, = 301.15 K Tbottom 323.15 K = = Figure 6.4 - Schematic of controlled slow cool processing conditions (cases 1 & 2). Figure 6.5 presents the first case examined for the CSC model. The latent heat associated with phase transformations is neglected in this case and the thermal diffusion model compares relatively well with the A B A Q U S F E M model. The cooling paths of all the nodes show similar behavior between the two models. Only a slight variation of several degrees is observed. 1500 600 -I 1 1 1 500 1 1 1 1 500 1000 1500 2000 0 " M T V > > ^ 1 — 2500 3000 Time (s) Figure 6.5 - Case 1: Comparison of C S C model to A B A Q U S F E M model - no latent heat evolution. Figure 6.6 demonstrates the C S C model's ability to account for latent heat evolution. The temperature profiles for the inner surface, mid-radius and outer surface nodes, produced Chapter 6: Model Validation 102 similar cooling paths. In addition, the evolution of latent heat seems to affect both models in a similar manner. 1500 - i 1 0 500 . . 1000 1500 1 1 2000 2500 . 3000 Time (s) Figure 6.6 - Case 2: Comparison of C S C model to A B A Q U S F E M model -latent heat evolution. Several additional cases were run (not presented) to further validate the model. Overall, both the A C M and the C S C models were found to exhibit favorable results, along with a robustness that enables them to handle a variety of different processing conditions. In addition, the PC run finite-difference models provide considerable improvements in terms of computational costs (i.e., computation times, storage memory, etc ...). In particular, runtimes on the PC were approximately 1/3 of those for the F E M model. This result is particularly promising if on-line control is being considered. Chapter 6: Model Validation 103 6.2 Experimental Analysis As described in Chapter 5 (Section 5.4), a series of six experimental tests were conducted on 5130 steel grade tubes to validate the CSC mathematical model. These tests were designed to simulate various controlled slow cooling conditions. Comparisons of the thermal fields were subsequently made between the simulated CSC model predictions and the experimental results. As with the numerical analysis, the C S C mathematical model was run on a Dell Dimension XPS T650r PHI machine with a 256 Mb of R A M and a 650 M H z processor. A simulation of 25 minutes of cooling (experimental time) required approximately 2.5 minutes of execution time. For the analyses, a mesh of 9 nodes in the radial direction and 17 nodes in the circumferential direction (producing 153 nodes and 153 elemental control volumes) was used. Tables 6.6 and 6.7 summarize the processing parameters that were used for the comparisons between the model predictions and the experimental data. The nodal locations for the thermocouple comparisons are depicted in Section 5.3 - Figure 5.5. Table 6.6 - As-rolled 5130 alloy steel test matrix. Test Tube 1 2 3 4 5 6 5130-1 5130-1 5130-2 5130-2 5130-3 5130-3 Preheat (°C) Hood 1145 1145 1145 1145 1145 1145 Yes Yes Yes No Yes No Doors No No Yes No Yes No Spacing (m) OD 0.419 0.203 0.203 0.203 0.419 0.419 0.08575 0.08575 0.08575 0.08575 0.08575 0.08575 (m) Wall Thickness (m) 0.0138 0.0138 0.01385 0.01385 0.01365 0.01365 Tube Length (m) 0.914 0.914 0.914 0.914 0.914 0.914 Table 6.7 - C S C model parameters. Test Initial Tube Temperature (°C) Austenite Grain Size Furnace Hood Temperature (°C) Floor Temperature (°C) Furnace Air Temperature (°C) Tube Velocity (m/s) Emissivity Thermocouple 1 2 3 4 5 6 1090 1120 883 1125 1133 1130 130 110 130 140 130 180 150 215 380 25 315 25 25 25 25 25 25 25 65 70 330 25 75 25 0.05 0.05 0.05 0.05 0.05 0.05 0.37 0.35 0.35 0.42 0.35 0.35 1 3 3 4 1 2 Indicates experimental thermocouple locations usedfor model validation (Section 5.3 - Figure 5 5). The material properties that have been used in the simulated model predictions are summarized in Tables 6.8 and 6.9. The density in these analyses was held constant at 7800 kg/m , to avoid altering the mass, since the volume of the computational domain does not vary 3 with temperature. Chapter 6: Model Validation 104 Table 6.8 - Material properties for 5130 alloy steel. Property Density (kg/m ) 3 Thermal conductivity (W/m/K) Symbol Value PA.F.P.B, M 7800 11.363+1.5586xl0" T 65.422-5.2176xl0" T+9.7673xl0" T 50.742-3.0567x 10" T+l. 1539x 10" T 2 k A ICF 2 kp.BM 6 2 2 7 674.73-0.2437 lT+2.0265x 10^T (T > 787°C): -10034.5+5.9668T+5.2xl0 T" (787°C < T< 769°C): 34754.5-31.9196T (769°C < T< 727°C): -11462.6+12.4346T (727°C < T< 527°C): -4704.5+4.568T+1.1057xl0 T" (527°C < T): 503.13-0.13068T-51.702xl0 T +4.47xlO^T 449.5+0.4501T (T < 720°C) 221656.4-864.4T+1.9795T -0.001478T (720°C < T< 780°C) -2.917xl0 +114590T-148.8T +0.06399T (780°C < T) 3277373.0-10575T+11.545T -0.00424T 2 CpA 9 Specific heat (J/kg/K) 2 CpF 2 9 2 5 CpP.B.M 2 2 Enthalpy of Formation (J/kg) H/F 7 3 2 3 2 2 , Hf 2 3 70651.0+225.23T-0.3469T +6.775x 10 T 92000 82600 The subscripts A, F, B, M refer to austenite, ferrite, bainite and martensite respectively. z _5 3 P Table 6.9 - Thermodynamic properties of dry air at atmospheric pressure. Temperature T (K) Density 273 293 313 333 353 373 473 573 673 773 1273 1.252 1.164 1.092 1.025 0.968 0.916 0.723 0.596 0.508 0.442 0.268 (kg/m ) 3 Specific Heat Absolute Viscosity fixlO (Ns/m ) Prandtl Number Pr (J/kg/K) Thermal Conductivity k (W/m/K) 1011 1012 1014 1017 1019 1022 1035 1047 1059 1076 1139 0.0237 0.0251 0.0265 0.0279 0.0293 0.0307 0.037 0.0429 0.0485 0.054 0.0762 17.456 18.24 19.123 19.907 20.790 21.673 25.693 39.322 32.754 35.794 48.445 0.71 0.71 0.71 0.71 0.71 0.71 0.71 0.71 0.72 0.72 0.74 c p 6 gpWxicr 8 (l/K/m ) 3 2 1.85 1.36 1.01 0.782 0.6 0.472 0.164 0.0709 0.035 0.0193 0.00236 The fit of the CSC model predictions to the experimental data was achieved largely through the adjustment of the emissivity value (~ 0.75) used in the model by means of an 'adjustment' parameter. This reduced emissivity corresponds to an 'effective emissivity' for the system rather than specifically for the 5130 steel. The adjustment parameter enables the 'tuning' or 'alignment' of the model to the particular experimental set-up, thereby providing a means of offsetting some of the simplifying assumptions postulated during model formulation. Chapter 6: Model Validation 105 For the purpose of validation, only two of the six tests are presented in the experimental analysis section, Tests 4 and 5 (plots for tests 1-6 may be found in Appendix C). These cases have been chosen to illustrate the effects of cooling with and without a furnace hood. In Test 4 the tubes cool without the cover of the furnace hood (i.e. cool in the ambient air). Figure 6.7 shows a comparison between the calculated cooling curve and the experimental data for Test 4. As can be seen, the agreement between the measurements and the predictions is good. However, the curves show a discrepancy between the predicted and measured data in the start of transformation kinetics and in the quantity of latent heat evolved. This would indicate that the microstructure subroutine needs further refinement to accurately predict the transformation kinetics. The predicted microstructure evolution of the various phases for Test 4 is shown in Figure 6.8. 1200 T 1 1 1 Time (s) Figure 6.7 - Predicted cooling curve and experimental data - 5130 Test 4 (without furnace hood). 106 Chapter 6: Model Validation 0 500 1000 1500 Time (s) Figure 6.8 - Microstructure evolution - 5130 test 4 (without furnace hood). Test 5 examines the controlled cooling of the tubes vvithin the furnace hood. The comparison between the calculated cooling curve and the experimental data for Test 5 is shown in Figure 6.9. Again, the agreement between the measurements and the predictions is good. A similar discrepancy in the start of the transformation kinetics between both curves is apparent. Figure 6.10 depicts the rnicrostructure evolution for Test 5. Chapter 6: Model Validation 107 1200 400 0 500 1000 1500 Time (s) Figure 6.9 - Predicted cooling curve and experimental data - 5130 Test 5 (with furnace hood). Chapter 6; Model Validation 0 108 500 1000 1500 Time (s) Figure 6.10 - Microstructure evolution - 5130 Test 5 (with furnace hood). For the six experimental test cases investigated, the model predictions agreed reasonably well with the measured cooling curves when effective emissivity values ranged in magnitude from 0.35 to 0.42 (corresponding to adjustment parameters of ~ 0.47-0.56 respectively). While these values may appear low, the feet that a relatively narrow range of emissivity values were capable of fitting a wide range of experimental conditions suggests that the model correctly captures the various factors effecting heat-transfer between the tube and its surroundings. The factors postulated to contribute to these low emissivity values are discussed in Section 7.1.3. The model will necessarily have to be 'retuned' if applied to a different system. Both the A C M and C S C models were numerically validated against the commercial F E M software package A B A Q U S for a variety of conditions. validated experimentally. Additionally, the C S C model was The results show that the models are capable of successfully predicting the thermal fields in tubes under industrial processing conditions. The microstructural models will need additional validation, which is outside the scope of this thesis project. Chapter 6: Model Validation 109 References 1. Y . Nagasaka, J.K. Brimacombe, E.B. Hawbolt, I.V. Samarasekera, B. Hernandez-Morales, and S. E . Chidiac, "Mathematical Model of Phase Transformations and Elasto-Plastic Stress in the Water Spray Quenching of Steel Bars", Metallurgical Transactions A , Vol. 24A, 1993, pp 795-808. Chapter 7: Discussion 110 CHAPTER 7 DISCUSSION The objective of this research was to develop a set of mathematical models for simulating "in-line" quenching and controlled slow cooling of tubular products for the Timken Company as part of their C T M P program. The long term goal is to optimize existing industrial processes and develop new processes leading to cost benefits and quality enhancements in tubular products. However, before these models can provide important insights into the tube manufacturing process, a sensitivity analysis of the model and process parameters must be performed. This analysis begins by determining model sensitivity to various solution parameters (i.e., initial time increment, mesh size, etc..) to ensure the models provide accurate and stable solutions that converge over a wide range of cooling conditions. The next step is to determine the effects of various model parameters that might have a significant influence on the cooling processes. The knowledge ascertained from these analyses can be subsequently used in the experimental validation process to fit the models with experimental data. Once the models have been calibrated to accurately and effectively simulate the cooling conditions in industrial processes, the optimization and design of processing conditions can proceed. An example of this application might be the use of the models, in conjunction with the knowledge of parameter influences on cooling, to simulate specific applications aimed at defining chemical compositions and cooling conditions necessary to develop different products (i.e., identify processing routes which provide the complex microstructures needed to give specialized mechanical properties). 7.1 Sensitivity Analysis The sensitivity of both the A C M and the CSC models was investigated. The sensitivity analyses were performed with three aims: 1) to determine the effects of certain parameters on the stability and convergence of the models, 2) to determine the effect of various parameters on the cooling of steel tubes and their influence on thermal fields and transformation product distributions, and 3) to determine the effect of the parameters on the temperature distribution in the tubes for use in fitting the model results to experimental data. The stability, accuracy and convergence of the transient heat-transfer solutions for the modeled cooling processes are essentially influenced by a number of solution parameters. Table 111 Chapter 7: Discussion 7.1 presents the solution parameters selected for the sensitivity analysis and postulated to have the greatest influence. Table 7.1 - Sensitivity analysis: solution parameters. Parameter Applicable Model Values Initial time increment, Ati i (sec) Maximum temperature difference, AT (K) Mesh size (nodes) A C M , CSC A C M , CSC A C M , CSC 0.01,0.1, 1, 10 5, 10, 20, 50 50, 90, 360, 720, 740 nitia max Two parameters associated with the dynamic time-stepping routine have been left out of the analysis: time increment and maximum volume fraction. The time increment parameter was neglected because it is not a user defined quantity as it is automatically calculated by the dynamic time-stepping routines. The maximum volume fraction parameter in the dynamic timestepping routine was neglected because the transformation kinetics routine was supplied by the Timken Company. Consequently, its analysis is outside the scope of the current study. In terms of processing, the cooling of seamless tubes is potentially influenced by a large number of parameters relating to the steel chemistry and the heat-transfer boundary conditions. The parameters chosen for this section of the sensitivity analysis were those postulated to have a significant influence on the cooling process, as well as those associated with the greatest uncertainty, or for which little data were available. These are summarized in Table 7.2. A number of parameters were obtained from literature; it was assumed that these parameters were sufficiently well established and that their inclusion in the sensitivity analysis was unnecessary. Table 7.2 - Sensitivity analysis: processing parameters. Parameter Wall thickness, wt (mm) Tube travel velocity, v (m/s) Latent heat Inter tube spacing, L (m) Quench location Tube end distance (mm) Furnace Temperature Emissivity, £ Applicable Model A C M , CSC A C M , CSC A C M , CSC CSC ACM ACM CSC A C M , CSC Values 5, 7.5, 10, 12.5 0.02, 0.05, 0.1,0.15 Trans, kinetics, No trans. Kinetics 0.124, 0.203, 0.419, 2.438 Int./Ext., Ext. 0, 12.5, 25, 37.5, 50, 62.5, 75 Hood, no hood 0.35, 0.45, 0.55, 0.65, 0.75 The results of the sensitivity analysis are shown in Figures 7.1 to 7.15. In the figures, the curves shown for each condition correspond to nodal locations that allow for the effects of the parameters on the stability and convergence criteria, and temperature fields to be readily assessed. Chapter 7: Discussion 112 7.1.1 Stability and Convergence In order to evaluate the stability and convergence of the A C M and C S C models, the effects of various parameters (initial time increment, maximum temperature difference, and mesh size) have been investigated. The issue of arriving at an acceptable solution is not a trivial one as the stability and convergence of the numerical scheme is dependent on specimen geometry, material properties and processing conditions as they can introduce non-linearities into the problem, making the effects of the aforementioned parameters process specific. Hence, both the A C M and CSC models have been examined using a variety of tube geometries and boundary conditions (i.e., radiation, convection, e t c . ) . Selected cases are presented in the forthcoming section. Initial Time Increment The sensitivity of the models to the initial time increment has been found to have a minimal effect on the stability and convergence of the solution owing to the adaptive nature of the dynamic time-stepping routine. Since the dynamic time-stepping routine automatically calculates the time increment for each time step based on temperature (AT ) and volume max fraction ( V f ) criterion, high values of the initial time increment automatically get reduced to the appropriate levels. Initial time increment values of 0.01 sec, 0.1 sec, 1 sec, and 10 sec were examined simulating identical quenching conditions over a period of 40 seconds. The boundary conditions for the simulated cooling regime are presented in Table 7.3. Table 7.3 - Sensitivity analysis: boundary conditions for Model Parameters OD boundary conditions ID boundary conditions AtjnMai comparisons. Type Length Value Convection 0.4 m /i=2000 W/m7K Natural Convection 1.0 m Convection Natural Convection Convection Natural Convection Adiabatic 0.3 1.0 1.0 2.0 0.3 m m m m m Convection 1.5 m Adiabatic 2.5 m — Boiling curve, Q=9.1 L/m /s 2 — Forced air, v =10 m/s ajr ... /?=2000 W/m /K 2 The predicted cooling curves from the A C M model for an outer surface node are presented in Figure 7.1. Chapter 7: Discussion 113 1300 860 858 856 1100 054 652 650~ 30 31 32 33 Time (s) 0) 900OD=76.2 mm ro i_ wt = 20 mm 0 CL E ,03 700 Tube length=1m f initial = ^ambient T 500- 300- 1273.15 K = Init. AT 0.01sec 0.1 sec - 1 sec 10 sec 298.15 K «tor = 293.15 K v = 0.1 m/s T 5 10 T 15 T 20 25 30 # Iter 2997 2978 3003 3032 35 Time (s) Figure 7.1 -Comparison of predicted cooling curves for an outer surface node on a tube experiencing a complex cooling regime. The results show that the initial time increment has a minimal effect on the stability and convergence of the solution. A maximum final temperature difference and percent difference in the total number of iterations of approximately 2°C and 1% respectively were realized between the initial time increment values. However, Figure 7.1 suggests that the accuracy of the solution is more sensitive to the initial time increment in regions where the overall convective heattransfer coefficient is high. A large temperature difference is observed between the solutions for Atmitaf=l0 seconds and the rest of Ati„uai values at solution times between t=14 to t=15 seconds. This is due to the fact that the non-linearity of the convective boundary condition is more severe (steep cooling rate between 14 and 15 seconds - spray water cooling). The radiation boundary condition has been found to be less sensitive to the variation in the initial time increment. Consequently, for highly non-linear transient heat-transfer problems, care must be taken when choosing the initial time increment. In addition, it should be noted that once a minimum initial Chapter 7: Discussion 114 time increment has been determined, the only significant effect of changing this value is in the increased computing times associated with a decreased initial time increment. Maximum Temperature Difference The maximum temperature difference parameter, AT , max in the dynamic time-stepping routine was found to have a significant effect on model stability, convergence and accuracy. AT max values of 5, 10, 20, and 50 K were investigated using the A C M model to simulate a tube experiencing an external quench (/j ,=2000 W/m /K) for a period of 40 seconds. The predicted 2 ex cooling curves for an outer surface node are presented in Figure 7.2. 0 5 10 15 20 25 30 35 40 Time (s) Figure 7.2 - Comparison of AT max values for an outer surface node experiencing an external quench. The results indicate a decrease in the solution accuracy with increasing ATmax where substantial solution differences (approximately 35°C) are observed in regions of high cooling rates. In addition, the ability of the model to form stable and converging solutions decreases as Chapter 7: Discussion the value of AT max for a AT max 115 increases. This is observed in Figure 7.2 where the simulated cooling curve value of 50 K is only stable for approximately 31 seconds, following which the solution begins to diverge and the program is unable to find a solution. Mesh Resolution The sensitivity of the models to mesh resolution was analyzed and found to have a significant effect. Mesh sizes of 50, 90, 360, 720, and 740 nodes were used in the CSC model to investigate mesh refinement on solution accuracy and execution time. The cooling of a tube spaced 0.419 m from the nearest adjacent tubes and experiencing natural convection on its outer surface for a period of 50 seconds was simulated. Figures 7.3 and 7.4 present the predicted cooling curves of an outer surface node. 1280' i 0 10 1 ' 20 1 30 ' 1 40 r 50 Time (s) Figure 7 3 - Simulated cooling curves of an outer surface node experiencing natural convection for several different mesh configurations. Chapter 7: Discussion 1158 OD = 85.75 mm wt = 6.9 mm T = 1273.15 K ^ a m Want 298.15 K v = 0,1 m/s taMM Rad. Nodes # Circ. Nodes # 10 5 = CD Exec. Time (s) 4 18 5 36 10 6 22 72 -*- 74 10 49 10 51 i _ 3 ro 1150 CD Q. E 1148 1146- 1144" 114245 46 47 48 49 50 Time (s) Figure 7.4 - Enlarged view of simulated cooling curves of an outer surface node experiencing natural convection for several different mesh configurations. A temperature difference of approximately 4°C was realized between the 50 and 740 node meshes. However, the execution time of the 740 node mesh was approximately 13 times longer than that of the 50 node mesh. In addition, a negligible difference in the final temperature of the 360 and 740 node meshes was observed, but the execution time of the former was approximately half of the latter. This indicates that there are inherent limits in achievable accuracy and in the optimum range of mesh refinement given a target accuracy. It should be mentioned that the optimum values are likely problem/process dependent and not distinct. Furthermore, the accuracy, stability and convergence were found to be sensitive to mesh geometries with excessive aspect ratios. Chapter 7: Discussion 117 Solution Parameters - Summary In general, it was observed from the models that accuracy may be improved if the mesh is refined for regions which experience large temperature gradient or if the time-step is reduced during periods when temperature changes quickly. Since the time-step is calculated within the models based on certain temperature difference and volume fraction criteria, solution stability and convergence are automatically realized: leading to significant computational time savings in solving real problems. The maximum temperature difference and mesh size were found to be the most important variables affecting stability and convergence. The initial time increment was found to have less of an influence owing to the adaptive nature of the dynamic time-stepping routine. The adaptive time-stepping routine is particularly important because the non-linearities (latent heat associated with phase transformations, boundary conditions, material properties, e t c . . ) occur at different times during processing for each problem. The computational time for each initial time increment selected was almost identical for every initial time increment owing to the adaptive nature of the dynamic time-stepping routine, however, it did rise markedly with mesh refinement. Thus, careful consideration must be given when refining the mesh for analysis of new tube geometries and process conditions. Likewise, the maximum temperature difference, AT , and the initial time increment, At max i, inita must also be chosen carefully to improve efficiency and eliminate inaccurate results. 7.1.2 Temperature Profiles The range of thermal processing routes covered by both models extends from controlled slow cooling through to intense "in-line" quenching. Under these cooling conditions, the effects of processing parameters on the thermal profiles vary. A n analysis has been carried out to investigate the influence of the main processing variables (Table 7.2) on the predicted thermal responses. The models enabled the quantification of the effects of the parameters on the uniformity of the thermal fields. A selection of simulations performed by the mathematical models are presented in the subsequent sections. Wall Thickness The wall thickness has a strong influence on the local cooling rate and microstructure of externally quenched tubes. Wall thickness values of 5 mm, 7.5 mm, 10 mm and 12.5 mm were examined by simulating identical quenching conditions over a period of 50 seconds using the Chapter 7: Discussion A C M model. 118 The heat-transfer coefficient (/z ,=2000 W/m /K) on the outer surface of the 2 ex simulated tube was applied over its entire length. A n adiabatic boundary condition was applied to the inner surface of the simulated tube. The predicted through-thickness cooling curves for the four wall thickness values are presented in Figure 7.5. 1300 3001 0 1 10 1 1 1 20 1 | 30 1 1 40 r—=l 50 Time (s) Figure 7.5 - Predicted through-thickness cooling curves for tubes of varying wall thickness experiencing an external quench. The simulations indicate a trend toward increased cooling rates as the tube thickness decreases, as would be expected for a given cooling condition. A difference of ar^oxirnately 300°C in the final through-thickness temperature was observed between the 5 mm and 12.5 mm wall thickness cases. The trend was such that the through-thickness thermal gradients (top and bottom curves of shaded thermal-field regions respectively) decrease as the wall thickness was decreased. Consequently, wall thickness plays a significant role in the development of the final tube microstructure. For example, varying the tube thickness for a given cooling intensity could result in a transition from a fully martensitic structure to partially martensitic structures within the tube. Furthermore, there is also the potential to change the thermal response at the surface to include a thermal rebound in the thicker tubes, thus allowing some capacity for auto tempering of the surface martensite to occur. Chapter 7: Discussion 119 Tube Travel Speed The model predictions confirm that the external surface temperature of the tubes after quenching is very sensitive to travel velocity. The A C M model was used to simulate the cooling of a 1 m long tube with a wall thickness of 20 mm and an outer diameter of 76.2 mm experiencing an external quench, 0.5 m in length, at different tube travel velocities. Tube travel velocities of 0.02 m/s, 0.05 m/s, 0.1 m/s, and 0.15 m/s were investigated over a period of 60 seconds. The through-thickness thermal profiles for each of the tube travel velocities are presented in Figures 7.6 to 7.9. 1300 Time (s) Figure 7.6 - Simulated through-thickness cooling curves for an externally quenched tube traveling at 0.02 m/s ( A C M model). Chapter 7: Discussion 120 1300 Figure 7.7 - Simulated through-thickness cooling curves for an externally quenched tube traveling at 0.05 m/s ( A C M model). 1300 OD - 76.2 mm Tube length -1 m h 2000 W/m /K v - 0.1 m/s 2 0 mm 20 mm Outer Diameter Inner Diameter Adiabatic • Convection 0.0 mm — 5.0 mm — - 10.0 mm Radiation CD — —- 15.0 mm 20.0 mm ID OD 0 10 20 30 30 40 50 50 60 Time (s) Figure 7.8 - Simulated through-thickness cooling curves for an externally quenched tube traveling at 0.1 m/s ( A C M model). Chapter 7: Discussion 121 1300OD = 76.2 mm Tube length = 1 m h „,» 2000 W/mVK Tsink - 293.15 K v - 0.15 m/s 800 700 m Radiation I=I Adiabatic • Convection — 0.0 mm — 5.0 mm — - 10.0 mm —- 15.0 mm - - 20.0 mm ID 600. |OD 10 10 0 20 30 30 40 40 50 50 60 Time (s) Figure 7.9 - Simulated through-thickness cooling curves for an externally quenched tube traveling at 0.15 m/s (ACM model). Generally, as the travel velocity increases, over a given external quench zone, the cooling rate decreases. Other simulations confirmed the importance of the travel velocity (i.e. quenching time) for ensuring a fully martensitic structure, even if high values of heat-transfer coefficient were employed. With increasing travel velocity a critical value can be achieved above which the temperatures of the internal part of the tube are higher than M„ and bainite can replace martensite. Internal/External Quenching The calculated thermal field of a 20 mm thick tube subjected to external quenching is compared to the predicted thermal field for the case of an internal and external quench in Figure 7.10. A n identical value of the heat-transfer coefficient at the outer and inner surface was adopted (/j=2000 W/m /K). 2 The calculations show that the introduction of the internal quenching increases the mean cooling rate at the inner surface from 12°C s" to 42°C s". The 1 cooling improvement is also evident at mid-thickness. 1 Chapter 7: Discussion 122 1300 3oo H 1 # 10 1 20 1 30 1 40 1 1 50 60 Time (s) 1 1 70 80 1 1 90 100 Figure 7.10 - Predicted cooling curves for internal/external quench. Consequently, higher cooling rates and reduced temperature gradients (low thermal stresses) can be achieved using a combination of internal/external quenching. Internal quenching may be achieved by the introduction of a lance incorporating a quenching head inside the tube or through tube immersion in a water tank with internal axial flow. It should be noted that the effects on microstructure should also be examined when considering the use of internal quenching units. A situation may arise where martensite is formed on both the O D and ID with a bainite core, thus the applicability of internal/external quenching units may ukimately be dictated by the evolved microstructure and associated mechanical properties. Latent Heat The evolution of the thermal field in steel during heat treating operations can be influenced significantly by latent evolution. As the steels pass through their transformation temperature range, the austenite decomposes into a mixture of different transformation products, which release heat. The amount of latent heat released varies depending on the type and quantity of the phases produced. Thus the evolution of this heat is linked to the phase transformation kinetics. A n illustration of latent heat effects on the thermal fields is presented in Figure 7.11. In Figure 7.11, the A C M model has been used to compare the results for a tube subject to an Chapter 7: Discussion 123 external quench with and without latent heat evolution included (with and without the phase transformation kinetics included). OD = 76.2 mm Time (s) Figure 7.11 - Comparison of predicted cooling curves of an external node submitted to an external quench with and without phase transformation kinetics. Figure 7.11 shows that the evolution of latent heat significantly alters the cooling of the steel at moderate cooling rates. In one study, Hodgson et a/. [1] modeled and measured cooling curves for low and high carbon steels under slow cooling conditions. They showed that the evolving microstructure significantly altered the cooling curve of the steel. In the low carbon steels investigated, relatively little heat of transformation was evolved and the distinct change in celling rate was due to the change in the thermal-physical properties. However, for the eutectoid steels, the change in cooling behavior was almost entirely due to the heat of transformation, with a smaller effect from the difference in thermal-physical properties between the various phases. At high quench rates the need to include latent heat effects is less certain as the effect is less obvious. The key point is to carefully establish the range of cooling rates and transformation products that are going to be investigated with the model. If equilibrium products and bainites are going to be examined, then the model must include latent heat effects. On the other hand, if Chapter 7: Discussion 124 the quench rates are to be severe and the transformation products are going to be exclusively martensitic then latent heat may be ignored. Tube Spacing Within the C S C model, inter tube spacing (center to center distance between adjacent tubes) was found to have an effect on the developing thermal profiles during slow cooling conditions. Inter tube spacing values of 0.124 m, 0.203 m, 0.419 m, and 2.438 m were investigated using tubes experiencing natural convection and radiation on the outer tube surface over a period of 1000 seconds. presented in Figure 7.12. The simulated cooling curves for an outer surface node are These cooling curves are the result of the heat-transfer between the body and its surroundings (i.e., the surrounding air, the floor, and its nearest neighboring tubes). 1300 100 200 300 400 i 500 r 600 700 800 900 1000 Time (s) Figure 7.12 - Predicted cooling curves of an outer surface node undergoing natural convection at various inter tube spacing values. The results show that as inter tube spacing increases, so does the mean cooling rate at the outer surface. This is because as the distance between adjacent tubes increases, the adjacent tubes provide less of a shielding effect to radiative heat-transfer to and from the surroundings. However, limits associated with tube spacing distances do exist for adjacent tube shielding Chapter 7: Discussion effects. 125 As the inter tube spacing distance approaches a distance of the outer diameter of the tube, minimum radiative heat losses occur. Conversely, as the inter tube spacing distance increases, a distance is reached (a function of outer diameter) where additional increases in inter tube distance cause a negligible change in the heat-transfer and the adjacent tube shielding effects approach zero. In essence, cooling is analogous to a single tube. Furnace Temperature During controlled slow cooling, the furnace temperature was found to alter the developing thermal fields within the tubes. In order to investigate the effects of furnace temperature, the CSC model was used to simulate two conditions: 1) slow cooling without a furnace hood (i.e., air cooling - tube experiencing natural convection and radiation at room temperature), and 2) slow cooling with a furnace hood (furnace hood temperature of 623.15 K, floor temperature of 298.15 K, and a furnace gas temperature of 473.15 K). 1300 0 1 1 1 1 1 1 1 r r 100 200 300 400 500 600 700 800 900 1000 Time (s) Figure 7.13 - Effects of furnace temperature on an outer surface node under slow cooling conditions. The results show that by increasing the furnace temperature consistent with the addition of a furnace hood, a significant decrease in the outer surface cooling rate may be achieved. Chapter 7: Discussion 126 Tube End Effects Tube end effects, or axial thermal gradients, largely occur due to the additional heat losses (convective or radiative) occurring at the tube ends resulting from the additional exposed surface area. The developing axial thermal gradients (at the tube ends) were investigated by examining the nodal temperatures at the mid wall thickness position along the axial length of the tube for a 20 mm thick, 1 m long tube, experiencing an external quench over a period of 50 seconds. 0 10 20 30 40 50 Time (s) Figure 7.14 - Comparison of tube end effects using the A C M model. The results show that the axial thermal gradients increase with distance and time. For the specific case depicted in Figure 7.14, after a period of 50 seconds, the axial thermal gradients are negligible at distances of approximately 37.5 mm from the tube ends. These effects are of some significance as they can alter the microstructure and mechanical properties in these regions. However, it is important to note that tube end effects are largely dependent on specimen geometry, material properties and boundary conditions. Hence, the importance of tube end effects is process specific. In the case of on-line process control, the minimization of tube end effects might be realized by introducing slight axial thermal gradients to offset the additional Chapter 7: Discussion 127 heat-transfer occurring at the tube ends. It should be mentioned that in 10-20 m long tubes these end effects represent a small portion of the overall tube length and are often simply removed in post processing. 7.1.3 Fitting Parameters Owing to the empirical nature of the relations used to quantify the heat-transfer coefficients, some adjustment or tuning of these correlations is often required when comparing the model predictions to experimental data. For the controlled slow cooling experimental runs conducted at the Timken Company, it was necessary to reduce the rate of heat transfer within the model to bring the model predictions in line with the measurements. This was achieved by multiplying the emissivity by an 'adjustment' parameter to lower effective emissivity. Figure 7.15 presents the effects of emissivity on the CSC model for a tube experiencing a combination of natural convection and radiation.. The cooling curves are of an outer surface node using five different emissivity values (0.35, 0.45, 0.55, 0.65, and 0.75). 0 300 600 900 1200 1500 Time (s) Figure 7.15 - Comparison of emissivity values for predicted cooling curves. The results show that lower cooling rates may be achieved as the effective emissivity value decreases. As stated in Section 6.2, the controlled slow cooling experimental tests were 128 Chapter 7: Discussion effectively 'fit' with the CSC model using a narrow range of effective emissivity values, 0.35 to 0.42 (Section 6.2 - Table 6.7). There are several different factors that might necessitate the need for utilizing a lower than expected effective emissivity. Some of these factors include: loose scale formation on the tubes; the assumptions with respect to the temperatures of the furnace hood, floor and furnace gas; the fact that the influence of the support rig on heat-transfer was neglected; and finally, the basic geometrical assumptions within the model used to approximate the radiation view factors. Particularly, in the case of the radiation view factor approximations, the angle of incident radiation was neglected and all incident radiation from a specific region (i.e. furnace top, etc..) was assumed to be equal. This exclusion of the incident radiation angle in the view factor approximations might have led to an overestimation in radiative heat-transfer within the C S C model and necessitated the need for a reduction in the emissivity values. Improvements with radiation view factor calculations could be made using view factor integration routines (incorporating ray tracing algorithms), however, these calculations are rather complex and would in all likelihood substantially increase the computation time of the models. Importantly, though, the relatively simple approach taken here does provide a reasonable representation of complex cooling systems. Scale Formation During the course of the experimental tests, it was noted that dark patches of scale were present on the instrumented tubes (Section 5.4 - Figure 5.6). It has been postulated that the scale layer on the steel surface affects the cooling conditions, with patches of coarse scale actually decreasing the cooling rate. This suggests that quite variable cooling curves might be obtained for identical set-up conditions depending on surface scale presence. This decrease in cooling rate is caused by the scale acting as a thermal insulator. If the steel is heavily worked and various descaling practices are employed to minimize the amount of scale, then obviously heat transfer will be improved and a different emissivity may be necessary. 7.2 Industrial Applications Once the models have been properly validated (i.e., comparing computed results with industrial results for various steel grades), the CSC and A C M models may be used: 1) to provide qualitative guidance in optimizing the design of a specific cooling process, 2) to determine the 129 Chapter 7: Discussion feasibility of processing various steel grades under actual industrial conditions, and 3) to provide the necessary data to formulate an on-line control solution (i.e., simulation data used to populate neural nets). As the preceding sensitivity analysis demonstrated, numerous parameters have a significant influence on the developing thermal fields during cooling processes. By utilizing the knowledge of how the process related parameters affect the cooling conditions experienced within a tube and the associated evolution of microstructure, conditions can be tailored to produce desired products. A n example of this might be in the development of a cooling process to produce a self-tempered steel tube. The production of self-tempered tubes requires the surface layers of the product to be cooled below the M , while leaving the core of the steel at a s sufficiently high temperature to temper the martensite by conduction of heat from the center to the surface, once the material has left the quenching zone. In the aforementioned example, tube travel velocity, water quench rate, steel chemistry, etc... might be examined using the A C M model to determine the optimum steel grade and processing conditions in order to reduce costs associated with conducting numerous experimental tests. Thus, the analysis of the cooling processes and the schemes proposed to model the process confirm the necessity of assessing the interactions among thermal and microstructural fields quantitatively when modeling the cooling of steel tubes. Successfully predicting the response enables the models to confidently optimize process parameters, at least in a comparative fashion and for their practical application in the following industrial applications: • Rationalization of quenching and controlled slow cooling conditions for steels currently in production; • Identification of new steel grades to achieve the desired structure and strength with the present plant constraints; • Prediction of improvements obtainable through plant modifications; • Design of new cooling systems (quench units, controlled slow cool furnaces, etc.); • Identification of simplified models, as well as the use of simulation results to populate neural nets for on-line process control. Chapter 7: Discussion 130 References 1. P.D. Hodgson, K . M . Browne, R.K. Gibbs, T.T. Pham, and D.C. Collinson, "The Mathematical Modeling of Temperature and Microstructure During Spray Cooling", Proceedings of the First International Conference on Quenching & Control of Distortion, Chicago, IL, 1992, pp 41-49. Chapter 8: Summary and Recommendations 131 CHAPTER 8 SUMMARY AND RECOMMENDATIONS 8.1 Summary and Conclusions Steel companies are using microstructural engineering models and the application of controlled thermo-mechanical processing (CTMP) to improve product quality, processing efficiency, and cost effectiveness. In the current study, two mathematical models were developed to predict the variation in temperature for seamless tubes during the water spray quenching and controlled slow cooling stages of the manufacturing process at the Timken Bearing and Steel Company of Canton, Ohio. These models identified as the accelerated cooling model (ACM) and the controlled slow cool model (CSC), have been coupled with a microstructure evolution model (incorporating recrystallization kinetics, grain growth kinetics and austenite decomposition kinetics for the ferrite, pearlite, bainite, and martensite reactions) provided by the Timken Company. The microstructural evolution model influences the thermal models through changes to the thermal-physical properties and the latent heat of transformation. A n evolutionary approach was employed for model development and model validation, with each new version being tested prior to further development. The models were coded in F O R T R A N to run on a personal computer and use finite-difference methods to solve the general transient heat conduction equation, with an alternating-difference-implicit scheme to solve the system of equations. Both the A C M and C S C models have been numerically validated against the commercial F E M software package A B A Q U S for a variety of conditions. Additionally, the CSC model has been validated through a series of experiments, conducted in the laboratory, to measure the thermal response of instrumented steel tubes under controlled slow cooling conditions. A summary of the key aspects of the A C M and CSC models is presented below. The A C M model employs a 2-D axisymmetric formulation (circumferential conduction is neglected), based on the premise that the tubes are subject to uniform cooling around their circumference and enter the cooling operation (i.e. water-spray stands) separately, with their axial lengths parallel to the direction of motion. The model has been designed as a flexible platform that allows for moving boundary conditions and multiple cooling zones, with on/off switches controlling the implementation of boundary conditions, applied to the outer diameter (OD) and/or the inner diameter (ID) of the tube. Chapter 8: Summary and Recommendations The CSC model employs 132 a 2-D circumferential geometry where radial and circumferential conduction have been included in the model and axial conduction has been ignored. This stems from the premise that the tubes enter the cooling process side-by-side (i.e., no vertical stacking) with the axial length of the tube perpendicular to the direction of motion. Conduction in the axial direction was neglected on the postulate that the tube is experiencing relatively uniform boundary conditions along the axial length, except at the tube ends. Since the axial length is relatively long in comparison to the tube diameter and wall thickness, end effects are assumed to be minimal and are, therefore, neglected. The model allows for multiple cooling zones that permit variable tube speeds, variable zone temperatures, and heat-transfer conditions. For the purpose of increasing robustness in the A C M and CSC models, a generalized heat balance approach (incorporating the latent heat evolved from any phase transformations that may occur in the process) was used. This approach allows for a variety of different boundary conditions to be applied to the control volumes where the boundary conditions describing the heat flow at the surface have been characterized by a Cauchy-type boundary condition. An overall heat-transfer coefficient, h , accounts for heat-transfer in cases where there is ov a combination of convection (natural, forced, and water) and radiation. Natural and forced-air convection have been characterized by a variety of empirical correlations, that are functions of Re, Pr, and Ra numbers, to describe the heat-transfer of air flowing over horizontal tubes. In the A C M model, water spray cooling (accelerated water cooling) has been modeled using boiling curve data (for a variety of different spray configurations) to characterize the heattransfer coefficients as a function of temperature. Radiation occurring during the process has been accounted for by "linearizing" the radiation equation so that the heat-transfer rate can be calculated using a conventional Cauchy-type boundary condition, with a "radiation-based" heattransfer coefficient. In the CSC model, the heat-transfer on the outer diameter surface due to radiation is more complex. Consequently, a geometric view factor approximation has been integrated into the model to account for the influence of adjacent tubes and the surrounding environment. A dynamic time-stepping subroutine has been added to the models to increase execution efficiency. The solution algorithm of this subroutine takes into account a maximum temperature 133 Chapter 8: Summary and Recommendations and volume fraction difference in order to determine the appropriate time increment needed for accurate solution integration and to ensure convergence. The incorporation of the material property data has been conducted using two different approaches: 1) a simple look-up table approach in which tabulated data, relating the thermalphysical properties to chemistry, structure and temperature, is entered into a model input data file; and 2) polynomial equations to characterize the thermal-physical property data. Linear interpolation techniques have been used to retrieve the appropriate data for the given temperature or phase. In those cases where there is a mixture of phases present, a simple rule of mixtures for any property P has been applied. The numerical validation comparing simulated results revealed that the thermal fields within the steel tubes can be predicted reasonably well using finite-difference modeling and that the latent heat effects of the phase transformations and the temperature dependence of material properties can be effectively incorporated in the model. The experimental and simulated results under controlled slow cooling conditions agreed well for the AISI 5130 steel alloy tubes. A n adjustment parameter was used to reduce the emissivity and yield an effective emissivity value that fit the CSC model results to the controlled slow cooling data obtained during the experimental runs performed at the Timken facilities. The relatively narrow range of the reduced emissivity values used to fit a wide range of experimental conditions suggests that the model correctly captures the various factors affecting hear-transfer between the tube and its surroundings. Loose scale formation on the tube surfaces, assumptions in surrounding temperatures, and the assumptions used to determine radiative view factors, were speculated to be the reason why it proved necessary to use a low effective emissivity in order to get the model to agree with the experimental results. A sensitivity analysis of the A C M and C S C models was performed. This involved investigating the influences of various solution and processing parameters and comparing the predicted temperature response and microstructural evolution with experimental temperature field measurements. The following sections summarize this sensitivity analysis. 8.1.1 Sensitivity Analysis - Solution Parameters The stability and convergence of the heat-transfer solutions for the modeled cooling processes are essentially influenced by the initial time increment, maximum temperature 134 Chapter 8: Summary and Recommendations difference and mesh size. The initial time increment was found to have minimal effect on the stability and convergence of the solution except in regions where the overall convective heattransfer coefficient is high. The maximum temperature difference had a significant effect on model stability, convergence and accuracy. temperature difference. Solution accuracy decreased with increasing The mesh size was found to have a significant effect on the solution accuracy especially for elements with excessive aspect ratio. Increasing refinement of the mesh for a given step criteria does not continuously improve accuracy, but may considerably raise the computing time. 8.1.2 Sensitivity Analysis - Process Parameters The range of thermal processing routes covered by both models extends from controlled slow cooling through to intense 'in-line' quenching. A n analysis was carried out to investigate the influence of the main processing variables on the predicted thermal responses and the results are presented below: • Wall thickness has a strong influence on the local cooling rate and microstructure of externally quenched tubes. The simulations indicate increased cooling rates on inner tube surfaces as tube thickness decreases. The influence of wall thickness on thermal profiles and resulting microstructure distributions plays a significant role in the development of conditions that result in fully martensitic structures or partially martensitic structures within the tube. • Generally, as tube travel velocity increases, over a given external quench zone, the cooling rate on the inner surface decreases. Controlling the travel velocity (i.e., quenching time) can influence percentage of evolved martensite, even if high values of the heat-transfer coefficient are employed. • The introduction of internal quenching may be used to increase the mean cooling rate at the inner surface and within the tube. The achievement of higher cooling rates and reduced temperature gradients within the tubes may be realized using a combination of internal/external quenching. • Transformation kinetics and associated latent heats of transformation were found to have a significant effect on the cooling rate experienced within the steel tubes. However, when the Chapter 8: Summary and Recommendations 135 latent heat corresponding to transformation is released and cooling rates are very high, the effect observed is relatively small. • Within the C S C model, inter tube spacing was found to have an effect on the developing thermal profiles during slow cooling conditions. The mean cooling rate at the outer surface was found to decrease as inter tube spacing increases. This is because as the distance between adjacent tubes increases, less of a shielding effect for radiative heat-transfer to and from the surroundings is provided. Lower and upper bounds were observed for this shielding effect. As the inter tube spacing distance approaches a distance approximately equal to the outer diameter of the tube, minimum radiative heat losses occurs. Conversely, as the inter tube spacing distance increases, a distance is reached where with additional increases in inter tube distance the effect becomes negligible and the adjacent tube shielding effects approach zero. • During controlled slow cooling, the furnace temperature was found to alter the developing thermal fields within the tubes. By increasing the furnace temperature, a significant decrease in the outer surface cooling rate may be achieved. • Tube end effects were found to be minimal in relation to overall tube length as well as being process specific. Where on-line process control is achievable, the minimization of tube end effects might be realized by introducing slight axial thermal gradients to offset the additional heat-transfer occurring at the tube ends. 8.2 Recommendations for Future Work The current work has shown that it is both necessary and possible to integrate thermal and microstructural models for the post deformation cooling of steel. However, there are still a number of areas that warrant further development before these models can be more widely applied. Suggestions for future work include: • The experimental validation of the A C M model. • The further development of material property data for industrially relevant steel grades to be included within the models. 136 Chapter 8: Summary and Recommendations • The further development of the microstructural evolution models to improve transformation kinetics routines and resulting microstructure distributions (involving the experimental analysis of tubes during cooling conditions - such data would provide a base for developing and validating heat flow versus microstructure correlations for seamless tube processing). • Boiling heat-transfer behavior of the steel tubes with spray cooling water is extremely complicated and a number of factors can influence this behavior. These include the effect of water velocity at the impingement point on the calculated boiling curve, as well as the cylindrical nature of the tubes. Thus, the development and further investigation of boiling water heat-transfer correlations that accurately describe the heat-transfer in tubes during water quenching is warranted. In addition, the quantification of the water spray quenching condition (i.e., water flow rate, nozzle type, and tube temperature) effects on boiling water heat-transfer for various steel alloys is needed. • The addition of different stacking configurations for the CSC model (i.e., tube stacking), as well as an investigation into improving the view factor approximations via view factor integration routines should be considered. Such an approach, while adding significant complexity, may provide more accurate predictions of the radiative heat-transfer in seamless tube cooling. 137 Appendix A: Heat-Transfer Finite-Difference Equations - ACM Model (Interior Node) APPENDIX A Heat-Transfer Finite-Difference Equations - A C M Model (Interior Node) The set of heat-transfer finite-difference equations for the seamless tube has been established, based on a heat balance at each node, by an energy analysis. The heat balance in both the radial and axial directions has been calculated by finite-difference approximation. An ADI method has been used to solve the parabolic equations (using tri-diagonal matrices and a time-step divided into two increments). In the first half time-step (f^ = f+ At/2) the implicit direction is in the radial direction. 112 In the second half time-step (t F+I =t p+1/2 + At/2) the implicit direction is in the axial direction. Figure A. 1 depicts the two implicit directions that have been used to determine the tri-diagonal equations. r 4- • (a) 1 direction (b) 2 st nd direction Figure A.1 - Schematic of ADI directions. The seamless tube has been spatially discretized, with a subsequent heat balance being applied to a control volume. The following derivations represent the numerical formulation of the heat-transferfinite-differenceequations for an interior node: Applying the 1 law of thermodynamics to the interior sub-volume depicted in Figure A.2 results st in: Aq + Aq + QAV = Aq r z r+Ar + Aq :+&z +E AV st (A.1) Where r and z are the spatial coordinates, Aq the heat energy acting on the control volume, Q the Appendix A: Heat-Transfer Finite-Difference Equations - ACM Model (Interior Node) 138 internal heat generation term due to phase transformations, AV is the volume of the control element, and E is the energy storage term. st control volume outer Diameter Az Az Figure A.2 - Schematic of A C M model interior node and associated control volume. 1 ADI Step: (f st =t +1/2 + At/2) - radial direction is implicit and axial direction is explicit: p (T* T -&Arr .AO Ar Ar ,1 / fpP+l/2 _ rp P+l/2 \ m+l,n m,n 1 1 AOAz 1 Ar --T T \ P m,n . + O Arr AOAz = m,n-l Az f jP+\2 ( j P _ j mji+l m,n - &Arr AO _ j P m,n + pc m.n t rn,n \ Arr_ .AOAz (A.2) Rearranging: 2 fcAOAzA/ f ^Arr AOAzt' •' n A r Y C ^ - C ^ V fcAOAzAr ( 2^ Ar J 2 ^ ^ ^ AOAz f"" , m>l + - kArr_ AOAr ^ , - 2 r f . + r r ^ 1 2pcArr • Arr_ .AOAzA/ A7 AOAz 2/x? Arr .AOAz * Substitute Fourier's numbers: Fo. = AAr p ^°J_(JP _2j p + + T l P i pc Az' 2r_ ArVp n,,n r + ^ 2 ) +1/2 m,n p ^ Fo ( r M 2r I m,n I, fc -C)(A.3) V2 pete' - T P + V p Fo. = Ar „ kAt - & Fo. = - and rearrange: 2pc + • 2 \ -T^ m.n D JP ArJT^ 1 2 2r_ m.n Ar^ ~2~,CIS + ^V A 2r. r p+i/2 (A.4) Appendix A: Heat-Transfer Finite-Difference Equations - ACM Model (Interior Node) 139 Finally yields: PC 2 2 P Fo. = (\ + Fo )T p+l/2 2 1 = f* ' +l ( Ar\ 2r_ ADI Step: (f° nd f Ar n,,n+ + r „ m 2. AOAz J (A.5) 2 ( rpP + \ _ *T./> + l ( rpP+\/2 _ rpP+\/2 \ V ••P+\/2 r + At/2) - axial direction is implicit and radial direction is explicit: 1 2 r / P + ] / 2 m-\,n m-\,n ' m,n •AArrAO Ar V A + 0 Arr Az r ( Ar^ AOAz r„ „ H v ^J V &Arr A O Ar + pc Az i (A.6) AOAz - ' jiP+\ _ jiP+\/2 ^ { '/ Arr A 2 AOAz J Rearranging: 2 f JiP+\/2 _ rpP+\/2 ^ Ar ArAOAzAr Pc Arr „AOAz p Ar mi AArr AOAf m,n+\ + 2pc • Arr AOAz m,n ' 1 Fo, P m,n ~ tj,p+\ 2T _ V m,n+l p + l i m,n T p + i • m,n-\ 2 P C A / \ —• A n A/ xi _ zpc 2r un V " 2) ' n 2r T P m , n ^ A kAt pe„Ar- 2 =(C-C ) / 2 (A-7) and rearrange: Ar jiP+V2 zr, Fo. m ,n rpP+\j2 _rpP+\/2 ^ m+\,n m,n z rpp+yi m,n ~ Ar m f kAt — - & Fo. z m,n Arr AOAzAr +0 J m,n-\ ) Ar" Z/x^Arr^AOAz Az Substitute Fourier's numbers: F o . = j +\ r AAOAzAr m—i.n 2r Ar n,r, r jP+\/2 m+\,n H (A.8) Finally yields: (i-^)C +2r ,/2 &r\p /2 +i 1 m-\,n , f + V Ar rpP+\/2 m+\,n At +Q 2 P , C (A.9) 140 Appendix B: Heat-Transfer Finite-Difference Equations - CSC Model (Interior Node) APPENDIX B Heat-Transfer Finite-Difference Equations - CSC Model (Interior Node) A similar set of heat-transfer finite-difference equations for the seamless tube has been established, based on a heat balance at each node, by an energy analysis for the CSC model. The heat balance in both the radial and circumferential directions has been calculated by finitedifference approximation. The time-step has been divided into two increments and an ADI method has been used to solve the parabolic equations. In the first half time-step (t p+1/2 = t + At/2) a tri-diagonal solver is used to solve the p implicit equations for conduction in the radial direction. For the second half time-step (t p+l t p+1/2 = + At/2) a full matrix solver is used to solve for conduction in the circumferential direction (Figure B.l). y (a) 1 direction -radial st (b) 2 nd direction - circumferential Figure B . l - Schematic of ADI directions. The seamless tube has been spatially discretized, with a subsequent heat balance being applied to a control volume. The following derivations represent the numerical formulation of the heat-transfer finite-difference equations for an interior node. 141 Appendix B: Heat-Transfer Finite-Difference Equations - CSC Model (Interior Node) Applying the 1 law of thermodynamics to the interior sub-volume depicted in Figure B.2 results st in: Aq +Aq r (t> + QAV = Aq r+Ar +A q ^ + E AV (B.l) sl where r is the radius, 0 is the angle between adjacent nodes, Aq the heat energy acting on the control volume, Q the internal heat generation term due to phase transformations, AV is the volume of the control element, and E, is the energy storage term. t F control volume Outer Diameter Inner Diameter Figure B.2 - Schematic of CSC model interior node and associated control volume. 1 ADI Step: (t st =f p+I/2 + At/2) - radial direction is implicit and circumferential direction is explicit: 2) -k\ r_.+ { Ar ) ( rpP+\/2 _ / T T P + 1 / 2 \ m+l.n m,n Ar AOAz 1 J 1 V 1 Ar m,n - kArAzl f -kArAzl m,n-\ r_.AO j P m,n+\ r J _ + QArr m A®Az = n J j> m.n _AO m,n At/2 J Arr_ .AOAz J (B.2) Rearranging: fcAOAzA/ ( 2pc ,A/r AOAzl ' •' r B ; (BIfl , _ 4 C - C 2 { Ar , kAOAzAt ( 2pc Arr A<&Az\ ' m n p mn ArHrj -T ^\ 2 | 2^ Ar P 142 Appendix B: Heat-Transfer Finite-Difference Equations - CSC Model (Interior Node) AArAzAr +• -7T +T ^ P JP f m,n+l P m,n • Arr AOAzAr m,n-\ + Q r AO / „ „, "bl „ \ { P+V2_ P = \ t T ( B 3 ) m,n Substitute Fourier's numbers: Fo^ - kAt /* (r ,„Ao) TP | **> (T Fo m,n^ Fo. + 2r • " m,n Ar\ ,n p P V m,n+\ ~ +T -OT P T ' m , n - l / _Fo P+l/2 P+ i 2 ~ m.n ~ 1 rp P+l/2 m-\,n 1 ,n Fo, -T V t Ar r 2r V \\* nAr 2pc A T V t 0 Ar r „ '" 2 r 2r~ r P+l/2 m Ar 2 >" „ m + — '" M 2r and rearrange: pc A r ' m p ± kAt & Fo, = 2 p P+l/2 m+l,n (B.4) Finally yields: A Ar 0 - ^ ) C + ^ ( C . + C-.)+e-, + Fo. Ar\ = (i+^)C -2r ,/2 2 ADI Step: (f nd m.n l J =f P+l Ar f P + y 2 - \ m-l,n P+I/2 rp P+l/2 m+l,n (B.5) V + At/2) - circumferential direction is implicit and radial direction is explicit: Ar -k AOAz f rp P+l/2 _ p P + l/2 \ ' m.n m-l.n ' Ar J AOAz f rpP+\ _ rpP+l \ m.n m,n-l ArArAzO ( pP+\ ( p P+l/2 _ J, P+\/2 \ m+l,n m.n Ar J _ rpP+l m.n+l - AArAz + 0Arr „AOAz = m r_ A O V \ m.n ( rpP+l _ rp P+l/2 \ m.n m.n r_ A O Arr At/2 AOAz (B.6) Rearranging: AAOAzA/ 2pcArr 2 nc Arr n I ( J, P+l/2 _ rp P+l/2 \ m—l,n m.n AOAz kArAzAt + • Ar m p m,n AOAz { ArAOAzAr +• 2pc Arr „AOAz Ar p ( rpP+l _ypP+l rpP + l \ m,n+l m.n m,n-l • Ar tn.n m A„,„A«>AzAf , rp P+l/2 _ rp P+l/2 ^ m+l.n m.n J{ Ar , r„ A O Substitute Fourier's numbers: Fo^ = £A/ & Fo. = kAt pc Ar' p and rearrange: J Appendix B: Heat-Transfer Finite-Difference Equations - CSC Model (Interior Node) JP+I _ Fo (PI T + P\ \ ' _-, P I T T + At + Fo ( P + l / 2 r Pp 2 Fo p i/2 r 2r_ + m,n Fo m,n \ 2r ArV/> i/2 f + r „ - 2r„ 2 C _Ar| m-l,n < Fo ( 2r 2 Ar r 143 P + l / 2 ) P+l/2 (B.8) Finally yields: (i-^)C +2r_ 1/2 r Ar — • T \l + P m 2 r mn + Ar rpP+\/2 m+l,n (B.9) 144 Appendix C: Experimental Validation Cooling Curves - CSC Model APPENDIX C Experimental Validation Cooling Curves - CSC Model The experimental validation cooling curves for the six experimental cases described in Chapters 5 and 6 are presented in Figures C l - C.6. Tables 6.6 - 6.9 (Chapter 6) summarize the processing parameters that have been used for the comparisons between the model predictions and the experimental data. The nodal locations for the thermocouple comparisons are depicted in Chapter 5, Section 5.3, Figure 5.5. 1200 ~i 1 400 - | 0 500 1 1 1 1000 , 1 1500 Time (s) Figure C l - Predicted cooling curve and experimental data - 5130 Test 1 (furnace hood). Appendix C: Experimental Validation Cooling Curves - CSC Model 145 1200 400 -I 0 1 1 500 1000 1 1500 Time (s) Figure C.2 - Predicted cooling curve and experimental data - 5130 Test 2 (furnace hood). Appendix C: Experimental Validation Cooling Curves - CSC Model 146 1200 1100 < Experimental — C S C Model 1000 O o 0) 3 2 o 900 800 Q. E .0) 700 600 500 400 0 500 1000 1500 Time (s) Figure C 3 - Predicted cooling curve and experimental data - 5130 Test 3 (furnace hood). Appendix C: Experimental Validation Cooling Curves - CSC Model 147 1200 0 500 1000 1500 Time (s) Figure C.4 - Predicted cooling curve and experimental data - 5130 Test 4 (without furnace hood). Appendix C: Experimental Validation Cooling Curves - CSC Model 148 1200 400 0 500 1000 1500 Time (s) Figure C.5 - Predicted cooling curve and experimental data - 5130 Test 5 (furnace hood). Appendix C: Experimental Validation Cooling Curves - CSC Model 149 1200 0 500 1000 1500 Time (s) Figure C.6 - Predicted cooling curve and experimental data - 5130 Test 6 (without furnace hood).
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Mathematical modeling of heat transfer during water...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Mathematical modeling of heat transfer during water spray cooling and controlled slow cooling of steel… Mwenifumbo, Steve 2003
pdf
Page Metadata
Item Metadata
Title | Mathematical modeling of heat transfer during water spray cooling and controlled slow cooling of steel tubes |
Creator |
Mwenifumbo, Steve |
Date Issued | 2003 |
Description | Mathematical models were developed which incorporate heat flow, phase transformation kinetics, and temperature and composition dependent thermal-physical properties to predict through-thickness temperature evolution during the manufacturing of steel tubular products under industrial cooling conditions. The two industrial cooling conditions investigated were an 'in-line' water quench, used to control austenite grain growth, and controlled slow cooling, applied at the end of the process to produce the desired microstructure and mechanical properties. This investigation combined computer modeling with laboratory measurements on industrial samples (i.e., AISI 5130 as-rolled steel tubes). A 2-D finite-difference axisymmetric accelerated cooling model (ACM) and a 2-D finite-difference circumferential controlled slow cooling model (CSC) were coded to run on a personal computer. The ACM model can handle multi-stage cooling with varying boundary conditions including radiation, natural convection, and forced convection. The model has the capacity to handle temperature and microstructural dependent thermal-physical properties associated with different phases that can form in the material. The CSC is capable of simulating the effects of tube spacing, forced air-cooling, and various ambient conditions. The basic heat diffusion parts of both the ACM and CSC models, including the routines associated with the liberation of latent heat, were validated against a benchmark commercial FEM software package, ABAQUS™. In addition, the results of the CSC model were verified experimentally by acquiring thermal data from a tube instrumented with thermocouples during the slow cooling of 5130 steel tubes. The CSC model predictions agreed reasonably well with the measured data, capturing the various factors affecting heat-transfer between tubes and their surroundings. A sensitivity analysis conducted with the CSC model indicated that the radiation exchange between different components (e.g., adjacent tubes, tubes and the fume hood, tubes and the floor) plays a significant role in heat-transfer. In addition, several factors influenced the effective emissivity within the model, including loose scale formation on the tubes, geometric assumptions in the model to approximate radiation view factors, and assumptions related to the temperature of the surrounding furnace environment. |
Extent | 25151926 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Date Available | 2009-11-04 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
IsShownAt | 10.14288/1.0078682 |
URI | http://hdl.handle.net/2429/14672 |
Degree |
Master of Applied Science - MASc |
Program |
Materials Engineering |
Affiliation |
Applied Science, Faculty of Materials Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 2003-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_2003-0706.pdf [ 23.99MB ]
- Metadata
- JSON: 831-1.0078682.json
- JSON-LD: 831-1.0078682-ld.json
- RDF/XML (Pretty): 831-1.0078682-rdf.xml
- RDF/JSON: 831-1.0078682-rdf.json
- Turtle: 831-1.0078682-turtle.txt
- N-Triples: 831-1.0078682-rdf-ntriples.txt
- Original Record: 831-1.0078682-source.json
- Full Text
- 831-1.0078682-fulltext.txt
- Citation
- 831-1.0078682.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0078682/manifest