A finite element calculation method for solidification heat transfer in continuous casting of profiled blanks

Abstract

The invention discloses a finite element calculation method for the solidification heat transfer of beam blank continuous casting, and relates to the technical field of steel continuous casting. The method comprises establishing a 1 / 2 beam blank two-dimensional geometric model according to the cross-sectional size of the beam blank; meshing the 1 / 2 beam blank two-dimensional geometric model to establish a 1 / 2 beam blank two-dimensional finite element model; establishing a 1 / 2 beam blank two-dimensional finite element solidification heat transfer model for the beam blank in view of the solidification heat transfer problem of the beam blank continuous casting; and calculating the temperature field change during the continuous casting process using the 1 / 2 beam blank two-dimensional finite element solidification heat transfer model. The invention fully considers the complex cross-sectional shape of the beam blank, and adopts triangular units and general quadrilateral units to better simulate the actual beam blank cross section; selects appropriate boundary conditions according to different cold zones, more accurately simulates the solidification heat transfer of the entire process of beam blank continuous casting production, and describes the temperature field change during the beam blank continuous casting process.

Description

Technical Field

[0001] The invention relates to the technical field of steel continuous casting, and in particular to a finite element calculation method for solidification heat transfer during continuous casting of profiled blanks. Background Art

[0002] The so-called special-shaped billet refers to the continuous casting billet with complex cross-sections other than square billets, slabs, round billets, and rectangular billets. Among them, the most widely used is H-shaped steel. H-shaped steel is an economical cross-section profile with the characteristics of thin wall, parallel inner and outer flanges, and right angles at the leg ends. It has excellent mechanical properties and performance. It is widely used in various fields of national economic construction, such as bridges, large ships, earthquake-resistant facilities and high-rise buildings. H-shaped steel is widely used. At present, there is a lot of room for the development of H-shaped steel in my country, especially large and medium specifications have a large gap, and the product quality is also a certain distance away from the advanced level of foreign countries. There are many reasons for the quality defects of special-shaped billets. Compared with conventional billets such as square billets and slabs with relatively simple cross-sections, special-shaped billets with more complex cross-sections will also produce some special quality defects during the production process. The quality defects of special-shaped billets can be summarized into two categories: surface defects and internal defects. Surface defects mainly include surface longitudinal cracks, corner longitudinal cracks, surface transverse cracks, star-shaped cracks, surface slag inclusions, pores and pits. Internal defects mainly include center cracks, straightening cracks, center segregation, center porosity and triangular area cracks.

[0003] Surface longitudinal cracks are manifested as cracks distributed longitudinally along the length of the billet, with a crack depth of less than 2mm. They may exist on both the web and the flange, but are more common on the web. They are one of the main quality defects of profiled billets. The main reason for the occurrence of surface longitudinal cracks is improper cooling water ratio, etc. The heat dissipation conditions at various points on the cross section of the profiled billet are quite different, the cooling is uneven, and the temperature difference at different positions on the surface is large, resulting in a higher incidence of surface cracks in profiled billets than in conventional billets. Surface transverse cracks are manifested as cracks distributed transversely along the length of the billet, with a crack depth of less than 2mm. Under the condition of high secondary cooling water ratio, the temperature of the billet entering the straightening zone is low, and the incidence of transverse cracks on the surface of the profiled billet is high. Straightening cracks are manifested as cracks distributed longitudinally along the thickness direction of the web at the tip of the cross section, and the crack depth is generally less than 2mm. When the hot billet is straightening, if it is in the low-temperature brittle zone, straightening cracks are prone to occur. Center segregation is manifested as a long strip shape distributed horizontally along the thickness direction at the center of the cross-section web, and the length is generally 3 to 4 mm; center porosity is manifested as fine pores distributed at the center line of the cross-section web. Center segregation and center porosity often exist together, which deteriorates the mechanical properties of steel, reduces the toughness and corrosion resistance of steel, and seriously affects product quality. When the secondary cooling zone adopts a strong cooling system, it is easy to cause center segregation and center porosity. Corner cracks are manifested as cracks distributed longitudinally along the thickness direction of the web at the cross-section flange, and the crack depth is generally greater than 2 mm. When the flange plate bulges, the cooling speed of the flange plate is too fast, and the cooling is asymmetric, cracks in the triangular area of ​​the profiled blank may be caused.

[0004] Therefore, uniform cooling of the billet in the crystallizer and the secondary cooling zone plays a vital role in reducing surface and internal defects of the billet. In order to achieve uniform cooling of the billet, it is necessary to first establish an online control system for the profiled billet, quickly and accurately predict the temperature change of the billet at each moment, and adjust the cooling water volume in real time. The key to the online control system for the profiled billet lies in the determination of the core algorithm for solidification heat transfer, that is, the establishment of a solidification heat transfer model for the profiled billet.

[0005] The boundary conditions of the continuous casting process are relatively complex. The current commercial software cannot calculate its solidification heat transfer process in a short time. The existing online solidification heat transfer models developed by individuals are mainly applied to conventional billets such as square billets and slabs, and there are not many for special-shaped billets. For example, the Chinese patent application with document number CN 103192048 B, "A method for simulating the solidification and cooling process of continuous casting billets based on accurate thermophysical property parameters", needs to rely on continuous casting cooling process simulation and optimization software, and is currently only for medium and thick slabs and cannot be extended to special-shaped billets; the Chinese patent application with document number CN 109446748 A, "A method for simulating the solidification process of continuous casting round billets", simulates the solidification process of continuous casting round billets by establishing a mathematical model of round billet solidification heat transfer, and is only applicable to continuous casting round billets. The numerical simulation method used is the finite difference method, which is cumbersome to handle irregular areas and has strict requirements on the continuity of the area; the Chinese patent application with document number CN 111199119 A, "A method for simulating the temperature of the head of continuous casting special-shaped billets", mainly targets the head of special-shaped billets, and does not consider the overall temperature change of the continuous casting billets. The numerical simulation method used is the finite volume method and divides the orthogonal grid units. The simulation of complex sections of special-shaped billets is poor, and the calculation accuracy is difficult to guarantee; the document number is CN110941889 A's Chinese patent application "Research Method for the Initiation and Extension of Micro and Macro Cracks in Continuously Casting Special-shaped Billets" simulates the temperature field and stress field of the special-shaped billets in the crystallizer through ANSYS to find out the positions where cracks in the billets are most likely to initiate and extend. However, the simulation of the temperature field is based on the ANSYS finite element software, and its principle is not explored, so it cannot be applied to the online control system of special-shaped billets.

[0006] In summary, establishing a solidification heat transfer model for profiled billets and calculating the temperature field at each moment of the continuous casting process are of great significance for real-time control of cooling water, achieving uniform cooling of the continuous casting billets, and improving the quality of the billets. However, the existing models do not yet have solidification heat transfer calculations for the entire process of profiled billet continuous casting. Summary of the invention

[0007] In view of the deficiencies in the above-mentioned prior art, the present invention provides a finite element calculation method for solidification heat transfer of profiled blanks during continuous casting, aiming to predict the temperature changes at each moment during the solidification process of profiled blanks, promote simulation and heat transfer analysis of profiled blanks, and provide a theoretical basis for establishing an online solidification heat transfer model of profiled blanks, optimizing water distribution in the secondary cooling zone, reducing surface defects of profiled blanks, and improving the quality of profiled blanks.

[0008] The technical solution of the present invention is:

[0009] A finite element calculation method for solidification heat transfer of beam blank continuous casting, the method comprising the following steps:

[0010] Step 1: According to the cross-sectional dimensions of the profiled blank, a 2D geometric model of the 1 / 2 profiled blank is established using finite element analysis software;

[0011] Step 2: Mesh the 1 / 2 profiled blank 2D geometric model, establish the 1 / 2 profiled blank 2D finite element model and save the mesh information;

[0012] Step 3: Aiming at the solidification heat transfer problem of beam blank continuous casting, a two-dimensional finite element solidification heat transfer model of 1 / 2 beam blank is established;

[0013] Step 4: Use the 1 / 2 profiled blank two-dimensional finite element solidification heat transfer model to calculate the temperature field changes during the continuous casting process.

[0014] Further, according to the finite element calculation method of the solidification heat transfer of the special-shaped billet continuous casting, the grid information includes unit information, node information, and boundary information; the unit information includes the unit number and the nodes contained in the unit; the node information includes the node number and the node coordinates; the boundary includes the upper boundary, i.e. the inner arc side, the lower boundary, i.e. the outer arc side, and the left boundary or the right boundary, i.e. the side arc side; the boundary information includes the unit information and node information on the boundary.

[0015] Further, according to the finite element calculation method for the heat transfer during the continuous casting of the beam blank, the step 3 comprises the following steps:

[0016] Step 3.1: Construct the basic equation of finite element heat transfer;

[0017] Step 3.2: Construct a finite element calculation model for transient heat transfer problems;

[0018] Step 3.3: Construct a planar three-node triangular element and perform heat transfer analysis based on the element for different heat transfer boundary conditions;

[0019] Step 3.4: Construct a planar four-node quadrilateral element and perform heat transfer analysis based on the element for different heat transfer boundary conditions.

[0020] Further, according to the finite element calculation method of the solidification heat transfer of the special-shaped billet continuous casting, the construction method of the plane three-node triangular unit is: encoding the three nodes as i, j, and m in a counterclockwise direction to construct a plane three-node triangular heat transfer unit; the construction method of the plane four-node quadrilateral unit is: selecting four vertices as interpolation nodes, and encoding the nodes as i, j, k and l in a counterclockwise direction to construct a typical plane four-node quadrilateral unit.

[0021] Further, according to the finite element calculation method for the solidification heat transfer of the profiled blank continuous casting, the content of the heat transfer analysis based on the plane four-node quadrilateral unit for different heat transfer boundary conditions is: selecting a square as a reference square unit in the parameter space; writing the geometric mapping function of the square to the quadrilateral in the physical space as a combination of node shape functions

[0022]

[0023] The square node shape function is

[0024]

[0025] The finite element method is used to transform the calculation method of the temperature stiffness matrix element of any quadrilateral unit in physical space to the global coordinate function into the calculation method of the square unit to the local coordinate in parameter space:

[0026]

[0027] The heat capacity matrix elements for any quadrilateral cell in physical space are also calculated in square cells in parameter space:

[0028]

[0029] Nodal loads due to heat sources in any quadrilateral element in physical space are also calculated in square elements in parameter space:

[0030]

[0031] In the above formula, x, y are Cartesian coordinate systems; ξ, η are coordinate systems in parameter space; x i 、x j 、x k 、x l are the horizontal coordinates of nodes i, j, k, and l respectively; y i ,y j ,y k ,y l are the ordinates of nodes i, j, k, and l respectively; k mn is the element of unit temperature stiffness matrix; k is the heat transfer coefficient; N m and Nn are all shape functions; a and b are the upper and lower limits of the integration in the parameter space, which are the coordinate variation ranges of the square respectively; J is the Jacobian matrix; |J| is the Jacobian matrix determinant; J * is the adjoint matrix of the Jacobian matrix; c mn represents the unit heat capacity matrix element; ρ is the material density; c T is the specific heat of the material; p m represents the node load; q v is the heat generated per unit volume inside the material, called the internal heat source density; l is the heat transfer boundary length of the unit.

[0032] Further, according to the finite element calculation method for the heat transfer during the continuous casting of the beam blank, the step 4 comprises the following steps:

[0033] Step 4.1: Read the steel grade information of the profiled billet, the structural parameters of the continuous casting machine and the continuous casting process parameters, and calculate the solidus temperature and the liquidus temperature;

[0034] Step 4.2: Import the mesh information of the 1 / 2 profiled blank 2D finite element model;

[0035] Step 4.3: Determine the time step; the time step is the time for calculating the movement of each slice; the slice represents the 1 / 2 profiled blank two-dimensional finite element model;

[0036] Step 4.4: Determine the unit phase region and calculate the physical property parameters of the unit;

[0037] Step 4.5: Calculate the unit heat transfer matrix, node load, unit heat capacity matrix and integrate the total heat transfer matrix, total temperature load, and total heat capacity matrix; the integrated total heat transfer matrix is ​​to decompose each unit heat transfer matrix into sub-matrices, and then "fill in" the block sub-matrices determined by each node sequence into the total heat transfer matrix according to the unit node number; the integrated total temperature load and total heat capacity matrix are the same as the integrated total heat transfer matrix;

[0038] Step 4.6: Solve the temperature field, that is, solve the finite element two-dimensional unsteady heat conduction control equation At the same time, determine whether the slice position exceeds the air cooling zone. If so, the calculation ends. If not, the slice continues to move and go to step 4.4.

[0039] Further, according to the finite element calculation method of the solidification heat transfer of the continuous casting of the special-shaped billet, the steel type information includes the steel grade and chemical composition; the continuous casting machine structural parameters include the crystallizer height and the secondary cooling zone structural parameters; the secondary cooling zone structural parameters include the number, length, inlet position and outlet position of the secondary cooling zone; the continuous casting process parameters include the casting temperature, pulling speed, crystallizer heat flux density, water flow and water temperature in each secondary cooling zone, ambient temperature and the cross-sectional dimensions of the special-shaped billet; the cross-sectional dimensions of the special-shaped billet include the wing tip thickness, wide face width, narrow face width, web length, web thickness, as well as the inner edge arc radius, flange radius and fillet radius.

[0040] Furthermore, according to the finite element calculation method of the heat transfer of the continuous casting solidification of the special-shaped billet, the phase region refers to the liquid phase region, the solid phase region and the solid-liquid two-phase region; the physical property parameters include unit density, unit specific heat, unit solid phase ratio and unit thermal conductivity coefficient.

[0041] Furthermore, according to the finite element calculation method of the solidification heat transfer of the special-shaped billet continuous casting, the method for judging the unit phase zone is: first, the center temperature of the plane four-node quadrilateral unit is calculated by the inverse distance weighted interpolation method, and the center temperature of the triangular unit is calculated by the two-dimensional linear interpolation method; finally, the phase zone of the unit is judged according to the unit center temperature.

[0042] Further, according to the finite element calculation method for the solidification heat transfer of the profiled blank continuous casting, the Gauss-Seidel iteration method is used to solve the matrix equation Calculate the temperature at each node.

[0043] In general, the above technical scheme conceived by the present invention has the following beneficial effects compared with the prior art: the finite element calculation method for the solidification heat transfer of the beam blank continuous casting provided by the present invention fully considers the complex cross-sectional shape of the beam blank, and adopts triangular units and general quadrilateral units to better simulate the actual cross-sectional shape of the beam blank; according to the different cold zones, appropriate boundary conditions are selected to more accurately simulate the solidification heat transfer of the entire process of beam blank continuous casting production. In short, the actual situation of beam blank continuous casting is fully considered, and the temperature field changes during the beam blank continuous casting process are described, thereby providing support and guidance for uniform cooling of the beam blank and improving the quality of the beam blank. BRIEF DESCRIPTION OF THE DRAWINGS

[0044] Figure 1 A schematic flow chart of a finite element calculation method for heat transfer during continuous casting of beam blanks provided in an embodiment of the present invention;

[0045] Figure 2 A schematic diagram of the cross-sectional dimensions of a profiled blank provided in an embodiment of the present invention;

[0046] Figure 3 A schematic diagram of a two-dimensional finite element model of a 1 / 2 profiled blank provided in an embodiment of the present invention;

[0047] Figure 4 Temperature distribution cloud diagrams of the profiled billets at different positions of the continuous casting machine provided by the embodiments of the present invention; wherein Figure (a) is a temperature distribution cloud diagram of the billet at the outlet of the crystallizer; Figure (b) is a temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 1; Figure (c) is a temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 2; Figure (d) is a temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 3; Figure (e) is a temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 4; Figure (f) is a temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 5;

[0048] Figure 5 A temperature variation curve diagram of each characteristic point of the profile blank provided by an embodiment of the present invention;

[0049] Figure 6 A comparison diagram of the calculated temperature and the measured temperature at the center of the flange surface of the profiled blank provided in an embodiment of the present invention. DETAILED DESCRIPTION

[0050] In order to facilitate the understanding of the present application, the present application will be described more fully below with reference to the relevant drawings. The preferred embodiments of the present application are given in the drawings. However, the present application can be implemented in many different forms and is not limited to the embodiments described herein. On the contrary, the purpose of providing these embodiments is to make the disclosure of the present application more thoroughly and comprehensively understood.

[0051] This embodiment takes the continuous casting process of 450mm×350mm×90mmQ235 steel profiled billets in a domestic steel plant as an example, and adopts the finite element calculation method of the profiled billet continuous casting solidification heat transfer of the present invention to calculate the temperature changes at various moments in the continuous casting process of the Q235 steel profiled billets.

[0052] The finite element calculation method for solidification heat transfer of beam blank continuous casting in this embodiment is as follows: Figure 1 As shown, the specific steps include:

[0053] Step 1: According to the cross-sectional dimensions of the profiled blank, a 2D geometric model of the 1 / 2 profiled blank is established using finite element analysis software;

[0054] According to the cross-sectional dimensions of the profiled blank, the coordinates of the key points are determined and imported into the ANSYS finite element analysis software, and then the key points are connected and the geometric surfaces are generated in ANSYS to establish a 1 / 2 profiled blank two-dimensional geometric model. The cross-sectional dimensions of the profiled blank include the wing tip thickness, wide face width, narrow face width, web length, web thickness, inner edge arc radius, flange radius, and fillet radius. The key points include the web endpoint, wing tip endpoint, inner edge arc endpoint, flange endpoint, and fillet endpoint.

[0055] In this embodiment, the width of the profiled blank is 0.435m, the width of the narrow face is 0.320m, the wing tip thickness is 0.045m, the web length is 0.229m, the web thickness is 0.090m, the flange radius is 0.020m, the inner edge arc radius is 0.100m, and the fillet radius is 0.010m. The schematic diagram of the cross-sectional dimensions of the profiled blank in this embodiment is shown in FIG. Figure 2 shown.

[0056] Step 2: Mesh the 1 / 2 profiled blank 2D geometric model, establish the 1 / 2 profiled blank 2D finite element model and save the mesh information;

[0057] Determine the unit type and unit shape in ANSYS, and determine the unit size, then mesh the 1 / 2 profiled blank 2D geometric model, and establish the 1 / 2 profiled blank 2D finite element model. Output the mesh information as a DAT format file and save it. However, due to the complex cross-sectional shape of the profiled blank, ANSYS cannot fully use quadrilateral units for meshing when performing regional segmentation, and some triangular units exist.

[0058] The mesh information includes unit information, node information, and boundary information. The unit information includes the unit number and which nodes are included in the unit; the node information includes the node number and the node coordinates; the boundary includes the upper boundary, i.e., the inner arc side, the lower boundary, i.e., the outer arc side, the left boundary or the right boundary, i.e., the side arc side. Since it is a 1 / 2 profiled blank, there is no right boundary or no left boundary; the boundary information includes the unit information and node information on the boundary.

[0059] In this embodiment, the unit type is selected as PLANE55, the unit shape is selected as a quadrilateral unit, the unit size is selected as 2mm, and the regional unit is automatically divided by the ANSYS finite element software. Figure 3 shown.

[0060] Step 3: Aiming at the solidification heat transfer problem of the profiled blank, a two-dimensional finite element solidification heat transfer model of 1 / 2 profiled blank is established;

[0061] Step 3.1: Derive the basic equations of finite element heat transfer;

[0062] The fundamental variable in a heat transfer process is temperature, which is a function of the geometric position in the object and time.

[0063] According to Fourier's heat transfer law and the law of conservation of energy, the governing equation for the heat conduction problem can be established, that is, the transient temperature field T(x, y, z, t) of the object should satisfy the following equation (1).

[0064]

[0065] Where ρ is the material density, kg / m 3 ;c T is the specific heat of the material, J / (kg·K); κ x , κ y , κ z The thermal conductivity coefficients along the x, y, and z directions are W / (m·K) respectively; Q(x, y, z, t) is the heat source intensity inside the object, W / kg.

[0066] There are three types of heat transfer boundary conditions, namely

[0067] The first type of boundary condition S1 is the Dirichlet condition, which gives a temperature value on the boundary

[0068]

[0069] The second type of boundary condition S2 is the Neumann condition for a given heat flux density.

[0070]

[0071] The third type of boundary condition S3 is the Neumann condition for given convective heat transfer

[0072]

[0073] Where n x ,n y ,n z is the direction cosine of the normal outside the boundary; is the temperature given on the heat transfer boundary that meets the first type of boundary condition S1, °C; is the given solidification heat flux density on the heat transfer boundary that meets the second type of boundary condition S2, W / m 2 ; is the convection heat transfer coefficient between the object and the surrounding medium, W / (m 2 ·K); T ∞ is the ambient temperature; t is the time, s; and the boundary of the object Ω is

[0074] If the initial condition IC of the problem is

[0075]

[0076] in, is the initial temperature, ℃.

[0077] The corresponding variational formulation is that in the permissible temperature field that satisfies the boundary conditions S1, S2, S3 and the initial condition IC, the real temperature field makes the following functional I take the minimum value, that is,

[0078]

[0079] In the process of dealing with practical problems, boundary conditions S2 and S3 are difficult to meet in advance. Therefore, these two conditions can be coupled into the functional, that is,

[0080]

[0081] Step 3.2: Establish a finite element calculation model for transient heat transfer problems;

[0082] In transient heat transfer problems, the temperature field of the element will change with time, that is,

[0083]

[0084] Where N(x,y,z) is the shape function, and the node temperature here is is time-varying, i.e.

[0085]

[0086] Substitute (8) into (7) and Finding the extreme value of variation, we can get

[0087]

[0088] in

[0089]

[0090]

[0091]

[0092]

[0093] Equation (10) is called the unit heat transfer equation; It is called the unit heat transfer matrix; is the unit node temperature array; is the unit heat capacity matrix; is the equivalent temperature load array of the unit node.

[0094] Step 3.3: Construct a planar three-node triangular element and perform heat transfer analysis based on the element for different heat transfer boundary conditions;

[0095] The three nodes are coded as i, j, and m in a counterclockwise direction to construct a planar three-node triangular heat transfer unit. The planar three-node triangular heat transfer unit composed of three nodes i, j, and m derives the unit matrices of the following three cases:

[0096] (1) No heat transfer boundary, that is, completely internal unit;

[0097] (2) If a certain edge of the unit, such as the jm edge, is a heat transfer boundary that meets the second type of heat transfer boundary condition S2: Constant to describe;

[0098] (3) If a certain edge of the unit, such as the jm edge, is a heat transfer boundary that meets the third type of heat transfer boundary condition S3: Constant to describe.

[0099] The nodal temperature array for this element is

[0100]

[0101] The interpolation relationship of the unit temperature field is taken as

[0102]

[0103] The shape function matrix N is

[0104] N=[N i N j N m ] (17)

[0105] in

[0106]

[0107]

[0108]

[0109] (i,j,m loop)

[0110]

[0111] Where A is the unit area, m 2 .

[0112] The following are several different heat transfer boundaries and the corresponding heat transfer matrices are calculated respectively. and nodal equivalent temperature load array

[0113] (1) Completely internal unit (no heat transfer boundary)

[0114] Substitute the shape function expression (17) into and In the calculation formulas (11) and (12) (note that only two-dimensional problems are considered here), we have

[0115]

[0116]

[0117] (2) When a certain edge of the unit, such as the jm edge, is a heat transfer boundary that meets the second type of heat transfer boundary condition S2

[0118] Substitute the shape function expression (17) and BC(S2) expression (3) into and In the calculation formulas (11) and (12) (note that only two-dimensional problems are considered here), we have

[0119]

[0120]

[0121] Where L is the length of the jm edge, m.

[0122] (2) When a certain edge of the unit, such as the jm edge, is a heat transfer boundary that meets the third type of heat transfer boundary condition S3

[0123] Similarly, substitute the shape function expression (17) and BC(S2) expression (4) into and In the calculation formulas (11) and (12), we have

[0124]

[0125]

[0126]

[0127] Step 3.4: Construct a planar four-node quadrilateral element and perform heat transfer analysis based on the element for different heat transfer boundary conditions.

[0128] Select four vertices as interpolation nodes, and encode the nodes as i, j, k and l in a counterclockwise direction to construct a typical planar four-node quadrilateral unit. The temperature field inside the unit is calculated using a bilinear interpolation function:

[0129] T(x,y)=a1+a2x+a3y+a4xy (25)

[0130] In the formula, a1, a2, a3, and a4 are unknown coefficients.

[0131] The calculation of the coefficients of the bilinear interpolation function in any plane four-node quadrilateral element is very difficult, and some integral terms in the heat transfer matrix, heat capacity matrix and load vector expressions cannot be obtained by analytical primitive functions. Therefore, the finite element method does not directly calculate the quadrilateral element model of arbitrary shape, but uses the isoparametric mapping method to indirectly calculate the quadrilateral element model.

[0132] The present invention selects the square [-1,1]×[-1,1] as the reference square unit in the parameter space. That is, the coordinates of the vertices of the square in the parameter space are i(-1,-1), j(1,-1), k(1,1) and l(-1,1). The geometric mapping function of the square [-1,1]×[-1,1] in the parameter space to the quadrilateral in the physical space is written as a combination of node shape functions, that is,

[0133]

[0134] The square node shape function is

[0135]

[0136] The finite element method is used to transform the calculation method of the temperature stiffness matrix elements of arbitrary quadrilateral units in physical space to the global coordinate function into the calculation method of the square unit to the local coordinate in the parameter space, that is,

[0137]

[0138] The upper and lower limits of the integration in the parameter space are the coordinate variation range of the square. For the [-1,1]×[-1,1] square, a=-1, b=1.

[0139] The heat capacity matrix elements of any quadrilateral unit in physical space can also be calculated in square units in parameter space, that is,

[0140]

[0141] The nodal loads generated by the heat source in any quadrilateral element in the physical space are also calculated in the square element in the parameter space, that is,

[0142]

[0143] Among them, q v The heat generated per unit volume inside the material is called the internal heat source density, W / m 2 ; J is called the Jacobian matrix

[0144]

[0145] |J| is called the Jacobian matrix determinant

[0146]

[0147] J * The adjoint matrix of the Jacobian matrix is ​​called

[0148]

[0149] Regarding the treatment of heat transfer boundaries, a four-node quadrilateral unit consisting of four nodes i, j, k, and l is used to perform heat transfer analysis for different heat transfer boundary conditions, and the unit matrices for the following three situations are derived:

[0150] (1) No heat transfer boundary, that is, completely internal unit;

[0151] (2) If a certain edge of the unit, for example, the ij edge, is a heat transfer boundary that meets the second type of heat transfer boundary condition S2: Constant to describe;

[0152] (3) If a certain edge of the unit, for example, the ij edge, is a heat transfer boundary that meets the third type of heat transfer boundary condition S3: Constant to describe.

[0153] (1) No heat transfer boundary, that is, completely internal unit;

[0154] When there is no heat transfer boundary and the unit is completely internal, the calculation formula for the unit temperature stiffness matrix element is formula (28), the calculation formula for the unit heat capacity matrix element is formula (29), and the calculation formula for the node load is formula (30).

[0155] (2) If a certain edge of the unit, such as the ij edge, is a boundary that meets the second type of heat transfer boundary condition S2: it is described by a constant;

[0156] The calculation formula of the unit temperature stiffness matrix element is formula (28), the calculation formula of the unit heat capacity matrix element is formula (29), and the calculation formula of the node load is as follows:

[0157]

[0158] (3) If one edge of the unit, for example, the ij edge, meets the third type of heat transfer boundary condition S3: Constant to describe.

[0159] The calculation formula of the unit heat capacity matrix element is formula (19), and the calculation formula of the unit temperature stiffness matrix element and the node load is as follows:

[0160]

[0161]

[0162]

[0163] In the above formula, k mn and k mm are all element temperature stiffness matrix elements; k is the heat transfer coefficient; N m and N nare all shape functions; a is the upper limit of the integration in the parameter space; b is the lower limit of the integration in the parameter space; J is the Jacobian matrix; |J| is the determinant of the Jacobian matrix; J * is the adjoint matrix of the Jacobian matrix; x, y are Cartesian coordinate systems; ξ, η are coordinate systems in the parameter space; c mn represents the unit heat capacity matrix element; ρ is the material density; c T is the specific heat of the material; p m represents the node load; q v is the heat generated in the unit volume of the material inside the unit volume, which is called the internal heat source density; L is the heat transfer boundary length of the unit; is the given solidification heat flux density on the boundary that meets the boundary condition S2; is the convective heat transfer coefficient between the object and the surrounding medium; T ∞ is the ambient temperature; t is the time.

[0164] Step 4: Use the 1 / 2 profiled blank two-dimensional finite element solidification heat transfer model to calculate the temperature field changes during the continuous casting process.

[0165] Step 4.1: Read the steel grade information of the profiled billet, the structural parameters of the continuous casting machine and the continuous casting process parameters, and calculate the solidus temperature and the liquidus temperature;

[0166] The steel grade information includes steel grade and chemical composition; the continuous casting machine structural parameters include mold height and secondary cooling zone structural parameters; the secondary cooling zone structural parameters include the number, length, inlet position and outlet position of the secondary cooling zone; the continuous casting process parameters include casting temperature, casting speed, heat flux density of the crystallizer, water flow and temperature of each secondary cooling zone, ambient temperature and cross-sectional dimensions of the profiled blank. The cross-sectional dimensions of the profiled blank include wing tip thickness, wide face width, narrow face width, web length, web thickness, inner edge arc radius, flange radius and fillet radius.

[0167] In this embodiment, the steel grade is Q235, and the chemical composition of the steel is shown in Table 1.

[0168] Table 1 Chemical composition of 450mm×350mm×90mm Q235 profiled blanks from a domestic steel plant

[0169] chemical composition C Si Mn P S Mass fraction (%) 0.19 0.21 0.47 0.025 0.027

[0170] In this embodiment, the crystallizer height of the profiled blank continuous casting machine is 0.7m; there are 5 secondary cooling zones, of which the length of secondary cooling zone 1 is 0.647m, and the entrance and exit are 0.7m and 1.347m away from the meniscus respectively; the length of zone 2 is 1.486m, and the entrance and exit are 1.347m and 2.833m away from the meniscus respectively; the length of zone 3 is 2.337m, and the entrance and exit are 2.833m and 5.170m away from the meniscus respectively; the length of zone 4 is 2.378m, and the entrance and exit are 5.170m and 7.548m away from the meniscus respectively; the length of zone 5 is 2.224m, and the entrance and exit are 7.548m and 9.772m away from the meniscus respectively. The casting temperature of Q235 profiled blank is 1545℃, and the pulling speed is 0.9m / min; the heat flux density of the wide surface of the crystallizer is 247w / m 2 The heat flux density of the narrow surface is 243.677w / m 2 The cooling conditions of the secondary cooling zone are shown in Table 2. The secondary cooling water temperature is 25°C, and the ambient temperature of the air cooling zone is 90°C.

[0171] Table 2 Secondary cooling conditions of 450mm×350mm×90mm Q235 profiled billet in a domestic steel plant

[0172] Second cooling zone 1 2 3 4 5 Inner arc water volume (L / min) 36 33 21 13 13 Outer arc water volume (L / min) 36 33 21 13 13 Side arc water volume (L / min) 40 37 19 8.4 8.4 Total water volume (L / min) 112 103 60 34..4 34.4

[0173] In this embodiment, the calculation formula for the corresponding solid-liquid phase temperature is selected according to the C content of Q235 steel:

[0174]

[0175]

[0176] Where T l is the liquidus temperature; T s is the solidus temperature; [C], [Si], [Mn], [P], [S], [Cr], [Al], [Ni], [V], and [Mo] are the mass fractions of the corresponding elements.

[0177] Step 4.2: Import the mesh information of the 1 / 2 profiled blank 2D finite element model;

[0178] Import the DAT format file exported in step 2 into the calculation program and store it. The grid information includes unit information, node information, and boundary information. The unit information includes the unit number and which nodes are included in the unit; the node information includes the node number and node coordinates; the boundary information includes which units and nodes are on the boundary.

[0179] Step 4.3: Determine the time step;

[0180] The time step is the time taken to calculate the movement of each slice (called a slice because the two-dimensional model has no thickness). The time step multiplied by the pulling speed is the distance the slice moves each time, and the accumulated moving distance is the position of the slice.

[0181] In this embodiment, the time step is 0.1 s and the pulling speed is 0.9 m / min.

[0182] Step 4.4: Determine the unit phase region and calculate the physical property parameters of the unit;

[0183] The center temperature of the quadrilateral element is calculated by the inverse distance weighted interpolation method, and the center temperature of the triangular element is calculated by two-dimensional linear interpolation. The phase region where the element is located is determined according to the center temperature of the element and the physical property parameters of the element are calculated.

[0184] The phase region refers to the liquid phase region, the solid phase region and the solid-liquid two-phase region;

[0185] The physical property parameters include unit density, unit specific heat, unit solid phase ratio and unit thermal conductivity.

[0186] In this embodiment, it can be known from step 4.2 that the solidus temperature is 1485° C. and the liquidus temperature is 1515° C. The unit phase region is determined according to the unit temperature and the physical property parameters of the unit are calculated using a suitable formula.

[0187] Step 4.5: Calculate the unit heat transfer matrix, node load, unit heat capacity matrix and integrate the total heat transfer matrix, total temperature load, and total heat capacity matrix;

[0188] According to the formulas in steps 3.2 and 3.3, the unit heat transfer matrix, node load, and unit heat capacity matrix of each unit are calculated respectively, and boundary conditions are added to the units on the boundary. The inner arc side, outer arc side, and side arc side are calculated according to the different positions of the slice in the continuous casting machine and the cooling conditions of the crystallizer and the secondary cooling zone read in step 4.1. It should be noted that the spray water is not evenly distributed in actual production, so the water volume of the boundary unit at different positions on the boundary needs to be considered when calculating the boundary conditions.

[0189] The integrated total heat transfer matrix is ​​to decompose each unit matrix into sub-matrices, and then fill the block sub-matrices determined by each node sequence into the total heat transfer matrix according to the unit node number. The same is true for the integrated total temperature load and total heat capacity matrices.

[0190] Step 4.6: Solve the temperature field, that is, solve the finite element two-dimensional unsteady heat conduction control equation At the same time, determine whether the slice position exceeds the air cooling zone. If yes, the calculation ends. If not, the slice continues to move and go to step 4.4;

[0191] Solve the matrix equation using Gauss-Seidel iteration method Calculate the temperature of each node. After the solution is completed, determine whether the slice position exceeds the air cooling zone. If it has exceeded the air cooling zone, the calculation ends and the solidification heat transfer calculation of the special-shaped billet continuous casting is completed. If it has not exceeded the air cooling zone, the slice moves to the next position and goes to step 4.4.

[0192] This embodiment uses the Gauss-Seidel iteration method to solve the matrix equation. The error is the square root of the sum of the squares of the differences between two iterations of all nodes. The error tolerance is 0.001. In order to improve the calculation efficiency, the parallel calculation function parallel_for in OpenCV is used for CPU parallel calculation.

[0193] Step 5: Visualize and post-process the calculation results of the solidification heat transfer of the beam blank continuous casting;

[0194] After the solidification heat transfer calculation of the beam blank continuous casting is completed, the calculation results at different positions of the continuous casting machine are extracted and imported into the Tecplot software for visualization to form a temperature cloud map; the temperature changes of the web center, inner edge center, narrow surface center and fillet center are extracted respectively, and imported into the Origin software for visualization to form a temperature change curve.

[0195] This embodiment provides a temperature distribution cloud diagram of a 450mm×350mm×90mm Q235 profiled blank at different positions of a continuous casting machine in a domestic steel plant, such as Figure 4 As shown, Figure (a) is the temperature distribution cloud diagram of the billet at the outlet of the crystallizer; Figure (b) is the temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 1; Figure (c) is the temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 2; Figure (d) is the temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 3; Figure (e) is the temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 4; Figure (f) is the temperature distribution cloud diagram of the billet at the outlet of the secondary cooling zone 5. The corresponding temperature change curves of each characteristic point are as follows, Figure 5 The calculated temperature at the center of the flange surface is compared with the measured temperature. Figure 6 As shown. Figure 6 It can be seen that the simulation results are basically consistent with the actual temperature. Due to the shape characteristics and cooling characteristics of the shaped blank, as the cooling progresses, the high temperature area in the center of the shaped blank gradually decreases and gradually moves to the inner edge. Since the inner edge is concave, the heat transfer effect is poor there, and the surface temperature of the inner edge is higher than the temperature at other locations on the surface. The corners of the shaped blank are in a two-dimensional heat transfer mode, and the heat dissipates faster, so the corner temperature is the lowest.

[0196] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or make equivalent replacements for some or all of the technical features therein. However, these modifications or replacements do not cause the essence of the corresponding technical solutions to deviate from the scope defined by the claims of the present invention.

Claims

1. A finite element calculation method for solidification heat transfer of beam blank continuous casting, characterized in that: The method comprises the following steps: Step 1: According to the cross-sectional dimensions of the profiled blank, a 2D geometric model of the 1 / 2 profiled blank is established using finite element analysis software; Step 2: Mesh the 1 / 2 profiled blank 2D geometric model, establish the 1 / 2 profiled blank 2D finite element model and save the mesh information; Step 3: Aiming at the solidification heat transfer problem of beam blank continuous casting, a two-dimensional finite element solidification heat transfer model of 1 / 2 beam blank is established; Step 4: Use the 1 / 2 profiled blank two-dimensional finite element solidification heat transfer model to calculate the temperature field changes during the continuous casting process; The step 3 comprises the following steps: Step 3.1: Construct the basic equation of finite element heat transfer; Step 3.2: Construct a finite element calculation model for transient heat transfer problems; Step 3.3: Construct a planar three-node triangular element and perform heat transfer analysis based on the element for different heat transfer boundary conditions; Step 3.4: Construct a planar four-node quadrilateral element and perform heat transfer analysis based on the element for different heat transfer boundary conditions; The step 4 comprises the following steps: Step 4.1: Read the steel grade information of the profiled billet, the structural parameters of the continuous casting machine and the continuous casting process parameters, and calculate the solidus temperature and the liquidus temperature; Step 4.2: Import the mesh information of the 1 / 2 profiled blank 2D finite element model; Step 4.3: Determine the time step; the time step is the time for calculating the slice movement each time; The slice represents the two-dimensional finite element model of the 1 / 2 profiled blank; Step 4.4: Determine the unit phase region and calculate the physical property parameters of the unit; Step 4.5: Calculate the unit heat transfer matrix, node load, unit heat capacity matrix and integrate the total heat transfer matrix, total temperature load, and total heat capacity matrix; the integrated total heat transfer matrix is ​​to decompose each unit heat transfer matrix into sub-matrices, and then "fill in" the block sub-matrices determined by each node sequence into the total heat transfer matrix according to the unit node number; the integrated total temperature load and total heat capacity matrix are the same as the integrated total heat transfer matrix; Step 4.6: Solve the temperature field, that is, solve the finite element two-dimensional unsteady heat conduction control equation At the same time, determine whether the slice position exceeds the air cooling zone. If so, the calculation ends. If not, the slice continues to move and go to step 4.

4.

2. The finite element calculation method for solidification heat transfer of beam blank continuous casting according to claim 1 is characterized in that: The grid information includes unit information, node information, and boundary information; the unit information includes the unit number and the nodes contained in the unit; the node information includes the node number and the node coordinates; the boundary includes the upper boundary, i.e. the inner arc side, the lower boundary, i.e. the outer arc side, and the left boundary or the right boundary, i.e. the side arc side; the boundary information includes the unit information and node information on the boundary.

3. The finite element calculation method for solidification heat transfer of beam blank continuous casting according to claim 1 is characterized in that: The construction method of the planar three-node triangular unit is: three nodes are encoded as i, j, and m in a counterclockwise direction to construct a planar three-node triangular heat transfer unit; the construction method of the planar four-node quadrilateral unit is: four vertices are selected as interpolation nodes, and the nodes are encoded as i, j, k, and l in a counterclockwise direction to construct a typical planar four-node quadrilateral unit.

4. The finite element calculation method for solidification heat transfer of beam blank continuous casting according to claim 3 is characterized in that: The heat transfer analysis based on the plane four-node quadrilateral unit for different heat transfer boundary conditions is as follows: a square is selected as a reference square unit in the parameter space; a geometric mapping function of the square to a quadrilateral in the physical space is written as a combination of node shape functions. The square node shape function is The finite element method is used to transform the calculation method of the temperature stiffness matrix element of any quadrilateral unit in physical space to the global coordinate function into the calculation method of the square unit to the local coordinate in parameter space: The heat capacity matrix elements for any quadrilateral cell in physical space are also calculated in square cells in parameter space: Nodal loads due to heat sources in any quadrilateral element in physical space are also calculated in square elements in parameter space: In the above formula, x, y are Cartesian coordinate systems; ξ, η are coordinate systems in parameter space; x i 、x j 、x k 、x l are the horizontal coordinates of nodes i, j, k, and l respectively; y i ,y j ,y k ,y l are the ordinates of nodes i, j, k, and l respectively; k mn is the element of unit temperature stiffness matrix; k is the heat transfer coefficient; N m and N n are all shape functions; a and b are the upper and lower limits of the integral in the parameter space, which are the coordinate variation ranges of the square respectively; J is the Jacobian matrix; |J is the Jacobian matrix determinant; J * is the adjoint matrix of the Jacobian matrix; c mn represents the unit heat capacity matrix element; ρ is the material density; c T is the specific heat of the material; p m represents the node load; q v is the heat generated per unit volume inside the material, called the internal heat source density; l is the heat transfer boundary length of the unit.

5. The finite element calculation method for solidification heat transfer of beam blank continuous casting according to claim 1 is characterized in that: The steel type information includes steel grade and chemical composition; the continuous casting machine structural parameters include crystallizer height and secondary cooling zone structural parameters; the secondary cooling zone structural parameters include the number, length, inlet position and outlet position of secondary cooling zones; the continuous casting process parameters include casting temperature, drawing speed, crystallizer heat flux density, water flow and water temperature in each secondary cooling zone, ambient temperature and cross-sectional dimensions of the profiled blank; the cross-sectional dimensions of the profiled blank include wing tip thickness, wide face width, narrow face width, web length, web thickness, inner edge arc radius, flange radius and fillet radius.

6. The finite element calculation method for solidification heat transfer of beam blank continuous casting according to claim 5 is characterized in that: The phase region refers to the liquid phase region, the solid phase region and the solid-liquid two-phase region; the physical property parameters include unit density, unit specific heat, unit solid phase ratio and unit thermal conductivity coefficient.

7. The finite element calculation method for solidification heat transfer of beam blank continuous casting according to claim 6 is characterized in that: The method for determining the unit phase zone is: firstly, the center temperature of the plane four-node quadrilateral unit is calculated by the inverse distance weighted interpolation method, and the center temperature of the triangle unit is calculated by the two-dimensional linear interpolation method; finally, the phase zone of the unit is determined according to the unit center temperature.

8. The finite element calculation method for solidification heat transfer of beam blank continuous casting according to claim 1, characterized in that: Solve the matrix equation using Gauss-Seidel iteration method Calculate the temperature at each node.

Images