Whirl flutter optimisation-based solution of twin turboprop aircraft using a full-span model

Whirl flutter is a specific type of flutter instability, relevant for turboprop aircraft, caused by the effect of rotating parts as a propeller or a gas-turbine engine rotor. The proposed optimisation-based analytical procedure is used to determine the critical values of the engine attachment stiffness parameters for the preselected flutter speed. For the half-span model, two design variables are used. The objective function is defined as the minimization of the engine vibration mode frequency sum. Design constraints keep the engine frequency ratio and the flutter stability at the selected velocity. However, application of a full-span model is necessary in some cases. In this case, special models of both symmetric and antisymmetric engine vibrations and four design variables must be used. Design constraints maintain the pitch mode frequency ratio, the yaw mode frequency ratio and the critical mode frequency ratio. Critical modes are dependent on the relation between the rotational direction of both propellers (identical or inverse). A flutter design constraint is applied as well. The described methodology is demonstrated on the application example of a twin-engine commuter aircraft. Demonstrated cases include symmetrical revolutions of propellers for both identical and inverse directions of rotation, cases of single engine failure and single propeller feathering, and finally, cases of unsymmetrical revolutions including the reduced and increased revolutions of a single propeller, for both identical and inverse directions of rotation. c © 2017 University of West Bohemia. All rights reserved.


Introduction
Airworthiness regulation standards require certification of turboprop aircraft structures with respect to whirl flutter stability.Whirl flutter is induced by the effect of rotating parts, such as a propeller or a gas turbine engine rotor.Rotating mass increases the number of degrees of freedom and induces additional forces and moments.A rotating propeller also causes complicated flowfield and aerodynamic interference effects with a nacelle and a wing structure.The instability is driven by motion-induced unsteady aerodynamic propeller forces and moments acting on the propeller plane.Whirl flutter may cause unstable vibration of a propeller mounting, or even a failure of an engine installation or a whole wing.It was a cause of several serious accidents.Apart from the nominal state of an engine mount, regulation standards also require taking into account the influence of stiffness and damping variation of propeller, engine, engine-mount and aeroplane structures on whirl-flutter stability (FAR/CS § 23.629(e)(1)(2)).
During the aircraft development process, reliable data regarding the engine attachment stiffness are not usually available until the ground vibration test of the prototype is performed and the updating of the analytical model stiffness characteristics is possible.However, it is worth performing whirl flutter calculations in the earlier phase to save time in the final development phase.For this purpose, an optimisation-based whirl flutter analytical procedure was prepared.It allows for the determination of whirl flutter stability margins for the certification speed, which is set by regulations.The solution employs a gradient-based algorithm and includes modal-based and flutter-based design responses.Design variables are represented by the engine attachment stiffness parameters.
The initial solution was prepared for the half-span models of aircraft.Such models, with either symmetric or antisymmetric boundary conditions, are ordinarily used for flutter analyses.However, application of a full-span model is necessary in some specific cases.The typical example is the whirl flutter of the usual twin wing-mounted tractor engine aircraft, for which whirl flutter stability characteristics are influenced also by the directions of rotation of both propellers.For this reason, the applicability of the optimisation-based solution was enlarged to include the full-span models.The solution is, compared to that of a half-span model, more complicated and the specific adjustment of a structural model is necessary.This paper describes the background and the methodology of the solution.The application is demonstrated on a reference model of a twin-engine turboprop aircraft.Whirl flutter stability boundaries are constructed for different conditions.Demonstrated cases include symmetrical revolutions of propellers for both identical and inverse directions of rotation, single engine failure and single propeller feathering and, finally, non-symmetrical revolutions including reduced and increased revolutions of a single propeller, again for both identical and inverse rotational directions.

Theoretical background
The principle of the whirl flutter phenomenon is outlined on a simple mechanical system with two degrees of freedom.The propeller and hub are considered as rigid.A flexible engine mounting is represented by two rotational springs of stiffnesses K Ψ and K Θ , as illustrated in Fig. 1.
Such a system has two independent mode shapes (yaw and pitch) with angular frequencies ω Ψ and ω Θ .For a propeller rotation with angular velocity Ω, the primary motion changes and the gyroscopic effect cause both independent mode shapes to merge into a whirl motion.The axis of rotation of the propeller exhibits an elliptical movement.The trajectory of this elliptical movement depends on both angular frequencies ω Ψ and ω Θ .The orientation of the gyroscopic movement is backward relative to the propeller rotation for the mode with the lower frequency Fig. 1.Gyroscopic system with propeller (backward whirl mode) and forward relative to the propeller rotation for the mode with the higher frequency (forward whirl mode).Because the yaw and pitch motions have a 90 • phase shift, the mode shapes in the presence of gyroscopic effects are complex.
The gyroscopic motion results in changes in the propeller blades' angles of attack, consequently leading to unsteady aerodynamic forces.These forces may, under specific conditions, induce whirl flutter instability.The flutter state is defined as the neutrally stable state with no damping of the system, and the corresponding airflow (V ∞ = V FL ) is called the critical flutter speed.If the air velocity is lower than flutter speed (V ∞ < V FL ), the system is stable and the gyroscopic motion is damped (see Fig. 2a).If the airspeed exceeds the flutter speed (V ∞ > V FL ), the system becomes unstable, and the gyroscopic motion is divergent (see Fig. 2b).An analytical solution is sought to determine the aerodynamic force caused by the gyroscopic motion on each of the propeller blades.The equations of motion were derived for the system shown in Fig. 1.The kinematical scheme is shown in Fig. 3.We select three angles (ϕ, Θ, Ψ) as the independent generalised coordinates.The rotating part is assumed to be cyclically symmetric with respect to both mass and aerodynamics (i.e., a propeller with a minimum of three blades).The propeller angular velocity is considered constant (ϕ = Ωt).Non-uniform mass moments of inertia of the engine with respect to the pitch and yaw axes (J Z = J Y ) are considered.The equations of motion become when small angles are considered, where a is the distance between the propeller plane and the vibration mode node point.Quantities γ Ψ and γ Θ are structural damping ratios of an engine yaw and pitch vibration modes.Neglecting the aerodynamic inertia terms ( Θ * ≈ Θ, Ψ * ≈ Ψ), the propeller aerodynamic forces and moments at the propeller plane (P Y , P Z , M Y,P , M Z,P ) are calculated as where ρ is the air density and R is the propeller radius.Aerodynamic derivatives (c-terms) are given by the propeller blade integrals.These integrals are usually calculated according to Houbolt and Reed [2] or according to Ribner [7,8].We use quasi-steady theory [3] and the effective angles (Θ * , Ψ * ) then become The effective angles are given as the effective static angles (terms in brackets).Optionally, the downwash (w 1 /V ∞ ) and sidewash (w 2 /V ∞ ) angles behind the propeller describing the interference between a propeller and nacelle may be added.The induced downwash and sidewash terms, which are dependent on the reduced frequency, can be obtained from the lift solution by partitioning the interference coefficients.The downwash effect influences the aerodynamic stiffness matrix; the influence on the aerodynamic damping matrix is neglected.These interference effects may be important, especially for the wing-mounted engine aircraft, as the effect is usually destabilizing.
Finally, seeking for the critical (flutter) state (assuming the harmonic motion) has a character of an eigenvalue problem.The final whirl flutter matrix equation can be expressed as where [M] is the mass matrix, [D] is the structural damping matrix, [K] is the structural stiffness matrix.[D A ] and [K A ] are the aerodynamic damping and aerodynamic stiffness matrices, respectively: These matrices include aerodynamic derivatives (c-terms) and geometry of the system (a, R).Finally, [G] is the gyroscopic matrix.The critical state emerges when the angular velocity ω is real.The critical state can be reached by increasing either V ∞ or Ω. Increasing the propeller advance ratio (V ∞ /(ΩR)) has a destabilizing effect.Another important parameter is the distance of the propeller plane from the node points of the engine vibration modes.Structural damping is a significant stabilization factor.By contrast, the influence of the propeller thrust is negligible.This small influence comes from the fact that the aerodynamic derivatives of the thrust propeller and windmilling propeller variance can be high in the low speed region, but at high velocities (where whirl flutter is expected), the variance is less than 5 % [6].For the rigid propeller blades, the whirl flutter inherently appears at the backward gyroscopic mode.
The most critical state is ω Θ = ω Ψ when the interaction between both independent motions is maximal and the trajectory of the gyroscopic motion is circular.
The described model that considers a rigid propeller is obviously applicable to conventional propellers, for which the propeller blade frequencies are much higher compared to the nacelle pitch and yaw frequencies.For the large multi-bladed propellers of heavy turboprop aircraft, the assumption of a rigid propeller appears to be too conservative and the blade flexibility must also be modelled.Obviously, the whirl flutter investigation of tilt-rotor aircraft must include even more complex analytical models [5].
Comprehensive information regarding the whirl flutter phenomenon can be found in [1].

Half-span model
For the purpose of aeroelastic analyses, simple dynamic beam structural models (stick models) are used.Structural stiffness is modelled using massless beam elements and structural inertia is modelled using concentrated mass elements, including the appropriate mass moments of inertia.Moreover, the model also includes spring elements, various conditions, multi-point constraints and auxiliary elements.The half-span model (see Fig. 4) includes a single side only with the half-values of stiffness and inertia at the plane of symmetry, as well as either symmetric or antisymmetric boundary conditions.Attachment of the engine to the wing is realized simply by means of two spring elements, which model engine pitch and yaw vibration modes.The stiffness constant of the spring determines the natural frequency of the appropriate mode.The spring element is stationed at the node point of the mode.

Full-span model
In the case of the full-span model (see Fig. 5), both symmetric and antisymmetric engine vibration modes must be modelled.These engine modes, in the typical order according to frequency, are 1) symmetric pitch, 2) antisymmetric pitch, 3) symmetric yaw and 4) antisymmetric yaw.The modes have diverse natural frequencies and diverse node points.Node points are typically stationed in the direction from the rear to the front (in the flight direction) in the order 1-3-2-4.The systems for modelling pitch and yaw engine vibration modes are separate.The appropriate rotational degree of freedom (i.e., around the lateral or vertical axis) is connected to the central system, which consists of the grounded spring element and two rod elements.The grounded spring is placed at the plane of symmetry at the station of the node of the symmetric mode, while the rod elements are placed at the station of the node point of the antisymmetric mode.Rod elements are oriented in the appropriate direction (i.e., laterally or vertically).The node point of the grounded spring is connected with the central node of the rod elements by means of a multi-point constraint.Apart from the appropriate degree of freedom, the other degrees of freedom are omitted from the analysis.The spring constant of the grounded spring element K δ and the torsional stiffness of the rod elements (GI k ) then determine the natural frequencies of both the symmetric and antisymmetric modes.The spring constant is decisive for the symmetric mode frequency, while the rod torsional stiffness is decisive for the antisymmetric frequency.However, there is also cross-influence; therefore, both parameters must be used jointly to set both frequencies.The examples of engine vibration modes are shown in Fig. 6.

Standard approach
In the standard approach, the input data are parameters of a structure and the outputs of the analysis are whirl flutter characteristics, i.e., V-g-f diagrams, and flutter speed and frequency.The solution is based on the strip aerodynamic theory [10] for the rigid propeller at the windmilling mode.For the residual structure, the unsteady doublet-lattice method is used in combination with the wing-body interference aerodynamic theory [4].Flutter stability analysis is performed using the p-k method [9].The basic flutter equation is expressed as where [M hh ], [B hh ] and [K hh ] are the modal mass, damping and stiffness matrices, respectively, which are functions of the Mach number (M) and the reduced frequency (k).Aerodynamic loads are included in the damping and stiffness matrices.[Q Re hh ] and [Q Im hh ] are the real and the imaginary parts of a complex aerodynamic matrix, which is also a function of parameters M and k.The parameter ρ is the air density, c is a reference length, and {u h } is a modal amplitude vector.The eigenvalue λ is given as and γ is a ransient decay rate coefficient.The structural damping coefficient (g) is expressed as In the case of the standard approach, the solution for the whirl flutter is performed for multiple velocities.The resulting quantities are V-g-f curves, i.e., the dependences of the damping and frequencies of the analysed modes on the flight velocity.The state with the zero damping represents the critical flutter state and the corresponding flight velocity is the critical flutter speed.

Optimisation-based approach
The optimisation-based approach employs gradient-based algorithms to obtain the whirl flutter solution.In this case, the flutter speed is set as an input parameter (certification speed), and the results are critical values of structural parameters.This solution enables us to obtain the stability margin for the specified certification speed using the calculated critical structural parameters.The analysed states are then compared with respect to the structural parameters and the relationship to the stability margin only.Such an approach can save large amounts of time because the number of whirl flutter analyses required by the regulations is dramatically reduced.
Two types of design responses (eigenvalue and flutter) are employed.The eigenvalue equation is: where λ n and φ n are the n-th eigenvalue and eigenvector, respectively; [K] is the structural stiffness; and [M] is the structural mass matrix.The equation is differentiated with respect to the i-th design variable x i : If equation ( 9) is premultiplied by φ T n , the first term becomes zero and it can then be solved for the eigenvalue derivatives: The solution method used for equation ( 10) is based on the semi-analytical approach in practice.The derivatives of the mass and stiffness matrices are approximated using the finite differences.The equation is solved for each retained eigenvalue referenced in the design model and for each design variable.
The aeroelastic flutter stability matrix equation is given by equation ( 5), which represents the p-k method of the flutter solution.This is the only method applicable for the purpose of the design optimisation.
Flutter design response computes the rates of change of the transient decay rate coefficient γ with respect to changes of the design variables.This equation is differentiated with respect to the design variables (∂γ/∂x i ).The solution is semi-analytical in nature, with derivatives approximated using either forward differences or central differences.
Contrary to the standard solution, the optimisation-based whirl flutter solution is performed for a single velocity.The resulting quantities are structural parameters, for which the flutter speed is equal to the specified certification speed.

Description of optimisation-based solution
The solution is demonstrated on an example involving the variation of engine pitch and yaw frequencies, which are the key parameters influencing whirl flutter.The inertial characteristics of the engine and propeller system, as well as the characteristics of the residual structure, are considered to have been reliably determined; thus, the engine attachment stiffness properties will be used as parameters for optimisation.

Solution for half-span model
For the half-span model, two design variables are defined: 1) effective stiffness of the engine attachment in pitch and 2) effective stiffness of the engine attachment in yaw.The design variables are directly related to the spring constants of the spring elements used to model the engine attachment (K Θ for the pitch and K Ψ for the yaw).
First, the preparatory step is performed.The purpose is to set the initial design variables for the main optimisation.The target frequency ratio (TFR), representing the ratio of the yaw and pitch frequencies, is set and both stiffness values are adjusted, so that the frequencies will reach the target ratio.
During the preparatory step, one of the stiffness values is selected as the design variable (while the other one is fixed) and the optimisation is performed.The objective function (OBJ) is defined as Note that ABS denotes the absolute value.The symbols f 1 and f 2 represent the pitch and yaw engine frequencies (f 2 is the higher frequency, whereas f 1 is the lower frequency).The solution may be obtained regardless of the frequency order; however, the yaw frequency is usually higher than the pitch frequency.The resulting quantities of this preparatory analysis are the input values of K Θ and K Ψ for the main optimisation.The ratio of the two frequencies (f Ψ for the yaw and f Θ for the pitch) is equal to the TFR.
After that, the second step (main optimisation) is performed.Both design variables (K Θ and K Ψ ) are used.Design constraints include the requirement to maintain the target frequency ratio, for which the ±2 % band is usually used in practical applications Another constraint requires the flutter stability (i.e., negative damping) at the selected certification speed V cert .This requirement is expressed as In the practical solution, the interval shift from the null value is given due to the numerical character of the solution, to prevent division by zero.The constraint is modified to The flutter constraint should also prevent another type of flutter instability below the certification speed that may be caused by the changes in the design variables.Typically, the modes within the frequency up to 100-120 Hz are included.The constraint should therefore be applied to all modes included in the solution.
The objective function of the main optimisation problem is defined simply as the minimization of the sum of the pitch and yaw frequencies, which is expressed as This simple formulation can be used, because, considering the fixed frequency ratio, there is a direct relation between the flutter speed and the frequencies of engine vibration modes.As the output, we will obtain the engine pitch and yaw stiffnesses, for which the flutter speed is equal to the specified certification speed and the yaw-to-pitch frequency ratio is equal to the specified target value.The described procedure is then repeated for several yaw-to-pitch frequency ratios, typically ranging from 1.05 to 2.0, to obtain enough points to construct a stability margin curve.The procedure is applicable regardless of whether the downwash effect (see Section 2) is included or not.If the downwash effects are to be included, the appropriate downwash terms must be calculated prior to the optimisation.
After any optimisation iteration, the cross-orthogonality correlation analysis of modes before and after the iteration is performed.The reason for this is the possibility of the switching of the mode order (with respect to frequency).If such a switch occurs, the modes must be re-ordered.The correlation analysis is performed using the Modal Assurance Criterion (MAC), which is expressed as where ψ 1 and ψ 2 are the correlated mode shapes.Note that only the engine pitch and yaw modes must be correlated with the residual modes during the verification.

Solution for full-span model
The optimisation solution for a full-span model includes four design variables: 1) effective stiffness of the engine attachment for symmetric pitch; 2) effective stiffness of the engine attachment for antisymmetric pitch; 3) effective stiffness of the engine attachment for symmetric yaw and 4) effective stiffness of the engine attachment for antisymmetric yaw.These design variables are related to the spring constants of the two grounded spring elements (K δ1 and K δ2 ) and to the torsional stiffnesses of the two pairs of rod elements [(GI k ) 1 and (GI k ) 2 ].However, the relation is not direct here due to the above mentioned cross-influences (see Section 3.2).Three frequency ratios are defined: 1) pitch frequency ratio (VFR = f AΘ /f SΘ ), 2) yaw frequency ratio (HFR = f AΨ /f SΨ ) and, finally, 3) critical frequency ratio (CFR).VFR and HFR are not changeable.The values are set according to the ground vibration test results or are guessed (by experience).The typical ratios for real aircraft structures range from 1.12 to 1.18.Subscripts of frequencies are denoted as follows: SΘ for symmetric pitch mode, AΘ for antisymmetric pitch mode, SΨ for symmetric yaw mode and, AΨ for antisymmetric yaw mode.
CFR is analogous to TFR, which was defined above.Critical modes are those modes whose combination causes flutter instability.The choice of the critical modes is mainly dependent on the mode order and on the relation between the directions of rotation of the two propellers.In the cases with the identical directions (i.e., CW-CW or CCW-CCW), the critical modes are symmetric pitch and antisymmetric yaw, or antisymmetric pitch and symmetric yaw.However, in cases with the inverse directions (i.e., CW-CCW or CCW-CW), the critical modes are symmetric pitch and symmetric yaw, or antisymmetric pitch and antisymmetric yaw.Note that CW denotes the clockwise direction and CCW denotes the counter-clockwise direction.
In the following description, we assume the typical frequency order (i.e., the pitch frequency is lower than the yaw frequency).
Similarly to the half-span model, the preparatory step is performed first to set the initial design variables for the main optimisation.The design constraint includes the requirement to maintain the highest frequency within the engine modes (SΘ, AΘ, SΨ and AΨ) which is typically f AΨ at the selected frequency value (f AΨT ).As usual, the ±2 % band is used for this purpose: The objective function is defined as the minimization of the frequency ratio error and is expressed as one of the following equations where SSQ denotes the sum of squares.
With respect to the character of whirl flutter mechanism, the first equation of ( 18) is applicable to the case of identical directions of propeller rotation, while the second is applicable to the case of inverse directions of propeller rotation.The resulting quantities represent the initial values of the design variables for the main optimisation.The frequency ratios are equal to the VFR, HFR and CFR values, respectively.
The main optimisation is performed similarly to that for the half-span model.The design constraints again include the requirement to maintain the frequency ratios, within the ±2 % band and Again, the first equation of ( 21) is applicable to the case of identical directions of propeller rotation, while the second is applicable to the case of inverse directions of propeller rotation.The constraint on flutter stability (i.e., negative damping) at the selected certification speed V cert is expressed in the same way as for the half-span model, i.e., by equation ( 17).The objective function is also defined similarly to that for the half-span model, i.e., as the minimization of the frequency sum, expressed here as The output quantities include the values of the design variables (K δ1 , K δ2 , (GI k ) 1 , (GI k ) 2 ), for which the flutter speed is equal to the specified certification speed and the three specified frequency ratios are equal to the specified target values.The described procedure is then repeated for several CFR values, typically ranging from 1.05 to 2.0, to obtain enough points to construct a stability margin curve.The notes regarding the downwash effect and the mode switches mentioned in Section 5.1 are also valid for the full-span model.

Reference model and parameters of analyses
The described methodology is demonstrated on the example of the new Czech twin wingmounted engine commuter aircraft for 19 passengers with the maximal take-off weight of 7 000 kg.For the purpose of demonstrating the method, the simplified model with four degrees of freedom, which are represented by both the symmetric and antisymmetric engine pitch and yaw vibrations, was used.The residual structure was considered to be rigid and the control surface and tab actuation drives were blocked.The aerodynamic model included only the wing, nacelles and wing-tip tanks, as the aerodynamics of the fuselage and tail surfaces have negligible effects on the whirl flutter phenomenon.The wings were modelled as Doublet-Lattice Panels, and the nacelles and wing-tip tanks were modelled as Slender and Interference Bodies.The aerodynamic model also included the correction factors for the propeller slipstream, applied to the appropriate aerodynamic elements of the wings and nacelles and the correction factors for the aerodynamic forces and moments at the nose parts of the control surfaces.Values of VFR and HFR were considered at three levels (1.00, 1.05, 1.10) and no downwash effect was considered.

Results -symmetrical revolutions, identical directions
For the case of identical directions of propeller revolution (CW-CW), two mechanisms of the whirl flutter appear: 1) a combination of symmetric pitch and antisymmetric yaw modes (SΘ/AΨ) and 2) a combination of antisymmetric pitch and symmetric yaw modes (AΘ/SΨ).Fig. 7 presents the resulting stability margins with respect to both mechanisms of whirl flutter.The required engine pitch and yaw frequencies are higher for the former mechanism (SΘ/AΨ).Thus, this mechanism of whirl flutter is more critical than the latter one.Comparing the curves for various VFRs and HFRs, we can see the negligible influence of VFR and HFR on the stability margin; thus, the influence of the remaining modes on the stability is negligible as well.Note that the labels "stability" and "flutter" denote the areas of stability (right from the margin curve) and flutter (left from the margin curve), respectively.Fig. 8 shows an example of a V-g-f (velocity -damping -frequency) diagram calculated by the standard approach.A flutter state is represented by the crossing of a V-g curve from negative to positive damping values.The mode No. 1 (SΘ-mode) crossing represents the primary flutter mechanism (SΘ/AΨ).The flutter speed is equal to the certification speed.The mode No. 2 (AΘ-mode) crossing represents the secondary flutter mechanism (AΘ/SΨ).The flutter speed is greater than the certification speed, i.e., above the stability margin, as this type of instability is less critical.

Results -symmetrical revolutions, inverse directions
In the case of inverse directions of propeller revolutions (CW-CCW), the character of the whirl flutter is different.The instability is caused by the combination of antisymmetric pitch and antisymmetric yaw modes (AΘ/AΨ).Stability margins are shown in Fig. 9. Compare to the previous case, the required engine pitch and engine yaw frequencies are considerably higher.Furthermore, this type of instability is influenced also by the remaining modes (SΘ and SΨ) as the influence of VFR and HFR on the stability margin is remarkable.Fig. 10 shows the example of a V-g-f diagram.The crossing of the zero damping line by the V-g curve of the mode No. 2 (AΘ mode) represents the flutter state.The flutter speed is here equal to the certification speed.

Results -a single engine failure and propeller feathering
A single engine failure and a single propeller feathering are modelled by omitting the dynamic and aerodynamic effects of a single rotor from the analysis.Fig. 11 shows the stability margins for both starboard engine (Ω L = 2 080 rpm, Ω R = 0 rpm) and port engine (Ω L = 0 rpm, Ω R = 2 080 rpm) omission.Flutter is caused by the combination of antisymmetric pitch and the symmetric yaw modes (AΘ/SΨ).In this case, the dynamic and aerodynamic effects of the rotating mass are considered on a single side only; therefore, flutter is caused by those pitch and yaw modes, which are closest in terms of frequency.Note that the order of the modes is assumed to be SΘ, AΘ, SΨ, and AΨ.Fig. 11 demonstrates the minimal difference between the cases of omission of the starboard rotor and omission of the port rotor.

Results -unsymmetrical revolutions, identical directions
The cases of unsymmetrical revolutions, including the reduction or increase in the revolutions of the starboard rotor by 15 % (Ω L = 2 080 rpm, Ω R = 1 768 rpm and Ω L = 2 080 rpm, Ω R = 2 392 rpm) for the identical directions of rotation of both propellers (CW-CW), are considered here.The flutter behaviour of the system is similar to that of the corresponding case of the symmetrical revolutions (see Section 6.2).A reduction in revolutions makes the system more stable, while an increase in revolutions makes the system less stable.Fig. 12 shows the stability margins for the reduced revolutions, while Fig. 13 shows the margins for the increased revolutions.Both figures include the primary flutter mechanism (SΘ/AΨ) only.

Results -unsymmetrical revolutions, inverse directions
The cases of unsymmetrical revolutions including the reduction or increase in the revolutions of the starboard rotor by 15 % (Ω L = 2 080 rpm, Ω R = 1 768 rpm and Ω L = 2 080 rpm, Ω R = 2 392 rpm), for the inverse directions of rotation of both propellers (CW-CCW) are considered here.Compared to the previous cases, the flutter behaviour is similar to the corresponding case of the symmetrical revolutions (see Section 6.3).Reduction of the revolutions makes the system more stable, while increase in revolutions makes the system less stable.Fig. 14 shows the stability margins for the reduced revolutions, while Fig. 15 shows the margins for the increased revolutions.

Conclusion
This paper presents an optimisation-based approach to whirl flutter analysis and the application of this method to the full-span model of the aircraft structure.The necessary adjustments of a stick computational model to make it applicable to the full-span case and the modification of the optimisation solution are described.The methodology is demonstrated on the reference model of a twin-engine turboprop commuter aircraft.Presented analyses include the cases of symmetric propeller revolution for both propeller rotation directions, i.e., CW-CW and CW-CCW, the cases with the absence of single propeller rotation representing a failure case of a single propeller feathering and the cases of unsymmetrical revolutions representing the failure cases of a single propeller, e.g., a single propeller overspeed.The typical flutter mechanisms are presented using stability margins.The most critical case is the case of inverse directions of propeller revolution (CW-CCW), for which the critical flutter modes are engine antisymmetric pitch and antisymmetric yaw.Future work will be focused on the assessment of the downwash effect, which represents the aerodynamic interference between propeller, nacelle and wing.

Fig. 2 .
Fig. 2. Stable (a) and unstable (b) states of gyroscopic vibrations for the backward flutter mode