Composite material contact stress calculation method and device, electronic equipment and storage medium

Abstract

The application discloses a composite material contact stress calculation method and device, electronic equipment and storage medium, and is used for solving the problems of low calculation efficiency and low solution accuracy in solving solid contact nonlinear problems in the prior art. The method comprises the following steps: obtaining a geometric model of an object to be analyzed, and performing mesh division on the geometric model based on a pre-established composite material contact pair to obtain a geometric model grid; according to the finite volume method and a preset boundary condition, the integral equation of the geometric model grid is discretized, and an iterative solver is selected to perform nonlinear iterative solution on the discrete result to obtain grid element displacement data of each interface in the geometric model grid; according to the grid element displacement data, the contact stress of the interface is calculated; when the friction force direction meets the preset judgment condition, the contact friction force of the interface is calculated according to the contact stress, so that the contact stress analysis of the object to be analyzed is quickly and accurately realized.

Description

Technical Field

[0001] This invention relates to the field of contact stress analysis technology, and in particular to a method, apparatus, electronic device, and storage medium for calculating contact stress of composite materials. Background Technology

[0002] Contact is prevalent in practical physical analysis engineering. Contact problems are also an important aspect of elasticity mechanics analysis, and their nonlinearity makes engineering analysis and research quite challenging. For nonlinear contact problems in solid mechanics, existing methods mainly rely on classical analytical methods or the finite element method. However, these methods may face challenges in terms of accuracy and stability when dealing with nonlinear contact problems.

[0003] Classical analytical methods based on Hertzian theory are mainly used to analyze contact forms such as cylindrical contact, spherical contact, and similar evolutionary forms. However, they cannot reflect the various complex contact patterns encountered in engineering practice, thus limiting the application range of Hertz's formula and resulting in inaccurate calculation results.

[0004] The finite element method (FEM) essentially transforms a continuum with infinitely many degrees of freedom into an aggregate with a finite number of degrees of freedom. It simplifies complex problems into simpler ones solvable using numerical methods. However, the FEM requires suitable type functions and has poor mesh adaptability. Furthermore, due to its discrete nature, it may exhibit numerical instability or low accuracy when dealing with nonlinear contact problems. Additionally, the FEM suffers from high computational complexity when handling large-scale contact problems, and its poor mesh adaptability for complex geometries leads to long computation times. Summary of the Invention

[0005] This invention provides a method, apparatus, electronic device, and storage medium for calculating contact stress in composite materials, which solves or partially solves the technical problems of low computational efficiency and low solution accuracy in solving nonlinear contact problems of solids in existing related technologies.

[0006] This invention provides a method for calculating contact stress in composite materials, the method comprising:

[0007] Obtain the geometric model of the object to be analyzed, and based on the pre-established composite material contact pairs, perform mesh generation on the geometric model to obtain the geometric model mesh;

[0008] Based on the finite volume method and preset boundary conditions, the geometric model mesh is discretized by integral equations, and an iterative solver is selected to perform nonlinear iterative solution on the discretized results to obtain the mesh element displacement data of each interface in the geometric model mesh.

[0009] Calculate the contact stress at the interface based on the displacement data of the mesh elements;

[0010] When the direction of the frictional force meets the preset judgment conditions, the contact frictional force at the interface is calculated based on the contact stress.

[0011] Optionally, the step of discretizing the geometric model mesh using integral equations based on the finite volume method and preset boundary conditions, and selecting an iterative solver to perform nonlinear iterative solutions on the discretized results to obtain the mesh element displacement data of each interface in the geometric model mesh, includes:

[0012] Based on the preset boundary conditions, the integral equation is discretized in each computational grid of the geometric model grid using the finite volume method to obtain the discretization coefficients. Then, based on the discretization coefficients, a system of discretized linear equations is constructed.

[0013] Select an iterative solver, input the matrix coefficients of each term of the discrete linear equation system into the iterative solver for nonlinear iterative solution, and stop iterating when the solution results meet the preset calculation requirements, and output the grid cell displacement data of each interface in the geometric model mesh.

[0014] Optionally, the mesh element displacement data includes the grid center displacement vectors and displacement gradients of the left and right elements of different materials on both sides of the interface. The step of calculating the contact stress at the interface based on the mesh element displacement data includes:

[0015] By iterating through all interfaces, and based on the displacement gradient, the grid center displacement vector gradient is interpolated onto the interface to obtain the relative displacement of the interface.

[0016] Calculate the face center displacement of the interface based on the relative displacement;

[0017] Based on the face center displacement, it is determined whether there is contact between the elements of the interface. If so, the contact stress of the interface is calculated using the penalty function method based on the face center displacement.

[0018] Optionally, the face-centered displacement includes the normal displacement vector of the interface, and determining whether there is contact between the elements of the interface based on the face-centered displacement includes:

[0019] Determine whether the normal displacement vector is greater than zero and exceeds a preset error value. If so, determine that there is contact between the units at the interface.

[0020] The step of calculating the contact stress at the interface using the penalty function method based on the face center displacement includes:

[0021] The contact stress at the interface is calculated using the penalty function method based on the normal displacement vector.

[0022] Optionally, the face-center displacement includes the tangential displacement vector of the interface, and the step of calculating the contact friction force of the interface based on the contact stress when the friction force direction meets the preset judgment condition includes:

[0023] If there is contact at the interface, determine whether the tangential displacement vector is not zero. If so, determine that there is friction between the units at the interface.

[0024] Based on the tangential displacement vector, calculate the theoretical traction force of the interface, and based on the contact stress, calculate the maximum static friction force of the interface.

[0025] If the theoretical traction force is less than the maximum static friction force, then the theoretical traction force is taken as the contact friction force of the interface.

[0026] If the theoretical traction force is greater than the maximum static friction force, then the interface dynamic friction force of the interface is calculated based on the contact stress, and the interface dynamic friction force is used as the contact friction force of the interface.

[0027] Optionally, the geometric model corresponds to multiple geometric sub-models, each geometric sub-model corresponding to a material distribution region. The step of meshing the geometric model based on pre-established composite material contact pairs to obtain a geometric model mesh includes:

[0028] Based on the pre-established composite material contact pairs, the meshing method of the geometric sub-model is determined according to the region where the material distribution area is located;

[0029] Based on the meshing method corresponding to each of the geometric sub-models, the computational mesh of the geometric model is divided to obtain the geometric model mesh of the object to be analyzed.

[0030] Optionally, the method further includes:

[0031] Obtain the object to be analyzed, which corresponds to multiple different material distribution regions;

[0032] A geometric sub-model is established for each of the aforementioned material distribution regions, and material property values ​​are assigned to each of the aforementioned geometric sub-models during the parameter initialization process based on the material properties of different material distribution regions.

[0033] Determine the interfaces between each pair of material distribution regions, and establish composite material contact pairs with penalty functions within the contact range that meet the preset contact conditions, based on the stress deformation of each geometric sub-model.

[0034] The individual geometric sub-models are assembled into the geometric model of the object to be analyzed.

[0035] The present invention also provides a composite material contact stress calculation device, comprising:

[0036] The geometric model meshing module is used to acquire the geometric model of the object to be analyzed, and to mesh the geometric model based on the pre-established composite material contact pairs to obtain the geometric model mesh.

[0037] The nonlinear iterative solution module is used to discretize the geometric model mesh using integral equations based on the finite volume method and preset boundary conditions, and select an iterative solver to perform nonlinear iterative solution on the discretized results to obtain the mesh element displacement data of each interface in the geometric model mesh.

[0038] The contact stress calculation module is used to calculate the contact stress at the interface based on the displacement data of the mesh elements.

[0039] The contact friction calculation module is used to calculate the contact friction force at the interface based on the contact stress when the direction of the friction force meets the preset judgment conditions.

[0040] The present invention also provides an electronic device, the device comprising a processor and a memory:

[0041] The memory is used to store program code and transmit the program code to the processor;

[0042] The processor is used to execute the composite material contact stress calculation method as described above, according to the instructions in the program code.

[0043] The present invention also provides a computer-readable storage medium for storing program code for executing the composite material contact stress calculation method as described in any of the preceding claims.

[0044] As can be seen from the above technical solutions, the present invention has the following advantages:

[0045] A method for calculating contact stress in composite materials is provided. First, a geometric model of the object to be analyzed is obtained. Based on pre-established composite material contact pairs, the geometric model is meshed to obtain a geometric model mesh. Next, the geometric model mesh is discretized using the finite volume method and preset boundary conditions through integral equations. An iterative solver is then selected to perform nonlinear iterative solutions on the discretized results to obtain the displacement data of each interface element in the geometric model mesh. Then, based on the mesh element displacement data, the contact stress at the interface is calculated. When the direction of the frictional force meets preset judgment conditions, the contact frictional force at the interface is calculated based on the contact stress, thereby achieving rapid and accurate contact stress analysis of the object to be analyzed. Attached Figure Description

[0046] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0047] Figure 1 A flowchart illustrating the steps involved in a method for calculating contact stress in composite materials;

[0048] Figure 2 This is a schematic diagram of grid numbering and generation matrix;

[0049] Figure 3 This is a schematic diagram of the process for calculating contact stress in a composite material.

[0050] Figure 4 This is a schematic diagram of a two-dimensional nozzle model;

[0051] Figure 5 This is a schematic diagram of the bottom of an ablation layer (inner wall of the nozzle);

[0052] Figure 6 This is a schematic diagram of the free boundary of a nozzle.

[0053] Figure 7 This is a schematic diagram of a fixed boundary of a fixed structural layer;

[0054] Figure 8 This is a schematic diagram of contact stress distribution;

[0055] Figure 9 This is a structural block diagram of a composite material contact stress calculation device. Detailed Implementation

[0056] This invention provides a method, apparatus, electronic device, and storage medium for calculating contact stress in composite materials, which solves or partially solves the technical problems of low computational efficiency and low solution accuracy in solving nonlinear contact problems of solids in existing related technologies.

[0057] To make the objectives, features, and advantages of this invention more apparent and understandable, the technical solutions of the embodiments of this invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the embodiments described below are only some embodiments of this invention, and not all embodiments. Based on the embodiments of this invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this invention.

[0058] As an example, existing methods for solving nonlinear contact problems in solid mechanics are mainly based on classical analytical methods or the finite element method. However, these methods may face challenges in terms of accuracy and stability when dealing with nonlinear contact problems.

[0059] Classical analytical methods based on Hertzian theory are mainly used to analyze contact forms such as cylindrical contact, spherical contact, and similar evolutionary forms. However, they cannot reflect the various complex contact patterns encountered in engineering practice, thus limiting the application range of Hertz's formula and resulting in inaccurate calculation results.

[0060] For complex multi-body contacts, the interactions between the entities are quite intricate. These contact problems are difficult to handle using Hertz's formula and can only be solved numerically using the finite element method (FEM). Currently, there are two types of FEM models: one is the penalized spring method, which simulates contact behavior by setting normal and tangential springs on the contact interface; the other is the Lagrange method, such as the Lagrange multiplier method and the modified Lagrange method, which mainly uses Lagrange multipliers to indirectly obtain the contact forces between objects to simulate contact behavior.

[0061] When using relevant models for finite element method (FEM) solutions, the FEM essentially transforms a continuum with infinitely many degrees of freedom into an aggregate with a finite number of degrees of freedom. It simplifies complex problems into simpler ones solvable using numerical methods. However, the FEM requires setting appropriate type functions and has poor mesh adaptability. Furthermore, due to its discrete nature, it may exhibit numerical instability or low accuracy when dealing with nonlinear contact problems. Additionally, the FEM suffers from high computational complexity when handling large-scale contact problems, and its poor adaptability to complex geometries leads to long computation times.

[0062] Previous research on composite materials has primarily focused on their thermodynamic properties, including extensive studies ranging from empirical to statistical constitutive models. However, detailed research and discussion on the stresses in different contact zones of the structural surfaces have been lacking.

[0063] Under normal stress, the actual contact area between structural surfaces often occupies only a small fraction of the nominal contact area. The actual contact stress within this area differs significantly from the nominal average contact stress of the structural surface. In practice, the actual contact stress in the contact area of ​​the structural surface is the direct cause of micro-protrusion failure during normal compression or shear slip. Therefore, assessing the actual contact stress on structural surfaces is fundamental to predicting their compressive or shear behavior.

[0064] Therefore, one of the core inventive points of this invention is as follows: to solve the current problem of contact nonlinear calculation with complex geometry, irregular contact area shape, uneven distribution of internal spatial properties of composite materials, and potential anisotropy and large differences in properties, the lattice-centered finite volume method is extended from fluid mechanics calculation to solid mechanics calculation. A contact nonlinear algorithm based on the lattice-centered finite volume method is proposed to calculate the contact stress between material interfaces. By fully considering the material nonlinearity, contact nonlinearity, and geometric nonlinearity of composite materials under actual engineering conditions, considering the influence of uneven local load distribution on the stress and strain distribution of composite materials, and performing contact analysis on the three-dimensional solid of the load transfer process, the problems of classical analytical methods, such as inability to reflect various complex contact forms in actual engineering, limited application scope, and inaccurate calculation results, can be avoided. By discretizing and solving the elastic equation using the finite volume method, the displacement increment of the contact interface is accurately calculated, and then the contact stress is calculated. The direction of friction force is also determined, and finally the contact friction force is calculated, thereby improving the efficiency and accuracy of contact stress calculation for composite materials. Compared with the finite element method, the finite volume method is more adaptable to unstructured meshes. By adopting the technical solution provided in the embodiments of the present invention, numerical methods can be provided for stress-strain analysis of composite materials, reducing development costs; at the same time, in terms of contact stress prediction, it can be better integrated with actual engineering, predicting areas of stress concentration caused by complex geometric structures and complex contact areas.

[0065] To facilitate a better understanding by those skilled in the art, the relevant technical features or means involved in the technical solution of this invention will be explained and described first.

[0066] The Finite Volume Method (FVM), also known as the Control Volume Method (CVM), is a numerical method for solving partial differential equations. Its working principle involves first dividing the computational domain into a series of non-overlapping control volumes (corresponding sub-regions, i.e., computational grids), with each grid point surrounded by a unique control volume. Then, the differential equation to be solved is integrated over each control volume to obtain a discrete equation. In deriving the discrete equation, certain assumptions need to be made regarding the distribution of the function being solved and its derivatives on the interface. Discrete equations derived using the finite volume method guarantee conservation properties, and the coefficients of the discrete equations have clear physical meanings, resulting in relatively low computational complexity. In practical applications, nodes can be considered as representatives of control volumes; during the discretization process, physical quantities on a control volume are defined and stored at that node.

[0067] After discretizing the region using the finite volume method, the following four geometric elements can be obtained:

[0068] Node: The geometric location of the unknown physical quantity that needs to be solved.

[0069] Control volume: The smallest geometric unit to which the governing equations or conservation laws are applied.

[0070] Interface (face): Defines the location of the interface of the control volume corresponding to each node. Among them, the interface refers to the boundary between two regions or entities, and these two regions must have at least one face that overlaps or partially overlaps.

[0071] Grid line: A cluster of curves that connect two adjacent nodes.

[0072] In current applications, the finite volume method is mainly used to solve fluid flow and heat transfer problems.

[0073] The lattice-centered finite volume method stores the variables to be solved at the center of discrete elements, which are the control volumes. The reason for applying the lattice-centered finite volume method to analyze solid and structural problems in this embodiment of the invention is threefold: first, the finite volume method strictly adheres to physical conservation laws; second, the finite volume method is applicable to complex geometries and has good mesh adaptability; and third, it can be developed within the same program framework as existing finite volume methods used for fluid dynamics calculations, and extended to multiphysics coupling. Through the lattice-centered finite volume method, the computational domain can be divided into discrete control volumes, and conservation equations can be established on each control volume. Compared to the finite element method, the lattice-centered finite volume method is more direct and accurate in handling conservation equations.

[0074] Reference Figure 1The diagram illustrates a flowchart of a method for calculating contact stress in composite materials according to an embodiment of the present invention, which may specifically include the following steps:

[0075] Step 101: Obtain the geometric model of the object to be analyzed, and based on the pre-established composite material contact pairs, perform mesh generation on the geometric model to obtain the geometric model mesh;

[0076] In practical applications, before conducting formal contact stress analysis based on the lattice-centered finite volume method, a corresponding geometric model of the object to be analyzed can be constructed first. Specifically, the construction of the geometric model of the object to be analyzed can be achieved through the following steps S01 to S04:

[0077] Step S01: Obtain the object to be analyzed, which corresponds to multiple different material distribution regions;

[0078] Step S02: Establish a geometric sub-model for each material distribution region, and assign material property values ​​to each geometric sub-model based on the material properties of different material distribution regions during the parameter initialization process;

[0079] CAD (Computer-Aided Design) software, such as AutoCAD (Autodesk Computer-Aided Design) and UG (Unigraphics NX, an interactive CAD / CAM system), can be used to create a geometric sub-model of each material distribution area based on the object to be analyzed.

[0080] Next, using general-purpose finite element analysis software such as ABAQUS and ANSYS, or self-developed software, a stress-strain analysis model is constructed based on multiple established geometric sub-models. Specifically, by fully considering the changes in the mechanical properties of composite materials with temperature and the non-uniform spatial distribution, characteristic curves of different materials are obtained, and material properties are set in the software. During the material parameter initialization process, material property values ​​are assigned to the geometric sub-models corresponding to different regions.

[0081] Different materials have different property curves, so it is necessary to determine the anisotropy of material properties and the variation of properties with temperature before calculation. For example, based on a known table of material property changes with temperature (e.g., by consulting a property library), a property-temperature relationship curve is plotted, and a formula for the relationship between property and temperature is fitted. The property distribution of the object to be analyzed is determined by the temperature distribution.

[0082] Step S03: Determine the interface between each pair of material distribution areas, and establish composite material contact pairs with penalty functions within the contact range that meets the preset contact conditions, based on the stress deformation of each geometric sub-model.

[0083] Based on the material distribution areas of the object to be analyzed, such as the material distribution within laminated materials (layers of different materials pressed together) or the material distribution in each region of a composite nozzle, the interfaces between different materials are clearly defined. For example, based on actual needs, partitions of different materials within the nozzle to be analyzed can be defined.

[0084] Based on the stress conditions of each geometric sub-model, the deformation trend (structural stress state) of the model is predicted, and composite material contact pairs are established for all possible contact ranges in combination with preset boundary conditions. The contact types of these contact pairs are either frictional contact pairs with penalty functions or frictionless contact pairs.

[0085] The principles for establishing contact pairs in composite materials (i.e., determining whether the contact conditions are met) are as follows: 1. After the materials are subjected to force, there may be several tendencies between different materials, such as mutual compression, separation, or maintaining their original state without interaction; 2. The boundaries that may have mutual compression tendencies are taken as contact pairs; 3. Contact can only occur between a certain material partition of the object to be analyzed and partitions of other different material types (such as between material A and material B of the composite material to be analyzed). Unconstrained boundaries will not make contact.

[0086] In real-world scenarios, consider two horizontally placed blocks. If they are stretched to the sides, they will not make contact; however, if they are compressed inward, they may make contact. Therefore, the "possible contact range" can be understood as the range of contact surfaces that tend to deform after the object being analyzed is subjected to force.

[0087] Based on the material properties of the object to be analyzed and the composite material contact pairs established in the preceding steps, the geometric nonlinear conditions limiting the stress solution process can be determined.

[0088] Geometric nonlinearity refers to the influence of force on the contact area, resulting in large strain due to geometric deformation, and whether the displacement of mesh elements is negligible compared to the mesh size itself. It can also be understood as whether the deformation of the object under analysis after being subjected to force is negligible compared to its own size. If it is not negligible (exceeding the range of elastic deformation), then geometric nonlinearity needs to be considered in the calculation. If the deformation is negligible compared to its own size, then geometric nonlinearity does not need to be considered in the solution calculation.

[0089] Step S04: Assemble the various geometric sub-models into a geometric model of the object to be analyzed.

[0090] Finally, the individual geometric sub-models can be assembled into a complete geometric model for computational analysis.

[0091] Next, the solid model of the object to be analyzed can be meshed. As discussed above, a geometric model can correspond to multiple geometric sub-models, and each geometric sub-model can correspond to a material distribution region. Therefore, in the specific implementation, the process of meshing the geometric model based on pre-established composite material contact pairs to obtain the geometric model mesh can include the following sub-steps S11 to S12:

[0092] Step S11: Based on the pre-established composite material contact pairs, determine the meshing method of the geometric sub-model according to the region where the material distribution area is located;

[0093] Import the geometric model into meshing software such as ICEMCD (The Integrated Computer Engineering and Manufacturing code for Computational Fluid Dynamics) or Hypermesh. Based on the location of different material distribution areas and the composite material contact pairs established in the previous steps, select the appropriate meshing method for each region (geometric sub-model) to ensure the quality of the contact pair mesh and the matching of the mesh face center and nodes.

[0094] Step S12: Divide the computational mesh of the geometric model according to the meshing method corresponding to each geometric sub-model, and obtain the geometric model mesh of the object to be analyzed.

[0095] Next, the computational grid of the geometric model corresponding to the object to be analyzed can be divided according to the meshing method. At the same time, during the meshing process, poor quality meshes can be deleted or modified according to actual needs to obtain the final geometric model mesh.

[0096] Step 102: Discretize the geometric model mesh using integral equations based on the finite volume method and preset boundary conditions, and select an iterative solver to perform nonlinear iterative solution on the discretized results to obtain the mesh element displacement data of each interface in the geometric model mesh.

[0097] Next, it is necessary to solve the displacement field and stress-strain field of the geometric model mesh of the object to be analyzed, and select an iterative solver to perform nonlinear iterative solution in order to obtain the displacement data of the mesh element at each interface in the geometric model mesh.

[0098] The solution principle is as follows: based on the given boundary conditions, the finite volume method is used to discretize the data and obtain the discretization coefficients of each computational grid element. Then, based on the discretization coefficients, a system of equations is constructed, and combined with nonlinear coupling effects (contact nonlinearity and material nonlinearity), an iterative solver is selected for iterative solution.

[0099] Material nonlinearity refers to the changes in material properties as loads are applied and damage occurs, depending on the calculation process and the load applied. Examples include stress exceeding the elastic range, irreversible damage to the material, and the temperature dependence of the material's properties.

[0100] Contact nonlinearity refers to the possibility that the magnitude and direction of contact stress in different contact areas may change over time and with the application of loads. For example, materials may debond and remain in a detached state without contact, or materials may transition from a detached state to a contact state, and the contact force distribution may be uneven.

[0101] Specifically, step 102 can be implemented through the following sub-steps S21 to S22:

[0102] Step S21: Based on the preset boundary conditions, the integral equation is discretized in each computational grid of the geometric model grid using the finite volume method to obtain the discretization coefficients. Then, based on the discretization coefficients, a system of discretized linear equations is constructed.

[0103] The preceding steps involved meshing the geometric model, resulting in a geometric model mesh. This process effectively involved meshing the computational domain of the object being analyzed into a finite volume mesh. To facilitate computation, each computational mesh cell can be numbered after meshing.

[0104] Next, based on the finite volume method, the integral equations are discretized within each computational grid cell. The discretization results are stored in a matrix with a dimension equal to the number of grid cells to construct a system of linear equations offline. For better understanding by those skilled in the art, exemplarily, Figure 2 A schematic diagram of the mesh numbering and generation matrix of a current geometric model mesh is shown.

[0105] For the generated matrix, the diagonal coefficients represent the influence on the cell numbered in the current row, and the off-diagonal coefficients (non-zero elements) represent the influence on adjacent cells. Combined with... Figure 2 Taking the third computational grid cell as an example, there is a correspondence between the computational grid cell numbered 3 and the corresponding position in the generator matrix. The influence on this computational grid cell is stored on the diagonal of the third row of the matrix, and the positions of the other non-zero elements in this row are located in the first, second, and fourth columns, respectively. This means that cells numbered 1, 2, and 4 are affected by cell number 3.

[0106] Step S22: Select the iterative solver, input the matrix coefficients of each term of the discretized linear equation system into the iterative solver for nonlinear iterative solution, and stop the iteration when the solution results meet the preset calculation requirements, and output the grid cell displacement data of each interface in the geometric model mesh.

[0107] Matrix solving involves inputting the coefficients of each term in the matrix into a matrix solver (iterative solver), calling the solver to perform repeated iterations until the final solution meets the requirements for stopping the calculation, such as reaching the maximum number of iterations, or the calculation accuracy meeting the preset accuracy requirements, such as being less than a certain preset error value, at which point the iteration stops and the displacement data of each interface in the geometric model mesh are output.

[0108] Step 103: Calculate the contact stress at the interface based on the displacement data of the mesh elements;

[0109] Then, based on the displacement data of the mesh elements, the contact stress at the interface can be calculated using the penalty function method.

[0110] The calculation principle of the penalty function method is as follows: At each step of the calculation, it is necessary to first check whether each slave node has penetrated the master surface. If no penetration has occurred, no action is taken on that slave node; if penetration has occurred, a large interfacial contact force is applied between the slave node and the penetrated master surface. The magnitude of the contact force is proportional to the penetration depth and the stiffness of the master surface. The value of this interfacial contact force is called the penalty function value. Its physical meaning can be understood as placing a normal spring between the slave node and the penetrated master surface to limit the penetration between the slave node and the master surface. The processing procedure for the master node is the same.

[0111] Specifically, after nonlinear iterative solution, the mesh element displacement data can include the grid center displacement vectors and displacement gradients of the left and right elements of different materials on both sides of the interface. The specific process of calculating the contact stress at the interface based on the mesh element displacement data can be achieved through the following sub-steps S31 to S33:

[0112] Step S31: Loop through all interfaces and, based on the displacement gradient, interpolate the grid center displacement vector gradient onto the interfaces to obtain the relative displacement of the interfaces;

[0113] The first step is to calculate the relative displacement of the interface.

[0114] Based on the grid center displacement vectors and displacement gradients of the left and right cells of different materials at the interface, all interfaces are iterated, and the displacement vectors of the left and right cell centers are interpolated onto the interface using gradient interpolation to obtain the relative displacement of the interface. The specific calculation formula is as follows:

[0115] U Lf =U L +(▽U) L ·d Lf

[0116] U Rf =U R +(▽U) R ·d Rf

[0117] U f =U Lf -U Rf

[0118] Where, d Lf and d Rf These represent the displacement vectors from the center of the left and right cells to the cell interface, respectively; (▽U) L and (▽U) R U represents the displacement gradient of the left and right elements, respectively; Lf The displacement of the center of the left cell to the contact interface is interpolated through the displacement gradient and distance vector of the left cell; U Rf This refers to the displacement amount interpolated from the displacement gradient and distance vector of the right-side element to the contact interface; U f The relative displacement of the interface refers to the displacement difference.

[0119] Step S32: Calculate the face center displacement of the interface based on the relative displacement;

[0120] Next, the face center displacement (including normal and tangential directions) of the interface can be calculated using the following formula:

[0121] U fn =(U f ·n f )·n f

[0122] U fs =U f -U fn

[0123] Among them, U fn U is the normal displacement vector of the interface; fs n is the tangential displacement vector of the interface; f Let be the unit surface normal vector of the element interface.

[0124] Step S33: Based on the face center displacement, determine whether there is contact between the elements at the interface. If so, calculate the contact stress at the interface using the penalty function method based on the face center displacement.

[0125] After calculating the interface displacement U f U fn U fs After that, the direction of displacement and the friction force can be determined, and then the contact stress and contact friction force can be calculated based on the determination results.

[0126] Based on step S32, it can be seen that the face-center displacement can include the normal displacement vector of the interface. Therefore, based on the face-center displacement, it can be determined whether there is contact between the elements of the interface. Specifically, it can be determined whether the normal displacement vector is greater than zero (U).fn If the error value is greater than 0 and exceeds a preset error value (such as a pre-defined error value), then it is determined that there is contact between the elements at the interface. Otherwise, there is no contact between the elements at the interface.

[0127] Furthermore, based on the face-center displacement, the contact stress at the interface is calculated using the penalty function method. Specifically, based on the normal displacement vector, the contact stress F (also known as normal stress) at the interface is calculated using the penalty function method.

[0128] F = kU fn

[0129]

[0130] Where k is the penalty function coefficient; f scale is the scaling factor for the penalty function, which is generally set to a number less than 0.1; K is the bulk modulus of the material, A is the interface area, and V is the volume of the computational mesh element.

[0131] Step 104: When the direction of frictional force meets the preset judgment conditions, calculate the contact frictional force of the interface based on the contact stress.

[0132] The face-centered displacement also includes the tangential displacement vector of the interface. Therefore, the calculation process in step 104 can be achieved through the following steps S41 to S43:

[0133] Step S41: If there is contact at the interface, determine whether the tangential displacement vector is not zero. If so, determine that there is friction between the elements at the interface.

[0134] If the interface is in contact and U fs If the value is not equal to 0, then there is a tendency for tangential relative displacement, meaning there is friction between the elements at the interface. In this case, contact friction can be further calculated. Conversely, if the value is not equal to 0, then there is no friction between the elements at the interface.

[0135] Step S42: Calculate the theoretical traction force of the interface based on the tangential displacement vector, and calculate the maximum static friction force of the interface based on the contact stress.

[0136] It is understandable that contact stress exists at the interface, which is a prerequisite for the existence of contact friction. Further analysis suggests that there may be a relative displacement tendency between the interfaces, but if the traction force of the displacement is less than the maximum static friction force, static friction occurs between the interfaces, meaning there is no tangential displacement between the elements, and the magnitude of the static friction force is proportional to the traction (relative displacement).

[0137] In practical calculations, a calculation coefficient ε can be given to determine the relationship between frictional force and tangential displacement. TCalculate the theoretical friction force and compare it with the maximum static friction force to determine whether sliding friction has occurred.

[0138] Specifically, the theoretical traction force f F It can be calculated using the following formula:

[0139]

[0140] f F =ε T |U fs |

[0141] Among them, f fric is the scaling factor for friction calculation, which is generally set to a number less than 0.1; G is the shear modulus of the material; A is the interface area; V is the volume of the computational mesh element.

[0142] Considering that the maximum static friction is greater than the sliding friction, a relative coefficient is set. This is used to calculate the relationship between sliding friction and maximum static friction. The coefficient varies depending on the material; for example, it can be set to 1.05. Maximum static friction f max The calculation formula is as follows:

[0143]

[0144] μ is the coefficient of kinetic friction between materials.

[0145] Step S43: If the theoretical traction force is less than the maximum static friction force, then the theoretical traction force is used as the contact friction force of the interface; if the theoretical traction force is greater than the maximum static friction force, then the interface dynamic friction force is calculated based on the contact stress, and the interface dynamic friction force is used as the contact friction force of the interface.

[0146] If f F <f max If the theoretical traction force is less than the maximum static friction force, static friction occurs, and the magnitude of the friction force is proportional to the tangential relative displacement. At this time, the contact friction force at the interface is f = f F .

[0147] If f F >f max If the theoretical traction force is greater than the maximum static friction force, then kinetic friction occurs. The magnitude of the friction force is proportional to the normal stress, and the interface friction force is f = f d .

[0148] Among them, the interfacial dynamic friction force f d Calculated using the following formula:

[0149] f d =μ|F|

[0150] The direction of the contact friction force at the interface is parallel to U. fs same.

[0151] To better illustrate, Figure 3 A schematic diagram of the general process for calculating contact stress in composite materials is shown.

[0152] First, the displacement data of the mesh elements at each interface in the geometric model mesh are obtained by combining the finite volume method with nonlinear iterative solution.

[0153] Mark the master and slave surfaces on both sides of the interface, and calculate the relative displacement of the interface using gradient interpolation based on the mesh element displacement data, and then calculate the displacement of the face center.

[0154] Based on the calculation results, determine whether the interface is debonded or penetrated. If it is debonded, the interface is not in contact; if it is penetrated, the penalty function method is used to calculate the contact stress of the interface.

[0155] Determine whether there is a tangential relative displacement trend among the elements at the interface. If so, calculate the contact friction force based on the contact stress; otherwise, there is no friction force at the interface.

[0156] Based on the contact analysis results, the interface of the geometric model of the object to be analyzed is marked.

[0157] Furthermore, different computational requirements and solution precision can be set according to varying computational capabilities. For example, it can be configured whether to consider the deformation of the structural surface matrix during the calculation process, i.e., whether to consider mesh deformation, or whether to incorporate frictional forces into the calculation. Specifically, different computational methods can be selected based on the number of meshes in the actual computational object (related to computation time) and the computation time itself, thereby improving work efficiency.

[0158] The actual engineering calculation accuracy requirements can be reflected in the following aspects: (1) Set the calculation convergence residual, analyze the calculation results of different residual settings, and select appropriate residual settings; (2) Set the scaling factor of contact stress and friction force. The larger the scaling factor, the faster the convergence, but the easier it is to cause the calculation residual to oscillate. The smaller the scaling factor, the slower the calculation convergence, but it will not cause the calculation residual to oscillate; (3) In the calculation process, if it is necessary to consider the mesh deformation, the deformation of each mesh needs to be calculated after each step of the calculation when using the traditional method. However, by using the method provided by the embodiment of the present invention, it is possible to fully adapt to the calculation of the deformed mesh, that is, to perform adaptive calculation on the deformed mesh without the need for additional manual modification.

[0159] In this invention, the lattice-centered finite volume method is extended from fluid mechanics calculations to solid mechanics calculations. A contact nonlinear algorithm based on the lattice-centered finite volume method is proposed for calculating contact stress between material interfaces. By fully considering the material nonlinearity, contact nonlinearity, and geometric nonlinearity of composite materials under actual engineering conditions, and considering the influence of uneven local load distribution on the stress-strain distribution of composite materials, and performing contact analysis on the three-dimensional solid during load transfer, the classical analytical method avoids the problems of being unable to reflect various complex contact patterns in actual engineering, having a limited application scope, and inaccurate calculation results. By discretizing and solving the elastic equation using the finite volume method, the displacement increment of the contact interface is accurately calculated, and then the contact stress is calculated. The direction of friction force is also determined, and finally the contact friction force is calculated, improving the efficiency and accuracy of contact stress calculation for composite materials. Compared with the finite element method, the finite volume method is more adaptable to unstructured meshes. By adopting the technical solution provided in this invention, a numerical method for stress-strain analysis of composite materials can be provided, reducing development costs. At the same time, in terms of contact stress prediction, it can be better integrated with actual engineering, predicting areas of stress concentration caused by complex geometric structures and complex contact regions.

[0160] For ease of understanding, the following description uses a specific example to illustrate an embodiment of the present invention.

[0161] Taking a nozzle made of a certain composite material as an example, the ablation layer of a typical nozzle is composed of carbon / phenolic (CP) and carbon / carbon (CC) materials. The insulation layer is mainly composed of high silica / phenolic (HSC-C) materials. The fixing structure layer is made of steel. For example, Table 1 gives the material properties of carbon / carbon, carbon / phenolic, and high silica / phenolic materials:

[0162] CC 3273 5.79 0.1 CP 1500 5.52 0.05 HSC-C 1700 5.38 0.05

[0163] Table 1: Material Performance Parameters

[0164] Where ρ represents density, with units of kg / m³. 3 E represents Young's modulus, in GPa; λ represents Poisson's ratio.

[0165] (I) Geometry Establishment and Mesh Generation

[0166] Figure 4 This is a simplified structural diagram of a two-dimensional model of a rocket nozzle used in this example.

[0167] Combination Figure 4The nozzle features a typical three-layer structure, consisting of an ablation-resistant layer (green grid in the color image, corresponding to a gray background in the grayscale image), an insulation layer (blue grid in the color image, corresponding to a black background in the grayscale image), and a fixing structure layer (red grid in the color image, corresponding to a white background in the grayscale image). The ablation-resistant layer withstands the pressure and erosion of the combustion gases; the insulation layer, with its low thermal conductivity, blocks heat transfer; and the fixing structure connects to external structures, secures the nozzle, prevents rigid displacement, and transmits thrust.

[0168] Regarding the establishment of contact pairs in composite materials, Figure 4 For example, this nozzle contains three materials. Based on the known internal material distribution of the nozzle (i.e., the object to be analyzed), the distribution area of ​​each material and its initial properties are determined. This nozzle can only have two contact pairs: green and blue, and blue and red. Furthermore, the prerequisite for contact between the material distribution areas is that these three materials have a tendency to compress each other, not a tendency to separate. If an upward force is applied to the top and a downward force to the bottom, causing the materials to separate, then contact is impossible.

[0169] (II) Structural Boundary Conditions

[0170] Figure 5 A schematic diagram of the bottom of the ablation layer (inner wall surface of the nozzle) is shown.

[0171] Among them, a pressure load of 4MPa is applied to the inner wall surface of the nozzle (e.g. Figure 5 The yellow border in the color map (corresponding to) Figure 5 The bolded boundary in the grayscale image is set as the second type of boundary condition for the statutory pressure.

[0172] The nozzle assemblies are designed to be in rigid contact. The contact stress calculation method proposed in this embodiment of the invention is used to discretize and calculate the contact stress between the materials, which are respectively... Figure 4 Contact between the medium ablation resistant layer (green in color image / grayscale image) and the insulation layer (blue in color image / black in grayscale image), and between the insulation layer (blue in color image / black in grayscale image) and the fixed structure layer (red in color image / white in grayscale image).

[0173] Figure 6 This is a schematic diagram of the free boundary of the nozzle. The longitudinal interface and other boundary conditions of the nozzle are set. Figure 6 The white boundaries in the color diagram are all considered frictionless constraints (corresponding to...) Figure 6 Apply free boundary conditions to the bolded boundary in the grayscale image.

[0174] Figure 7This is a schematic diagram of the fixed boundary of the nozzle's fixed structural layer. The top fixed structural layer of the nozzle is typically connected to the combustion chamber shell via bolts and flanges. For the fixed or displacement constraint boundary of the nozzle flange, this example uses [method / method] during simulation. Figure 7 The light blue boundary in the color image is subject to a fixed constraint (corresponding to...). Figure 7 (The bolded boundary in the grayscale image) is used as a fixed constraint at the connection between the flange and the combustion chamber shell, with a displacement of 0.

[0175] (III) After the calculation is completed, output the contact stress at the interface and statistically analyze the average displacement and load.

[0176] Figure 8 This is a schematic diagram of the contact stress distribution of the nozzle. Figure 8 In the calculation results, ppface represents the contact stress, and the numerical value represents the magnitude of the contact stress.

[0177] To facilitate comparison, contact stress calculations were performed on the same object to be analyzed (i.e., the nozzle) using a commercial software, an open-source software, and the method provided in this embodiment of the invention. The comparison of the calculation results is shown in Tables 2 and 3.

[0178] ablation layer and insulation layer 7.2 7.0 7.0 Insulation layer and fixing layer 0.32 0.38 0.36

[0179] Table 2: Comparison of Contact Stress Results

[0180]

[0181] Table 3: Calculation Time by Different Methods

[0182] The calculation method provided in this invention can directly handle nonlinear problems such as geometric nonlinearity and contact nonlinearity without additional modifications, since the solution process has already iterated on the equation system. During the calculation, no inter-surface interpolation is required; only the normal contact force and tangential friction force need to be calculated locally on the two contact surfaces. Comparison with commercial software calculation results further verifies the accuracy of the method provided in this invention. Simultaneously, the maximum contact pressure shows good consistency, with small differences related to the Hertzian assumption of the loading force and the loading displacement. As shown in Table 3, although the single-core calculation time of this invention differs from that of commercial software, the method provided in this invention is adapted to multi-core parallel computing. By increasing the number of parallel cores, its computational efficiency is comparable to that of commercial software.

[0183] Based on this example, it can be concluded that the method for calculating the contact force of composite materials provided in this invention fully considers material nonlinearity, contact nonlinearity, and geometric nonlinearity. By performing contact analysis on the three-dimensional solid of the composite material during the load transfer process, the stress-strain distribution and contact stress distribution at the interface can be accurately calculated, thereby realistically simulating material properties and accurately obtaining the load magnitude under real-world conditions. Furthermore, using the technical solution of this invention for composite material contact stress analysis also offers advantages such as high accuracy, shortened development cycle, and reduced costs. It can be developed within the same program framework for fluid computation based on the finite volume method and extended to multiphysics coupling.

[0184] Reference Figure 9 The diagram illustrates a structural block diagram of a composite material contact stress calculation device provided in an embodiment of the present invention, which may specifically include:

[0185] The geometric model mesh generation module 901 is used to acquire the geometric model of the object to be analyzed, and to perform mesh generation on the geometric model based on the pre-established composite material contact pairs to obtain the geometric model mesh.

[0186] The nonlinear iterative solution module 902 is used to discretize the geometric model mesh by integral equations according to the finite volume method and preset boundary conditions, and select an iterative solver to perform nonlinear iterative solution on the discretized results to obtain the mesh element displacement data of each interface in the geometric model mesh.

[0187] The contact stress calculation module 903 is used to calculate the contact stress at the interface based on the displacement data of the mesh elements.

[0188] The contact friction force calculation module 904 is used to calculate the contact friction force of the interface based on the contact stress when the direction of the friction force meets the preset judgment conditions.

[0189] In one optional embodiment, the nonlinear iterative solution module 902 includes:

[0190] The finite volume method discretization module is used to discretize the integral equation in each computational grid of the geometric model grid based on preset boundary conditions, obtain the discretization coefficients, and then construct the discretized linear equation system based on the discretization coefficients.

[0191] The nonlinear iterative solution submodule is used to select an iterative solver, input the matrix coefficients of each term of the discrete linear equation system into the iterative solver for nonlinear iterative solution, and stop iterating when the solution result meets the preset calculation requirements, and output the grid cell displacement data of each interface in the geometric model mesh.

[0192] In one optional embodiment, the mesh element displacement data includes the grid center displacement vectors and displacement gradients of the left and right elements of different materials on both sides of the interface, and the contact stress calculation module 903 includes:

[0193] The gradient interpolation module is used to iterate through all interfaces and interpolate the grid center displacement vector gradient onto the interface based on the displacement gradient to obtain the relative displacement of the interface.

[0194] The face-center displacement calculation module is used to calculate the face-center displacement of the interface based on the relative displacement.

[0195] The contact stress calculation submodule is used to determine whether there is contact between the elements of the interface based on the face center displacement. If there is contact, the contact stress of the interface is calculated using the penalty function method based on the face center displacement.

[0196] In one optional embodiment, the face-centered displacement includes the normal displacement vector of the interface, and the contact stress calculation submodule is specifically used for:

[0197] Determine whether the normal displacement vector is greater than zero and exceeds a preset error value. If so, determine that there is contact between the units at the interface.

[0198] The contact stress at the interface is calculated using the penalty function method based on the normal displacement vector.

[0199] In one optional embodiment, the face-centered displacement includes the tangential displacement vector of the interface, and the contact friction force calculation module 904 includes:

[0200] The friction force determination module is used to determine whether the tangential displacement vector is not zero when there is contact at the interface. If so, it is determined that there is friction force between the units of the interface.

[0201] The friction force related calculation module is used to calculate the theoretical traction force of the interface based on the tangential displacement vector, and to calculate the maximum static friction force of the interface based on the contact stress.

[0202] The contact friction force determination module is used to take the theoretical traction force as the contact friction force of the interface when the theoretical traction force is less than the maximum static friction force.

[0203] The interface dynamic friction calculation module is used to calculate the interface dynamic friction force of the interface based on the contact stress when the theoretical traction force is greater than the maximum static friction force, and to use the interface dynamic friction force as the contact friction force of the interface.

[0204] In one optional embodiment, the geometric model corresponds to multiple geometric sub-models, each geometric sub-model corresponds to a material distribution region, and the geometric model mesh generation module 901 includes:

[0205] The mesh generation method determination module is used to determine the mesh generation method of the geometric sub-model based on the pre-established composite material contact pairs and the region where the material distribution area is located.

[0206] The computational mesh generation module is used to divide the computational mesh of the geometric model according to the mesh generation method corresponding to each of the geometric sub-models, so as to obtain the geometric model mesh of the object to be analyzed.

[0207] In one alternative embodiment, the device further includes:

[0208] The object to be analyzed module is used to acquire the object to be analyzed, which corresponds to multiple different material distribution regions.

[0209] The material property assignment module is used to establish a geometric sub-model for each of the material distribution regions, and to assign material property values ​​to each of the geometric sub-models during the parameter initialization process based on the material properties of different material distribution regions.

[0210] The composite material contact pair establishment module is used to determine the interface between each pair of material distribution areas, and establish composite material contact pairs with penalty functions within the contact range that meets the preset contact conditions based on the stress deformation of each geometric sub-model.

[0211] The geometric model assembly module is used to assemble the various geometric sub-models into a geometric model of the object to be analyzed.

[0212] As the device embodiment is basically similar to the method embodiment, it is described in a relatively simple way. For relevant details, please refer to the description of the method embodiment above.

[0213] This invention also provides an electronic device, which includes a processor and a memory:

[0214] The memory is used to store program code and transfer the program code to the processor;

[0215] The processor is used to execute the composite material contact stress calculation method of any embodiment of the present invention according to the instructions in the program code.

[0216] This invention also provides a computer-readable storage medium for storing program code for executing the composite material contact stress calculation method of any embodiment of the invention.

[0217] Those skilled in the art will clearly understand that, for the sake of convenience and brevity, the specific working processes of the systems, devices, and units described above can be referred to the corresponding processes in the foregoing method embodiments, and will not be repeated here.

[0218] In the several embodiments provided in this application, it should be understood that the disclosed systems, apparatuses, and methods can be implemented in other ways. For example, the apparatus embodiments described above are merely illustrative; for instance, the division of units is only a logical functional division, and in actual implementation, there may be other division methods. For example, multiple units or components may be combined or integrated into another system, or some features may be ignored or not executed. Furthermore, the coupling or direct coupling or communication connection shown or discussed may be an indirect coupling or communication connection between apparatuses or units through some interfaces, and may be electrical, mechanical, or other forms.

[0219] The units described as separate components may or may not be physically separate. The components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the units can be selected to achieve the purpose of this embodiment according to actual needs.

[0220] Furthermore, the functional units in the various embodiments of the present invention can be integrated into one processing unit, or each unit can exist physically separately, or two or more units can be integrated into one unit. The integrated unit can be implemented in hardware or as a software functional unit.

[0221] If the integrated unit is implemented as a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the present invention, in essence, or the part that contributes to the prior art, or all or part of the technical solution, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute all or part of the steps of the methods described in the various embodiments of the present invention. The aforementioned storage medium includes various media capable of storing program code, such as USB flash drives, portable hard drives, read-only memory (ROM), random access memory (RAM), magnetic disks, or optical disks.

[0222] The above-described embodiments are only used to illustrate the technical solutions of the present invention, and are not intended to limit it. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A method for calculating contact stress in composite materials, characterized in that, include: Obtain the geometric model of the object to be analyzed, and based on the pre-established composite material contact pairs, perform mesh generation on the geometric model to obtain the geometric model mesh; Based on the finite volume method and preset boundary conditions, the geometric model mesh is discretized by integral equations, and an iterative solver is selected to perform nonlinear iterative solution on the discretized results to obtain the mesh element displacement data of each interface in the geometric model mesh. Calculate the contact stress at the interface based on the displacement data of the mesh elements; When the direction of the frictional force meets the preset judgment condition, the contact frictional force of the interface is calculated based on the contact stress. The grid cell displacement data includes the grid center displacement vectors and displacement gradients of the left and right cells of different materials on both sides of the interface; the calculation of the contact stress at the interface based on the grid cell displacement data includes: By iterating through all interfaces, and based on the displacement gradient, the grid center displacement vector gradient is interpolated onto the interface to obtain the relative displacement of the interface. Calculate the face center displacement of the interface based on the relative displacement; Based on the face center displacement, determine whether there is contact between the elements of the interface. If there is, calculate the contact stress of the interface using the penalty function method based on the face center displacement. The face-center displacement includes the tangential displacement vector of the interface; when the frictional force direction meets the preset judgment condition, the contact frictional force of the interface is calculated based on the contact stress, including: If there is contact at the interface, determine whether the tangential displacement vector is not zero. If so, determine that there is friction between the units at the interface. Based on the tangential displacement vector, calculate the theoretical traction force of the interface, and based on the contact stress, calculate the maximum static friction force of the interface. If the theoretical traction force is less than the maximum static friction force, then the theoretical traction force is taken as the contact friction force of the interface. If the theoretical traction force is greater than the maximum static friction force, then the interface dynamic friction force of the interface is calculated based on the contact stress, and the interface dynamic friction force is used as the contact friction force of the interface.

2. The method for calculating contact stress of composite materials according to claim 1, characterized in that, The geometric model mesh is discretized using integral equations based on the finite volume method and preset boundary conditions. An iterative solver is then selected to perform nonlinear iterative solutions on the discretized results to obtain the displacement data of the mesh elements at each interface in the geometric model mesh, including: Based on the preset boundary conditions, the integral equation is discretized in each computational grid of the geometric model grid using the finite volume method to obtain the discretization coefficients. Then, based on the discretization coefficients, a system of discretized linear equations is constructed. Select an iterative solver, input the matrix coefficients of each term of the discrete linear equation system into the iterative solver for nonlinear iterative solution, and stop iterating when the solution results meet the preset calculation requirements, and output the grid cell displacement data of each interface in the geometric model mesh.

3. The method for calculating contact stress of composite materials according to claim 1, characterized in that, The face-centered displacement includes the normal displacement vector of the interface. Determining whether there is contact between the elements of the interface based on the face-centered displacement includes: Determine whether the normal displacement vector is greater than zero and exceeds a preset error value. If so, determine that there is contact between the units at the interface. The step of calculating the contact stress at the interface using the penalty function method based on the face center displacement includes: The contact stress at the interface is calculated using the penalty function method based on the normal displacement vector.

4. The method for calculating contact stress of composite materials according to claim 1, characterized in that, The geometric model corresponds to multiple geometric sub-models, and each geometric sub-model corresponds to a material distribution region. The process of meshing the geometric model based on pre-established composite material contact pairs to obtain a geometric model mesh includes: Based on the pre-established composite material contact pairs, the meshing method of the geometric sub-model is determined according to the region where the material distribution area is located; Based on the meshing method corresponding to each of the geometric sub-models, the computational mesh of the geometric model is divided to obtain the geometric model mesh of the object to be analyzed.

5. The method for calculating contact stress of composite materials according to any one of claims 1 to 4, characterized in that, Also includes: Obtain the object to be analyzed, which corresponds to multiple different material distribution regions; A geometric sub-model is established for each of the aforementioned material distribution regions, and material property values ​​are assigned to each of the aforementioned geometric sub-models during the parameter initialization process based on the material properties of different material distribution regions. Determine the interfaces between each pair of material distribution regions, and establish composite material contact pairs with penalty functions within the contact range that meet the preset contact conditions, based on the stress deformation of each geometric sub-model. The individual geometric sub-models are assembled into the geometric model of the object to be analyzed.

6. A composite material contact stress calculation device, characterized in that, include: The geometric model meshing module is used to obtain the geometric model of the object to be analyzed, and to perform meshing on the geometric model based on the pre-established composite material contact pairs to obtain the geometric model mesh. The nonlinear iterative solution module is used to discretize the geometric model mesh using integral equations based on the finite volume method and preset boundary conditions, and select an iterative solver to perform nonlinear iterative solution on the discretized results to obtain the mesh element displacement data of each interface in the geometric model mesh. The contact stress calculation module is used to calculate the contact stress at the interface based on the displacement data of the mesh elements. The contact friction calculation module is used to calculate the contact friction force at the interface based on the contact stress when the direction of the friction force meets the preset judgment conditions. The grid cell displacement data includes the grid center displacement vectors and displacement gradients of the left and right cells of different materials on both sides of the interface; the contact stress calculation module includes: The gradient interpolation module is used to iterate through all interfaces and interpolate the grid center displacement vector gradient onto the interface based on the displacement gradient to obtain the relative displacement of the interface. The face-center displacement calculation module is used to calculate the face-center displacement of the interface based on the relative displacement. The contact stress calculation submodule is used to determine whether there is contact between the elements of the interface based on the face center displacement. If there is contact, the contact stress of the interface is calculated using the penalty function method based on the face center displacement. The face-center displacement includes the tangential displacement vector of the interface; the contact friction force calculation module includes: The friction force determination module is used to determine whether the tangential displacement vector is not zero when there is contact at the interface. If so, it is determined that there is friction force between the units of the interface. The friction force related calculation module is used to calculate the theoretical traction force of the interface based on the tangential displacement vector, and to calculate the maximum static friction force of the interface based on the contact stress. The contact friction force determination module is used to take the theoretical traction force as the contact friction force of the interface when the theoretical traction force is less than the maximum static friction force. The interface dynamic friction calculation module is used to calculate the interface dynamic friction force of the interface based on the contact stress when the theoretical traction force is greater than the maximum static friction force, and to use the interface dynamic friction force as the contact friction force of the interface.

7. An electronic device, characterized in that, The device includes a processor and a memory: The memory is used to store program code and transmit the program code to the processor; The processor is used to execute the composite material contact stress calculation method according to any one of claims 1-5 according to the instructions in the program code.

8. A computer-readable storage medium, characterized in that, The computer-readable storage medium is used to store program code for executing the composite material contact stress calculation method according to any one of claims 1-5.

Images