Research Article Fixture Design in Flexible Tooling of Aircraft Panel Based on Thin Plate Theory

Flexible and conformal positioning of the skin is of great significance for the efficient and low deformation assembly of aircraft panels. As the main conformal positioning part of the skin, the fixture positioning surface is complicated to design when taking into account the flexible positioning. In this study, taking the low deformation and flexible clamping of the skin under complex stress conditions as the target, a fixture positioning profile design method is proposed to predict the deformation of the skin and provide moderate support for the skin based on thin plate theory (TPT). In the face of complex stress and boundary conditions of the skin, the differential method is used to solve the stress and deformation of the skin under different support conditions, which is correctly verified by finite element simulation. The proposed method reveals the influence of the fixture positioning profile on the deflection of the aircraft skin under concentrated load (drilling force), which provides a good design theoretical basis for the fixture design for the low deformation flexible clamping of the skin. relation in large deflection theory.


Introduction
In the fuselage, wing, and other structures of aircraft, thinwalled parts are widely used to reduce the weight of the whole aircraft. e clamping deformation and assembly deviation of this kind of thin-walled parts in the assembly stage will accumulate and spread along the assembly dimension chain of products [1] and ultimately a ect the quality of aircraft manufacturing [2]. So how to design a suitable xture for thin-walled parts is of great signi cance in the eld of aircraft assembly.
Aiming at the problem of xture design for weak sti ness parts, Cai et al. [3] came up with the "N-2-1" positioning principle in 1996. e "N-2-1" positioning principle is more suitable for the clamping of thin-walled parts than the "3-2-1" positioning principle by increasing the number of positioning elements on the main positioning surface and enhancing the overall sti ness of thin-walled parts. e increase in the number of positioning elements causes trouble in the xture arrangement. e intelligent optimization algorithm has the characteristics of high search e ciency and convenient mathematical description, which is commonly used in xture positioning point layout optimization. For xture layout optimization, di erent xtures develop di erent solutions. Minh et al. [4] presented a single-axis rod type exible xture system for thin-walled components in machining processes utilizing the N-2-1 location principle and proposed a geometry-based method to nd the optimal location of the workpiece relative to the xture system, which could nd the optimal location within 30 s. Chen et al. [5] proposed a "N-M" positioning principle based on the "N-2-1" principle and used the genetic algorithm to search the optimal solution of xture layout, which reduced the machining deformation of thin-walled parts by 56.5%. Abolfazl et al. [6] used the genetic algorithm to optimize the unconstrained problem of car body parts clamping error propagation and realized the optimal fixture configuration. Ma et al. [7] proposed a hybrid optimization algorithm GAOT, which generated the global optimal fixture layout, and obtained a better optimal solution than the existing optimization algorithm. Zhou et al. [8] designed the fixture layout scheme based on the hybrid particle swarm optimization algorithm with the objective of minimizing the normal deformation of metal sheet in the flexible tooling system. Wan et al. [9] applied the fast nondominated multiobjective optimization algorithm (NSGA-II) to optimize the four defined objective functions to minimize the position error of the workpiece and obtain the optimal positioner layout. Yang et al. [10] proposed a fixture location layout optimization method based on the support vector regression model and NSGA-II algorithm and realized the optimal design of support points for aircraft skin parts. Zeshan et al. [11] proposed a "N-3-2-1" positioning principle for sheet metal fixture layout optimization and combined with finite element and genetic algorithm to realize fixture layout optimization in the multipoint spot welding process. Li et al. [12] proposed a method of fixture location layout optimization design for thin-walled parts based on the Kriging agent model and flower pollination algorithm and completed the "4-2-1" fixture layout optimization design for curved thin-walled parts. Aimed at aircraft weak stiffness parts, Wang et al. [13] calculated the maximum deformation under different positioning schemes by finite element simulation and combined with the firefly algorithm to achieve the iterative optimization of positioning point layout. Li et al. [14] proposed an optimization method for the clamping scheme of aircraft thin-walled parts based on the genetic algorithm. Combined with finite element simulation, the synchronous optimization of clamping layout and clamping sequence was realized. To sum up, in the design of thin-walled parts fixture, the application of "N-2-1" positioning principle, and the positioning point layout optimization guided by the intelligent optimization algorithm, can significantly reduce the deformation of thin-walled parts in the assembly and manufacturing process and improve the final quality of products.
In the processing of thin-walled parts, the workpiece and cutting tool will have elastic deformation under the action of cutting force, which will lead to the change of process parameters and affect the processing quality.
in plate problem is a basic problem in elastic mechanics. e key to study this kind of problem is how to deal with the boundary conditions and how to get the exact analytical or numerical solution [15][16][17]. Based on the thin plate theory of elasticity, the machining deformation of thin-walled parts can be analyzed. e solution of thin plate problems usually includes the finite element method, meshless method, finite difference method, boundary element method, Rayleigh-Ritz method, and Galerkin method [18]. Slark et al. [19] studied the large deformation problem of simply supported thin plates by the meshless method and solved the nonlinear equations. Based on Kirchhoff plate theory, Shabana et al. [20] proposed a 12-node plate element model, in which the displacement continuity can be ensured by rigid connection between nodes, and the motion and deformation of flexible plate can be accurately described. Chen et al. [21] proposed a method of machining error compensation for thin-walled parts based on bicubic B-spline interpolation considering the influence of elastic deformation. In order to reduce the machining deformation of thin-walled titanium alloy parts, Li et al. [22] proposed a nonuniform allowance design strategy of discrete allowance volume element based on the Rayleigh Ritz method, which makes better use of the workpiece stiffness and reduces the machining error. Tang et al. [23] solved the bending problems of Kirchhoff and Winkler thin plates by using the generalized finite difference method. Numerical experiments show that the method can effectively solve the bending problems of two kinds of thin plates under different transverse loads.
In addition, as a large aviation thin-walled part, the assembly sequence of aircraft panels has a crucial impact on the assembly efficiency and assembly quality of the panels. In order to improve assembly efficiency and quality, many scholars home and abroad began to study the mechanical assembly sequence as early as 4 decades ago. e optimization of the assembly sequence in the early stage was mainly realized by manual calculation, and its defects were low efficiency and single results. With the rapid development of computer computing power and the introduction of intelligent algorithms, the optimization of assembly sequence has become efficient, flexible, and rich in calculation results [24]. Yang et al. [25] developed a novel PSOBC algorithm based on the conventional PSO algorithm for complex products; this method avoids a large number of matrix calculation by establishing synchronized assembly Petri net (SAPN) to describe the precedence relationships, and the provided PSOBC not only prevents premature convergence to a high degree but also keeps a more rapid convergence rate than the standard PSO algorithm. Bahubalendruni et al. [26] proposed a novel and efficient method to obtain all valid assembly sequences and optimized assembly sequence for a given assembled product, which considers four basic predicates, namely, "liaison predicate, geometrical feasibility, mechanical feasibility, and stability" to validate each sequence, and the validity of the method is validated by different example products. Bahubalendruni et al. [27] proposed an efficient computational method to find the optimal sequence from a huge set of all assembly sequences, and the method is proven in generating optimal solutions for any given product effectively.
Based on the above references, it can be found that both scholars domestic and abroad have made a lot of achievements in the research of thin-walled parts fixture design. As a kind of large thin-walled parts, how to reduce the load deformation of aircraft skin and achieve reliable clamping is the paramount problem of aircraft panel flexible tooling development. Scholars domestic and abroad have made research on the fixture layout and deformation analysis of thin-walled parts, but the research on the fixture layout mainly focuses on the point supporting flexible tooling, and the research on the typical fixture type tooling like surface supporting is too little. Few scholars established the mechanical model to guide fixture design through elastic deformation theory of aircraft skins. In this study, in order to minimize the deformation of aircraft skins under load, the research on the design of exible xture is carried out.

Analysis of Fixture Design
e skin of large aircraft is usually designed with variable thickness. e thickness of the skin is reduced in the area with less stress, which forms the inner surface of the skin as concave, so as to reduce the weight of the whole skin structure. Figure 1 shows the skin of the front fuselage panel of a certain civil aircraft, and the blue areas are the chemical milling area of the skin. When it comes to the exible tooling xture, the chemical milling area of di erent skin may be di erent, which makes the xture surface di cult to coordinate or even share the same xture surface. In an ideal situation, the longer the length of the xture surface is, the longer the length of the e ective bonding section with the skin is, and the better the xture sti ness can be achieved. However, due to the location of the truss axis, the length of the clamping surface is limited within a certain range, as shown in Figure 2 erefore, it is necessary to study the quantitative relationship between the length of the surface and the e ect of skin support sti ness to guide the design of the surface.
In this study, aiming at the design problem of the xture, the design goal is to ensure the e ect of skin support, and the solution of thin plate problem in the elasticity is used to guide the xture surface design according to the suppression e ect of the elastic deformation of the skin.

Mechanical Model Based on Thin Plate Theory
In the assembly process, the skin is a ected by gravity, strap pressure, board support force, frictional force, and drilling force, so it is often necessary to establish a complex mechanical analysis model, which increase the di culty of solving, so the establishment of a simpli ed mechanical model is an important part of solving engineering problems. For the case where the skin is thin and large and only the stress-strain conditions near small load area are analyzed, the classical thin plate deformation theory is adopted, where the skin radian is ignored, the skin is considered to be an ideal elastic body, and the theoretical model is based on elasticity of a rectangular thin plate. Hence, several assumptions should be given in advance [28]: (1) Zero radian assumption: since the stress area is small when drilling, the surface radian of the skin can be ignored; (2) Continuity assumption: it is assumed that the object is a continuous medium and there is no gap between particles; (3) Complete elasticity assumption: the object completely obeys Hooke's law, and the strain of the object is proportional to the stress; (4) Isotropic assumption: the physical properties of any point in the object are the same in all directions; (5) Uniformity assumption: the object is composed of the same type of uniform material, so that the elastic coe cient of the object does not change with the change of position coordinates.
e thin plate will be bent under the action of transverse load, and the stress, strain, and displacement problems caused by it are calculated according to the bending problem of thin plate [29]. e bending problems of thin plates are generally divided into small de ection bending problems and large de ection bending problems. After the skin is positioned and clamped on the frame, the preconnection holes between the skin and truss, and the skin and corner piece need to be drilled. In this process, the xed e ect of the strap and the clamp on the skin is regarded as the xed boundary condition of the skin, and the skin is only a ected by the drilling force. e drilling force is a transverse force and perpendicular to the neutral plane. erefore, the stress and deformation of the skin under the action of drilling force can be solved as the bending problem of thin plate.

Mathematical Problems in Engineering
When dealing with the problem of thin plate bending, it can be considered as the problem of small de ection bending of thin plate when the sti ness of thin plate is large and the de ection is far less than the thickness under the action of external load.
(1) On any normal line perpendicular to the neutral plane, the displacement of each point in the thickness of the thin plate is the same, that is, the de ection ω, which is only related to the x and y coordinates. Let the neutral plane of the thin plate be the xy plane ( Figure 3), so that (2) e deformation caused by the normal stress σ z is negligible, so the equation can be simpli ed to (3) ere is only vertical displacement and no displacement parallel to the neutral plane in the thin plate, that means e bending deformation problem of the skin is simpli ed to the small de ection bending problem of the thin plate, and a di erential equation for ω(x, y) needs to be established. As shown in Figure 3, take out a parallelepiped microunit, whose length and width are dx and dy, respectively, from the rectangular thin plate with thickness δ. e four sides of the unit are subjected to the internal force of the plate, the upper surface is subjected to the external load qdxdy, and the lower surface has no load. On the side where x is a constant, stress components s x , τ xy , and τ xz act; on the side where y is a constant, stress components σy, τ yx (equal to τ xy ), and τ yz act. ese stresses can be synthesized into bending moments M x , M y , torque M xy , and transverse shear forces F sx , F sy , respectively.
Combine the stress components per unit width of the ank: Substituting (5) into (4), the bending moment and torque of the microunit can be obtained as In order to keep balance of the microunit, the sum of the moments around x-axis and y-axis is 0, and the force in zdirection is balanced, that means the following balance conditions must be satisfied: After omitting the high-order trace, it can be obtained that and zF sx zx Regarding Mxy�Myx, substitute (7) and (9) into equation (10) to get Arrange (11) as follows, which is the equilibrium differential equation for bending deformation of the skin with small deflection.
where ∇ 2 � z 2 /zx 2 + z 2 /zy 2 is the Laplace operator. However, the assumption for small deflection theory is that the longitudinal displacement of thin plate, which is perpendicular to the neutral plane is so small than the thickness of thin plate that the membrane force can be ignored. For a thin metal plate such as the aircraft skin, the longitudinal displacement of each point in the plane is not necessarily much less than the thickness of the plate. If so, it can be solved as a large deflection bending problem of thin plate. When establishing the equilibrium differential equation of large deflection bending problem of thin plate, it is necessary to take into account the neutral plane strain and membrane force caused by the longitudinal displacement of each point in the neutral plane. e bending differential equation of the elastic plate is e deflection ω and the film forces (caused by lateral load q) F tx , F ty , and F txy in this equation are unknown. In order to simplify the equation and reduce the number of unknowns, the stress function Φ is introduced, so that At this point, the elastic plate differential equation becomes an equation containing two unknown ω and Φ, which still cannot be solved. A compatibility equation needs to be constructed for the relationship between film force and deflection. Considering the geometric relationship between the normal strain, shear strain, and three-dimensional displacement in the microunit of the thin plate, the following geometric relationship equation can be obtained: By eliminating the lateral displacements u and v in the equation, the following compatibility equation can be obtained: According to (2) and (14), the compatibility equation can be simplified; then, the differential equations of elastic plate can be simultaneously combined to obtain the following differential equation system based on the large deflection theory: Mathematical Problems in Engineering

Theoretical Calculation Based on the Difference Method
In the process of using small deflection theory or large deflection theory to solve the skin deformation problem, the same boundary of the skin may contain both fixed boundary and free boundary or other types of boundaries due to the panel is partiality supported by the board. Moreover, during the drilling or milling, the distribution of the load and size of the load is irregular; thus, the description of the load is also complex. To sum up, there are mathematical difficulties in solving the differential equation of skin deformation because the mathematical models are often not that accurate to include the different boundary conditions of different skin segments and the randomly distributed load conditions. Difference method is a common method to solve partial differential equations. e basic method is as follows: the domain is discretized by difference grid, the independent variable is represented by the function of the discrete variable at the point of difference grid, and the derivative or partial derivative is represented by the difference quotient of the function between the points of difference grid. e original differential equation and boundary condition are approximately expressed by difference equation in the form of algebraic equation, and then, the differential equation is replaced by the algebraic equation. It is considered that the solution of the algebraic equation is the approximate solution of the original problem.

Difference Calculation Based on Small Deflection eory.
For the convenience of calculation, the skin area is divided by square mesh, and the side length of square mesh is represented by Δ d. As shown in Figure 4, the skin is divided into (i − 1) × (j − 1) mesh regions, including i × j nodes.
Using the square grid as shown in Figure 4, from column x � 1 to column x � i, write the difference equations for the unknown nodes in each column and list the difference equations. By using the matrix form, we can get the following results: where K is the coefficient matrix of order (i×j)×(i×j) of the difference equations, which can be expressed as where A, B, C, are the block matrices, all of which are square matrices of order j × j. If the number of rows is j and the number of columns is i in this matrix, A, B, and C are the coefficient matrices corresponding to the internal i-4 column nodes, A 2 , B 2 , C 2 , D 2 , A 1 ′ , B 1 ′ , C 1 ′ , A 2 ′ , B 2 ′ , C 2 ′ , D 2 ′ are the coefficient matrices corresponding to a column of nodes in the boundary, and A 3 , B 3 , C 3 , A 3 ′ , B 3 ′ , and C 3 ′ are the coefficient matrices corresponding to the nodes of the boundary column. In the construction, the above matrix needs to be adjusted according to the actual boundary conditions. W is a (i × j)×1-order column matrix composed of unknown deflection w i,j at each node and is expressed as follows: Q is the value on the right side of the difference equation and is a known condition in actual calculation, which can be expressed as where q(i, j) is the load value per unit area at the node (i,j).
When the boundary condition is free, the virtual nodes of one line outside the boundary can be represented by two lines of nodes close to the boundary and the virtual nodes of two lines outside the boundary can be represented by three lines of nodes close to the boundary. erefore, when constructing the coe cient matrices A, B, and C corresponding to the internal i-4 column nodes, the j-4 between them acts as the conventional coe cient, which is derived from the original fourth-order di erence equation. e rst two lines and the second two lines of coe cient matrices A, B, and C are still missing the virtual nodes of the line outside the boundary and the two lines of virtual nodes outside the boundary. erefore, it is necessary to rewrite the correlation coe cient in the matrix, as shown in t 1 -t 9 and t 1 ′ -t 9 ′ in the following formula. Since the coe cients to be rewritten in dealing with simply supported boundary and xed boundary have been included in t 1 -t 9 and t 1 ′ -t 9 ′ , the coe cient matrices A, B, and C can be expressed in the following general form: After the construction of coe cient matrix K and matrix Q, the numerical solution of skin de ection under di erent boundary conditions can be obtained by using the MATLAB software, and then, the real solution of the original small de ection bending di erential equation can be approximated. We can use the left division command in MATLAB, which has better numerical stability than calculating the inverse matrix [30]. e calculation formula is as follows:

Di erence Calculation Based on Large De ection eory.
According to (17), the di erence equations based on large de ection theory are composed of di erential equations and compatibility equations of elastic surface. To solve this problem, we need to solve the unknown functions in the equations, de ection function W(x, y) and stress function ϕ(x, y), according to the actual boundary condition of the mechanical model. Since the exact solutions of multiple nonlinear di erential equations is di cult to obtain, based on the di erence calculation method, the equilibrium differential equations of large de ection theory are transformed into the following di erence equations: ere are two unknown matrices in the di erence equations, which are de ection matrix W and stress function matrix ϕ. In the original system of di erential equations, the unknown matrix ϕ is on the right side of the di erential equation of elastic surface. In order to unify the coe cients, it is moved to the left side, that is, the unknown matrix K 1 is included in ϕ. G is a matrix containing an unknown matrix W. K 1 and K 2 are the coe cient matrices to be constructed. e equations cannot be solved directly by the di erence method. e approximate solution is obtained by using the idea of gradual approximation. First, the initial value of ϕ is assumed and brought into the di erential equation of elastic surface to solve W by the di erence method (usually, the initial value of ϕ is 0, which is solved according to the small de ection theory). en, the solved W is brought into the compatible equation to solve the ϕ value by the di erence method, and then, ϕ is brought into the di erential equation of elastic surface to compare the solved W with the previous W. It is considered that the solution is completed by iterating in cycles according to this law until the W values of the previous and the following two times approach. Mathematical Problems in Engineering e approximation conditions of the two W values (W t and W t+1 ) are expressed by the Frobenius norm of the difference between the two matrices: e differential equation of elastic surface is solved. According to the difference idea, the difference relation of each basic node w i,j is represented by the difference diagram as shown in Figure 5.
It can be seen that the partial differential expression of the stress function ϕ(x, y) is needed when constructing the coefficient matrix of the difference equations. e second-order partial differential of stress function ϕ(x, y) can be expressed by calculating the second-order degree of matrix ϕ twice with gradient function ϕ in MATLAB software. At the same time, it should be noted that in the difference equation of the basic node ω i,j , the partial differential of ϕ in its coefficient must be the second step degree of one-to-one corresponding node ϕ i,j in the stress function matrix ϕ. In the coefficient matrix, the coefficients z 2 ϕ/zx 2 , z 2 ϕ/zy 2 , and z 2 ϕ/zx zy used by ω i,j in the corresponding rows of the coefficient matrix are the twostep values of ϕ i,j in x, y , and xy directions, respectively. erefore, the difference between the large deflection theory and the small deflection theory is that when the large deflection theory is used, the values of each line in the coefficient matrix are affected by ϕ and change with ϕ. In the process of solving the differential equation of elastic surface in large deflection theory, the treatment of boundary conditions is the same as that in solving the differential equation of small deflection theory. e coefficient matrix K j of elastic surface differential equation can be expressed as follows: In the matrix, A i , B i , C i , D i , and H i are the square matrices of order j, which can be given by the following equation (27).

Mathematical Problems in Engineering
Similar to the boundary treatment method based on the theory of small de ection, the corresponding coe cients of the rst two and the latter two lines of A i , B i , C i , D i , and H i should be adjusted according to the change of boundary conditions, and the parameters with b subscripts are used in the matrix. Other parameters are given in Table 1. e compatibility equation is solved. e boundary condition of stress function ϕ is treated as follows: the ϕ value of the boundary node is 0, and the ϕ value of the virtual node of a row (or column) outside the boundary is the ϕ value of the symmetric node of a row (or column) inside the boundary. e construction form of coe cient matrix can refer to the construction method of coe cient matrix in the calculation of small de ection theory. e right term of the compatibility equation also deals with the de ection matrix W by using gradient function.

Calculation Error Analysis
In this section, the "4-2-1" support layout of skin is taken as an example to evaluate the accuracy of the elastic mechanics model by using the di erence calculation method based on the small de ection theory and the one based on the large de ection theory.
Taking the thickness of the skin as 2 mm, two xtures are used, and each xture has two sections of surface to t with the skin. e concentrated force is applied on the middle part of the skin to calculate its de ection. As shown in Figure 6, the loads and constraints of the skin are shown, AB and CD are the free boundaries and AD and BC are the freexed composite boundaries. In the case of free-xed composite boundary conditions, it is di cult to solve the analytical solution of the di erential equation of thin plate bending deformation by using conventional methods. erefore, using the di erence method to deal with the complex boundary problems can fully re ect the convenience of the di erence method to deal with the complex boundary conditions. In order to evaluate the accuracy of the theoretical calculation, ABAQUS nite element simulation is used for comparative calculation. Aluminum alloy 2024-T3 and SC8R shell elements are selected as the skin material in the simulation.
First, the deformations of two kinds of the skin with dimensions of 200 × 200 mm, xed boundary length of 80 mm, and dimensions of 400 × 400 mm, xed boundary length of 100 mm under 65 N concentrated force, are compared. e di erence method based on the large deection theory and the di erence method based on the small de ection theory are used for theoretical calculation. e di erence method is 40 × 40 units. e simulation process of ABAQUS is calculated by linear analysis and nonlinear analysis with large deformation. As shown in Figures 7 and  8, the curves of the maximum de ection of two kinds of skin  Figure 6: Skin loads and constraints. Mathematical Problems in Engineering changing with the concentrated force are shown. It can be seen that when the maximum de ection is less than 1 mm (that is, half of the skin thickness), the error of the four calculation results can be controlled within 0.10 mm. With the increase of de ection, the calculation results based on the small de ection theory and ABAQUS linear analysis are approximate, and the maximum error within 2 mm deection is 5.76%. However, based on the small de ection theory and ABAQUS linear analysis, the maximum error within 2 mm de ection is 28.43%. is is because the strain and internal force caused by the longitudinal displacement of each point in the neutral plane are not taken into account in the calculation based on the small de ection theory. With the increase of deformation, the neutral plane stress and strain caused by longitudinal displacement increase, which makes the calculation based on small de ection theory deviate greatly from the actual situation. e calculated results based on the large de ection theory are similar to the results of ABAQUS nonlinear analysis, and the variation trend is close. e maximum error within 2 mm de ection is 12.70%. If it is considered that the nite element nonlinear analysis results with large deformation are closer to the real situation, the calculation result based on large de ection theory is more accurate in the case of large deformation of the skin. e de ection of ABAQUS nonlinear analysis results, calculation results based on small de ection theory, and calculation results based on large de ection theory at the red dot position shown in Figures 9-11 are compared. e skin with 400 mm size and 100 mm xed boundary length is selected to analyze the results under the action of 50 N concentrated force. Figure 9 shows the nite element simulation results, Figure 10 shows the theoretical di erence calculation results of large de ection, and Figure 11 shows the theoretical di erence calculation results of small deection. Figure 12 shows the comparison of the three calculation results. It can be seen that the de ection change trend of the three calculation methods is consistent at different node positions, which indicates that the di erence method based on thin plate theory not only has more accurate calculation results at the maximum de ection but also can accurately re ect the de ection at di erent positions of the skin. When the de ection exceeds 1 mm, the calculated results based on the large de ection theory are closer to the results of ABAQUS nonlinear analysis, which is the same as the previous conclusion. By comparison, when the actual deformation is large, the calculation results based on the small de ection theory exaggerate the deformation to a certain extent. At this time, the design and calculation according to the small de ection is safe and wasteful, which is consistent with the conclusion mentioned in the literature [12]. e calculation results based on the large de ection theory are compared under the di erent grid densities of 20×20, 30×30, 40×40, and 50×50, and the simulation results are shown in Figure 13. In addition to the large deviation of the calculation results when using the di erence grid density of 20×20, the de ection of the other three cases keeps a high consistency with the change of load, the maximum deviation is 1.3%, and the deviation is 0.03 mm, which can be ignored. However, with the increase of grid density, the computing time will be doubled. For example, under the same computing conditions, the computing time of 20×20, 30×30, 40×40, and 50×50 di erential grid density is 1.5 s, 2.9 s, 8.8 s, and 20.0 s, respectively. is is due to the increase of the unknowns in the di erence equations and the capacity of the coe cient matrix with the increase of the grid density. Moreover, under the same conditions, with the increase of de ection, the calculation time will also increase, which is caused by the gradual approximation calculation method of large de ection theory. Under the assumption of the initial value, the initial error of large and small de ection theory is large, so that the number of iterations increases. erefore, it is necessary to choose the suitable grid density in di erent cases. e following conclusions can be obtained by the error evaluation of the di erence method above.   analysis, and the results based on the large de ection theory are also similar to the results of the nite element nonlinear analysis. It is proved that the model based on thin plate theory is accurate, and the approximate result of di erence calculation is reliable. If it is considered that the results of nonlinear nite element analysis are closer to the real situation, the large de ection theory is more accurate.
(2) When the de ection is less than half of the skin thickness, the calculation results based on the large/ small de ection theory are both accurate; When the de ection is more than half of the skin thickness, the calculation results based on the large de ection theory can still approach the accurate value of the simulation, and the calculation results based on the small de ection theory will exaggerate the deformation to a certain extent. At this time, if the design and calculation are carried out according to the small de ection, it is too safe and even wasteful.  (i) e di erence method is derived by using the Taylor series expansion of the function and omitted the higher-order terms above the second order, and its error cannot be ignored. erefore, the more reserved terms and the smaller the mesh size is, the more accurate results can be obtained, but the amount of calculation will increase. (ii) In the calculation process based on the large de ection theory, the treatment of the boundary conditions of the stress function is relatively rough. In practice, the xed boundary will be subject to certain longitudinal constraints, and the stress function of the boundary is di erent from the approximate solution, which leads to errors. (iii) When using the di erence method, the di erence of grid density will also a ect the calculation results of de ection. Generally speaking, the denser the grid, the higher the calculation accuracy.

Application Examples
e application veri cation is carried out by taking the skin corner position of a certain type of aircraft as an example. In Figure 14, a and b represent the edge size of the area and c represents the length of the xture surface. Study the deection after applying drilling force at the red point F1, as shown in Figure 14. Figure 15 shows the relationship between the length of the xture surface and the maximum de ection of the skin when the two adjacent surfaces of the xture are e ective support segments (b 280 mm) under the action of 65 N drilling force. e curve of a 420 mm is added as the control group. It can be seen from Figure 15 that the maximum de ection of the skin decreases with the increase of the length of the board surface. Under the condition of a 420 mm, the maximum de ection is 1.12 mm and the minimum de ection is 0.70 mm in the range of 20-110 mm. Under the condition of a 840 mm, the maximum de ection is 1.29 mm and the minimum de ection is 0.74 mm in the range of 20-110 mm. When the length of the surface is more than 60 mm, the deection di erence of the skin is basically maintained within 0.05 mm when the distance between the two plates is di erent. At this time, the in uence of the distance between the two plates on the de ection is very small.
is shows that the skin de ection under concentrated load is mainly a ected by the boundary condition of the surrounding area. When the area is large enough, the in uence of the change of the far boundary on the skin de ection can be almost ignored. Figure 16 shows the di erence calculation results based on the large de ection theory when a 840 mm and the support length is 100 mm, and Figure 17 shows the nite element analysis results under the same conditions. Comparing the two gures, it can be seen that the de ection error of the theoretical calculation results and the nite element analysis results at the load point is 0.05 mm, and the error of the maximum de ection at the edge is 0.04 mm, and the overall de ection di usion trend is the same, indicating that the theoretical calculation results are accurate.

Conclusions
In this study, aiming at the problem of xture design in exible tooling, the load deformation problem of the skin supported by the xture is transformed into the thin plate bending problem in elastic mechanics, and the mechanical models are established based on the small de ection theory and the large de ection theory, respectively. Considering the complexity of the boundary conditions of the skin under the support of the xture, the di erence method is used to solve the theoretical model. e coe cients of the di erence equation are transformed into the form of coe cient matrix, and the approximate values of the de ection of each node of the skin are obtained by matrix operation. e di erence results based on small de ection theory and large de ection theory are compared with the nite element simulation results. Finally, the skin corner position of a certain type of aircraft is taken as an example to verify. e support area and support e ect of the xture on the skin can be calculated and optimized by the numerical calculation method proposed in the manuscript, which guarantees the skin deformation within the allowable range when reducing the e ective support section of the xture to the skin. It has important theoretical guiding signi cance for the exible design of xture.
Data Availability e data, models, and code generated or used to support this study are included within the article.

Conflicts of Interest
e authors declare that they have no con icts of interest.