A finite element modeling method for simulating explosive welding of nickel-copper steel

Abstract

The present invention provides a finite element modeling method for simulating nickel-copper steel explosion welding, and relates to the technical field of computer simulation. The present invention establishes a three-dimensional model and sets material properties and initial conditions, utilizes the birth and death unit function in combination with welding parameters to simulate the explosion welding process, refines the unit based on an adaptive refinement algorithm in combination with strain gradients and temperature gradients, maps grid control point variables to the next time period based on a Lagrange interpolation method, and outputs curves of temperature, pressure and strain changes over time as well as temperature fields, pressure fields and strain fields at each moment, thereby achieving modeling.

Description

Technical Field

[0001] The invention relates to the technical field of computer simulation, in particular to a finite element modeling method for simulating nickel-copper-steel explosion welding. Background Art

[0002] Explosive welding is an efficient welding method that uses the huge impact energy generated by the explosion to firmly combine two or more materials with different properties in a short time. The explosion process involves rapid changes in multiple physical properties including temperature, pressure and strain. Traditional experimental methods are difficult to capture and analyze these rapidly changing physical quantities in a timely and accurate manner. With the development of finite element technology, a new idea has been provided for analyzing the process of explosive welding. Finite element modeling is used to numerically simulate and predict the explosive welding process, so as to carefully analyze the response behavior of different materials under such extreme conditions, and then optimize the process parameters and improve the welding quality. However, the difficulty in the current finite element modeling and analysis lies in how to reproduce the transient temperature field, pressure field and strain field in the explosive welding process and ensure the accuracy of the mesh in the simulation process. There are still many deficiencies in dealing with high-temperature, high-speed dynamic processes and complex multi-material combinations. Therefore, it is particularly important to develop a sophisticated and efficient finite element modeling method for explosive welding.

[0003] In the prior art, publication number CN111428416A discloses a finite element modeling method for simulating high-energy beam welding. A geometric image represented by NURBS is used as a welding structure model to obtain a Bezier unit, welding process parameters and material property parameters are input, a weld life and death unit set of the welding structure model is established, and a local refined Bezier unit is generated based on a quadtree or octree grid subdivision strategy. The Euclidean norm is used to map the control variables of two adjacent grid points in two specific time periods, and the welding temperature thermal cycle curve of each sampling point in the entire welding time history and the temperature field at each moment are output.

[0004] The main problems with the above method are: modeling based on geometric figures represented by NURBS, the calculation of the model is more complicated when dealing with mesh refinement and unit life and death state operations, which increases the difficulty of solving the model; the mesh is subdivided through a quadtree or octree strategy, while the quadtree and octree strategies are suitable for static or gradually changing detail requirements, and are insufficient to capture details for the simulation process of fast changes and dynamic special effects, resulting in excessive refinement and increased calculation; the final output of the welding temperature thermal cycle curve and temperature field data fails to take into account physical quantities such as pressure field and strain field that have a greater impact on the explosive welding process, resulting in an incomplete analysis of the welding process.

[0005] The above information disclosed in this Background section is only for enhancement of understanding of the background of the present disclosure and therefore it may contain information that does not constitute the prior art that is already known to one of ordinary skill in the art. Summary of the invention

[0006] The object of the present invention is to provide a finite element modeling method for simulating nickel-copper steel explosion welding to solve the problems raised in the above background technology.

[0007] To achieve the above object, the present invention provides the following technical solutions:

[0008] A finite element modeling method for simulating nickel-copper steel explosion welding, the specific steps comprising:

[0009] Step 1: Establish a three-dimensional model of the explosive welding structure through finite element software, select Bezier unit for meshing, collect material properties of nickel, copper and steel materials, set initial temperature and initial stress, set the degree of freedom of all nodes at the bottom of the model to 0, and calculate the explosion pressure load at the top of the model as the boundary condition. The material properties include elastic modulus, Poisson's ratio, density, thermal conductivity and specific heat capacity;

[0010] Step 2: Use the built-in birth and death unit function of the finite element software, combine the welding speed and energy input, define the birth and death unit set in the welding interface area, set the time window, update the birth and death unit state in each time window, and simulate the welding cladding process in discrete time periods;

[0011] Step 3: Based on the grid adaptive refinement algorithm, a refinement factor is generated according to the strain gradient and temperature gradient calculated in real time to determine whether unit refinement is needed. For the units that need to be refined, a local refined Bezier unit is generated for the bonding interface area of ​​each time window, and the temperature field, pressure field and strain field corresponding to the discrete time series are generated;

[0012] Step 4: Mapping the grid control point variables of the current time period to the next time period by Lagrange interpolation method, wherein the grid control point variables include temperature, pressure and strain;

[0013] Step 5: Output the explosion welding temperature curve, pressure curve and strain curve of each observation point in the entire explosion welding time history, as well as the temperature field, pressure field and strain field at each moment.

[0014] Furthermore, the formula for calculating the explosion pressure load on the top of the model is:

[0015] P(t)=P0·e -αt

[0016] Wherein, P(t) represents the explosion pressure at time t, t represents the time elapsed from the moment of the explosion, P0 represents the initial explosion pressure, and α represents the attenuation constant.

[0017] Furthermore, the principle for defining the life and death unit set and updating the life and death unit status in each time window is as follows:

[0018] The welding interface area is selected as the definition area of ​​the dead and alive unit set, including the interface of nickel, copper and steel contact. The initial state is set as a dead unit. The formula for generating energy input is:

[0019]

[0020] Among them, Q represents the energy input per unit area during welding, P represents the total input energy, and A represents the area of ​​the welding interface;

[0021] In each time window, the units in the weld interface area determine the transition between life and death units according to the changes in welding speed and energy input, based on the formula:

[0022]

[0023] Where i, j, and k represent the coordinate indexes of the Bezier unit in the x, y, and z directions, respectively. D(i,j,k) represents the life and death unit state of the Bezier unit with coordinates (i,j,k). represents the position vector of the Bezier cell with coordinates (i,j,k) in three-dimensional space, represents the welding speed vector, t m represents the mth time window, and m represents the index of the time window.

[0024] Furthermore, the principle for generating the local refined Bezier elements for the bonding interface region for each time window is based on:

[0025] The formulas for generating strain gradient and temperature gradient are:

[0026]

[0027] in, represents the strain gradient of the mth time window, ε xx (t m ) represents the strain tensor in the x direction caused by the displacement in the x direction in the mth time window, represents the strain tensor in the y direction caused by the displacement in the y direction in the mth time window, represents the strain tensor in the z direction caused by the displacement in the z direction in the mth time window, ε xy (tm ) represents the strain tensor in the x direction caused by the displacement in the y direction in the mth time window, ε yz (t m ) represents the strain tensor in the y direction caused by the displacement in the z direction in the mth time window, ε zx (t m ) represents the strain tensor in the z direction caused by the displacement in the x direction in the mth time window, where x, y, and z represent directions in three-dimensional space, respectively. represents the temperature gradient of the mth time window, T(t m ) represents the temperature of the mth time window;

[0028] Generate a refinement factor and determine whether unit refinement is needed based on the refinement factor and the refinement threshold. The formula is:

[0029]

[0030] Among them, f(t m ) represents the refinement factor of the mth time window, ε thresh represents the strain gradient threshold, T thresh represents the temperature gradient threshold;

[0031] When f(t m )>f thresh When the unit is refined, f thresh represents the refinement threshold;

[0032] The control points are adjusted according to the Bezier difference function, and the Bezier element is refined based on the adjusted control points. The formula is:

[0033]

[0034] Among them, B(u,v,w) represents the Bezier body, u,v,w respectively represent the parameters describing the three directions of constructing the Bezier body, p,q,r respectively represent the index of the control point in the u,v,w direction, I,J,K respectively represent the upper limit of the index in the u,v,w direction, P pqr represents the control points of the Bezier surface, represent the Bezier basis functions in the u, v, and w directions respectively;

[0035] When the refinement factor is greater than the refinement threshold, the upper limit of the index in the u, v, and w directions is multiplied by the refinement factor, and the Bezier basis function is updated to generate the refined control points, and the Bezier unit is refined based on the control points.

[0036] Furthermore, the principle for generating the temperature field, pressure field and strain field corresponding to the discrete time series is as follows:

[0037] The principle of generating temperature field is based on:

[0038] The discretized form of the heat conduction equation is:

[0039]

[0040] Where C represents the heat capacity matrix, represents the temperature vector, represents the temperature change rate, R represents the thermal conductivity matrix, represents the internal heat source vector;

[0041] The formula for generating the temperature field in discrete time based on the implicit backward Euler method is:

[0042]

[0043] Among them, T m+1 and T m They represent the temperature of the m+1th time window and the mth time window respectively, and △t represents the length of the time window;

[0044] The principle of generating pressure field is based on:

[0045] The discretized form of the dynamic equation is:

[0046]

[0047] Where M represents the mass matrix, represents the acceleration vector, N represents the damping matrix, represents the velocity vector, G represents the stiffness matrix, represents the displacement vector, represents the external load vector;

[0048] The formula for generating the pressure field in discrete time based on the central difference method is:

[0049]

[0050] in, Respectively represent the displacement vectors of the m+1th time window, the mth time window, and the m-1th time window, represents the external load vector of the mth time window;

[0051] The formula for generating strain fields is:

[0052]

[0053] Among them, ε m+1 represents the strain vector of the m+1th time window, and B represents the strain-displacement matrix.

[0054] Furthermore, the principle for mapping the grid control point variables of the current time period to the next time period by Lagrange interpolation is:

[0055] The Lagrange interpolation basis function is constructed based on the formula:

[0056]

[0057] S(X′,Y′,Z′)=S1·L1(X′,Y′,Z′)+S2·L2(X′,Y′,Z′)·S3·L3(X′,Y′,Z′)

[0058] Among them, (X1, Y1, Z1), (X2, Y2, Z2), and (X3, Y3, Z3) represent the control points of the variables in the current time period, S1, S2, and S3 represent the variable values ​​corresponding to the three control points, respectively. L1(X′, Y′, Z′), L2(X′, Y′, Z′), and L3(X′, Y′, Z′) are the Lagrangian basis functions of the interpolation points of the next time period corresponding to the control points of the three variables, respectively, and S(X′, Y′, Z′) represents the variable value of the interpolation point.

[0059] Compared with the prior art, the present invention has the following beneficial effects:

[0060] The present invention constructs a three-dimensional model according to the detailed material properties of nickel-copper steel to ensure that the material properties are accurately reflected in the simulation, providing precise initial conditions for the subsequent use of birth and death unit technology and grid adaptive refinement; through the birth and death unit function, combined with the welding speed and energy input, the birth and death unit state is dynamically updated to simulate the addition and removal of materials during the welding process, ensuring an accurate description of the thermal and mechanical behaviors during the welding process, obtaining more refined temperature fields and stress fields, and making the overall finite element simulation more realistic and efficient.

[0061] The present invention also uses a grid adaptive algorithm to adjust the grid density according to the real-time strain gradient and temperature gradient, performs necessary grid refinement in local stress concentration and heat-affected key areas, captures key features such as material deformation and heat conduction, and compared with global refinement, local refinement is both accurate and reduces the amount of calculation, improves resource utilization and the ability to respond to real-time changes. The generated discrete sequence can record in detail the changes in temperature, pressure and strain during the entire explosion time history, and combines the interpolation method to map the grid control point variables to the next time period, which can accurately simulate the subtle changes in the welding interface area, ensure smooth transition of data, effectively reduce time and space errors, avoid numerical oscillations, and improve the reliability of simulation results; at the same time, it outputs physical quantities closely related to the explosion welding process, including temperature, pressure and strain. The combination of the three variables can more accurately analyze the response characteristics of materials and structures under dynamic loading, monitor and adjust temperature, pressure and strain in real time, ensure that process parameters are within the optimal range, and effectively improve welding quality. BRIEF DESCRIPTION OF THE DRAWINGS

[0062] Figure 1 This is a schematic diagram of a method flow of an embodiment of the present invention;

[0063] Figure 2 This is a pressure cloud diagram of an embodiment of the present invention;

[0064] Figure 3 It is a temperature cloud diagram of an embodiment of the present invention;

[0065] Figure 4 It is a plastic strain cloud diagram of an embodiment of the present invention;

[0066] Figure 5 This is a schematic diagram of the metal welding structure before explosion according to an embodiment of the present invention;

[0067] Figure 6 This is a schematic diagram of the welding structure when the explosion starts in an embodiment of the present invention;

[0068] Figure 7 , Figure 8 The pressure-time curves of two observation points in the embodiment of the present invention;

[0069] Fig. 9 , Fig.10 The temperature-time curves of two observation points in the embodiment of the present invention;

[0070] Fig.11 , Fig.12 1 is a plastic strain-time curve of two observation points in an embodiment of the present invention. DETAILED DESCRIPTION

[0071] In order to make the objectives, technical solutions and advantages of the present invention more clearly understood, the present invention is further described in detail below in conjunction with specific embodiments.

[0072] It should be noted that, unless otherwise defined, the technical terms or scientific terms used in the present invention should be understood by people with ordinary skills in the field to which the present invention belongs. The words "first", "second" and similar words used in the present invention do not indicate any order, quantity or importance, but are only used to distinguish different components. "Include" or "comprise" and similar words mean that the elements or objects appearing before the word include the elements or objects listed after the word and their equivalents, without excluding other elements or objects. "Connect" or "connected" and similar words are not limited to physical or mechanical connections, but may include electrical connections, whether direct or indirect. "Up", "down", "left", "right" and the like are only used to indicate relative positional relationships. When the absolute position of the described object changes, the relative positional relationship may also change accordingly.

[0073] Example:

[0074] See also Figures 1 to 12 , the present invention provides a technical solution:

[0075] A finite element modeling method for simulating nickel-copper steel explosion welding, the specific steps comprising:

[0076] Step 1: Establish a three-dimensional model of the explosive welding structure through finite element software, select Bezier unit for meshing, collect material properties of nickel, copper and steel materials, set initial temperature and initial stress, set the degree of freedom of all nodes at the bottom of the model to 0, and calculate the explosion pressure load at the top of the model as the boundary condition. The material properties include elastic modulus, Poisson's ratio, density, thermal conductivity and specific heat capacity;

[0077] In this embodiment, explosion welding is carried out in three-dimensional space, but in order to facilitate analysis and display, various variables in the explosion process are analyzed through two-dimensional cloud maps, including pressure cloud maps, temperature cloud maps and plastic strain cloud maps. The two-dimensional cloud maps can more clearly show the changes in pressure, temperature and plastic strain during the explosion process, and it is easier to make horizontal comparisons between data, thereby analyzing the relationship between multiple variables in the explosion welding process.

[0078] In this embodiment, the formula for calculating the explosion pressure load on the top of the model is:

[0079] P(t)=P0·e -αt

[0080] Wherein, P(t) represents the explosion pressure at time t, t represents the time elapsed from the moment of the explosion, P0 represents the initial explosion pressure, and α represents the attenuation constant.

[0081] The explosion pressure reflects the explosion pressure applied to the top of the model at any time. It is a dynamic load that changes with time. The initial explosion pressure P0 reflects the highest pressure applied at the moment the explosion begins. It is obtained through experiments or calibration data. The attenuation constant is obtained by experiments. The larger the attenuation constant, the faster the explosion pressure decays. The explosion pressure at time t is proportional to the initial explosion pressure and inversely proportional to the time since the explosion began.

[0082] The degrees of freedom of all nodes at the bottom of the model are set to 0, which means that the bottom nodes of the model cannot move or rotate. The purpose is to simulate the real fixed support to accurately reflect the mechanical behavior in explosive welding.

[0083] Step 2: Use the built-in birth and death unit function of the finite element software, combine the welding speed and energy input, define the birth and death unit set in the welding interface area, set the time window, update the birth and death unit state in each time window, and simulate the welding cladding process in discrete time periods;

[0084] In this embodiment, the welding interface area refers to the area where welding materials are added and the materials are in contact during the welding process. The time for activating and deactivating the unit, the welding speed and the energy input are set by the birth and death unit function of the finite element software ANSYS.

[0085] In this embodiment, the principle for defining the birth and death unit set and updating the birth and death unit status in each time window is as follows:

[0086] The welding interface area is selected as the definition area of ​​the dead and alive unit set, including the interface of nickel, copper and steel contact. The initial state is set as a dead unit. The formula for generating energy input is:

[0087]

[0088] Among them, Q represents the energy input per unit area during welding, P represents the total input energy, and A represents the area of ​​the welding interface;

[0089] Measure the historical data of the energy generated by the explosion in the explosive welding area during the historical explosive welding process, and record the type and mass of explosives. Compare the type and mass of explosives used in the current explosive welding with those used in the historical data, find a group of historical data with the closest type and mass of explosives, and use the explosion energy measured by this group of historical data as the total energy currently input; within each time window, the unit in the welding interface area determines the transition of the state of the life and death unit according to the changes in the welding speed and energy input, and the formula based on this is:

[0090]

[0091] Where i, j, and k represent the coordinate indexes of the Bezier unit in the x, y, and z directions, respectively. D(i,j,k) represents the life and death unit state of the Bezier unit with coordinates (i,j,k). represents the position vector of the Bezier cell with coordinates (i,j,k) in three-dimensional space, represents the welding speed vector, t m represents the mth time window, and m represents the index of the time window.

[0092] Welding speed and energy input are key parameters in simulating explosive welding. The welding speed determines the movement rate of the welding interface, the energy input determines the energy distribution of the welding, and the length of the time window determines the time resolution of the simulation. A smaller time window can capture rapidly changing processes, while a larger time window may result in loss of details. Within each time window, the unit in the welding interface area decides whether to change from a "dead" state to a "live" state based on the changes in welding speed and energy input. The unit reflecting the welding interface area in the mth time window is The energy coverage range moved by the influence. When the coverage range is smaller than the energy input, it means that the unit has been affected by the welding process in the current time window. Therefore, the state changes from the "dead" state to the "alive" state, which means that the unit is changed from a unit not participating in the calculation to a unit participating in the calculation, and the temperature field, pressure field and strain field of the unit will be updated with the time window.

[0093] Step 3: Based on the grid adaptive refinement algorithm, a refinement factor is generated according to the strain gradient and temperature gradient calculated in real time to determine whether unit refinement is needed. For the units that need to be refined, a local refined Bezier unit is generated for the bonding interface area of ​​each time window, and the temperature field, pressure field and strain field corresponding to the discrete time series are generated;

[0094] In this embodiment, the principle for generating the local refined Bezier element for the bonding interface region in each time window is:

[0095] The formulas for generating strain gradient and temperature gradient are:

[0096]

[0097] in, represents the strain gradient of the mth time window, ε xx (t m ) represents the strain tensor in the x direction caused by the displacement in the x direction in the mth time window, represents the strain tensor in the y direction caused by the displacement in the y direction in the mth time window, represents the strain tensor in the z direction caused by the displacement in the z direction in the mth time window, ε xy (t m ) represents the strain tensor in the x direction caused by the displacement in the y direction in the mth time window, ε yz (t m ) represents the strain tensor in the y direction caused by the displacement in the z direction in the mth time window, ε zx (t m ) represents the strain tensor in the z direction caused by the displacement in the x direction in the mth time window, where x, y, and z represent directions in three-dimensional space, respectively. represents the temperature gradient of the mth time window, T(t m ) represents the temperature of the mth time window;

[0098] The strain gradient represents the partial derivative of the strain tensor with respect to position, and reflects the rate at which the internal strain of the material changes with spatial position. The higher the strain gradient, the greater the local deformation change of the material, and there may be stress concentration. During the explosive welding process, the junction of different materials will produce extremely high strain gradients under the action of the explosion impact, and a refined mesh is required to simulate these local changes; the temperature gradient reflects the rate at which the temperature changes with spatial position, and the heat flux density is proportional to the temperature gradient. Therefore, the size of the temperature gradient directly affects the heat conduction rate and direction. During the explosive welding process, the temperature of the material changes very rapidly and extremely high temperature gradients may occur locally, and a refined mesh is required to simulate these local changes.

[0099] Generate a refinement factor and determine whether unit refinement is needed based on the refinement factor and the refinement threshold. The formula is:

[0100]

[0101] Among them, f(t m ) represents the refinement factor of the mth time window, ε thresh represents the strain gradient threshold, T thresh represents the temperature gradient threshold;

[0102] The purpose of the refinement factor is to normalize the strain gradient and temperature gradient so that different dimensions can be compared under the same standard. Both the strain gradient and the temperature gradient affect the mesh refinement, so the maximum value between the two is taken as the indicator to determine the refinement.

[0103] When f(t m )>f thresh When the unit is refined, f thresh represents the refinement threshold;

[0104] Establish initial mesh models with different degrees of refinement, conduct preliminary simulations to obtain mesh independence, compare the calculation results of temperature, pressure and strain at different degrees of refinement, and when the degree of refinement increases but the calculation results change very little, it indicates convergence. Select the finest mesh as the reference mesh, set multiple refinement thresholds from large to small, perform numerical simulations on each set refinement threshold, record the convergence of the physical quantities in the simulation results, and select the refinement threshold whose convergence is closest to the reference mesh as the refinement threshold actually used;

[0105] The control points are adjusted according to the Bezier difference function, and the Bezier element is refined based on the adjusted control points. The formula is:

[0106]

[0107] Among them, B(u,v,w) represents the Bezier body, u,v,w respectively represent the parameters describing the three directions of constructing the Bezier body, p,q,r respectively represent the index of the control point in the u,v,w direction, I,J,K respectively represent the upper limit of the index in the u,v,w direction, P pqr represents the control points of the Bezier surface, represent the Bezier basis functions in the u, v, and w directions respectively;

[0108] The Bezier basis functions have the form:

[0109]

[0110] in, represents the binomial coefficient;

[0111] When the refinement factor is greater than the refinement threshold, the upper limit of the index in the u, v, and w directions is multiplied by the refinement factor, and the Bezier basis function is updated to generate the refined control points, and the Bezier unit is refined based on the control points.

[0112] The Bezier body B(u,v,w) defines a Bezier body composed of multiple Bezier surfaces in the u, v, and w directions in three-dimensional space. When generating a refined Bezier unit, new control points are generated between the vertices of the Bezier body for refinement according to the size of the refinement factor.

[0113] The principle for generating the temperature field, pressure field and strain field corresponding to the discrete time series is:

[0114] The principle of generating temperature field is based on:

[0115] The discretized form of the heat conduction equation is:

[0116]

[0117] Where C represents the heat capacity matrix, represents the temperature vector, represents the temperature change rate, R represents the thermal conductivity matrix, represents the internal heat source vector;

[0118] The heat capacity matrix is ​​composed of the volume of the Bezier unit, the density of the material and the specific heat capacity. Each element of the heat capacity matrix reflects the heat stored in each node of each Bezier unit under unit temperature change. The basic form of the heat capacity matrix is:

[0119] C=ρ·c·V e ∫Z T ZGar e

[0120] Where C represents the heat capacity matrix, ρ represents the density of the material, c represents the specific heat capacity of the material, V e represents the volume of the e-th finite element unit, e represents the e-th finite element unit currently being processed, Z represents the shape function, and Z T The transposed matrix representing the shape function;

[0121] The basic form of the thermal conductivity matrix is:

[0122]

[0123] Where H represents the thermal conductivity matrix, d represents the thermal conductivity of the material, and Z d represents the derivative matrix of the shape function, The transposed matrix representing the derivative matrix of the shape function;

[0124] The formula for generating the temperature field in discrete time based on the implicit backward Euler method is:

[0125]

[0126] Among them, T m+1 and T m They represent the temperature of the m+1th time window and the mth time window respectively, and △t represents the length of the time window;

[0127] The equation of the temperature field is based on the heat conduction equation and is obtained in the finite element method by the weighted residual method and time discretization;

[0128] The principle of generating pressure field is based on:

[0129] The discretized form of the dynamic equation is:

[0130]

[0131] Where m represents the mass matrix, represents the acceleration vector, N represents the damping matrix, represents the velocity vector, G represents the stiffness matrix, represents the displacement vector, Represents the external load vector; the mass matrix represents the distribution of matter in the system and its influence on the motion, and all characteristic groups in the mass matrix are positive numbers; the damping matrix reflects the energy loss inside and outside the system, that is, the energy consumed by the system during the motion through friction, internal dissipation of materials and other mechanisms; the stiffness matrix reflects the resistance of the material or structure to external forces, quantifies the reaction theory of unit displacement hardness, reflects the rigidity and deformation characteristics of the material, and determines the displacement response when the load is applied;

[0132] The formula for generating the pressure field in discrete time based on the central difference method is:

[0133]

[0134] in, Respectively represent the displacement vectors of the m+1th time window, the mth time window, and the m-1th time window, represents the external load vector of the mth time window;

[0135] The equations of the pressure field are based on the linear elastic dynamics equations, which discretize the equations of motion of the continuum;

[0136] The formula for generating strain fields is:

[0137]

[0138] Among them, ε m+1 represents the strain vector of the m+1th time window, and B represents the strain-displacement matrix.

[0139] The relationship between strain and displacement is described by the small deformation theory. When performing finite element discretization, the strain at the node is calculated by the shape function and the node displacement. The strain-displacement matrix represents the spatial derivative of the shape function. The specific form is:

[0140]

[0141] Among them, H h represents the h-th shape function, h represents the index of the shape function, represents the strain of the node displacement in the x direction, represents the strain of the node displacement in the y direction, represents the strain of the node displacement in the z direction;

[0142] Step 4: Mapping the grid control point variables of the current time period to the next time period by Lagrange interpolation method, wherein the grid control point variables include temperature, pressure and strain;

[0143] In this embodiment, the principle of mapping the grid control point variables of the current time period to the next time period by Lagrange interpolation method is as follows:

[0144] The Lagrange interpolation basis function is constructed based on the formula:

[0145]

[0146]

[0147] S(X′,Y′,Z′)=S1·L1(X′,Y′,Z′)+S2·L2(X′,Y′,Z′)·S3·L3(X′,Y′ r ,Z′)

[0148] Among them, (X1, Y1, Z1), (X2, Y2, Z2), and (X3, Y3, Z3) represent the control points of the variables in the current time period, S1, S2, and S3 represent the variable values ​​corresponding to the three control points, respectively. L1(X′, Y′, Z′), L2(X′, Y′, Z′), and L3(X′, Y′, Z′) are the Lagrangian basis functions of the interpolation points of the next time period corresponding to the control points of the three variables, respectively, and S(X′, Y′, Z′) represents the variable value of the interpolation point.

[0149] For each of the three grid control point variables, temperature, pressure and strain, a Lagrangian interpolation basis function is constructed. Three control points are selected for each variable to construct the basis function. All basis functions together constitute the Lagrangian interpolation polynomial. Based on the interpolation basis function of the current time period, the interpolation point variable value of the next time period is generated.

[0150] Step 5: Output the explosion welding temperature curve, pressure curve and strain curve of each observation point in the entire explosion welding time history, as well as the temperature field, pressure field and strain field at each moment.

[0151] In this embodiment, the observation point represents the welding interface area, and the temperature value, pressure value and strain value of the observation point in each time window are extracted, and curves of temperature, pressure and strain changing with time are generated in chronological order, and the temperature field, pressure field and strain field corresponding to the discrete time series are generated according to step 3. The temperature field, pressure field and strain field at the corresponding moment in the set time window are generated.

[0152] The above formulas are all dimensionless and numerical calculations. The formula is a formula for the most recent real situation obtained by collecting a large amount of data and performing software simulation. The preset parameters in the formula are set by technicians in this field according to actual conditions.

[0153] The above embodiments may be implemented in whole or in part by software, hardware, firmware or any other combination thereof. When implemented by software, the above embodiments may be implemented in whole or in part in the form of a computer program product. Those skilled in the art may appreciate that the units and algorithm steps of each example described in conjunction with the embodiments disclosed herein may be implemented by electronic hardware, or a combination of computer software and electronic hardware. Whether these functions are performed by hardware or software methods depends on the specific application and design constraints of the technical solution.

[0154] The units described as separate components may or may not be physically separated, and the components shown as units may or may not be physical units, and may be located in one place or distributed on multiple network units. Some or all of the units may be selected according to actual needs to achieve the purpose of the solution of this embodiment.

[0155] The above description is only a specific implementation manner of the present application, but the protection scope of the present application is not limited thereto. Any technician familiar with the technical field can easily think of changes or substitutions within the technical scope disclosed in the present application, which should be included in the protection scope of the present application.

Claims

1. A finite element modeling method for simulating nickel-copper steel explosion welding, characterized in that: The specific steps include: Step 1: Establish a three-dimensional model of the explosive welding structure through finite element software, select Bezier unit for meshing, collect material properties of nickel, copper and steel materials, set initial temperature and initial stress, set the degree of freedom of all nodes at the bottom of the model to 0, and calculate the explosion pressure load at the top of the model as the boundary condition. The material properties include elastic modulus, Poisson's ratio, density, thermal conductivity and specific heat capacity; Step 2: Use the built-in birth and death unit function of the finite element software, combine the welding speed and energy input, define the birth and death unit set in the welding interface area, set the time window, update the birth and death unit state in each time window, and simulate the welding cladding process in discrete time periods; Step 3: Based on the grid adaptive refinement algorithm, a refinement factor is generated according to the strain gradient and temperature gradient calculated in real time to determine whether unit refinement is required. For units that need to be refined, a local refined Bezier unit is generated for the bonding interface area of ​​each time window, and the temperature field, pressure field, and strain field corresponding to the discrete time series are generated; Step 4: Mapping the grid control point variables of the current time period to the next time period by Lagrange interpolation method, wherein the grid control point variables include temperature, pressure and strain; Step 5: Output the explosion welding temperature curve, pressure curve and strain curve of each observation point in the entire explosion welding time history, as well as the temperature field, pressure field and strain field at each moment.

2. A finite element modeling method for simulating nickel-copper steel explosion welding according to claim 1, characterized in that: The formula for calculating the explosion pressure load on the top of the model in step 1 is: P(t)=P0·e -αt Wherein, P(t) represents the explosion pressure at time t, t represents the time elapsed from the moment of the explosion, P0 represents the initial explosion pressure, and α represents the attenuation constant.

3. The finite element modeling method for simulating nickel-copper steel explosion welding according to claim 1, characterized in that: The principle for defining the birth and death unit set in step 2 and updating the birth and death unit status in each time window is: The welding interface area is selected as the definition area of ​​the dead and alive unit set, including the interface of nickel, copper and steel contact. The initial state is set as a dead unit. The formula for generating energy input is: Among them, Q represents the energy input per unit area during welding, P represents the total input energy, and A represents the area of ​​the welding interface; In each time window, the units in the weld interface area determine the transition between life and death units according to the changes in welding speed and energy input, based on the formula: Where i, j, and k represent the coordinate indexes of the Bezier unit in the x, y, and z directions, respectively. D(i, j, k) represents the life and death unit state of the Bezier unit with coordinates (i, j, k). represents the position vector of the Bezier cell with coordinates (i, j, k) in three-dimensional space, represents the welding speed vector, t m represents the mth time window, and m represents the index of the time window.

4. The finite element modeling method for simulating nickel-copper steel explosion welding according to claim 3 is characterized in that: The principle behind the generation of the local refined Bezier elements for the bond interface region for each time window is: The formulas for generating strain gradient and temperature gradient are: in, represents the strain gradient of the mth time window, ε xx (t m ) represents the strain tensor in the x direction caused by the displacement in the x direction in the mth time window, represents the strain tensor in the y direction caused by the displacement in the y direction in the mth time window, represents the strain tensor in the z direction caused by the displacement in the z direction in the mth time window, ε xy (t m ) represents the strain tensor in the x direction caused by the displacement in the y direction in the mth time window, ε yz (t m ) represents the strain tensor in the y direction caused by the displacement in the z direction in the mth time window, ε zx (t m ) represents the strain tensor in the z direction caused by the displacement in the x direction in the mth time window, where x, y, and z represent directions in three-dimensional space, respectively. represents the temperature gradient of the mth time window, T(t m ) represents the temperature of the mth time window; Generate a refinement factor and determine whether unit refinement is needed based on the refinement factor and the refinement threshold. The formula is: Among them, f(t m ) represents the refinement factor of the mth time window, ε thresh represents the strain gradient threshold, T thresh represents the temperature gradient threshold; When f(t m )>f thresh When the unit is refined, f thresh represents the refinement threshold; The control points are adjusted according to the Bezier difference function, and the Bezier element is refined based on the adjusted control points. The formula is: Among them, B(u, v, w) represents the Bezier body, u, v, w represent the parameters describing the three directions of constructing the Bezier body, p, q, r represent the index of the control point in the u, v, w direction, I, J, K represent the upper limit of the index in the u, v, w direction, P pqr represents the control points of the Bezier surface, represent the Bezier basis functions in the u, v, and w directions respectively; When the refinement factor is greater than the refinement threshold, the upper limit of the index in the u, v, and w directions is multiplied by the refinement factor, and the Bezier basis function is updated to generate the refined control points, and the Bezier unit is refined based on the control points.

5. The finite element modeling method for simulating nickel-copper steel explosion welding according to claim 1, characterized in that: The principle for generating the temperature field, pressure field and strain field corresponding to the discrete time series in step 3 is: The principle of generating temperature field is based on: The discretized form of the heat conduction equation is: Where C represents the heat capacity matrix, represents the temperature vector, represents the temperature change rate, R represents the thermal conductivity matrix, represents the internal heat source vector; The formula for generating the temperature field in discrete time based on the implicit backward Euler method is: Among them, T m+1 and T m They represent the temperature of the m+1th time window and the mth time window respectively, and Δt represents the length of the time window; The principle of generating pressure field is based on: The discretized form of the dynamic equation is: Where M represents the mass matrix, represents the acceleration vector, N represents the damping matrix, represents the velocity vector, G represents the stiffness matrix, represents the displacement vector, represents the external load vector; The formula for generating the pressure field in discrete time based on the central difference method is: in, Respectively represent the displacement vectors of the m+1th time window, the mth time window, and the m-1th time window, represents the external load vector of the mth time window; The formula for generating strain fields is: Among them, ε m+1 represents the strain vector of the m+1th time window, and B represents the strain-displacement matrix.

6. The finite element modeling method for simulating nickel-copper steel explosion welding according to claim 1, characterized in that: The principle for mapping the grid control point variables of the current time period to the next time period by Lagrange interpolation method in step 4 is: The Lagrange interpolation basis function is constructed based on the formula: S(X′,Y′,Z′)=S1·L1(X′,Y′,Z′)+S2·L2(X′,Y′,Z′)·S3·L3(X′,Y′,Z′) Among them, (X1, Y1, Z1), (X2, Y2, Z2), and (X3, Y3, Z3) represent the control points of the variables in the current time period, S1, S2, and S3 respectively represent the variable values ​​corresponding to the three control points, L1(X′, Y′, Z′), L2(X′, Y′, Z′), and L3(X′, Y′, Z′) are the Lagrangian basis functions of the interpolation points of the next time period corresponding to the control points of the three variables, and S(X′, Y′, Z′) represents the variable value of the interpolation point.

Images