Active Vibration Control of Cantilever Structures by Integrating the Closed Loop Control Action into Transient Solution of Finite Element Model and an Application to Aircraft Wing

1. Introduction

Flexible cantilever structures are widely used in engineering applications such as civil engineering structures, construction machines, robot manipulators, aircraft wings, satellite antennas, solar panels, and marine risers. Flexible structures undergo higher vibrations and thus vibration control attracted the attention of many researchers [1]. Wang et al. [1] studied a cantilever beam with end mass and designed adaptive boundary controller to suppress vibration. They used partial differential equation of the beam. Lyu et al. [2] proposed a vibration control for cantilever beams which uses a combined macro-fiber composite (MFC) and rotary motor actuators. They modeled the system using Euler-Bernoulli beam theory using orthogonal polynomials as mode shapes and obtained superior vibration control with combined actuation as compared with using actuators alone. Li et al. [3] used a hybrid PID-FxLMS algorithm for active vibration control (AVC) of a cantilever beam with a piezoelectric stack actuator. They modeled the system in ABAQUS to obtain state space model for simulation and verified the simulation results with experiments. They observed that vibration suppression performance of the hybrid control is much better than that of the FxLMS algorithm or traditional PID controller alone. Amer et al. [4] considered the linear velocity feedback control, the cubic velocity feedback control, the non-linear saturation controller and the positive position feedback (PPF) controllers to compare and they observed that PPF controller is the best controller for AVC of cantilever beams. Saif et al. [5] explored the direct application of deep learning technique and the approach which integrates deep learning into a PID controller for AVC of cantilever beams. They compared the two techniques to evaluate effectiveness and strengths and offered their potential applications. Huang et al. [6] studied AVC of composite cantilever beams with the finite element (FE) method and Golla–Hughes–McTavish model. They optimized the linear quadratic feedback gain using particle swarm optimization optimization. Djokoto et al. [7] studied a micro-cantilever beam (MCB) with electrorheological fluids (ERF) support. They presented results showing the vibration controlling of MCB by activating ERF with an electric field. Cui et al. [8] proposed a finite element model of a piezoelectric cantilever beam using Hamilton variation principle. They employed the independent modal space method based on linear quadratic regulator and verified the effectiveness of the proposed model and AVC with numerical solutions and experiments. Umar et al. [9] employed integrating a non-symmetric Bouc–Wen model with an Auto-Regressive model with eXogenous inputs to enhance AVC of a cantilever beam using MFC. Their approach shows potential for practical engineering applications. Teoh et al. [10] used the receptance method for the robust control of DigiTwin Operating Platform. They proposed an algorithm which involves estimating the transfer function approximately by Rational Fraction Polynomial curve fitting method from the frequency response function matrix found experimentally, and then obtaining the control gains. Karagülle et al. [11] realized the integration of the control action into ANSYS solution for the AVC of a cantilever beam with piezoelectric actuation.

There are studies on AVC of cantilever structures in civil engineering and aerospace. Ito et al. [12] proposed using cameras and an image processing for the displacement for feedback masurements for AVC in the low frequency band (0-1.5 Hz), which is applicable in high rise buildings. Ramírez-Neria et al. [13] presented active disturbance rejection control for seismically excited building structures. They used five-storey building structure, active mass damper as controlling actuator and a shaking table as the disturbance. Gehni et al. [14] explored deep reinforcement learning (DRL) for controlling vibrations in building structures. They proposed a DRL-based control framework that learns to adjust the control signal using adaptive linear quadratic regulator.

Zhang et al. [15] presented an approach to mitigate vibrations induced by external shock on composite structures, which are used in the fuselage of airplanes. They used deep learning network to enhance control efficiency. He et al. [16] studied AVC of a Blended-Wing-Body airplane wing by using a smooth-switching linear-parameter-varying dynamic output-feedback control. Li et al. [17] proposed an improved parallel multi-input multi-output multi-frequency filtered-x least mean square control algorithm for AVC of aircraft panel structure using piezoelectric stack actuators. Dong et al. [18] presented phase compensation active disturbance rejection control for shimmy vibration using magnetorheological damper. Sahin et al. [19] studied smart beams and fins structures with bonded piezoelectric patches to suppress free and forced vibrations using H∞ and sliding mode controllers. They presented simulation and experimental results. Prakash et al. [20] designed a Field Programmable Gate Array based multichannel feedforward adaptive Filtered-X Least Mean Square algorithm based AVC. They deployed the developed system for AVC testing of the full-scale wing of an all composite two-seater transport aircraft and obtained encouraging results.

In this study, it is aimed to study AVC of cantilever structures with ANSYS. Commercial computer aided engineering (CAE) programs are widely used in the simulation age with computers for analyzing real complex engineering systems. The user defines the system and inputs only, and CAE software develops mathematical model and gives the solutions for static, harmonic and transient analyses. There is no direct modeling and solutions for closed loop control systems in ANSYS. A procedure which integrates control action into the transient solution in ANSYS APDL is given in this work. The validity of the procedure is tested by comparing the results with the solutions obtained with the Laplace transform and Newmark methods for a cantilever beam with end mass. Experimental results are also given. Then the procedure is used for AVC of an aircraft wing.

This paper consists of following sections. Simulations of AVC of cantilever beam with end mass with FE modeling, Laplace transform method, Newmark method and ANSYS APDL are given in Section 2. The experimental system is described in Section 3. Simulation of AVC of aircraft wing using ANSYS APDL is explained in Section 4. Results and discussions are presented in Section 5. Finally, the conclusions of the study are given in Section 6.

2. Simulation Of Active Vibration Control of Cantilever Beam with End Mass

2.1. Finite Element Model for Laplace Transform and Newmark Method Solution

The mathematical model of the system under study is developed for the analysis with Laplace Transform and Newmark methods. The cantilever beam considered is shown in Figure 1.

The nodes are numbered in normal parantheses as (1) to (5). The FE assembly is defined in Table 1.

The local origin is at Node-4 and the end point is at Node-5; the local x axis, x₂, is in the direction from Node-4 to Node-5; the length is 282 mm; the orientation angle is 90^o for FE-2 for example. The nodal displacement indices of the assembled system in global coordinates are shown in curly parantheses in Figure 1. For example, the indices for Node-5, {10,11,12}, correspond to [d₁₀, d₁₁, d₁₂]^’, where d₁₀, d₁₁, and d₁₂ are the displacement in the x direction, the displacement in the y direction, and the bending rotation of the cross section about z axis at Node-5, respectively. The prime symbol indicates the transpose of the matrix. The known nodal displacements are numbered at last following the symbol asterix (*). For FE-2, the displacement vector with the indices for the assembled system is [d₇,d₈,d₉,d₁₀,d₁₁,d₁₂]’.

The end mass is located at Node-3. The degree of freedom of the system, n_dof, is 12. The impulse force, f_d, applied in the x direction at Node-3 (in the direction of d₄) is considered as the disturbance. The feedback signal in AVC, d_m, is considered as the displacement in the x direction of the end mass. Thus, d_m=d₄. The active vibration controlling velocity actuation, v_s, is applied to the support in the x direction at Node-1. Thus, v_s=$\dot{\mathrm{d}}$_s=$\dot{\mathrm{d}}$₁₃.

Each finite element (FE) has stiffness and mass matrices in local coordinates as [K_FE-L]_n, and [M_FE-L]_n, respectively. Here, n is the index of FE. The FE stiffness and mass matrices in local coordinates for the Bernoulli-Euler beam are given below [21].

(1)

Here E, ρ, h, A, and I are the elasticity modulus, density, the length, cross sectional area, and bending moment of inertia of the cross section. The transformed FE stiffness and mass matrices in global coordinates of FE-n are [K_FE]_n, and, [M_FE]_n, can be found from the following equations [21].

(2)

where the transformation matrix is given as

(3)

Here, $\mathsf{\phi }$ is the orientation angle of the FE.

The mathematical model of the system can be written as [M]{$\ddot{\mathrm{d}}$}+[C]{$\dot{\mathrm{d}}$}+[K]{d}={f}, where [M], [C], and [K] are mass, damping, and stiffness matrices of the system having a size of (15x15); {d} is the displacement vector of the system with 15 rows; and {f} is the nodal forces vector of the system with 15 rows. The assembled system matrices can be obtained using Table 1 as given below as examples:

9th row, 7th column of [M]=6th row, 4th column of [M_FE]₁+3rd row, 1st column of [M_FE]₂

10th row, 10th column of [K]=4th row, 4th column of [K_FE]₂+1st row, 1st column of [K_FE]₃

The system damping matrix of the system can be written as [C]=α[M]+β[K], where α and β are Rayleigh damping coefficients [21].

The system displacement vector can be written as

{d}=[d₁, d₂, d₃, d₄, d₅, d₆, d₇, d₈, d₉, d₁₀, d₁₁, d₁₂, d₁₃, 0, 0]’

Here, d₁₃=0 for without AVC (Control-Off case). $\dot{\mathrm{d}}$₁₃=$\dot{\mathrm{d}}$_s=v_s is the controlling input for AVC (Control-On case).

The system nodal force vector can be written as

{f}=[0, 0, 0, f_d, 0, 0, 0, 0, 0, 0, 0, 0, F_Ax, F_Ay,M_A]’

where F_Ax, F_Ay, and M_A are the reaction forces at the Node-1.

The system mass matrix of the system, [M], is sub-divided to four matrices; [M_a], [M_b], [M_c], and [M_d]. [M_a] is obtained by taking 1-12 rows, and 1-12 columns of [M]. [M_b] is obtained by taking 1-12 rows, and 13-15 columns of [M]. [M_c] is obtained by taking 13-15 rows, and 1-12 columns of [M]. [M_d] is obtained by taking 13-15 rows, and 13-15 columns of [M]. [K_a], [K_b], [K_c], [K_d], [C_a], [C_b], [C_c], [C_d] matrices are obtained similarly.

The assembled system displacement vector is sub-divided to two vectors as

{d_a}=[d₁, d₂, d₃, d₄, d₅, d₆, d₇, d₈, d₉, d₁₀, d₁₁, d₁₂]’ and {d_b}=[d₁₃, 0,0]’

The assembled system nodal force vector is sub-divided to two vectors as

{f_a}=[0, 0, 0, f_d, 0, 0, 0, 0, 0, 0, 0, 0]’ and {f_b}=[F_Ax, F_Ay,M_A]’

Here {d_b} and {f_a} are time dependent inputs, and {d_a} and {f_b} are to be calculated as outputs.

The following equations can be written:

(4)

(5)

Equation (4) is the mathematical FE model of the system, which is a set of linear set of differential equations and can be solved analytically by the Laplace transform method, or by numerical methods like Newmark method [21]. Equation (4) is solved to find {d_a}, and then Equation (5) is solved to find {f_b].

A MatLAB program is developed to obtain [M_a], [M_b], [M_c], [M_d], [C_a], [C_b], [C_c], [C_d], [K_a], [K_b], [K_c], [K_d] from the inputs in Table 1, E, ρ, I, A.

2.2. Simulation by Laplace Transform Method

The following equation can be found by taking the Laplace transform of Equation (4) [22]

(6)

where

(7)

Here, {D_a}, {F_a}, and {D_b} are the Laplace transforms of {d_a}, {f_a}, and {d_b}, respectively.

The displacement in the x direction at Node-3, d₄=d_m, is considered as the feedback for the closed loop control. The velocity in the x direction at Node-1, $\dot{\mathrm{d}}$₁₃= $\dot{\mathrm{d}}$_s=v_s, is considered as the controlling actuator input. Let the transfer function from d_s to d_m be G_a. Then,

(8)

The disturbance force, f_d, is applied in the x direction at Node-3 (in the direction of d_m). Let the transfer function from f_d to d_m be G_f. Then,

G_f(s)= The 4th row, 4th column of [H_o]

(9)

Let G_a(s)=N_a(s)/D_a(s), and G_f(s)=N_f(s)/D_f(s); where N_a(s), D_a(s), N_f(s), and D_f(s) are polynomial functions of s. The open loop block diagram of the control system is shown in Figure 2.

D_s, F_d, and D_m are the Laplace transforms of d_s, f_d, and d_m, respectively.

A MatLAB program is developed to find the polynomials N_f, D_f, N_a, and D_a using the FE model developed in Section 2.1, and the Equations (8) and (9). N_f and D_f is divided to a constant (C_f), and N_a and D_a are divided to a constant (C_a), so that D_f=D_a is obtained. The commands conv, expand, and sym2poly are used.

The closed loop block diagram of AVC is shown in Figure 3.

D_mr=D_m-D_s is the relative displacement is multiplied by the sensor constant (K₂) to obtain the output voltage, V_o. V_o is the feedback signal and V_r is the reference voltage. The reference voltage is taken as zero for AVC. N_c(s) and D_c(s) are numerator and denominator of the control block. K₁ is the actuator constant, v_s the time derivative of d_s, which is generated by the actuator (servo motor). The transfer function from F_d to V_o, H_f(s), can be obtained as

(10)

Only proportional control is considered, so N_c=K_p and D_c=1 where K_p is the proportional control constant.

The impulse response is found from the transfer function, H_f(s), using the commands tf and impulse in MatLAB. The roots of the denominator of H_f are eigenvalues. The time step is taken as Δt<T_max/20 where T_max is the period for the highest natural frequency considered [23].

Note that K₁ is multiplied by -1. This can be explained as the actutator applies velocity and the resulting controlling inertia force is in the opposite direction of the acceleration.

2.3. Simulation by Newmark Method

Newmark method [21] is used to solve Equation (4) numerically. The solution starts with initial values at t=0. A time step, Δt, is used. The solution at the time t_n+1=(n+1) Δt is found using the solution at the previous time t_n=nΔt by numerical integration.

A MatLAB program is developed to study AVC of the system in Figure 1. The flow chart of the program is given in Figure 4.

The mathematical model of the system is established in B1 using the procedure given in Section 2.1 The matrices shown in square brackets are square matrices with size of 12x12. K_p, A_d, Δt, t_end, K₁ and K₂ are defined in B2. A_d is the amplitude of the impulse disturbance, t_end is the time at which the transient analysis is ended. The Newmark parameters are defined in B3. The initial values are given in B4. n_s is the number of steps. u(1:n_dof,1:n_s) is the array for nodal displacements. The size of the two dimensional array is (n_dof, n_s). u(5,18) corresponds to the value of the displacement, d₅, at the time step t₁₈, for example. u₁(1:n_dof,1:n_s) is the array for nodal velocities. u₂(1:n_dof,1:n_s) is the array for nodal accelerations. f_a(1:n_dof,1:n_s) is the array for nodal forces. The nodal force is applied to the degree of freedom, 4, only as an impulse to create the disturbance. The impulse is modeled by assigning the non-zero value of the nodal force as A_d at nodal degree of freedom 4 at the time step, t₂, only, by defining f(4,2)=A_d. The unit impulse can be applied if A_d is taken as 1/ Δt [23].

The Newmark method is applied in B5 to calculate nodal displacements, velocities, and accelerations at the time steps. The control action is integrated to the solution by taking error as

errv=0-V_o, where 0 is the reference value, and V_o is the sensor signal at the time step n. The controlling velocity actuation is v_s=-K₁K_perrv. The velocity is applied in the direction of d_s. The displacement at the time step n is calculated by numerical integral, d_s=d_s+v_sΔt. The acceleration at the time step is calculated by numerical differentiation, a_s=(v_s-v_sp)/Δt, where v_sp is the velocity at the previous time step. {f_sum} is the total force (the right side of Equation (4)) at the time step n, which is a vector of 12 rows. The output signal is found as V_o=K₂[u(4,n)-d_s]. Vo is stored to the array {v_out}, where v_out is a vector with n_s rows. The calculations in B5 is continued as a loop until n>n_s, and {v_out} is stored to a file.

2.4. Simulation by ANSYS APDL

Ansys Parametric Design Language (APDL) can be used to study active vibration control of structures [11]. The system in Figure 1 is modeled using the Preprocessor of APDL in the first part of the program starting with /prep command. Material properties are defined using the mp command. The Rayleigh damping coefficients are defined by the alphad and betad commands. FE types are defined by the et command. The BEAM188 and MASS21 elements are used. The section properties of the beam are defined using sectype and secdata commands. The mass property (171 gr) is defined using the r command. The nodes (1)-(5) are defined by the n command using the locations of the nodes with the given values in Table 1. The beam FE’s are defined by e command using Table 1. The mass FE is defined at Node-3. All the nodes are constrained in the y and z directions, and rotations about x and y axes using the d command. Node-1 is constrained in the x direction and rotation about z axis also because it is clamped. So, all the nodes except Node-1 have the degree of freedom in the x direction and rotation about z axis.

The transient analysis with the integration of the control action is programmed in the second part of the APDL code starting with /solu command. The values of K_p, A_d, Δt, t_end, K₁ and K₂ are defined. A new transient solution is started using antype command. The time step, Δt, is assigned by the deltim command. A file is opened to write the results. The transient solution and the integration of the control action are achieved using do-enddo loop as given below:

*do,n,2,ns,1

 t=(n-1)*dt

 fd=0

 *if,t,eq,dt,then

  fd=Ad

 *endif

 f,3,fx,fd

 time,t

 solve

 *get,dm,node,3,u,x

 *get,ds,node,1,u,x

 Vo=K2*(dm-ds)

 errv=0-Vo

 vs=-K1*Kp*errv

 d,1,velx,vs

 *vwrite,t,Vo

 (E15.8,',',E15.8)

*enddo

The variable ns is the samples of the transient analysis and it is defined before the loop as ns=tend/dt+1. The value of the disturbance, fd, is zero for all the samples except the sample n=2, to model the impulse as given in Section 2.3. The disturbance force is defined to Node-3 in the x direction by the command “f,3,fx,fd”. The displacements of Node-1 and 3 in the x direction are assigned to the variables dm and ds, respectively. The output value for the sample n is as calculated as Vo=K2(dm-ds). The error value is calculated as errv=0-Vo. The control actuation is calculated as vs=-K1*Kp*errv. Note that the sign of vs is negative as shown in Figure 3 and explained in Section 2.2. The velocity actuation is applied to Node-1 by the command d,1,velx,vs. The time and output signal, Vo, is written to an output file in the loop. The output file is closed after exiting the loop, and the program is ended.

3. Experimental System of Active Vibration Control of Cantilever Beam

The experimental system is shown in Figure 5.

The cantilever beam is mounted on the car of a linear motion system driven by a servo motor. The ball screw of the linear motion system has a lead of 5 mm per motor revolution. The motor has maximum speed as 3000 rpm, and it is driven in the velocity mode by applying ±10 V range. So the relation between the control voltage applied to the motor driver, V_motor, and the speed of the car, v_s, can be written K₁=v_s/V_motor, and it can be calculated as K₁=0.025 (m/s)/Volt. The 200 watt servo motor is Mitsubishi-Model HC-KFS23B and the driver is Mitsubishi-Model MR-J2S20A. The geometry given for the simulation in Figure 1 is used in the experiments.

Four strain-gauges located near the bottom point (Node-1 in Figure 1) connected with Wheatstone bridge, and the output voltage of a load-cell amplifier (LCA) is used as the feedback signal in experiments. Two strain gauges are mounted at the surface at x=1 mm, and the other two strain gauges are mounted at the surface at x=-1 mm. The direction of the strain gauges are in the y direction. The strain gauges located at x=1 and the strain gauges located at x=-1 give resistance changes in opposite signs and a Wheatstone bridge is constructed as in load cells [24]. The strain gauges are Tokyo Measuring Instruments FLAB-5-11-1LJC-F having a length of 5 mm, a gauge resistance of 120 Ω, and a gauge factor of 2.1. The Weighing Transmitters Amplifier Weight Sensor Voltage Current Converter DC 0-5V/0-10V 4-20MA (Model JY-S60) is used for amplification. The output voltage of LCA amplifier, V_lca, is used as the strain feedback signal for AVC in the experiments.

The laser displacement sensor was used to compare displacement signals with V_lca. Keyence LK-G157, and the base station Keyence LK-3001 was used for the displacement measurements. The base station gives the output voltage, V_lds, in the ±10 V range for the displacement set values ±40 mm. The relation between V_lds, and the displacement, d_mr (relative displacement in the x direction at Node-3 in Figure 1) can be written K₂= V_lds / d_mr, and it can be calculated as K₂=250 Volt/m.

Adlink-9222 PCI Multi-Function DAQ card is used as the controller. It has technical specifications as 16 bit 250 kHz ±10 V 16 analog inputs, 16 bit 1 MHz ±10 V 2 analog outputs. Ubuntu 20.04 with Kernel 5.4.0-47-generic is used as the operating system in the desktop computer. Adlink provides drivers for C/C++ with the “dask.h” library.

The flow chart of the C program for AVC of the cantilever beam with Adlink-9222 is shown in Figure 6.

The libraries are included, the values K_p, t_delay,t_end, and α are defined. V_lca is read from the analog input (AI) of Adlink-9222 card before the loop starts and assigned to V_lca0. V_lca0 is the voltage when there is no vibration before the disturbance impulse. t_delay is the delay time to change the sampling frequency, f_s. The disturbance impulse is applied manually to Node-3 in Figure 1 in the x direction in an approximately 1 s after the loop starts. The strain feedback voltage, V_lca read from AI at a time step in the loop is replaced by V_lca-V_lca0. The resulting V_lca is filtered by using Infinite Impulse Response (IIR) filter [25], to eliminate the noise in the experiments. The filtered voltage, V_filtered, is calculated as

V_filtered=αV_lca-(1-α)V_filtered. The cut-off frequency, f_c, of the low-pass filter is calculated as

f_c=αf_s/[2π(1-α)]. The error value is calculated as errv=0-V_filtered. The voltage, V_motor=-Kp*errv, is applied to the servo motor driver through the analog output (AO) of the card. The time is read, and a delay is applied. The loop continues until t<_tend. The time history of V_lca is saved to a file, and the program is ended. The C program for the AVC of the cantilever beam is given in Appendix A.

4. Simulation of Active Vibration Control of Aircraft Wing by ANSYS APDL

The simulation of AVC of an aircraft wing by ANSYS is considered in this section. Solidworks model of and ANSYS beam-shell model the aircraft wing is shown in Figure 7.

The geometry of NACA 4412 airfoil [26] is used to model the ribs. The length of the foil is 300 mm. The thickness of the ribs is 4.2 mm. The wing has 13 ribs. The mid-plane of Rib-1 is at z=0. The z coordinates of the mid-axes of the 13 ribs are 0, 50.8, 130.8, 210.8, 290.8, 370.8, 450.8, 530.8, 611.4, 661.4, 711.4, 760.8, 810.8 mm, respectively. There are 5 beams along z axis. The origins of the beams located at z=0 have the center locations (x=5.5 mm, y=0.6 mm), (x=81.2 mm, y=10.7 mm), (x=100.6 mm, y=6 mm), (x=198 mm, y=9.8 mm), (x=211 mm, y=5.6 mm), respectively. The beam cross sections have shapes and dimensions as circular pipe with D=5 mm and t_k=0.5 mm; circular pipe with D=16 mm and t_k=1 mm; rectangular with w=3 mm, h=25 mm; circular pipe with D=14 mm and t_k=1 mm; rectangular with w=3 mm, h=15 mm, respectively. Here, D is the outer diameter, t_k is the thickness, w is the width and h is the height. Beam-6 is used to connect the beams 2 and 4 (Beams with D=16mm and D=14) with the dimensions 10 mm x 30 mm x 157.5 mm. Node-1 on Beam-6 is fixed for the Control-Off case and the controlling velocity input is applied to Node-1 for the Control-On case. The disturbance impulse is applied to Node-3 on the rib at z=0 and the displacement feedback is gathered from the node. Node-1 and 2 are located at the centroid of the airfoil area.

In the first part of the ANSYS-APDL program, the ribs are modeled as shell FE’s (SHELL181) and the beams are modeled as beams (BEAM188) in ANSYS. The element types are defined by et command. Material properties are defined by mp command. Rayleigh damping coefficients are defined by alphad and betad commands. Beam and shell sections are defined using sectype and secdata commands. Keypoints are created using k command at the locations of the rib foil sample coordinates and at the locations of beam centers for all the ribs. Nodes 1 and 3 are defined using n command. Areas are created to model the ribs using a command and the keys created. Sub-areas are used to model a rib so that the centers of the beams coincide with the corners of sub-areas. aglue command is used to connect the sub-areas. Lines are created to model beams using l command and the keypoints created. Areas are meshed using amesh command, and lines are beams are meshed using lmesh command. nummrg command is used to merge the nodes created at the same locations. The keypoints and areas created for Rib-1 at z=0 are repeated for subsequent 12 ribs with the same x-y coordinates and changing z coordinates and APDL commands for this are created using VisualBASIC program. The second part of the APDL program is the same as given in Section 2.4.

5. Results and Discussions

5.1. Comparison of Experimental Strain and Displacement Sensor Output Signals

The amplitude, A_d, of the impulse in the simulations is selected so that the output signal has a maximum value of 2.5 V in the solution with Laplace transform method for Control-Off case (Kp=0). The value, 2.5 V, corresponds to 10 mm maximum mass displacement. Let the value of the maximum of the output voltage is V_p for the gain K_p for Control-On case in the solution with Laplace transform method. The impulse disturbance in the experiments are applied manually, and thus the amplitude of the impulse is random. The output signal in the experiments, V_o(t), is scaled with C_p=(V_p/V_max)V_o, where V_max is the maximum of V_o, for comparisons. The output signals are shifted so that the time at which the first peak occurs is 0.5 s.

Let rms₀ and rms_p are the rms values [24] for the output signals for Control-Off case and Control-On case for K_p. The suppression of the vibration with AVC is quantified with the relative value of r_p (%) and given as r_p=100(rms₀-rms_p)/rms₀.

The scaled strain output signal, V_lca, obtained from the experimental system for the free vibration of the cantilever beam after an impulse disturbance is shown in Figure 8a. Windowed 2V_lca and V_lsd signals are shown in Figure 8b. It is observed that 2V_lca≈V_lds. Displacement output signals are used in simulations because it can be implemented in analyses. ANSYS gives displacement values but does not give strain values in the *do *enddo loop in the transient solutions. Strain values are obtained using etable command in the postprocessing stage after the solution stage. On the other hand, displacement sensors require additional mounting for installation and are not practical in applications. Strain gauges can be easily installed and give a cost-effective solution for feedback signals. The observation in Figure 8b shows that the results obtained in simulations using the displacement feedback can be used in applications using the strain feedback.

5.2. Determining Damping Constant Experimentally

Rayleigh damping [21] is considered in simulations. The damping matrix is defined as [C]=α[M]+β[K] where α and β are the proportional damping constants. Modal damping ratios are given as ξ_i=α/(2ω_i)+βω_i/2, where i is the i-th mode, and ω_i is the i-th natural frequency. The modal damping ratio can be found using the logarithmic decrement formulation [22]. The stiffness damping only is considered in this study. The amplitude versus frequency plot obtained from the Fast Fourier Transform [25] of the signal in Figure 8a is shown in Figure 9.

The peak values can be read on the signal in Figure 8a and it is found that the amplitude of the 3rd peak, and the amplitude of the 8th peak are 2.3209 and 2.2363, respectively. Then the damping ratio for the first natural frequency (7 Hz) is found from ξ₁=δ/$\sqrt{4{\mathsf{\pi }}^2+{\mathsf{\delta }}^2}$ where

δ=[1/(8-3)] Ln(2.2363/2.3209). Then Rayleigh stiffness damping coefficient, β, is found from ξ₁=βω₁/2 where ω₁=43.9823 rad/s as β=5.3702 (10^-5). This value found from experimental data is used for simulations.

5.3. AVC of Cantilever Beam Results Obtained by Laplace Transform Method, Newmark Method, ANSYS and Experiments

The computer processor and t_delay determines the sampling frequency, f_s, in the experimental system. t_delay is taken as 1000 µs and f_s is measured as 790.1420 Hz. The FFT gives the spectrum 0-790.1420/2 Hz and only the peak value at f=7 Hz is observed in the spectrum as shown in Figure 9. Thus, the first natural frequency is considered in simulations. The first un-damped natural frequency obtained from the FE model given in Section 2.1 is found as 7.2222 Hz, using eig command in MatLAB. The first natural frequency from ANSYS model (Section 2.4) is found as 7.22 Hz, using modal analysis in ANSYS. It is observed that simulation and experimental results are in good agreement for the first natural frequency.

The time step, Δt, is taken as 0.005 s, which is smaller than 1/7.2222/20=0.0069 s, in simulations. The amplitude of the unit impulse is taken as A_d Δt in Laplace transform method. A_d is taken as the value of the disturbance force at the time Δt and zero other time steps in Newmark method and ANSYS. A_d value is found as 21.989 so that the output signal has a maximum value of 2.5 V in the solution with Laplace transform method for Control-Off case (Kp=0). V_p values for scaling are found as 2.5, 2.4698, 2.4401, 2.3822, 2.3539, 2.2750, 2.2301, and 2.0237 for K_p values 0, 1, 2, 4, 5, 8, 10, and 20, respectively.

The values of K₁ and K₂ are taken as 0.0025 and 250; Rayleigh damping coefficients α and β values are taken as 0 and 5.3702(10^-5) in simulations. The constants C_f and C_a defined in Section 2.2 are 10²²⁴ and 128(10²²⁴), respectively, for the selected damping constants. The output signals for selected K_p values are shown in Figure 10, Figure 11, Figure 12 and Figure 13. The video for experiments is given in Supplementary Materials.

It is observed that the free vibrations after the impulse disturbance are suppressed with AVC (Blue lines: Control-Off, Red lines: Control-On). The rms values of the output signals are shown in Table 2.

The suppression level of vibration for different K_p values are shown in Figure 14.

It is observed that simulation results with Laplace Transform method, Newmark method and ANSYS are in good agreement. Experimental results are in acceptable agreement with simulation results. r_p values are higher for experimental results than simulation results. r_p values approach to an asymptotic value with increasing K_p.

5.4. AVC of Aircraft Wing by ANSYS APDL

The APDL program developed for the AVC of the cantilever beam and verified with the results of Laplace transform and Newmark methods is extended for the aircraft wing structure given in Section 4. The smallest three natural frequencies found with ANSYS are found as f₁=10.1272, f₂=24.2032, and f₃=41.4421 Hz. Mode shapes are shown in Figure 15. Mode-1 is bending about x axis, Mode-2 is bending about y axis, and Mode-3 is torsion about z axis.

The impulse disturbance at Node-3 shown in Figure 7 is in the y direction and excites the wing mainly as bending about x axis, and thus Mode-1 is considered for AVC.

Δt is taken as 0.004 s, which is less than 1/10.1272/20=0.0049 s. The Rayleigh damping constants are taken as α=0 and β=0.0002 by inspection. The amplitude, A_d, of the impulse is selected so that the output signal has a maximum value of 2.5 V. A_d is found as 71.1269. The output signals for selected K_p values are shown in Figure 16.

The rms values of the output signals for AVC of aircraft wing are shown in Table 3.

The suppression level of vibration for AVC of aircraft wing for different K_p values are shown in Figure 17.

It is observed that the APDL program developed for the AVC of the cantilever beam and extended for the aircraft wing can be used successively. r_p values approach to an asymptotic value with increasing K_p for AVC of cantilever aircraft wing as observed for AVC of cantilever beam with end mass.

6. Conclusions

In this study, it is aimed to study active vibration control (AVC) of cantilever structures with ANSYS. Commercial computer aided engineering (CAE) programs are widely used in the simulation age with computers for analyzing real complex engineering systems. Developing mathematical model and using Laplace transform method was the major approach to analyze closed loop control systems before the computer age. Graphical methods were very important and analyses of complex real engineering systems were limited. Numerical methods like Newmark were started to be used by developing computer codes at the early age of computers. Today, the commercial computer aided engineering (CAE) programs like ANSYS enables one to model systems and define inputs only, and CAE generates the mathematical model and solves the equations for static, harmonic, and transient analyses. Developing engineering systems followed produce, test, and revise procedures before the computer age changed to simulate, revise, and produce procedures today. There is no direct modeling and solutions for closed loop control systems in ANSYS. A procedure which integrates control action into the transient solution in ANSYS APDL is proposed in this work. It should be noted that one must be very careful using the simulation methods because a minor mistake in the process can give erroneous results. Testing the validity of the CAE procedures by Laplace transform and numerical methods for a basic structure is important before applying it to complex structures.

In this study, a cantilever beam with an end mass is considered first. The disturbance is considered as an impulse to the end mass. The controlling action is considered as the velocity actutation to the base, which can be generated with servo motors and servo valves. The strain signals amplified by a load cell amplifier are considered as the feedback in experiments. Strain feedback cost effective and can easily be implemented in real engineering applications. Displacement feedback is used in simulations because of the difficulty to obtain strain signals in the transient solution stage in ANSYS. Strain and displacement feedback signals compared experimentally and it is observed that there is a linear relation between them approximately. Finite element model of cantilever beam with end mass is developed for AVC analysis with Laplace transform and Newmark methods. The control action is integrated to the transient solution in Newmark method. There is no need to develop mathematical model after defining the systems in ANSYS. The control action is integrated to the transient solution like in Newmark method. An experimental system, which uses Adlink-9222 PCI DAQ card and a C program is introduced. It is observed that the controlling velocity sign must be in the oppsite direction of error sign.

It is concluded that the developed ANSYS procedure can be used for AVC analysis of cantilever structures. Finally, the ANSYS procedure is applied for AVC of a cantilever aircraft wing successfully. The initial selections of sensor, actuator and control parameters can be determined with iterations inANSYS, and then the control parameters can be improved on produced systems using like the receptance method.

Only proportional control action is used in this study. The sudy can be extended for various control methods. AVC of multiple input and output systems, torsional vibrations in aircraft wings, different structures can be considered for further studies.