FEM-Validated Optimal Design of Laminate Process Parameters Based on Improved Genetic Algorithm

: In tape placement process, the laying angle and laying sequence of laminates have proven their signiﬁcant effects on the mechanical properties of carbon ﬁbre reinforced composite material, speciﬁcally, laminates. In order to optimise these process parameters, an optimisation algorithm is developed based on the principles of genetic algorithms for improving the precision of traditional genetic algorithms and resolving the premature phenomenon in the optimisation process. Taking multi-layer symmetrically laid carbon ﬁbre laminates as the research object, this algorithm adopts binary coding to conduct the optimisation of process parameters and mechanical analysis with the laying angle as the design variable and the strength ratio R as the response variable. A case study was conducted and its results were validated by the ﬁnite element analyses. The results show that the stresses before and after optimisation are 116.0 MPa and 100.9 MPa, respectively, with a decrease of strength ratio by 13.02%. The results comparison indicates that, in the iterative process, the search range is reduced by determining the code and location of important genes, thereby reducing the computational workload by 21.03% in terms of time consumed. Through multiple calculations, it validates that “gene mutation” is an indispensable part of the genetic algorithm in the iterative process.


Introduction
Carbon fibre reinforced composite materials (CFRP) are widely used in the fields of aviation equipment structures [1], fuselages and wings manufacturing due to their high corrosion resistance, high strength and lightweight [2,3]. As a typical representative of the practical application of composite materials, laminates have a direct impact on the mechanical properties of CFRPs, which depend highly on their laying angle and laying sequence. Hence, many investigations have been conducted in-depth on exploring the influence of these two variables and searching for their optimal combinations.
A genetic algorithm (GA) is a random search method abstracted by simulating gene selection and genetic mechanisms in the natural environment. It can not only solve linear problems but also can be used to solve nonlinear problems; and genetic algorithm is sensitive to the problem of discretization, so it has attracted the attention of scholars [4], and it has been applied to many occasions. At present, scholars worldwide favor the use of genetic algorithms to optimise the design of the mechanical properties of laminates, and certain research results have been obtained. Wang et al. [5] applied a multi-island genetic algorithm to optimise the damage of different loads to the laminate during low-velocity impact and optimised the laminate sequence. Strain energy was taken as the optimisation objective and the laying angle as the design variable. Their results show that the impacted area was reduced by 42.1%, and the residual compression strength increased by 10.79%. It 2 of 24 is notable that the laying angle was optimised as a continuous variable rather than discrete as in the other investigations mentioned above. Chen et al. [6] aimed at optimising the fibre laying angle to achieve the minimum stress of the laminate, and used a genetic algorithm to obtain the best laying angle and minimum stress of the laminate. Their method realised the optimisation of the laying angle of two-layer laminates and three-layer laminates by taking the design variable as a continuous variable.
In order to improve the reliability, some optimisation strategies were proposed by Xiu and Cui [7], which took the buckling load of laminate as the optimisation objective, and the laying angle of laminate as constraint condition, and used the neural network model and genetic algorithm to optimise the ratio of laying angle to the buckling load of the laminate. The optimisation method proposed in their investigation is capable of producing a better laying sequence with the same material weight under the condition that the laying angles are given and the number of layers is a constant; nevertheless, the critical process parameter, laying angle, cannot be optimised by this method and is just taken as a given constant. Feng et al. [8] took the strength of laminate as the optimisation objective, and the laying angle of laminate as the constraint condition. In order to obtain a reasonable laying angle and ply number, Jin et al. [9] proposed a two-step optimisation strategy and developed a three-level optimisation model with embedded genetic algorithm and Tsai-Wu failure criterion. The concept of sublayer was introduced into the optimisation method for the optimisation of the location of layers, geometric dimensions, number of layers, and laying sequence. The authors claimed that the introduction of this concept improved the manufacturability of laminates, however, the laying angle was not optimised in this work. By embedding the Tsai-Hill failure criterion, Park et al. [10] used genetic algorithms to apply different loads to the symmetrical composite laminates under different boundary conditions to perform mechanical optimisation design. Tournament selection and uniform crossover methods were adopted in the optimisation process. The laying angle was taken as the design variable as a continuous variable, and a low failure index was obtained. The results showed that the modified genetic algorithm was capable of identifying a global optimal solution for the parameter's optimisation of laid laminates.
Furthermore, new encoding approaches and operators have been developed for genetic algorithms to accelerate the convergence speed and computational accuracy. A genetic algorithm was used to calculate the fibre laying angle and validated that the calculation results can effectively improve the strength of the laminate. This work improved the genetic algorithm by modifying the encoding approach and the integer encoding was embedded in the iteration process for the laminate parameters optimisation. The advantage of this modification is that integer encoding has the potential to avoid unnecessary iterations of floating number, thereby improving the computational efficiency theoretically. The laying angle in this optimisation was specified to 0 • , ±45 • , and 90 • , leading to weakening the designability of laminates to some extent. By designing a new mutation operator and crossover operator, Wang et al. [11] adopted a symbol encoding method for the development of a self-adaptive genetic algorithm with the bending stiffness parameter as the fitness function. They deeply explored the influence of the adaptively changing crossover operator and the mutation operator on the genetic algorithm and further optimised the laying order. Compared with the traditional genetic algorithm, their modification enhanced the reliability, convergence speed, and computational efficiency of the optimisation process. Another self-adaptive genetic algorithm was proposed by Yang et al. [12] through the automatic change of the mutation rate and crossover rate according to fitness value. This algorithm took the maximum deformation of laminate as optimisation objective, the laminate angle as constraint condition, and called ABAQUS to analyse the stress of the laminate. The results show that this algorithm is easier to get convergence than the traditional genetic algorithm. Hwang et al. [13] designed an elite comparison operator for the genetic algorithm to accelerate the iteration process for the optimisation of laminates. This elite comparison operator was used to identify the differences in the design variables from two different generations, thereby maintaining the similar ones and changing the different ones to find their possible values.
Two optimisation problems were solved by this modified genetic algorithm and the results indicate that the designed operator has the potential to speed up the convergence of the optimisation process.
The above investigations mainly focus on the change of optimisation parameters, design of encoding methods, and improvement of operators in the application of traditional genetic algorithms. However, in order to improve the precision of traditional genetic algorithms and eliminating the premature phenomenon in the optimisation process, this article improves the genetic algorithm (improved genetic algorithm, referred as IGA below) on the basis of the principles of traditional genetic algorithms through the application of binary coding to explore important gene codes and their locations. In this IGA, the laying angle is set as the constraint condition, and the minimum stress as the optimisation objective, which customised this algorithm for the optimisation of process parameters of laminates tape placement. Based on the Tsai-Wu failure criterion [9], the laying angle and laying sequence were optimised; and the obtained stress was compared with that of the traditional genetic algorithm to validate the IGA. In addition, the optimal solutions were further validated by a finite element analysis of a laminate with the same geometry, loading, and boundary conditions.

On-and off-Axis Stress Conversion of Laminates
Generally speaking, a multi-layer laminate can be formed by stacking multiple monolayer boards at different laying angles. According to the designability of composite materials and the theoretical basis of elastic mechanics, the on-axis direction of the monolayer boards can be laid at different angles and in different orders according to the loading conditions to meet the requirements of strength and stiffness design of the structural parts under the load conditions.
Research for the classical laminate theory (CLT), the establishment and derivation of formulae are mainly based on the following assumptions:

•
The layers of the laminate need to be firmly bonded so that the deformation between the layers is consistent, and relative slippage between the layers can be avoided; • Laminates are thin plates with an invariable thickness, and their strength in the thickness direction is relatively small to negligible in respect to that in other directions; • The bending deformation of the laminate needs to be in a small deflection range, and the straight line perpendicular to the midplane must be maintained before and after the bending deformation while the length of the straight line remains unchanged (viz. γ xz = γ yz = 0; ε z = 0); • It is considered that each monolayer of the laminate is in a plane stress state, that is, the theory of plane stress is suitable for the analysis of the laminate structure of thin planes, curved surfaces, or shells.
For simplification of the formula derivation process, the median plane of the laminate is taken as the reference plane, viz. the plane perpendicular to the z-direction at the middle point of the thickness of the laminate. According to the above assumptions, let the stress function x and y be ϕ, then the corresponding stress components σ x , σ y , τ xy , are given in Equation (1) as below: Since the classic laminate theory only considers plane stress and ignores the existence of interlayer stress in the iterative process, the stiffness coefficient of each layer of the laminate is different so that the internal force of the laminate can only be integrated layer by layer. Assuming that the stress of the kth layer of the laminate is σ k , the internal force of the laminate can be expressed as: The stress-strain relationship of the kth layer of the laminate is as follows: Substituting Equation (2) into Equation (1), and through integration, the stress-strain relationship of the laminate can be obtained as follows: where, i, j = 1, 2, 6; A ij is the in-plane stiffness coefficient, B ij is the coupling stiffness coefficient, and D ij is the bending stiffness coefficient. The flexibility matrix of the laminate is where, As the presence of the matrix [B] indicated, the laminate is prone to the coupling of bending deformation and tension and compression deformation when subjected to force; in addition, in-plane deformation is relatively easy to occur when subject to bending stress, and further leads to the warpage of the laminate after curing (that is, B ij = 0, which means it is hard for the calculation of laminates). For this reason, the composite material laminates generally used in engineering should be symmetrically laid (B ij = 0) to avoid the calculation difficulties caused by the coupling stiffness [14]. Note that symmetry is assumed in this study to make sure the laminate only subject to tensile strain but no flexure [14], though it is difficult to always keep symmetry exactly in industry. Therefore, the flexibility matrix S of the laminate at this time can be simplified as: At the same time, through the stress conversion matrix [T σ ], the off-axis stress of each layer can be converted into the on-axis stress. The specific conversion formula is shown in Equation (8): where, m = cos θ, n = sin θ, and θ is the laying angle of each layer of fibres.

Failure Criterion for Laminates
After comprehensively considering the advantages and disadvantages of multiple strength criteria, Tsai and Wu [15] proposed a new strength criterion based on the tensor form. That is, the matrix breaks based on Tsai-Wu criterion, as below: where, σ 1 is the longitudinal stress, σ 2 is the transverse stress, σ 6 is the shear stress, X T is the longitudinal tensile strength, X C is the longitudinal compressive strength, Y T is the transverse tensile strength, Y C is the transverse compressive strength, S is the shear strength; F 1 , F 11 , F 2 , F 22 , and F 66 are strength parameters as shown below in Equation (10).
The value of F 12 here has little effect on the calculation of the strength criterion, and this is validated by the investigation cited here [16]. Then let F 12 = 0, the calculation result of Equation (9) and that of the Hoffman strength criterion are relatively close with an error of 10%. Hence, Equation (9) is improved as below in Equation (11).
The value of F 16 , F 26 here is obtained by a semi-empirical formula [17], and the value of F 12 is the value when the Tsai-Wu strength criterion and the Hoffman criterion are the same. See Equation (12) for details.
Note that the value of F 12 here is not −1/2 √ F 11 F 22 because: • The absolute value of Equation (12) is relatively not very large, and it is a negative value, which is reasonable; • When part of the curves in Figure 1 [15] is under the 0 MPa line, there is a situation where the value of σ 1 is greater than 2X C ; and this situation does not exist in reality because when greater than, fracture has already appeared for composite materials.

Failure Criterion for Laminates
After comprehensively considering the advantages and disadvantages of multiple strength criteria, Tsai and Wu [15] proposed a new strength criterion based on the tensor form. That is, the matrix breaks based on Tsai-Wu criterion, as below: where, is the longitudinal stress, is the transverse stress, is the shear stress, is the longitudinal tensile strength, is the longitudinal compressive strength, is the transverse tensile strength, is the transverse compressive strength, S is the shear strength; , , , , and are strength parameters as shown below in Equation (10).
The value of here has little effect on the calculation of the strength criterion, and this is validated by the investigation cited here [16]. Then let = 0, the calculation result of Equation (9) and that of the Hoffman strength criterion are relatively close with an error of 10%. Hence, Equation (9) is improved as below in Equation (11).
The value of , here is obtained by a semi-empirical formula [17], and the value of is the value when the Tsai-Wu strength criterion and the Hoffman criterion are the same. See Equation (12) for details.
Note that the value of here is not −1/2 because: • The absolute value of Equation (12) is relatively not very large, and it is a negative value, which is reasonable; • When part of the curves in Figure 1 [15] is under the 0 MPa line, there is a situation where the value of is greater than 2 ; and this situation does not exist in reality because when greater than, fracture has already appeared for composite materials. From Equation (11), when the calculated value on the left side of the equation equates 1, it indicates that the material has just reached the critical failure state; when the left side is greater than 1, it indicates that damage has begun; when the left side is less than 1, it indicates that the material is in a safe state, and there is still a certain safety margin.
In order to obtain better calculation speed and expected values in genetic algorithms, a suitable. fitness value is generally specified. Here, the strength ratio is applied to the From Equation (11), when the calculated value on the left side of the equation equates 1, it indicates that the material has just reached the critical failure state; when the left side is greater than 1, it indicates that damage has begun; when the left side is less than 1, it indicates that the material is in a safe state, and there is still a certain safety margin.
In order to obtain better calculation speed and expected values in genetic algorithms, a suitable. fitness value is generally specified. Here, the strength ratio is applied to the laminate, namely where R is the ratio of a certain component of the ultimate stress to the corresponding component of loading stress, and i = 1, 2, 6. The meaning of "correspondence" here [18] assumes that the loading is linear, that is, the stress components increase simultaneously in a certain proportion, and this situation is in line with the actual application in industrial production. Transforming Equation (13) to get Rσ = σ i , substitute this into Equation (11) and after simplification to produce Equation (14) below.
This is a quadratic equation in one variable. It is easy to get Equation (15) as below. where: Equation (15) is called the strength ratio equation. From the definition of the strength ratio, it can be seen that when R > 1, the material has not failed; when the applied stress is increased to (R − 1) times, the material fails. R cannot be less than 1 because this is mathematically meaningless; on contrast, R < 1 in engineering applications still has a certain reference value, that is, it indicates that the applied stress must be reduced or the size of the structure must be increased to ensure that the structure is not damaged.

Genetic Algorithm and Its Improvement
Genetic algorithm is based on the natural law of survival of the fittest by simulating the theory of genetics and biological evolution proposed by Darwin, thereby retaining the excellent genes and individuals [19].
The core content of genetic algorithm mainly includes four parts, namely selection, crossover, mutation, and fitness evaluation. The optimisation process is shown in Figure 2.
where R is the ratio of a certain component of the ultimate stress to the corresponding component of loading stress, and = 1, 2, 6. The meaning of "correspondence" here [18] assumes that the loading is linear, that is, the stress components increase simultaneously in a certain proportion, and this situation is in line with the actual application in industrial production. Transforming Equation (13) to get = , substitute this into Equation (11) and after simplification to produce Equation (14) below.
This is a quadratic equation in one variable. It is easy to get Equation (15) as below. where: Equation (15) is called the strength ratio equation. From the definition of the strength ratio, it can be seen that when R > 1, the material has not failed; when the applied stress is increased to (R − 1) times, the material fails. R cannot be less than 1 because this is mathematically meaningless; on contrast, R < 1 in engineering applications still has a certain reference value, that is, it indicates that the applied stress must be reduced or the size of the structure must be increased to ensure that the structure is not damaged.

Genetic Algorithm and Its Improvement
Genetic algorithm is based on the natural law of survival of the fittest by simulating the theory of genetics and biological evolution proposed by Darwin, thereby retaining the excellent genes and individuals [19].
The core content of genetic algorithm mainly includes four parts, namely selection, crossover, mutation, and fitness evaluation. The optimisation process is shown in Figure  2.  At present, the improvement of genetic algorithm mainly focuses on self-adaptive genetic operators, that is, generating a certain number of individuals, improving the crossover rate and genetic probability between individuals during genetic operations, and using average fitness to improve calculation efficiency. The method proposed in this paper is to start from the population; when generating the population, only one individual is generated. Then the genes of this individual mutate later in the evolution process, and the gene mutation means gene transfer from one vertex to another according to the definition of the hyperplane. The calculation is performed by gene mutation; the better genes and positions are selected through fitness, and negative feedback adjustment is used to generate individuals with important genes, and then identify the corresponding laying angle through decoding. The basic mathematical theory mentioned above is: For the full set space composed of individuals whose gene length is l [20], that is, the individual space is given in Equation (17) where we regard this as the spatial vertex set with the dimension l. When l = 3, it becomes a three-dimensional cube. As Figure 3 indicated, for three-dimensional individual space, the genes at each vertex can be connected by 12 straight lines. It can also be seen that in a straight line, we only need to change the value of one of the genes, then the gene can be converted into the gene at the other end of the same line; in any plane, as long as the values of two genes are changed, the genes at the four vertices can be switched to each other. Therefore, in a three-dimensional individual space, the straight line and the plane are called the "hyperplane", and the interconnection between them constitutes the "subspace" in the space. At present, the improvement of genetic algorithm mainly focuses on self-adaptive genetic operators, that is, generating a certain number of individuals, improving the cross over rate and genetic probability between individuals during genetic operations, and us ing average fitness to improve calculation efficiency. The method proposed in this pape is to start from the population; when generating the population, only one individual is generated. Then the genes of this individual mutate later in the evolution process, and the gene mutation means gene transfer from one vertex to another according to the definition of the hyperplane. The calculation is performed by gene mutation; the better genes and positions are selected through fitness, and negative feedback adjustment is used to gener ate individuals with important genes, and then identify the corresponding laying angle through decoding. The basic mathematical theory mentioned above is: For the full set space composed of individuals whose gene length is l [20], that is, the individual space is given in Equation (17) = 0,1 (17 where we regard this as the spatial vertex set with the dimension l. When l = 3, it becomes a three-dimensional cube. As Figure 3 indicated, for three-dimensional individual space, the genes at each ver tex can be connected by 12 straight lines. It can also be seen that in a straight line, we only need to change the value of one of the genes, then the gene can be converted into the gene at the other end of the same line; in any plane, as long as the values of two genes are changed, the genes at the four vertices can be switched to each other. Therefore, in a three dimensional individual space, the straight line and the plane are called the "hyperplane" and the interconnection between them constitutes the "subspace" in the space. Once the "hyperplane" is determined, it is possible to locate the gene value of the corresponding position after confirming that the gene value of some gene position is valid thereby reducing the search range and further reducing the computational workload. The specific process of the IGA is shown in Figure 4. Therefore, according to the above definition of "hyperplane", for the i-th gene, if for any x k = 0, 1(k = i) holds, or for any x k = 0, 1(k = i) holds, the gene i is taken as an important gene. When Equation (18) is established, important gene 1 is the better choice; while when Equation (19) is established, important gene 0 is the better choice.
Once the "hyperplane" is determined, it is possible to locate the gene value of the corresponding position after confirming that the gene value of some gene position is valid, thereby reducing the search range and further reducing the computational workload. The specific process of the IGA is shown in Figure 4.

Ply Design Requirements for Laminates
For the ply design of laminates, the coupling effect caused by stretching and bending during the curing process should be avoided since it causes resin cracking, warping deformation and even delamination for laminates. Therefore, a design generally complies with the following process design rules [21].

•
Standard laid layers should be adopted, which is supposed to include four angle values of 0°, ±45°, and 90°, respectively; • Adjacent 4 monolayers cannot be formed at the same laying angle so as to prevent

Ply Design Requirements for Laminates
For the ply design of laminates, the coupling effect caused by stretching and bending during the curing process should be avoided since it causes resin cracking, warping deformation and even delamination for laminates. Therefore, a design generally complies with the following process design rules [21].

•
Standard laid layers should be adopted, which is supposed to include four angle values of 0 • , ±45 • , and 90 • , respectively; • Adjacent 4 monolayers cannot be formed at the same laying angle so as to prevent the substrate from cracking, and simultaneously a full 90 • ply should be avoided; • The proportion of 0 • laid layers is between 20% and 40%, ±45 • layers between 40% and 60%, and 90 • layers must be between 10% and 30%.

Genetic Algorithm Coding
In the process of searching optimal solutions, the genetic algorithm generally uses a finite length string code to define the solution to be solved, such as a binary code. The layup of a laminate is a matter of discrete variables, so the choice of variables affects the length of the code string to some extent. Since the other encoding schemes, such as integer encoding and floating-point number encoding, are not unique in the subsequent mutations, they are not conducive to genetic operations such as crossover and mutation. Therefore, traditional binary encoding is selected hereof. The concrete corresponding method is: Since the initial population generates randomly, there must be a certain number of initial individuals that do not meet the requirements of the laminate layup design; these individuals must be eliminated, and the others that meet the requirements should remain. The subsequent operations conduct all the above operations for each subsequent generation, thereby ensuring that the population individuals are suitable. The operation is like a "human intervention" in nature.

Fitness Function
In the genetic algorithm, fitness is used to evaluate the goodness of individuals adapting to the "environment", thereby further weighing the degree to which individuals can reach or approach the optimal solution in the process of genetic evolution [22]. Therefore, the fitness function selected for this calculation is where, n is the total number of laminates, R(i) is the strength ratio of each layer of fibre, n 1 is the number of adjacent four layers with the same laying angle, n 2 is the number of layers that do not meet the requirements of the layup ratio, 5n 1 + 10n 2 is the penalty factor. Because the initial individuals generated may not all meet the design requirements of the laminate, as mentioned above, a certain high penalty is required for the fitness factor in the searching process. Hence, through reduction of the fitness function value of unsuitable. Individuals, the probability of being selected in selection replication is further reduced. In actual operations, there may be a phenomenon that a high penalty makes R s a negative value, and such a value is not conducive to the calculation of the subsequent results. Therefore, when it is a negative value, set R s = 0.05, thereby forcing it to become a small positive value by "human interference" to facilitate subsequent calculations. According to the equations above, the optimisation process of the strength of laminate by the IGA algorithm is the process of analysing and calculating the best strength ratio of the single layer in the population.

Relevant Parameters Selection
In this paper, T300/5228 epoxy carbon fibre with good specific strength is selected, and the material parameters are shown in Table 1. In an optimisation process, the communication of individuals between populations is not a complete fusion, but on the basis of absorbing each other's excellent genes and simultaneously ensuring that their own excellent genes are preserved. Therefore, individuals need to use different crossover rates P C and mutation rates P m . Furthermore, the values of P C and P m also affect the global and local search capabilities and lead to different optimisation results; empirically, researchers recommend to adopt a large value of P C and small value of P m . Here in this article, the two probability factors are specified as 0.95 and 0.005, respectively. The specific parameter settings related to the genetic algorithm are detailed in Table 2. At the same time, in order to speed up the optimisation process, the evaluation parameters of the adaptive genetic operator are introduced [23], namely where, R max is the maximum fitness value in each generation of individuals; R avg is the average fitness value in each generation of the population; R min is the minimum fitness value in each generation; ε is an infinitesimal with the main purpose to prevent the extreme case of 0 in the denominator. Through Equation (21), the value of P C and P m in the two cases of U > 1 and U < 1 is further derived so as to speed up the calculation speed.

Relevant Parameters Selection
Given that a rectangular composite laminate beam with four sides simply supported is designed, the size of the laminate is, length b = 2l with the value of 1300 mm, width a = 132 mm, thickness h = 2.25 mm, and the bearing pressure perpendicular to the surface q = 0.5 MPa. The design adopts a symmetrical 18-layer layup scheme. The geometric model of the laminate is illustrated in Figure 5.
J. Compos. Sci. 2022, 6, x FOR PEER REVIEW optimisation results; empirically, researchers recommend to adopt a larg and small value of . Here in this article, the two probability factors are sp and 0.005, respectively. The specific parameter settings related to the gen are detailed in Table 2. At the same time, in order to speed up the optimisation process, the rameters of the adaptive genetic operator are introduced [23], namely where, is the maximum fitness value in each generation of individua average fitness value in each generation of the population; is the mi value in each generation; is an infinitesimal with the main purpose to p treme case of 0 in the denominator.
Through Equation (21), the value of and in the two cases of U > further derived so as to speed up the calculation speed.

Relevant Parameters Selection
Given that a rectangular composite laminate beam with four sides sim is designed, the size of the laminate is, length b = 2l with the value of 1300 132 mm, thickness h = 2.25 mm, and the bearing pressure perpendicular to 0.5 MPa. The design adopts a symmetrical 18-layer layup scheme. The ge of the laminate is illustrated in Figure 5.  The schematic diagram of the simply supported laminate beam is shown in Figure 6, as below.   Since σ y is caused by the lateral load q, and q is a fixed value, which means that σ y is a function of y, here we can set (See the detailed derivation process of the equations in this section at the Appendix A.) Then integrating Equation (22) it becomes Substituting Equation (23) into compatibility Equation (26) it becomes where, the balance equation and compatibility equation are shown in Equations (25) and (26), and the stress function ϕ is shown in Equation (1) above If for any x, Equation (23) holds, then all the coefficients of Equation (24) need to be 0, that is Hence, it can set f (y) as According to the principle of Saint-Venant it gives At the same time, the component force according to the shear stress is Q, namely Substituting Equation (28) into boundary condition Equations (29) and (30), and according to the axial component force in the axial direction is 0, the stress expression of σ x , σ y , τ xy can be derived as [24] In order to conduct optimisation, R is first obtained by combining Equations (31) and (15) according to the above requirements, and then it is substituted into Equation (20) to produce the penalty fitness function; and the MATLAB script was written. For speeding up the convergence process, here the adaptive operator of Equation (21) was introduced, then the size of the crossover rate P C was calculated and so was the mutation rate P m between individuals, thereby further obtaining the target value. The compiled calculation program was conducted on a PC with the CPUs (AMD Ryzen 71700 Eight-Core Processor, 2.99 GHz); and the important gene was 0 with 2 important gene positions, which were the 1st and 5th genes, respectively. In order to avoid contingency caused by one or two cycles, the program for finding important genes was run repeatedly and the important genes generated after each cycle were counted and sorted, and then through counting the number of occurrences and position numbers, finally it got the important gene positions as the 1st and 5th genes. Then injecting the obtained important genes into new individuals randomly generated, the strength ratio R was obtained during the evolution process. It needs to run multiple times to prevent local convergence that may occur accidentally and fail to reach the global convergence value.

Results Discussion
The optimisation diagram of final strength ratio R is shown in Figure 7. It suggests that after multiple runs, almost all optimisation curves had fluctuated in the fitness value in the early stage. This is because it is necessary to determine whether the population individual meets the requirements of the layup design after each iteration; if not, eliminate them and introduce new individuals. Therefore, in the continuous genetic evolution process, individuals that did not meet the requirements were eliminated, and new individuals that meet the requirements were added in thereby ensuring that the original design layup criteria were always met.
This value always changed in the optimisation process in the later stage. As indicated by the first, fourth, fifth, and seventh optimisation curves in Figure 7, the best strength ratio appeared at a certain genetic stage, and the genetic generation number was kept short. The "best value" is produced by "gene mutation" and is not the global optimal solution. Simultaneously, the 4th and 7th optimisation curves gradually increase with the number of iterations. Although the optimal strength ratio was found due to the "gene mutation", the final optimisation curve converges locally, which also indicates that the genetic algorithm has the disadvantage of local convergence, and thus, unable to find the global optimal solution smoothly.
reach the global convergence value.

Results Discussion
The optimisation diagram of final strength ratio R is shown in Figure 7. It suggests that after multiple runs, almost all optimisation curves had fluctuated in the fitness value in the early stage. This is because it is necessary to determine whether the population individual meets the requirements of the layup design after each iteration; if not, eliminate them and introduce new individuals. Therefore, in the continuous genetic evolution process, individuals that did not meet the requirements were eliminated, and new individuals that meet the requirements were added in thereby ensuring that the original design layup criteria were always met.  In order to determine whether the curve is reasonable and whether the calculation result meets the initial design requirements, the curve was processed through polynomial data fitting. The seven strength ratios were averaged in the fitting process, and the final fitting result is shown in Figure 8. This value always changed in the optimisation process in the later stage. As indicated by the first, fourth, fifth, and seventh optimisation curves in Figure 7, the best strength ratio appeared at a certain genetic stage, and the genetic generation number was kept short. The "best value" is produced by "gene mutation" and is not the global optimal solution. Simultaneously, the 4th and 7th optimisation curves gradually increase with the number of iterations. Although the optimal strength ratio was found due to the "gene mutation", the final optimisation curve converges locally, which also indicates that the genetic algorithm has the disadvantage of local convergence, and thus, unable to find the global optimal solution smoothly.
In order to determine whether the curve is reasonable and whether the calculation result meets the initial design requirements, the curve was processed through polynomial data fitting. The seven strength ratios were averaged in the fitting process, and the final fitting result is shown in Figure 8. The expression of the red curve in Figure 8 is given below as, where The expression of the red curve in Figure 8 is given below as, According to statistics, Goodness of fitting T 2 is the indicator to evaluate the degree of fit [25]. It is expressed in Equation (34) below.
where,R i is the fitted value; R is the average value of the initial data; R i is the i-th initial data. Putting the average of the strength ratios of the IGA algorithm into Equation (34), it becomes T 2 = 0.943. The value range of T 2 is 0 to 1. The closer to 1, the better the fitting effect is. In order to further determine the reliability of this improved algorithm, the number of iterations was increased up to 550 generations, and the optimisation iteration was performed three times with their results averaged. At the same time, the calculation result of the fitted curve was substituted for comparison, as shown in Figure 9.
J. Compos. Sci. 2022, 6, x FOR PEER REVIEW 14 o Figure 9. The curve of reliability of the IGA. Figure 9 shows that after 550 generations with the maximum difference ratio of 0. the minimum difference ratio of −0.18%, the maximum difference range of 0.208%, a of 0.943 at the same time; hence, it can be considered that the IGA algorithm and fitted curve are reliable.
It can be seen from Tables 3 and 4 and Figure 7 that when the number of iterati evolves to 52 generations at the earliest stage, the strength ratio R reaches the maxim value of 3.825, suggesting a significant improvement with the increase rate of 48.83%, a further, in the iterative process, the laying angle ±45° gradually approaches the outerm layer, which can also effectively improve the impact resistance and stability of the la nate. At the same time, through decoding the optimised layering scheme inversely, it s gests that the number of genes 0 is gradually increasing. This gradual increase also sho that 0 is an important gene, which is consistent with the initial calculation results. Co  Figure 9 shows that after 550 generations with the maximum difference ratio of 0.1%, the minimum difference ratio of −0.18%, the maximum difference range of 0.208%, and T 2 of 0.943 at the same time; hence, it can be considered that the IGA algorithm and the fitted curve are reliable.
It can be seen from Tables 3 and 4 and Figure 7 that when the number of iterations evolves to 52 generations at the earliest stage, the strength ratio R reaches the maximum value of 3.825, suggesting a significant improvement with the increase rate of 48.83%, and further, in the iterative process, the laying angle ±45 • gradually approaches the outermost layer, which can also effectively improve the impact resistance and stability of the laminate. At the same time, through decoding the optimised layering scheme inversely, it suggests that the number of genes 0 is gradually increasing. This gradual increase also shows that 0 is an important gene, which is consistent with the initial calculation results. Compared with the traditional genetic algorithm, in order to ensure a single variable, the population initially formed by the IGA was retained, and simultaneously ensuring that other initial settings remained unchanged; then optimisation iterations were performed. The calculation results are shown in Figure 10 and Table 5.    Figure 10 presents that the strength ratio R obtained by the GA algorithm mainly   [9]" denotes there are nine layers with a laying angle of 0 • . Figure 10 presents that the strength ratio R obtained by the GA algorithm mainly shows an increasing trend at the beginning, and the maximum strength ratio R is 8.819. At the same time, Table 5 shows that the laying angles are all 0 • at this stage, and then the strength ratio R fluctuates. This can be attributed to the "gene mutation" due to mutation during the evolution process. From Table 4, it is known that the "gene mutation" has produced Gene 1; however, as indicated by the layup scheme in Table 4, the laying angle does not meet the original design criteria, which led to the phenomenon of "prematurity".
At the same time, through the comparison between Figures 7 and 10 it results that: • Due to the existence of "gene mutation", the IGA algorithm identifies the optimal strength ratio, but due to the existence of local convergence, it is necessary to rely on "gene mutation" to find the optimal strength ratio. Hence, "gene mutation" cannot be removed in actual operation. At the same time, in order to avoid contingency, it is necessary to run multiple times to ensure the unity of the results, thereby preventing local convergence produced by a single run. The GA algorithm has an "unreasonable" maximum strength ratio in the 149th generation, and this reversely indicates the occurrence of the "premature" phenomenon.

•
The average running time of the IGA algorithm is 217.41 s, while the traditional genetic algorithm is 275.293 s. It is obvious that the IGA algorithm has a certain improvement in its iterative operation efficiency, with an increase of 21.03%.

FEM Validation
In order to validate the proposed IGA, a finite element model of a laminate was developed with a length of 1200 mm, a width of 130 mm, and thickness of 2.25 mm. The laminate has a total of 18 monolayers with a thickness of each layer of 0.125 mm. The Laminate plate was constrained by fixed the x, y, and z directions on all four sides (U1 = U2 = U3 = 0). That is to say, the laminate plate was set as a simply supported beam. Then the plate was meshed with hexagonal shell elements with eight nodes (SC8R) and reduced integration. A uniform load of 0.5 MPa was applied perpendicular to the x-y plane of the laminate plane. Finally, the finite element analysis was conducted through the commercial software ABAQUS and the simulation results before and after optimisation were obtained, as shown in Figure 11.
As indicated by Figure 11a,b, the maximum stress in the fibre direction before optimisation is 116.0 MPa, and the maximum stress after optimisation is 100.9 MPa. The stress value is greatly reduced, and the reduction ratio reaches 13.02%, indicating that the optimisation results meet the initial requirements. The stresses at the other two directions, namely s22 and s12 changed from 8.57 MPa and 8.80 MPa to 5.464 MPa and 11.11 MPa, respectively. These two directions are not the directions mainly subject to applied forces and the corresponding stresses are also relatively small to s11. At the same time, it can be found from the figures that the maximum stress is mainly concentrated on the two long sides of the laminate. This is because the use of four-sided simple support leads to stress concentration. U2 = U3 = 0). That is to say, the laminate plate was set as a simply supported beam. Th the plate was meshed with hexagonal shell elements with eight nodes (SC8R) and reduc integration. A uniform load of 0.5 MPa was applied perpendicular to the x-y plane of laminate plane. Finally, the finite element analysis was conducted through the commerc software ABAQUS and the simulation results before and after optimisation were o tained, as shown in Figure 11. As indicated by Figure 11a,b, the maximum stress in the fibre direction before op misation is 116.0 MPa, and the maximum stress after optimisation is 100.9 MPa. The str Figure 11. Stress contour of the laminate before and after optimisation: (a) Stress contour before optimisation; (b) Stress contour after optimisation.

Conclusions
Through the development of an IGA, the laying angle and laying sequence of laminates were optimised and the optimal results were compared with those of traditional genetic algorithms. The reliability of the optimisation was further validated by comparison with the finite element analysis results and the conclusions are drawn as below.

•
Based on the genetic algorithm, this paper optimised the layering angle of laminate, and adopted the identification of important genes to determine the important genes and their positions and reduced the query range so that the calculation cost in terms of time was reduced by 21.03%.

•
In the process of genetic evolution, the operation "gene mutations" is essential and indispensable, and it proves that the "gene mutation" has the potential to facilitate the identification of global optimal solutions. • Through the development of MATLAB script, it found the important genes were the first and fifth genes. The optimisation results show that the improved strength ratio is 3.825, and the optimal laying angle is sequentially [−45 The stresses before and after the optimisation were 116.0 MPa and 100.9 MPa, respectively, with a decrease of 13.02%. This comparison validates the IGA and the optimal results can provide a reference for engineering design. coupling flexibility matrix B ij coupling stiffness of row i and column j in matrix D bending flexibility matrix D ij bending stiffness of row i and column j in matrix F i , F ij strength parameter, i, j = 1, 2, 6 respresnting x, y, and x-y direction, respectively R the average value of the initial data R avg mean fitness values of individuals in each generation R i the i-th initial data R max maximum fitness in each generation of individuals R min minimum fitness value in each generation of individuals R s fitness function S flexibility matrix S 11 flexibility coefficients S 16 flexibility coefficients S 66 flexibility coefficients S ij flexibility coefficient of row i and column j in matrix S sp individual space transverse tensile strength n 1 number of same layering angles for four consecutive layers n 2 number of layers that do not meet the ratio requirement x i the kth gene P C crossover probability P m mutation probability length of laminate E 1 longitudinal modulus of elasticity E 2 transverse modulus of elasticity G 12 shear modulus of elasticity in the plane of the lamina h thickness of laminates k kth layers M x the bending moment of x-direction M xy the torque of x-y plane M y the bending moment of y-direction n total number of layers N x the force of x-direction N xy the shear force of x-y plane N y the force of y-direction Q shear stress q uniform loading pressure R the ratio of a certain component of the ultimate stress to the corresponding component of loading stress. S the shear strength in x-y plane T the indicator to evaluate the degree of fitting x the fiber direction y the vertical fiber direction Z k the z coordinate of kth layer Z k−1 the z coordinate of k − 1th layer µ 1 poisson's ratio in the x-direction ρ density ϕ stress function U genetic operator evaluates parameters l half-length of laminate ε an infinitesimal value

Appendix A
Considering a composite beam with simply supported ends, it has σ z = 0 τ yz = τ zy = τ zx = τ xz = 0 (A1) Here, only σ x , σ y , τ xy are not zero. Then according to statics, it obtains the following expressions. Let stress ∅ be a function of x, y, the following relationship holds: • Geometric compatibility equations The expressions of σ x , σ y , τ xy from (A7) as below 1 Since σ y is caused by the constant shear force q, and its value does not depend on x, it is a function of y. Therefore, we can set σ y = f (y) Through integral formula (A8), it can produce: 2 Therefore, substituting (A9) into the compatible equation, each item of (A7) can be obtained as follows: 2S 12 + S 66 S 11 y 4 + C 1 6 y 3 + C 2 2 y 2 +C 3 y + C 4 (A15) Because (C 3 y + C 4 ) is not used in next equation, it can be omitted. Finally, through combining the Equations (A13)-(A15), it obtains stress function as below, ∅ = 1 2 x 2 f (y) + x f 1 (y) + f 2 (y) = 1 2 x 2 A 1 y 3 6 + A 2 y 2 2 + A 3 y + A 4 + 1.
When the upper and lower boundaries y = ±h/2, τ xy = 0 (this holds for any x), through the above equation, the relational expressions of A 2 , A 3 , B 1 , B 2 , B 3 , and A 1 can be derived.