The Application Research of Inverse Finite Element Method for Frame Deformation Estimation

A frame deformation estimation algorithm is investigated for the purpose of real-time control and health monitoring of flexible lightweight aerospace structures. The inverse finite element method (iFEM) for beam deformation estimation was recently proposed by Gherlone and his collaborators. The methodology uses a least squares principle involving section strains of Timoshenko theory for stretching, torsion, bending, and transverse shearing. The proposed methodology is based on staindisplacement relations only, without invoking force equilibrium. Thus, the displacement fields can be reconstructed without the knowledge of structural mode shapes, material properties, and applied loading. In this paper, the number of the locations where the section strains are evaluated in the iFEM is discussed firstly, and the algorithm is subsequently investigated through a simple supplied beam and an experimental aluminum wing-like frame model in the loading case of end-node force. The estimation results from the iFEM are compared with reference displacements from optical measurement and computational analysis, and the accuracy of the algorithm estimation is quantified by the root-mean-square error and percentage difference error.


Introduction
Aircraft flexible wings with embedded conformal antennas and large frame structures that carry antennas require accurate real-time deformation estimation to provide feedback for their actuation and control systems [1][2][3].Using the measured strain to reconstruct the shape of the structure is a key technology in the accurate real-time deformation estimation, which has been studied by many researchers.
The computation of the displacement field of the deformed structure is commonly performed on the basis of strain data measured in real time by a network of strain gauges [4][5][6].For example, fiber Bragg grating (FBG) sensors have been extensively researched for deformation estimation due to their lightness, accuracy, and easy embedment.The strategies that reconstruct the deformed shape displacement field of the structure with in situ strain data can be divided into two kinds.One of the two kinds trains the mapping relation between the displacement field and the measured strain by model learning algorithms, for example, neural network model and fuzzy network algorithm [7,8].When the training system has enough measured strains and displacement field data of the structure, a certain relation matrix can be determined and the stable relationship between the measured strain and the displacement field can be obtained.But the strategy requires a large number of training data, and the mapping relation is easy to fail when the actual loading is beyond the range of the training cases.
The other establishes the mapping relation between the displacement field and in situ strain data without the model learning.In the literature [9][10][11], the global or piecewise continuous basis function methods were employed to fit the surface-measured strain into the structure strain field, and then, the structure deformation displacement was obtained from the strain-displacement relationship.These methods are easy to implement, but the reconstruction accuracy of deformation estimation depends on the appropriate selection of basis function and weight coefficients.Mode shapes have been used as basis function in [12,13].The deformation displacements are reconstructed from measured strains by using the modal transformation method.However, there exist the following disadvantages in this method.(1) The detailed material elasticity and inertial parameters are needed to precisely construct mode shapes.(2) The accuracy of deformation reconstruction is severely limited to the modeling precision of the structure, and it is quite difficult to precisely model the complex structure.Maincon [14] developed a finite element-based methodology involving an inverse interpolation formulation that employs the surfacemeasured strain to determine the loads and structural response of aerospace vehicles, while this algorithm needs an appropriate quality function to adapt the different loading cases and this function is constructed based on a mass of computer simulation and experimental statistics.
Based on the Euler-Bernoulli beam equation, Jute et al. evaluated the deflection of beam by the integration of discretely measured strains directly [15].The algorithm applied the classical beam equation and piecewise continuous polynomials to approximate beam curvature through integration.Derkevorkian et al. [16] compared the algorithm with the modal transformation method and demonstrated the additional benefits of the algorithm to achieve a robust method for monitoring ultralightweight flying wings or next-generation commercial airplanes.Though this onedimensional scheme has displayed high accuracy in predicting deflection, it fails to estimate the element deformation under multidimensional complex loads.
Tessler and Spangler [17] proposed the inverse finite methodology (iFEM), which can be used to reconstruct the displacement field of shear-deformable structures, not only beam structure but also plate and shell structures.The main idea is reconstructing a three-dimensional displacement field of beam structure from the surface-measured strains according to a least squares approach.Due to the fact that only the displacement-strain relationship is used, the deformation reconstruction can be accomplished by the methodology without the prior knowledge of loads, materials, and inertial and damping properties.To model arbitrary plate and shell structures [18,19], Tessler and his partners developed the iFEM algorithm using the first-order shear deformation theory and a three-node inverse shell element.FBG sensors were applied to measure the surface strains on the slender beams, and then, the deformed displacement was reconstructed by using an iFEM shell model.
The beam-deformed displacement and cross-sectional torsion were reconstructed by Gherlone et al. who employed the inverse finite element formulation to achieve high reconstruction accuracy of deformed displacement [20][21].The authors first used the Timoshenko beam theory to model the beam kinematic accurately, then used the C 0 or C 1 inverse frame elements and least squares formulation to establish the relationship between the measured strain and the arbitrary node displacement field of the beam element with no prior knowledge of the finite element model and loads.The reconstruction equation is nonsingular when the boundary conditions are applied; that is, the status of one end node of the beam element must be known.As the displacement field of the beam element in iFEM is constructed based on an isotropic straight beam structure, whose cross section is invariable along the whole beam, the deformation of the tapered beam cannot be estimated.Meanwhile, in view of different load cases, the authors discussed different types of shape functions that are used to interpolate the kinematic variables in the beam element; then, the number of surface strain measurements and the deformation reconstruction equation are confirmed.Unfortunately, the minimum number of locations where the section strains are evaluated, which is critical for obtaining a correct solution to the algorithm, is not discussed clearly in different load cases.If the minimum number of the locations where the section strains are evaluated is inappropriately set, the solution will be nonunique; that is, the reconstructed displacement may be incorrect.
The contribution of this paper is twofold.(i) The relationship between the loading case and the number of the section strain locations is discussed in detail, and the following conclusions are verified; on the one hand, the minimum number of locations should be 2 in the loading case of endnode forces; on the other hand, the minimum number should be 3 under the uniform distributed loading.(ii) In order to examine the minimum number which is 2, estimation studies are carried on a beam structure and a wing-like threedimensional frame structure, in the loading case of endnode force.Six fiber optic strain sensors are placed at two nodes of each beam to capture in situ strain dates.To assess the estimation error of the iFEM algorithm, the RMS and difference percent parameters are employed.

iFEM Algorithms for Beam Estimation
A beam deflection estimation algorithm was developed by Gherlone and his collaborators [20].In the algorithm, the Timoshenko beam theory is first used to analyze the expression of the displacement field in a straight beam element (Figure 1).
where u x , v y , and w z are the point displacements along the x, y, and z axes, respectively.u x , v x , and w x denote the displacements at y = z = 0; θ x x , θ y x , and θ z x are the rotations about the three coordinate axes.The six kinematic variables in the middle axes can be grouped in vector form as follows: T can be obtained by (1).
The six kinematic variables u can be interpolated by the right shape functions where N x and u e denote the shape function and nodal degrees of freedom, respectively.Substituting ( 4) into (3) gives arbitrary section strains in terms of the nodal degrees of freedom as follows: where the matrix x contains the derivatives of the shape functions N x .Once the section strains e u are obtained, the nodal displacement u e is determined; then, the kinematic variables u can be acquired by ( 4).However, the section strains e u are derived from the kinematic variables u theoretically, rather than the strain measurements.So, iFEM uses in situ section strains e ε computed from the measured strains to replace e u when the least squares error function φ u reaches the minimum.In view of the effect of the axial stretching, bending, twisting, and transverse shearing, the improved least squares error functional φ e u is obtained by the dot product of the weighting coefficient vector W and the original vector where w 0 k k = 1, 2, …, 6 denote dimensionless weighting coefficients whose initial values are identically set as 1; A e , I e y and I e z , and I e p are, respectively, the cross-sectional area, second moments of the area according to the y-and z-axes, and polar moment of the area of the beam element.L is the length of the beam element; x i 0 ≤ x i ≤ L and n are, respectively, the axial coordinate of the locations where the section strains are evaluated and the number of locations, that is, the axial coordinate of sections where section strains are distributed in, and the number of sections.
For a straight beam member of constant circle cross section, the in situ section strains at x = x i , e ε x i = e ε 1 x i , e ε 2 x i , …, e ε 6 x i T , could be derived from the measured strains of the beam surface with appropriate straintensor transformations from the θ, x, r to x X , x Y , x Z coordinates [22].
with c β ≡ cos β, s β ≡ sin β, c θ ≡ cos θ, and s θ ≡ sin θ In (8), v is the Poisson ratio, ε * 2 x i , θ, β is the in situ strains that are obtained from strain sensors.R denotes the radius of the beam cross section (Figure 2).Substituting ( 5) into (7) results in the following quadratic form: Herein, c e is a constant, and k e and f e are indicated as follows: Finally, the relationship between the deformation and in situ section strains, shown in (11), can be confirmed when the minimization of functional φ e u is performed.
Once the appropriate shape functions and the problemdependent displacement boundary conditions (e.g., setting the displacement of one end point to zero, which means that one of end nodes of the beam is fixed) are given, u e can be derived from a nonsingular system, and the vector f e depends on the measured strain values that change during 3 International Journal of Aerospace Engineering deformation.Once the nodal degrees of freedom u e are confirmed, the displacements and rotations of every node along the centroid axis of the beam element are obtained by ( 4) and the deformed shape of the whole beam can be reconstructed by (1).

Minimum Number of the Section Strain Location Selection
Although the displacement shape function and the number of surface strain measurements in different loading cases have been discussed in [20,21], n, the number of the sections where the section strains are distributed in, is not shown clearly.In (11), once the problem-dependent displacement boundary conditions are determined, the nonsingularity of the equation will depend on the selection of n.In our investigation, it is found that n can be determined by a specific loading case.More specifically, the orders of the section forces and the moments can be obtained by the loading case, which determines the order of the section strains.With the value of the section strain order, n is determined immediately.With the equilibrium equations ( 12), the section forces (N, Q y , and Q z ) and moments (M x , M y , and M z ) can be obtained by the load concentrations q x x , q y x , and q z x along the x, y, and z directions (Figure 3).
The relationship between section forces and moments and the section strains e i i = 1, 2, …, 6 can be interpreted as the following constitutive equations: where A x ≡ EA is the axial rigidity; G y ≡ k 2 y GA and G z ≡ k 2 z GA are the shear rigidities, with k 2 y and k 2 z denoting the shear correction factors; and J x ≡ GI p and D y ≡ EI y and D z ≡ EI z are, respectively, the torsional rigidity and the bending rigidities.For the uniform section beam element, parameters mentioned are constant.Then, the order of the section strains e i will be identical to that of the section forces and moments.
Substituting ( 13) into (12) gives the relationship between the section strains and the loads as follows: In most cases, the load in the x direction is zero, that is, q x x = 0; then, The section strains e 1 and e 6 are constant in ( 14) and (15), and the order of the residual section strains, e 2 , e 3 , e 4 , and e 5 , will be discussed as the following two loading cases.
Consider a beam element loaded by the end-node forces, q y x and q z x are zero along the x-axis.Then, e 2 , e 3 , e 4 , and e 5 are determined by (14).
where a 1 , a 2 , a 3 , a 4 , b 1 , and b 2 are unknown constant parameters.
As in (15), the highest order of the section strains e u is linear e 2 , e 3 , which means that the distribution of the bending moments is a skew line, that the corresponding errors in (7) can be obtained from the shaded area between two skew lines, and that each line can be confirmed by two different nodes, that is, two section strains in different sections 4 International Journal of Aerospace Engineering (Figure 4(a)).Thus, n = 2 for e 2 and e 3 in ( 7).Meanwhile, the other section strains, e 1 , e 4 , e 5 , and e 6 , are constant, which means that their distributions are lines parallel to the x-axis, that the associative errors in ( 7) can be obtained from the shaded area between two parallel lines, and that every line can be confirmed by one node, that is, one section strain (Figure 4(b)).Thus, the number of the sections where section strains are distributed is n = 1.Finally, the number of the sections, n, is 2 in ( 7) and ( 10) in the loading case of end-node force and moments.
Another case is the beam element loaded by uniformly distributed forces, q y x and q z x are constant along the x-axis, and the order of the section strains e 2 , e 3 , e 4 , and e 5 can be deduced as follows. Assumption: Then, Herein, d i i = 1, …, 8 are unknown constant parameters.The highest order of the section strains e u is quadratic e 2 , e 3 , and then, the distribution of the bending moments is a parabola.The corresponding errors in (7) be obtained from the shaded area between two parabolas which are confirmed by three different nodes (Figure 4(c)), that is, three section strains in three different sections.Thus, n = 3 for e 2 and e 3 in (7).Meanwhile, e 4 and e 5 are linear, where the distribution is a skew line, and e 1 and e 6 are constant.Similar to the loading case of end-node force, the corresponding numbers of sections are 1 and 2, respectively.Finally, the number of the sections, n, is 3 in ( 7) and (10) in the loading case of uniformly distributed forces.

Verification
A simple cantilevered solid beam and a three-dimensional frame structure were subjected to the end-node static loads to assess the iFEM potential for the flexible wing deformation estimation.The beam and frame structures were made of 6061-T6 aluminum alloy.The Young's modulus is E = 73000 Mpa, the Poisson ratio is v = 0 3, and its density is ρ = 2712 63kg/m 3 .The frame was composed of three solid beams and middle plates; the span of each beam L = 640 mm and the radius R = 10 mm (see Figure 5).For the solid circular cross section, the shear correction factors are k 2 y = k 2 z = 0 887 [21].In our verification, every whole solid beam is regarded as one beam element.Accordingly, the principle axis of frame structure is divided into three elements.
The experimentally measured surface strains are obtained by fiber optic strain sensors.Displacement measurements are captured by 3D optical measurement instruction (NDI, Optotrak Certus), which determines the position of the identification point by using three CCD cameras to capture the infrared lights emitted by position sensors (Figure 6).The instruction is also used to assess the iFEMrecovered deflections, where the accuracy of the NDI is 0.1 mm.A position sensor was placed as close as possible to   For the frame structure, each beam is regarded as an element and the optical fiber sensor location scheme used for each element is identical.For the end-node static loads, the displacement field of the whole beam element is interpolated by C 0 continuous shape function and the number of required strain sensors is 6 (see [21]).As the radius of the beam element is small (R = 10 mm) and the grid length of every fiber grating sensor is 10 mm, it is difficult to stick the six strain gauges on one section; six fiber optic strain gauges are placed on two different sections along the beam (Table 1).
The accuracy of the reconstitution is assessed by rootmean-square (RMS) and percentage difference (%Diff).The superscript "iFEM" refers to the reconstitution by iFEM while "NDI" refers to the experimental measure from 3D optical measurement instruction; x i is the ith location along the axis where the displacement u is measured.
At first, five different vertical force cases (Table 2) were applied at the free end of the beam (Figure 5(a)).Figure 7 shows the comparison between iFEM reconstruction and NDI measurement along the beam axis in force direction, for the loading case 68 N. Symbol x/L indicates the location of identification point along the beam surface.Zero means identification point at the clamped end node, and 1 means identification point at the free end node.v(y) indicates the corresponding displacement in force direction.
For the frame structure, five different vertical force cases (Table 3) were applied at the free end of the frame (Figure 5(c)).The two deformation contour plots of the frame are plotted by high-fidelity direct FE analysis (ANSYS 12.0), which shows the comparison between iFEM reconstruction and NDI measurement along the beam axis in force direction, for the loading case 45.5 N (Figure 8).The deformation unit in the two figures is meter.
It is seen from the preceding results that the iFEM methodology shows a good potential as a reliable estimation technique.The estimation results gained from the iFEM algorithm have a better approximation to the NDI measurement results.For the beam test, the accuracy of the reconstitution assessed by percentage difference does not exceed 6.4% and the error brings down with the loading increase.For the frame test in the loading case of 18 N, the percentage difference is 1.3%, but the difference grows up with the loading increase; especially from the loading case 23 N to 35 N, the difference increases 4%.The cause for this phenomenon is that the gripper of the frame structure is not very stable, which leads to the fact that the gap spacing between the frame and its carrier grows up with the loading increase that has a great influence on the boundary conditions.

Conclusion
This study investigates the application potential of the iFEM algorithm for the flexible wing and other frame structures.This method employs the in situ strain measurements and proper displacement shape function to estimate the deformation shape of the beam and frame structure, without the need to know the applied loads and material properties and to use modal shapes.When the conditions of the reconstruction equation and the shape function are determined, the deformation shape will be accurately estimated in a nonsingular system.Although the boundary conditions are known in hypothesis, the number of location where the section strains are evaluated is not shown clearly; thereby, the singularity of the equation will be influenced.
This paper discusses the minimum number of location where the section strains are evaluated with the loading case, using equilibrium equations in detail, and verifies the result that for the loading case of end-node force, the number is 2. The results of the static loads in the beam and the frame model show that the iFEM algorithm has a good potential for the flexible frame structure estimation; especially for a simply supported beam, the accuracy of the estimation by

Figure 2 :
Figure 2: Location and coordinates of a strain gauge placed on the beam external surface.

Figure 3 :
Figure 3: Beam section forces and moments.

Figure 5 :
Figure 5: Loading case and the distributions of NDI identification points.

Figure 4 :
Figure4: Distribution and error of section strains in different orders.e ε is a section strain in e ε computed from surface strain measurements by(8).e is a section strain in e u deduced by(3).

2 2 ,
%Dif f τ = 100 * RMS max τ NDI x i , 18 where τ = v, w is the deformation displacement along the y-or z-axis and m is the number of deflection shape displacement measurement identification points of NDI (Figure 5(b)A) in the structure.For the beam, m = 6 (Figure 5(a)), and for the frame, m = 14 (Figure 5(c)).

Figure 7 :
Figure 7: Comparison of iFEM reconstruction for a beam to NDI measurement in Y.

Table 2 :
Loading case for the clamped-free beam.Max deflections are captured from NDI.

Table 1 :
Optical fiber sensor location in one element.

Table 3 :
Loading case for the frame structure.Max deflections are captured from NDI.