A modeling method for unidirectional fiber periodic reinforced composite materials including fiber waviness

Abstract

The present invention proposes a modeling method for unidirectional fiber periodic reinforced composite materials with fiber corrugation, which includes the creation of a soft-core fiber system and the correction of a hard-core fiber system; the soft-core fiber system does not consider the overlap between the fibers, allows overlap between the fibers, and sets the geometric periodicity of the entire fiber system; the hard-core fiber system is further adjusted on the basis of the soft-core system to eliminate the overlap between the fibers without destroying the fiber structure. The adjustment process introduces a force bias algorithm and applies it to the created soft-core fiber system. Through the present invention, a periodic three-dimensional representative volume element model of a unidirectional fiber reinforced composite material with the desired fiber corrugation can be generated, which helps to understand the relationship between the component properties / microstructure and macroscopic properties of the material at a deeper level.

Description

Technical Field

[0001] The invention relates to a high-fidelity modeling method of a unidirectional fiber-reinforced composite material containing fiber corrugation and periodic conditions, and belongs to the field of composite material manufacturing. Background Art

[0002] Fiber reinforced composites (FRP) are widely used in high-end equipment manufacturing fields such as aerospace due to their light weight and high strength. Among them, unidirectional fiber reinforced composites (UD-FRP) have fibers arranged in the same direction and are the most basic and typical form of FRP. Other types of FRP can be regarded as a superposition of various UD-FRPs.

[0003] Since composite materials themselves have the characteristics of multiphase structure, the factors affecting the mechanical properties of composite materials include component properties (matrix, fiber and interface, etc.) and microstructural characteristics (shape, size and distribution characteristics, etc.). In the manufacturing process of composite materials, due to the influence of the material properties and thermal residual stress generated during the curing process, the fiber configuration in the actual composite material usually inevitably has a certain degree of corrugation. The fiber corrugation will cause large fluctuations in the mechanical properties of the composite material such as axial compression and transverse shear. In order to accurately characterize the mechanical behavior of UD-FRP, it is necessary to establish a high-fidelity representative volume element (RVE) model of the microstructure to analyze its homogenization characteristics, so as to understand the relationship between the component properties / microstructure and macroscopic properties of the material at a deeper level.

[0004] In order to describe the influence of fiber waviness characteristics on the axial damage characteristics of UD-FRP, Catalanotti, Sebaey and others developed a semi-random algorithm to characterize fiber waviness in a statistical way. The algorithm first discretizes continuous fibers into multiple groups of spheres, determines the deflection angle of each group of fiber spheres inside UD-FRP by electronic computed tomography (CT), statistically fits the CT data according to the von-Mises distribution, and characterizes the size of fiber waviness with concentration parameters. Finally, a random perturbation program and an optimization program are used to ensure that the RVE has the same concentration parameters as the actual specimen, while ensuring the geometric periodicity of the model itself. Although the algorithm constructs an RVE model with the same fiber waviness as the actual UD-FRP in a statistical sense, the random perturbation program only acts on each fiber sphere individually, without considering the interaction between adjacent spheres on the fiber, resulting in excessive distortion of the fiber locally, which does not exist in the actual fiber structure. Summary of the invention

[0005] In order to accurately characterize the fiber waviness and geometric periodicity characteristics of UD-FRP, this paper proposes a UD-FRP three-dimensional periodic corrugated fiber RVE modeling method based on fiber random walk algorithm, which can be generally divided into two parts: fiber soft core system creation and hard core system correction, such as Figure 1 shown.

[0006] The construction process of the soft-core system first creates two-dimensional random fiber position points of the fiber based on the requirements of the fiber radius and fiber volume fraction in the RVE (Representative Volume Element). The creation method has been reported in many literatures and will not be repeated here. Then, along the main direction of the fiber (set as the global Z coordinate direction in this method), the random walk algorithm is used to create three-dimensional fiber configurations at the existing two-dimensional fiber position points in sequence, and finally the entire fiber system is set to geometric periodicity. This process allows fibers to overlap with each other, so it is called a soft-core fiber system.

[0007] The hard-core system adjustment process is to further adjust the soft-core fiber system that has been constructed, and eliminate the overlap between fibers while keeping the fiber configuration basically unchanged. The adjustment process introduces a force bias algorithm, which defines the relationship between the repulsion / attraction between fibers and the relative distance and the bending force-curvature of the fibers themselves, and applies it to the soft-core fiber system that has been created above. In addition, in order to speed up the adjustment efficiency, a three-dimensional neighbor list algorithm is introduced and applied to each fiber.

[0008] In order to accurately characterize the fiber waviness and geometric periodicity characteristics of UD-FRP, the present invention first creates each fiber individually based on the fiber random walk algorithm. Each fiber is characterized as a Bézier curve, and the direction of each control point of the curve obeys the Von-Mises Fischer distribution. Given the fiber waviness parameters and RVE size, the Bézier curves of individual fibers are generated in turn. In this algorithm, μ represents the axial vector of UD-FRP, u represents the direction vector of each control point of the Bézier curve, and k is the corresponding reliability parameter. To ensure that the average axial direction of each fiber follows a uniform global orientation, a rotation matrix is ​​created.

[0009] D(ω,e,α)=(e·ω)e+cosα((e×ω)×e)+sinα(e×ω)

[0010] Rotate the fiber curve globally so that the average fiber orientation is aligned with the global principal direction μ of the RVE 1 =(0,0,1) to keep it consistent. Assign the coordinate values ​​of the fiber Bézier curve to a set of equally spaced r tspheres, and the sphere spacing along the z direction is Divide the fibers into r t Plane P i (i=1,,,r t ), set the periodic conditions for each ball chain point on each plane separately in To avoid overlap between fibers, a three-dimensional neighbor list is established for each sphere in the fiber system to speed up the overlap judgment. Repulsion is applied to overlapping spheres between different fibers. (K r is the repulsive force constant of the ball chain, Δl is the overlap of the two balls), and in order to ensure that the fiber structure is not destroyed in the process, a spring force, i.e., gravitational force, is applied between the balls in the same fiber. and bending force (K T and K B are the tension and compression spring constants and bending spring constants of the ball chain, both given by actual material parameters, The final fiber digital model is imported into the modeling software, and a UD-FRP three-dimensional RVE solid model with fiber waviness and periodicity characteristics can be created through simple curve fitting and sweeping operations.

[0011] The method specifically includes the following process:

[0012] Step 1: According to the size of the representative volume element RVE and the requirements for the fiber radius and fiber volume fraction, a random perturbation method is used to create random two-dimensional fiber position points; a random walk algorithm is used at each position point to create each fiber in the system, and each fiber is characterized by a Bézier curve;

[0013] Step 2: The average intercept of each fiber Bézier curve along the global principal direction (i.e., the RVE axial direction) is r t coordinates, respectively assigned to r t sphere (r t The value of is selected as appropriate. The larger the value, the better the smoothness of the fiber generated, but the greater the amount of calculation). Therefore, the ball chain is used to characterize the curve of each fiber, in order to facilitate the subsequent adjustment of the local coordinates of the fiber.

[0014] For the fibers generated in the previous step, the orientation of the fibers Not necessarily equal to the global principal direction μ of RVE 1 =(0,0,1). Therefore, in order to ensure that all created fibers follow the global orientation distribution, a rotation matrix is ​​created to rotate each fiber while keeping the fiber structure unchanged. and the global main direction μ of RVE1 =(0,0,1) to keep consistent, such as Figure 3 As shown. During this process, the structure of the fiber (length, radius, curvature) remains unchanged;

[0015] Step 3: If Figure 4 As shown, each ball chain is connected along the global main direction with a ball chain spacing L 3 / r t -1 is divided into r t Plane P i ; At this point, it can be found that because the ball chain and the plane are divided at equal intervals along the Z direction at the same distance, there are exactly n fiber spheres on each plane. At this time, the periodic conditions are set for each ball chain point on each plane separately;

[0016] Step 4: After the periodic conditions are set, because the fibers are previously generated separately, there may be overlaps between the fibers. To avoid this, a three-dimensional neighbor list is established for all fiber spheres in the fiber system. The neighbor list of each sphere center contains all other sphere centers that may overlap with it (the distance is less than the sum of the radii). In the neighbor list, the neighbor sphere centers of all sphere centers on the target fiber sphere chain are searched, so that the overlapping fiber sphere chains and the overlapping spheres can be quickly found, avoiding the need to traverse all fiber spheres for each judgment, which greatly improves the efficiency of the algorithm;

[0017] Step 5: For each fiber sphere, search the neighbor list for the fiber sphere chains that overlap with it and the fiber spheres that generate the overlap;

[0018] Apply repulsive force to the overlapping spheres between different fibers found

[0019]

[0020] K r is the repulsive force constant of the ball chain, which is calculated from the material parameters of the fiber, Δl is the overlap between the two balls, q 1 ,q 2 are the coordinates of the centers of the two balls;

[0021] At the same time, in order to ensure that the fiber structure is not destroyed during the process, a spring force is applied between the spheres in the same fiber, including a tension and compression spring restoring force.

[0022]

[0023] and the bending spring restoring force Where K T and K B are the tension and compression spring constants and bending spring constants of the ball chain, K B =Ef πr 3 , where E f is the elastic modulus of the fiber, l 0 is the initial distance between the balls, which is usually taken as is the friction factor, q i is the coordinate of the center of the i-th sphere on the fiber, and r is the radius of the sphere. This process can eliminate the overlap between fibers.

[0024] Step 6: Import the final fiber digital model into the modeling software, and create a UD-FRP three-dimensional RVE solid model with fiber waviness and periodicity characteristics through curve fitting and sweeping operations.

[0025] Furthermore, in step 1, according to the given fiber waviness parameter and RVE size, the Bézier curve of a single fiber is

[0026]

[0027] Where s is the number of control points of the curve, which is defined by the user, 0≤t≤1, x i (i=0,1,2,...,s) represents the coordinates of each control point.

[0028] Furthermore, in step 1, the tangent direction at each control point of the Bézier curve obeys the given orientation distribution function von Mises-Fischer distribution Where μ represents the axial vector of the unidirectional fiber reinforced composite material, u represents the direction vector at the tangent of each control point of the Bézier curve, and k is the reliability parameter corresponding to u relative to μ, which determines the correlation between u and μ.

[0029] Furthermore, in step 1, the tangent direction at each control point of the Bézier curve follows the multivariate von Mises-Fisher distribution

[0030]

[0031] in μ 1 =(0,0,1) represents the axial direction of unidirectional fiber reinforced composite materials, μ 2 is the direction vector representing the tangent direction of each control point, and its modulus is 1; k 1 With k 2 The tangent directions of each control point of the fiber are With the global principal direction μ 1 =(0,0,1) and the tangent direction of the previous control point on the same fiber The corresponding correlation parameter, k 1 The value determines the degree of deviation between the tangent direction of the fiber control point and the global main direction.2 The value determines the waviness characteristics of the fiber itself.

[0032] Furthermore, in step 2, the process of rotating each fiber is as follows:

[0033] Define e as and μ 1 =(0,0,1) normal vector, α is With μ 1 The deviation angle of the fiber before rotation is k (k = 1, ..., r t ) spheres Relative to the first sphere on the fiber Direction; Create a rotation matrix

[0034] D(ω,e,α)=(e·ω)e+cosα((e×ω)×e)+sinα(e×ω)

[0035] The coordinates of the kth sphere on the fiber after rotation are the sum of its rotation matrix and the coordinates of the first sphere on the fiber.

[0036] Furthermore, in step 3, the process of setting the periodic condition is:

[0037] For a sphere whose center is beyond the geometric range of RVE, a new sphere is generated at its symmetric position to ensure the geometric periodicity of RVE; take the periodicity parameter

[0038]

[0039] Then the coordinates of the newly generated periodic sphere at the boundary are

[0040] For the spheres at the four corner points of the plane, spheres with corresponding periodic conditions are simultaneously created at the corresponding positions of the other three corner points.

[0041] Furthermore, in step 5, the friction factor

[0042]

[0043] δ represents the change in length of the spring due to extension or compression.

[0044] Furthermore, in step 5, the corresponding friction factor δ in the restoring force of the tension and compression spring is min and δ max Taken as 1% and 2% respectively, the corresponding friction factor δ in the bending spring restoring force min With δ max Take 0.28% and 0.56% respectively.

[0045] Beneficial Effects

[0046] The modeling method of unidirectional fiber periodic composite materials with fiber waviness proposed in the present invention includes the creation of a soft-core fiber system and the correction of a hard-core fiber system. The soft-core fiber system refers to a fiber system with only waviness characteristics, without considering the overlap between the fibers, that is, allowing overlap between the fibers. The process explains how to incorporate the structural parameters of the fiber system into the random walk algorithm and how to control the bending of the fibers. Subsequently, the geometric periodicity of the entire fiber system is set. The hard-core fiber system is further adjusted on the basis of the soft-core system to eliminate the overlap between the fibers without destroying the fiber structure. The adjustment process introduces a force bias algorithm, that is, the repulsion, attraction and bending force between the fibers, and applies it to the created soft-core fiber system. In addition, a three-dimensional neighbor list is introduced to accelerate the force bias algorithm in the adjustment process of the hard-core system, and the implementation details are given. Through the present invention, a periodic three-dimensional representative volume element model of a unidirectional fiber-reinforced composite material with the required fiber waviness can be generated, which helps to understand the relationship between the component properties / microstructure and macroscopic properties of the material at a deeper level.

[0047] Additional aspects and advantages of the present invention will be given in part in the following description and in part will be obvious from the following description, or will be learned through practice of the present invention. BRIEF DESCRIPTION OF THE DRAWINGS

[0048] The above and / or additional aspects and advantages of the present invention will become apparent and easily understood from the description of the embodiments in conjunction with the following drawings, in which:

[0049] Figure 1 Algorithm framework diagram

[0050] Figure 2 Control point space parameters of a single fiber curve

[0051] Figure 3 Schematic diagram of the creation of a single fiber

[0052] Figure 4 RVE periodic boundary condition settings

[0053] Figure 5 Schematic diagram of 3D neighbor list

[0054] Figure 6 Force bias algorithm diagram

[0055] Figure 7 Spring friction factor

[0056] Figure 8 RVE model for different fiber waviness DETAILED DESCRIPTION

[0057] According to the required size of RVE, fiber radius and fiber volume fraction, a random perturbation method is used to create random two-dimensional fiber position points. A random walk algorithm is used at each position point to create each fiber in the system, and the fiber is characterized by a Bézier curve. In the random walk algorithm, in order to ensure the structural authenticity between the fibers, that is, to avoid local excessive buckling, the tangent direction at each control point is made to obey the given orientation distribution function von Mises-Fischer (vMF) distribution. (In this algorithm, μ represents the axial vector of UD-FRP, u represents the direction vector at the tangent of each control point of the Bézier curve, and k is the reliability parameter corresponding to u relative to μ, which determines the correlation between u and μ).

[0058] Furthermore, the tangent direction at each control point can also follow the multivariate von Mises-Fisher (vMF) distribution in The multivariate von Mises-Fisher distribution is a simple extension of the previous von Mises-Fisher distribution, μ 1 =(0,0,1) represents the axial direction of UD-FRP, μ 2 is the direction vector representing the tangent direction of each control point, and its modulus is 1, that is, |μ 2 |=1. k 1 With k 2 The tangent directions of each control point of the fiber are With the global principal direction μ 1 =(0,0,1) and the tangent direction of the previous control point of the same fiber The corresponding correlation parameter, k 1 The value determines the degree of deviation between the tangent direction of the fiber control point and the global main direction. 2 The value will determine the waviness characteristics of the fiber itself. The specific process of generating the periodic three-dimensional RVE model of UD-FR P with different fiber waviness is as follows:

[0059] 1. Generate tangent direction vectors of each control point of the fiber curve like Figure 2 As shown, r is the spacing between control points, which is taken as the radius of the fiber. The RVE contains n fibers in total. represents the axial direction of the jth fiber. Define w i (i=1,2,3) is μ 2 The size of the projection component on each axis of the coordinate system, that is, μ 2 =(w 1 ,w 2 ,w 3 ), so there is For each fiber curve, to ensure And obey the multivariate vMF distribution, define w 3 =M (M and N are both independent random variables, N is a random number on the unit circle ), M follows the density formula of the von Mises-Fisher distribution (where k = |k 1 μ 1 +k 2 μ 2 |). Calculate the cumulative distribution function of the vMF distribution To generate random variables that follow the vMF distribution, calculate the inverse of the cumulative distribution function Take y as a random number on (0,1) and substitute it into the above formula to make M = F -1 (y). At this point, we get M and N, and we can calculate them one by one to get the direction vector of the tangent line of each control point of the fiber curve.

[0060] According to the given geometric dimensions L of the RVE model i (i=1,2,3), the number of control points on each fiber in this algorithm is defined as s=(L 3 / r)+1, the coordinates of each control point on the fiber can be obtained through the above random walk algorithm, and the Bézier curve, that is, each fiber curve, can be further generated.

[0061] 2. Take the average intercept r of each fiber Bézier curve along the global Z direction (i.e., the RVE axial direction) t coordinates, respectively assigned to r t sphere (r t The value of is selected as appropriate. The larger the value, the better the smoothness of the final generated fiber, but the greater the amount of algorithm calculation). That is, the ball chain is used to specifically characterize the curve of each fiber. 2 The parameters in the final fiber orientation have a certain randomness. The axial direction μ of UD-FRP 1 =(0,0,1) produces a certain deviation. Therefore, in order to ensure the global orientation distribution of UD-FRP, it is necessary to rotate the fiber as a whole to generate a new fiber curve while keeping the structure (length, radius, curvature) of the fiber generated in the previous step unchanged.

[0062] like Figure 3 As shown, define e as and μ 1 =(0,0,1) normal vector, α is With μ 1 The deviation angle of the fiber before rotation is k (k = 1, ..., r t ) spheres Relative to the first sphere on the fiber direction, we can calculate Then create the rotation matrix

[0063] D(ω,e,α)=(e·ω)e+cosα((e×ω)×e)+sinα(e×ω)

[0064] Therefore, the coordinates of the kth sphere on the fiber after rotation are the sum of its rotation matrix and the coordinates of the first sphere on the fiber. The above operation rotates each fiber in the RVE model generated in the previous step so that these fibers have a uniform axial direction.

[0065] 3. Set periodic boundary conditions for the RVE model. First, determine the coordinates of the fiber ball chain points built in the previous step. Place all the spheres of the RVE in the z direction with a ball chain spacing L. 3 / (r t -1) divided into r t Plane P i (i=1,,,r t ), set the periodic conditions for each ball-chain point on each plane separately, such as Figure 4 As shown. For the sphere in each plane, we need to consider n sphere centers ( Figure 4 The coordinates of the medium-light sphere (i=1,2, i.e. global X and Y directions respectively; k=1,…,n) exceeds the geometric range L of RVE i When , a new sphere should be generated at its symmetric position to ensure the geometric periodicity of RVE, and the coordinates of other spheres completely located within the RVE range do not need to be changed. Then the coordinates of the periodic sphere corresponding to the newly generated fiber sphere at the boundary are Spheres located at the four corners of the plane ( Figure 4 The dark corner point sphere) should correspond to the spheres at the other three corner points and ensure their periodicity, that is, spheres corresponding to the periodic conditions are created at the corresponding positions of the other three corner points. t The above-mentioned periodic arrangement of the spheres on the plane can ensure the overall geometric periodicity of the soft-core fiber system.

[0066] 4. In the RVE model with periodic boundary conditions created in the previous step, overlapping is allowed between fibers. This step will create a three-dimensional neighbor list for each fiber sphere, in preparation for the next step of applying a force bias algorithm between overlapping fiber spheres in the model to separate the overlapping fibers.

[0067] Divide the space into several small cubic areas, such as Figure 5As shown. The side length of the small cubic area is not less than the fiber radius, and each vertex is set as a node. The coordinates of all the sphere center points in the area are q i 'Round down and assign to the corresponding node q i Place (q i '∈q i ), which is convenient for narrowing the search range. Take the center point of each sphere as the origin of its three-dimensional neighbor list (the center point of the neighbor list), then the neighbor list includes all neighbor spheres that may overlap with the origin. At the same time, the neighbor lists of the outer fiber ball chain and other internal fiber ball chains that determine the geometric periodicity of RVE can be stored separately. First, the hard core system creation process is performed on the outer fiber ball chain, and then its center coordinates are fixed, and the process is performed on other internal fiber ball chains. In this way, the neighbor list of each fiber sphere is created without destroying the geometric periodicity of RVE.

[0068] 5. After generating the three-dimensional neighbor list, the neighboring sphere centers of all sphere centers on the target fiber sphere chain can be searched in each sphere neighbor list to determine whether there is overlap. Repulsive force is applied to the overlapped fiber spheres to restore their rigid contact. The position points of the two fiber spheres are q 1 ,q 2 , the radius is r, the repulsive force between the two spheres is related to their overlap Δl, Δl = |2r-||q 1 -q 2 |||, then the ball q 1 The repulsive force is Where K r is the repulsive force constant of the ball chain, which can be defined by relevant material parameters and set according to actual conditions. 2 The total repulsive force on each sphere is the vector sum of the repulsive forces exerted on it by other spheres in its neighbor list. At the same time, in order to prevent the fiber from destroying its own structure due to external forces during the creation of the hard core system, tension and compression spring restoring forces and bending spring restoring forces are applied between the spheres in the same fiber ball chain, such as Figure 6 shown.

[0069] The tension and compression spring limits the axial extension and compression between the ball chains, and its elastic constant is K T In order to ensure the compactness of the ball chain, the ball q i Not only affected by the two adjacent balls q i-1 and q i+1 The spring force is also affected by the two balls q i-2 and q i+2 The spring force, adding the friction factor Limit the spring restoring force of the fiber and stabilize the movement of the ball chain. When the spring changes length δ due to elongation or compression, minWhen the spring length change gradually increases to δ max When the spring restoring force reaches the normal level, Figure 7 As shown, its mathematical expression is The magnitude of the spring force is linearly related to the change in the length of the spring. i The next ball on the ball chain q i+1 The restoring force of the tension and compression spring is The friction factor δ in this model is min and δ max They are taken as 1% and 2% respectively (can be adjusted according to actual needs). i At the same time, it is also affected by the sphere q i-2 ,q i-1 and q i+2 By analogy, we can get the spring restoring force of tension and compression on all balls in the whole ball chain system.

[0070] In order to ensure the smoothness of the fiber structure, a bending spring restoring force is applied at the same time. Assume that the initial angle between the three balls is π and the angle after the change is γ, where the ball q i To the ball i-2 , ball q i+2 The foot of the perpendicular line is P, and the angle force and the sphere q are defined. i The distance from the initial coordinate is linearly related, and considering the influence of the friction factor, the ball q i By the ball i-2 and ball q i+2 The restoring force of the bending spring is Ball Q i-2 and ball q i+2 At the same time, the ball i The angle force acts in the direction of On the contrary, the size is half of it, and similarly the ball q i Also affected by the ball i-1 and ball q i+ 1 angle force. For the model in this paper, the bending spring friction factor δ min With δ max Take 0.28% and 0.56% respectively, the bending spring constant K B Definition and K T The same is given by the fiber parameters required for practical applications.

[0071] In this step, the repulsive force, spring force and angle force between the spheres in the neighbor list are defined. The total force on a sphere is the vector sum of these forces applied to it by other spheres in its neighbor list. By applying the above forces to the external fiber ball chain and the internal fiber ball chain respectively, the coordinates of each ball chain with periodic geometric conditions are generated, and then a UD-FRP high-fidelity RVE model containing fiber waviness characteristics is created.

[0072] 6. In order to quickly convert the digital model of the fiber generated previously (the coordinates of the center of the ball on the ball chain and the radius of each ball chain) into a solid model, this paper builds a framework based on the Microsoft Visual Studio carrier and adapts the modeling software development and integration. By reading the digital model of the fiber and importing it into the modeling software, the fiber is quickly solid modeled through curve fitting and sweeping operations, and finally a UD-FRP periodic three-dimensional solid model with different fiber corrugation characteristics is generated, such as Figure 8 shown.

[0073] Although the embodiments of the present invention have been shown and described above, it is to be understood that the above embodiments are exemplary and are not to be construed as limitations on the present invention. A person skilled in the art may change, modify, replace and modify the above embodiments within the scope of the present invention without departing from the principles and intent of the present invention.

Claims

1. A modeling method for unidirectional fiber periodic reinforced composite materials with fiber corrugation, Features: The following steps are involved: Step 1: According to the size of the representative volume element RVE and the requirements for the fiber radius and fiber volume fraction, a random perturbation method is used to create random two-dimensional fiber position points; a random walk algorithm is used at each position point to create each fiber in the system, and each fiber is characterized by a Bézier curve; Step 2: The average intercept r of each fiber Bézier curve along the global principal direction t coordinates, respectively assigned to r t spheres, so that each fiber curve can be represented by a ball chain; Under the premise of keeping the fiber structure unchanged, each fiber is rotated to make the fiber orientation μ 1 j and the global principal direction μ of the representative volume element 1 =(0,0,1) to be consistent, where j = 1,2,...,n, n is the number of two-dimensional fiber position points; Step 3: Arrange each ball chain in the global main direction with a ball chain spacing L 3 / r t -1 is divided into r t Plane P i ; There are n fiber spheres on each plane, and the periodic conditions are set for the ball-chain points on each plane separately; Step 4: Establish a three-dimensional neighbor list for all fiber spheres; Step 5: For each fiber sphere, search the neighbor list for the fiber sphere chains that overlap with it and the fiber spheres that generate the overlap; Apply repulsive force to the overlapping spheres between different fibers found K r is the repulsive force constant of the ball chain, which is calculated from the material parameters of the fiber, Δl is the overlap between the two balls, q 1 ,q 2 are the coordinates of the centers of the two balls; At the same time, spring force is applied between the balls in the same fiber, including tension and compression spring restoring force. and the restoring force of the bending spring where K T and K B are the tensile and compressive spring constant and the bending spring constant of the ball chain, K B = E f πr 3 where E f is the elastic modulus of the fiber, l 0 is the initial distance between the small balls, is the friction factor, q i is the coordinate of the center of the i-th ball on the fiber, and r is the radius of the sphere; Step 6: Import the final fiber digital model into the modeling software, and create a UD-FRP three-dimensional RVE solid model with fiber waviness and periodicity characteristics through curve fitting and sweeping operations.

2. A modeling method for a unidirectional fiber periodically reinforced composite material containing fiber corrugation according to claim 1, Features: In step 1, according to the given fiber waviness parameter and RVE size, the Bézier curve of a single fiber is Where s is the number of control points of the curve, which is defined by the user, 0≤t≤1, x i (i=0,1,2,...,s) represents the coordinates of each control point.

3. A modeling method for a unidirectional fiber periodically reinforced composite material containing fiber corrugation according to claim 1, Features: In step 1, the tangent direction at each control point of the Bézier curve obeys the given orientation distribution function von Mises-Fischer distribution Where μ represents the axial vector of the unidirectional fiber reinforced composite material, u represents the direction vector at the tangent of each control point of the Bézier curve, and k is the reliability parameter corresponding to u relative to μ, which determines the correlation between u and μ.

4. A modeling method for a unidirectional fiber periodically reinforced composite material containing fiber corrugation according to claim 1, Features: In step 1, the tangent direction at each control point of the Bézier curve follows the multivariate von Mises-Fisher distribution in μ 1 =(0,0,1) represents the axial direction of unidirectional fiber reinforced composite materials, μ 2 is the direction vector representing the tangent direction of each control point, and its modulus is 1; k 1 With k 2 The tangent directions of each control point of the fiber are With the global principal direction μ 1 =(0,0,1) and the tangent direction of the previous control point on the same fiber The corresponding correlation parameter, k 1 The value determines the degree of deviation between the tangent direction of the fiber control point and the global main direction. 2 The value determines the waviness characteristics of the fiber itself.

5. A method for modeling a unidirectional fiber periodically reinforced composite material containing fiber corrugation according to claim 1, Features: In step 2, the process of rotating each fiber is as follows: Define e as and μ 1 =(0,0,1) normal vector, α is With μ 1 The deviation angle of the fiber before rotation is k (k = 1, ..., r t ) spheres Relative to the first sphere on the fiber direction; Creating a rotation matrix D(ω,e,α)=(e·ω)e+cosα((e×ω)×e)+sinα(e×ω) The coordinates of the kth sphere on the fiber after rotation are the sum of its rotation matrix and the coordinates of the first sphere on the fiber.

6. A modeling method for a unidirectional fiber periodically reinforced composite material containing fiber corrugation according to claim 1, Features: In step 3, the process of setting the periodic condition is: For a sphere whose center is beyond the geometric range of RVE, a new sphere is generated at its symmetric position to ensure the geometric periodicity of RVE; take the periodicity parameter Then the coordinates of the newly generated periodic sphere at the boundary are For the spheres at the four corner points of the plane, spheres with corresponding periodic conditions are simultaneously created at the corresponding positions of the other three corner points.

7. A modeling method for a unidirectional fiber periodically reinforced composite material containing fiber corrugation according to claim 1, Features: In step 5, the friction factor δ represents the change in length of the spring due to extension or compression.

8. A method for modeling a unidirectional fiber periodically reinforced composite material containing fiber corrugation according to claim 7, Features: In step 5, the friction factor δ in the restoring force of the tension and compression spring is min and δ max Taken as 1% and 2% respectively, the corresponding friction factor δ in the bending spring restoring force min With δ max Take 0.28% and 0.56% respectively.

Images