WO2021139130A1 - 复合材料多尺度模型的动态渐进失效分析方法 - Google Patents
复合材料多尺度模型的动态渐进失效分析方法 Download PDFInfo
- Publication number
- WO2021139130A1 WO2021139130A1 PCT/CN2020/102477 CN2020102477W WO2021139130A1 WO 2021139130 A1 WO2021139130 A1 WO 2021139130A1 CN 2020102477 W CN2020102477 W CN 2020102477W WO 2021139130 A1 WO2021139130 A1 WO 2021139130A1
- Authority
- WO
- WIPO (PCT)
- Prior art keywords
- damage
- stress
- fiber
- matrix
- macroscopic
- Prior art date
Links
- 239000002131 composite material Substances 0.000 title claims abstract description 67
- 230000000750 progressive effect Effects 0.000 title claims abstract description 15
- 238000004458 analytical method Methods 0.000 title claims description 51
- 239000000835 fiber Substances 0.000 claims abstract description 60
- 239000011159 matrix material Substances 0.000 claims abstract description 54
- 230000003321 amplification Effects 0.000 claims abstract description 14
- 238000003199 nucleic acid amplification method Methods 0.000 claims abstract description 14
- 238000012217 deletion Methods 0.000 claims abstract description 10
- 230000037430 deletion Effects 0.000 claims abstract description 10
- 230000009471 action Effects 0.000 claims abstract description 6
- 238000005520 cutting process Methods 0.000 claims abstract description 3
- 239000000463 material Substances 0.000 claims description 28
- 238000000034 method Methods 0.000 claims description 25
- 230000014509 gene expression Effects 0.000 claims description 17
- 230000008569 process Effects 0.000 claims description 13
- 230000010354 integration Effects 0.000 claims description 10
- 238000009826 distribution Methods 0.000 claims description 9
- 238000004364 calculation method Methods 0.000 claims description 8
- 230000006835 compression Effects 0.000 claims description 7
- 238000007906 compression Methods 0.000 claims description 7
- 230000015556 catabolic process Effects 0.000 claims description 6
- 238000006731 degradation reaction Methods 0.000 claims description 6
- 238000004088 simulation Methods 0.000 claims description 6
- 238000011156 evaluation Methods 0.000 claims description 5
- 239000000758 substrate Substances 0.000 claims description 4
- 238000004422 calculation algorithm Methods 0.000 claims description 3
- 230000002427 irreversible effect Effects 0.000 claims description 3
- 238000002360 preparation method Methods 0.000 claims description 3
- 238000000354 decomposition reaction Methods 0.000 claims description 2
- 238000013461 design Methods 0.000 claims description 2
- 238000012546 transfer Methods 0.000 claims description 2
- 230000009466 transformation Effects 0.000 claims 1
- 230000004044 response Effects 0.000 abstract description 3
- 239000011347 resin Substances 0.000 abstract description 2
- 229920005989 resin Polymers 0.000 abstract description 2
- 230000006399 behavior Effects 0.000 abstract 3
- 230000007850 degeneration Effects 0.000 abstract 1
- 230000035882 stress Effects 0.000 description 54
- 238000005553 drilling Methods 0.000 description 24
- 238000010586 diagram Methods 0.000 description 18
- 230000032798 delamination Effects 0.000 description 15
- 239000004918 carbon fiber reinforced polymer Substances 0.000 description 11
- 238000006073 displacement reaction Methods 0.000 description 9
- 238000011068 loading method Methods 0.000 description 8
- 238000012360 testing method Methods 0.000 description 8
- 239000010410 layer Substances 0.000 description 6
- 238000012545 processing Methods 0.000 description 5
- 230000008859 change Effects 0.000 description 4
- 238000003754 machining Methods 0.000 description 4
- 230000000737 periodic effect Effects 0.000 description 4
- 238000000926 separation method Methods 0.000 description 4
- 238000006243 chemical reaction Methods 0.000 description 3
- 238000010835 comparative analysis Methods 0.000 description 3
- 238000002474 experimental method Methods 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- 229920000049 Carbon (fiber) Polymers 0.000 description 2
- 239000004917 carbon fiber Substances 0.000 description 2
- 239000011248 coating agent Substances 0.000 description 2
- 238000000576 coating method Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 230000008676 import Effects 0.000 description 2
- 230000000977 initiatory effect Effects 0.000 description 2
- 239000011229 interlayer Substances 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- VNWKTOKETHGBQD-UHFFFAOYSA-N methane Chemical compound C VNWKTOKETHGBQD-UHFFFAOYSA-N 0.000 description 2
- 239000011148 porous material Substances 0.000 description 2
- 239000011208 reinforced composite material Substances 0.000 description 2
- 239000002356 single layer Substances 0.000 description 2
- 239000007787 solid Substances 0.000 description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 1
- 229920001410 Microfiber Polymers 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 239000003822 epoxy resin Substances 0.000 description 1
- 230000006355 external stress Effects 0.000 description 1
- 238000012482 interaction analysis Methods 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000000691 measurement method Methods 0.000 description 1
- 239000003658 microfiber Substances 0.000 description 1
- 238000002156 mixing Methods 0.000 description 1
- 238000012821 model calculation Methods 0.000 description 1
- 238000013433 optimization analysis Methods 0.000 description 1
- 239000011120 plywood Substances 0.000 description 1
- 229920000647 polyepoxide Polymers 0.000 description 1
- 238000011002 quantification Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000002441 reversible effect Effects 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 238000012916 structural analysis Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/17—Mechanical parametric or variational design
-
- B—PERFORMING OPERATIONS; TRANSPORTING
- B32—LAYERED PRODUCTS
- B32B—LAYERED PRODUCTS, i.e. PRODUCTS BUILT-UP OF STRATA OF FLAT OR NON-FLAT, e.g. CELLULAR OR HONEYCOMB, FORM
- B32B15/00—Layered products comprising a layer of metal
- B32B15/04—Layered products comprising a layer of metal comprising metal as the main or only constituent of a layer, which is next to another layer of the same or of a different material
- B32B15/08—Layered products comprising a layer of metal comprising metal as the main or only constituent of a layer, which is next to another layer of the same or of a different material of synthetic resin
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N3/00—Investigating strength properties of solid materials by application of mechanical stress
- G01N3/02—Details
- G01N3/06—Special adaptations of indicating or recording means
- G01N3/066—Special adaptations of indicating or recording means with electrical indicating or recording means
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N2203/00—Investigating strength properties of solid materials by application of mechanical stress
- G01N2203/003—Generation of the force
- G01N2203/0053—Cutting or drilling tools
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N2203/00—Investigating strength properties of solid materials by application of mechanical stress
- G01N2203/0058—Kind of property studied
- G01N2203/0096—Fibre-matrix interaction in composites
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N2203/00—Investigating strength properties of solid materials by application of mechanical stress
- G01N2203/02—Details not specific for a particular testing method
- G01N2203/0202—Control of the test
- G01N2203/0212—Theories, calculations
- G01N2203/0218—Calculations based on experimental data
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N3/00—Investigating strength properties of solid materials by application of mechanical stress
- G01N3/58—Investigating machinability by cutting tools; Investigating the cutting ability of tools
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/26—Composites
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/02—Reliability analysis or reliability optimisation; Failure analysis, e.g. worst case scenario performance, failure mode and effects analysis [FMEA]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Definitions
- the invention relates to the field of mechanical response and progressive damage generated under dynamic load conditions of a composite material laminated plate structure.
- the traditional and simple and direct research method relies on experimental analysis, but the traditional method still has a lot of shortcomings, such as long preparation time, large material consumption, The cost is higher, and more importantly, it is impossible to observe and understand the mechanism of mechanical processing deformation and material failure from a microscopic perspective.
- the present invention provides a three-dimensional multi-scale dynamic progressive damage evolution model based on micro-failure criteria.
- the relationship between macro-stress and micro-stress is established by using the stress amplification factor, and the structural analysis of composite laminates is converted from macro-scale to micro-scale. .
- the three-dimensional damage constitutive relationship and damage evolution analysis of composite materials are established.
- the model was verified under dynamic load conditions and process parameters, and the composite (hole wall surface morphology, import and export delamination damage, etc.) was composited according to the mechanical processing quality damage evaluation criteria.
- the machining behavior of materials is evaluated to realize precision machining.
- the multi-scale analysis method proposed in this application is established based on the macro-micro mechanical failure theory of composite materials.
- the core of the method is to establish the macro-micro structure relationship between the components in the composite material through the macro-micro mechanical interaction analysis theory.
- the damage and failure behavior of each component structure under the action of external load is mainly divided into the following three steps:
- the first step is to ignore the initial damages such as bubbles and fiber breakage generated during the preparation of the composite material laminate.
- the composite material is idealized as a unidirectional laminate by superimposing the preset layer sequence, and the composite material is macroscopically structured. Carry out analysis, consider the influence of macro-structure design, material properties, layering method and overall structure load form on the stress and strain distribution inside the structure.
- the second step is to rationally simplify the fiber distribution of the unidirectional laminate based on the mechanical analysis of the unidirectional laminate.
- Based on the fiber volume fraction establish a representative volume unit of micromechanics for the apparent mechanical properties of the composite structure, and adopt a multi-scale analysis method. Transfer the load of the macro-structure unit or integration point to the micro-model, simulate the stress distribution of each unit integration point of the fiber and the matrix component material under this load, and use the micro-mechanical theory to judge the failure of the matrix and the fiber separately.
- the stiffness of the damaged element is reduced, the failure is deleted, and then the entire unit cell model is homogenized to obtain the macro-mechanical properties of the damaged material.
- the third step is to assign the macroscopic properties of the damaged material to the corresponding unit in the macroscopic structure, and then perform the next iteration of the macroscopic analysis of the overall structure, and repeat iteratively until the set time analysis step is completed, and finally realize that the composite material is assembled by the Multi-scale numerical simulation analysis of structural parts under the cutting action of tools.
- Figure 1 is a schematic diagram of a drilling multi-scale analysis method
- Figure 2 is a schematic diagram of representative volume unit models and reference points
- Figure 3 is a schematic diagram of the application of load conditions of a representative volume element model
- Figure 4 is a schematic diagram of unit characteristic length calculation and fiber failure state
- Figure 5 is a schematic diagram of the failure state of the substrate
- Figure 6 is a schematic diagram of a drilling multi-scale explicit finite element algorithm program
- Figure 7 is a schematic diagram of a layered failure model of cohesive unit
- Figure 8 is a flow chart of micromechanical analysis of representative volume element models
- Figure 9 is a schematic diagram of loading periodic boundary conditions of a representative volume element model
- Figure 10 is a schematic diagram of the analysis results of a representative volume element model
- Figure 11 is a schematic diagram of the stress magnification factor of key points on the fiber under different working conditions
- Figure 12 is a schematic diagram of the stress magnification factor of key points on the substrate under different working conditions
- Figure 13 is a schematic diagram of a dagger drilling test for T700S-12K/YP-H26 composite material.
- Figure 14 is a schematic diagram of the multi-scale finite element model of dagger drilling and drilling composite materials
- Figure 14 is a schematic diagram of comparative analysis of axial force and torque in drilling
- Figure 15 is a schematic diagram of the analysis and comparison of the surface morphology of the composite material hole wall
- Figure 16 is a schematic diagram of damage defects at the entrances and exits of composite material prefabricated holes
- Figure 17 is a schematic diagram of composite material export damage coefficient and error under different processing parameters
- the macroscopic stress of the composite laminate structure and the microscopic stress of each component structure are bridged by the stress amplification factor, and the stress amplification factor can be obtained by the finite element analysis result of the representative volume element containing the fiber and the matrix.
- the macro-micro conversion equation is used to realize the conversion of macro-stress to micro-stress.
- the expression of the macro-micro conversion equation is shown in formulas (1), (2), (3):
- the stress amplification factor method is used to extract the microscopic stress of the fiber and the matrix.
- Figure 2 according to the observation of the interface of the unidirectional composite laminate by the ultra-depth-of-field electron microscope, the staggered arrangement of the representative volume unit model of the rectangular structure is closer to the random distribution of fibers, and the rectangular model (Figure 2c) is similar. It is more accurate than the square representative volume unit model ( Figure 2b), so this application will use the rectangular representative volume unit model for analysis.
- the micro stress at the key points is extracted according to the model calculation results.
- the stress amplification factor obtained at the key reference points set in the representative volume element model is Based on this principle, six different numerical macro-loads are applied to the representative volume element model to obtain the stress amplification coefficients.
- the expression of each intensity tensor is as follows:
- the damage evolution law of the fiber in the X direction mainly under tension and compression can be expressed as:
- the constitutive relationship of the fiber is defined as the damage at the microscopic scale:
- ⁇ f , C f , and ⁇ f respectively represent the stress, stiffness and strain of the fiber on the microscopic scale; with Respectively represent the longitudinal and transverse damage modulus of the fiber on the micro scale; E f1 and E f2 represent the fiber’s modulus of elasticity without damage in the longitudinal and transverse directions on the micro scale, respectively.
- T mi and C mi respectively represent the initial tensile strength and compressive strength of the matrix
- ⁇ VM and I 1 respectively represent the ⁇ Von Mises equivalent stress and the first invariant of stress, which can be expressed as:
- the Stassi equivalent stress ⁇ eq is used to define the damage evolution of the matrix, and the expression is:
- the damage of the matrix is determined according to the stress state, and the damage value is between 0 and 1. After the initial damage is judged, the damage evolution law And the equivalent stress can be expressed as:
- ⁇ m , C m , and ⁇ m respectively represent the stress, stiffness and strain of the matrix at the microscopic scale;
- E m respectively represent the damaged elastic modulus and the undamaged elastic modulus of the matrix at the microscopic scale; Represents the maximum damage coefficient of all reference points on the substrate in the microscopic scale in the representative volume element model.
- d L and d T represent the longitudinal and transverse macroscopic damage variables of the composite material, respectively; Indicates the undamaged stiffness coefficient in the macroscopic constitutive of the composite material, and its expression is:
- V f represents the fiber volume fraction
- V m represents the matrix volume fraction
- the laminate model of the composite material in this application adopts the hexahedral element model, which is correct Solving the characteristic roots of the matrix can obtain the principal stretch ratios ⁇ 1 , ⁇ 2 , and ⁇ 3 in the three directions of the element. According to the definition of Biot strain (nominal strain), three principal strain values can be obtained:
- the auxiliary element deletion criterion can be expressed as:
- the multi-scale explicit finite element algorithm is shown in Figure 6.
- the realization of the entire multi-scale explicit finite element algorithm can be divided into three steps:
- Step (1) At the beginning of the global time analysis step n, pass the previous time analysis step Step with strain increment Sum to calculate the global total strain Through the previous composite element stiffness matrix Solve and calculate the macroscopic stress of each integration point
- Step (2) Multiply the key reference points respectively set for the component structural units according to the macro stress with Calculate the microscopic stress of the fiber and matrix separately with Based on the established failure criteria of the micromechanical failure theory, determine the failure coefficient of the fiber and the matrix, and obtain the microscopic damage variables of the fiber and the matrix through the established progressive damage model; determine the fiber damage coefficient and the matrix damage coefficient in the representative volume element model All reference points with The maximum value.
- Step (3) By establishing a new damage evolution criterion, the macroscopic damage variable at each integration point can be calculated from the maximum value of d L and d T to achieve element stiffness degradation, performance attenuation, etc., and based on the element deletion criterion, The deletion of the unit is realized. Since the macroscopic damage variable is irreversible at the n-1th time analysis step, the macroscopic damage variable should be selected with a larger value, and the composite material stiffness matrix will be realized at the n time analysis step. If the set time analysis step is reached, the next time analysis step will be entered. If the time analysis step is terminated, the program will be terminated.
- ⁇ i represents the nominal strain in the corresponding direction, which The value is the initial thickness of the cohesive element after the separation position at the interface is removed;
- K ii represents the constitutive stiffness of the cohesive element;
- the damage state variable D s in the mixed mode can be expressed as
- ⁇ > represents Macaulay brackets
- ⁇ represents the BK criterion parameter in the hybrid mode.
- the research object of this application is carbon fiber reinforced composite materials with a value of 1.45
- each unit stress load under six working conditions is respectively applied to the representative volume element model to obtain the microscopic stress distribution, as shown in Figure 10.
- 10 red solid circles are used to indicate the key reference points of micro-fibers and 15 purple solid circles are used to indicate the stress magnification factor calculated from the reference points of the micro-matrix of the matrix, as shown in Figure 2.
- the F10 reference point and the M10 ⁇ M15 points are located at the center of the fiber and the center of the matrix, respectively.
- F1 ⁇ F9 reference points and M1 ⁇ M9 consider the key reference points of the position where the fiber and the matrix join.
- the microscopic stress distribution of each key reference point is the stress amplification factor.
- the stress of each key reference point on the fiber and matrix under different working conditions is shown in Figure 10 and Figure 11. Show. Since the inherent properties of the material will not change, the representative volume element model only needs to be solved once and the key reference point results are extracted. After obtaining the stress amplification factor, store it in a separate file in the Fortran code. In the macro composite drilling finite element model, use the “#include” function to call the file based on the ABAQUS VUMAT code to realize the dynamics of the composite material. Multi-scale analysis of progressive failure.
- This application uses T700S-12K/YP-H26 high temperature resistant epoxy resin composite material as the research object, and the fiber volume fraction content is about 59%, using [(0°/90°/45°/-45°) s ] 4 Laying sequence, the overall size of the composite laminate structure is 15mm ⁇ 15mm ⁇ 5.76mm, a total of 32 layers, under the preset global coordinates, the single-layer composite material is given orthotropic anisotropy along the fiber direction
- the parameters of the single-layer composite materials are shown in Table 2, and the parameters of the cohesion unit are shown in Table 3.
- a cemented carbide dagger drill with a diameter of 6mm (AlTiN coating) is set as a discrete rigid body model.
- the total mass and moment of inertia of the dagger drill are constrained At the top reference point position, and based on the reference point, impose boundary conditions such as the feed speed and rotation speed of the drill bit.
- the mesh refinement process is carried out on the drilling contact area (a circular area with a diameter of 8mm) involving composite laminates, and coarser meshes are used for the area near the edge of the laminate and the dagger drill. grid.
- the composite material uses C3D8R elements (8 nodes, three-dimensional reduced integration elements), the layered area uses COH3D8 viscous elements (0mm thickness), and the dagger drill uses C3D10M discrete rigid body elements.
- the mesh size of the composite material unit model is guaranteed to be a meso scale. According to the optimization analysis result of the unit mesh sensitivity model, there are a total of 4,203438 in the entire drilling finite element model.
- the composite material model contains 3,382,122 hexahedral units, the layered area contains 819,108 cohesion units, and the dagger drill contains 8,208 tetrahedral units.
- the minimum size of the composite unit is about 45 ⁇ 45 ⁇ 40 ⁇ m.
- boundary conditions such as speed and feed rate are applied to the entire drilling finite element model. Since the motion state of the dagger drill model is constrained at the top reference point, the reference point is The displacement in the X and Y directions is restricted, and the feed speed is applied in the Z direction. Similarly, the rotation speed in the X and Y directions is restricted, and the clockwise rotation speed is applied in the Z direction, and the composite material layer The four vertical surfaces of the plywood are fixed. After the simulation calculation is completed, this application verifies the model through an experimental platform built.
- the experimental platform mainly includes three systems of CFRP drilling, data acquisition, and experimental observation. The schematic diagram of the experimental scheme and the connection of the test device are shown in Figure 13.
- the axial force and torque of the drilling multi-scale finite element model and experiment have the same value and change trend at each stage of drilling, only in the T IV and T V In the phase, the reduction speed of axial force and torque is slower than that of the experiment.
- the main reason is that the fiber adopts the degradation method based on fracture toughness in the multi-scale finite element model, and the set dissipation energy value is based on empirical evaluation.
- the bonding area CEs in the simulation of the interlayer delamination phenomenon some macro-units in the model are not deleted, resulting in too long contact time between the unit and the drill bit, resulting in some smaller forces and torques.
- the CFRP pore wall is mainly composed of the resin-coated surface and the fiber fracture surface, and the pore wall coating surface near the entrance is more obvious.
- the damage of CFRP from the entrance to the exit of drilling is mainly manifested as splitting at the entrance, tearing at the exit, burrs, separation between layers, radial crushing, microcracks and other damages.
- the damage position of each layer in Figure 15 (a) and (b) it is found that the pits and the hole wall damage positions are almost the same. Therefore, the multi-scale finite element model established in this paper can simulate the damage state of the hole wall of CFRP drilling more realistically.
- CFRP is prone to burrs, tearing, delamination and other processing damages at the entrance and exit during the process of making holes.
- the dagger drill is a special CFRP drill with integrated drilling and reaming, only the fiber fracture at the entrance of the hole when drilling the CFRP laminate It is relatively flat, the entrance tearing phenomenon is not obvious and there are fewer burrs, but through the finite element analysis model and test results, there are more obvious burrs, tears, and delamination damage at the exit, as shown in Figure 16(a) .
- the multi-scale finite element model is based on the deletion of macroscopic elements to realize the simulation of material removal, the scale effect of the element can only simulate the damage phenomenon at the macro level, and cannot achieve the quantification of burr and tear damage.
- the delamination damage system is based on the tool diameter.
- the area calculated by the number of deleted elements around the hole is divided by the area of the reference hole to calculate the delamination damage coefficient.
- the measurement method is to export the picture in the ABAQUS software, and then import it into the CAD software.
- the belt tool adopts the equal proportion method to realize the measurement respectively. Because the drill bit rotates too fast under actual working conditions, and the deletion of the multi-scale finite element model unit is irrecoverable, the area of the damaged unit in the center of the layered area is used to calculate the layered damage coefficient, as shown in Figure 24(b) Shown. Therefore, the delamination damage coefficient can be expressed as:
- D d represents the delamination damage coefficient of the outlet
- a nom represents the area of the central aperture area
- a d represents the area of the damage unit caused by the delamination.
- the multi-scale finite element model established in this paper can truly simulate the CFRP damage state under actual working conditions.
- the maximum errors of the axial force, torque and outlet delamination damage coefficient are 3.37%, 7.69%, and 4.28%, respectively.
- it can also simulate the phenomenon of hole wall surface damage and exit delamination more realistically.
Landscapes
- Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Geometry (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Hardware Design (AREA)
- General Engineering & Computer Science (AREA)
- Evolutionary Computation (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Analysis (AREA)
- Computational Mathematics (AREA)
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Chemical & Material Sciences (AREA)
- Analytical Chemistry (AREA)
- Biochemistry (AREA)
- General Health & Medical Sciences (AREA)
- Immunology (AREA)
- Pathology (AREA)
- Investigating Strength Of Materials By Application Of Mechanical Stress (AREA)
Abstract
一种复合材料层合板结构在动态工况下产生的力学响应及渐进损伤行为，首先，考虑纤维和基体在微观状态下的不同力学行为，提出了一种基于微观失效理论的三维多尺度动态渐进损伤演化规律模型，基于典型代表性体积单元模型中微观组分退化弹性参数，提出了新型的纤维和树脂基体的损伤演化规律模型及失效单元的辅助删除准则；其次，采用应力放大系数建立复合材料模型中的宏观应力与代表性体积单元微观应力关系模型，结合双线性粘聚力单元模型，模拟了复合材料在匕首钻切削作用下的层内及层间损伤行为。
Description
本发明涉及一种复合材料层合板结构动态载荷工况下产生的力学响应及渐进损伤领域。
为减少对复合材料结构在动态载荷工况下出现的损伤，传统且简单直接的研究方法是依赖于实验分析，但是传统的方法依旧存在大量的不足，如前期准备时间长、材料消耗量大、成本较高等，更重要的是无法从微观上观察和理解机械加工变形、材料失效等机理。
发明内容
本发明提供了一种基于微观失效准则的三维多尺度动态渐进损伤演化模型，利用应力放大系数建立了宏观应力与微观应力之间的关系，将复合材料层合板结构分析从宏观尺度转换到微观尺度。然后，建立复合材料三维损伤本构关系和损伤演化分析，通过引入组分失效的损伤变量，渐进退化材料刚度，模拟复合材料在动态载荷工况下各相材料的失效，分层等行为。基于相同铺层序列方式的复合材料试样，在动态载荷工况工艺参数下对模型进行了验证，并根据机械加工质量损伤评价标准对(孔壁表面形貌、进出口分层损伤等)复合材料的机械加工行为进行评价，从而实现精密机械加工。
本申请提出的多尺度分析方法是基于复合材料宏-微观力学失效理论而建立的，该方法的核心是通过宏-微观力学交互分析理论建立复合材料中组成成分之间的宏-微观结构联系及在受到外载荷作用下的各组分结构的损伤和失效行为。根据图1所示，多尺度分析方法主要分为以下三个步骤：
第一步，忽略复合材料层合板在工艺制备中产生的气泡、纤维丝断裂等初始损伤，将复合材料理想化为单向层合板预设的铺层序列叠加而成，对复合材料进行宏观结构进行分析， 考虑宏观结构设计，材料属性，铺层方式以及整体结构受载形式对结构内部的应力和应变分布的影响。
第二步，根据单向层合板的力学分析，对单向板纤维分布进行合理简化，基于纤维体积分数建立面向复合材料结构力学性能表观的微观力学代表性体积单元，采用多尺度分析方法，将宏观结构单元或积分点的受载情况传递到微观模型中，在该载荷下模拟纤维和基体组分材料各单元积分点的应力分布，采用微观力学理论分别对基体，纤维进行失效判断，假设各对应单元发生损伤，则对损伤单元进行刚度折减、失效删除，再对整个单胞模型进行均匀化处理，得到损伤后的材料宏观力学性能。
第三步，将损伤后的材料宏观性能赋予宏观结构中的对应单元中，对整体结构再进行下一迭代步的宏观分析，如此反复迭代直至设定时间分析步完成，最终实现复合材料由组分到结构件的在刀具切削作用下的多尺度数值模拟分析。
图1是钻削多尺度分析方法示意图；
图2是代表性体积单元模型和参考点示意图；
图3是代表性体积单元模型载荷工况的施加示意图；
图4是单元特征长度计算与纤维失效状态示意图；
图5是基体失效状态示意图；
图6是钻削多尺度显式有限元算法程序示意图；
图7是内聚力单元分层失效模型示意图；
图8是代表性体积单元模型微观力学分析流程图；
图9是代表性体积单元模型周期性边界条件加载示意图；
图10是代表性体积单元模型分析结果示意图；
图11是不同工况下纤维上关键点应力放大系数示意图；
图12是不同工况下基体上关键点应力放大系数示意图；
图13是匕首钻钻削T700S-12K/YP-H26复合材料试验示意图。
图14是匕首钻钻削复合材料多尺度有限元模型示意图
图14是钻削轴向力与扭矩对比分析示意图；
图15是复合材料孔壁表面形貌分析对比示意图；
图16是复合材料预制孔出入口位置损伤缺陷示意图；
图17是不同加工参数下的复合材料出口损伤系数及误差示意图；
复合材料层合板结构的宏观应力与各组分结构微观应力通过应力放大系数进行桥联，而应力放大系数的获取可以通过包含纤维、基体的代表性体积单元的有限元分析结果得到。得到应力放大系数后，通过宏-微观转化方程实现宏观应力向微观应力的转化，其中，宏-微观转化方程的表达式如式(1)、(2)、(3)所示：
式(1)～(3)中：σ
i表示纤维或基体n点处的应力，n为受力分析点，i＝1～6；
表示复合材料层合板的宏观应力分量，j＝1～6；
表示为n点的机械载荷作用下的应力放大系数，p＝1～6，k＝1～6。
基于微观力学失效理论，通过对代表性体积单元模型设置一些参考点，用应力放大系数 法提取纤维和基体的微观应力。如图2所示，根据超景深电子显微镜对单向复合材料层合板的界面观测，采用长方形结构代表性体积单元模型的交错型排列更接近于纤维随机分布情况，且长方形模型(图2c)相比于正方形的代表性体积单元模型(图2b)更准确，故本申请将采用长方形代表性体积单元模型进行分析。
对代表性体积单元模型加载周期性边界条件，其中代表性体积单元模型中相对面上每对节点的位移关系满足表达式(4)：
式(4)中：
和
分别表示代表性体积单元模型中x，y，z平面上节点沿j方向法向量边界面上的i方向位移；上标“+”和“-”表示相对边界面上的节点对；
是相对边界面上节点对在i方向上位移的差值，主要由代表性体积单元的平均应变确定；假定代表性体积单元的平均应变已确定，则
变为常数。
代表性体积单元模型中各组分结构选择的代表性关键点如图2所示，其中纤维选取10个参考点(F1～F10)，基体选择15个参考点(M1-M15)。
式(6)中：ζ
i(i＝1,2,3,4,5,6)表示代表性体积单元模型中每个关键参考点在六种宏观载荷工况下产生的微观应力矩阵；
表示代表性体积单元模型中六种宏观均匀应力载荷工况。其中ζ
i和
的计算过程如下式所示：
式(7)中：σ
i(i＝11,22,33,12,23,13)表示代表性体积单元模型在施加六种宏观应力载荷下的分析后得到微观应力，
表示代表性体积单元模型施加的六种宏观应力载荷，其中11,22,33方向表示三种沿x，y，z方向宏观拉伸，12,23,13方向表示沿z，x，y平面的三种宏观剪切，如图3所示。以x方向拉伸为例，对代表性体积单元模型施加单位宏观拉伸应力，宏观应力
矩阵等于[1 0 0 0 0 0]
T，根据模型计算结果提取关键点处的微观应力，在代表性体积单元模型中设置的关键参考点获取的应力放大系数为
基于此原理，分别对代表性体积单元模型施加六种不同的数值宏观载荷即可分别获取应力放大系数。
由于建立的钻削多尺度有限元模型是三维状态，涉及六个应力分量：
式(8)中：f
i,f
ij表示在各个i方向和平面ij的强度张量i,j＝1,2,3(1表示纤维纵向)；σ
i表示纤维在微观尺度上法向应力，其中各强度张量的表达式如下：
式(9)中：X
f，Y
f，Z
f表示在纤维在微观尺度上各方向的失效强度(上标T表示拉伸状态，C表示压缩状态)，
表示纤维在微观尺度上的各方向上的剪切强度(i,j＝X,Y,Z)。
针对纤维在X方向主要在拉压状态下的损伤演化规律可表示为：
式(10)中：
表示在代表性体积单元模型中纤维上所有参考点定义的损伤系数(上标T表示拉伸状态，C表示压缩状态)；
表示当损伤变量达1时(即纤维失效)的最大应变；ε
11表示实际应变；
可分别表示初始损伤应变；其中获取各类型损伤应变的表达式如下：
式(11)，(12)中：
表示纤维断裂应变能释放率临界值；
和
表示纤维在微观尺度上在纵向拉伸强度和压缩强度。L表示纤维细观单元特征长度，主要由单元的尺寸以及构造形式相关，复合材料模型中纤维所有单元均采用六面体单元类型，其可能发生破坏的区域均具有相同的面内长度，因而单元特征长度计算的方程可按如下表示：
式(13)中：L
initial表示有限元软件中获取的初始单元特征长度，在材料子程序中通过VUMAT内嵌“CharLength”函数控制；l
z表示单元体在全局坐标系下沿坐标系z方向的尺寸，如图4(a)所示。
对于拉伸模型，当
时，材料在1方向上彻底失效。对于压缩模型，纤维压缩是由微屈曲机制引起的，如扭折带，当纤维单元断裂后依旧存有部分残余承载能力
具体失效模式如图4(b)所示。
纤维的本构关系在微观尺度下损伤定义为：
式(14)中：σ
f、C
f、ε
f分别表示纤维在微观尺度下的应力、刚度和应变；
和
分别表示纤维在微观尺度下纵向和横向损伤弹性模量；E
f1和E
f2分别表示纤维在微观尺度下纵向和横向没有损伤弹性模量。
基体屈服失效准则表达式具体如下：
式(15)中：T
mi和C
mi分别表示基体的初始拉伸强度和压缩强度；σ
VM和I
1分别表示·Von Mises等效应力和第一应力不变量，可表示为：
式(16)中：σ
mi为基体的在微观尺度下的主应力(i＝1～3)和剪切应力(i＝4～6)。
考虑基体在拉压状态下强度不一致，基于微观力学失效理论，采用Stassi等效应力σ
eq对基体的损伤演化进行定义，其表达式为：
式(17)中：β表示基体的初始压缩强度与拉伸强度比，即β＝C
mi/T
mi。
在代表性体积单元模型中基体的损伤根据应力状态进行确定，损伤值的大小在0～1之间，在初始损伤判定后，其损伤演化规律
和等效应力可表达为：
f(σ
eq,κ
m)＝σ
eq-κ
m (19)
式(18)，(19)中：
表示在代表性体积单元模型中基体上所有参考点定义的损伤系数(上标T表示拉伸状态，C表示压缩状态)；γ表示基于测试数据下的校准单轴应力-应变曲线损伤状态参数，其值是基于测试数据，然后根据宏观材料属性反向求解得到；κ
m表示在一个压缩强度T
mi标准值的历史变量，其初始损伤和更新需要根据Kuhn-Tucker加载-卸载条件下进行确定，即：
然后，根据式(18)，基体在压缩和拉伸状态下的损伤系数
可表示为：
假设损伤初始和演化规律都是各向同性的，基体本构关系在微观尺度下考虑刚度退化的损伤如图5所示，其表达式定义为：
式(22)中：σ
m、C
m、ε
m分别表示基体在微观尺度下的应力，刚度和应变；
和E
m分别表示基体在微观尺度下损伤弹性模量和没有损伤的弹性模量；
表示代表性体积单元模型中在微观尺度下所有基体上参考点的最大损伤系数。
当获取每个单元的损伤变量，并满足失效准则后。复合材料宏观本构-损伤模型的刚度矩阵
更新表达式如下：
式(23)中：d
L和d
T分别表示复合材料的纵向和横向宏观损伤变量；
表示复合材料 的宏观本构中未损伤刚度系数,其表达式为：
式(24)中：E
i，v
ij，G
ij(i,j＝1,2,3)分别代表复合材料的宏观弹性模量、泊松比、宏观剪切模量。
对于各损伤变量的评估计算表达式如下：
式(25)，(26)中：V
f代表纤维体积分数；V
m代表基体体积分数。
为了预防在有限元软件中因出现严重畸变单元而导致在计算时出现不收敛情况，本申请引入最大主应变及最小主应变准则实现过度畸变单元的辅助删除，基于极分解定理，任意物体的运动都可以将其分解为沿三个方向的纯拉伸和纯刚体旋转，即：
式(27)中：
表示应变梯度矩阵；
表征材料的纯拉伸变形；
表示材料的纯刚体旋转；
和
分别表示微元的初始位置和变形后的位置；
本申请复合材料的层合板模型均采用六面体单元模型，对
矩阵求解特征根，可以获取单元的3个方向主拉伸比λ
1，λ
2，λ
3，根据Biot应变(名义应变)定义分别求出3个主应变 值：
ε
i＝λ
i-1 (28)
式(28)中：ε
i(i＝1,2,3)表示单元的主应变；λ
i(i＝1,2,3)表示单元的主拉伸比。
因而，根据求出的主应变，辅助单元删除准则可表示为：
式(29)中：
和
分别表示单元最大主应变和最小主应变；
和
分别表示手动设定的最大主应变和最小主应变临界值。
另外，
和
的设定值需要结合对比真实状态下的钻削试验进行设定，
的不能太小，
取值不能太大，否则会将尚为失效的单元强制删除，最终无法实现钻削过程准确模拟，本申请采用的是碳纤维增强复合材料，将其值分别设定为
基于建立的微观力学失效理论失效准则和动态渐进损伤模型，多尺度显式有限元算法如图6所示。整个多尺度显式有限元算法的实现主要可分为三个步骤：
步骤(1)：在开始的全局时间分析步n时，通过前一个时间分析步
与应变增量步
求和计算全局总应变
通过前一个复合材料单元刚度矩阵
求解计算出每个积分点的宏观应力
步骤(2)：根据宏观应力对组分结构单元分别设置的关键参考点乘以
和
分别计算出纤维和基体的微观应力
和
基于建立的微观力学失效理论失效准则，确定纤维和基体的失效系数，通过建立的渐进损伤模型，得到纤维和基体的微观损伤变量；对代表性体积单元模型中的纤维损伤系数和基体损伤系数确定所有参考点
和
的最大值。
步骤(3)：通过建立新的损伤演化准则，每个积分点上的宏观损伤变量可以通过d
L和d
T 的最大值计算，实现单元刚度的退化，性能衰减等，并基于单元删除准则，实现单元的删除。由于在第n-1时间分析步宏观损伤变量是不可逆的，宏观损伤变量的选择较大值，并在n时间分析步时实现复合材料刚度矩阵
的更新，如果达到设定的时间分析步，将进入下一个时间分析步，如果时间分析步终止，则程序终止。
裂纹萌生和扩展由表面牵引力-分离位移之间的本构关系模型控制，其失效损伤演化规律表达式为：
式(30)中：T
i(i＝n,s,t)表示各方向的名义应力，n，s，t分别表示法线方向和局部剪切方向；δ
i表示对应方向的名义应变，其值为界面处分离位移除以内聚力单元的初始厚度；K
ii表示内聚力单元本构的刚度；
表示层间损伤变量；
表示在加载过程中的初始失效时的等效位移；
表示在加载过程中等效位移的最大值，并且在历程分析中是不可逆转的；
表示损伤萌生的等效应力；其中
的分别可表示为：
式(31)，(33)，(34)中：
表示各方向的最大名义应力；
表示层间断裂韧性，即断裂能，具体如图7(c)所示。
牵引力-分离失效损伤演化规律图7(b)所示，其损伤失效演化规律表达式为：
其中混合模式下的损伤状态变量D
s可表示为：
式(36)中：
表示在加载过程中等效位移的最大值；
表示混合模型在加载过程中的等效位移；
表示混合模式下损伤萌生的等效应力；其中，
的表达式分别为：
式(37)、(38)、(39)中：<·>表示麦考利括号；η表示混合模式下的B-K准则参数，本申请的研究对象是碳纤维增强复合材料，取值为1.45；β表示混合比，其表达式为：
实际应力和最大应力的关系可以表示为：
在有限元软件的计算过程中，如果当某个内聚力单元的截面积分点的损伤变量D
s达到
该单元将被删除，ABAQUS软件中无量纲的刚度退化系数SDEG和模型中的损伤变量D
s具有相同的含义，即当SDEG＝1时，判断单元失效并删除。
纤维和基体的各项材料参数如表1所示。
表1 纤维丝和基体的材料参数
对代表性体积单元模型进行网格划分时，在各对称面上设置的单元节点在位置和数量上需保持完全一致，并采用映射的方式生成网格，如图9(a)所示。为避免代表性体积单元模型中顶点和边上的网格节点出现过约束的情况，在对节点集施加周期性边界条件时，必须将顶点和边上的节点从面中去除，同时必须将顶点上的节点从边上去除，一个顶点仅包含一个节点，但所有的面和边上的每个节点集合可以同时包含多个节点。
在全局坐标系下，将建立的代表性体积单元模型的各边长定义为W
x,W
y,和h，坐标原点位于D，在六种宏观应力载荷
下，根据式(4)的建立的周期性边界条件下通过以下线性约束方程组实现顶点、棱边和平面的响应载荷工况的加载。
在ABAQUS/Standard分析步下对代表性体积单元模型分别施加六种工况(三轴宏观拉伸及剪切)下的各单位应力载荷，得到微观应力分布，如图10所示。采用10红色实心圆表示 微观纤维关键参考点和15个紫色实心圆表示基体的微观基体的参考点计算获取的应力放大系数，如图2所示。F10参考点和M10～M15点分别位于纤维内部中心位置和基体内部中心位置。F1～F9参考点和M1～M9考虑纤维和基体结合处的位置关键参考点。由于施加的是外在为1MPa的应力载荷，故获取各关键参考点的微观应力分布即为应力放大系数，不同工况下的纤维和基体上各关键参考点的应力如图10、图11所示。由于材料的固有属性不会发生变化，代表性体积单元模型只需进行求解一次并提取关键参考点结果。在获取应力放大系数后将其以一个独立的文件存储在Fortran代码中，在宏观复合材料钻削有限元模型中，基于ABAQUS VUMAT代码使用“#include”函数调用该文件，实现复合材料材料的动态渐进失效多尺度分析。
本申请采用是T700S-12K/YP-H26的耐高温环氧树脂复合材料为研究对象，纤维体积分数含量约为59％，采用[(0°/90°/45°/-45°)
s]
4铺层顺序，复合材料层合板结构总体尺寸为15mm×15mm×5.76mm，总计32层，在预设的全局坐标下，对单层复合材料复合材料赋予沿纤维方向上的正交各向异性材料，其单层复合材料的参数如表2所示，内聚力单元的参数如表3所示。将直径为6mm的硬质合金匕首钻(AlTiN涂层)设为离散刚体模型，为完整且准确的预测钻削过程中的轴向力和扭矩等参数，将匕首钻的总质量和转动惯量约束在顶部参考点位置，并基于该参考点施加钻头的进给速度和转速等边界条件。同时为提升计算效率及节约计算机资源，对涉及复合材料层合板钻削接触区域(直径为8mm圆形区域)进行网格细化处理，对接近层合板边缘区域和匕首钻均采用较粗糙的网格。复合材料采用C3D8R单元(8节点、三维减缩积分单元)，分层区域使用COH3D8粘性单元(0mm厚度)，匕首钻采用C3D10M离散刚体单元。为了模拟出钻削过程中产生的毛刺现象，其复合材料单元的模型网格尺寸保证为细观尺度，根据对单元网格灵敏性模型的优化分析结果，整个钻削有限元模型中总共有4203438个单元，其中复合材料模型包含3382122个六面体单元，分层区域包含819108个内聚力单元和匕首钻包含8208个四面弹单元，复合材料单元最小尺寸约为 45×45×40μm。
表2 多向铺层结构下UD-CFRP的材料参数
表3 CEs材料属性
根据钻头沿轴线方向进给复合材料的实际工况，对整个钻削有限元模型中施加转速和进给速度等边界条件，由于将匕首钻模型的运动状态约束在顶部参考点处，对参考点在X和Y方向上的位移进行限制，在Z方向上施加进给速度，同理，在对X和Y方向上的转速进行限制，在Z方向施加顺时针方向转动速度，并对复合材料层合板的四个垂直表面进行固定。在仿真计算结束后，本申请通过搭建的实验平台对模型进行验证，实验平台主要包含CFRP钻削、数据采集、实验观测三个系统，实验方案示意图及测试装置连接如图13所示。
根据图14所示的实验与仿真对比分析结果，钻削多尺度有限元模型与试验的轴向力和扭矩在钻削的各个阶段值及变化趋势趋近一致，仅在第T
IV及T
V阶段时，轴向力和扭矩的减少速度相比于实验较为迟滞，主要原因是在多尺度有限元模型中纤维采用基于断裂韧性的退化方式，其设定的耗散能值是基于经验评估得到，加上由于模拟层间分层现象中粘结区域CEs的作用，使得模型中的有些宏观单元没有删除，导致单元与钻头接触时间过长，进而产生部 分较小的力与扭矩。
从图15(a)中可以清晰观测出，CFRP孔壁主要由树脂涂覆表面和纤维断口表面组成，靠近入口的孔壁涂覆表面比较明显。CFRP从钻削的入口到出口处的损伤主要表现为入口处劈裂、出口处撕裂、毛边、层间分离、径向挤伤、微裂纹等损伤。根据图15(a)与(b)中每层损伤的位置对比分析，发现出现凹坑和孔壁损伤位置几乎一致。因而，本文建立的多尺度有限元模型可以较真实的模拟出CFRP钻削的孔壁损伤状态。
CFRP在制孔过程中在出入口容易产生毛刺，撕裂、分层等加工损伤，由于匕首钻是一种钻扩孔一体的CFRP专用钻头，仅对CFRP层合板钻削时孔的入口处纤维断口较为平整，入口撕裂现象不明显及产生的毛刺较少，但是通过有限元分析模型和试验结果观测，在出口处出现较明显毛刺、撕裂，分层损伤，如图16(a)所示。但由于在多尺度有限元模型中是基于宏观单元的删除，实现材料去除的模拟，在单元的尺度效应，仅能够模拟出其宏观层面的损伤现象，无法实现毛刺及撕裂损伤的量化。CFRP制孔出现的分层损伤一般仅存在出口侧面的几层材料之间，其内部的损伤大小需要采用超声扫描等观测设备进行观察，但对于出口处的分层会沿着匕首钻进给方向往材料外部扩展，在损伤区域内有一部分隆起的高度，仅需采用超景深显微镜观测即可，如图16(a)所示。
分层损伤系统是以刀具直径为基准，通过孔周边删除单元数量计算面积除以基准孔的面积计算出分层损伤系数，测量方法是在ABAQUS软件里面导出图片，然后导入CAD软件，基于软件自带工具采用等比例法分别实现测量。由于钻头在实际工况下转速太快，而多尺度有限元模型单元的删除是不可恢复的，故取分层区域中心区域产生损伤单元面积用于计算分层损伤系数，如图24(b)所示。因而，分层损伤系数可表示为:
式(47)中，D
d表示出口分层损伤系数，A
nom表示中心孔径区域面积，A
d表示分层产生损伤单元面积。根据计算出多尺度有限元模型和试验工况中不同加工参数下的钻削CFRP出口处的 分层损伤系数。从图17中可以看出，有限元模型和试验结果在分析损伤系数的结果较为接近，同时上机随加工参数变化，其变化趋势也趋于一致。
综上所述，本文建立的多尺度有限元模型能够真实的模拟出实际工况下的CFRP损伤状态，轴向力、扭矩及出口分层损伤系数最大误差分别为3.37％，7.69％、4.28％，同时也可以较真实的模拟出孔壁表面损伤及出口分层等现象。
Claims (6)
- 一种复合材料多尺度模型的动态渐进失效分析方法，其特征在于包含以下步骤：第一步，忽略复合材料层合板在工艺制备中产生的初始损伤，将复合材料理想化为单向层合板预设的铺层序列叠加而成，对复合材料的宏观结构进行分析，考虑宏观结构设计，材料属性，铺层方式以及整体结构受载形式对结构内部的应力和应变分布的影响；第二步，根据单向层合板的力学分析，对单向板纤维分布进行合理简化，基于纤维体积分数建立面向复合材料结构力学性能表观的微观力学代表性体积单元，并采用多尺度分析方法，将宏观结构单元或积分点的载荷情况传递到微观模型中，在该载荷下模拟纤维和基体组分材料各单元积分点的应力分布；另外，采用微观失效理论分别对基体，纤维进行失效判断，假设各对应单元发生损伤，则对损伤单元进行刚度折减、失效删除，直至删除，如果单元没有删除，则对整个单胞模型进行均匀化处理，得到损伤后的材料宏观力学性能；第三步，将损伤后的材料宏观性能赋予宏观结构中的对应单元中，对整体结构再进行下一迭代步的宏观分析，如此反复迭代直至设定时间分析步完成，最终实现复合材料由组分到结构件的在刀具切削作用下的多尺度数值模拟分析。
- 根据权利要求1所述的复合材料多尺度模型的动态渐进失效分析方法，其特征在于复合材料层合板结构的宏观应力与各组分结构微观应力通过应力放大系数进行桥联，而应力放大系数的获取可以通过包含纤维、基体的代表性体积单元的有限元分析结果得到。
- 根据权利要求1所述的复合材料多尺度模型的动态渐进失效分析方法，其特征在于包含以下有限元算法：步骤(1)：在开始的全局时间分析步n时，通过前一个时间分析步与应变增量步
求和计算全局总应变
通过前一个复合材料单元刚度矩阵
求解计算出每个积分点的宏观应力
步骤(2)：根据宏观应力对组分结构单元分别设置的关键参考点乘以纤维的应力方法系数和基体的应力方法系数
分别计算出纤维和基体的微观应力
和
基于建立 的微观力学失效理论失效准则，确定纤维和基体的失效系数，通过建立的渐进损伤模型，得到纤维和基体的微观损伤变量；对代表性体积单元模型中的纤维损伤系数和基体损伤系数确定所有参考点
和
的最大值；
步骤(3)：通过建立新的损伤演化准则，每个积分点上的宏观损伤变量可以通过d L和d T的最大值计算，实现单元刚度的退化，性能衰减，并基于单元删除准则，实现单元的删除；由于在第n-1时间分析步宏观损伤变量是不可逆的，宏观损伤变量的选择较大值，并在n时间分析步时实现复合材料刚度矩阵的更新，如果达到设定的时间分析步，将进入下一个时间分析步，如果时间分析步终止，则程序终止。
- 根据权利要求1所述的复合材料多尺度模型的动态渐进失效分析方法，其特征在于宏-微观转化方程的表达式如式(1)、(2)、(3)所示：
式(1)～(3)中：σ i表示纤维或基体n点处的应力，n为受力分析点，i＝1～6；表示复合材料层合板的宏观应力分量，j＝1～6；
表示为n点的机械载荷作用下的应力放大系数，p＝1～6，k＝1～6。
- 根据权利要求1所述的复合材料多尺度模型的动态渐进失效分析方法，其特征在于单元特征长度计算的方程可按如下表示：
式(4))中：L initial表示有限元软件中获取的初始单元特征长度，在材料子程序中通过VUMAT内嵌“CharLength”函数控制；l z表示单元体在全局坐标系下沿坐标系z方向的尺寸。 - 根据权利要求1所述的复合材料多尺度模型的动态渐进失效分析方法，其特征在于为了预防在有限元软件中因出现严重畸变单元而导致在计算时出现不收敛情况，本申请引入最大主应变及最小主应变准则实现过度畸变单元的辅助删除，基于极分解定理，任意物体的运动都可以将其分解为沿三个方向的纯拉伸和纯刚体旋转，即：
式(5)中：表示应变梯度矩阵；
表征材料的纯拉伸变形；
表示材料的纯刚体旋转；
和
分别表示微元的初始位置和变形后的位置；
当获取每个单元的损伤变量，并满足失效准则后，复合材料宏观本构-损伤模型的刚度矩阵更新表达式如下：
式(26)中：d L和d T分别表示复合材料的纵向和横向宏观损伤变量；表示复合材料的宏观本构中未损伤刚度系数,其表达式为：
式(27)中：E i，v ij，G ij(i,j＝1,2,3)分别代表复合材料的宏观弹性模量、泊松比、宏观剪切模量；复合材料宏观本构-损伤模型演化规律中的损伤变量直接由代表性体积单元模型中确定，通过损伤组分失效准则计算微观组分退化弹性参数实现损伤评估，其中对于各损伤变量的评估计算表达式如下：
式(28)，(29)中：V f代表纤维体积分数；V m代表基体体积分数，E m表示基体在微观尺度下没有损伤的弹性模量，和
分别表示纤维在微观尺度下纵向和横向损伤弹性模量；E f1和E f2分别表示纤维在微观尺度下纵向和横向没有损伤弹性模量，上标T代表拉伸状态，C代表压缩状态，
表示基体在拉伸和压缩状态下的损伤状态变量。
Priority Applications (2)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| US17/595,100 US20220284150A1 (en) | 2020-07-16 | 2020-07-16 | Dynamic Progressive Failure Analysis Method For Composite Multi-Scale Model |
| PCT/CN2020/102477 WO2021139130A1 (zh) | 2020-07-16 | 2020-07-16 | 复合材料多尺度模型的动态渐进失效分析方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| PCT/CN2020/102477 WO2021139130A1 (zh) | 2020-07-16 | 2020-07-16 | 复合材料多尺度模型的动态渐进失效分析方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| WO2021139130A1 true WO2021139130A1 (zh) | 2021-07-15 |
Family
ID=76787423
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| PCT/CN2020/102477 WO2021139130A1 (zh) | 2020-07-16 | 2020-07-16 | 复合材料多尺度模型的动态渐进失效分析方法 |
Country Status (2)
| Country | Link |
|---|---|
| US (1) | US20220284150A1 (zh) |
| WO (1) | WO2021139130A1 (zh) |
Cited By (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114861494A (zh) * | 2022-04-25 | 2022-08-05 | 西北工业大学 | 考虑热解碳各向异性特征的c/c复合材料弹性性能预测方法 |
| CN114970170A (zh) * | 2022-05-31 | 2022-08-30 | 西北工业大学 | 一种考虑固化缺陷的Z-pin增韧复合材料强度预测方法 |
| CN116644619A (zh) * | 2023-07-27 | 2023-08-25 | 北京理工大学 | 纤维预制体成型工艺的压力优化方法、装置、设备及介质 |
| CN116757128A (zh) * | 2023-08-22 | 2023-09-15 | 中国人民解放军火箭军工程大学 | 强度极限线的获取和推进剂的失效判断的方法及相关设备 |
| CN117371272A (zh) * | 2023-09-22 | 2024-01-09 | 天津大学 | 适用于不同各向异性材料及尺寸的夹持式单边缺口拉伸试样裂纹长度及断裂性能计算方法 |
| CN117875139A (zh) * | 2024-03-13 | 2024-04-12 | 中南大学 | 无砟轨道多尺度损伤演变分析方法及系统 |
Families Citing this family (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115290447A (zh) * | 2022-08-10 | 2022-11-04 | 香港理工大学深圳研究院 | 一种超螺旋状复合纤维驱动器的模量分析测试方法 |
| CN115906570B (zh) * | 2022-11-22 | 2023-07-11 | 大连理工大学 | 一种基于数据驱动的复合材料开孔层合板渐进损伤预测方法 |
| CN116644617B (zh) * | 2023-07-27 | 2023-12-22 | 北京航空航天大学 | 基于纤维计算晶粒的航空纤维复合材料细观力学分析方法 |
| CN117610344A (zh) * | 2023-11-10 | 2024-02-27 | 中山大学 | 一种飞行汽车海上迫降结构损伤评价方法 |
| CN117727409B (zh) * | 2024-02-18 | 2024-05-10 | 西北工业大学 | 孔隙随机分布的c/c复合材料的强度预测方法及装置 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106777769A (zh) * | 2017-01-08 | 2017-05-31 | 浙江大学 | 预测低速冲击下复合材料多层厚板渐进失效的有限元方法 |
| CN109241650A (zh) * | 2018-09-25 | 2019-01-18 | 南京航空航天大学 | 基于跨尺度仿真的碳纤维增强复合材料力学性能预测方法 |
| US20190278760A1 (en) * | 2008-11-14 | 2019-09-12 | Georgetown University | Process and Framework For Facilitating Information Sharing Using a Distributed Hypergraph |
| CN110688789A (zh) * | 2019-08-29 | 2020-01-14 | 西安建筑科技大学 | 一种碳纤维复合单向板偏载冲击的多尺度响应分析方法 |
| CN111832209A (zh) * | 2020-07-16 | 2020-10-27 | 南京航空航天大学 | 复合材料多尺度模型的动态渐进失效分析方法 |
-
2020
- 2020-07-16 US US17/595,100 patent/US20220284150A1/en active Pending
- 2020-07-16 WO PCT/CN2020/102477 patent/WO2021139130A1/zh active Application Filing
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20190278760A1 (en) * | 2008-11-14 | 2019-09-12 | Georgetown University | Process and Framework For Facilitating Information Sharing Using a Distributed Hypergraph |
| CN106777769A (zh) * | 2017-01-08 | 2017-05-31 | 浙江大学 | 预测低速冲击下复合材料多层厚板渐进失效的有限元方法 |
| CN109241650A (zh) * | 2018-09-25 | 2019-01-18 | 南京航空航天大学 | 基于跨尺度仿真的碳纤维增强复合材料力学性能预测方法 |
| CN110688789A (zh) * | 2019-08-29 | 2020-01-14 | 西安建筑科技大学 | 一种碳纤维复合单向板偏载冲击的多尺度响应分析方法 |
| CN111832209A (zh) * | 2020-07-16 | 2020-10-27 | 南京航空航天大学 | 复合材料多尺度模型的动态渐进失效分析方法 |
Cited By (12)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114861494A (zh) * | 2022-04-25 | 2022-08-05 | 西北工业大学 | 考虑热解碳各向异性特征的c/c复合材料弹性性能预测方法 |
| CN114861494B (zh) * | 2022-04-25 | 2024-04-02 | 西北工业大学 | 考虑热解碳各向异性特征的c/c复合材料弹性性能预测方法 |
| CN114970170A (zh) * | 2022-05-31 | 2022-08-30 | 西北工业大学 | 一种考虑固化缺陷的Z-pin增韧复合材料强度预测方法 |
| CN114970170B (zh) * | 2022-05-31 | 2024-02-23 | 西北工业大学 | 一种考虑固化缺陷的Z-pin增韧复合材料强度预测方法 |
| CN116644619A (zh) * | 2023-07-27 | 2023-08-25 | 北京理工大学 | 纤维预制体成型工艺的压力优化方法、装置、设备及介质 |
| CN116644619B (zh) * | 2023-07-27 | 2023-09-22 | 北京理工大学 | 纤维预制体成型工艺的压力优化方法、装置、设备及介质 |
| CN116757128A (zh) * | 2023-08-22 | 2023-09-15 | 中国人民解放军火箭军工程大学 | 强度极限线的获取和推进剂的失效判断的方法及相关设备 |
| CN116757128B (zh) * | 2023-08-22 | 2023-11-03 | 中国人民解放军火箭军工程大学 | 强度极限线的获取和推进剂的失效判断的方法及相关设备 |
| CN117371272A (zh) * | 2023-09-22 | 2024-01-09 | 天津大学 | 适用于不同各向异性材料及尺寸的夹持式单边缺口拉伸试样裂纹长度及断裂性能计算方法 |
| CN117371272B (zh) * | 2023-09-22 | 2024-04-19 | 天津大学 | 适用于不同各向异性材料及尺寸的夹持式单边缺口拉伸试样裂纹长度及断裂性能计算方法 |
| CN117875139A (zh) * | 2024-03-13 | 2024-04-12 | 中南大学 | 无砟轨道多尺度损伤演变分析方法及系统 |
| CN117875139B (zh) * | 2024-03-13 | 2024-05-24 | 中南大学 | 无砟轨道多尺度损伤演变分析方法及系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| US20220284150A1 (en) | 2022-09-08 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| WO2021139130A1 (zh) | 复合材料多尺度模型的动态渐进失效分析方法 | |
| CN111832209B (zh) | 复合材料多尺度模型的动态渐进失效分析方法 | |
| Mei et al. | A novel fabrication method and mechanical behavior of all-composite tetrahedral truss core sandwich panel | |
| Du et al. | Open-hole tensile progressive damage and failure prediction of carbon fiber-reinforced PEEK–titanium laminates | |
| Heimbs | Foldcore sandwich structures and their impact behaviour: an overview | |
| Liu et al. | Recent developments on damage modeling and finite element analysis for composite laminates: A review | |
| He et al. | Low-velocity impact behavior of X-Frame core sandwich structures–experimental and numerical investigation | |
| Jiang et al. | Flexural performances of fiber face-sheets/corrugated core sandwich composite structures reinforced by horizontal stiffeners | |
| CN109101692B (zh) | 基于最大应力准则的复合材料层合板极限载荷计算方法 | |
| CN113420376B (zh) | 基于多尺度的碳纤维复合材料抗冲击力学性能仿真方法 | |
| CN112528415B (zh) | 一种复合材料轴结构宏-细观失效模式分析方法 | |
| McElroy | Use of an enriched shell finite element to simulate delamination-migration in a composite laminate | |
| Mishra et al. | An experimental approach to free vibration response of woven fiber composite plates under free-free boundary condition | |
| Al Rifaie | Resilience and toughness behavior of 3D-printed polymer lattice structures: Testing and modeling | |
| Liu et al. | Scale-span modelling of dynamic progressive failure in drilling CFRPs using a tapered drill-reamer | |
| Zhu et al. | Experimental study and modeling analysis of planar compression of composite corrugated, lattice and honeycomb sandwich plates | |
| CN113158508A (zh) | 一种复合材料层合板贴补修复后极限载荷确定方法 | |
| CN112926244A (zh) | 一种复合材料层合板开孔件极限载荷确定方法 | |
| Kang et al. | Parametric study on the buckling load after micro-bolt repair of a composite laminate with delamination | |
| Jäger et al. | A discrete model for simulation of composites plate impact including coupled intra-and inter-ply failure | |
| Mishnaevsky Jr | Micromechanical modelling of wind turbine blade materials | |
| CN112906263B (zh) | 含制孔分层损伤的复合材料层合板强度预测方法 | |
| McElroy | An enriched shell element for delamination simulation in composite laminates | |
| Kress et al. | Honeycomb sandwich residual stress deformation pattern | |
| Rammerstorfer et al. | Buckling in thin walled micro and meso structures of lightweight materials and material compounds |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| 121 | Ep: the epo has been informed by wipo that ep was designated in this application |
Ref document number: 20912552 Country of ref document: EP Kind code of ref document: A1 |
|
| NENP | Non-entry into the national phase |
Ref country code: DE |
|
| 122 | Ep: pct application non-entry in european phase |
Ref document number: 20912552 Country of ref document: EP Kind code of ref document: A1 |
|
| 122 | Ep: pct application non-entry in european phase |
Ref document number: 20912552 Country of ref document: EP Kind code of ref document: A1 |