Initial Investigation of Energy Finite Element Validation on High-Frequency Flexural Vibration of Stiffened Thin Orthotropic Plates

Energy finite element analysis (EFEA) has unique advantages in solving high-frequency dynamic responses of orthotropic structures, due to its ability to obtain detailed local response information. In order to accurately predict high-frequency vibration response of the stiffened orthotropic plate, EFEA theory on the propagation of bending wave in the orthotropic structure and the energy transfer coefficient which express the energy transfer at the stiffener was investigated. Based on the EFEA theory presented, high-frequency dynamic responses of a stiffened orthotropic plates were predicted. Furthermore, tests were done for the same problem, and differences between prediction and test were discussed. Finally, the future works were pointed out.


Introduction
In high-speed flight of aircraft, high-frequency vibration of fuselage structures can influence the reliability of airborne electronic equipment.Considering this, the prediction of high-frequency dynamic responses of aircraft fuselage structure is one of important contents in the process of designing aircraft.
EFEA is a developing method for predicting highfrequency dynamic responses of a structure.In the EFEA, the governing differential equations are formulated for an energy variable that is spatially averaged over a wavelength and time averaged over a period.Differential equations are derived for all wave bearing domains within a system, and it represents a power balance over a control volume [1].The geometric characteristics and nonuniform damping characteristics of a structure can be completely expressed and the characteristics of distributed loads are strictly described.In this way, the energy response information of all interested points in the structure could be obtained.Compared with conventional finite element methods for predicting dynamic responses of structures, EFEA shows the advantages including simple model and less calculation in predicting high-frequency responses.In addition, compared with the statistical energy analysis, EFEA can obtain the dynamic response at every interested point in the structures.Therefore, this method is suitable for predicting the high-frequency dynamic responses of composite materials with orthotropic characteristics.
EFEA is originated from the study on energy flows.In the early studies of energy flows, it is assumed that the energy within different wave forms is irrelevant and can be superposed.In the light of this assumption, the energy flow equation for bars was deduced by using heat conduction equation.Based on this equation, Nefske and Sung proposed to solve energy flow equation by utilizing finite element method in 1989 [2], and, since then, the EFEA has been known in the world.Early researches on EFEA mainly focus on the energy propagation in bars and beams.In 1992, Bouthier and Bernhard introduced the approximation assumption of lossy plane waves in far fields and small damping assumption of structures [3] and studied the governing differential equation of energy density of plates and membranes.In this way, the application of EFEA was extended from bars and beams to plates and membranes.Afterwards, many scholars pay attention to EFEA and its theoretical and application fields are researched and developed constantly.Since 2008, some scholars have attempted to apply the EFEA to orthotropic structures.Yan developed EFEA to study the energy propagation in laminated composite materials [4].Lee applied EFEA to analyze the vibration responses of aircraft rotor made of multilayered laminated composite plates [5].In addition, Cai et al. studied the application of the EFEA in coupled composite laminated beam structures [6].
Since the energy density, variable of governing equation of EFEA, is discontinuous in the coupled joints, the overall matrix equation needs to be specially processed if the structures with coupled joints are considered.The coupled joints exist in location of a change of geometry, a change in material properties, multiple components being connected together, or different media interfacing with each other.The energy transfer relationship in the coupled joints is needed in EFEA.In 1993, Cho put forward an idea of employing the EFEA to a coupled vibration system including plates, membranes, and acoustic fields and studied the characteristics at coupled joints of various coupled structures [7].In 1999, Bernhard and Huff used the EFEA to research the couple problems of vibration and acoustic fields in chamber structures [8].In addition, by using the EFEA, Yan explored the propagation of energy in linearly coupled composite laminated structures [4].However, stiffened plate structures widely used in engineering are rarely studied deeply.
Orthotropic composite materials are the simplest and most typical composite materials.Therefore, a stiffened orthotropic plate was used as the research object in this investigation; the problems in predicting high-frequency vibration responses of the stiffened thin orthotropic plate by using the EFEA are theoretically analyzed.Based on wave theory, the energy transfer coefficient at the stiffener of the stiffened orthotropic plate is deduced.Furthermore, the energy transfer coefficient is integrated into the energy finite element matrix.The equation is used for programming and calculation, and the predicted results are discussed and compared with the experimental results.

EFEA Theory for Predicting the Response of High-Frequency Flexural Vibration of the Stiffened Plate
There are two typical processes in EFEA: deriving governing differential equation in terms of energy density variables and employing finite element approach to solve those equations numerically.
As to the thin orthotropic plate, the plane in the thickness direction is the isotropic plane.As shown in Figure 1, the  plane with  = 0 in the  coordinate system is used as the middle plane and the displacements in , , and  directions are , V, and , respectively.
To derive the EFEA governing equation of an orthotropic plate, the assumptions have been applied.The assumptions of thin plate theory and the invariant straight-line method are as follows.
(1) The bending deformation of the thin plate occurs within a small deflection range.The straight line that is initially vertical to the middle plane of the plate before deformation remains a straight line and vertical to the middle plane, and its length does not change, after the plate is stretched and bent.
(2) The thickness of the thin plate does not change and the stress   in  direction is ignored, which is generally acceptable in-plane stress state analysis.
On the premise of aforementioned assumptions, the general equation of motion for the bending vibrations of an orthotropic plate is as follows [4]: where When a plane wave is considered, the general far-field solution for the movement equation is where  is a constant related to the travelling wave amplitude and   and   are complex wave numbers related to the natural frequency of vibration damping in the and directions, respectively.
where  is hysteresis damping loss factor.Assuming that  is the total wave number in the plate and   ,   are the and -direction components of total wave number, far-field solution is substituted into the equation of motion; then the following scattering relationship is obtained.
where  is the wave propagation angle.Then, ( 5) and ( 6) are substituted into scattering relationship formula (4) to obtain the expression for the relationship between the total wave number and wave deflection angle, as follows: Assuming that there is a scattering wave field in the plate, the integral of the wave number from 0 to 2 is calculated and the average wave number can be represented as Then, the average complex wave numbers  *  and  *  can be described by the averaged wave number, as follows: The energy density  is defined as the energy per unit volume in a given system.The intensities  are defined as the energy that passes through a given area of the medium in unit time.The total energy density is the sum of the kinetic energy density and potential energy density.The time averaged total energy density ⟨⟩ can be represented via displacement as follows: where superscript * indicates the complex conjugate.
Then, and -components of the time averaged energy intensities ⟨  ⟩ and ⟨  ⟩ are represented via the shear forces   ,   , bending moments   ,   , and twisting moments   ,   , as follows: Bending moments, twisting moments, and shear forces are represented via displacement as follows: ) . ( Given that structural damping loss factor  is much less than 1, the second-order items and high-order items of the damping loss factor can be ignored.The far-field displacement solution is substituted into the time averaged energy density and intensities expressions, and a spatial averaging operation is applied.After simplifying the arrangement, the expressions of time and spatial averaged energy density and intensities are as follows: The time and spatial averaged energy density and intensities expressions show that the -component and component of the two parts of the average intensities are proportional to the derivatives of the two parts of the average energy density over  and , respectively. where Considering the energy conservation relationship at a stable state, the energy balance equation is The relationship between the dissipation power and energy density is Then the governing differential equation for the bending wave of an orthotropic composite material plate is To obtain the numerical solution to EEEA governing differential equation (20), the matrix form of the differential equation should be obtained via the variation method.The linear equation is where Γ is the element boundary,  is the area domain of the element, and   is the shape function.The first term of above equation is essentially ∫ Γ    ⋅  Γ. is the normal direction of element boundary.
For the coupled structure, due to the fact that the energy density is discontinuous in the joint, the power transmission relationship is needed to calculate the first term in (21).When the bending wave only occurred in plates, power transmission relationship can be expressed by reflection coefficient  and transmission coefficient .Therefore, reflection coefficient  and transmission coefficient  will be derived in the next section.

Calculation of Energy Transfer Coefficient in Joint.
For the coupled structure, the core step of calculating responses using EFEA is to solve the energy transmission coefficient and reflection coefficient in the joint.
Figure 2 illustrates a stiffened orthotropic plate.The plate is divided into two parts, namely, plates 1 and 2, by the stiffener.The axial direction of the stiffener is defined as axes, and the height direction of the stiffener is defined as axes.
According to the balance of force and moment in the stiffener and joints of the plates, the relationships of force and moment in the stiffener with force and moment in the joints of the two plates are where  = 2,  = (  ,   )  represents the force and moment per unit length in the stiffener and  indicates the tension  = (    )  acted on the joint of the th plate, respectively.
1 = ( 1 0 0 1 ),  2 = ( −1 0 0 1 )  is related to the displacement of the stiffener.Considering the shear deformation and moment of inertia, the governing equation of motion of the stiffener can be expressed as follows [9]: where   ,   , and   represent the shear elastic modulus, the mass density, and the cross-sectional area of the stiffener, respectively;   and   denote the second-order moment of and the shape factor of the cross section, respectively;   and   indicate the force and moment per unit length on the stiffener caused by the connected plates separately.
Thus, we can get the relationship of  and the displacement  = (, )  of the stiffener: From ( 23), the component of  is shown as follows: 12 =  21 = −     .
According to local coordinate system, the dynamic governing equation of the deformation of the th plate is where  and   represent the displacement in  direction and the mass density per unit area, respectively.The relationship between the displacement and tension in the joint ( = 0) of plates is where  and  denote the shear force in  direction and the bending moment in the joint of every plate, respectively.
In the following, it is considered that the vibration waves are plane waves which propagate from a semi-infinite plate to the plate joint, and some are reflected and the others are transmitted in the joint.Supposing that the incident wave is in the form of exp(−++), the responses in all plates should be the same dependence relationship exp(− + ) on , while the dependence on  is determined by the equation of motion of plates according to Snell's law.Assume that the displacement in the form of exp(− +    + ),   can be expressed as where where  1 and  2 represent two effective roots of (27), while  1 and  2 indicate the complex amplitudes relating to the two roots, respectively.The rotation angle   can be expressed as According to the above two equations, the displacement and rotation in the joints of plate  can be calculated through  = 0 and they can be expressed in accordance with  1 and  2 as well as  1 and  2 : The tensions   and   in the boundary can be expressed according to  1 and  2 : Based on (30),  1 and  2 can be expressed according to   and   as By substituting formula (32) into formula (31), the relationship between displacement and tension in the joint of plates is obtained: The above equation is rewritten in the matrix form as follows: where   = (    )  and   = (    )  .
If the incident waves propagate in plate 1 and amplitude is  in , the relationship between force and displacement of plate 1 needs to be modified as follows: 1 and   1 represent the boundary displacement and force induced by incident waves, respectively.The vector   1 can be expressed as follows: where  =   sin  and  =   cos .
By substituting the displacement expression into the expressions of bending moment and rotation angle,   1 is obtained, that is, Substitute ( 24), (34), and (36) into (22); then The compatibility condition between two plates and the stiffener is that the displacement   = (  ,   )  of the joint in the two plates equals the displacement  = (, )  of the stiffener: Substitute (40) to (39); then (39) is changed as follows: By solving the above formula, the displacement  of the stiffener can be obtained; then the joint displacement of each plate can be obtained.
For a wave with an amplitude, radial frequency, and the angle between incident direction and joint being   , , and , the energy transferred to the joint per unit length can be written as 2  sin , where   and  represent the group velocity and the average kinetic energy per unit area caused by fluctuation.Then The calculation formula for transmission coefficient is The transfer coefficient  denotes the ratio of power transferred to plate 2 to incident power, and  denotes the ratio of power reflected to plate 1 to incident power.

Prediction of High-Frequency Flexural Vibration Response of a Stiffened Plate
With the aid of EFEA, the high-frequency vibration response of a square orthotropic stiffened thin plate is calculated.The plate is prepared by pasting an aluminum stiffener with rectangular sections on a square orthotropic thin plate.As shown in Figure 3, the dimensions of composite plate and stiffener are 2 mm × 1000 mm × 1000 mm and 2.5 mm × 50 mm × 1000 mm, respectively.The bonding position and excitation points are demonstrated in 4.
Young's moduli of the square plate in  and  directions are   = 25 Gpa and   = 0.625 Gpa, while shear modulus and primary Poisson's ratio are   = 0.375 Gpa and ]  = 0.25, respectively.The plate density and thickness are  = 1780 kg/m 3 and 0.002 m, respectively.As for the boundary condition, this square stiffened plate has four free edges.The grid model of the square stiffened plate is shown in Figure 5. Owing to the fact that the energy density value is discontinuous in the stiffener, double nodes are set in the stiffener.The entire grid model consists of 100 quadrilateral elements and 132 nodes.
The hysteresis damping factor of bending waves is 0.08 and the excitation of 2,000 Hz is applied at point .The boundary condition is set that there are no energy flows to outflow in the four sides of the square plate.Figure 6 shows the distribution of energy density on the square plate predicted by EFEA when point  is excited.
As shown in Figure 6, due to the existence of the stiffener, the vibration energy is mainly concentrated in the area of the vicinity of the excitation point.Particularly, if the excitation point is far away from the stiffener, very little energy passes through the stiffener; thus the stiffener exerts little influence on energy density of the whole stiffened plate.

Test Investigation
In order to verify the correctness of predicting high-frequency dynamic responses of the stiffened thin plate using EFEA, the same stiffened thin plate was analyzed through experiments.The stiffened thin plate is achieved by gluing the stiffener to the orthotropic plate using epoxy resin adhesive, as shown in Figure 7.The orthotropic plate is made of carbon fiber reinforced composite; the weaved fiber direction is 0 ∘ and 5 ∘ , as shown in Figure 3.The stiffener is made of aluminum.The cross section of the stiffener is 0.25 cm and 5 cm in width and height, respectively.In the experiments, the plate is hanged up by four elastic ropes attached to four corners of the plate, to make the plate perpendicular to the horizon.The experimental scene is shown in Figures 8 and   9.The single-frequency sinusoidal excitation was applied in excitation point b by utilizing the exciter.By measuring the vibration response of the stiffened plate using the laser vibrometer, the square of amplitude of vibration velocity at every measurement point in plate can be gotten; then the spatial distribution of vibration energy on the plate can be observed.
When the excitation points is point b, the energy distribution in the plate is obtained through data measurement and analysis; comparison of test results and prediction results were demonstrated in Figure 10.
From Figure 10, it can be seen that vibration energy mainly focuses on the areas near to the excitation point.Because the variant of EFEA is spatial averaged energy density, there are differences between the results of EFEA and test.So we can compare the range of energy distribution to conclude if EFEA are in agreement with test.After a comparison, it can be seen that the vibration energy mostly is distributed on the right side of stiffener; the excitation point b is just in the right side of stiffener.The test results are found to be basically consistent with EFEA.
For the reason of comparison, experiment is also conducted by breaking the symmetry of point  in vertical direction and using a plate without stiffener.The predicted  energy density distribution is shown in Figure 11.It can been seen from the figure that the energy can be transmitted into whole plate when the stiffener is removed and effect of the excitation location change in vertical direction could be clearly identified.
In this investigation, difference between prediction and test may be induced by the following factors: (1) in deriving the formulation of EFEA, shear, in-plane force, torsion of the stiffener, and ambient boundary condition of the plate have been ignored; (2) in EFEA prediction, diffusion wave filed assumption is used which only shows omnidirectional wave, while the edge reflection is presence in the test; (3) the variable of EFEA is time averaged and spatially averaged energy density, while result from the test is that the energy distribution on the plate is excited by single frequency; (4) the thin aluminum plate stiffener is not rigidly and accurately attached to the orthotropic plate with only epoxy glue along its side edges.

Conclusions
In this paper, the high-frequency dynamic response of a stiffened thin orthotropic plate is predicted by EFEA.In order to get the energy density distribution of the stiffened plate, the EFEA governing differential equation has been derived, and the transmission and reflection coefficients have been derived to indicate how wave energy is distributed when the bending wave propagates through the stiffener.Then the matrix of energy transmission coefficients is integrated into the governing equation matrix of energy density of the whole thin plate to obtain the distribution of flexural vibration energy density on the stiffened thin plate.In order to verify the correctness of the results, the vibration response of the same plate in single-frequency excitation is measured by using the laser vibrometer to obtain the distribution of vibration energy on the thin plate.The main conclusions are as follows: (1) the prediction results of the EFEA are roughly consistent with the experimental results.The differences in response between test and EFEA are discussed and the possible reasons responsible for the discrepancy are analyzed.
(2) Stiffener greatly influences the distribution of energy density on the plate.When the excitation points are far away from the stiffener, little energy passes through the stiffener.The existence of the stiffener blocks the transmission of the vibration energy and energy is mainly concentrated in the areas where the excitation points are located.
This paper is only an initial investigation on the EFEA of stiffened thin orthotropic plates.There are still many works that should or could be conducted in the future.For example, comparison of predicted and measured results should be made at various frequencies and, as a function of frequency, the problem considered could be simplified by using an isotropic plate with different sized stiffeners to investigate

Figure 2 :
Figure 2: Schematics of the stiffener and plate.

Figure 5 :
Figure 5: Grid model of the square stiffened plate.

Figure 6 :
Figure 6: Energy density distribution predicted by EFEA while the plate is excited at point .

Figure 7 :
Figure 7: The detailed picture of stiffened orthotropic plate.

Figure
Figure 8: Position relationship between the exciter and the stiffened plate.

8 :
Figure 8: Position relationship between the exciter and the stiffened plate.

Figure 9 :
Figure 9: Position relationship between the laser vibrometer and the stiffened plate.
Figure 10: Comparison of test results and prediction results.

Figure 11 :
Figure 11: The energy density distribution of a plate without stiffener.
is displacement in the  direction;  is exciting force;  is radian frequency; ( −  0 )( −  0 ) is delta function;   is mass density per unit area;  is time;   is bending rigidity which describes the relationship among bending curvature, torsional curvature, and internal flexural force;   and   are the longitudinal and horizontal elastic moduli, respectively;   is the shear modulus;   and   are the primary Poisson ratio and secondary Poisson ratio, respectively; and ℎ is thickness of the plate.