CN113609735A - Geometric-mechanical model integration method for weaving composite material - Google Patents

Geometric-mechanical model integration method for weaving composite material Download PDF

Info

Publication number
CN113609735A
CN113609735A CN202110906645.1A CN202110906645A CN113609735A CN 113609735 A CN113609735 A CN 113609735A CN 202110906645 A CN202110906645 A CN 202110906645A CN 113609735 A CN113609735 A CN 113609735A
Authority
CN
China
Prior art keywords
warp
unit
yarn
section
model
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110906645.1A
Other languages
Chinese (zh)
Other versions
CN113609735B (en
Inventor
高希光
张煦
宋迎东
张盛
于国强
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing University of Aeronautics and Astronautics
Original Assignee
Nanjing University of Aeronautics and Astronautics
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing University of Aeronautics and Astronautics filed Critical Nanjing University of Aeronautics and Astronautics
Priority to CN202110906645.1A priority Critical patent/CN113609735B/en
Publication of CN113609735A publication Critical patent/CN113609735A/en
Application granted granted Critical
Publication of CN113609735B publication Critical patent/CN113609735B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

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
    • 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/14—Force analysis or force optimisation, e.g. static or dynamic forces
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Woven Fabrics (AREA)

Abstract

The invention discloses a geometric-mechanical model integration method for weaving a composite material, which comprises the following steps: on the basis of a ply yarn path model, setting every 4 adjacent interweaving points to form a three-dimensional rhombic unit cell, dividing each unit cell into two Delaunay triangles, and constructing to obtain a Delaunay triangular grid; establishing a Thiessen polygonal mesh based on the Delaunay triangular mesh; obtaining a main direction of the unit according to a finite element model of the Thiessen polygonal grid and the whole area; selecting a complex geometric area of the ply yarn path model for local refining to generate a local refining finite element model; and judging the material properties of each unit in the local refined finite element model, wherein the material properties comprise warps, wefts and a matrix. The method can be used for quickly judging the material properties of each unit in the finite element model.

Description

Geometric-mechanical model integration method for weaving composite material
Technical Field
The invention relates to the technical field of pretreatment of mechanical property finite element analysis of plain weave composite materials, in particular to a geometric mechanical model integration method for weaving composite materials.
Background
The Ceramic Matrix Composites (CMC) reinforced by fiber weaving has the excellent characteristics of high specific strength, high specific modulus, high reliability, strong high temperature resistance and the like, has great potential to replace metal materials to become core materials of hot end parts of future high-performance aviation gas turbine engines, and can improve the efficiency of the turbine engines revolutionarily. The woven CMC not only can maintain good high-temperature performance of the ceramic material, but also has better fracture toughness due to the reinforcing effect of the fiber. However, the introduction of the new material is greatly hindered by the high manufacturing cost, and the mechanical property of the woven CMC structure needs to be accurately calculated by an effective performance prediction method, so as to achieve the purposes of reducing the number of iterative tests and reducing the manufacturing cost and the development period.
In previous work, researchers have generally used a periodic method to predict the mechanical properties of CMC structures. However, CMC components in engineering applications usually have complex preform structures, and it is not suitable to perform mechanical simulation on complex CMC structures with high aperiodicity by using a periodic method, and before predicting the mechanical properties of a woven CMC structure, a mechanical solving method capable of considering the non-periodic mesoscopic structure of the whole CMC component needs to be established. In mesomechanics, both the fibres and the matrix of the composite are well established, which effectively reflects the non-periodic nature of the composite structure. However, at present, mesomechanics methods are only applied on the unit cell level. Due to the limitations at the computer level today, it can take a lot of time to perform mesomechanical analysis of composite materials at the component level, especially for geometrically complex composite components.
Therefore, it is necessary to design a method that combines the advantages of the computation speed of the homogenization method and the computation accuracy of the mesomechanics method, and the determination of the material properties of each element in the mechanical finite element model is a key step of the method.
Disclosure of Invention
Aiming at the defects in the prior art, the invention provides a geometric and mechanical model integration method for weaving composite materials, which is used for the pretreatment of mechanical property simulation of components of the woven composite materials. The method takes a composite material prefabricated body geometric model and a finite element mechanical model of a component as input, and realizes the integration of geometry and mechanics aiming at a structural uniform region and a non-uniform region of the component respectively to obtain the material property of each unit in the finite element model of the component.
In order to achieve the purpose, the invention adopts the following technical scheme:
a geometric-mechanical model integration method for woven composites, the material property simulation method comprising the steps of:
s1, inputting a composite material prefabricated body geometric model and a finite element mechanical model of the component; the geometric model of the composite material preform comprises a layering yarn path model and a preform mesoscale model; the finite element mechanical model comprises a finite element model of the whole area and a finite element model of the local non-uniform area;
s2, on the basis of the ply yarn path model, setting every 4 adjacent interweaving points to form a three-dimensional rhombic unit cell, dividing each unit cell into two Delaunay triangles, and constructing to obtain a Delaunay triangular grid;
s3, establishing a Thiessen polygonal mesh based on the Delaunay triangular mesh;
s4, obtaining a main direction of the unit according to a finite element model of the Thiessen polygonal grid and the whole area;
s5, selecting a complex geometric area of the ply yarn path model for local refining to generate a local refining finite element model;
and S6, judging the material properties of each unit in the local refined finite element model, wherein the material properties comprise warp yarns, weft yarns and a matrix.
In order to optimize the technical scheme, the specific measures adopted further comprise:
further, in step S2, the process of constructing the Delaunay triangular mesh based on the ply yarn path model includes the following steps:
s21, settingThe warp base line and the weft base line divide the yarn path model into intersection points
Figure BDA0003201913970000021
4 quadrants as the center, each 4 adjacent interweaving points form a three-dimensional rhombic unit cell, and 4 vertexes of the three-dimensional rhombic unit cell are set as
Figure BDA0003201913970000022
And
Figure BDA0003201913970000023
i is the coordinate value of the warp yarn constituting the three-dimensional rhombohedral unit cell, and j is the coordinate value of the weft yarn constituting the three-dimensional rhombohedral unit cell;
s22, selecting a three-dimensional rhomboid unit cell, and calculating the corresponding side
Figure BDA0003201913970000024
And edge
Figure BDA0003201913970000025
Angle therebetween
Figure BDA0003201913970000026
S23, judging the angle
Figure BDA0003201913970000027
Whether it is an obtuse angle: if the angle is obtuse, the diagonal line is used
Figure BDA0003201913970000028
Divide the unit cell into two Delaunay triangles
Figure BDA0003201913970000029
And
Figure BDA00032019139700000210
if not obtuse, by diagonal
Figure BDA00032019139700000211
Dividing unit cellIs two Delaunay triangles
Figure BDA00032019139700000212
S24, repeating steps S22 to S23, dividing the entire yarn path model into Delaunay triangular meshes.
Further, in step S3, the process of building Thiessen polygon meshes based on the Delaunay triangulation mesh includes the following steps:
s31, selecting any one warp and weft yarn intersection point
Figure BDA00032019139700000213
Arranged counterclockwise, with
Figure BDA00032019139700000214
Denotes that k is 0,1,2, …, kTWherein k isTAs the crossing point of warp and weft yarns
Figure BDA00032019139700000215
The number of the neighborhood triangles is Delaunay triangle containing the crossing point of the warp and weft; i is the coordinate value of the warp yarn, j is the coordinate value of the weft yarn;
s32, calculating the crossing point of the warp and weft yarns
Figure BDA00032019139700000216
The center of the circumscribed circle of each neighborhood triangle is recorded as
Figure BDA00032019139700000217
S33, connected in a counterclockwise order
Figure BDA00032019139700000218
Obtaining the intersection point of the warp and weft
Figure BDA00032019139700000219
A corresponding Thiessen polygon;
and S34, repeating the steps S31 to S33 until obtaining Thiessen polygons corresponding to all warp and weft yarn intersections.
Further, in step S4, the process of obtaining the main directions of the elements according to the Thiessen polygon meshes and the finite element model of the whole area includes the following steps:
s41, assuming that the node number of any unit A in the finite element model of the whole area is n1The unit node is denoted as Pl，l＝1,2,…,n1And calculating the coordinates of the center point of the unit A according to the following formula:
Figure BDA0003201913970000031
s42, traversing the polygons in the Thiessen polygon mesh;
s43, the number of the traversed Thiessen polygon vertexes is m1The vertex is
Figure BDA0003201913970000032
m＝1,2,…m1The crossing point of the corresponding warp and weft yarns is
Figure BDA0003201913970000033
The ply is at the point
Figure BDA0003201913970000034
Normal vector at point is
Figure BDA0003201913970000035
i is the coordinate value of the warp yarn, j is the coordinate value of the weft yarn; calculating the direction vector of each side of the Thiessen polygon by adopting the following formula:
Figure BDA0003201913970000036
s44, calculating the center point P of the cell A by the following formulacenterAt the passing point
Figure BDA0003201913970000037
And by vectors
Figure BDA0003201913970000038
Projected points on a plane as a normal vector:
Figure BDA0003201913970000039
s45, calculating the point by the following formula
Figure BDA00032019139700000322
Direction vectors to vertices of Thiessen polygons:
Figure BDA00032019139700000310
s46, by judging
Figure BDA00032019139700000311
And
Figure BDA00032019139700000312
cross direction of and
Figure BDA00032019139700000313
whether the included angle is an acute angle or not is judged, whether the central point of the unit is in the Thiessen polygon or not is judged, and matching parameters between the unit A and the Thiessen polygon are defined:
Figure BDA00032019139700000314
if m is for all positive integers m ∈ [1, m1]All are provided with
Figure BDA00032019139700000315
The cell a belongs to the traversed Thiessen polygon, and the process goes to step S47, otherwise, the process goes back to step S42;
s47, assuming crossing points of warp and weft yarns
Figure BDA00032019139700000316
In the warp and weft directions of
Figure BDA00032019139700000317
And
Figure BDA00032019139700000318
local coordinate system X of unitR-YR-ZRThe coordinate axis in (1) is represented by unit vectors in the world coordinate system X-Y-Z as:
Figure BDA00032019139700000319
the main directions of the materials in the unit are all represented by the three vectors;
and S48, repeating the steps S41 to S47 until the main material directions of all the units in the finite element model of the whole area are obtained.
Further, the complex geometric region includes a torsional geometric feature region.
Further, in step S6, the warp yarn cross section and the weft yarn cross section are respectively set as
Figure BDA00032019139700000320
And
Figure BDA00032019139700000321
respectively representing the x-th section of the ith warp yarn and the y-th section of the jth weft yarn, and the process of judging the material property of each unit in the local refined finite element model comprises the following steps:
s61, assuming that the number of nodes of a certain unit B in the local refined finite element model is n2Each node of the cell is denoted as Eg，g＝1,2,…n2The center node coordinates of cell B are calculated according to the following formula:
Figure BDA0003201913970000041
s62, traversing all the warp yarn sections;
s63, assuming the traversed warp yarn section
Figure BDA0003201913970000042
The number of contour control points of (1) is m2Where x is 1,2, …, c-1, c is the number of sections of the ith warp yarn, and the contour point is
Figure BDA0003201913970000043
h＝1,2,…m2And calculating the coordinates of the center point of the section by adopting the following formula:
Figure BDA0003201913970000044
in the formula, point FcenterAt the same time, is a sectional point of the yarn central line, and the direction vector of the sectional point is
Figure BDA0003201913970000045
The direction vector and the cross section
Figure BDA0003201913970000046
Vertically;
s64, judging whether the unit B is positioned on the cross section of the warp yarn
Figure BDA0003201913970000047
Cross section of warp
Figure BDA0003201913970000048
And defining and calculating a parameter rho between the planes by adopting the following formula:
Figure BDA0003201913970000049
if rho is less than or equal to 0, the unit B is positioned on the cross section of the warp yarn
Figure BDA00032019139700000410
Cross section of warp
Figure BDA00032019139700000411
Between the planes, go to step S65, otherwise go through the next warp yarn section and return to step S63;
s65, calculating E by the following formulacenterAt the cross section of warp
Figure BDA00032019139700000412
Projected point on:
Figure BDA00032019139700000413
s66, calculating the cross section of the warp yarn by adopting the following formula
Figure BDA00032019139700000414
Direction vector of contour edge:
Figure BDA00032019139700000415
s67, calculating the point by the following formula
Figure BDA00032019139700000416
To warp cross section
Figure BDA00032019139700000417
Direction vector of contour control point of (1):
Figure BDA00032019139700000418
s68, judging whether the unit B is positioned on the cross section of the warp yarn
Figure BDA00032019139700000419
Cross section of warp
Figure BDA00032019139700000420
In the yarn section between, fixDefining parameters:
Figure BDA00032019139700000421
if h is an arbitrary positive integer e [1, m2]All are provided with
Figure BDA00032019139700000422
Then it means that unit B is in warp cross-section
Figure BDA00032019139700000423
Cross section of warp
Figure BDA00032019139700000424
In the yarn segment between the two, the unit B is judged to belong to the warp, the warp direction of the segment is taken as the material main direction of the unit B, and the step S610 is carried out; otherwise, traversing the next warp yarn section and returning to the step S63, and if all the warp yarn interfaces are traversed, entering the step S69;
s69, according to the method from S61 to S68, judging whether the unit B is in a weft yarn subsection, if the unit B is not in any weft yarn subsection, judging that the unit B belongs to the matrix, and if the unit B is in one of the weft yarn subsections, judging that the unit B belongs to the weft yarn;
s610, repeating the steps S61 to S69 until the material properties of all the units in the local subdivision finite element model are judged.
The invention has the beneficial effects that:
the invention can use the composite material prefabricated body geometric model and the finite element mechanical model of the component as input to realize the rapid judgment of the material attribute of each unit in the finite element model. The invention provides a practical pretreatment tool for predicting the mechanical property of the woven composite material. The method can be applied to various mechanical property prediction methods such as a homogenization method, a submodel method, a multi-scale method and the like.
Drawings
FIG. 1 is a schematic view of a characteristic ply model of an embodiment of the invention.
FIG. 2 is a schematic view of a characteristic ply yarn path model of an embodiment of the invention.
FIG. 3 is a schematic view of a detailed model of a featured ply preform of an embodiment of the invention.
FIG. 4 is a schematic view of a finite element model of a feature layup ensemble according to an embodiment of the present invention.
FIG. 5 is a schematic diagram of a finite element model for local refining of a feature mat according to an embodiment of the present invention.
FIG. 6 is a schematic diagram illustrating a solution of a layered Delaunay triangulation network according to an embodiment of the present invention.
FIG. 7 is a schematic diagram of a solution for a layered Thiessen polygon mesh according to an embodiment of the present invention.
FIG. 8 is a schematic view of the integrated process of the woven composite material geometric-mechanical model according to an embodiment of the present invention.
FIG. 9 is a flow chart of a geometric-mechanical-model integration method for braiding composite materials, in accordance with an embodiment of the present invention.
Detailed Description
The present invention will now be described in further detail with reference to the accompanying drawings.
It should be noted that the terms "upper", "lower", "left", "right", "front", "back", etc. used in the present invention are for clarity of description only, and are not intended to limit the scope of the present invention, and the relative relationship between the terms and the terms is not limited by the technical contents of the essential changes.
According to the embodiment of the invention, the composite material prefabricated body geometric model and the finite element mechanical model of the component are used as input, and the rapid judgment of the material properties of each unit in the finite element model is realized. As shown in FIG. 1, the detailed steps of this embodiment will be described by taking the example of a feature overlay consisting of arched geometric features and twisted geometric features. Wherein the geometric model of the composite preform comprises a ply yarn path model as shown in FIG. 2 and a preform mesoscale model as shown in FIG. 3; the finite element mechanical model includes a global region finite element model as shown in FIG. 4 and a local non-uniform region finite element model as shown in FIG. 5.
Referring to fig. 9, the material property simulation method for determining the material property of each element in the finite element model through the coordinate relationship between the geometric model of the composite material preform and the finite element mechanical model comprises the following steps:
s1, inputting a composite material prefabricated body geometric model and a finite element mechanical model of the component; the geometric model of the composite material preform comprises a layering yarn path model and a preform mesoscale model; the finite element mechanical model comprises a finite element model of the whole area and a finite element model of the local non-uniform area.
S2, on the basis of the ply yarn path model, setting every 4 adjacent interweaving points to form a three-dimensional rhombic unit cell, dividing each unit cell into two Delaunay triangles, and constructing to obtain a Delaunay triangular grid.
S3, establishing a Thiessen polygonal mesh based on the Delaunay triangular mesh.
And S4, obtaining the main directions of the units according to the Thiessen polygonal meshes and the finite element model of the whole area.
And S5, selecting a complex geometric area of the ply yarn path model for local refining, and generating a local refining finite element model.
And S6, judging the material properties of each unit in the local refined finite element model, wherein the material properties comprise warp yarns, weft yarns and a matrix.
The specific implementation steps are illustrated by taking the implementation flow of fig. 8 as an example (the method can be applied to other models): 1) the yarn path model shown in FIG. 2 and the characteristic ply integral finite element model shown in FIG. 4 are input. 2) And (3) determining the main material direction of each unit in the finite element model of the whole paving layer according to the coordinate relation between the two models input in the step 1), namely the processes of (a) to (c) of the drawing 8. 3) The torsion area is selected as a local complex structure in FIG. 1, and a feature ply preform microscopic model as shown in FIG. 3 and a local refined finite element model as shown in FIG. 5 are input. 4) Determining the material property of each unit in the local refined finite element model according to the coordinate relation between the two models input in the step 3), and judging that each unit belongs to the warp or weft or the matrix, namely the processes of (b) - (d) in the step 8). The method comprises the following specific steps:
the method comprises the following steps: the Delaunay triangular mesh is constructed based on the yarn path model shown in fig. 2.
Step 1.1: as shown in FIG. 6, the warp and weft base lines divide the yarn path model into intersections
Figure BDA0003201913970000061
4 quadrants as the center, every 4 adjacent interweaving points can form a three-dimensional rhombic unit cell, and 4 vertexes of a rhombic unit cell are set as
Figure BDA0003201913970000062
And
Figure BDA0003201913970000063
step 1.2: computing edges
Figure BDA0003201913970000064
And edge
Figure BDA0003201913970000065
Angle therebetween
Figure BDA0003201913970000066
Step 1.3: judging angle
Figure BDA0003201913970000067
Whether it is an obtuse angle: if the angle is obtuse, the diagonal line is used
Figure BDA0003201913970000068
Divide the unit cell into two Delaunay triangles
Figure BDA0003201913970000069
And
Figure BDA00032019139700000610
if not obtuse, by diagonal
Figure BDA00032019139700000611
Will unit cellDivided into two Delaunay triangles
Figure BDA00032019139700000612
Step 1.3: the entire yarn path model is divided into Delaunay triangular meshes according to step 1.2.
Step two: and (3) establishing a Thiessen polygonal mesh based on the Delaunay triangular mesh obtained in the first step, wherein a schematic diagram is shown in FIG. 7.
Step 2.1: a warp-weft intersection point
Figure BDA00032019139700000613
Arranged according to a counterclockwise line, with
Figure BDA00032019139700000614
Figure BDA00032019139700000615
Is represented by (a) wherein kTAs the crossing point of warp and weft yarns
Figure BDA00032019139700000616
The neighborhood triangle here is the Delaunay triangle that contains the point.
Step 2.2: calculating the center of a circumscribed circle of each adjacent triangle of the intersection point of the warp and weft yarns and recording the center of the circumscribed circle
Figure BDA00032019139700000617
Step 2.3: are connected in a counterclockwise sequence
Figure BDA00032019139700000618
Obtaining the intersection point of the warp and weft
Figure BDA00032019139700000619
Corresponding Thiessen polygons.
Step 2.4: and (4) obtaining Thiessen polygons corresponding to the intersection points of the other warps and wefts according to the method of the step 2.1 to the step 2.3.
Step three: and obtaining the main direction of the unit according to the Thiessen polygonal grid established in the step two and the integral finite element model.
Step 3.1: suppose that the node number of a certain element A in the whole finite element model shown in FIG. 4 is n1The unit node is denoted as Pl(l＝1,2,…,n1) Calculating the coordinates of the center point of the unit:
Figure BDA0003201913970000071
step 3.2: polygons in the Thiessen polygon mesh are traversed.
Step 3.3: let the number of traversed Thiessen polygon vertices be m1The vertex is
Figure BDA0003201913970000072
The corresponding warp and weft yarn intersection point is
Figure BDA0003201913970000073
The ply is at the point
Figure BDA0003201913970000074
Normal vector at point is
Figure BDA0003201913970000075
Calculating the direction vector of each side of the Thiessen polygon:
Figure BDA0003201913970000076
step 3.4: calculating the center point P of the unit AcenterAt the passing point
Figure BDA0003201913970000077
And by vectors
Figure BDA0003201913970000078
Projected points on a plane as a normal vector:
Figure BDA0003201913970000079
step 3.5: calculating points
Figure BDA00032019139700000720
Direction vectors to vertices of Thiessen polygons:
Figure BDA00032019139700000710
step 3.6: by making a judgment
Figure BDA00032019139700000711
And
Figure BDA00032019139700000712
cross direction of and
Figure BDA00032019139700000713
whether the included angle is an acute angle or not is judged, whether the central point of the unit is in the Thiessen polygon or not is judged, and matching parameters between the unit A and the Thiessen polygon are defined:
Figure BDA00032019139700000714
if m is for all positive integers m ∈ [1, m1]All are provided with
Figure BDA00032019139700000715
Then cell a belongs to the traversed Thiessen polygon, otherwise go back to step 3.2.
Step 3.7: assuming crossing of warp and weft yarns
Figure BDA00032019139700000716
In the warp and weft directions of
Figure BDA00032019139700000717
And
Figure BDA00032019139700000718
the unit local coordinate system XR-YR-ZRThe coordinate axis in (1) is represented by a unit vector in the world coordinate system X-Y-Z:
Figure BDA00032019139700000719
the main direction of the material within the cell can be represented by the three vectors described above.
Step 3.8: the material main direction of all units is determined according to the method from step 3.1 to step 3.7.
Step four: selecting a complex geometric area of the layer model for local thinning, selecting a torsion area as a local thinning area for the characteristic layer shown in figure 1, and obtaining a thinned finite element grid model as shown in figure 5.
Step five: judging the material property of each unit in the local refined finite element model, namely whether a certain unit belongs to warp yarns, weft yarns or a matrix, wherein the cross section of each warp yarn and the cross section of each weft yarn are respectively assumed to be
Figure BDA0003201913970000081
And
Figure BDA0003201913970000082
respectively showing the x-th cross section of the ith warp yarn and the y-th cross section of the jth weft yarn.
Step 5.1: suppose the number of nodes of a certain unit B in the local refined finite element model is n2Each node of the cell is denoted as Eg(g＝1,2,…n2) Calculating the center node coordinates of the unit B:
Figure BDA0003201913970000083
step 5.2: all warp yarn sections are traversed.
Step 5.3: assuming traversed cross-section
Figure BDA0003201913970000084
The number of contour control points of (1) is m2C is the number of sections of the ith warp yarn and the contour point is
Figure BDA0003201913970000085
Calculating the coordinates of the center point of the cross section:
Figure BDA0003201913970000086
point FcenterAt the same time, is a sectional point of the yarn central line, and the direction vector of the sectional point is
Figure BDA0003201913970000087
The vector and the cross section
Figure BDA0003201913970000088
And is vertical.
Step 5.4: determine whether the unit B is located on the cross section
Figure BDA0003201913970000089
And cross section
Figure BDA00032019139700000810
Between the planes, the parameters are defined and calculated:
Figure BDA00032019139700000811
if rho is less than or equal to 0, the cell B is positioned on the section
Figure BDA00032019139700000812
And cross section
Figure BDA00032019139700000813
Between the planes, otherwise, traverse the next section and return to step 5.3.
Step 5.5: calculation of EcenterIn cross section
Figure BDA00032019139700000814
Projected point on:
Figure BDA00032019139700000815
step 5.6: calculating cross sections
Figure BDA00032019139700000816
Direction vector of contour edge:
Figure BDA00032019139700000817
step 5.7: calculating points
Figure BDA00032019139700000818
To cross section
Figure BDA00032019139700000819
Direction vector of contour control point of (1):
Figure BDA00032019139700000820
step 5.8: determine whether the unit B is located on the cross section
Figure BDA00032019139700000821
And cross section
Figure BDA00032019139700000822
In the yarn section in between, the parameters are defined:
Figure BDA00032019139700000823
if h is an arbitrary positive integer e [1, m2]All are provided with
Figure BDA00032019139700000824
It means that the cell B is in cross section
Figure BDA00032019139700000825
And cross section
Figure BDA00032019139700000826
In the yarn section between, and the warp direction of this section is taken as the material main direction of the unit B. Otherwise, traversing the next section and returning to the step 5.3, and if all the warp yarn sections are not traversed, entering the next step.
Step 5.9: and judging whether the unit B is positioned in a certain weft yarn section according to the methods from the step 5.1 to the step 5.8, and if the unit B is not positioned in any weft yarn section, the unit B belongs to the base body.
Step 5.10: and judging the material properties of all the units in the local subdivision finite element model according to the methods of the step 5.1 to the step 5.9.
The above is only a preferred embodiment of the present invention, and the protection scope of the present invention is not limited to the above-mentioned embodiments, and all technical solutions belonging to the idea of the present invention belong to the protection scope of the present invention. It should be noted that modifications and embellishments within the scope of the invention may be made by those skilled in the art without departing from the principle of the invention.

Claims (6)

1. A geometric-mechanical model integration method for weaving composite materials is characterized in that the material property simulation method comprises the following steps:
s1, inputting a composite material prefabricated body geometric model and a finite element mechanical model of the component; the geometric model of the composite material preform comprises a layering yarn path model and a preform mesoscale model; the finite element mechanical model comprises a finite element model of the whole area and a finite element model of the local non-uniform area;
s2, on the basis of the ply yarn path model, setting every 4 adjacent interweaving points to form a three-dimensional rhombic unit cell, dividing each unit cell into two Delaunay triangles, and constructing to obtain a Delaunay triangular grid;
s3, establishing a Thiessen polygonal mesh based on the Delaunay triangular mesh;
s4, obtaining a main direction of the unit according to a finite element model of the Thiessen polygonal grid and the whole area;
s5, selecting a complex geometric area of the ply yarn path model for local refining to generate a local refining finite element model;
and S6, judging the material properties of each unit in the local refined finite element model, wherein the material properties comprise warp yarns, weft yarns and a matrix.
2. The integrated geometric-mechanical model method for weaving composite material according to claim 1, wherein in step S2, the process of constructing Delaunay triangular mesh based on the ply yarn path model comprises the following steps:
s21, dividing the yarn path model into crossing points by the warp base lines and the weft base lines
Figure FDA0003201913960000011
4 quadrants as the center, each 4 adjacent interweaving points form a three-dimensional rhombic unit cell, and 4 vertexes of the three-dimensional rhombic unit cell are set as
Figure FDA0003201913960000012
And
Figure FDA0003201913960000013
i is the coordinate value of the warp yarn constituting the three-dimensional rhombohedral unit cell, and j is the coordinate value of the weft yarn constituting the three-dimensional rhombohedral unit cell;
s22, selecting a three-dimensional rhomboid unit cell, and calculating the corresponding side
Figure FDA0003201913960000014
And edge
Figure FDA0003201913960000015
Angle therebetween
Figure FDA0003201913960000016
S23, judging the angle
Figure FDA0003201913960000017
Whether it is an obtuse angle: if the angle is obtuse, the diagonal line is used
Figure FDA0003201913960000018
Divide the unit cell into two Delaunay triangles
Figure FDA0003201913960000019
And
Figure FDA00032019139600000110
if not obtuse, by diagonal
Figure FDA00032019139600000111
Divide the unit cell into two Delaunay triangles
Figure FDA00032019139600000112
S24, repeating steps S22 to S23, dividing the entire yarn path model into Delaunay triangular meshes.
3. The method for integrating geometric-mechanical models of woven composite materials according to claim 1, wherein in step S3, the process of building Thiessen polygonal mesh based on Delaunay triangular mesh comprises the following steps:
s31, selecting any one warp and weft yarn intersection point
Figure FDA00032019139600000113
Arranged counterclockwise, with
Figure FDA00032019139600000114
Denotes that k is 0,1,2, …, kTWherein k isTAs the crossing point of warp and weft yarns
Figure FDA00032019139600000115
The number of the neighborhood triangles is Delaunay triangle containing the crossing point of the warp and weft; i is the coordinate value of the warp yarn, j is the coordinate value of the weft yarn;
s32, calculating the crossing point of the warp and weft yarns
Figure FDA00032019139600000116
The center of the circumscribed circle of each neighborhood triangle is recorded as
Figure FDA00032019139600000117
S33, connected in a counterclockwise order
Figure FDA0003201913960000021
Obtaining the intersection point of the warp and weft
Figure FDA0003201913960000022
A corresponding Thiessen polygon;
and S34, repeating the steps S31 to S33 until obtaining Thiessen polygons corresponding to all warp and weft yarn intersections.
4. The method for integrating geometric-mechanical models of woven composite materials according to claim 1, wherein in step S4, the step of obtaining main directions of elements according to the finite element model of Thiessen polygonal meshes and integral regions comprises the following steps:
s41, assuming that the node number of any unit A in the finite element model of the whole area is n1The unit node is denoted as Pl，l＝1,2,…,n1And calculating the coordinates of the center point of the unit A according to the following formula:
Figure FDA0003201913960000023
s42, traversing the polygons in the Thiessen polygon mesh;
s43, the number of the traversed Thiessen polygon vertexes is m1The vertex is
Figure FDA0003201913960000024
The corresponding warp and weft yarn intersection point is
Figure FDA0003201913960000025
The ply is at the point
Figure FDA0003201913960000026
Normal vector at point is
Figure FDA0003201913960000027
i is the coordinate value of the warp yarn, j is the coordinate value of the weft yarn; calculating the direction vector of each side of the Thiessen polygon by adopting the following formula:
Figure FDA0003201913960000028
s44, calculating the center point P of the cell A by the following formulacenterAt the passing point
Figure FDA0003201913960000029
And by vectors
Figure FDA00032019139600000210
Projected points on a plane as a normal vector:
Figure FDA00032019139600000211
s45, calculating to obtain point P 'by the following formula'centerDirection vectors to vertices of Thiessen polygons:
Figure FDA00032019139600000212
s46, by judging
Figure FDA00032019139600000213
And
Figure FDA00032019139600000214
cross direction of and
Figure FDA00032019139600000215
whether the included angle is an acute angle or not is judged, whether the central point of the unit is in the Thiessen polygon or not is judged, and matching parameters between the unit A and the Thiessen polygon are defined:
Figure FDA00032019139600000216
if m is for all positive integers m ∈ [1, m1]All are provided with
Figure FDA00032019139600000217
The cell a belongs to the traversed Thiessen polygon, and the process goes to step S47, otherwise, the process goes back to step S42;
s47, assuming crossing points of warp and weft yarns
Figure FDA00032019139600000218
In the warp and weft directions of
Figure FDA00032019139600000219
And
Figure FDA00032019139600000220
local coordinate system X of unitR-YR-ZRThe coordinate axis in (1) is represented by unit vectors in the world coordinate system X-Y-Z as:
Figure FDA0003201913960000031
the main directions of the materials in the unit are all represented by the three vectors;
and S48, repeating the steps S41 to S47 until the main material directions of all the units in the finite element model of the whole area are obtained.
5. The geometric-mechanical-model integration method for braided composite material of claim 1, wherein said complex geometric regions comprise torsional geometric feature regions.
6. The method of claim 1, wherein in step S6, the warp yarn cross-section and the weft yarn cross-section are defined as
Figure FDA0003201913960000032
And
Figure FDA0003201913960000033
respectively representing the x-th section of the ith warp yarn and the y-th section of the jth weft yarn, and the process of judging the material property of each unit in the local refined finite element model comprises the following steps:
s61, assuming that the number of nodes of a certain unit B in the local refined finite element model is n2Each node of the cell is denoted as Eg，g＝1,2,…n2The center node coordinates of cell B are calculated according to the following formula:
Figure FDA0003201913960000034
s62, traversing all the warp yarn sections;
s63, assuming the traversed warp yarn section
Figure FDA0003201913960000035
The number of contour control points of (1) is m2Where x is 1,2, …, c-1, c is the number of sections of the ith warp yarn, and the contour point is
Figure FDA0003201913960000036
Calculating the coordinates of the center point of the cross section by adopting the following formula:
Figure FDA0003201913960000037
in the formula, point FcenterAt the same time, is a sectional point of the yarn central line, and the direction vector of the sectional point is
Figure FDA0003201913960000038
The direction vector and the cross section
Figure FDA0003201913960000039
Vertically;
s64, judging whether the unit B is positioned on the cross section of the warp yarn
Figure FDA00032019139600000310
Cross section of warp
Figure FDA00032019139600000311
And defining and calculating a parameter rho between the planes by adopting the following formula:
Figure FDA00032019139600000312
if rho is less than or equal to 0, the unit B is positioned on the cross section of the warp yarn
Figure FDA00032019139600000313
Cross section of warp
Figure FDA00032019139600000314
Between the planes, go to step S65, otherwise go through the next warp yarn section and return to step S63;
s65, calculating E by the following formulacenterAt the cross section of warp
Figure FDA00032019139600000315
Projected point on:
Figure FDA00032019139600000316
s66, calculating the cross section of the warp yarn by adopting the following formula
Figure FDA00032019139600000317
Direction vector of contour edge:
Figure FDA00032019139600000318
s67, calculating the point by the following formula
Figure FDA0003201913960000041
To warp cross section
Figure FDA0003201913960000042
Direction vector of contour control point of (1):
Figure FDA0003201913960000043
s68, judging whether the unit B is positioned on the cross section of the warp yarn
Figure FDA0003201913960000044
Cross section of warp
Figure FDA0003201913960000045
In the yarn section in between, the parameters are defined:
Figure FDA0003201913960000046
if h is an arbitrary positive integer e [1, m2]All are provided with
Figure FDA0003201913960000047
Then it means that unit B is in warp cross-section
Figure FDA0003201913960000048
Cross section of warp
Figure FDA0003201913960000049
In the yarn segment between the two, the unit B is judged to belong to the warp, the warp direction of the segment is taken as the material main direction of the unit B, and the step S610 is carried out; otherwise, traversing the next warp yarn section and returning to the step S63, and if all the warp yarn interfaces are traversed, entering the step S69;
s69, according to the method from S61 to S68, judging whether the unit B is in a weft yarn subsection, if the unit B is not in any weft yarn subsection, judging that the unit B belongs to the matrix, and if the unit B is in one of the weft yarn subsections, judging that the unit B belongs to the weft yarn;
s610, repeating the steps S61 to S69 until the material properties of all the units in the local subdivision finite element model are judged.
CN202110906645.1A 2021-08-09 2021-08-09 Geometric mechanical model integration method for braiding composite material Active CN113609735B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110906645.1A CN113609735B (en) 2021-08-09 2021-08-09 Geometric mechanical model integration method for braiding composite material

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110906645.1A CN113609735B (en) 2021-08-09 2021-08-09 Geometric mechanical model integration method for braiding composite material

Publications (2)

Publication Number Publication Date
CN113609735A true CN113609735A (en) 2021-11-05
CN113609735B CN113609735B (en) 2024-02-20

Family

ID=78339926

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110906645.1A Active CN113609735B (en) 2021-08-09 2021-08-09 Geometric mechanical model integration method for braiding composite material

Country Status (1)

Country Link
CN (1) CN113609735B (en)

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110536307A (en) * 2019-09-03 2019-12-03 南京华苏科技有限公司 A kind of method and device of the adjacent section planning based on Thiessen polygon

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110536307A (en) * 2019-09-03 2019-12-03 南京华苏科技有限公司 A kind of method and device of the adjacent section planning based on Thiessen polygon

Also Published As

Publication number Publication date
CN113609735B (en) 2024-02-20

Similar Documents

Publication Publication Date Title
Winslow et al. Multi-objective optimization of free-form grid structures
Glaessgen et al. Geometrical and finite element modelling of textile composites
Ullah et al. A unified framework for the multi-scale computational homogenisation of 3D-textile composites
CN112182929A (en) Size control-considered cross-scale reliability topological optimization method for porous material
Changizi et al. Stress-based topology optimization of steel-frame structures using members with standard cross sections: Gradient-based approach
CN112836411A (en) Method and device for optimizing structure of stiffened plate shell, computer equipment and storage medium
Tellier et al. Linear Weingarten surfaces for conceptual design
CN111310366B (en) Method for establishing micro model of plain weave composite material preform
Gowayed et al. Modification and application of a unit cell continuum model to predict the elastic properties of textile composites
CN113609735A (en) Geometric-mechanical model integration method for weaving composite material
Wu et al. Design and optimization of the variable-density lattice structure based on load paths
Geleta et al. Prediction of effective material properties for triaxially braided textile composite
Zhao et al. Triaxial woven fabric (TWF) composites with open holes (part I): finite element models for analysis
Haasemann et al. Multi-scale modelling and simulation of textile reinforced materials
Pastore et al. Applications of a meso-volume-based analysis for textile composite structures
CN113987882B (en) Digital modeling method for microstructure yarn structure of woven composite material
CN114692468A (en) Method for predicting cross-scale damage of continuous fiber reinforced ceramic matrix composite
Winslow et al. Mapping two-way grids onto free-form surfaces
CN116187074A (en) Multi-scale topological optimization method of anisotropic periodic structure material based on isogeometry
CN114912306A (en) Modeling method of 2.5D woven composite material flat plate structure multi-scale model
Georgiou Interactive Structural Analysis and Form-Finding
CN114254541B (en) Arch dam full dam face equivalent stress calculation method and system
CN112818576B (en) Multi-level optimization method for curve fiber composite structure design
Dzwierzynska Shaping of Curvilinear Steel Bar Structures for Variable Environmental Conditions Using Genetic Algorithms—Moving towards Sustainability. Materials 2021, 14, 1167
Jialiang et al. Research on the surface flattening model for carbon fiber plain-woven composite preforming processing

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant