A simulation method and system for progressive damage and failure removal of composite materials

Abstract

This invention discloses a method and system for simulating progressive damage and failure removal in composite materials, belonging to the field of computer-aided engineering. The method includes the following steps: S1, coupled solution: solving the partial differential governing equations describing heat transfer and solid mechanics, considering the coupling effects caused by thermal strain and the influence of temperature on material properties to obtain the temperature field and displacement field distribution of the current time step or iteration step; S2, material state calculation: based on the temperature and displacement fields obtained in the coupled solution step, calculating the material state for each integral point of the model using a predefined composite material constitutive model; S3, element processing and mesh updating: performing element processing and mesh updating after the material state calculation is completed. This invention can more comprehensively and accurately predict the final failure behavior of the structure.

Description

Technical Field

[0001] The present invention relates to the field of computer-aided engineering, and more particularly to a method and system for simulating progressive damage and failure removal of composite materials. Background Art

[0002] Composite materials are widely used in key sectors such as aerospace, energy, and transportation due to their high specific strength, high specific modulus, and excellent environmental resistance. However, in real-world operating environments, these materials are often subjected to complex loading conditions, particularly the combined effects of high temperatures and mechanical loads, known as coupled thermal-mechanical conditions. Accurately predicting the structural response, damage evolution, and ultimate failure of composite materials in such environments is crucial to ensuring the safety and reliability of structures.

[0003] Although numerical techniques, especially the finite element method (FEM), have become important tools for analyzing the behavior of composite materials, existing technologies still face many challenges in simulating progressive damage and failure under coupled thermal-mechanical conditions:

[0004] 1) Limitations of Traditional Methods in Simulating Material Removal: While current methods can reasonably predict the evolution of stress and strain distributions and damage variables (such as stiffness degradation factors) within a structure, they typically cease analysis after reaching a certain level of damage or only provide a "cloud map" of the damage state. However, in many practical scenarios, such as high-temperature ablation, high-velocity impact erosion, or when a structure is loaded until it fractures and breaks, the accumulation of damage ultimately leads to the physical removal of material. Traditional methods struggle to simulate this critical process and are unable to naturally account for the dramatic geometric changes caused by material removal, the dynamic redistribution of loads across the remaining structure, and the emergence of new boundary conditions on the newly created free surfaces (e.g., exposed internal material beginning to experience heat flux or pressure). Predicting the damage state without simulating material removal can significantly bias predictions of the structure's ultimate load-bearing capacity and failure modes.

[0005] 2) Difficulty integrating refined constitutive models with dynamic mesh modification: To accurately describe damage mechanisms, researchers have developed numerous refined material constitutive models, such as progressive damage models based on criteria such as Hashin, Puck, or LaRC04. These models are able to distinguish between different failure modes and account for their evolution. However, driving element deletion (a dynamic mesh modification technique) based on the physical states predicted by these models and tightly and efficiently integrating the two within a unified thermal-mechanical coupling numerical framework is a challenge in existing technologies. This involves complex data transfer (transferring the failure state calculated by the constitutive model to the mesh modification module), numerical stability (element deletion may lead to singular or drastic changes in the stiffness matrix), and computational efficiency (frequent mesh modifications can be very time-consuming).

[0006] 3) Lack of effective methods to trigger element removal based on physical damage states and automatically handle the consequences: Some current methods for simulating material failure may rely on simplified criteria, such as directly deleting elements when stress or strain reaches a single threshold. This approach may not accurately reflect the actual material failure process caused by the accumulation of multiple damage modes. A more ideal approach would allow a sophisticated material damage model to determine whether a material point has completely lost its load-bearing capacity and should be removed based on its complex internal physical state (for example, whether multiple damage variables reach critical values, whether energy dissipation meets specific conditions, etc.). In addition, when an element is removed, the original internal surfaces of its neighboring elements are exposed as new free surfaces. The lack of a general and robust method that can automatically identify these new surfaces and apply appropriate physical boundary conditions based on the physical scenario (for example, the exposed surface after ablation is subjected to new heat flux) is also a significant shortcoming of the existing technology. Summary of the Invention

[0007] The purpose of the present invention is to overcome the shortcomings of the existing technology and provide a simulation method and system for the progressive damage and failure removal of composite materials. The method can not only accurately simulate the complex progressive damage process of composite materials in a thermal-mechanical coupling environment, but also dynamically remove failed units according to the actual physical damage state inside the material, and automatically handle the resulting changes in geometry and boundary conditions, thereby more comprehensively and accurately predicting the final failure behavior of the structure.

[0008] The object of the present invention is achieved through the following solutions:

[0009] A method for simulating progressive damage and failure removal of composite materials, comprising the following steps:

[0010] S1, coupled solution: Solve the partial differential governing equations describing heat transfer and solid mechanics, taking into account the coupling effects caused by thermal strain and the influence of temperature on material properties, to obtain the temperature field and displacement field distribution of the current time step or iteration step;

[0011] S2, material state calculation: Based on the temperature and displacement fields obtained in the coupled solution step, the material state is calculated for each integration point of the model using the predefined composite material constitutive model;

[0012] S3, element processing and mesh updating: After the material state calculation is completed, perform element processing and mesh updating.

[0013] Furthermore, in step S2, the material state is calculated using a predefined composite material constitutive model, specifically including the following sub-steps:

[0014] calculating at least one or more internal damage state variables;

[0015] and calculation of effective material properties taking into account the current damage state;

[0016] And calculate and store an "element deletion flag" material property that indicates whether the material at the corresponding integration point has reached a preset complete failure state.

[0017] Furthermore, in step S3, the execution unit processes and updates the grid, specifically including the following sub-steps:

[0018] S31, access the calculated "element deletion mark" material properties;

[0019] S32, for each cell in the computational domain, determine it as a cell to be deleted based on whether the value of the "cell deletion flag" of one or more integration points satisfies a preset deletion threshold;

[0020] S33, executing a cell deletion operation to remove all cells determined to be to be deleted from the computational grid;

[0021] S34, performing a boundary update operation to identify newly exposed cell surfaces due to cell deletion, and assigning predefined boundary labels to these newly exposed surfaces to the original eroded boundary;

[0022] S35, updating the topology and related data of the mesh to reflect the results of the cell deletion and boundary update.

[0023] Furthermore, in step S32, the preset deletion threshold includes a flag value indicating failure, and the criteria for the "unit deletion flag" reaching failure include:

[0024] at least one of the one or more internal damage state variables reaches or exceeds a predefined critical damage threshold;

[0025] and / or a critical energy dissipation threshold value calculated based on the fracture energy of the material and the characteristic size of the element is reached or exceeded by the accumulated energy dissipation density of the corresponding integration points.

[0026] Furthermore, in step S34, the execution of the boundary update operation specifically includes the following sub-steps:

[0027] After a unit E_del is determined to be to be deleted, an adjacent unit E_adj that is adjacent to E_del and not determined to be to be deleted is identified;

[0028] For each neighboring cell E_adj, determine the face S_adj that it shares with E_del;

[0029] A predefined boundary identifier is assigned to the surface S_adj so that it can assume the corresponding boundary conditions in subsequent calculations.

[0030] Furthermore, it also includes the steps of detecting and deleting isolated units:

[0031] defining at least one or more stability boundaries;

[0032] Starting from the cells adjacent to the stable boundary, identifying all active cells connected to the stable boundary by a mesh connectivity traversal algorithm;

[0033] All active cells that are not identified as connected to the stable boundary are judged as island cells;

[0034] All cells that are identified as island cells are removed from the computational grid.

[0035] Furthermore, the grid connectivity traversal algorithm includes a BFS or DFS algorithm.

[0036] Furthermore, in step S2, the composite material constitutive model includes a progressive damage model based on the Hashin failure criterion.

[0037] Furthermore, the energy dissipation density is calculated specifically based on stress, strain and damage increment.

[0038] A simulation system for progressive damage and failure removal of composite materials includes a computer device, wherein the computer device includes a processor and a memory, wherein a computer program is stored in the memory, and when the computer program is loaded by the processor, any of the above methods is executed.

[0039] The beneficial effects of the present invention include:

[0040] (1) Improving simulation physical fidelity: By introducing a "unit failure marker" calculated based on the material's internal physical state (fine damage variables or energy state) to trigger unit removal, rather than a simple stress / strain threshold, the proposed method can more accurately reflect the material's true process from damage accumulation to ultimate loss of load-bearing capacity and physical separation. This is crucial for predicting the ultimate failure mode and structural life of composite materials under complex loads.

[0041] (2) Simulating the dynamic evolution of structural topology: The proposed method can dynamically change the geometric topology of the computational model by physically deleting failed elements. This makes it possible to simulate changes in structural shape caused by material removal (such as ablation, erosion, and shedding of fractured fragments), overcoming the limitations of traditional fixed-grid methods and realizing the non-uniform retreat of materials.

[0042] (3) Accurately capture load redistribution and stress concentration: Element removal disrupts the original load transfer path, leading to complex stress redistribution across the remaining structure and potentially generating stress concentrations at the edges of the new geometry. The proposed method can naturally capture these effects by dynamically updating the mesh and boundaries, thereby more accurately predicting the further expansion of damage and the ultimate failure of the structure.

[0043] (4) Simulating the physical interactions of newly generated free surfaces: Element removal generates new free surfaces. By identifying and marking these new boundaries and allowing the original physical boundary conditions (such as heat flow, pressure, chemical reaction, etc.) to be applied to them, the proposed method can simulate the continuous interaction between these new surfaces and the external environment, which is crucial for problems such as ablation and heat transfer analysis.

[0044] (4) Enhanced modularity and versatility of the simulation method: This paper separates the computation of complex material constitutive behavior (including failure judgment) from the execution logic of element removal and mesh processing into separate computational modules. This decoupling design improves the modularity of the method, facilitating users to replace or improve material models, failure criteria, or element deletion strategies as needed, thereby enhancing the versatility and scalability of the method.

[0045] (5) Providing a powerful tool for analyzing structural responses under extreme working conditions: Combining the above advantages, the present invention provides a powerful numerical tool with higher physical fidelity for analyzing the behavior of composite materials under thermomechanical coupling conditions, such as ablation of aircraft thermal protection systems, anti-erosion of composite materials structures and other complex scenarios. BRIEF DESCRIPTION OF THE DRAWINGS

[0046] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following briefly introduces the drawings required for use in the embodiments or the description of the prior art. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying any creative work.

[0047] Figure 1 Delete the thermal-mechanical coupling numerical calculation flow chart for hashin damage-based elements;

[0048] Figure 2 It is the energy triangle diagram;

[0049] Figure 3 A diagram comparing the breadth-first search (BFS) and depth-first search (DFS) algorithms.

[0050] Figure 4 is the calculation domain distribution diagram;

[0051] Figure 5a This is the unit failure distribution diagram at time t = 0.07s;

[0052] Figure 5b This is the unit failure distribution diagram at time t = 0.68s;

[0053] Figure 5c This is the unit failure distribution diagram at time t = 1.68s;

[0054] Figure 6a The stress distribution diagram at time t = 0.07s;

[0055] Figure 6b is the stress distribution diagram at time t = 0.68s;

[0056] Figure 6c This is the stress distribution diagram at time t = 1.68s. DETAILED DESCRIPTION

[0057] All features disclosed in all embodiments in this specification, or steps in all methods or processes implicitly disclosed, except for mutually exclusive features and / or steps, can be combined and / or expanded or replaced in any manner.

[0058] The present invention aims to address several limitations of existing technologies in simulating the complex failure behavior of advanced materials such as composites in a thermomechanically coupled environment. It tightly couples continuum mechanics solutions, material state assessment based on a bilinear damage model, and dynamic mesh topology modification within a unified computational framework.

[0059] Specifically, a method for simulating progressive damage and failure removal of composite materials is provided, specifically a numerical method for simulating progressive damage and unit failure removal of composite materials under thermal-mechanical coupling. Figure 1 As shown, the method is executed in a numerical calculation framework (such as the finite element method) in a time step sequence, and each time step or the iterative step therein includes at least:

[0060] S1, coupled solution step: solve the partial differential control equations describing heat transfer and solid mechanics, taking into account the coupling effect caused by thermal strain and the influence of temperature on material properties, to obtain the temperature field and displacement field distribution of the current time step or iteration step;

[0061] S2, material state calculation step: Based on the temperature and displacement fields obtained in the coupled solution step, a predefined composite material constitutive model (such as the progressive damage model based on the Hashin failure criterion) is used for the integration point of each element in the model:

[0062] calculating at least one or more internal damage state variables (e.g., fiber tension / compression damage, matrix tension / compression damage);

[0063] Compute effective material properties (e.g., degraded elastic stiffness matrix) that take into account the current damage state;

[0064] Calculates and stores an "element deletion marker" material property that indicates whether the material at that integration point has reached a predetermined state of complete failure.

[0065] Figure 1 The following is a flow chart for the thermomechanical coupling numerical calculation of element deletion based on Hashin damage. It consists of two main parts: a material failure property calculation module and an element deletion module. The former provides the latter with the material properties of the integration points at the current time step, which is used to determine whether element deletion is necessary.

[0066] Among them, material properties and status updates (corresponding to Figure 1 The material property calculation module in the calculation domain performs the following calculation process for each integration point of each unit:

[0067] Calculate mechanical strain (A1): Calculate pure mechanical strain based on total strain (from displacement field) and thermal strain (from temperature field).

[0068] Get historical damage (A2): Read and use the damage state variables stored in the previous time step or iteration step.

[0069] Calculate effective stress (A3): Calculate the effective stress of the material taking into account the current accumulated damage.

[0070] Calculate damage initiation (A4): Based on the effective stress state and the material's damage initiation criterion (such as the Hashin criterion), determine whether a new damage mode is activated or whether the existing damage mode needs to be further developed.

[0071] Calculate damage evolution (A5): If the damage initiation criterion is met, the state variables of each damage mode in the current step are calculated according to the damage evolution law. If the damage initiation criterion is not met, the stiffness matrix is ​​directly calculated to calculate the lossless stress.

[0072] Update damage variables (A6): Compare the calculated damage state variable with the damage state variable in the previous step, take the maximum value to ensure that the damage is irreversible, and obtain the new damage state variable of the current integration point (such as fiber damage state variable df, matrix damage state variable dm).

[0073] Determine Element Failure (A7): Based on the updated damage state variable, determine whether the current integration point meets the preset complete failure criteria (corresponding to "Does the element meet the failure criteria?" in Figure 1, such as df > 0.95, dm > 0.95, or energy dissipation meeting specific conditions). Set Element Failure Index (A7): Based on the result of the above complete failure criteria, explicitly set a Boolean or real number "Element Failure Index" (element_deleted): set to 1 (or true) if the failure criteria are met; set to 0 (or false) if not.

[0074] Calculate damage stiffness (A8): Calculate the effective stiffness matrix of the material in the current damage state based on the new damage state variables.

[0075] Calculate true stress (A9): Use the stiffness matrix and mechanical strain after damage to calculate the true stress of the current integration point and complete the material property calculation of the integration point.

[0076] Furthermore, for material property calculations, the damage model can consider the Hashin failure model to calculate material properties, mainly including the material damage initiation judgment and damage evolution process.

[0077] The material loss is initially judged using a two-dimensional hashin illustration as an example, and subsequently a three-dimensional hashin failure model can be considered for numerical calculations in higher dimensions.

[0078] Determination of material damage onset:

[0079] Fiber stretching :

[0080] (1)

[0081] Fiber compression :

[0082] (2)

[0083] Matrix stretching :

[0084] (3)

[0085] Matrix compression :

[0086] (4)

[0087] Where: : Transverse tensile strength; : longitudinal compressive strength;

[0088] : longitudinal tensile strength; : longitudinal compressive strength;

[0089] : longitudinal shear strength; : transverse shear strength;

[0090] : The influence coefficient of shear stress in the fiber stretching onset criterion, which can usually be considered to be 1;

[0091] 、 、 is the effective stress tensor The weight. 、 、 、 are the starting judgment criteria for fiber stretching, fiber compression, matrix stretching, and matrix compression, respectively. The subscript f represents fiber, the subscript m represents matrix, the superscript t represents stretching, and the superscript c represents compression.

[0092] It is worth noting that if the calculation is done in three dimensions, the damage initiation criterion will change if the interlayer effect is further considered. For example, the damage initiation criterion for the tensile and compressive failure of the matrix is ​​converted to and In addition, the interlayer failure Fi failure mode is added, and the initial criterion for tensile and compressive failure is and , the corresponding interlaminar tensile and compressive failure onset criteria are formulas (5) and (6), respectively.

[0093] Matrix stretching :

[0094] (5)

[0095] Matrix compression :

[0096] (6)

[0097] in, : interlaminar tensile strength; : interlaminar compressive strength; : Interlaminar shear strength. It represents the interlayer failure criterion, where the subscript i means interlayer.

[0098] In addition, the conversion relationship between the effective stress tensor and the true stress tensor in the above formula is: , taking the two-dimensional hashin failure model as an example, the parameters As shown in formula (7).

[0099] (7)

[0100] Where, : fiber damage variable; : matrix damage variable; : shear damage variable;

[0101] Three internal injury variables of appeal 、 、 There are four different injury patterns. 、 、 、 Calculation yields:

[0102] (8)

[0103] (9)

[0104] (10)

[0105] At this time, the constitutive equation with damage is:

[0106] (11)

[0107] Damage evolution of materials: In the ideal bilinear damage model, the stress-displacement relationship can be expressed as:

[0108] (12)

[0109] in: is stress, is the equivalent displacement, is the peak stress intensity, is the equivalent displacement corresponding to the peak stress, is the equivalent displacement of complete damage. No damage occurs, Indicates that damage has begun to occur, and the corresponding failure criterion is just greater than 1. ~ This is the damage process.

[0110] The bilinear damage evolution model is based on the energy release rate principle in fracture mechanics and describes the energy required for crack propagation: fracture energy It can be expressed as the area under the stress-displacement curve, that is:

[0111] (13)

[0112] Substituting the bilinear stress-displacement relationship into the above equation, we can obtain:

[0113] (14)

[0114] Compute the first integral:

[0115] (15)

[0116] Compute the second integral:

[0117] (16)

[0118] Add the two integral results:

[0119] (17)

[0120] Therefore, we can get:

[0121] (18)

[0122] Solving this equation yields the fully damaged equivalent displacement :

[0123] (19)

[0124] Finally, we get the "energy triangle", and the area of ​​the energy triangle is ,like Figure 2 shown.

[0125] In order to reduce the dependence of material damage evolution on the grid, the characteristic length .

[0126] The equivalent stress intensity and equivalent displacement of the four damage modes are defined as follows:

[0127] (1) Fiber stretching :

[0128] (20)

[0129] (twenty one)

[0130] (2) Fiber compression :

[0131] (twenty two)

[0132] (twenty three)

[0133] (3) Matrix stretching :

[0134] (twenty four)

[0135] (25)

[0136] (4) Matrix compression :

[0137] (26)

[0138] (27)

[0139] characteristic length It depends on the element geometry and form, for example it is the characteristic length of a straight line in a first-order element. represents the Macaulay bracket operator, for any have .

[0140] Behavior after injury, determined by specific pattern injury variables The relationship between and equivalent displacement is given by:

[0141] (28)

[0142] in, is the damage variable, ranging from 0 to 1. =0 means no damage, =1 indicates complete damage corresponding to the fiber damage variable , matrix damage variable , shear damage variable .

[0143] For the treatment of damage variables, we can also consider the viscosity regularization method to ensure the stability of the value and limit the damage increment to prevent excessive changes in a single step. Viscous regularization adds a viscosity parameter To smooth the damage evolution:

[0144] (29)

[0145] S3, element processing and mesh update step: After the material state calculation is completed (for example, after the nonlinear iterations converge or at the end of the time step):

[0146] Access calculated "element deletion marker" material properties;

[0147] For each element in the computational domain, the element is determined to be a to-be-deleted element based on whether the value of the "element deletion flag" of one or more of its integration points meets the preset deletion threshold (for example, the flag value indicates failure).

[0148] Execute the cell deletion operation to remove all cells that are determined to be deleted from the computational grid;

[0149] Performing a boundary update operation to identify newly exposed cell surfaces due to cell deletion and assigning predefined boundary labels to these newly exposed surfaces to the original eroded boundary;

[0150] Updates the mesh's topology and related data to reflect the results of element deletions and boundary updates.

[0151] In a further embodiment, the triggering basis of the deletion mark is designed as follows: the criteria for the "unit deletion mark" to reach a failure state include: at least one of at least one or more internal damage state variables reaches or exceeds a predefined critical damage threshold; and / or a critical energy dissipation threshold calculated based on the material fracture energy and the unit characteristic size is reached or exceeded by the energy dissipation density accumulated at the integration point (calculated based on stress, strain and damage increment).

[0152] In a further embodiment, the boundary update mechanism is designed as follows: the boundary update operation specifically includes: after a cell (denoted as E_del) is determined to be to be deleted, identifying neighboring cells (denoted as E_adj) adjacent to E_del that have not been determined to be to be deleted; for each neighboring cell E_adj, determining the surface (denoted as S_adj) previously shared with E_del; and applying a predetermined boundary identifier (e.g., the original eroded surface) to the surface S_adj so that the surface can bear the corresponding boundary conditions in subsequent calculations.

[0153] The island unit processing mechanism is designed as follows: define at least one or more stable boundaries; starting from the cells adjacent to the stable boundaries, identify all active cells connected to the stable boundaries through a grid connectivity traversal algorithm (such as BFS or DFS); determine all active cells that are not identified as connected to the stable boundaries as island cells; and remove all cells determined to be island cells from the computational grid.

[0154] More specifically, the element deletion process (corresponding to the element deletion module in Figure 1) executes the following element deletion logic after the material properties of all integration points are calculated:

[0155] Get failure index (B1): Traverse the elements in the calculation domain and read the "element failure index" (element_deleted) set by the material property calculation module.

[0156] Determine whether to delete (B2): Check whether the read "unit failure index" meets the preset deletion trigger condition (for example, the index value is equal to 1 or greater than or equal to a certain threshold, such as 0.99). If so, the unit that meets the deletion condition is added to the list to be deleted.

[0157] Update neighbor boundaries (B3): Before actual deletion, identify the neighbor cells of the cell to be deleted and mark the faces shared by these neighbor cells and the cell to be deleted as new free surfaces (for example, add them to the eroded boundary ID to ensure that the newly exposed boundary continues to be affected by erosion).

[0158] Delete Elements (B4): Removes the marked elements from the finite element mesh.

[0159] Update grid topology (B5): The grid system is updated to ensure data structure consistency.

[0160] Island processing (B6): Detect whether there are "island" units or unit clusters that are no longer connected to the main structure or the specified stable boundary due to unit deletion, and delete them.

[0161] Final grid update (B7): Perform necessary global grid data synchronization and prepare for numerical calculations in the next time step.

[0162] Time step advance (B8): Checks whether the simulation end time has been reached. If not, it proceeds to the next time step and repeats the entire process of solving the physics field, updating the material state, and deleting the elements.

[0163] The process of simulating material loss by element deletion is further detailed as follows:

[0164] (1) Unit removal decision and marking step: After the coupled field solver converges (for example, at the end of the time step), a unit processing module is started. This module traverses all units in the computational domain and reads the unit failure mark value of each unit (or its integration point). If the mark value is greater than a certain threshold, the unit is added to a "to be removed list". Here, the unit failure mark value can be the damage variable described above, or the energy judgment index, or both. Only when both are met can the unit be marked for removal.

[0165] (2) Neighboring Boundary Update Step: Traverse each cell in the "To Be Removed List" and find all its neighboring cells. If a neighboring cell is not in the "To Be Removed List," identify the shared interface (face or edge) between the two cells. Mark the shared interface as a new boundary and assign it a preset boundary identifier (meaning that after the material loss, the new material is affected by the boundary). This new boundary identifier can be used to apply heat flow, pressure, or other boundary conditions in subsequent calculations.

[0166] (3) Cell removal execution steps: Traverse the "to be removed list" and accurately remove the specified cells and no longer needed nodes from the existing mesh data structure. Update the connection relationship and physically delete these cells and their associated nodes from the computational mesh data structure without changing the coordinates of the remaining nodes. Unlike existing mesh reconstruction techniques that regenerate the mesh, this method directly modifies the mesh topology data structure.

[0167] (4) Island unit detection and removal step: After executing the unit removal, the connectivity analysis is performed. Starting from the pre-defined "stable boundary" (such as the fixed constraint boundary of the model, which can be considered as the boundary that will not be eroded in practice), a graph search algorithm (such as breadth-first search BFS or depth-first search DFS, such as Figure 3 (as shown in the figure) traverses all active cells reachable through cell adjacency relationships. Any unreachable active cells are identified as "island cells" and added to the "to-be-removed list," and then the removal step is performed. DFS, after reaching the end of a path, requires returning to the previous step to search for a second path. BFS, on the other hand, searches a larger area, first visiting all adjacent nodes of the current vertex, and then visiting all nodes in the next layer.

[0168] In other embodiments, based on the above method, the following specific calculation example is provided:

[0169] Computational domain: 2D rectangular entity with a width of 10 units and a height of 5 units. Figure 4 shown.

[0170] Mesh: A structured mesh is used for spatial discretization, consisting of 5000 four-node quadrilateral elements (QUAD4) arranged in a 100x50 grid.

[0171] Solving equations: mainly including transient heat conduction equation, force balance equation and constitutive equation with damage (Equation (11))

[0172] Transient heat conduction equation:

[0173] (30)

[0174] in, is the density, is the specific heat capacity, k is the thermal conductivity, and internal heat sources are not considered.

[0175] Force balance equation:

[0176] (31)

[0177] in, is the Cauchy stress tensor.

[0178] Calculation method: Due to material damage and possible geometric nonlinearity, an iterative method is used to solve a nonlinear algebraic equation system. For example, the Newton-Raphson method solves a linear equation system in each time step or load step. To iteratively update the solution vector Where u is a global vector containing all nodal unknowns (temperature and displacement components), R(u) is the global residual vector (representing the imbalance between the discretized governing equations and boundary conditions), J = ∂R / ∂u is the Jacobian matrix (or tangent stiffness matrix), and k is the number of iterations. Furthermore, adaptive time step control algorithms can be considered to accelerate convergence.

[0179] Boundary conditions: heat flux and pressure boundary conditions (1 kW / m², 1e6╳cos(30╳pi / 180)Pa (normal pressure) and -1e8╳sin(30╳pi / 180)Pa (shear force)) are given on the upper boundary, adiabatic boundaries are applied on the left and right sides, a constant temperature boundary (300K) at the bottom, and a fixed boundary (displacement in the x and y directions is 0).

[0180] The main material parameters are shown in Table 1 and can also be changed according to different composite material properties:

[0181] Table 1 Main material parameters of Hashin failure model

[0182]

[0183] Calculation results: Calculation results at different times (t = 0.07s, 0.68s, 1.68s), including the material failure attribute element_deleted such as Figure 5a 、 Figure 5b and Figure 5c As shown in Figure 2, the von Mises stress vonmises_stress is Figure 6a 、 Figure 6b and Figure 6c As shown in the figure, the material removal is not a uniform layered peeling, but rather forms a rough, jagged surface. This is consistent with the expectations of the physics-based damage model: material removal is determined by the local stress / strain state. Stress is redistributed as the geometry changes, and high-stress areas tend to be concentrated at the tips, protrusions, and roots of depressions of the remaining material. These areas become potential locations for subsequent damage extension and element removal. Some areas preferentially reach the damage threshold, resulting in the element_deleted being marked, which triggers the removal of local elements and causes non-uniform material retreat.

[0184] In addition, the calculation results show the strong coupling interaction between thermal loads, mechanical loads, material damage and geometric changes. Material removal is triggered by the internal material state (Hashin damage reaches a critical value) rather than based on simple geometric rules or externally defined removal rates, which makes the simulation closer to physical reality.

[0185] Computational results demonstrate that the Hashin damage model (other loss models, such as the Tsai-Wu failure criterion, are also considered) simulates progressive damage in composite materials under coupled thermal-mechanical loading, and dynamically simulates the material erosion process by removing elements using the element_deleted tag. This includes key physical characteristics such as surface recession, roughening, and stress redistribution. The innovative simulation of non-uniform recession and surface roughness provides a powerful tool for studying material failure in complex environments.

[0186] The units involved in the embodiments of the present invention may be implemented in software or hardware, and the units described may also be provided in a processor. In some cases, the names of these units do not limit the units themselves.

[0187] According to one aspect of an embodiment of the present invention, a computer program product or computer program is provided, comprising computer instructions stored in a computer-readable storage medium. A processor of a computer device reads the computer instructions from the computer-readable storage medium and executes the computer instructions, causing the computer device to perform the methods provided in the various optional implementations described above.

[0188] As another aspect, embodiments of the present invention further provide a computer-readable medium, which may be included in the electronic device described in the above embodiments, or may exist independently and not incorporated into the electronic device. The computer-readable medium carries one or more programs, and when executed by the electronic device, the electronic device implements the methods described in the above embodiments.

Claims

1. A method for simulating progressive damage and failure removal of composite materials, characterized in that: The following steps are involved: S1, coupled solution: Solve the partial differential governing equations describing heat transfer and solid mechanics, taking into account the coupling effects caused by thermal strain and the influence of temperature on material properties, to obtain the temperature field and displacement field distribution of the current time step or iteration step; S2, material state calculation: Based on the temperature and displacement fields obtained in the coupled solution step, the material state is calculated for each integration point of the model using the predefined composite material constitutive model; S3, unit processing and mesh updating: after the material state calculation is completed, unit processing and mesh updating are performed; In step S2, the material state is calculated using a predefined composite material constitutive model, specifically including the following sub-steps: calculating at least one or more internal damage state variables; and calculation of effective material properties taking into account the current damage state; and calculate and store an "element deletion flag" material property that indicates whether the material at the corresponding integration point has reached a preset complete failure state; In step S3, the execution unit processes and updates the grid, specifically including the following sub-steps: S31, access the calculated "element deletion mark" material properties; S32, for each cell in the computational domain, determine it as a cell to be deleted based on whether the value of the "cell deletion flag" of one or more integration points meets a preset deletion threshold; S33, executing a cell deletion operation to remove all cells determined to be to be deleted from the computational grid; S34, performing a boundary update operation to identify newly exposed cell surfaces due to cell deletion, and assigning predefined boundary labels to these newly exposed surfaces to the original eroded boundary; S35, updating the topology and related data of the grid to reflect the results of the cell deletion and boundary update; In step S34, the execution of the boundary update operation specifically includes the following sub-steps: After a unit E_del is determined to be to be deleted, an adjacent unit E_adj that is adjacent to E_del and not determined to be to be deleted is identified; For each neighboring cell E_adj, determine the face S_adj that it shares with E_del; A predefined boundary identifier is assigned to the surface S_adj so that it can assume the corresponding boundary conditions in subsequent calculations.

2. The method for simulating progressive damage and failure removal of composite materials according to claim 1, characterized in that: In step S32, the preset deletion threshold includes a flag value indicating failure, and the criteria for the "unit deletion flag" reaching failure include: at least one of the one or more internal damage state variables reaches or exceeds a predefined critical damage threshold; and / or a critical energy dissipation threshold value calculated based on the fracture energy of the material and the characteristic size of the element is reached or exceeded by the accumulated energy dissipation density of the corresponding integration points.

3. The method for simulating progressive damage and failure removal of composite materials according to any one of claims 1 to 2, characterized in that: It also includes the steps of island unit detection and deletion: defining at least one or more stability boundaries; Starting from the cells adjacent to the stable boundary, identifying all active cells connected to the stable boundary by a mesh connectivity traversal algorithm; All active cells that are not identified as connected to the stable boundary are judged as island cells; All cells that are identified as island cells are removed from the computational grid.

4. The method for simulating progressive damage and failure removal of composite materials according to claim 3, characterized in that: The grid connectivity traversal algorithm includes a BFS or DFS algorithm.

5. The method for simulating progressive damage and failure removal of composite materials according to claim 1, characterized in that: In step S2, the composite material constitutive model includes a progressive damage model based on the Hashin failure criterion.

6. The method for simulating progressive damage and failure removal of composite materials according to claim 2, characterized in that: The energy dissipation density is calculated based on stress, strain and damage increment.

7. A simulation system for progressive damage and failure removal of composite materials, characterized in that: The method comprises a computer device, wherein the computer device comprises a processor and a memory, wherein a computer program is stored in the memory, and when the computer program is loaded by the processor, the method according to any one of claims 1 to 2 is executed.

Images