Topological optimization method of fiber reinforced composite structure based on double-layer orthogonal unit

Abstract

The invention belongs to the technical field related to structure topological optimization, and discloses a topological optimization method of a fiber reinforced composite structure based on a double-layer orthogonal unit, which comprises the following steps: (1) carrying out grid discretization on a design domain of a to-be-optimized structure to obtain the double-layer orthogonal unit; (2) interpolating the elastic tensor of each double-layer orthogonal unit to obtain a unit stiffness matrix; (3) filtering and numerical relaxation are carried out on the density variable, and then updating of fiber orientation is carried out; (4) solving the sensitivity of the target function and the constraint function about the density variable, updating the density variable based on the obtained sensitivity information, judging whether the density variable converges or not, if so, outputting the topological configuration and the corresponding discrete fiber orientation, and turning to the step (5); otherwise, turning to the step (2); and (5) fitting the obtained fiber orientation based on Runge-Kutta integration to obtain a spatial continuous additive manufacturing fiber path. According to the invention, the stability and applicability are improved.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of structural topology optimization, and more particularly relates to a topology optimization method for a fiber-reinforced composite structure based on double-layer orthogonal units. BACKGROUND

[0002] A fiber-reinforced composite structure is composed of an isotropic matrix material and a fiber material with high modulus. Compared with traditional metal materials, the fiber-reinforced composite structure has excellent corrosion resistance, fatigue resistance, and high specific strength / stiffness, and has become a research hotspot in the industrial field. However, due to the constraints of design methods and manufacturing methods, most of the fiber-reinforced composite structures used in the industry are 0°, 45°, 90°, 135°, etc. Fiber layout, which can only change the orientation of the fibers in different composite plates, cannot fully exploit the high design flexibility and mechanical potential of the fiber-reinforced composite structure.

[0003] As a kind of efficient structural optimization method developed rapidly in recent years, structural topology optimization can fully exert the design freedom of the structure, and give the optimal structure topology according to the input conditions such as boundary conditions, material properties, constraints, etc. through scientific calculation, which is highly complementary to the high designability of fiber morphology, and becomes an effective design method for obtaining fiber-reinforced composite structures with excellent performance.

[0004] Currently, the topology optimization method of the fiber-reinforced composite material structure mainly includes two categories: one is the design method based on discrete fiber orientation, and the other is the topology optimization design method based on curve fiber modeling. In the design method based on discrete fiber orientation, the design domain is divided into multiple discrete units, and the fiber orientation of each unit becomes the object of optimization. However, due to the introduction of the rotation matrix, the optimization problem often presents strong non-convexity, which not only increases the complexity of the problem, but also makes the algorithm convergence extremely difficult. In addition, the discrete fiber orientation cannot be used for subsequent additive manufacturing. On the other hand, the topology optimization design method based on curve fiber modeling, although it can improve the continuity and smoothness of the design to a certain extent through the fitting of fiber path by using level set function or Bezier curve tools in the design stage, this method has certain limitations in design freedom. Since it relies too much on the pre-set curve model, this method may not be able to fully explore and utilize all possible design spaces, thereby limiting the possibility of obtaining higher performance fiber-reinforced composite structures. Finally, the performance of the fiber-reinforced composite structure is very sensitive to direction, and the direction of the load in engineering will have a certain deviation. When there is a slight deviation in the load or structure cooperation, it may cause a significant decrease in the mechanical properties of the fiber-reinforced composite structure or even failure. SUMMARY

[0005] In view of the above defects or improvement needs of the prior art, the present application provides a topology optimization method of a fiber-reinforced composite structure based on double-layer orthogonal units, which aims to solve the problems of low stability and applicability of the parallel design of the existing fiber orientation and structure topology.

[0006] To achieve the above-mentioned object, according to one aspect of the present application, a topology optimization method of a fiber-reinforced composite structure based on double-layer orthogonal units is provided, which comprises the following steps: Step one, the design domain of the structure to be optimized is discretized by mesh, and each mesh is divided into double-layer orthogonal units, the double-layer orthogonal units include first-layer orthogonal units and second-layer orthogonal units, the fiber orientation of the first-layer orthogonal units is consistent with the direction of the local principal stress, and the fiber orientation of the second-layer orthogonal units is consistent with the direction of the local secondary principal stress; Step two, the elastic tensor of each double-layer orthogonal unit is interpolated based on the density interpolation method of orthotropic anisotropy to obtain the element stiffness matrix of the double-layer orthogonal unit; the interpolated elastic tensor should be the density of the fiber material and the fiber orientation of two design variables; Step three, the density variable is filtered and numerically relaxed first, and then the fiber orientation is updated by using the method based on local stress analysis; or, the fiber orientation is updated by using the method based on local stress analysis first, and then the density variable is filtered and numerically relaxed; Step four, the sensitivity of the objective function and the constraint function corresponding to the topology optimization problem of the fiber-reinforced composite structure with respect to the density variable is solved, and then the density variable is updated based on the obtained sensitivity information, and it is judged whether it converges or not, if it converges, the topology configuration and the corresponding discrete fiber orientation are output, and step five is turned to; otherwise, step two is turned to; Step five, the obtained fiber orientation is fitted based on the Runge-Kutta integral to obtain the spatial continuous fiber path which can be additively manufactured.

[0007] Further, the elastic tensor of the fiber composite material after interpolation is: In the formula, D 0 is the elastic tensor of the fiber composite material along the horizontal direction, and T is the rotation matrix of the elastic tensor; since the principal stress and the secondary principal stress are always orthogonal, the element stiffness matrix of the double-layer orthogonal unit is established: Further, the density variable is numerically relaxed: The above two equations are respectively the relaxation operation of the design variable and the projection operation of the filtered density field , where represents the centroid coordinate of the i-th element, R is the filtering radius, is the projection threshold value; β is gradually doubled from 1 until reaching the maximum value 128.

[0008] The elastic tensor D of the unrotated fiber composite material 0 and the rotation matrix are written as: In the formula, E1 and E2 are respectively the elastic modulus of the fiber composite material in the horizontal direction and the vertical direction, μ 12 and μ 21 are respectively the Poisson's ratio in the corresponding direction, G 12 is the shear modulus of the fiber composite material; ρ e represents the pseudo-density corresponding to the element e; the penalty coefficient p=3; is the included angle between the fiber orientation and the horizontal direction.

[0009] Further, the stress tensor matrix is given by four values: In the formula, and are respectively the normal stress in the x direction and the y direction, and are equal in value, and respectively represent the shear force in the x direction and the y direction; The eigenvalue decomposition of the stress tensor matrix is obtained as: In the formula, and are respectively the principal stress and the secondary principal stress, and v1 and v2 are respectively the direction of the principal stress and the secondary principal stress.

[0010] Further, v1 and v2 are orthogonal and are solved by the following formula: .

[0011] Further, the update formula of the fiber orientation is: .

[0012] Further, the objective function, i.e., the structural flexibility of the continuous fiber reinforced composite structure based on the double-layer orthogonal sub-element, is specifically: Where u is the nodal displacement vector, k is the element stiffness matrix, and U and K are the corresponding global displacement matrix and stiffness matrix; The expression for the constraint function: Where N represents the total number of elements in the finite element model. Let frac represent the pseudo-density variable of the i-th cell, and frac is the upper limit of the allowed volume fraction.

[0013] Furthermore, the steps for fitting the obtained fiber orientation based on the Runge-Kutta integral are as follows: First, set the relative density of the seed points. Every -1 unit adds a seed unit, i.e., numbered 1, 1+ 1+2 ...until the Nth cell is identified as the seed cell, the middle seed cell is selected as the starting point of the first iteration, and the fiber orientation of any point in the design domain is obtained by inverse distance interpolation. Next, a fixed integration step size is set, the fiber orientation of the seed unit is selected as the initial input for fitting, and the coordinates of the next integration point are calculated using a curve fitting strategy based on the second-order Runge-Kutta method. This process is repeated until the structural boundary is reached, and a complete continuous fiber path is obtained. Finally, a merging threshold is given. After generating the first continuous fiber path, calculate the minimum distance from all seed points to that fiber path. If the minimum distance is less than the merging threshold, then... The seed points are merged with the fiber paths. The merged seed points will no longer participate in the generation of subsequent fiber paths. This process is repeated until all seed points are merged.

[0014] Furthermore, define It is a point For the fiber orientation at a given location, and for the fiber orientation of every four adjacent element centers, the fiber orientation at any point within this square domain can be obtained by inverse distance interpolation: The weighting factors satisfy: .

[0015] Furthermore, the integral coordinates are calculated using the following formula: Where k1 and k2 are based on coordinates The two estimates of fiber orientation are specifically written as follows: .

[0016] In summary, compared with the prior art, the topology optimization method for fiber-reinforced composite structures based on double-layer orthogonal units provided by this invention has the following advantages: 1. Each grid is subdivided into two layers of orthogonal units, which include a first layer of orthogonal units and a second layer of orthogonal units. The fiber orientation of the first layer of orthogonal units is consistent with the direction of the local principal stress, while the fiber orientation of the second layer of orthogonal units is consistent with the direction of the local secondary principal stress. In this way, the fibers can be orthogonally distributed according to the magnitude of the local principal stress and the secondary principal stress. The orthogonally distributed fibers are distributed in different layers, which not only further enhances the mechanical properties of the fiber-reinforced composite structure, but also avoids the difficulties brought about by fiber distribution in subsequent additive manufacturing, and improves stability and applicability.

[0017] 2. The fiber orientation obtained is fitted based on Runge-Kutta integral to obtain a spatially continuous additively manufactured fiber path. This method can automatically adjust the fiber content according to the relative density of the input fiber path and the merging threshold without any other human intervention, and has a high degree of automation.

[0018] 3. Determine the sensitivity of the objective function and constraint function to the density variable for topology optimization of the fiber-reinforced composite structure, and then update the density variable based on the obtained sensitivity information. Fiber orientation was performed using a method based on local stress analysis. This hybrid variable update method not only ensures the accuracy of density variable updates, but also improves the efficiency and stability of parallel optimization design of structural topology and fiber paths. Attached Figure Description

[0019] Figure 1 (a), (b), and (c) in the diagram represent the local element stress state, the principal stress and secondary principal stress of the element, and the orthogonal fiber orientation modeling diagram under the element, respectively. Figure 2 This is an interpolation diagram illustrating global fiber orientation; Figure 3 This is a schematic diagram of a curve fitting method based on the second-order Runge-Kutta method; Figure 4 This is a schematic diagram of an automatic seed point merging strategy; Figure 5 This is a flowchart of a topology optimization method for a fiber-reinforced composite structure based on a double-layer orthogonal unit provided by the present invention; Figure 6 This is a schematic diagram of the boundary conditions and design domain of the Michel beam; Figure 7 In the diagram, (a), (b), (c), (d), (e), and (f) represent the macroscopic topology and discrete fiber orientation design of the first layer, the macroscopic topology and fitted fiber path of the first layer, the macroscopic topology and discrete fiber orientation design of the second layer, the macroscopic topology and fitted fiber path of the second layer, the macroscopic topology and discrete fiber orientation design of the first and second layers, and the macroscopic topology and fitted fiber path of the first and second layers, respectively. Figure 8 In the figure, (a), (b), and (c) represent the macroscopic topology and fitted fiber path of the first layer, the macroscopic topology and fitted fiber path of the second layer, and the macroscopic topology and fitted fiber path of the first and second layers, respectively; the relative fiber density is 20, and the object is a Michel beam. Detailed Implementation

[0020] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention. Furthermore, the technical features involved in the various embodiments of this invention described below can be combined with each other as long as they do not conflict with each other.

[0021] Please see Figure 1 This invention improves a topology optimization method for fiber-reinforced composite structures based on double-layer orthogonal units. The topology optimization method mainly includes the following steps: Step 1: Discretize the design domain of the structure to be optimized into a mesh. Each mesh is subdivided into two layers of orthogonal elements. The two layers of orthogonal elements include a first layer of orthogonal elements and a second layer of orthogonal elements. The fiber orientation of the first layer of orthogonal elements is consistent with the direction of the local principal stress, while the fiber orientation of the second layer of orthogonal elements is consistent with the direction of the local secondary principal stress.

[0022] Specifically, the design domain is discretized into a mesh, and each mesh is subdivided into two layers of orthogonal elements. The fiber orientation of the first layer of orthogonal elements is aligned with the direction of the local principal stress, while the fiber orientation of the second layer of orthogonal elements is aligned with the direction of the local secondary principal stress. Thus, the fiber orientations of the two layers of orthogonal elements remain orthogonal. Figure 1 As shown.

[0023] Fiber composite laminate structures can generally be considered as 2.5D structures, so only the macroscopic topology and fiber paths of the fiber composite laminate structure need to be designed on a plane. Unlike traditional designs, in this embodiment, every two layers of fiber composite laminate are considered as a novel orthogonal composite laminate. The first fiber layer has its fibers arranged according to the principal stress direction, and the second fiber layer has its fibers arranged according to the secondary principal stress direction. The fibers of the two fiber composite laminates always remain orthogonal. The design domain is divided into two layers of orthogonal elements using Q4 elements.

[0024] Step two: Based on the density interpolation method of orthogonal anisotropy, interpolate the elastic tensor of each bilayer orthogonal element to obtain the element stiffness matrix of the bilayer orthogonal element; the interpolated elastic tensor should be the density of the fiber material. and fiber orientation A function of two design variables.

[0025] Specifically, in the material modeling stage, each bilayer orthogonal unit is considered as being filled with an orthogonally anisotropic fiber composite material, and the material's pseudo-density... and fiber orientation The design variables are considered as optimization variables; the elastic tensor of each bilayer orthogonal element is interpolated using the density method based on orthogonal anisotropy, and the interpolated elastic tensor should be the material density. and fiber orientation For a general two-dimensional problem, the interpolated elastic tensor of the fiber composite material, a function of two design variables, can be written in the following form: In the formula, D 0 Let be the elastic tensor of the fiber composite material with fibers distributed along the horizontal direction, and T be the rotation matrix of the elastic tensor. Due to the symmetry of the fiber material, Covering only an interval of length π is sufficient to satisfy all fiber orientations. To further enhance the mechanical properties and robustness of the fiber-reinforced composite structure, the fibers can be distributed along the principal stress and secondary principal stress directions, respectively. Since the principal stress and secondary principal stress are always orthogonal, the element stiffness matrix of a two-layer orthogonal element can be established: The elastic tensor D of unrotated fiber composites 0 The rotation matrix can be written as: In the formula, E1 and E2 are the elastic moduli of the fiber composite material in the horizontal and vertical directions, respectively, μ 12 and μ 21 These are the Poisson's ratios in the corresponding directions, and G. 12 This is the shear modulus of the fiber composite material. Due to the presence of μ...12 :μ 21 = E1: E2, therefore D 0 It is still a symmetric matrix. To avoid the singularity of the matrix, E is generally chosen. min =10 -9 . ρ e This represents the pseudo-density corresponding to unit e. In order to make the pseudo-density closer to 0 and 1, the penalty coefficient p is generally taken as 3. The angle between the fiber orientation and the horizontal direction is due to the symmetry of the fiber material. All fiber orientations can be satisfied by covering an interval of length π.

[0026] Step 3: Filter and numerically relax the density variable.

[0027] Generally speaking, It should be a set of Boolean variables taking values ​​of 0 and 1. However, such large-scale integer optimization problems cannot be solved using gradient-based methods, making numerical solutions difficult. Therefore, to make the optimization problem solvable and improve its efficiency and stability, numerical relaxation is applied to the density variables in this implementation: The two formulas above represent the relaxation operation on the design variables and the operation on the filtered density field, respectively. Perform projection operation, where This represents the centroid coordinates of the i-th unit, and R is the filtering radius. This is the projection threshold. In this embodiment, it is set to 0.5, and β gradually doubles from 1 until it reaches the maximum value of 128.

[0028] Step four: Fiber orientation is performed using a method based on local stress analysis. Update.

[0029] For fiber orientation If a gradient-based solver is used for updating, the introduction of a rotation matrix T containing many trigonometric function variables increases the non-convexity of the solution, leading to instability and difficulty in achieving convergence during iteration. To address these issues, this implementation method employs a local stress analysis-based approach for fiber orientation. For the stress state of a general two-dimensional problem, the stress tensor matrix can be characterized by four values: In the formula, and These are the normal stresses in the x and y directions, respectively. and Numerically they are equal, representing the shear forces in the x and y directions, respectively.

[0030] Eigenvalue decomposition of the stress tensor matrix yields: In the formula, and v1 and v2 represent the magnitudes of the principal and secondary principal stresses, respectively, and their directions, respectively. It's important to note that after eigenvalue decomposition, v1 and v2 are orthogonal, and can be solved using the following formula: Therefore, fiber orientation The update can be performed using the following formula: In other methods, the order of steps three and four can be reversed.

[0031] Step 5: Calculate the sensitivity of the objective function and constraint function to the density variable in the topology optimization problem of the fiber-reinforced composite structure, and then update the density variable based on the obtained sensitivity information. Then determine whether it converges. If it converges, output the topology and the corresponding discrete fiber orientation, and go to step six; otherwise, go to step two.

[0032] For density variable Solving and updating requires providing the objective function and constraint functions with respect to the density variables in advance. The expression for the objective function, i.e., the structural compliance of the continuous fiber-reinforced composite structure based on bilayer orthogonal subunits, can be written as: Where u is the nodal displacement vector, k is the element stiffness matrix, and U and K are the corresponding global displacement matrix and stiffness matrix.

[0033] The specific assembly of the element stiffness matrix can be performed according to the following formula: In the formula, B is the strain-displacement matrix in finite element analysis, which is constructed based on the shape function of the element and the derivative of the displacement field. It can be calculated, and thus the structural flexibility of the entire fiber-reinforced composite structure can be obtained. Generally, in topology optimization design, the optimization objective is to minimize the structural flexibility.

[0034] Generally, the constraint in a topology optimization problem is the amount of all materials used, and the expression for the constraint function can be given as follows: Where N represents the total number of elements in the finite element model. Let frac represent the pseudo-density variable of the i-th cell, and frac is the upper limit of the allowed volume fraction.

[0035] This optimization problem employs gradient-based MMA to calculate the density variable. The update, the sensitivity of the objective function and constraint function to the density variable, can be written as: The sensitivity of the intermediate variable can be obtained by the following formula: The sensitivity of the constraint function to the density variable can be written as: Based on the sensitivity information of the obtained objective function and constraint function, the density variable is updated using the MMA method. After each update, check if convergence has occurred. If convergence has occurred, directly output the optimal topology and the corresponding discrete fiber orientation, and proceed to step six; if convergence has not occurred, return to step two, as shown below. Figure 5 As shown.

[0036] Step 6: Fit the obtained fiber orientation based on the Runge-Kutta integral to obtain a spatially continuous additively manufactured fiber path.

[0037] The steps for fitting the obtained fiber orientation based on Runge-Kutta integral are as follows: First, set the relative density of the seed points. Every -1 unit adds a seed unit, i.e., numbered 1, 1+ 1+2 ...until the Nth cell is identified as the seed cell, the middle seed cell is selected as the starting point of the first iteration, and the fiber orientation of any point in the design domain is obtained by inverse distance interpolation.

[0038] Generally, for a rectangular design domain, the seed cell for the initial iteration is located at the center of the rectangle. Before fitting the discrete fiber orientation, it is necessary to interpolate the fiber orientation of the discrete cells to obtain the fiber orientation at any position globally. Definition It is a point The fiber orientation at each location, for the fiber orientation of every four adjacent element centers, such as Figure 2 As shown, the fiber orientation at any point within this square domain can be obtained through inverse distance interpolation: The weighting factors satisfy: Next, a fixed integration step size is set, the fiber orientation of the seed unit is selected as the initial input for fitting, and the coordinates of the next integration point are calculated using a curve fitting strategy based on the second-order Runge-Kutta method. This process is repeated until the structural boundary is reached, resulting in a complete continuous fiber path.

[0039] Specifically, the integral coordinates can be calculated using the following formula: Where k1 and k2 are based on coordinates The two estimates of fiber orientation can be specifically written as: The continuous curve integral process based on the second-order Runge-Kutta method is as follows: Figure 3 As shown, at the integration point The fiber orientation at point k i-1 By searching in this direction with a step size of h, we can obtain... Next, calculate. place and The arithmetic mean of the fiber orientation at each location is used as the search direction from... Restart the search with a step size of h to obtain the next integration point. The coordinate information is obtained. This process is repeated until the structural boundary is reached. Similarly, at the integration point... By searching in the opposite direction of fiber orientation and repeating this operation until the structural boundary, a complete continuous fiber path can be obtained.

[0040] Finally, a merging threshold is given. After generating the first continuous fiber path, calculate the minimum distance from all seed points to that fiber path. If the minimum distance is less than the merging threshold, then... The seed points are merged with the fiber paths. The merged seed points (marked in yellow in the diagram) will no longer participate in the generation of subsequent fiber paths. This process is repeated until all seed points have been merged. Figure 4 As shown.

[0041] Among them, the merging threshold The value is equal to the smaller of the design domain's length and width divided by the relative fiber density. Seed points near existing fiber paths are merged according to the merging threshold, and their seed point characteristics are removed, meaning these seed points are no longer used as initial integration points for fiber path fitting. Then, the seed point closest to an existing fiber path is selected from the remaining seed points to generate a new fiber path, until all seed points are merged.

[0042] This invention proposes a topology optimization method for continuous fiber-reinforced composite structures based on two-layer orthogonal subunits. Orthogonal fiber paths are generated based on local principal and secondary principal stresses. Then, at the same seed point, the fiber orientation is rotated by 90°, and continuous fiber paths are generated again using the method described above, and the merging threshold is also applied. Merging seed points near fibers ensures a relatively uniform fiber path globally, facilitating subsequent additive manufacturing. The above two operations are considered a fitting operation. After a fitting operation is completed, the minimum distance from each remaining seed point to an existing fiber path is recalculated. The seed point with the smallest minimum distance is selected as the starting seed point for the next fitting operation, and the fitting operation is repeated until no new seed points remain in the entire design domain.

[0043] The present invention will be further described in detail below with reference to specific embodiments.

[0044] Michel beams have long been considered a classic example of topology optimization for fiber-reinforced composite structures, demonstrating the effectiveness of the proposed algorithm. To simplify the problem, all numerical values ​​used in the calculations are dimensionless. For example... Figure 6 As shown, the rectangular design domain has an aspect ratio of 4:1. A downward vertical load of F=1 is applied to the upper center point of the rectangular domain. The lower left node is completely fixed, while the lower right node is subject to vertical constraints. The design domain is subdivided into 200×50 Q4 elements, each of which is further divided into two layers of orthogonal sub-elements for arranging orthogonally distributed fiber reinforcement materials. The Young's modulus of the matrix material is 1, and the Young's modulus of the fiber reinforcement is 10. The upper limit of the overall structure volume is constrained to 0.7, and the optimization objective is to minimize the overall structural flexibility.

[0045] With the relative density of the fibers set to 10, the optimized macroscopic topology and fiber path of the Michel beam are as follows: Figure 6 As shown, the topology and fiber paths of the structure are optimized in parallel. Specifically, Figure 7 In the diagram, (a), (c), and (e) represent the topology-optimized structure, and it can be seen that each unit corresponds to a fiber orientation. Figure 7(b), (d), and (f) in the figure represent the structures after fitting the corresponding discrete fibers. The fitted results are uniformly and continuously distributed throughout the space, without any intersections or abrupt breaks. To further verify the stability of this fiber fitting method and its scalability regarding fiber volume fraction, the relative density of the fibers was set to 20, and continuous fitting of continuous fibers was performed. The fitting results are shown below. Figure 8 As shown, it can be observed that the number of fibers is relatively lower than that of... Figure 7 The fiber density is significantly increased, meeting diverse engineering needs. The entire process is fully automated, requiring no additional manual intervention. Furthermore, the fibers remain relatively uniformly distributed in space, without being excessively dense or sparse. The fiber distribution largely conforms to the directions of principal and secondary principal stresses, and there are no fiber intersections within each layer, providing traversal for subsequent additive manufacturing. Figure 8 As can be seen from (c), the fibers always remain orthogonal at the interlayer intersections, which meets the design requirements.

[0046] This invention provides a topology optimization method for continuous fiber-reinforced composite structures based on bilayer orthogonal elements. The traditional Q4 element is further divided into bilayer orthogonal sub-elements, ensuring that fiber orientation satisfies the directions of local principal and secondary principal stresses between different layers, thereby improving the mechanical properties and robustness of the fiber-reinforced composite structure. Furthermore, this invention employs a second-order Runge-Kutta method for fiber path fitting, allowing for the generation of continuous fiber paths with varying densities by inputting relative fiber density. This enables subsequent additive manufacturing, and the entire process is highly automated, meeting diverse engineering requirements.

[0047] Those skilled in the art will readily understand that the above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units, characterized in that, The method includes the following steps: Step 1: Discretize the design domain of the structure to be optimized into a mesh. Each mesh is subdivided into two layers of orthogonal elements. The two layers of orthogonal elements include a first layer of orthogonal elements and a second layer of orthogonal elements. The fiber orientation of the first layer of orthogonal elements is consistent with the direction of the local principal stress, while the fiber orientation of the second layer of orthogonal elements is consistent with the direction of the local secondary principal stress. Step two: Based on the density interpolation method of orthogonal anisotropy, interpolate the elastic tensor of each bilayer orthogonal element to obtain the element stiffness matrix of the bilayer orthogonal element; the interpolated elastic tensor should be the density of the fiber material. and fiber orientation A function of two design variables; Step 3: First, filter and numerically relax the density variable, then use a method based on local stress analysis to determine fiber orientation. Update; or, first use a method based on local stress analysis for fiber orientation. The density variable is then updated, and then filtered and numerically relaxed. Step four: Calculate the sensitivity of the objective function and constraint function to the density variable in the topology optimization problem of the fiber-reinforced composite structure, and then update the density variable based on the obtained sensitivity information. Then determine whether convergence has occurred. If convergence has occurred, output the topology and the corresponding discrete fiber orientation, and proceed to step five; otherwise, proceed to step two. Step 5: Fit the obtained fiber orientation based on the Runge-Kutta integral to obtain a spatially continuous additively manufactured fiber path.

2. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in claim 1, characterized in that: The elastic tensor of the interpolated fiber composite material is: ； In the formula, D 0 Let T be the elastic tensor of the fiber composite material with fibers distributed along the horizontal direction, and let T be the rotation matrix of the elastic tensor. Since the principal stresses and secondary principal stresses are always orthogonal, the element stiffness matrix of the bilayer orthogonal element is established: ； The elastic tensor D of unrotated fiber composites 0 The rotation matrix can be written as: ; In the formula, E1 and E2 are the elastic moduli of the fiber composite material in the horizontal and vertical directions, respectively, μ 12 and μ 21 These are the Poisson's ratios in the corresponding directions, and G. 12 It is the shear modulus of fiber composite materials; ρ e This represents the pseudo-density corresponding to element e; the penalty coefficient p = 3. It is the angle between the fiber orientation and the horizontal direction.

3. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in claim 1, characterized in that: Numerical relaxation was performed on the density variable: ； ； The two formulas above represent the relaxation operation on the design variables and the operation on the filtered density field, respectively. Perform projection operation, where This represents the centroid coordinates of the i-th unit, and R is the filtering radius. The projection threshold is β; β starts from 1 and gradually doubles until it reaches the maximum value of 128.

4. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in claim 1, characterized in that: The stress tensor matrix is ​​characterized by four values: ； In the formula, and These are the normal stresses in the x and y directions, respectively. and Numerically they are equal, representing the shear forces in the x and y directions, respectively; Eigenvalue decomposition of the stress tensor matrix yields: ； In the formula, and v1 and v2 are the magnitudes of the principal stress and secondary principal stress, respectively, and v1 and v2 are the directions of the principal stress and secondary principal stress, respectively.

5. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in claim 4, characterized in that: v1 and v2 are orthogonal, and can be solved using the following formula: ； 。 6. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in claim 5, characterized in that: Fiber orientation The update formula is: 。 7. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in any one of claims 1-6, characterized in that: The objective function, namely the structural compliance of the continuous fiber-reinforced composite structure based on bilayer orthogonal subunits, is as follows: ； Where u is the nodal displacement vector, k is the element stiffness matrix, and U and K are the corresponding global displacement matrix and stiffness matrix; The expression for the constraint function: ； Where N represents the total number of elements in the finite element model. Let frac represent the pseudo-density variable of the i-th cell, and frac is the upper limit of the allowed volume fraction.

8. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in any one of claims 1-6, characterized in that: The steps for fitting the obtained fiber orientation based on Runge-Kutta integral are as follows: First, set the relative density of the seed points. Every -1 unit adds a seed unit, i.e., numbered 1, 1+ 1+2 ...until the Nth cell is identified as the seed cell, the middle seed cell is selected as the starting point of the first iteration, and the fiber orientation of any point in the design domain is obtained by inverse distance interpolation. Next, a fixed integration step size is set, the fiber orientation of the seed unit is selected as the initial input for fitting, and the coordinates of the next integration point are calculated using a curve fitting strategy based on the second-order Runge-Kutta method. This process is repeated until the structural boundary is reached, and a complete continuous fiber path is obtained. Finally, a merging threshold is given. After generating the first continuous fiber path, calculate the minimum distance from all seed points to that fiber path. If the minimum distance is less than the merging threshold, then... The seed points are merged with the fiber paths. The merged seed points will no longer participate in the generation of subsequent fiber paths. This process is repeated until all seed points are merged.

9. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in claim 8, characterized in that: definition It is a point For the fiber orientation at a given location, and for the fiber orientation of every four adjacent element centers, the fiber orientation at any point within this square domain can be obtained by inverse distance interpolation: ； The weighting factors satisfy: 。 10. The topology optimization method for fiber-reinforced composite structures based on bilayer orthogonal units as described in claim 8, characterized in that: The integral coordinates are calculated using the following formula: ； Where k1 and k2 are based on coordinates The two estimates of fiber orientation are specifically written as follows: ； 。

Images