A Level-Set-Based Density Method for Buckling Optimization of Structure with Curved Grid Stiffeners

: Curved grid stiffeners, compared to straight stiffeners, offer greater flexibility in adjusting the force transmission paths and give better structural performance. In this paper, a level-set-based density method is employed to generate layouts of curved grid stiffeners so that the critical buckling load factor (BLF) of the stiffened structures is improved. During the optimization process, volume constraint is incorporated to control material utilization, and gradient constraints are employed to maintain uniformity in the width of the stiffeners. Finally, the proposed method is demonstrated through several numerical examples


Introduction
Stiffeners can effectively enhance structures in bending stiffness, strength.Therefore, they have been widely utilized in aerospace, automotive, architectural, and other fields [1][2][3].Curved grid stiffeners typically consist of many curved and crossed ribs, and they can adjust the load paths according to the different requirements of different local regions of the structure.Therefore, they have greater potential in reducing structural weight and improving structural performance [4][5][6][7][8][9][10], and there has been growing interests in the layout optimization of curved grid stiffeners [6,8,[11][12][13][14][15].
The stability of a structure is an important indicator of its performance.When subjected to axial compression and the load reaches a certain value, structures often experience buckling failure [16].Therefore, the buckling performance of structures becomes a key factor in design and optimization.In related studies, Cai et al. [17] employed generalized beam theory to analyze buckling in thin-walled structures.Ferrari et al. [18] presented a topology optimization of structures under plane stress condition for linear buckling with 250 lines of code.Xu and others [19] presented an innovative linear interpolation approach to improve the buckling resistance of structures in topology optimization.In related experimental research of buckling response, Falkowicz [20] performed both linear and nonlinear buckling analyses, as well as experimental tests, on a thin-walled composite plate with a central cutout under axial compression.In latest developments [21,22], asymmetric systems and their mechanical couplings have been utilized to create elements that can function as elastic components, and comprehensive experimental and numerical studies have been conducted using multiple detection techniques and interdisciplinary approaches.
For stiffened structures, the introduction of stiffeners significantly improves the buckling resistance.In addition, optimizing the layout of stiffeners can further amplify this advantage [4,7,8,23].The design optimization of stiffened structures considering buckling characteristic has become an important research topic.Kapania et al. [4] optimized stiffened plates primarily governed by buckling phenomena and demonstrated that, compared to straight stiffeners, curved stiffeners provide enhanced design space, resulting in better design outcomes.A buckling analysis of composite plates was proposed [8], indicating that arranging curved stiffeners for composite plates can improve buckling load compared to straight stiffeners.Shi and co-workers [24] explored the static, vibration, and buckling behaviors of stiffened plates under in-plane compressive and tensile stresses.Hao et al. [25] and Wang et al. [26] utilized the smeared stiffener method (SSM) and the Rayleigh-Ritz method to forecast and determine the buckling load of hierarchical stiffened plates.Paschero and Hyer [5] significantly improved the buckling performance for elliptical lattice cylinders by adjusting the distribution of circumferential stiffeners.An effective optimization strategy for cylindrical structure with reinforced cutouts was proposed by Hao et al. [7], employing curved stiffeners to enhance the load distribution and local stiffness of cylindrical stiffened shells.Luo et al. [27] proposed an approach for optimizing the arrangement of stiffeners in slender structure, considering the maximum critical buckling load under geometric uncertainties that vary spatially.Wang et al. [13,28] calculated the global and local buckling loads of periodically grid-stiffened composite plates and shells and optimized the layout of curved stiffeners.A data-driven optimization approach was presented to maximize the critical BLF [29].These research works on the design optimization considering buckling of stiffened structures further promote the progress and development of structural optimization, especially in the field of structures with curved stiffeners.
Theoretically, structures have many linear buckling modes linked to relevant critical load values.But in engineering applications, the smallest critical load value, namely, the critical buckling load factor (BLF), is under consideration, as exceeding this load leads to structural instability and failure.Therefore, in optimization, it is essential to maximize the critical BLF of the structure to enhance its stability [16,28,30,31].However, it should also be considered that in practice, structural instability often occurs before reaching the critical buckling load.This is due to the potential presence of imperfections in the structure or the material reaching its yield point.
In this paper, a level-set-based density method proposed in our previous study [15] is used to optimize the critical BLF of structures with curved grid stiffeners.This provides new insights for the buckling research of stiffened structures.Two level set functions are converted by cosine function, combined using the maximum function, and projected to real densities using an approximate Heaviside function.The Mindlin plate theory is employed for finite element analysis of buckling.After that, the sensitivity analysis is performed, and the design variables are updated.Ultimately, the optimized results are obtained.The effectiveness of this approach in enhancing structural stability is proven by several numerical examples.
The main structure of the paper is outlined as follows.Section 2 presents the definition of the grid stiffener structure.Section 3 discusses the relevant constraints for stiffener optimization.The optimization framework and sensitivity analysis are detailed in Section 4.
Several numerical examples and result analyses are presented in Section 5. Section 6 gives the conclusions.

Definition of Stiffener Height
The total height h e at the center of the e-th element of the grid-stiffened structure is obtained by adding the base height h b and the stiffener height h s .h e is calculated from the function h(x) at the element center x e = (x e , y e ), and h(x) is given by where ρ(x) is regarded as the real density of stiffeners, describing their layout.When ρ(x) is 1, it is considered as stiffeners being present, while values close to 0 indicate the absence of stiffeners.
The real density ρ(x) is derived through a threshold projection of the function ρ(x).The threshold projection employed in this paper is an approximate Heaviside function [32] ρ where the value of parameter τ affects the approximation of Equation ( 2) to the ideal Heaviside function; η controls the location of the step point of the function.As shown in Figure 1, as τ increases, this function approaches the ideal Heaviside function more closely.In this paper, τ is set to 7. In addition, as depicted in Figure 1, when η = 0.5, the Step Point 1 is positioned at x = 0.5.When η = 0.The function ρ(x) is obtained by taking the maximum value from two functions, Γ(x) and Θ(x), i.e., ρ(x) = max{Γ(x), Θ(x)} where the function max{•, where abs ε (a) = √ a 2 + ε; and ε is set to 0.01.The functions Γ(x) and Θ(x) are derived from two functions γ(x) and θ(x) through cosine transformations [33], i.e., where γ(x) and θ(x) represent fundamental level set functions; µ 1 and µ 2 correspond to the wavelengths of two cosine functions.It can be easily observed that the ranges of the functions Γ(x) and Θ(x) are both between 0 and 1.When τ in Equation ( 2) is sufficiently large and the norms of gradient vector of the fundamental level set functions γ(x) and θ(x) equal 1, the width w s of the stiffeners can be expressed as [15] where µ can be µ 1 or µ 2 .From Equation ( 7), it can be observed that the width w s varies in proportion to the wavelength µ of the cosine functions, and there is a relationship between the parameters w s and η.Therefore, η can be used as a design variable to change the width of the stiffeners.
Based on the aforementioned theories, an example process for generating grid stiffeners is provided.Figure 2 shows the images of the fundamental level set functions γ(x) and θ(x).And Figure 3 shows their corresponding contours.It can be observed that the contours of the level set functions are at 0 • and 90 • .The colors transition from blue to green from bottom to top, indicating a gradual change from smaller to larger values in Figures 2 and 4   Figure 4 shows the images of functions Γ(x) and Θ(x), which are obtained by applying the cosine transformation in Equations ( 5) and (6) to the level set functions γ(x) and θ(x), respectively.When projected onto a horizontal plane, it contains width information.The grayscale images of Γ(x) and Θ(x) are shown in Figure 5.
After applying the maximum combination from Equations ( 3) and ( 4), Figure 6a is obtained.Subsequently, the threshold projection in Equation ( 2) is used to generate the final real density distribution map, as seen in Figure 6b.It can be observed that the threshold projection effectively reduces the amount of intermediate density, making the grid stiffener structure clearer and more visually intuitive.

Fundamental Level Set Function Based on RBFs
The radial basis functions (RBFs) [34][35][36][37][38] are utilized here to construct the fundamental level set functions γ(x) and θ(x), i.e., where α i and β i represent the coefficients of the i-th RBF; q i is the knot of the i-th RBF; m stands for the total number of RBFs.The RBF at knot q i is denoted as where the ∥ • ∥ represents the Euclidean norm.
In this study, the compactly supported RBF (CS-RBF) with C 2 continuity is utilized.It is expressed as [37] φ where (•) + represents max{•, 0}; t is formulated as where d s is the support radius of CS-RBF; κ is set as 10 −4 to prevent division by zero.

Constraint of Volume
In order to ensure that the optimization meets practical volume or mass requirements, material volume constraint for stiffeners is defined as where V denotes the total volume of stiffeners; h s is the thickness of stiffener; A e denotes the area of the e-th element; and V stands for the maximum allowable volume of stiffeners.

Constraint of Uniform Width
The uniformity of stiffener width holds substantial significance for practical engineering manufacturing because in some applications, stiffeners with uniform widths are easier to manufacture.One effective approach to achieve this goal is to ensure that the norm of the gradient vector of the fundamental level set functions is equal to 1 across nearly all points.The norm of the gradient vector of γ(x) is provided by where ∂φ/∂x and ∂φ/∂y are, respectively, the partial derivatives of φ with respect to x and y.According to Equation (11), they are defined as [37] ∂φ ∂x (x, The ∥∇θ(x)∥ of the function θ(x) can be obtained similarly.
The subsequent two constraints can ensure that at the center point x e of each element, ∥∇γ(x e )∥ and ∥∇θ(x e )∥ are close to 1. Let d ϕ e and d ψ e denote the constraints of the e-th element on γ(x) and θ(x).Following the methods in [39][40][41], equality constraints can be transformed into inequality constraints as The p-norm approach [42] is utilized to consolidate the constraints as where ξ is the upper limit of the two constraints; n is the amount of finite elements; and p 2 > 0 is the parameter of the p-norm.

Optimization Framework and Sensitivity Analysis 4.1. Definition of Buckling Load Optimization Problem
The buckling optimization problem is defined as where λ 1 represents the fundamental BLF; φ 1 denotes the eigenvector corresponding to the eigenvalue λ 1 ; K stands for the global stiffness matrix, while G represents the global geometric stiffness matrix; δ min and δ max , respectively, stand for the lower and upper limits of the design variables α i and β i ; η min and η max indicate the minimum and maximum values of the design parameter η; g γ and g θ are the gradient constraint functions for the two fundamental level set functions, and ξ represents the maximum value of the gradient constraints; V denotes the total volume of stiffeners, while V indicates the volume constraint of stiffeners.

Finite Element Analysis
The Mindlin plate theory [43][44][45][46] is used in the finite element analysis.Additionally, we assumed the substrate and stiffeners are made of the same material.Therefore, the stiffness matrix k e for each element is defined as where k b , k s and k m represent the bending, shear, and axial stiffness matrices of the element, respectively, and they can be specified as where B b , B s and B m are the strain-displacement matrices for bending, shearing, and axial deformation, respectively; D b , D s and D m are the corresponding elastic matrices.Ω e represents the element domain.The stiffness matrix of finite element in Equation ( 22) can be modified to The geometric stiffness matrix is defined as where g represents the partial derivative of the shape functions, and σ x , σ y , and τ xy are the plane axial and shear stresses.The expression for g e in Equation ( 27) can be modified to where g 1 depicts the geometric stiffness matrix of unit thickness.

Sensitivity Analysis
The sensitivity of the optimization objective J with respect to the design variables α i is given by From Equation ( 21), the partial derivative of the critical BLF λ 1 with respect to α i can be acquired as (likewise, the partial derivative of λ 1 regarding the design variables β i can be derived) From Equations ( 26) and ( 28), we have Based on Equations ( 1 where Next, the sensitivity of the optimization objective J with respect to the design variable η is formulated as Similarly, based on Equations ( 26) and ( 28), we have ∂g e ∂η = ∂h e ∂η g 1 (42) According to Equations ( 1) and (2), we have where where In accordance with Equation ( 13), the sensitivity of the volume constraint to α i , β i , and η are Based on Equation ( 19), the sensitivity of the gradient constraint g γ to α i is According to Equations ( 14)-( 17), we have Under the same conditions, the sensitivity of the gradient constraint g θ concerning the design variables β i can be acquired.All design variables are updated by using the method of move asymptote (MMA) [47].

Numerical Examples
In this part, several numerical examples of optimizing the fundamental BLF for the structures with grid stiffeners are provided.For all examples, the Young's modulus is set to 200 GPa, the Poisson's ratio is set to 0.3, the thickness of the base plate h b is 0.5 mm, and the thickness of the stiffener h s is 1 mm.The parameter τ in Equation ( 2) is set to 7, the parameter p 2 in Equations ( 19) and ( 20) is 12, the support radius parameter d s of CS-RBF is set to 8, and the maximum and minimum limits for the design variables in optimization problem are set as follows: The convergence condition for optimization is characterized as where E err represents the error between the relevant objective functions during the optimization process; k indicates the current iteration step; δ E is the upper limit of the error, and it is set to 0.2%.Additionally, if the iteration step k reaches 600, the optimization procedure will be terminated.

Example 1
Within this example, a square plate with dimensions of 10 mm × 10 mm, supported on both sides with simply supported boundary conditions, is considered.Uniform compressive loads of magnitude F = 1 N/mm are applied simultaneously on both sides.The boundary conditions are illustrated in Figure 7: the bottom-left corner of the plate is restricted in the x and y-direction translational degrees of freedom (i.e., u = 0 and v = 0), the translational degree of freedom in the y-direction at the bottom-right corner is constrained, and the z-direction degrees of freedom on all four sides of the plate are fixed (i.e., w = 0).The wavelength parameters in Equations ( 5) and ( 6) are set to µ 1 = µ 2 = 1.4.The design domain consists of 200 × 200 finite element elements, with 20 × 20 CS-RBF knots distributed uniformly within the design domain.To maintain uniformity in the width of the stiffeners post-optimization, the parameter ξ in Equations ( 19) and ( 20) is set to 0.3.
Several upper bounds of volume are applied in this example, and the stiffeners of the initial design are set in two directions: horizontal (0 • ) and vertical (90 • ), as shown in Figure 8a.In Equation ( 13), the volume constraints for stiffeners are defined as V = 0.8V 0 , V = V 0 , and V = 1.2V 0 , respectively, where V 0 is the initial volume of the stiffeners.The optimization results are illustrated in Figure 8b-d  It can be seen that despite varying volume constraints on the stiffeners, the layout of the stiffeners remains essentially consistent.The effect of volume constraint manifests in the width of the stiffeners, where a larger volume leads to wider stiffeners.Regarding the layout, it is notable that the distance between stiffeners increases in the central region.This enhances the structural strength in the central region.The path of the stiffeners forms outwardly curved arcs, aiming to enhance the load-bearing capacity against axial compression loads on both sides.Figure 9a-c depict the corresponding first buckling modes of the three optimized structures, and it should be noted that the buckling mode in this example, as well as in all subsequent examples, represents the structure's normalized displacement in the w-direction.
Figure 10a-c depict the first two BLFs and volume variation curves under three different volume constraint scenarios.It is apparent that there is no intersection between the curves of the first and second BLFs, indicating the absence of repeated first BLFs.Additionally, the volume changes for each case follow the prescribed volume constraint values, eventually reaching their respective maximum allowable volume.Furthermore, it can be seen that the optimization of the critical BLF converges under all three constraint conditions, demonstrating the validity of the optimization.
Table 1 presents the corresponding outcomes of the optimization, indicating a clear trend where increasing volume leads to higher buckling load.Under the maximum allowable volumes of 1 times and 1.2 times the initial volume, the increment in critical buckling load reaches 13.3% and 44.5%, respectively.When V = 0.8V 0 , the buckling load decreases on account of the reduced volume of the stiffeners.

Table 1. Comparison of results under divergent volume constraints (λ 0
1 , V 0 : initial critical BLF and total volume of stiffeners; λ # 1 , V # : final critical BLF and total volume of stiffeners;

Example 2
In this case, the size of the design domain is the same as that in Example 1, with the left end fixed and a uniform compressive load of magnitude F = 1 N/mm applied on the right side, as shown in Figure 11.The wavelength parameters in Equations ( 5) and (6), as well as the finite element mesh division and the distribution of CS-RBF nodes, are the same as those in Example 1.The maximum allowable volume V is configured as V 0 .The initial layout of the stiffeners is in two directions, 45 • and −45 • , as shown in Figure 12a.Both scenarios with and without gradient constraints are considered.In the case of with gradient constraints, the gradient constraint values in Equations (19) and (20) are set to ξ = 0.06 and ξ = 0.3, respectively.The optimization results are illustrated in Figure 12b-d.It is clear that the optimization results with gradient constraints have more uniform width of the stiffeners.Moreover, the stricter the gradient constraint, the more pronounced this effect.Regarding the stiffener layout, it is evident that the stiffeners on the left side are optimized to align with the direction of force transmission (horizontal direction in the figure), with increased width, thereby enhancing the ability to withstand axial compressive loads more effectively.Figure 13a-c, respectively, depict the first buckling modes under gradient constraints of ξ = 0.06, ξ = 0.3, and without gradient constraints.
Figure 14 illustrates the variation of gradients of the level set functions during the optimization process.Figure 14a-c, respectively, represent the cases with gradient constraints ξ = 0.06, ξ = 0.3, and without gradient constraints.One can see that under the gradient constraint conditions, the gradients eventually reach the set values, satisfying the gradient constraint requirements.In the absence of gradient constraints, the gradient values are ultimately optimized to 0.54.
In Figure 15, the history of the objective function and volume are displayed under three different cases.It can be observed that the volume remains close to a constant value throughout, validating the effectiveness of the volume constraint.The optimization magnitude of the critical BLF increases as the gradient constraint weakens.Particularly, the optimization magnitude is greatest when there is no gradient constraint.Hence, relaxing the gradient constraint can lead to better performance.In Table 2, it is clear to see the optimization magnitude is 26.1% with a gradient constraint value of ξ = 0.06, 45.2% with ξ = 0.3, and 47.1% without gradient constraints.

Example 3
The design domain of this case study is a rectangular shape with dimensions of 10 mm × 20 mm.The bottom edge is completely fixed, and a vertically downward distributed load is applied on top with a length of 2 mm and a magnitude of F = 1 N/mm.Boundary conditions are depicted in Figure 16.The wavelength parameters in Equations ( 5) and (6) are set to µ 1 = µ 2 = 2.1; 200 × 400 finite elements is used; and 20 × 40 CS-RBF knots are uniformly distributed.The gradient constraint parameter ξ in Equations ( 19) and ( 20) is configured as 0.3.The maximum allowable volume V is set to be V 0 .In this example, the initial layout of the stiffeners is oriented in two directions, at 45 • and −45 • , as illustrated in Figure 17a.Figure 17b exhibits the final outcome, where the stiffeners distribution at the bottom, near the center region, resembles that of Example 2, optimized in the direction of the load path, thereby enhancing the load-bearing capacity.The stiffeners in the upper region follow an arch-shaped path, similarly aimed at improving the load-bearing capacity. Figure 18a-c represents the first-, second-, and third-order buckling modes corresponding to the optimization results, respectively, demonstrating the absence of duplicate first-order BLFs.
The plots depicting the evolution of the critical BLF and the volume of stiffeners are illustrated in Figure 19.The optimization procedure converged eventually, with the critical BLF λ 1 increasing by 29.8%.

Conclusions
In this study, we optimized the critical buckling load of curved grid stiffener structures based on a level-set-based density method.By transforming two fundamental level set functions, two clusters of stiffeners were obtained, and the real density of each element was obtained using density-based methods.Then, the overall structural element height distribution was obtained by adding the heights of the matrix and stiffeners.The Mindlin plate theory was used for finite element analysis of buckling, and the stiffness matrix and geometric stiffness matrix of each element were determined and assembled into the overall matrix.After that, the eigenvalue equation was solved to obtain the corresponding eigenvalues and eigenvectors.Combining this with the optimization objective, relevant sensitivity analysis calculations were performed.Simultaneously, volume and gradient constraints were added to improve material utilization and maintain uniform stiffeners width.Finally, the MMA algorithm was used to update the design variables.The effectiveness of the proposed method was validated through several numerical examples.
In practical applications, this research method improves the stability of grid stiffener structures under specific load conditions in engineering applications.Optimizing the stiffener layout instead of increasing the stiffener thickness or base wall thickness helps the structure move towards lightweight development.In addition, introducing gradient constraints to control the minimum spacing between stiffeners enhances the manufacturability of the structure.In the future, the layout and shape of the obtained grid stiffener structure can be subjected to certain post-processing to meet further requirements, and specific experimental analysis can be added for verification and explanation.In addition, the method can be combined with topology optimization to change the shape of base structure while optimizing the layout of curved grid stiffeners. .

Figure 2 .
Figure 2. Examples of the fundamental level set functions γ(x) and θ(x).

Figure 3 .
Figure 3. Contour lines of the level set functions.

Figure 4 .
Figure 4.The functions after the cosine transformation.

Figure 6 .
Figure 6.The images of the function after maximization combination and threshold projection. .

Figure 7 .
Figure 7.The design domain and boundary conditions for Example 1.

Figure 11 .
Figure 11.The design domain and boundary conditions for Example 2.

Figure 16 .
Figure 16.The design domain and boundary conditions for Example 3.

Figure 17 .
Figure 17.Initial design and optimized result.