Nonlinear efficient iterative solving method for thermal strain load

Abstract

The invention relates to the technical field of computer-aided engineering, in particular to a nonlinear efficient iterative solving method for thermal strain load, which is used for solving an obtained finite element equation by using a modified Newton-Raphson method and comprises the following steps of: firstly, establishing a geometric model to be analyzed and performing grid division; the method comprises the following steps: firstly, setting a model to be analyzed, then setting boundary conditions and temperature field change conditions of the model to be analyzed, finally, inputting the conditions into an algorithm solver, solving a nonlinear equation set, iteratively inputting a load continuously according to a modified Newton-Raphson method in the solving process, calculating node displacement until convergence, and solving stress-strain conditions in the welding process according to a material constitutive relation. In the solving process, by judging the size of the two norms of the displacement increment obtained in the adjacent increment step heuristic stage, the convergence of nonlinear solving of the structure is enhanced through the calculation mode, and the solving speed of the nonlinear problem in the thermal elastic-plastic constitutive structure is increased.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of computer aided engineering, in particular to a nonlinear efficient iterative solution method of thermal strain load. BACKGROUND

[0002] In the welding process, complex physical and metallurgical phenomena are accompanied by high-density energy input, which can form a significant temperature gradient in the welding area. The significant temperature gradient can directly induce severe plastic strain. At the same time, the welding material can undergo solid phase transformation and molten flow. The above influences ultimately lead to significant welding deformation and residual stress in the local structure, and directly affect the structural integrity and long-term reliability of the key components under normal operating conditions and extreme service environments.

[0003] In the numerical simulation analysis of thermal-mechanical coupled engineering structures, the Newton-Raphson iteration method is often used to solve nonlinear problems, that is, the nonlinear problem is linearized into multiple linear equation systems by linearization, so as to realize nonlinear solution. In the welding simulation problem, the displacement increment needs to be solved according to the temperature field change at each time step, and the current stress is obtained according to the displacement increment, and the cycle is repeated until the welding process is completed. However, the constitutive relationship of the material is a nonlinear term with thermal physical curve, and the solving process will cause significant nonlinear effect. The temperature field and displacement field need to be recalculated through multiple iterations to solve the thermal load term, which not only leads to slow overall calculation process, but also easily causes numerical instability problems such as stiffness matrix singularity due to load mutation and material property nonlinear mutation, which affects the effectiveness and result accuracy of the iterative solution process. Therefore, a nonlinear efficient iterative solution method of thermal strain load is proposed. SUMMARY

[0004] In order to solve the above technical problems in the prior art, the present application provides a nonlinear efficient iterative solution method of thermal strain load.

[0005] To solve the above technical problems, the present application provides the following technical scheme: a nonlinear efficient iterative solution method of thermal strain load, the nonlinear efficient iterative solution method of thermal strain load comprising the following steps:

[0006] S1, using Hypermesh to discretize the three-dimensional geometric model of the workpiece to be welded into tetrahedral mesh elements, obtaining mesh data and node data, and importing the calculation program;

[0007] S2, setting the heat source model, welding parameters, boundary conditions and material parameters required for welding;

[0008] S3, constructing a heat transfer matrix through a three-dimensional heat conduction equation, and then calculating the temperature field change result with time through the heat transfer matrix;

[0009] S4, converting the unit temperature difference into the equivalent load of the node based on the virtual displacement principle;

[0010] S5, calculating the residual of the node internal stress and the external load in each iteration process by the modified Newton-Raphson method, constantly correcting the load term of the linear equation system, constructing the linear equation system, and solving the matrix equation to obtain the node displacement, using the radial return algorithm, first calculating the trial stress and bringing it into the yield function of the workpiece to be welded in the state, judging the state of the workpiece to be welded, and updating the stress increment, and finally obtaining the accurate stress and strain state;

[0011] S6, calculating the residual and determining whether the result meets the convergence condition, if not, repeating the previous step, otherwise outputting the stress and strain result.

[0012] Preferably, in the step S2, the heat source model includes the heat source type, heat source size and heat source power required for welding.

[0013] The welding parameters include the welding speed.

[0014] The boundary conditions include relevant mechanical boundaries and heat exchange boundaries.

[0015] The material parameters are the corresponding material types of the workpiece to be welded in the material library, and the corresponding material properties are assigned to the geometric model.

[0016] Preferably, the step S3 specifically includes the following steps:

[0017] S31, inputting the heat source load changing with time to the nodes of the discrete tetrahedral grid elements, and constructing the heat conduction matrix and the balance equation.

[0018] S32, combining the temperature field distribution of the tetrahedral grid model at each time, and calculating the node temperature field data difference between adjacent time steps.

[0019] Preferably, the balance equation in the step S31 is:

[0020] ;

[0021] In the formula, is the specific heat capacity of the material, is the material density, is the interpolation function matrix form, is the unit temperature column array, is the derivative with respect to time, 、 、 is the heat conduction coefficient of the material in three directions, is the transpose of the shape function matrix, is the partial derivative of the shape function in x, y, z direction, is the heat exchange coefficient, is the area integral symbol, is the ambient temperature, is the heat source acting on the inside of the workpiece;

[0022] that is,

[0023] ;

[0024] in the formula, is the unit heat conduction matrix, is the unit temperature of the calculation step, is the unit heat capacity matrix, is the unit temperature of the previous calculation step, is the interval time between two calculation steps, is the unit load array;

[0025] The unit heat conduction matrix is assembled to obtain the overall heat transfer matrix and the heat capacity matrix, and the linear equation set to be solved after assembly and arrangement is:

[0026] ;

[0027] The final solution of the linear equation set obtains the temperature field distribution in the discrete domain at this moment, and the results are stored in a file.

[0028] Preferably, in the step S4, the formula for converting the unit temperature difference into the node equivalent load is:

[0029] ;

[0030] in the formula, is the thermal expansion coefficient of the material, which is related to the material, The components in three directions of the isotropic material are consistent, D is the constitutive matrix, is the node temperature difference array of the unit, which needs to be expanded according to the matrix operation, and the unit is .

[0031] Preferably, the step S5 specifically includes the following steps:

[0032] S51, according to the elastic modulus and Poisson's ratio of the material, the unit stiffness matrix is constructed, and then the unit stiffness matrix is assembled into the overall tangent stiffness matrix through the corresponding node number, specifically:

[0033] ;

[0034] in the formula, is the strain-displacement matrix, D is the constitutive matrix, and

[0035] S52, the initial equivalent node load vector array is brought into the global tangent stiffness matrix, and the linear equation group to be solved is constructed as:

[0036] ;

[0037] In the formula, is the node displacement array, and the node displacement increment is obtained by solving the linear equation group;

[0038] S53, the residual term at the beginning of iteration is calculated by the difference between the initial equivalent node load and the node internal stress and external force load, and in the subsequent iteration process, the difference between the residual of the previous step and the node internal stress increment caused by the strain increment in the iteration step is obtained, and the unit load formula in the iteration process is:

[0039] ;

[0040] ;

[0041] In the formula, is the true strain calculated by the displacement and strain-displacement matrix, is the strain increment of the last iteration step, is the thermal strain of the free expansion of the material in the unconstrained state, is the external force load on the node, is the total load term at the nth iteration, is the thermal expansion coefficient of the material, is the volume integral of the unit, is the constitutive matrix after plastic correction, and the expression is:

[0042] ;

[0043] In the formula, D is the constitutive matrix in the elastic state, T is the node temperature, F is the yield function under the Mises yield criterion, is the normal direction along the yield surface F=0 in the stress space;

[0044] By continuously correcting the load term, the linear equation group is solved iteratively, so that the calculated node displacement approaches the correct value;

[0045] S54, whether the convergence condition is met is determined by the two norm of displacement, if the two norm is less than 1e-5, it is determined that the convergence is met, the iteration is ended, the stress and strain results are output, and the iteration of the next time step is carried out; if it does not converge, it is repeated, and the convergence judgment value is calculated as follows:

[0046] ;

[0047] In the formula, This represents the total number of nodes in the discrete structure. , , Let Re be the components of the node displacement in the three directions, and Re be the convergence calculation value. The convergence calculation value is compared with the convergence criterion value 1e-5. If it is less than 1e-5, then convergence is determined.

[0048] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0049] This invention utilizes a modified Newton-Raphson method to solve strongly nonlinear problems arising from thermophysical property curves and material constitutive models. It reduces the minimum number of iterations required for convergence, improves the solution speed for nonlinear problems using thermo-elastic-plastic constitutive models, and reduces the overall time consumption for structural welding simulation. Furthermore, it is adaptable to different material constitutive equations and is widely applicable to nonlinear solutions in thermo-mechanical coupled field simulations. Attached Figure Description

[0050] Figure 1 This is a schematic representation of the number of iterations in the improved NR method of the present invention.

[0051] Figure 2 This is a schematic diagram of the geometric model of Embodiment 2 of the present invention;

[0052] Figure 3 This is a schematic diagram of the discretized mesh model according to Embodiment 2 of the present invention;

[0053] Figure 4 This is a schematic diagram of stress distribution in the calculation process of Embodiment 2 of the present invention;

[0054] Figure 5 This is a schematic diagram of the traditional NR iteration process of the present invention;

[0055] Figure 6 This is a schematic diagram of the improved NR iteration process of the present invention. Detailed Implementation

[0056] The present invention will be further described below with reference to the accompanying drawings and embodiments, which illustrate the above and other technical features and advantages of the present invention. However, the following embodiments are merely preferred embodiments of the present invention and are not exhaustive.

[0057] Example 1:

[0058] like Figure 1 , Figure 5 and Figure 6 As shown, this invention provides a nonlinear efficient iterative solution method for thermal strain loads, which includes the following steps:

[0059] S1 uses Hypermesh to discretize the three-dimensional geometric model of the workpiece to be welded into tetrahedral mesh elements according to the set meshing parameters, realizes the spatial discretization of geometric models of different shapes, obtains mesh data and node data, and imports them into the calculation program.

[0060] S2, set the heat source model, welding parameters, boundary conditions and material parameters required for welding;

[0061] S3, a heat transfer matrix is ​​constructed through the three-dimensional heat conduction equation, and then the temperature field changes with time is calculated through the heat transfer matrix;

[0062] S4, based on the principle of virtual displacement, converts the temperature difference of the elements (all elements are tetrahedral mesh elements) into equivalent loads at the nodes;

[0063] S5. The residuals of internal stress and external load at each iteration node are calculated by the modified Newton-Raphson method. The load terms of the linear equation system are continuously modified to form a linear equation system. The matrix equation is solved to obtain the node displacement. The radial return algorithm is used to first obtain the trial stress and substitute it into the yield function of the workpiece to be welded under the temperature state to determine the state of the workpiece to be welded and update the stress increment. Finally, the accurate stress-strain state is obtained.

[0064] S6, calculate the residuals and determine whether the results meet the convergence conditions. If they do not converge, repeat the previous step; otherwise, output the stress and strain results.

[0065] In this embodiment, in step S2, the heat source model includes the type of heat source required for welding, the size of the heat source, and the power of the heat source;

[0066] Welding parameters include welding speed;

[0067] Boundary conditions include relevant mechanical boundaries and heat transfer boundaries;

[0068] Material parameters refer to the material type of the workpiece to be welded in the material library, and the corresponding material properties are assigned to the geometric model.

[0069] In this embodiment, step S3 specifically includes the following steps:

[0070] S31, the time-varying heat source load is input to the nodes of the discrete tetrahedral mesh element to construct the heat conduction matrix and equilibrium equation;

[0071] S32, combining the temperature field distribution of the tetrahedral mesh model at each time step, calculates the difference in node temperature field data between adjacent time steps.

[0072] In this embodiment, the equilibrium equation in step S31 is:

[0073] ;

[0074] In the formula, For the specific heat capacity of the material, For material density, It is in the form of an interpolation function matrix. For unit temperature array, For the time derivative, , , represents the thermal conductivity coefficient of the material in three directions. For the transpose of the shape function matrix, Let be the partial derivatives of the shape function in the x, y, and z directions. The heat transfer coefficient, The symbol for area fractions. For ambient temperature, It is a heat source acting inside the workpiece;

[0075] Right now:

[0076] ;

[0077] In the formula, For the unit heat conduction matrix, The unit temperature for this calculation step. For the element heat capacity matrix, The unit temperature from the previous calculation step. The time interval between two calculation steps. For element load array;

[0078] The unit heat conduction matrices are then assembled to obtain the overall heat transfer matrix and heat capacity matrix. The assembled and rearranged linear equation set is as follows:

[0079] ;

[0080] Finally, the linear equations are solved to obtain the temperature field distribution in the discrete domain at that moment, and the results are saved to a file.

[0081] In this embodiment, the formula for converting the element temperature difference into the nodal equivalent load in step S4 is as follows:

[0082] ;

[0083] In the formula, This is the coefficient of thermal expansion of the material, which is material-dependent. For isotropic materials, the components in all three directions remain consistent, and D is the constitutive matrix. This is the nodal temperature difference array of the element, which is expanded according to the needs of matrix operations. The element contains... .

[0084] In this embodiment, step S5 specifically includes the following steps:

[0085] S51, based on the material's elastic modulus and Poisson's ratio, construct the element stiffness matrix, and then assemble the element stiffness matrix into the global tangent stiffness matrix using the corresponding node numbers, specifically:

[0086] ;

[0087] In the formula, The strain-displacement matrix is... Let D be the global tangent stiffness matrix and D be the constitutive matrix.

[0088] S52, substituting the initial equivalent nodal load vector column matrix into the global tangent stiffness matrix, the resulting system of linear equations to be solved is as follows:

[0089] ;

[0090] In the formula, The nodal displacement matrix is ​​used, and the nodal displacement increments are obtained by solving a system of linear equations.

[0091] S53, the residual term at the beginning of the iteration is calculated from the difference between the initial equivalent nodal load and the nodal internal stress and external load. The residual term in subsequent iterations is calculated from the difference between the residual of the previous step and the nodal internal stress increment caused by the strain increment in the current iteration step. The element load formula during the iteration is:

[0092] ;

[0093] ;

[0094] In the formula, To calculate the true strain from the displacement and strain-displacement matrix, This is the strain increment compared to the previous iteration. The thermal strain represents the free expansion of the material under unconstrained conditions. External force loads on the nodes, This is the total load term in the nth iteration. The coefficient of thermal expansion of the material. For the unit volume integral, The constitutive matrix after plasticity correction is expressed as:

[0095] ;

[0096] In the formula, D is the constitutive matrix under elastic conditions, T is the nodal temperature, and F is the yield function under the Mises yield criterion. This is along the normal direction of the yield surface F=0 in stress space;

[0097] By continuously correcting the load terms and iteratively solving the linear equation system, the calculated nodal displacements can be made to approximate the correct values.

[0098] S54 determines whether the convergence condition is met by using the L2 norm of the displacement. If the L2 norm is less than 1e-5, convergence is determined, the iteration ends, the stress-strain results are output, and the iteration proceeds to the next time step. If convergence is not achieved, the iteration is repeated. The convergence criterion value is calculated as follows:

[0099] ;

[0100] In the formula, This represents the total number of nodes in the discrete structure. , , Let Re be the components of the node displacement in the three directions, and Re be the convergence calculation value. The convergence calculation value is compared with the convergence criterion value 1e-5. If it is less than 1e-5, then convergence is determined.

[0101] Example 2:

[0102] like Figure 2 As shown, a T-shaped weld was performed on sample 1 (100mm long and wide, 10mm thick) and sample 2 (200mm long, 100mm wide, 10mm thick). The blue area represents the solder filling area. Both samples are made of 1050A aluminum alloy with an elastic modulus of 6.97-3.68e11 N / m², varying with temperature. A Poisson's ratio of 0.3 was applied using a double ellipsoidal heat source with a loading power of 1000W and a velocity of 10mm / s. Fixed constraints were applied to the upper and lower surfaces. The discretized tetrahedral element mesh number was 28838. Figure 3 As shown.

[0103] Based on this model, under the same temperature field input, the proposed method in this embodiment and the traditional NR iterative method are used to solve the problem and obtain the stress and strain conditions.

[0104] The above description is merely a preferred embodiment of the present invention and is illustrative rather than restrictive. Those skilled in the art will understand that many changes, modifications, and even equivalents can be made within the spirit and scope defined by the claims of the present invention, all of which will fall within the protection scope of the present invention.

Claims

1. A nonlinear, efficient iterative solution method for thermal strain loads, characterized in that, The nonlinear efficient iterative solution method includes the following steps: S1. Use Hypermesh to discretize the three-dimensional geometric model of the workpiece to be welded into tetrahedral mesh elements, obtain mesh data and node data, and import them into the calculation program. S2, set the heat source model, welding parameters, boundary conditions and material parameters required for welding; S3, a heat transfer matrix is ​​constructed through the three-dimensional heat conduction equation, and then the change of temperature field with time is calculated through the heat transfer matrix; S4, based on the principle of virtual displacement, converts the element temperature difference into an equivalent nodal load; S5. The residuals of internal stress and external load at each iteration process are calculated by the modified Newton-Raphson method. The load terms of the linear equation system are continuously modified to form a linear equation system. The matrix equation is solved to obtain the nodal displacement. The radial return algorithm is used to first obtain the trial stress and substitute it into the yield function of the workpiece to be welded under this state to determine the state of the workpiece to be welded and update the stress increment. Finally, the accurate stress-strain state is obtained. S6, calculate the residuals and determine whether the results meet the convergence conditions. If they do not converge, repeat the previous step; otherwise, output the stress and strain results.

2. The nonlinear efficient iterative solution method for thermal strain load as described in claim 1, characterized in that, In step S2, the heat source model includes the type of heat source, the size of the heat source, and the power of the heat source required for welding. The welding parameters include welding speed; The boundary conditions include relevant mechanical boundaries and heat transfer boundaries; The material parameters refer to the material type of the workpiece to be welded in the material library, and the corresponding material properties are assigned to the geometric model.

3. The nonlinear efficient iterative solution method for thermal strain load as described in claim 1, characterized in that, Step S3 specifically includes the following steps: S31, the time-varying heat source load is input to the nodes of the discrete tetrahedral mesh element to construct the heat conduction matrix and equilibrium equation; S32, combining the temperature field distribution of the tetrahedral mesh model at each time step, calculates the difference in node temperature field data between adjacent time steps.

4. The nonlinear efficient iterative solution method for thermal strain load as described in claim 1, characterized in that, The equilibrium equation in step S31 is: ； In the formula, For the specific heat capacity of the material, For material density, It is in the form of an interpolation function matrix. For unit temperature array, For the time derivative, , , The thermal conductivity coefficients of the material in three directions are . For the transpose of the shape function matrix, Let be the partial derivatives of the shape function in the x, y, and z directions. The heat transfer coefficient, The symbol for area fractions. For ambient temperature, It is a heat source acting inside the workpiece; Right now: ； In the formula, For the unit heat conduction matrix, The unit temperature for this calculation step. For the element heat capacity matrix, The unit temperature from the previous calculation step. The time interval between two calculation steps. For element load array; The unit heat conduction matrices are then assembled to obtain the overall heat transfer matrix and heat capacity matrix. The assembled and rearranged linear equation set is as follows: ； Finally, the linear equations are solved to obtain the temperature field distribution in the discrete domain at that moment, and the results are saved to a file.

5. The nonlinear efficient iterative solution method for thermal strain load as described in claim 1, characterized in that, In step S4, the formula for converting the unit temperature difference into the nodal equivalent load is as follows: ； In the formula, This is the coefficient of thermal expansion of the material, which is material-dependent. For isotropic materials, the components in all three directions remain consistent, and D is the constitutive matrix. This is the nodal temperature difference array of the element, which is expanded according to the needs of matrix operations. The element contains... .

6. The nonlinear efficient iterative solution method for thermal strain load as described in claim 1, characterized in that, Step S5 specifically includes the following steps: S51, based on the material's elastic modulus and Poisson's ratio, construct the element stiffness matrix, and then assemble the element stiffness matrix into the global tangent stiffness matrix using corresponding node numbers, specifically: ； In the formula, The strain-displacement matrix, This is the global tangent stiffness matrix; S52, substituting the initial equivalent nodal load vector column matrix into the global tangent stiffness matrix, the resulting system of linear equations to be solved is as follows: ； In the formula, The nodal displacement matrix is ​​used, and the nodal displacement increments are obtained by solving a system of linear equations. S53, the residual term at the beginning of the iteration is calculated from the difference between the initial equivalent nodal load and the nodal internal stress and external load. The residual term in subsequent iterations is obtained from the difference between the residual of the previous step and the nodal internal stress increment caused by the strain increment in that iteration step. The element load formula during the iteration is: ； ； In the formula, To calculate the true strain from the displacement and strain-displacement matrix, This is the strain increment compared to the previous iteration. The thermal strain represents the free expansion of the material under unconstrained conditions. External force loads on the nodes, This is the total load term in the nth iteration. The coefficient of thermal expansion of the material. For the unit volume integral, The constitutive matrix after plasticity correction is expressed as: ； In the formula, D is the constitutive matrix under elastic conditions, T is the nodal temperature, and F is the yield function under the Mises yield criterion. This is along the normal direction of the yield surface F=0 in stress space; By continuously correcting the load terms and iteratively solving the linear equation system, the calculated nodal displacements can be made to approximate the correct values. S54 determines whether the convergence condition is met by using the L2 norm of the displacement. If the L2 norm is less than 1e-5, convergence is determined, the iteration ends, the stress-strain results are output, and the iteration proceeds to the next time step. If convergence is not achieved, the iteration is repeated. The convergence criterion value is calculated as follows: ； In the formula, This represents the total number of nodes in the discrete structure. , , Let Re be the components of the node displacement in the three directions, and Re be the convergence calculation value. The convergence calculation value is compared with the convergence criterion value 1e-5. If it is less than 1e-5, then convergence is determined.

Images