Topological optimization method and device for GPU acceleration fiber reinforced composite material structure

Abstract

The invention relates to the technical field of structure optimization, in particular to a topological optimization method and device for a GPU acceleration fiber reinforced composite material structure, and the method comprises the steps: obtaining parameters such as material attribute parameters and design domain size; 21 basic stiffness matrixes are calculated; respectively filtering the density variables and the angle variables; a geometric multi-grid preprocessing conjugate gradient method iteration solver is adopted to calculate the current displacement of each node of the densest grid level; and calculating an objective function according to the current design variable and the node displacement, performing sensitivity analysis on the objective function and the constraint condition to update the density and angle variables, judging whether a topological optimization iteration convergence condition is met or not, and outputting an optimization result if the topological optimization iteration convergence condition is met. Therefore, the problems that the stiffness matrix of each element in the topological optimization problem of the large-scale fiber reinforced composite material structure cannot be stored due to the limited memory of the global memory of the GPU, and the calculation efficiency of the stiffness matrix of the dynamic calculation element in the finite element solving process is low are solved.

Description

Technical Field

[0001] The present invention relates to the technical field of structural optimization, and particularly relates to a topological optimization method and device for a GPU-accelerated fiber-reinforced composite material structure. Background Art

[0002] In recent years, fiber-reinforced composite materials have become an important part of high-tech aerospace products such as aircraft, launch vehicles, and satellites. The application parts and usage amounts thereof have become one of the important indicators for measuring the advancement of the structure of aerospace vehicles. With the rapid development of 3D printing technology, it has become possible to process composite material components with complex structures and continuously varying fiber directions in space.

[0003] For the optimization design problems of large aerospace components such as wings, hatches, and shells, there are several orders of magnitude differences between the size of the overall structure and the current processing accuracy. To obtain a refined design result with smaller characteristic dimensions, it is necessary to perform high-resolution mesh division on the overall structure, which makes the scale of the optimization problem reach tens of millions or even billions of element numbers. During the process of using the topological optimization method to determine the material distribution in the design domain, the huge number of design variables and the solution of large linear systems in the finite element analysis pose huge challenges to the storage and computing capabilities of ordinary computers.

[0004] As a high-performance processor, GPU is widely used in the accelerated solution of large-scale problems. Compared with the CPU, although the GPU is designed to have more powerful parallel computing capabilities and has achieved remarkable results in the solution of finite element analysis problems of large-scale isotropic materials, the global memory of the GPU is extremely limited, which is an undeniable shortcoming for the topological optimization problem of fiber-reinforced composite materials. Due to the continuously varying fiber direction in space, the element stiffness matrices in each finite element model are different, and storing each element stiffness matrix consumes a huge amount of the GPU's global memory; if not stored but using the traditional Gaussian integration method to dynamically calculate each element stiffness matrix, there will be a large number of repeated accumulation calculations, resulting in a significant reduction in the performance of the optimization algorithm. That is to say, currently, it is difficult to solve the high-resolution three-dimensional fiber-reinforced composite material topological optimization problem through the method of GPU parallel acceleration. Summary of the Invention

[0005] The present invention provides a topological optimization method and device for a GPU-accelerated fiber-reinforced composite material structure to solve the problems that due to the limited memory of the GPU global memory, the element stiffness matrices in the topological optimization problem of large-scale fiber-reinforced composite material structures cannot be stored, and the calculation efficiency of dynamically calculating the element stiffness matrix using the traditional Gaussian integration method during the finite element solution process is low, etc.

[0006] An embodiment of the first aspect of the present invention provides a topology optimization method for a GPU-accelerated fiber-reinforced composite material structure, including the following steps: obtaining the property parameters of a preset fiber-reinforced composite material and the element size of a target finite element model; establishing a topology optimization model, initializing the design variables of the topology optimization model, and defining the constraint conditions of the topology optimization model, where the design variables include density design variables and angle design variables; calculating 21 basic stiffness matrices according to the element size of the target finite element model, and storing the 21 basic stiffness matrices in a preset GPU global memory; performing an optimization iteration process based on the topology optimization model, filtering the density design variables and the angle design variables to obtain filtered density design variables and filtered angle design variables; calculating the constitutive relation matrix after rotation of each element according to the property parameters of the preset fiber-reinforced composite material, the filtered density design variables, and the filtered angle design variables; adopting an inter-element parallel manner, calculating the template stiffness matrix of each node at each grid level according to the template stiffness matrix solution format and the template stiffness matrix transfer relationship between each grid level, and storing the template stiffness matrices of all other levels except the densest grid level in the preset GPU global memory, where when calculating the stiffness matrix of each element, calculating the stiffness matrix of each element according to the 21 basic stiffness matrices and the constitutive relation matrix after rotation of each element; calculating the preprocessed smoothed residual of the densest grid level by using the "V" cycle algorithm of geometric multigrid, and processing the preprocessed smoothed residual of the densest grid level based on the preconditioned conjugate gradient method to obtain the current displacements of each node at the densest grid level; solving the objective function according to the current design variables and the current displacements of each node at the densest grid level, and performing a sensitivity analysis on the objective function value and the constraint conditions to obtain sensitivity values; updating the current design variables according to the objective function value and the sensitivity values to obtain a new optimized design, determining whether the current optimization iteration times reach a preset threshold, if so, ending the iteration and outputting the new optimized design, otherwise continuing the iteration.

[0007] Optionally, the property parameters of the preset fiber-reinforced composite material include at least one of longitudinal elastic modulus, transverse elastic modulus, Poisson's ratio, or in-plane shear modulus.

[0008] Optionally, the 21 basic stiffness matrices are:

[0009]

[0010] Wherein, is the basic stiffness matrix, n is n = i(i - 1) / 2 + j, B e is a shape function matrix composed of 8 submatrices, Dij is a 6×6 matrix with the element in the i-th row and j-th column replaced by 1 and elements in other positions replaced by 0, Ω e is the integration region represented by each unit, is D ji is a 6×6 matrix with the element in the j-th row and i-th column replaced by 1 and elements in other positions replaced by 0.

[0011] Optionally, filtering the density design variable and the angle design variable to obtain a filtered density design variable and a filtered angle design variable includes:

[0012] Filtering the density design variable by using a preset filtering method to obtain the filtered density design variable, where the preset filtering method is:

[0013]

[0014] where is the filtered density design variable, ω is a weight function, ω(x k ) is ω(x k ) = max(R - ||x k - x e ||, 0), Φ e is a set of all units whose distance from the center point to the center point of the e-th unit is not greater than the filtering radius, V k is the volume of the k-th unit whose distance from the center point to the center point of the e-th unit is not greater than the filtering radius, ρ k is the density design variable of the k-th unit whose distance from the center point to the center point of the e-th unit is not greater than the filtering radius, N e is the number of units;

[0015] According to trigonometric functions, converting the angle design variables α e and θ e in the spherical coordinate system into a vector in the Cartesian coordinate system

[0016] Filtering the vector in the Cartesian coordinate system respectively according to the preset fiber direction filtering radius to obtain a filtered direction vector

[0017] Normalizing the filtered direction vector and transforming the normalized direction vector according to the inverse trigonometric function to obtain the filtered angle design variable.

[0018] Optionally, the relationship for transferring the template stiffness matrix between each mesh level is:

[0019]

[0020] where, ω j =(2 - |j1|)(2 - |j2|)(2 - |j3|) / 8, ω k =(2 - |k1|)(2 - |k2|)(2 - |k3|) / 8, are both the weights for the transfer of the template stiffness matrix between different mesh levels, is the template stiffness matrix of each node on a coarser mesh level, is the template stiffness matrix of each node on a denser mesh level.

[0021] Optionally, calculating the smoothed residual after preprocessing at the densest mesh level by using the "V" - cycle algorithm of geometric multigrid, and processing the smoothed residual after preprocessing at the densest mesh level based on the preconditioned conjugate gradient method to obtain the current displacements of each node at the densest mesh level, including:

[0022] Based on the parallel - between - elements method, calculating the initial residual at the densest mesh level according to the template stiffness matrix at the densest mesh level, the preset initial displacement solution, and the preset load of each node;

[0023] Using the restriction method to transfer the initial residual at the densest mesh level to the coarsest mesh level to obtain the restricted residual at the coarsest mesh level, and using the direct method to solve the restricted displacements of each node at the coarsest mesh level according to the global stiffness matrix at the coarsest mesh level and the restricted residual at the coarsest mesh level;

[0024] By using the prolongation method to transfer the restricted displacements of each node at the coarsest mesh level to the densest mesh level to obtain the prolonged displacements of each node at the densest mesh level, and calculating the prolonged residual at the densest mesh level according to the template stiffness matrix at the densest mesh level and the prolonged displacements of each node at the densest mesh level;

[0025] Using the post - smoothing method to process the prolonged residual at the densest mesh level to obtain the smoothed residual after preprocessing at the densest mesh level;

[0026] Using the preconditioned conjugate gradient method to obtain the current displacements of each node at the densest mesh level according to the smoothed residual at the densest mesh level.

[0027] Optionally, updating the current design variables according to the objective function value and the sensitivity value to obtain a new optimized design, including:

[0028] Based on the optimization criterion method, updating the density design variables in the current design variables according to the objective function value and the sensitivity value to obtain new density design variables;

[0029] Based on the moving asymptote method, update the angular design variables in the current design variables according to the objective function value and the sensitivity value to obtain the angular design variables, where the new optimal design includes the new density design variables and the angular design variables.

[0030] An embodiment of the second aspect of the present invention provides a GPU-accelerated topology optimization device for fiber-reinforced composite material structures, including: an acquisition module for acquiring the property parameters of a preset fiber-reinforced composite material and the element size of a target finite element model; an initialization module for establishing a topology optimization model, initializing the design variables of the topology optimization model, and defining the constraint conditions of the topology optimization model, where the design variables include density design variables and angular design variables; a first calculation module for calculating 21 basic stiffness matrices according to the element size of the target finite element model and storing the 21 basic stiffness matrices in a preset GPU global memory; a filtering module for performing an optimization iteration process based on the topology optimization model to filter the density design variables and the angular design variables to obtain filtered density design variables and filtered angular design variables; a second calculation module for calculating the constitutive relation matrix after rotation of each element according to the property parameters of the preset fiber-reinforced composite material, the filtered density design variables, and the filtered angular design variables; a third calculation module for using an element-by-element parallel method to solve the template stiffness matrix of each node at each grid level according to the template stiffness matrix solution format and storing the template stiffness matrices of all other levels except the densest grid level in the preset GPU global memory, where when the stiffness matrix of each element is required during the calculation process, calculate the stiffness matrix of each element according to the 21 basic stiffness matrices and the constitutive relation matrix after rotation of each element; a fourth calculation module for using the "V" cycle algorithm of geometric multigrid to calculate the preprocessed smooth residual of the densest grid level and processing the preprocessed smooth residual of the densest grid level based on the preconditioned conjugate gradient method to obtain the current displacements of each node at the densest grid level; a solution and analysis module for solving the objective function according to the current design variables and the current displacements of each node at the densest grid level, and performing sensitivity analysis on the objective function value and the constraint conditions to obtain sensitivity values; an iterative optimization module for updating the current design variables according to the objective function value and the sensitivity values to obtain a new optimal design, determining whether the current optimization iteration times reach a preset threshold, if so, ending the iteration and outputting the new optimal design, otherwise continuing the iteration.

[0031] In a third aspect embodiment of the present invention, an electronic device is provided, including: a memory, a processor, and a computer program stored on the memory and executable on the processor, where the processor executes the program to implement the GPU-accelerated topological optimization method for fiber-reinforced composite material structures as described in the above embodiments.

[0032] In a fourth aspect embodiment of the present invention, a computer-readable storage medium is provided, where the computer-readable storage medium stores a computer program, and when the program is executed by a processor, it implements the above GPU-accelerated topological optimization method for fiber-reinforced composite material structures.

[0033] The GPU-accelerated topological optimization method and device for fiber-reinforced composite material structures proposed in the embodiments of the present invention can effectively solve problems such as the inability to store the stiffness matrices of each unit in the topological optimization of large-scale fiber-reinforced composite material structures due to limited memory in the GPU global memory, and the low calculation efficiency of dynamically calculating the unit stiffness matrix using the traditional Gaussian integration method in the finite element solution process. This enables the topological optimization problem of three-dimensional large-scale fiber-reinforced composite material structures to be accelerated and solved through the method of GPU parallel computing, greatly saving time costs, and obtaining a clear topological structure and continuously varying fiber directions in the optimization results.

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

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

[0036] Figure 1 is a flowchart of a GPU-accelerated topological optimization method for fiber-reinforced composite material structures provided by an embodiment of the present invention;

[0037] Figure 2 is a topological schematic diagram of a GPU-accelerated topological optimization method for fiber-reinforced composite material structures provided by an embodiment of the present invention;

[0038] Figure 3 is a schematic diagram of the definition method of angular variables provided by an embodiment of the present invention;

[0039] Figure 4 is a flowchart of the "V" cycle algorithm of the geometric multigrid method preprocessing process provided by an embodiment of the present invention;

[0040] Figure 5 is a schematic diagram of the working conditions of the topological optimization problem provided by an embodiment of the present invention;

[0041] Figure 6 The iterative history change diagram provided by the embodiment of the present invention;

[0042] Figure 7 The schematic diagram of the optimization result provided by the embodiment of the present invention;

[0043] Figure 8 The schematic diagram of the time comparison between the embodiment of the present invention and the traditional Gaussian integration method;

[0044] Figure 9 The block schematic diagram of the GPU-accelerated topology optimization device for fiber-reinforced composite material structures provided by the embodiment of the present invention;

[0045] Figure 10 The structural schematic diagram of the electronic device provided by the embodiment of the present invention. Detailed implementation manners

[0046] The embodiments of the present invention will be described in detail below. Examples of the embodiments are shown in the accompanying drawings, where the same or similar reference numerals denote the same or similar elements or elements having the same or similar functions throughout. The embodiments described below with reference to the accompanying drawings are exemplary and are intended to explain the present invention, and should not be construed as limiting the present invention.

[0047] The GPU-accelerated topology optimization method and device for fiber-reinforced composite material structures according to the embodiments of the present invention will be described below with reference to the accompanying drawings.

[0048] Specifically, Figure 1 The flowchart of a GPU-accelerated topology optimization method for fiber-reinforced composite material structures provided by the embodiment of the present invention.

[0049] As Figure 1 shown, the GPU-accelerated topology optimization method for fiber-reinforced composite material structures includes the following steps:

[0050] In step S101, obtain the property parameters of the preset fiber-reinforced composite material and the element size of the target finite element model.

[0051] Among them, the property parameters of the preset fiber-reinforced composite material include at least one of the longitudinal elastic modulus, the transverse elastic modulus, the Poisson's ratio, or the in-plane shear modulus.

[0052] Specifically, set parameters such as the property parameters of the fiber-reinforced composite material, the element size of the target finite element model, the Cartesian grid model resolution, and the geometric multigrid method level on the CPU side, store them in the CPU, and at the same time copy them to the GPU global memory.

[0053] In step S102, a topology optimization model is established, the design variables of the topology optimization model are initialized, and the constraint conditions of the topology optimization model are defined, where the design variables include density design variables and angle design variables.

[0054] Specifically, the topology optimization model is established on the GPU side as follows:

[0055] find: ρ, α, θ

[0056] min: c(ρ, α, θ) = F T U

[0057] s.t.: K(ρ, α, θ)U = F

[0058] : g k (ρ, α, θ) ≤ 0 (k = 1,..., N c )

[0059] : 0 ≤ ρ e ≤ 1 (e = 1,..., N e )

[0060] : -π ≤ α e ≤ π (e = 1,..., N e )

[0061] : -π / 2 ≤ θ e ≤ π / 2 (e = 1,..., N e )

[0062] where c is the structural strain energy, K is the global stiffness matrix, U is the nodal displacement vector, F is the global structural load vector, Ne is the number of elements, ρ is the density design variable, α and θ are the two angle design variables of each element, and g k is the volume constraint function.

[0063] Define the boundary conditions and loads of the topology optimization model, declare and initialize the density design variables and angle design variables in the topology optimization model, and their definition methods are as Figure 3 shown, as well as the field variables used to store the template stiffness matrix, preprocessing matrix, etc. corresponding to all nodes in each grid level.

[0064] In step S103, 21 basic stiffness matrices are calculated according to the element size of the target finite element model, and the 21 basic stiffness matrices are stored in the preset GPU global memory.

[0065] Among them, the 21 basic stiffness matrices are:

[0066]

[0067] Among them, is the basic stiffness matrix, and n = i(i - 1) / 2 + j, B e is the shape function matrix composed of 8 sub - matrices, D ij is a 6×6 matrix with the element in the i - th row and j - th column replaced by 1 and other elements replaced by 0, Ω e is the integration region represented by each element, is D ji is a 6×6 matrix with the element in the j - th row and i - th column replaced by 1 and other elements replaced by 0.

[0068] Among them, the matrix B ei can be expressed as:

[0069]

[0070] Due to the adoption of the SIMP method, D(ρ e ,α e ,θ e ) can be expressed as D(ρ e ,α e ,θ e ) = g(ρ e )D(α e ,θ e ), where usually take ε = 10 -6 and p = 3.

[0071] In step S104, based on the topology optimization model, an optimization iteration process is carried out to filter the density design variables and the angle design variables, and the filtered density design variables and the filtered angle design variables are obtained.

[0072] Specifically, at this time, the iteration process of topology optimization starts. In order to avoid numerical problems such as checkerboard phenomenon in the present invention, the density variables and angle variables are filtered first. Among them,

[0073] the density design variables are filtered by using a preset filtering method to obtain the filtered density design variables. The preset filtering method is:

[0074]

[0075] Among them, is the filtered density design variable, ω is the weight function, ω(x k ) is ω(x k ) = max(R - ||x k - x e ||, 0), Φ e is the set of all elements whose distance from the center point to the center point of the e - th element is not greater than the filtering radius, Vk is the volume of the k-th element whose distance from the center point to the center point of the e-th element is not greater than the filtering radius, ρ k is the density design variable of the k-th element whose distance from the center point to the center point of the e-th element is not greater than the filtering radius, N e is the number of elements;

[0076] According to trigonometric functions, transform the angular design variables α e and θ e in the spherical coordinate system into vectors in the Cartesian coordinate system where the trigonometric functions are: and

[0077] According to the preset fiber direction filtering radius, filter the vectors in the Cartesian coordinate system respectively to obtain the filtered direction vectors

[0078] Normalize the filtered direction vectors and transform the normalized direction vectors according to the inverse trigonometric functions to obtain the filtered angular design variables, where the inverse trigonometric functions are: and

[0079] In step S105, calculate the rotated constitutive relation matrix of each element according to the property parameters of the preset fiber-reinforced composite material, the filtered density design variable, and the filtered angular design variable.

[0080] Specifically, in each optimization iteration, first calculate the rotated constitutive relation matrix of each element according to the angular variable and store it in the global memory. The calculation method is as follows:

[0081]

[0082] where C is the constitutive relation matrix of the orthotropic material, which is expressed by the above material property parameters (the material parameters include longitudinal and transverse elastic moduli, Poisson's ratio, and in-plane shear modulus) as:

[0083]

[0084] The rotation matrix of the angular variable can be expressed as:

[0085]

[0086] where

[0087] In step S106, in a parallel manner between units, according to the solution format of the template stiffness matrix and the transfer relationship of the template stiffness matrix between each grid level, calculate the template stiffness matrix of each node at each grid level, and store the template stiffness matrices of all levels except the densest grid level in a preset GPU global memory. When calculating the stiffness matrix of each unit, calculate the stiffness matrix of each unit according to 21 basic stiffness matrices and the constitutive relationship matrix after rotation of each unit.

[0088] Among them, the stiffness matrix of each unit is:

[0089]

[0090] In the formula, n = i(i - 1) / 2 + j, is the elastic coefficient matrix of a certain unit on the densest grid level The element in the i-th row and j-th column.

[0091] Specifically, the calculation process of the template stiffness matrices of each level required by the assembly-free geometric multigrid preconditioning method is as follows: In a parallel manner between units, each node center at each level needs to store 27 block matrices of size 3×3 corresponding to surrounding nodes. The calculation process of the template stiffness matrices corresponding to 27 adjacent nodes of each node is respectively assigned to a separate GPU thread to complete. The 27 nodes respectively belong to 8 units adjacent to the central node. When the stiffness matrix of an adjacent unit is used to calculate its contribution to the template stiffness matrix of a certain node, it is obtained by multiplying and accumulating the 21 basic stiffness matrices with the elements in the corresponding position of the elastic coefficient matrix.

[0092] Among them, the transfer relationship of the template stiffness matrix between each grid level is:

[0093]

[0094] Among them, ω j =(2 - |j1|)(2 - |j2|)(2 - |j3|) / 8, ω k =(2 - |k1|)(2 - |k2|)(2 - |k3|) / 8, both are the weights for the transfer of the template stiffness matrix between different grid levels, is the template stiffness matrix of each node on a coarser grid level, is the template stiffness matrix of each node on a denser grid level.

[0095] In step S107, use the "V" cycle algorithm of geometric multigrid to calculate the smoothed residual after preprocessing of the densest grid level, and process the smoothed residual after preprocessing of the densest grid level based on the preconditioned conjugate gradient method to obtain the current displacements of each node on the densest grid level.

[0096] Specifically, as Figure 3 shown, based on the parallel mode between units, the initial residual of the densest grid level is calculated according to the stiffness matrix of the densest grid level template, the preset initial displacement solution, and the preset load of each node.

[0097] The initial residual of the densest grid level is transferred to the coarsest grid level by using a restriction method to obtain the restricted residual of the coarsest grid level, and a direct method is used to solve the restricted displacement of each node of the coarsest grid level according to the overall stiffness matrix of the coarsest grid level and the restricted residual of the coarsest grid level.

[0098] The restricted displacement of each node of the coarsest grid level is transferred to the densest grid level by using an extension method to obtain the extended displacement of each node of the densest grid level, and the extended residual of the densest grid level is calculated according to the stiffness matrix of the densest grid level template and the extended displacement of each node of the densest grid level.

[0099] The extended residual of the densest grid level is processed by using a post-smoothing method to obtain the smoothed residual after preprocessing of the densest grid level.

[0100] The preconditioned conjugate gradient method is used to process according to the smoothed residual of the densest grid level to obtain the current displacement of each node of the densest grid level.

[0101] In order to save video memory space, the stiffness matrix of the template of each node of the densest grid level is calculated by using a dynamic method, that is, when the stiffness matrix of the corresponding element of each unit is required, the method of multiplying and accumulating the elements of the basic stiffness matrix and the corresponding constitutive relation matrix proposed by the present invention is used to improve the calculation efficiency. The Jacobi method or the Gauss-Seidel method is used in the two smoothing processes.

[0102] In step S108, the objective function is solved according to the current design variables and the current displacement of each node of the densest grid level, and the sensitivity analysis is performed on the objective function value and the constraint conditions to obtain the sensitivity value.

[0103] Specifically, according to the current design variables and the current displacement of each node of the densest grid level, the objective function value of the optimized structure, that is, the structural strain energy, is solved; and the sensitivity analysis of the objective function and each constraint condition is performed to obtain the sensitivity value of each design variable.

[0104] In step S109, the current design variables are updated according to the objective function value and the sensitivity value to obtain a new optimized design, and it is judged whether the current optimization iteration times reach the preset threshold. If so, the iteration is ended and the new optimized design is output; otherwise, the iteration continues.

[0105] Specifically, based on the optimization criterion method, the density design variables are updated according to the objective function value and the sensitivity value to obtain new density design variables. Based on the moving asymptote method, the angle design variables are updated according to the objective function value and the sensitivity value to obtain angle design variables, and it is judged whether the current optimization iteration times reach a preset threshold. If so, the iteration is terminated and the new optimized design is output. Otherwise, the iteration continues and steps S104 - S109 are executed.

[0106] The following takes the three - dimensional cantilever beam structure obtained after topology optimization as an example to illustrate the finite - element analysis of large - scale fiber - reinforced composites. In this structure, the areas near both sides with a width of 1 / 4 each on the left - end face of the cubic initial design domain are fixed constraints. The height is 1.0 m, the length is 2.0 m, and the width is 1.0 m. The number of elements N after discretizing the design domain into a finite - element model e can be expressed as: N e = 128×256×128 = 4.1943×10 6 . The longitudinal elastic modulus E2 of the fiber - reinforced composite is 132×10 9 Pa, the transverse elastic modulus E1 is 9.6×10 9 Pa, the Poisson's ratio v 21 is 0.258, v 13 is 0.35, and the in - plane shear modulus G 21 is 5.8×10 9 . The structure is subjected to a vertical 20 - kN concentrated load at the mid - point of the lower end of the right - end face. With minimizing the structural strain energy as the optimization objective, a volume constraint is imposed, and the volume fraction is 10%.

[0107] Step 1, establish the mathematical model of the topology optimization model as follows:

[0108] find: ρ, α, θ

[0109] min: c(ρ, α, θ) = F T U

[0110] s.t.: K(ρ, α, θ)U = F

[0111] : g k (ρ, α, θ) ≤ 0 (k = 1, …, N c )

[0112] : 0 ≤ ρ e ≤ 1 (e = 1, …, N e )

[0113] : -π ≤ α e ≤ π (e = 1, …, N e )

[0114] : -π / 2 ≤ θe ≤π / 2 (e = 1, …, N e )

[0115] Where c is the structural strain energy, K is the global stiffness matrix, U is the nodal displacement vector, F is the global structural load vector, N e is the number of elements, ρ is the density design variable, α and θ are the two angular design variables of each element, and the volume constraint function g k can be expressed as follows:

[0116]

[0117] Where v e is the volume of each element, and ξ is the volume fraction.

[0118] Step 2, Initialize the design variables Both the density filtering radius and the angular filtering radius are set to 1.5 times the element size.

[0119] Step 3, Calculate 21 basic stiffness matrices, and the calculation method is as follows:

[0120]

[0121] Where is the basic stiffness matrix, n = i(i - 1) / 2 + j, k e is the element stiffness matrix, B e is the shape function matrix composed of 8 sub - matrices, B e = [B e0 B e1 B e2 B e3 B e4 B e5 B e6 B e7 , D ij is a 6×6 matrix with the element in the i - th row and j - th column replaced by 1 and other elements replaced by 0, Ω e is the integration region represented by each element, is D ji is a 6×6 matrix with the element in the j - th row and i - th column replaced by 1 and other elements replaced by 0.

[0122] Where the sub - matrix B ei can be expressed as:

[0123]

[0124] Due to the SIMP method, D(ρ e ,αe , θ e ) can be expressed as D(ρ e , α e , θ e ) = g(ρ e )D(α e , θ e ), where usually take ε = 10 -6 and p = 3.

[0125] Step 4, to avoid numerical problems such as checkerboard phenomenon, filter the density design variables and angle design variables. Among them, the density filtering is expressed as follows:

[0126]

[0127] Among them, is the filtered density design variable, ω is the weight function, ω(x k ) is ω(x k ) = max(R - ||x k - x e ||, 0), Φ e is the set of all elements whose distance from the center point to the center point of the e-th element is not greater than the filtering radius, V k is the volume of the k-th element whose distance from the center point to the center point of the e-th element is not greater than the filtering radius, ρ k is the density design variable of the k-th element whose distance from the center point to the center point of the e-th element is not greater than the filtering radius, N e is the number of elements.

[0128] The operation process of angle variable filtering is as follows:

[0129] (1) According to the following trigonometric function transformation: and Convert the angle design variables α e and θ e in spherical coordinates to a vector in Cartesian coordinates

[0130] (2) According to the preset fiber direction filtering radius, respectively filter to obtain the filtered direction vector

[0131] (3) Normalize the filtered direction vector to obtain

[0132] (4) According to the inverse trigonometric function relationship: and Calculate the filtered angle variables. ​

[0133] Step 5: In each optimization iteration, first calculate the constitutive relation matrix after rotation of each element according to the angular variables and store it in the global memory. The calculation method is as follows:

[0134]

[0135] Among them, C is the constitutive relation matrix of the orthotropic material, which is expressed by the above material property parameters as:

[0136]

[0137] The rotation matrix of the angular variables can be expressed as:

[0138]

[0139] Among them,

[0140] Step 6: Calculate the template stiffness matrices of each level required for the assembly-free geometric multigrid preconditioning method. The calculation process adopts the parallel method between elements. The calculation process of the template stiffness matrix corresponding to 27 nodes adjacent to each node is separately assigned to a single GPU thread to complete. The 27 nodes respectively belong to 8 elements adjacent to the central node. When the stiffness matrix of a certain adjacent element is used to calculate its contribution to the template matrix of a certain node, it is obtained by multiplying and accumulating the elements in 21 basic stiffness matrices and the elastic coefficient matrix at the corresponding positions. Among them,

[0141] The calculation formula of the element stiffness matrix is:

[0142]

[0143] The transfer relationship of the template stiffness matrix between each mesh level is:

[0144]

[0145] Among them, ω j =(2 - |j1|)(2 - |j2|)(2 - |j3|) / 8, ω k =(2 - |k1|)(2 - |k2|)(2 - |k3|) / 8, both are the weights for the transfer of the template stiffness matrix between different mesh levels. is the template stiffness matrix of each node on a coarser mesh level. is the template stiffness matrix of each node on a denser mesh level.

[0146] Step 7: Adopt the "V-cycle" algorithm shown in Figure 3 to perform the geometric multigrid preconditioning process. Its process includes:

[0147] (1) Calculate the template matrices and preconditioning matrices for each level except the densest level required by the assembly-free geometric multigrid method and store them in the GPU global memory;

[0148] (2) Smooth the densest level;

[0149] (3) Calculate the residual of the static equilibrium equation for the densest level grid;

[0150] (4) Transfer the residual to the next densest level as the right-hand side term of the static equilibrium equation through the restriction method, and loop the above process until the coarsest grid level;

[0151] (5) On the coarsest level, assemble the global stiffness matrix and solve the static equilibrium equation using the direct method based on Cholesky decomposition;

[0152] (6) Correct the displacement values of the next coarsest level through the prolongation method, and loop this step until the displacement solution of the densest level is obtained;

[0153] (7) Smooth the densest level;

[0154] In step seven, the calculation process of process (5) is completed using CPU calculation, and the displacement data of each node is copied to the GPU global memory. All other processes are completed using GPU parallel calculation. When the template stiffness matrix corresponding to each node is required in processes (1), (2), and (7), the template element stiffness matrix needs to be used to save storage space and improve calculation efficiency.

[0155] Step eight, solve the static equilibrium equation using the preconditioned conjugate gradient method to calculate the displacement of each node. The convergence tolerance of the conjugate gradient method iterative solver is set to 1×10 -5 .

[0156] Step nine, according to the element stiffness matrices calculated in step five and the displacements of each node calculated in step eight, solve the objective function of the optimized structure, that is, the structural strain energy; and perform sensitivity analysis of the objective function and each constraint condition;

[0157] Step ten, optimize and update the design variables: update the density design variables using the optimization criterion method, and update the angle design variables using the moving asymptotes method to obtain a new optimized design.

[0158] Step eleven, determine whether the iteration converges: determine whether the number of optimization iterations reaches 150 times. If it reaches 150 times, the convergence condition is satisfied, the iteration ends, and the optimized result is output; if the convergence condition is not satisfied, jump to step four to continue the iteration.

[0159] Next, the optimization design of the cantilever beam structure is carried out by using the method of the embodiment of the present invention and the large-scale topology optimization method of the traditional Gaussian integral online calculation unit stiffness matrix respectively. The embodiment of the present invention runs on a high-performance computing server equipped with Intel(R) Xeon(R) Gold 6240 and NVIDIA Tesla V100 GPU (32G video memory and 5120 CUDA cores). The topology optimization results are as Figure 8 shown, and the running time statistics of the optimization algorithm are shown in Table 1.

[0160] Table 1 Comparison of the results of the method of the embodiment of the present invention and the traditional method

[0161]

[0162]

[0163] The results show that the embodiment of the present invention can simultaneously optimize the density and fiber direction of large-scale fiber-reinforced composites, and obtain a clear topological structure and a continuously varying fiber direction design. At the same time, the calculation efficiency of the optimization algorithm proposed by the embodiment of the present invention has been greatly improved compared with the traditional Gaussian integral method, which greatly saves the time of the large-scale fiber-reinforced composite topology optimization problem, making the space cost of storing the anisotropic material unit stiffness matrix and the time cost of calculation reasonably balanced, which has important significance in practical applications.

[0164] In summary, according to the GPU-accelerated topology optimization method for fiber-reinforced composite structures proposed by the embodiment of the present invention, on the one hand, it greatly reduces the consumption of the GPU global memory for storing the fiber-reinforced composite unit stiffness matrix, and on the other hand, it dynamically calculates the stiffness matrix of each unit by multiplying and accumulating the elements in the elastic coefficient matrix with the basic stiffness matrix, and adopts the storage method of each grid level node by node and the displacement solution format without assembly, which greatly improves the solution efficiency of the finite element large-scale linear equations, and can be used in the fields of three-dimensional large-scale structure static analysis and structure topology optimization.

[0165] Next, a GPU-accelerated topology optimization device for fiber-reinforced composite structures proposed according to an embodiment of the present invention is described with reference to the accompanying drawings.

[0166] Figure 9 is a block diagram of a GPU-accelerated topology optimization device for fiber-reinforced composite structures provided by an embodiment of the present invention.

[0167] As Figure 9As shown, the topological optimization device 90 for a GPU-accelerated fiber-reinforced composite material structure includes: an acquisition module 901, an initialization module 902, a first calculation module 903, a filtering module 904, a second calculation module 905, a third calculation module 906, a fourth calculation module 907, a solution and analysis module 908, and an iterative optimization module 909. This device only stores the pre-computed basic stiffness matrices in the GPU global memory, solves the problem of insufficient GPU video memory for storing all element stiffness matrices, and can achieve the parallel topological optimization problem solving of fiber-reinforced composite material structures with low video memory consumption and high calculation efficiency.

[0168] Among them, the acquisition module 901 is used to acquire the property parameters of the preset fiber-reinforced composite material and the element size of the target finite element model. The initialization module 902 is used to establish a topological optimization model, initialize the design variables of the topological optimization model, and define the constraint conditions of the topological optimization model. Among them, the design variables include density design variables and angle design variables. The first calculation module 903 is used to calculate 21 basic stiffness matrices according to the element size of the target finite element model and store the 21 basic stiffness matrices in the preset GPU global memory. The filtering module 904 is used to perform an optimization iteration process based on the topological optimization model, filter the density design variables and angle design variables, and obtain the filtered density design variables and filtered angle design variables. The second calculation module 905 is used to calculate the rotated constitutive relation matrix of each element according to the property parameters of the preset fiber-reinforced composite material, the filtered density design variables, and the filtered angle design variables. The third calculation module 906 is used to solve the template stiffness matrix of each node at each grid level according to the template stiffness matrix solution format and calculate it in a parallel manner between elements, and store the template stiffness matrices of all other levels except the densest grid level in the preset GPU global memory. Among them, when each element stiffness matrix is required during the calculation process, each element stiffness matrix is calculated according to the 21 basic stiffness matrices and the rotated constitutive relation matrix of each element. The fourth calculation module 907 is used to calculate the preprocessed smoothed residual of the densest grid level using the "V" cycle algorithm of geometric multigrid and process the preprocessed smoothed residual of the densest grid level based on the preconditioned conjugate gradient method to obtain the current displacements of each node at the densest grid level. The solution and analysis module 908 is used to solve the objective function according to the current design variables and the current displacements of each node at the densest grid level, and perform sensitivity analysis on the objective function value and the constraint conditions to obtain the sensitivity value. The iterative optimization module 909 is used to update the current design variables according to the objective function value and the sensitivity value to obtain a new optimized design, determine whether the current optimization iteration times reach the preset threshold. If so, end the iteration and output the new optimized design. Otherwise, continue the iteration.

[0169] It should be noted that the foregoing explanation of the embodiments of the topological optimization method for GPU-accelerated fiber-reinforced composite material structures also applies to the GPU-accelerated topological optimization device for fiber-reinforced composite material structures in this embodiment, and will not be elaborated here.

[0170] The topological optimization device for GPU-accelerated fiber-reinforced composite material structures proposed according to the embodiments of the present invention, on the one hand, greatly reduces the consumption of GPU global memory for storing the stiffness matrix of fiber-reinforced composite material units. On the other hand, it dynamically calculates the stiffness matrix of each unit by multiplying and accumulating the elements in the elastic coefficient matrix with the basic stiffness matrix, and adopts a node-by-node storage method for each grid level and a displacement solution format without assembly, which greatly improves the solution efficiency of large-scale finite element linear equations and can be used in fields such as three-dimensional large-scale structural static analysis and structural topological optimization.

[0171] Figure 10 FIG. is a schematic structural diagram of an electronic device provided by an embodiment of the present invention. The electronic device may include:

[0172] A memory 1001, a processor 1002, and a computer program stored on the memory 1001 and executable on the processor 1002.

[0173] When the processor 1002 executes the program, it implements the topological optimization method for GPU-accelerated fiber-reinforced composite material structures provided in the foregoing embodiments.

[0174] Furthermore, the electronic device further includes:

[0175] A communication interface 1003 for communication between the memory 1001 and the processor 1002.

[0176] The memory 1001 is used to store a computer program executable on the processor 1002.

[0177] The memory 1001 may include a high-speed RAM memory and may also include a non-volatile memory, such as at least one disk memory.

[0178] If the memory 1001, the processor 1002, and the communication interface 1003 are implemented independently, the communication interface 1003, the memory 1001, and the processor 1002 can be interconnected via a bus and communicate with each other. The bus can be an Industry Standard Architecture (ISA) bus, a Peripheral Component Interconnect (PCI) bus, an Extended Industry Standard Architecture (EISA) bus, or the like. The bus can be divided into an address bus, a data bus, a control bus, etc. For ease of representation, Figure 10 it is represented by only one thick line in Figure 10 , but this does not mean that there is only one bus or one type of bus.

[0179] Optionally, in a specific implementation, if the memory 1001, the processor 1002, and the communication interface 1003 are integrated on a single chip, the memory 1001, the processor 1002, and the communication interface 1003 can communicate with each other through an internal interface.

[0180] The processor 1002 may be a Central Processing Unit (CPU), or an Application Specific Integrated Circuit (ASIC), or one or more integrated circuits configured to implement the embodiments of the present invention.

[0181] The embodiments of the present invention also provide a computer-readable storage medium, on which a computer program is stored, and when the program is executed by a processor, the above-mentioned topological optimization method of a GPU-accelerated fiber-reinforced composite material structure is implemented.

[0182] In the description of this specification, the descriptions with reference to the terms "one embodiment", "some embodiments", "example", "specific example", or "some examples", etc. mean that the specific features, structures, materials, or characteristics described in connection with the embodiment or example are included in at least one embodiment or example of the present invention. In this specification, the schematic representations of the above terms do not necessarily refer to the same embodiment or example. Moreover, the specific features, structures, materials, or characteristics described can be combined in any one or N embodiments or examples in a suitable manner. In addition, without contradiction, those skilled in the art can combine and combine the different embodiments or examples described in this specification and the features of different embodiments or examples.

[0183] In addition, the terms "first" and "second" are used for descriptive purposes only and should not be construed as indicating or implying relative importance or implicitly specifying the quantity of the technical features indicated. Thus, features defined with "first" and "second" may explicitly or implicitly include at least one such feature. In the description of the present invention, the meaning of "N" is at least two, such as two, three, etc., unless otherwise specifically defined.

[0184] Any process or method description represented in a flowchart or otherwise described herein may be understood to represent a module, segment, or portion of code including one or N executable instructions for implementing a customized logical function or process, and the scope of the preferred embodiments of the present invention includes additional implementations, where the functions may be executed in a substantially simultaneous manner or in a reverse order according to the functions involved, rather than in the order shown or discussed, which should be understood by those skilled in the art to which the embodiments of the present invention pertain.

[0185] The logic and / or steps represented in a flowchart or otherwise described herein, for example, may be considered as a sequenced list of executable instructions for implementing a logical function, and may be specifically implemented in any computer-readable medium for use by or in connection with an instruction execution system, apparatus, or device, such as a computer-based system, a system including a processor, or other systems that can fetch and execute instructions from the instruction execution system, apparatus, or device. For the purposes of this specification, a "computer-readable medium" may be any device that can contain, store, communicate, propagate, or transport a program for use by or in connection with an instruction execution system, apparatus, or device. More specific examples (a non-exhaustive list) of the computer-readable medium include the following: an electrical connection portion with one or N wirings (electronic device), a portable computer diskette (magnetic device), a random access memory (RAM), a read-only memory (ROM), an erasable programmable read-only memory (EPROM or flash memory), an optical fiber device, and a portable compact disc read-only memory (CDROM). Additionally, the computer-readable medium may even be paper or other suitable medium on which the program can be printed, as the program can be obtained electronically by optically scanning the paper or other medium, followed by editing, interpretation, or other appropriate processing as necessary, and then stored in a computer memory.

[0186] It should be understood that various parts of the present invention can be implemented by hardware, software, firmware or a combination thereof. In the above embodiments, the N steps or methods can be implemented by software or firmware stored in a memory and executed by a suitable instruction execution system. If implemented in hardware, as in another embodiment, any one or a combination of the following techniques well known in the art can be used: discrete logic circuits having logic gate circuits for implementing logical functions on data signals, application specific integrated circuits having appropriate combinational logic gate circuits, programmable gate arrays (PGAs), field programmable gate arrays (FPGAs), etc.

[0187] Those of ordinary skill in the art can understand that all or part of the steps carried by the method of the above embodiments can be completed by instructing relevant hardware through a program, and the program can be stored in a computer-readable storage medium. When the program is executed, it includes one or a combination of the steps of the method embodiments.

[0188] In addition, in each embodiment of the present invention, the functional units can be integrated into a processing module, or each unit can exist physically alone, or two or more units can be integrated into one module. The above integrated module can be implemented in the form of hardware or in the form of a software functional module. When the above integrated module is implemented in the form of a software functional module and sold or used as an independent product, it can also be stored in a computer-readable storage medium.

[0189] The above-mentioned storage medium can be a read-only memory, a magnetic disk or an optical disc, etc. Although the embodiments of the present invention have been shown and described above, it can be understood that the above embodiments are exemplary and should not be construed as limiting the present invention. Those of ordinary skill in the art can make changes, modifications, substitutions and variations to the above embodiments within the scope of the present invention.

Claims

1. A topological optimization method for a GPU-accelerated fiber-reinforced composite material structure, characterized in that, It includes the following steps: Obtain the property parameters of the preset fiber-reinforced composite material and the element size of the target finite element model; Establish a topology optimization model, initialize the design variables of the topology optimization model, and define the constraint conditions of the topology optimization model, where the design variables include density design variables and angle design variables; Calculate 21 basic stiffness matrices according to the element size of the target finite element model, and store the 21 basic stiffness matrices in the preset GPU global memory; Based on the topology optimization model, perform an optimization iteration process, filter the density design variables and the angle design variables to obtain the filtered density design variables and the filtered angle design variables; Calculate the rotated constitutive relation matrix of each element according to the property parameters of the preset fiber-reinforced composite material, the filtered density design variables, and the filtered angle design variables; Adopt the parallel mode between elements, calculate the template stiffness matrix of each node at each grid level according to the template stiffness matrix solution format and the template stiffness matrix transfer relationship between each grid level, and store the template stiffness matrices of all other levels except the densest grid level in the preset GPU global memory. When the stiffness matrix of each element is required in the calculation process, calculate the stiffness matrix of each element according to the 21 basic stiffness matrices and the rotated constitutive relation matrix of each element; Adopt the "V" cycle algorithm of geometric multigrid to calculate the preprocessed smoothed residual of the densest grid level, and process the preprocessed smoothed residual of the densest grid level based on the preconditioned conjugate gradient method to obtain the current displacements of each node at the densest grid level; Solve the objective function according to the current design variables and the current displacements of each node at the densest grid level, and perform sensitivity analysis on the objective function value and the constraint conditions to obtain the sensitivity value; Update the current design variables according to the objective function value and the sensitivity value to obtain a new optimized design, and determine whether the current optimization iteration times reach the preset threshold. If so, end the iteration and output the new optimized design; otherwise, continue the iteration.

2. The topological optimization method of the GPU-accelerated fiber-reinforced composite material structure according to claim 1, wherein The property parameters of the preset fiber-reinforced composite material include at least one of longitudinal elastic modulus, transverse elastic modulus, Poisson's ratio, or in-plane shear modulus.

3. The topology optimization method of the GPU-accelerated fiber-reinforced composite material structure according to claim 1, characterized in that The 21 basic stiffness matrices are: Among them, is the basic stiffness matrix, and n is n = i(i - 1) / 2 + j, B e is the shape function matrix composed of 8 sub - matrices, D ij is a 6×6 matrix with the element in the i - th row and j - th column replaced by 1 and other elements replaced by 0, Ω e is the integration region represented by each element, is D ji is a 6×6 matrix with the element in the j - th row and i - th column replaced by 1 and other elements replaced by 0.

4. The topology optimization method of the GPU-accelerated fiber-reinforced composite material structure according to claim 1, characterized in that, The filtering of the density design variables and the angle design variables to obtain the filtered density design variables and the filtered angle design variables includes: Filter the density design variables by using a preset filtering method to obtain the filtered density design variables, where the preset filtering method is: Among them, is the filtered density design variable, ω is the weight function, ω(x k ) is ω(x k ) = max(R - ||x k - x e ||, 0), Φ e is the set of all elements whose distance from the center point to the center point of the e-th element is not greater than the filtering radius, V k is the volume of the k-th element whose distance from the center point to the center point of the e-th element is not greater than the filtering radius, ρ k is the density design variable of the k-th element whose distance from the center point to the center point of the e-th element is not greater than the filtering radius, N e is the number of elements; Design the angular design variable α in the spherical coordinate system according to the trigonometric function e and θ e and convert them into vectors in the Cartesian coordinate system Filter the vectors in the Cartesian coordinate system according to the preset fiber direction filtering radii respectively to obtain the filtered direction vectors Normalize the filtered direction vector and obtain the filtered angular design variable according to the inverse trigonometric function transformation of the normalized direction vector ​ 5. The topological optimization method of the GPU-accelerated fiber-reinforced composite material structure according to claim 1, characterized in that The template stiffness matrix transfer relationship between each grid level is: where, ω j = (2 - |j1|)(2 - |j2|)(2 - |j3|) / 8, ω k = (2 - |k1|)(2 - |k2|)(2 - |k3|) / 8, both are the weights for the transfer of the template stiffness matrix between different mesh levels, is the template stiffness matrix of each node on a coarser mesh level, is the template stiffness matrix of each node on a denser mesh level.

6. The topological optimization method of the GPU-accelerated fiber-reinforced composite material structure according to claim 1, wherein The adoption of the "V" cycle algorithm of geometric multigrid to calculate the preprocessed smoothed residual of the densest grid level and the processing of the preprocessed smoothed residual of the densest grid level based on the preconditioned conjugate gradient method to obtain the current displacements of each node at the densest grid level includes: Based on the parallel manner between elements, calculate the initial residual of the densest grid level according to the stiffness matrix of the densest grid level template, the preset initial displacement solution, and the preset load of each node. Use the restriction method to transfer the initial residual of the densest grid level to the coarsest grid level to obtain the restricted residual of the coarsest grid level, and use the direct method to solve the restricted displacements of each node of the coarsest grid level according to the global stiffness matrix of the coarsest grid level and the restricted residual of the coarsest grid level. Transfer the restricted displacements of each node of the coarsest grid level to the densest grid level through the prolongation method to obtain the prolonged displacements of each node of the densest grid level, and calculate the prolonged residual of the densest grid level according to the stiffness matrix of the densest grid level template and the prolonged displacements of each node of the densest grid level. Use the post-smoothing processing method to process the prolonged residual of the densest grid level to obtain the smoothed residual after preprocessing of the densest grid level. Use the preconditioned conjugate gradient method to obtain the current displacements of each node of the densest grid level according to the smoothed residual of the densest grid level.

7. The topology optimization method of the GPU-accelerated fiber-reinforced composite material structure according to claim 1, wherein The updating of the current design variables according to the objective function value and the sensitivity value to obtain a new optimal design includes: Based on the optimization criterion method, update the density design variables among the current design variables according to the objective function value and the sensitivity value to obtain new density design variables. Based on the moving asymptotes method, update the angle design variables among the current design variables according to the objective function value and the sensitivity value to obtain new angle design variables.

8. A topological optimization device for a GPU-accelerated fiber-reinforced composite material structure, characterized in that, Includes: An acquisition module for acquiring the property parameters of the preset fiber-reinforced composite material and the element size of the target finite element model. An initialization module for establishing a topology optimization model, initializing the design variables of the topology optimization model, and defining the constraint conditions of the topology optimization model, where the design variables include density design variables and angle design variables. A first calculation module for calculating 21 basic stiffness matrices according to the element size of the target finite element model and storing the 21 basic stiffness matrices in the preset GPU global memory. A filtering module for performing an optimization iteration process based on the topology optimization model, filtering the density design variables and the angle design variables to obtain the filtered density design variables and the filtered angle design variables. A second calculation module for calculating the constitutive relation matrix after rotation of each element according to the property parameters of the preset fiber-reinforced composite material, the filtered density design variables, and the filtered angle design variables. A third calculation module for using the parallel manner between elements to solve the stiffness matrix of each node of each grid level according to the solution format of the template stiffness matrix, and storing the stiffness matrices of all other levels except the densest grid level in the preset GPU global memory, where when the stiffness matrix of each element is required in the calculation process, calculate the stiffness matrix of each element according to the 21 basic stiffness matrices and the constitutive relation matrix after rotation of each element. The fourth calculation module is used to calculate the smoothed residual after preprocessing of the densest grid level by using the "V" cycle algorithm of geometric multigrid, and process the smoothed residual after preprocessing of the densest grid level based on the preconditioned conjugate gradient method to obtain the current displacements of each node at the densest grid level; The solution and analysis module is used to solve the objective function according to the current design variables and the current displacements of each node at the densest grid level, and perform sensitivity analysis on the objective function value and the constraint conditions to obtain sensitivity values; The iterative optimization module is used to update the current design variables according to the objective function value and the sensitivity values to obtain a new optimized design, determine whether the current optimization iteration times reach a preset threshold, if so, end the iteration and output the new optimized design, otherwise continue the iteration.

9. An electronic device, characterized in that, Comprising: A memory, a processor, and a computer program stored on the memory and executable on the processor, where the processor executes the program to implement the GPU-accelerated topology optimization method for fiber-reinforced composite material structures as described in any one of claims 1-7.

10. A computer-readable storage medium having a computer program stored thereon, characterized in that, The program is executed by the processor to be used for implementing the GPU-accelerated topology optimization method for fiber-reinforced composite material structures as described in any one of claims 1-7.

Images