CN105260492A - Synthetic modeling method of rotor and airframe coupling dynamics modality - Google Patents

Abstract

The invention discloses a synthetic modeling method of rotor and airframe coupling dynamics modality, and belongs to the dynamic design technology of helicopters. The synthetic modeling method of the rotor and airframe coupling dynamics modality is suitable for the dynamic stability and dynamic response analysis modeling of single-rotor and multi-rotor coupling. The synthetic modeling method analyzes the air resonance instability of a rotor blade by aiming at the rotor blade operation pneumatic environment of a flight state stipulated by a flight envelope rule which needs to considers to cover and properly exceed the performance requirements of the helicopter, researches the influence of high-speed high-overload flight state parameters, looks for the structural design parameters of rotors, blades and a steering control system which cause or affect the generation of the air resonance, and determines the instability boundary of the air resonance. Through the sensitivity analysis of the air resonance by the structural design parameters, a structural design parameter range which affects the air resonance is determined, and a design measure for expanding the security boundary of the air resonance is put forward.

Description

A kind of rotor and body Coupled Dynamics modal synthesis modeling method

Technical field

The invention belongs to Helicopter Dynamics designing technique, particularly a kind of rotor and body Coupled Dynamics modal synthesis modeling method.

Background technology

The rotor stability Design analysis that is coupled with body occupies very important status in Helicopter Dynamics designs, its dynamics directly has influence on the safety of helicopter, performance and quality, be be related to helicopter model to be developed into the gordian technique lost, new configuration rotor (as bearingless rotor) and the research and development of multirotor helicopter model in recent years become the new direction of China's Helicopter Development.Bearingless rotor uses flexible beam to replace the level hinge of traditional articulated rotor propeller hub, vertical hinge and axial hinge, because flexible beam distortion is complicated, blade flapping, shimmy and displacement sports coupling are comparatively strong, make such lifting airscrew and body analysis of coupled system's stability must carry out accurate modeling to the power transmission of multichannel propeller hub and strong coupling yielding flexibility beam; More complicated than single-rotor helicopter of the instability that exists and many rotors are coupled with body.Therefore, the new configuration rotor Stability Modeling analytical technology that is coupled with body is one of key breaking through new configuration rotor and multirotor helicopter designing technique.

In prior art, the articulated rotor lifting airscrew Stability Modeling analysis and design technology that to be coupled with body is grasped comprehensively, but did not carry out bearingless rotor and the research of multirotor helicopter gruond and air resonance.Domesticly at present grind and be all radial type grinding helicopter model, blade and body are all worked as rigid body process by analytical model, do not study multirotor helicopter rotor to be coupled with body the modeling and analysis methods of stability and the body Dynamic Characteristics Test method required by analytical approach.Therefore, for adapting to new configuration heligyro research and development demand, carry out rotor to be coupled with body stability modal synthesis Modeling Method, the various new configuration rotor hub of accurate simulation, blade, helicopter body rigid motion and elastic deformation, grasp new configuration or multirotor helicopter coupling stability Design analytical technology, meet China's Helicopter Development new demand.

Summary of the invention

For the new configuration of modern helicopters rotor and many rotors development trend, need to solve rotor and body Coupled Dynamics design analysis technical matters, propose a kind of Modal Synthesis Technique that adopts and the method for modal synthesis modeling is carried out to rotor, this many-body dynamics coupled system of body, application Hamilton's principle derivation rotor and body the coupled dynamical equation, set up rotor and body Coupled Dynamics modal synthesis analytical model, mainly comprise:

(1) build complete description body, rotor hub, blade motion coordinate system and body, rotor blade and aerodynamic model;

(2) body and rotor blade and control structure system finite element method kinetic model is built, for the structural dynamic characteristics of the various rotor of accurate simulation and blade Configuration Design;

(3) set up rotor and body modal coupling kinetic model, adopt modular modeling conception, be suitable for assembling and the choice of many rotors, many blades and control system Coupling Dynamic Model;

(4) what derivation body foundation motion produced involves inertia and the load such as pneumatic, and being loaded on for rotor body derives body convected motion and express in blade inertia and the load such as pneumatic with the remaining nodal displacement of minor structure;

(5) unsteady aerodynamic model that single order and higher order dynamics become a mandarin is built, for stability analysis and dynamic response analysis are selected;

(6) the ONERA aerodynamic model of profile lift, profile drag and the moment coefficient response characteristic in unsteady flow environment is calculated based on aerofoil profile aerodynamic characteristics tests data, for calculating the aerodynamic characteristic of aerofoil profile in unsteady flow environment;

(7) build blade trailing edge winglet deflection Controlling model, to be coupled with body the model that stability controls and vibratory response or load control for lifting airscrew;

(8) build and realize lifting airscrew and be coupled with body the pitch control system model that stability controls and vibratory response or load control, for analyzing fly control and ACFS system to blade pitch control effect;

(9) set up complicated high-order nonlinear kinetics equation in matrix representation mode, set up the wing and body Coupled Dynamics modal synthesis model.

Preferably, rotor and the modeling of body Coupling Dynamics Analysis comprise and carry out Structural Dynamics modeling, Unsteady Aerodynamic Modeling and modal synthesis to rotor, body.This power system is assembled by multiple structure, adopts many-body dynamics method, set up coordinate system and transformational relation, describe the motion of each structure, adopt modeling method respectively to set up the model via dynamical response of rotor, body.First, isolated rotor blade and housing construction Dynamics Finite Element Model is set up respectively; They are minor structure kinetic models of whole analytical model, according to be concerned about rank number of mode and scope, choose lower mode and carry out comprehensively from substructure mode, then apply Modal Synthesis Technique and set up rotor and body coupled mode Comprehensive Analysis Model of Unit.

In above-mentioned either a program preferably, rotor and body the coupled dynamical equation are derived and are adopted Hamilton's principle derivation rotor and body the coupled dynamical equation.First body mode motion kinetic energy, potential energy variation and external force virtual work is derived, derive the deformation energy of blade, the variation of motion energy and external force (aerodynamic force) virtual work again, the variation of cumulative all rotor blade deformation energys, kinetic energy and external force virtual work, substitute into variation equation:

${\∫}_{t_1}^{t_2}(\mathrm{\δU}-\mathrm{\δT}-\mathrm{\δW})\mathrm{dt}=0---\left(1\right)$

By separating housing construction finite element kinetics equation, obtain body eigenfrequncies and vibration models, with selected body Mode Shape matrix, kinetics equation is transformed to Modal Space, obtain the kinetic energy of body mode motion and potential energy variation and external force virtual work, substitute into equation (1) and derive NF independently body mode and the rotor blade coupled vibrations differential equation.

By the variation of all rotor blade deformation energys, kinetic energy and external force virtual work, substitute into variation equation (1), derive rotor blade and body modal coupling kinetics equation.With independently blade movement node displacement { q _iediscrete blade motion { X _rsin space and temporal correlativity, obtain the nonlinear dynamical equation of corresponding blade joint movements.Again through blade trim and linearization process, obtain the rotor blade based on trim displacement and body modal coupling linear dynamics equation.By analyzing the eigenfrequncies and vibration models of paddle blade structure, selected NP rank blade mode, rotor blade and body modal coupling kinetics equation are transformed to the Modal Space of body and blade, decoupling zero is NP rank blade mode and body modal coupling oscillatory differential equation.

NF independently body mode and the rotor blade coupled vibrations differential equation, and the NP rank blade mode of N sheet blade and body modal coupling oscillatory differential equation, form rotor and to be coupled with body stability modal synthesis kinetics equation.

In above-mentioned either a program preferably, define the description body of complete set, rotor hub, blade motion coordinate system and body, rotor blade and aerodynamic model.Adopt coordinate system based on body, rotor hub rotation, blade motion deformation are the multi-body dynamics modeling method of relative coordinate system, and the coordinate transformation relation of foundation makes the coupled motions of the multiple rotor of description and blade and body simply clear and definite.

In above-mentioned either a program preferably, finite element method is adopted to set up body and rotor blade and control structure system dynamics model, can the structural dynamic characteristics of the various rotor of accurate simulation and blade Configuration Design.Eliminate conventional rotor to be coupled with body in Stability Model and to suppose body and blade rigid body, not only consider rigid motion, and consideration elastic movement, improve stability analysis precision, be applicable to radial type, partly cut with scissors formula, hingeless formula bearing-free formula and conformational rigidity rotor.Kinetic model was both applicable to stability analysis, was applicable to again dynamic response analysis,

In above-mentioned either a program preferably, adopt Modal Synthesis Technique, on the basis independently setting up and analyze body and rotor blade Structural Dynamics characteristics of mode, propose the method setting up rotor and body modal coupling kinetic model, reduce analytical scale, priorly there is provided modular modeling conception, be particularly suitable for assembling and the choice of many rotors, many blades and control system Coupling Dynamic Model.

In above-mentioned either a program preferably, what derivation body foundation motion produced involves inertia and the bearing method such as pneumatic.The feature of single-point on body is installed on according to a width rotor, describing the body foundation motion convected motion to rotor blade, the inertia that produces and the load such as pneumatic of deriving time, first adopt the remaining nodal displacement of minor structure to derive, be finally transformed into body Modal Space, simplify equation and deduce complexity.

In above-mentioned either a program preferably, single order and higher order dynamics is adopted to become a mandarin unsteady aerodynamic model, for stability analysis and dynamic response analysis supply a model.Higher order dynamics enters the introducing of flow model, and the expectation not only making rotor induce to become a mandarin is accurate, more simplifies (as free wake) calculating and induces the complicacy become a mandarin.

In above-mentioned either a program preferably, have employed the ONERA aerodynamic model calculating profile lift, profile drag and the moment coefficient response characteristic in unsteady flow environment based on aerofoil profile aerodynamic characteristics tests data, consider comperssibility influence, include dynamic stall, and dynamic inflow model combines aerodynamics evaluation can be made more accurate.

In above-mentioned either a program preferably, establish and to be coupled with body the model that stability controls and vibratory response or load control with the blade trailing edge winglet yaw motion lifting airscrew that is control variable, can by design analysis trailing edge winglet the exhibition of every sheet blade to position, length and width change, to provide initiatively or Passive Control rotor is coupled with body the Research approach of stability and vibratory response or load.

In above-mentioned either a program preferably, establish and utilize pitch control system to realize lifting airscrew to be coupled with body the model that stability controls and vibratory response or load control, simultaneously can considering that ACFS system is to blade pitch control effect, providing technological approaches for analyzing under various state of flight the rotor full bridge aero-elastic model (Air Resonance) that to be coupled with body.

In above-mentioned either a program preferably, the method for complicated high-order nonlinear kinetics equation of deriving in matrix representation mode is established.Do not need the expression formula of complexity to launch in derivation, only need derived relation formula: be added, phase multiplication and division or invert, differentiate, the operation relation such as integration, effectively utilize computer programming to realize formation and the calculating of the whole matrix of coefficients of kinetics equation.

Key point of the present invention is:

For the new configuration of lifting airscrew and many rotors development trend, propose a kind of Modal Synthesis Technique that adopts to rotor, this many-body dynamics coupled system of body carries out the method for comprehensive modeling, application Hamilton's principle derives rotor and body the coupled dynamical equation, the rotor set up and body Coupled Dynamics modal synthesis analytical model, introduce winglet and pitch control design parameter, and AFCS system, not only can be directly used in rotor and body analysis of coupled system's stability, the Helicopter Dynamics design analysiss such as helicopter vibration control and load expectation, rotor can also be used for be coupled with body stability Active Control Design and the research of rotor vibratory load ACTIVE CONTROL.

Described body modal synthesis kinetic model.

Described rotor blade modal synthesis kinetic model.

Described rotor and body modal coupling comprehensive dynamic model

Described winglet controls Aerodynamic Model

Described pitch control (or ACFS) model

Described rotor and body coupled mode comprehensive dynamic equation are derived.

Beneficial effect of the present invention: this modeling technique adopts Finite Element Method to consider body and blade rigid motion and elastic deformation, can the various rotor configuration of accurate simulation and power transmission, adapts to following and new configuration rotor research and development dynamics Design analysis demand.Utilization Modal Synthesis Technique sets up rotor and body coupled mode analytical model had both reduced kinetics equation scale, highlights the stability of paid close attention to coupled mode; Consider that elastic deformation eliminates rigidity hypothesis, make analytical model more accurately, more meet engineering reality, reduce design analysis error.Meanwhile, can need to select rotor blade and body mode according to dynamics Design analysis, rotor and the body coupled mode analytical model of foundation can be directly used in the Helicopter Dynamics design analysis such as helicopter vibration control, load expectation.And introduce winglet and pitch control input parameter, establish the rotor and body analysis of coupled system's stability model of considering AFCS systematic influence, be coupled with body stability Active Control Design and rotor vibratory load of rotor can be used for and control to study.Therefore, this modeling method defines a set of various configuration heligyro rotor and body Coupled Dynamics modeling technique, is the technological means that following new model, new configuration rotor aerodynamics design analysis provide reliably, are suitable for.And be used successfully to in-service and developing model dynamics Design analysis, obtained model checking, achieve substantial economics.

Accompanying drawing explanation

Fig. 1 is body according to a preferred embodiment of multichannel power transmission rotor hub Structural Dynamics modeling method of the present invention and rotor hub coordinate system model schematic.

Fig. 2 is the present invention's rotor hub embodiment illustrated in fig. 1 and blade coordinate system schematic diagram.

Fig. 3 is the present invention's blade non-deformation coordinate system model schematic diagram embodiment illustrated in fig. 1.

Fig. 4 is the present invention's blade distortion coordinate system model schematic embodiment illustrated in fig. 1.

Fig. 5 is the present invention's propeller-blade section coordinate system model schematic embodiment illustrated in fig. 1.

Fig. 6 is the present invention's blade embodiment illustrated in fig. 1 and winglet Aerodynamic Model schematic diagram.

Fig. 7 is the present invention's pitch control (or ACFS) model schematic embodiment illustrated in fig. 1.

Embodiment

Below in conjunction with accompanying drawing, multichannel power transmission rotor hub Structural Dynamics modeling method involved in the present invention is described in further details.

The first step: definition coordinate system and transformational relation

(1) body and rotor hub coordinate system is defined shown in Fig. 1.{ Og, Xg, Yg, Zg} are earth-fixed axis systems, and in order to use, { ig, jg, kg} represent coordinate vector, and aircraft gravity is along-kg direction.{ O _f, X _f, Y _f, Z _fbody axis system, true origin Of at full machine center of gravity place, coordinate vector { i _f, j _f, k _frepresent, relative to coordinate system, { 6 freedoms of motion of ig, jg, kg} are X to body _f, Y _f, Z _f, x _fjust be forward, Z _fjust be upwards, Y _fforward is determined by the right-hand rule.

{ O _h, X _h, Y _h, Z _hrotor hub coordinate system, be arranged on body { X for describing any one _f, Y _f, Z _fthe motion of rotor at coordinate place, its true origin is at propeller hub center, and coordinate vector is { i _h, j _h, k _h.When having multiple rotor, { O is numbered to it _h, X _hi, Y _hi, Z _hi(i=1,2 ..., IR), rotor shaft k _hrelative axis k _flean forward an angle α _h, it is for simulating rotor shaft forward tilted positions, given by overall design.The motion X at rotor hub center _h, Y _h, Z _h, describe, it is the local coordinate at rotor hub place on body, is convenient to the motion describing propeller hub.For tiltrotor, rotor hub is connected with oar axle and adopts universal hinge form, and propeller hub is also connected with body with damping element (non rotating) by flexible member, to realize the manipulation to rotor, is similar to the auto-bank unit of helicopter.Propeller hub these motions in housing construction are considered in finite element model, choose body mode when carrying out comprehensive, consider propeller hub motion these local modes in full machine vibration.

(2) definition rotor hub and blade coordinate system shown in Fig. 2, illustrate the position of a slice blade in propeller hub coordinate and rotating coordinate system { ik, jk, the kk} of blade.Blade rotates to (overlooking) counterclockwise with the rotating speed of Ω, and position angle ψ k is the angle of kth sheet blade apart from propeller hub coordinate axis iH negative sense.When blade clockwise direction (overlooking) rotates, position angle is-ψ k, and the force and moment R direction acted on blade is included in Conversion Matrix of Coordinate [Φ k].

(3) coordinate system before the distortion of the blade of definition shown in Fig. 3 and the coordinate system after distortion.Coordinate system { rotating coordinate system { ik, the jk of the relative blade of is, js, ks} before blade distortion, kk} has a pre-cone angle beta C, and blade coordinate origin is selected in propeller shank apart from EH place, propeller hub center, and on blade, any point is at coordinate system { is, js, with coordinate, { Xs, Ys, Zs} describe in position in ks}.After blade distortion the coordinate of arbitrary section be ib, jb, kb}, and the point in blade resilience axes at the deformation displacement u at this section place, v, w, represent, relative to the coordinate system before blade distortion, { is, js, ks} have rotated three angles these three angles give geometric meaning in figs. 4-6.

(4) propeller-blade section coordinate system is defined shown in Fig. 4.η and ζ is the coordinate axis of propeller-blade section.When η and ζ represents the elastic shaft of propeller-blade section, θ just represents the corner of blade elastic shaft relative coordinate axle iS; When η and ζ represents the axes of inertia of propeller-blade section, θ just represents the corner of propeller-blade section axes of inertia relative coordinate axle iS; When η and ζ represents propeller-blade section aerofoil profile axle, θ just represents the established angle of propeller-blade section relative coordinate axle iS.Blade elastic torsion angle is added from different θ for the calculating of blade elastic load, inertial load and aerodynamic loading.

(5) the transformational relation battle array between coordinate system

Suppose the rolling of body, pitching and crab angle be low-angle, then { ig, jg, kg} are to body coordinate vector { i for ground stationary coordinate vector _f, j _f, k _ftransition matrix [φ g] not relevant with the precedence of corner, can simply be expressed as unique form:

$\begin{array}{cc}\left\{\begin{array}{c}{i}_f\\ j_f\\ k_f\end{array}\right\}=\left[{\mathrm{\φ}}_g\right]\left\{\begin{array}{c}{i}_g\\ j_g\\ k_g\end{array}\right\}& \left[{\mathrm{\φ}}_g\right]=\left[\begin{array}{ccc}1& {\mathrm{\φ}}_{\mathrm{zf}}& -{\mathrm{\φ}}_{\mathrm{yf}}\\ -{\mathrm{\φ}}_{\mathrm{zf}}& 1& {\mathrm{\φ}}_{\mathrm{xf}}\\ {\mathrm{\φ}}_{\mathrm{yf}}& -{\mathrm{\φ}}_{\mathrm{xf}}& 1\end{array}\right]\end{array}---\left(2\right)$

Body coordinate vector { i _f, j _f, k _fto rotor hub coordinate vector { i _h, j _h, k _htransformational relation and transition matrix [α _h] be:

$\begin{array}{cc}\left\{\begin{array}{c}{i}_H\\ j_H\\ k_H\end{array}\right\}=\left[{\mathrm{\α}}_H\right]\left\{\begin{array}{c}{i}_g\\ j_g\\ k_g\end{array}\right\}& \left[{\mathrm{\α}}_H\right]=\left[\begin{array}{ccc}\cos {\mathrm{\α}}_H& 0& -\sin {\mathrm{\α}}_H\\ 0& 1& 0\\ \sin {\mathrm{\α}}_H& 0& \cos {\mathrm{\α}}_H\end{array}\right]\end{array}---\left(3\right)$

Rotor hub vector { i _h, j _h, k _hto blade rotational coordinates vector { transformational relation of ik, jk, kk} and transition matrix for:

$\begin{array}{cc}\left\{\begin{array}{c}{i}_k\\ j_k\\ k_k\end{array}\right\}=\left[{\mathrm{\φ}}_k\right]\left\{\begin{array}{c}{i}_H\\ j_H\\ k_H\end{array}\right\}& \left[{\mathrm{\φ}}_k\right]=\left[\begin{array}{ccc}-\cos {\mathrm{\ψ}}_k& -\sin {\mathrm{\ψ}}_k& 0\\ \sin {\mathrm{\ψ}}_k& -\cos {\mathrm{\ψ}}_k& 0\\ 0& 0& 1\end{array}\right]\end{array}---\left(4\right)$

Blade rotational coordinates vector ik, jk, kk} to blade distortion before coordinate vector the transformational relation of is, js, ks} and transition matrix [β pc] (if plunder angle in advance in addition, can revise formula (5), will plunder angle ξ C in advance) are:

$\begin{array}{cc}\left\{\begin{array}{c}{i}_s\\ j_s\\ k_s\end{array}\right\}=\left[{\mathrm{\β}}_{\mathrm{pC}}\right]\left\{\begin{array}{c}{i}_k\\ j_k\\ k_k\end{array}\right\}& \left[{\mathrm{\β}}_{\mathrm{pC}}\right]=\left[\begin{array}{ccc}\cos {\mathrm{\β}}_C& 0& \sin {\mathrm{\β}}_C\\ 0& 1& 0\\ -\sin {\mathrm{\β}}_C& 0& \cos {\mathrm{\β}}_C\end{array}\right]\end{array}---\left(5\right)$

Blade distortion before coordinate vector is, js, ks} to blade be out of shape after coordinate vector the transformational relation of ib, jb, kb} and transition matrix [T] are:

$\begin{array}{cc}\left\{\begin{array}{c}{i}_b\\ j_b\\ k_b\end{array}\right\}=\left[T\right]\left\{\begin{array}{c}{i}_s\\ j_s\\ k_s\end{array}\right\}& \left[T\right]=\begin{array}{c}\left[\begin{array}{ccc}{c}_{\stackrel{\‾}{\mathrm{\β}}}{c}_{\stackrel{\‾}{\mathrm{\ζ}}}& c_{\stackrel{\‾}{\mathrm{\β}}}{s}_{\stackrel{\‾}{\mathrm{\ζ}}}& s_{\stackrel{\‾}{\mathrm{\β}}}\\ -c_{\stackrel{\‾}{\mathrm{\θ}}}{s}_{\stackrel{\‾}{\mathrm{\ζ}}}-s_{\stackrel{\‾}{\mathrm{\θ}}}{c}_{\stackrel{\‾}{\mathrm{\ζ}}}{s}_{\stackrel{\‾}{\mathrm{\β}}}& c_{\stackrel{\‾}{\mathrm{\θ}}}{c}_{\stackrel{\‾}{\mathrm{\ζ}}}-s_{\stackrel{\‾}{\mathrm{\ζ}}}{s}_{\stackrel{\‾}{\mathrm{\β}}}{s}_{\stackrel{\‾}{\mathrm{\θ}}}& c_{\stackrel{\‾}{\mathrm{\θ}}}{s}_{\stackrel{\‾}{\mathrm{\θ}}}\\ s_{\stackrel{\‾}{\mathrm{\θ}}}{s}_{\stackrel{\‾}{\mathrm{\ζ}}}-c_{\stackrel{\‾}{\mathrm{\θ}}}{s}_{\stackrel{\‾}{\mathrm{\β}}}{c}_{\stackrel{\‾}{\mathrm{\ζ}}}& -s_{\stackrel{\‾}{\mathrm{\θ}}}{c}_{\stackrel{\‾}{\mathrm{\ζ}}}-s_{\stackrel{\‾}{\mathrm{\ζ}}}{s}_{\stackrel{\‾}{\mathrm{\β}}}{c}_{\stackrel{\‾}{\mathrm{\θ}}}& c_{\stackrel{\‾}{\mathrm{\β}}}{c}_{\stackrel{\‾}{\mathrm{\θ}}}\end{array}\right]\end{array}\end{array}---\left(6\right)$

In formula, $\stackrel{\‾}{\mathrm{\ζ}}={\tan }^{-1}\left(\frac{v^{\′}}{\sqrt{1-v^{\′2}-w^{\′2}}}\right),\stackrel{\‾}{\mathrm{\β}}={\tan }^{-1}\left(\frac{w^{\′}}{\sqrt{1-w^{\′2}}}\right),\stackrel{\‾}{\mathrm{\θ}}=\mathrm{\θ}+\mathrm{\φ}-{\∫}_0^x{v}^{\′\′}{w}^{\′}\mathrm{dx}.$

Second step, housing construction kinetic model

Housing construction kinetic model adopts Finite Element Method to set up, business software PATRAN is utilized directly to carry out FEM meshing to each parts of housing construction of design, according to the stressed and power transmission of structural elements, parts, from the cell library that NASTRAN software provides, select suitable unit simulation all kinds structure, set up housing construction Dynamics Finite Element Model.

When aircraft is in ground landing, landing gear provides rigidity and damping, vital to the rigid-body vibration mode of aircraft on landing gear and body thereof the stability that to be coupled with rotor power, therefore, need to carry out modeling to landing gear, simulate rigidity and the damping of X, Y that single landing gear provides and Z-direction with three springs, damping element respectively.

Rotor hub is the interface be connected with body, and propeller hub central point is a node in housing construction finite element model.Similarly, in housing construction finite element model using horizontal tail, vertical fin and tail-rotor Aerodynamic force action point also as a node.For the rotor hub of tiltrotor, its pitch motion describes in body kinetic model, and the local mode that namely propeller hub is formed under controlled mechanism and flexible member constraint puts into body kinetic model.Housing construction is discrete turn to finite element model after form a many-degrees of freedom system, the form of its undamped oscillation differential equation is:

$\left[M\right]\left\{\stackrel{\·\·}{X}\right\}+\left[K\right]\left\{X\right\}=0---\left(7\right)$

Natural frequency [f is calculated with NASTRAN software ₁, f ₂, f ₃..., f _np] and the vibration shape [X _fPM], by the eigenfrequncies and vibration models of analytical engine body structure, determine NF body vibration mode, wherein must comprise the rigid body mode of 6 bodies.With selected body Mode Shape matrix, equation (7) is transformed to Modal Space, decoupling zero is NF the independently modal vibration differential equation, and its form is:

$M_{\mathrm{FPi}}{\stackrel{\·\·}{X}}_{\mathrm{FPi}}+C_{\mathrm{FPi}}{\stackrel{\·}{X}}_{\mathrm{FPi}}+K_{\mathrm{FPi}}{X}_{\mathrm{FPi}}=0,(i=\mathrm{1,2},...,N_{\mathrm{FP}})---\left(8\right)$

Modal transformation relational expression is used: { X}=[X in formula _fPM] { X _fP, modal damping battle array is made up of structural damping and damping element, obtains through modal coordinate conversion.

3rd step, rotor blade kinetic model

Blade kinetic model describes.As shown in Figure 3 and Figure 4, blade is reduced to an elongated elastic beam to rotor blade kinetic model, supposes that blade elastic shaft passes through rotation center; Propeller-blade section barycenter, center of tension and aerodynamic center can not overlap, quality, the tensible rigidity distribution of section are not in relation to section chordwise axes η symmetry, propeller-blade section relative resilient axle along exhibition to pre-torsional angle be θ (Fig. 5), the relative surfaces of revolution of blade has a pre-cone angle beta _c, consider that blade deformation geometry is non-linear.Propeller-blade section mass inertia axle relative to coordinate axis corner be θ _i.Blade comprises along exhibition to the distortion of arbitrary section r place elastic shaft: axial displacement u, shimmy to displacement v, wave to displacement w and torsional deflection four motion and structure and inertia coupling, these displacements are relative propeller hub coordinate systems.The connection of propeller shank and propeller hub is by boundary condition and adopt boundary element to simulate its hinged, clamped or hinged band spring constraint and multichannel force-transmitting relation.

4th step, the existing body Coupled Dynamics modal synthesis model of rotor

Export body vibration modal parameter with body dynamics model analysis, these modal parameters are the amounts being transformed into Modal Space through modal coordinate, according to hyperelement definition, body modal coordinate { X _fPbe defined as the supernode coordinate of body minor structure.The linkage interface of body and rotor blade is propeller hub, and the freedom of motion of propeller hub is the nodal displacement of body minor structure and rotor blade minor structure linkage interface, includes the convected motion of propeller hub in rotor blade kinetic model.Therefore, rotor and body coupled mode Comprehensive Analysis Model of Unit are defined as the supernode of body and rotor blade minor structure respectively the body of Modal Space and isolated rotor blade mode of oscillation, linkage interface propeller hub is defined as residual texture, and the displacement of propeller hub is defined as the displacement of remaining node.

Setting up in rotor and body coupled mode Comprehensive Analysis Model of Unit process, first using blade as residual texture, linkage interface propeller hub node is as remaining node, set up the Dynamics Coupling analytical model that Modal Space mixes with physical space, this model is convenient to derivation rotor and body the coupled dynamical equation.On this basis, by finite element method discretize blade dynamics partial differential equation, calculate the vibration characteristics of isolated blade, choose N _pindividual blade mode, then blade motion transform to blade Modal Space, the nodal displacement modal coordinate { X of propeller hub _fPrepresent, obtain rotor and body coupled mode Comprehensive Analysis Model of Unit.

5th step, Aerodynamic Model

In order to be applicable to rotor/body coupling power stability and the Dynamic Response, on the basis adopting quasi-steady aerodynamic force model, the impact on stability and response of unsteady aerodynamic force and dynamic stall considered by rotor blade aerodynamic power model.The non-homogeneous single order dynamic inflow model of people's researchs such as D.M.Pitt and D.A.Peters is used for the dynamic inflow model of analysis of coupled system's stability model and the aero-elastic response analytical model that is coupled employing high-order.New configuration rotor blade is along exhibition to the aerofoil profile all adopting different aerodynamic characteristic, and the aerodynamic characteristic data of aerofoil profile are desirable from drying test.Fly enter dynamic stall with lift state blade upper part aerofoil profile front in addition, therefore, adopt with the ONERA model of the aerodynamic characteristics tests data of aerofoil profile, the dynamic response of rotor power system can be calculated more exactly, especially transient response.Lifting airscrew calculates the aerodynamic loading of fuselage (comprising wing, horizontal tail, vertical fin and tail-rotor) each parts and adopts test blowing data.

(1) blade quasi-steady aerodynamic force model.Fig. 6 is the aerodynamic force of propeller-blade section effect, considers the motion of trailing edge winglet.The aerodynamic force element acting on section is:

$\begin{array}{c}T=\frac{\mathrm{\ρc}}{2}V(C_L{U}_T-C_D{U}_P)-\frac{\mathrm{\ρ}{a}_{\∞}{c}^2}{8}{\stackrel{\·}{U}}_P\\ \≈\frac{\mathrm{\ρc}}{2}[C_{L0}{U}_T^2-a_{\∞}{U}_T{U}_P-C_{D0}{U}_T{U}_P+C_{D1}{U}_P^2]-\frac{\mathrm{\ρ}{a}_{\∞}{c}^2}{8}{\stackrel{\·}{U}}_P\\ =\frac{\mathrm{\ρc}}{2}[C_{L0}{U}_T^2-(a_{\∞}+C_{D0})U_T{U}_P+C_{D1}{U}_P^2]-\frac{\mathrm{\ρ}{a}_{\∞}{c}^2}{8}{\stackrel{\·}{U}}_P\end{array}---\left(9\right)$

$\begin{array}{c}S=-\frac{\mathrm{\ρc}}{2}V(C_L{U}_P+C_D{U}_T)\\ \≈-\frac{\mathrm{\ρc}}{2}[C_{L0}{U}_T{U}_P-a_{\∞}{U}_P^2+C_{D0}{U}_T^2{-C}_{D1}{U}_T{U}_P+C_{D2}{U}_P^2]\\ =-\frac{\mathrm{\ρc}}{2}[C_{L0}{U}_T^2+(C_{L0}-C_{D1})U_T{U}_P-(a_{\∞}-C_{D2})U_P^2]\end{array}---\left(10\right)$

When blade opens up the airfoil trailing edge band winglet to r place section, the aerodynamic loading on this section will add the contribution of winglet, and the Additional pneumatic load produced relative to the deflection angle δ of aerofoil profile chordwise axes Yb (η) by winglet is:

$T_{\mathrm{flap}}\≈\mathrm{\ρ}{\mathrm{cU}}_T(U_T{f}_1\mathrm{\δ}+\frac{c}{4}{f}_2\stackrel{\·}{\mathrm{\δ}})-\frac{\mathrm{\ρ}{c}^2}{4}{f}_3{U}_T\stackrel{\·}{\mathrm{\δ}}---\left(12\right)$

$\begin{array}{c}{M}_{\mathrm{flap}}\≈-\mathrm{\ρ}{c}^2{U}_T(U_T{f}_1\mathrm{\δ}+\frac{c}{4}{f}_2\stackrel{\·}{\mathrm{\δ}}){\stackrel{\‾}{e}}_{\mathrm{ac}}-\frac{\mathrm{\ρ}{c}^2}{4}{U}_T^2(f_1+f_3)\mathrm{\δ}\\ -\frac{\mathrm{\ρ}{c}^3}{4}{U}_T[T_1-T_8-(1.5+2{\stackrel{\‾}{e}}_{\mathrm{ac}}-2{\stackrel{\‾}{c}}_f)f_3+\frac{f_2}{3}]\stackrel{\·}{\mathrm{\δ}}\end{array}---\left(13\right)$

(2) blade unsteady aerodynamic model---single order dynamic inflow.Rotor shows for following form at a linear order harmonics distribution table of the induced velocity of oar Pan Chu:

${\stackrel{\‾}{v}}_i(r,t)=v_{\mathrm{ie}}\left(t\right)+{\mathrm{rv}}_{\mathrm{ic}}\left(t\right)\cos {\mathrm{\ψ}}_k+{\mathrm{rv}}_{\mathrm{is}}\left(t\right)\sin {\mathrm{\ψ}}_k---\left(14\right)$

In formula, dimensionless induced velocity, v _ie, v _ic, v _isdetermined by following kinetics equation:

$\left[\mathrm{\τ}\right]{\left[\begin{array}{ccc}{\stackrel{\·}{v}}_{\mathrm{ie}}& {\stackrel{\·}{v}}_{\mathrm{ic}}& {\stackrel{\·}{v}}_{\mathrm{is}}\end{array}\right]}^T+{\left[\begin{array}{ccc}{v}_{\mathrm{ie}}& v_{\mathrm{ic}}& v_{\mathrm{is}}\end{array}\right]}^T=\left[L\right]{\left[\begin{array}{ccc}{C}_T& C_{\mathrm{M\φx}}& C_{\mathrm{M\φy}}\end{array}\right]}^T---\left(15\right)$

(3) blade unsteady aerodynamic model---higher order dynamics becomes a mandarin.

${\stackrel{\‾}{v}}_i(r,t)=\underset{j=0}{\overset{\mathrm{JN}}{\mathrm{\Σ}}}\underset{l=j+1,j+3,...}{\overset{\mathrm{LN}}{\mathrm{\Σ}}}{P}_l^j\left(v\right)[v_{\mathrm{lc}}^j\left(t\right)\cos \left(j{\mathrm{\ψ}}_k\right)+v_{\mathrm{ls}}^j\left(t\right)\sin \left(j{\mathrm{\ψ}}_k\right)]---\left(16\right)$

$\left[M\right]\left\{\stackrel{\·}{v_{\mathrm{lc}}^j}\right\}+\left[L^C\right]\left\{v_{\mathrm{lc}}^j\right\}=\frac{1}{2}\left\{{\mathrm{\τ}}_n^{\mathrm{mc}}\right\},\left[M\right]\left\{\stackrel{\·}{v_{\mathrm{ls}}^j}\right\}+\left[L^S\right]\left\{v_{\mathrm{ls}}^j\right\}=\frac{1}{2}\left\{{\mathrm{\τ}}_n^{\mathrm{ms}}\right\}---\left(17\right)$

$\begin{array}{c}{L}_{\ln }^{0\mathrm{mc}}=\frac{1}{2\mathrm{\π}}{\∫}_0^{2\mathrm{\π}}{\∫}_0^1{P}_l^0\left(v_0\right)\frac{\∂}{\∂z}\left[P_n^m\left(v\right)Q_n^m\left(\mathrm{i\η}\right)\right]\cos \left(\mathrm{m\ψ}\right)\mathrm{d\ξd}{v}_0\mathrm{d\ψ}\\ L_{\ln }^{\mathrm{jmc}}=\frac{1}{2\mathrm{\π}}{\∫}_0^{2\mathrm{\π}}{\∫}_0^1{P}_l^j\left(v_0\right)\cos \left(r{\mathrm{\ψ}}_k\right){\∫}_0^{\∞}\frac{\∂}{\∂z}\left[P_n^m\left(v\right)Q_n^m\left(\mathrm{i\η}\right)\right]\cos \left(\mathrm{m\ψ}\right)\mathrm{d\ξd}{v}_0\mathrm{d\ψ}\\ L_{\ln }^{\mathrm{jms}}=\frac{1}{2\mathrm{\π}}{\∫}_0^{2\mathrm{\π}}{\∫}_0^1{P}_l^j\left(v_0\right)\cos \left(r{\mathrm{\ψ}}_k\right){\∫}_0^{\∞}\frac{\∂}{\∂z}\left[P_n^m\left(v\right)Q_n^m\left(\mathrm{i\η}\right)\right]\sin \left(\mathrm{m\ψ}\right)\mathrm{d\ξd}{v}_0\mathrm{d\ψ}\end{array}---\left(18\right)$

$\begin{array}{c}{\mathrm{\τ}}_n^{0c}=\frac{1}{\mathrm{\π}}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\∫}_0^1\frac{F_w^k}{\mathrm{\ρ}{\mathrm{\Ω}}^2{R}^3}{\mathrm{\φ}}_n^0\left(r\right)\mathrm{dr}\\ {\mathrm{\τ}}_n^{\mathrm{mc}}=\frac{1}{\mathrm{\π}}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\∫}_0^1\frac{F_w^k}{\mathrm{\ρ}{\mathrm{\Ω}}^2{R}^3}{\mathrm{\φ}}_n^m\left(r\right)\cos \left(m{\mathrm{\ψ}}_k\right)\mathrm{dr}\\ {\mathrm{\τ}}_n^{\mathrm{ms}}=\frac{1}{\mathrm{\π}}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{\∫}_0^1\frac{F_w^k}{\mathrm{\ρ}{\mathrm{\Ω}}^2{R}^3}{\mathrm{\φ}}_n^m\left(r\right)\cos \left(m{\mathrm{\ψ}}_k\right)\mathrm{dr}\end{array}---\left(19\right)$

(4) ONERA model.ONERA model is the aerodynamic model calculating profile lift, profile drag and the moment coefficient response characteristic in unsteady flow environment based on aerofoil profile aerodynamic characteristics tests data, it can consider comperssibility influence, and include dynamic stall, and dynamic inflow model combines aerodynamics evaluation can be made more accurate.Also do not have article report dynamic inflow to combine with ONERA model at present and calculate the research of unsteady aerodynamic force.This built-up pattern can process non-linear, for rotor/body coupling transient response, dynamic stability analysis.

The profile lift differential equation utilizing ONERA model to set up is:

$\begin{array}{c}{\stackrel{\·}{\mathrm{\Γ}}}_1+\frac{2U}{c}{R}_{\mathrm{\λ}}{\mathrm{\Γ}}_1=\frac{2U}{c}{R}_{\mathrm{\λ}}(C_{\mathrm{LSL}}U+\frac{c}{2}\stackrel{\·}{\mathrm{\α}})\\ \stackrel{\·\·}{{\mathrm{\Γ}}_2}+A_{\mathrm{LL}}\frac{2U}{c}{\stackrel{\·}{\mathrm{\Γ}}}_2+R_w^2{\left(\frac{2U}{c}\right)}^2{\mathrm{\Γ}}_2=-R_w^2{\left(\frac{2U}{c}\right)}^2\mathrm{U\Δ}{C}_L\end{array}---\left(20\right)$

The aerofoil profile aerodynamic drag differential equation is:

${\stackrel{\·\·}{{\mathrm{\Γ}}_D}+A}_{\mathrm{LD}}\frac{2U}{c}{\stackrel{\·}{\mathrm{\Γ}}}_D+R_P^2{\left(\frac{2U}{c}\right)}^2{\mathrm{\Γ}}_D=-R_P^2{\left(\frac{2U}{c}\right)}^2\mathrm{U\Δ}{C}_D---(21$

The aerofoil profile aerodynamic moment differential equation is:

${\stackrel{\·\·}{{\mathrm{\Γ}}_m}+A}_{\mathrm{Lm}}\frac{2U}{c}{\stackrel{\·}{\mathrm{\Γ}}}_m+R_P^2{\left(\frac{2U}{c}\right)}^2{\mathrm{\Γ}}_m=-R_P^2{\left(\frac{2U}{c}\right)}^2\mathrm{U\Δ}{C}_m---\left(22\right)$

In formula,

$\begin{array}{cc}{A}_{\mathrm{LL}}=a_{L0}+a_{L2}\mathrm{\Δ}{C}_L^2,& R_w=w_{\mathrm{ds}0}+w_{\mathrm{ds}2}\mathrm{\Δ}{C}_L^2\\ A_{\mathrm{LD}}=a_{D0}+a_{D2}\mathrm{\Δ}{C}_L^2,& A_{\mathrm{Lm}}=a_{m0}+a_{m2}\mathrm{\Δ}{C}_L^2\\ R_p=R_{p0}+R_{p2}\mathrm{\Δ}{C}_L^2,& \end{array}---\left(23\right)$

Δ C _l, Δ C _d, Δ C _mprofile lift, the resistance trial value of force and moment coefficient and the difference of its linear value, R _λ(23) in formula, other coefficient is the parameter of ONERA model, is determined by Airfoil Testing data.Equation (20) the 2nd formula is by Δ C to (22) _l, Δ C _d, Δ C _mthe aerofoil profile aerodynamic characteristic changing unit caused, Δ C _l, Δ C _d, Δ C _mbecome with the angle of attack, the aerofoil profile aerodynamic lift evoked, the forced response of resistance force and moment reflect the dynamic change of aerofoil profile aerodynamic characteristic, and equation (20) the 1st formula is the circular rector under linear not stall conditions.This prescription journey for calculating the pneumatic dynamic perfromance of aerofoil profile in linear zone and inelastic region, as time lag and dynamic stall etc.

6th step, pitch control (or ACFS) model

Conventional configuration helicopter flight control is realized by auto-bank unit, ACFS is automatic flight control system, it is also realize pitch control by controlling auto-bank unit motion, therefore, the pitch input θ by controlling auto-bank unit motion considered by pitch control (or ACFS) model _cC, θ _1Cand θ _1S(total apart from, in length and breadth to feathering).If the pitch at blade 0.75R place is expressed as: θ=θ _cC+ θ _1Ccos ψ _k+ θ _1Ssin ψ _k, pitch control (or ACFS) is input as: φ=φ _cC+ ι _1Ccos ψ _k+ φ _1Ssin ψ _k, the transmission characteristic from control inputs to the pitch of blade 0.75R is: θ _cc(s)=H _cc(s) φ _cc(s), θ _1c(s)=H _1c(s) φ _1c(s), θ _1s(s)=H _1s(s) φ _1s(s).S is Laplace variable, and pitch control input is given or ACFS system output by driver.Fig. 7 shows Controlling model.

7th step, the derivation of rotor blade kinetics equation

(1) a slice blade potential energy of strain.Definition U _brepresent the potential energy of strain of a slice blade, the Nb sheet blade potential energy of strain is NU _b, U _bcan following formula be expressed as:

$\begin{array}{c}{U}_b=\frac{1}{2}{\∫}_0^R{\∫\∫}_A({\mathrm{\σ}}_{x_s{x}_s}{\mathrm{\ϵ}}_{x_s{x}_s}+{\mathrm{\σ}}_{x_s\mathrm{\η}}{\mathrm{\ϵ}}_{x_s\mathrm{\η}}+{\mathrm{\σ}}_{x_s{\mathrm{\ζ}}_s}{\mathrm{\ϵ}}_{x_s\mathrm{\ζ}})\mathrm{d\ηd\ζdr}\\ =\frac{1}{2}{\∫}_0^R{\∫\∫}_A({\mathrm{E\ϵ}}_{x_s{x}_s}{\mathrm{\ϵ}}_{x_s{x}_s}+{\mathrm{G\ϵ}}_{x_s\mathrm{\η}}{\mathrm{\ϵ}}_{x_s\mathrm{\η}}+{\mathrm{G\ϵ}}_{x_s{\mathrm{\ζ}}_s}{\mathrm{\ϵ}}_{x_s\mathrm{\ζ}})\mathrm{d\ηd\ζdr}\end{array}---\left(24\right)$

Will $\begin{array}{ccc}{\mathrm{\σ}}_{x_s{x}_s}={\mathrm{E\ϵ}}_{x_s{x}_s}& {\mathrm{\σ}}_{x_s\mathrm{\η}}=G{\mathrm{\ϵ}}_{x_s\mathrm{\η}}& {\mathrm{\σ}}_{x_s\mathrm{\ζ}}={\mathrm{G\ϵ}}_{x_s\mathrm{\ζ}}\end{array}$ Substitute into above formula and obtain U _bvariation be:

$\mathrm{\δ}{U}_b={\∫}_0^R{\∫\∫}_A({\mathrm{E\ϵ}}_{x_s{x}_s}\mathrm{\δ}{\mathrm{\ϵ}}_{x_s{x}_s}{+\mathrm{G\ϵ}}_{x_s\mathrm{\η}}\mathrm{\δ}{\mathrm{\ϵ}}_{x_s\mathrm{\η}}+G{\mathrm{\ϵ}}_{x_s{\mathrm{\ζ}}_s}\mathrm{\δ}{\mathrm{\ϵ}}_{x_s\mathrm{\ζ}})\mathrm{d\ηd\ζdr}---\left(25\right)$

(2) rotor blade is out of shape total potential energy variation.The total potential energy of distortion of rotor blade is added by whole blade deformation energy, namely obtains rotor blade be out of shape total potential energy variation to its variation.

$\mathrm{\δ}{U}_R=\underset{\mathrm{ir}=1}{\overset{\mathrm{IR}}{\mathrm{\Σ}}}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}\mathrm{\δ}{U}_{\mathrm{bk}}---\left(26\right)$

(3) kinetic energy of a slice blade and variation

On blade along exhibition to any movement velocity matrix representation in Arbitrary Airfoils section:

$\begin{array}{c}\left\{V_s\right\}=\left(\right[\mathrm{Cvho}]+C_{\mathrm{\ψk}}[\mathrm{Cvhc}]+S_{\mathrm{\ψk}}[\mathrm{Cvhs}\left]\right)\left\{{\stackrel{\·}{X}}_{\mathrm{he}}\right\}+\left[\mathrm{Crvs}\right]\left\{{\stackrel{\·}{X}}_{\mathrm{rs}}\right\}+\\ +\mathrm{\Ω}\left[\mathrm{Crds}\right]\left\{X_{\mathrm{rs}}\right\}+\mathrm{\Ω}\left\{V_{50}\right\}\\ \left\{V_s\right\}={\left[\begin{array}{ccc}{V}_{\mathrm{is}}& V_{\mathrm{js}}& V_{\mathrm{ks}}\end{array}\right]}^T,\left\{V_{s0}\right\}={\left[\begin{array}{ccc}{V}_{\mathrm{is}0}& V_{\mathrm{js}0}& V_{\mathrm{ks}0}\end{array}\right]}^T\\ \left\{X_{\mathrm{he}}\right\}={\left[\begin{array}{cccccc}{X}_H& Y_H& Z_H& {\mathrm{\φ}}_{\mathrm{XH}}& {\mathrm{\φ}}_{\mathrm{YH}}& {\mathrm{\φ}}_{\mathrm{ZH}}\end{array}\right]}^T\end{array}---\left(27\right)$

In formula, { X _rs}=[uvw φ v'w'] ^t

Derive a slice blade kinetic energy again:

Blade kinetic energy variation is:

Speed and speed variation expression formula are substituted into:

$\begin{array}{c}{\mathrm{\&delta;T}}_b={\&Integral;}_0^R{\&Integral;\&Integral;}_A\mathrm{\&rho;}<{\left\{\mathrm{\&delta;}{\stackrel{\&CenterDot;}{X}}_{\mathrm{he}}\right\}}^T({\left[\mathrm{Cvho}\right]}^T+C_{\mathrm{\&psi;k}}{\left[\mathrm{Cvhc}\right]}^T+S_{\mathrm{\&psi;k}}{\left[\mathrm{Cvhs}\right]}^T)+{\left\{\mathrm{\&delta;}{\stackrel{\&CenterDot;}{X}}_{\mathrm{rs}}\right\}}^T{\left[\mathrm{Crvs}\right]}^T+\\ +{\left\{\mathrm{\&delta;}{X}_{\mathrm{rs}}\right\}}^T({\left[\mathrm{Crds}\right]}^T+{\left[\mathrm{Dcrds}\right]}^T+[{\mathrm{Dcrvs}}^T]+{\left[\mathrm{Dcvho}\right]}^T+C_{\mathrm{\&psi;k}}{\left[\mathrm{Dcvhc}\right]}^T+S_{\mathrm{\&psi;k}}{\left[\mathrm{Dcvhs}\right]}^T>\&CenterDot;\\ \&CenterDot;<\left(\right[\mathrm{Cvho}]+C_{\mathrm{\&psi;k}}[\mathrm{Cvhc}]+S_{\mathrm{\&psi;k}}[\mathrm{Cvhs}\left]\right)\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{he}}\right\}+\left[\mathrm{Crvs}\right]\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{rs}}\right\}+\left[\mathrm{Crds}\right]\left\{X_{\mathrm{rs}}\right\}>\mathrm{d\&eta;d\&zeta;dr}\end{array}---\left(30\right)$

(4) rotor total kinetic energy and variation

The total kinetic energy variation of a rotor is N number of (28) sum, and the variation of the total kinetic energy of IR rotor can be expressed as:

$\mathrm{\&delta;}{T}_R=\underset{r=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}\mathrm{\&delta;}{T}_H^r+\underset{r=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}\underset{k=1}{\overset{N_b^r}{\mathrm{\&Sigma;}}}\mathrm{\&delta;}{T}_{\mathrm{bk}}^r---\left(31\right)$

(5) propeller hub inertial load

By rotor blade total kinetic energy and variation, derive propeller hub inertial load.A corresponding rotor hub 6 freedom of motion X _heinertia item comprise N sheet blade, it is each synthesis of blade inertial load at propeller hub place.With δ X in formula (30) _hethe item be multiplied is all act on inertial load on propeller hub, and it is the inertial load of a slice (being defined as kth sheet) blade.

(6) blade inertial load

By rotor blade total kinetic energy and variation, derive propeller hub inertial load and blade inertial load.Corresponding a slice blade freedom of motion X _rsinertial load be in (30) formula with δ X _rsthe item be multiplied, is defined as the inertial load of kth sheet blade, and their product is the variation of corresponding a slice blade motion energy.

8th step, the aerodynamic force in dynamic inflow equation is derived

The rotor lift of dynamic inflow equation (15) right-hand member, rolling and pitching moment coefficient are the functions of body mode motion and blade motion, also be the function of induced velocity itself, the dynamic inflow establishing equation dynamic relationship of rotor induced velocity and body and blade sports coupling.The definition of rotor lift, rolling and pitching moment coefficient:

$\begin{array}{cc}{C}_T=\frac{T_R}{\mathrm{\&rho;\&pi;}{R}^2{\left(\mathrm{\&Omega;R}\right)}^2},& T_R=\underset{k=1}{\overset{N_b}{\mathrm{\&Sigma;}}}{\&Integral;}_{r_0}^R{F}_w^k\mathrm{dr}\\ C_{\mathrm{M\&phi;x}}=\frac{M_{\mathrm{\&phi;x}}}{\mathrm{\&rho;\&pi;}{R}^3{\left(\mathrm{\&Omega;R}\right)}^2},& M_{\mathrm{\&phi;x}}=-\underset{k=1}{\overset{N_b}{\mathrm{\&Sigma;}}}{S}_{\mathrm{\&psi;k}}{\&Integral;}_{r_0}^R{F}_w^k{r}_x\mathrm{dr}\\ C_{\mathrm{M\&phi;y}}=\frac{M_{\mathrm{\&phi;y}}}{\mathrm{\&rho;\&pi;}{R}^3{\left(\mathrm{\&Omega;R}\right)}^2},& M_{\mathrm{\&phi;y}}=\underset{k=1}{\overset{N_b}{\mathrm{\&Sigma;}}}{C}_{\mathrm{\&psi;k}}{\&Integral;}_{r_0}^R{F}_w^k{r}_x\mathrm{dr}\end{array},r_x=e_H{C}_{\mathrm{\&beta;c}}+r---\left(32\right)$

9th step, adopts the aerodynamic force of ONERA model to derive

When adopting ONERA model to calculate blade aerodynamic load, will solve the numerical solution of rotor and body power coupling equation, transient state and steady-state response, at this moment do not do linearization process.Therefore, aerodynamic force can represent by brief expression formula propeller hub motion and the blade virtual work done of moving.

Aerodynamic force to the blade virtual work done of moving is:

$\mathrm{\&delta;}{W}_{\mathrm{br}}={\&Integral;}_{r_0}^R{\left\{\begin{array}{c}\mathrm{\&delta;u}\\ \mathrm{\&delta;v}\\ \mathrm{\&delta;w}\\ \mathrm{\&delta;\&phi;}\end{array}\right\}}^T\left\{\begin{array}{c}{F}_u\\ F_v\\ F_w\\ M_{\mathrm{\&phi;a}}\end{array}\right\}\mathrm{dr}={\&Integral;}_{r_0}^R{\left\{\begin{array}{c}\mathrm{\&delta;u}\\ \mathrm{\&delta;v}\\ \mathrm{\&delta;w}\\ \mathrm{\&delta;\&phi;}\end{array}\right\}}^T\left[\begin{array}{ccc}{T}_{21}& T_{31}& 0\\ T_{22}& T_{32}& 0\\ T_{23}& T_{33}& 0\\ 0& 0& 1\end{array}\right]\left\{\begin{array}{c}S\\ T\\ M_{\mathrm{\&phi;a}}\end{array}\right\}\mathrm{dr}---\left(33\right)$

$\left\{\begin{array}{c}S\\ T\\ M_{\mathrm{\&phi;a}}\end{array}\right\}=\frac{\mathrm{\&rho;c}}{2}\left[\begin{array}{ccc}-C_{\mathrm{DD}}& -C_{\mathrm{YD}}& 0\\ C_{\mathrm{YD}}& -C_{\mathrm{DD}}& 0\\ c(C_{\mathrm{MD}}-C_{\mathrm{YD}}{\stackrel{\&OverBar;}{e}}_{\mathrm{ac}})& c{C}_{\mathrm{DD}}{\stackrel{\&OverBar;}{e}}_{\mathrm{ac}}& c{C}_{\mathrm{MD}}\end{array}\right]\left\{\begin{array}{c}{U}_T^2\\ U_T{U}_P\\ U_P^2\end{array}\right\}=\left[\mathrm{CYDM}\right]\left\{\begin{array}{c}{U}_T^2\\ U_T{U}_P\\ U_P^2\end{array}\right\}---\left(34\right)$

Aerodynamic force to the propeller hub virtual work done of moving is:

$W_{\mathrm{br}}={\&Integral;}_{r_0}^R\mathrm{\&delta;}{\left\{X_{\mathrm{he}}\right\}}^T\left[\mathrm{FABP}\right]\left[\begin{array}{ccc}{T}_{21}& T_{31}& 0\\ T_{22}& T_{32}& 0\\ T_{23}& T_{33}& 0\\ 0& 0& 1\end{array}\right]\left\{\begin{array}{c}S\\ T\\ M_{\mathrm{\&phi;a}}\end{array}\right\}\mathrm{dr}---\left(35\right)$

Tenth step, rotor, body aerodynamic force virtual work are derived

(1) fuselage aerodynamic force virtual work

The aerodynamic force that fuselage is subject to, moment act on the center of gravity of airplane, and the virtual displacement of the corresponding center of gravity of airplane represents by body rigid body mode broad sense virtual displacement, and does not consider body Elastic mode, namely only considers fuselage aerodynamic force, virtual work that moment is done body rigid body mode.Therefore fuselage aerodynamic force virtual work is:

${\mathrm{\&delta;W}}_{\mathrm{AF}}=\underset{i=1}{\overset{6}{\mathrm{\&Sigma;}}}{F}_{\mathrm{fai}}\mathrm{\&delta;}{X}_{\mathrm{oc}}=\underset{i=1}{\overset{6}{\mathrm{\&Sigma;}}}{F}_{\mathrm{fai}}\mathrm{\&delta;}{X}_{\mathrm{fpi}}={\left[\mathrm{\&delta;}{X}_{\mathrm{fp}}\right]}^T{\left\{F_{\mathrm{fa}}\right\}}_{6\&times;1}---\left(36\right)$

In formula, F _fabe three power and three moments that act on body center of gravity place, the vibration shape matrix of corresponding 6 rigid body modes in body center of gravity place is unit matrix, so the virtual displacement in equation (36) is exactly the variation δ X of 6 rigid body mode coordinates of body _fp.

(2) blade aerodynamic power to be moved the virtual work done to blade

Blade aerodynamic power is projected to coordinate system in, premultiplication blade virtual displacement array, and along blade exhibition to integration, obtaining a slice blade aerodynamic power to this blade virtual work done of moving is:

$\mathrm{\&delta;}{W}_{\mathrm{br}}={\&Integral;}_{r_0}^R{\left\{\begin{array}{c}\mathrm{\&delta;u}\\ \mathrm{\&delta;v}\\ \mathrm{\&delta;w}\\ \mathrm{\&delta;\&phi;}\end{array}\right\}}^T\left\{\begin{array}{c}{F}_u\\ F_v\\ F_w\\ M_{\mathrm{\&phi;a}}\end{array}\right\}\mathrm{dr}={\&Integral;}_{r_0}^R{\left\{\begin{array}{c}\mathrm{\&delta;u}\\ \mathrm{\&delta;v}\\ \mathrm{\&delta;w}\\ \mathrm{\&delta;\&phi;}\end{array}\right\}}^T\left[\begin{array}{ccc}{T}_{21}& T_{31}& 0\\ T_{22}& T_{32}& 0\\ T_{23}& T_{33}& 0\\ 0& 0& 1\end{array}\right]\left\{\begin{array}{c}S\\ T\\ M_{\mathrm{\&phi;a}}\end{array}\right\}\mathrm{dr}---\left(37\right)$

(3) rotor aerodynamic force to move the virtual work done to propeller hub

Rotor aerodynamics is divided into and rotates virtual displacement works to the move virtual work done of propeller hub to 3 translations of propeller hub and 3.Rotor aerodynamics can be expressed as to the propeller hub virtual work done of moving the form being similar to (38).

$\begin{array}{c}{\mathrm{\&delta;W}}_{\mathrm{Xhe}}=\underset{k=1}{\overset{N_b}{\mathrm{\&Sigma;}}}{\&Integral;}_{r_0}^R{\left\{{\mathrm{\&delta;X}}_{\mathrm{he}}\right\}}^T\left[\mathrm{FABP}\right]\left[\mathrm{STM}\right]\left\{U2S0\right\}\mathrm{dr}+\underset{k=1}{\overset{N_b}{\mathrm{\&Sigma;}}}{\&Integral;}_{r_0}^R{\left\{{\mathrm{\&delta;X}}_{\mathrm{he}}\right\}}^T\left[\mathrm{FABP}\right]\left[\mathrm{STM}\right]\left[U2\mathrm{ST}\right]*\\ *<\left\{v_i\right\}+\left(\right[\mathrm{Cavho}]+C_{\mathrm{\&psi;k}}[\mathrm{Cavhc}]+S_{\mathrm{\&psi;k}}[\mathrm{Cavhs}\left]\right)\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{he}}\right\}+\left[\mathrm{Carvs}\right]\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{rs}}\right\}>\mathrm{dr}+\\ +\underset{k=1}{\overset{N_b}{\mathrm{\&Sigma;}}}{\&Integral;}_{r_0}^R{\left\{\mathrm{\&delta;}{X}_{\mathrm{he}}\right\}}^T\left[\mathrm{FABP}\right]\left[\mathrm{TSB}\right]\left\{\mathrm{UPDT}\right\}\begin{array}{c}\left[\begin{array}{ccc}{T}_{31}& T_{32}& T_{33}\end{array}\right]<\left\{{\stackrel{\&CenterDot;}{V}}_{\mathrm{sa}0}\right\}+\left[\mathrm{Dacrds}\right]\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{rs}}\right\}+\end{array}\\ +\left\{{\stackrel{\&CenterDot;}{v}}_i\right\}+\mathrm{\&Omega;}\left(C_{\mathrm{\&psi;k}}\right[\mathrm{Cavhs}]-S_{\mathrm{\&psi;k}}[\mathrm{Cavhc}\left]\right)\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{he}}\right\}>\mathrm{dr}+\\ +\underset{k=1}{\overset{N_b}{\mathrm{\&Sigma;}}}{\&Integral;}_{r_0}^R{\left\{\mathrm{\&delta;}{X}_{\mathrm{he}}\right\}}^T\left[\mathrm{FABP}\right]\left[\mathrm{TSB}\right]\left\{\mathrm{UPDT}\right\}\left[\mathrm{TSDV}\right]{\left[\begin{array}{ccc}\stackrel{\&CenterDot;}{\mathrm{\&phi;}}& \stackrel{\&CenterDot;}{v^{\&prime;}}& \stackrel{\&CenterDot;}{w^{\&prime;}}\end{array}\right]}^T\mathrm{dr}+\\ +\underset{k=1}{\overset{N_b}{\mathrm{\&Sigma;}}}{\&Integral;}_{r_0}^R{\left\{\mathrm{\&delta;}{X}_{\mathrm{he}}\right\}}^T\left[\mathrm{FABP}\right]\left[\mathrm{TSB}\right]\left\{\mathrm{UDTk}\right\}{\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_k\mathrm{dr}\end{array}---\left(38\right)$

First the aerodynamic force of every sheet blade is projected to coordinate system in, reprojection is to coordinate system premultiplication propeller hub virtual displacement array { δ X _hδ Y _hδ Z _h} ^t, along blade exhibition to integration, to the summation of Nb sheet blade, obtain rotor aerodynamics to propeller hub motion { δ X _hδ Y _hδ Z _h} ^tthe virtual work done; Calculate rotor aerodynamic force to propeller hub motion { δ φ _xhδ φ _yhδ φ _zh} ^tthe virtual work done also needs M _{Φ a}project to coordinate system in, and by { F _uf _vf _w} ^tsquare is asked to propeller hub center, then these moments are projected to coordinate system in, premultiplication propeller hub virtual displacement array { δ φ _xhδ φ _yhδ φ _zh} ^t, along blade exhibition to integration, Nb sheet blade is sued for peace.

11 step, the derivation of rotor blade dynamics modal equations.

(1) blade kinetics equation and space-time discretize

Formula (26), (31), (37) and (38) are substituted into variation equation (1) derivation blade kinetics equation.Adopt finite element method to blade motion { X _rscarry out discretize, blade is divided into several beam elements, with the column joints deformation displacement { q of beam element _ie} _k,idescribe the elastic vibration of blade, obtain the nonlinear dynamical equation of all degrees of freedom on a node basis of corresponding blade, can be expressed as:

$\left[M_{\mathrm{xrsxrs}}\left(q_e\right)\right]\left\{{\stackrel{\&CenterDot;\&CenterDot;}{q}}_e\right\}+\left[C_{\mathrm{xrsxrs}}\left(q_e\right)\right]\left\{\stackrel{\&CenterDot;}{q}\right\}+\left[K_{\mathrm{xrsxrs}}\left(q_e\right)\right]\left\{q_e\right\}+\left\{F_{\mathrm{xrs}0}\left(q_e\right)\right\}=0---\left(39\right)$

(2) blade trim

Trim equation is: [K _xrsxrs(q _e)] { q _e}+{ F _xrs0(q _e)=0 (40) solve blade equilibrium state by equation (40) under the trim displacement of each node, if it is { q _e0.

(3) linearization of blade kinetics equation

The displacement of each for blade node is expressed as: { q _e}={ q _e0}+{ Δ q _e(41) by (41) formula substitute into equation (40) carry out linearization process, remain into Δ { q _eonce item, equation (39) available linearization is:

$\left[M_{\mathrm{qe}}\left(q_{e0}\right)\right]\left\{\mathrm{\&Delta;}{\stackrel{\&CenterDot;\&CenterDot;}{q}}_e\right\}+\left[C_{\mathrm{qe}}\left(q_{e0}\right)\right]\left\{\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{q}}_e\right\}+\left[K_{\mathrm{qe}}\left(q_{e0}\right)\right]\left\{\mathrm{\&Delta;}{q}_e\right\}=0---\left(42\right)$

(4) blade model analysis and mode contract suddenly

Obtain the blade eigenfrequncies and vibration models of (42) undamped state, choose NP blade flapping, shimmy and twisted coupling lower mode and body by model analysis and kinematic train carries out modal synthesis, namely contracting is carried out to power system coupled wave equation rapid.If the vibration shape matrix corresponding to this NP lower mode is:

{Δq _e}＝[X _MB]{X _bp}(43)

Formula (43) is substituted into this linearizing oscillatory differential equation (39), and premultiplication modal matrix [X _mB] transposition, draw rotor blade dynamics modal equations.

12 step, the derivation of rotor and body coupled mode comprehensive dynamic equation.

The rotor kinetic energy of deriving, potential energy and variation and external force virtual work expression formula are substituted into equation variation equation (1), according to variation { the δ X describing rotor blade, body movement variable _rs} ^t{ δ X _fp} ^tequation is divided into two: corresponding variation { the δ X of Part I _fp} ^tbody and the equation of rotor blade sports coupling, the second corresponding variation { δ X _rs} ^trotor blade to move the equation be coupled with body movement

(1) body rotor blade and modal coupling kinetics equation

By body mode variation { δ X corresponding in equation variation equation (1) _fp} ^tthe part be multiplied takes out the kinetics equation being body mode and rotor blade coupled motions, and in formula, each matrix of coefficients contains cosn ψ _k, sinn ψ _kthe factor, i.e. periodic coefficient, and { X _rs} ^tby modal transformation formula (43) by blade node coordinate variable in body and rotor blade the coupled dynamical equation transform to blade Modal Space, obtain body mode and rotor blade modal coupling kinetics equation.

$\begin{array}{c}{\&Integral;}_0^t{\left\{{\mathrm{\&delta;X}}_{\mathrm{fp}}\right\}}^T\left(\right[M_{\mathrm{fpfp}}\left]\right\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{fp}}\}+[C_{\mathrm{fpfp}}\left]\right\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{fp}}\}+[K_{\mathrm{fpfp}}\left]\right\{X_{\mathrm{fp}}\}+\{F_{\mathrm{fp}}\}+\underset{i=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}\underset{k=1}{\overset{N_b^i}{\mathrm{\&Sigma;}}}{\&Integral;}_0^{R_i}[X_{\mathrm{hfp}}^i\left]\right\{F_{\mathrm{He}}\}\mathrm{dr}+)\\ +\underset{i=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}{\&Integral;}_0^{R_i}\left[X_{\mathrm{hfp}}^i\right]\left[M_{\mathrm{HB}}\left(X_{\mathrm{rs}}\right)\right]\left[X_{\mathrm{hfp}}^i\right]\left\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{fp}}\right\}+\underset{i=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}\underset{k=1}{\overset{N_b^i}{\mathrm{\&Sigma;}}}{\&Integral;}_0^{R_i}\left[M_{\mathrm{fpxrs}}(X_{\mathrm{rs}},\mathrm{\&psi;})\right]\left\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{rs}}\right\}\mathrm{dr}+\\ +\underset{i=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}\underset{k=1}{\overset{N_b^i}{\mathrm{\&Sigma;}}}{\&Integral;}_0^{R_i}\left[C_{\mathrm{fpxrs}}(X_{\mathrm{rs}},\mathrm{\&psi;})\right]\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{rs}}\right\}\mathrm{dt}+\underset{i=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}\underset{k=1}{\overset{N_b^i}{\mathrm{\&Sigma;}}}{\&Integral;}_0^{R_i}\left[K_{\mathrm{fpxrs}}(X_{\mathrm{rs}},\mathrm{\&psi;})\right]\left\{X_{\mathrm{rs}}\right\}\mathrm{dr})\mathrm{dt}=0\end{array}---\left(44\right)$

(2) rotor blade mode and body modal coupling kinetics equation

By blade variation { δ X corresponding in equation variation equation (1) _rs} ^tthe part that is multiplied is taken out, and obtains rotor blade and to move the equation be coupled with body mode motion:

$\begin{array}{c}{\&Integral;}_0^t{\left\{{\mathrm{\&delta;X}}_{\mathrm{rs}}\right\}}_{k,i}^T{\&Integral;}_0^{R_i}\left({\left[M_{\mathrm{xrsfp}}(X_{\mathrm{rs}},\mathrm{\&psi;})\right]}_{k,i}\right\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{fp}}\}+{\left[C_{\mathrm{xrsfp}}(X_{\mathrm{rs}},\mathrm{\&psi;})\right]}_{k,i}\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{fp}}\}+\\ +{\left[K_{\mathrm{xrsfp}}(X_{\mathrm{rs}},\mathrm{\&psi;})\right]}_{k,i}\left\{X_{\mathrm{fp}}\right\}+\left\{\right[M_{\mathrm{xrsxrs}}(X_{\mathrm{rs}},\mathrm{\&psi;})\left]\right\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{rs}}\}+[C_{\mathrm{xrsxrs}}(X_{\mathrm{rs}},\mathrm{\&psi;})\left]\right\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{rs}}\}+\\ +\left[K_{\mathrm{xrsxrs}}(X_{\mathrm{rs}},\mathrm{\&psi;})\right]\left\{X_{\mathrm{rs}}\right\}+{\left\{F_{\mathrm{xrs}0}(X_{\mathrm{rs}},\mathrm{\&psi;})\right\}\}}_{k,i})\mathrm{drdt}=0,(i=1,...,\mathrm{IR},k=1,...,N_b^i)\end{array}---\left(45\right)$

Same modal transformation formula (43) of passing through is by blade node coordinate variable in rotor blade and body modal coupling kinetics equation transform to blade Modal Space, obtain blade mode and body modal coupling kinetics equation.

(3) expression-form of rotor and body coupled mode comprehensive dynamic equation

Through arranging equation (44) and (45), final derivation body and rotor coupled mode comprehensive dynamic equation:

$\begin{array}{c}\left[M_{\mathrm{fpfp}}\right]\left\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{fp}}\right\}+\left[C_{\mathrm{fpfp}}\right]\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{fp}}\right\}+\left[K_{\mathrm{fpfp}}\right]\left\{X_{\mathrm{fp}}\right\}+\left\{F_{\mathrm{fp}}\right\}+\underset{i=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}\left\{F_{\mathrm{Hefp}}\left(X_{\mathrm{bpl}}\right)\right\}+\\ +\underset{i=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}{\left[{\mathrm{MHB}}_{\mathrm{fpfp}}\left(X_{\mathrm{bpl}}\right)\right]}^i\left\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{fp}}\right\}+\underset{i=1}{\overset{\mathrm{IR}}{\mathrm{\&Sigma;}}}\underset{k=1}{\overset{N_b^i}{\mathrm{\&Sigma;}}}\left(\right[M_{\mathrm{fpxbp}}(X_{\mathrm{bpl}},\mathrm{\&psi;})\}\left\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{bpl}}\right\}+\\ +\left[C_{\mathrm{fpxbp}}(X_{\mathrm{bpl}},\mathrm{\&psi;})\right]\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{bp}}\right\}+\left[K_{\mathrm{fpxbp}}(X_{\mathrm{bpl}},\mathrm{\&psi;})\right]\left\{X_{\mathrm{bpl}}\right\}{)}^{k,i}=\left\{0\right\}\end{array}---\left(46\right)$

The modal synthesis kinetics equation that the K sheet blade of I rotor is coupled with body:

$\begin{array}{c}{\left[M_{\mathrm{bpfp}}(X_{\mathrm{bpl}},\mathrm{\&psi;})\right]}^{k,i}\left\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{fp}}\right\}+{\left[C_{\mathrm{hpfp}}(X_{\mathrm{bpl}},\mathrm{\&psi;})\right]}^{k,i}\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{fp}}\right\}+\left(\right[M_{\mathrm{bp}}(X_{\mathrm{bpl}},\mathrm{\&psi;})\left]\right\{{\stackrel{\&CenterDot;\&CenterDot;}{X}}_{\mathrm{bpl}}\}+\\ +\left[C_{\mathrm{bp}}\right(X_{\mathrm{bpl}},\mathrm{\&psi;}\left)\right]\left\{{\stackrel{\&CenterDot;}{X}}_{\mathrm{bpl}}\right\}+\left[K_{\mathrm{bp}}(X_{\mathrm{bpl}},\mathrm{\&psi;})\right]\left\{X_{\mathrm{bpl}}\right\}+\left\{F_{\mathrm{bp}}(X_{\mathrm{bpl}},\mathrm{\&psi;})\right\}{)}^{k,i}=0\\ (i=1,...,\mathrm{IR};k=1,...,N_b^i)\end{array}---\left(47\right)$

(46) { F and in (47) _fpand { F _bp(X _bpl, ψ) }) ^k,icontaining the relevant aerodynamic loading item that winglet yaw motion and pitch input control cause.Be homogeneous equation when nothing control, solving obtained dynamics is stability characteristic under rotor and body are coupling in basic flight reference, and when having control state, then can study winglet and pitch input control or AFCS system to be coupled with body on rotor the impact of stability, also can be used for rotor and to be coupled with body the research of stability Active Control Design.

Equation (46) and (47) form multiple rotor and body coupled mode comprehensive dynamic system of equations, to be coupled stability and dynamic response for analyzing rotor with body.And solve this system of equations need adopt concrete computing method according to stability or dynamic response analysis requirement, be wherein coupled steady-state response calculate (or trim) is an iterative process.The solution obtained is the steady-state response of rotor and body modal coupling.On this basis, then with trim solution lienarized equation (46) and (47), the dynamic stability analysis carried out is then the dynamic response under stable state state of flight and stability analysis.For maneuvering flight state, body and blade can not trims, need direct solution equation (46) and (47) to obtain transient response, then according to instantaneous response analysis dynamic stability.

It should be noted that; rotor of the present invention and body Coupled Dynamics modal synthesis modeling method comprise any one and combination in any thereof in above-described embodiment; but embodiment recited above is only be described the preferred embodiment of the present invention; not the scope of the invention is limited; do not departing under the present invention designs spiritual prerequisite; the various distortion that the common engineering technical personnel in this area make technical scheme of the present invention and improvement, all should fall in protection domain that claims of the present invention determine.

Claims (7)

1. rotor and a body Coupled Dynamics modal synthesis modeling method, is characterized in that:

(1) build for describe body, rotor hub, blade motion coordinate system and body, rotor blade and aerodynamic model;

(2) build for simulating the body of rotor and blade Configuration Design structural dynamic characteristics and rotor blade and control structure system finite element method kinetic model;

(3) for assembling and the choice of many rotors, many blades and control system Coupling Dynamic Model, adopt modular idea about modeling, set up rotor and body modal coupling kinetic model;

(4) be loaded on for rotor and body derives body convected motion with the remaining nodal displacement of minor structure express in blade inertia and the load such as pneumatic, what derivation body foundation motion produced involve inertia and the load such as pneumatic;

(5) unsteady aerodynamic model that single order and higher order dynamics become a mandarin is built;

(6) the ONERA aerodynamic model of profile lift, profile drag and the moment coefficient response characteristic in unsteady flow environment is calculated based on aerofoil profile aerodynamic characteristics tests data, for calculating the aerodynamic characteristic of aerofoil profile in unsteady flow environment;

(7) build blade trailing edge winglet deflection Controlling model, to be coupled with body the model that stability controls and vibratory response or load control for lifting airscrew;

(8) build and realize lifting airscrew and be coupled with body the pitch control system model that stability controls and vibratory response or load control, for analyzing fly control and ACFS system to blade pitch control effect;

(9) set up complicated high-order nonlinear kinetics equation in matrix representation mode, set up the wing and body Coupled Dynamics modal synthesis model.

2. rotor according to claim 1 and body Coupled Dynamics modal synthesis modeling method, it is characterized in that: the coordinate system of the description body that described structure is complete, rotor hub, blade motion and body, rotor blade and aerodynamic model, adopt coordinate system based on body, the multi-body dynamics modeling method that rotor hub rotates, blade motion deformation is relative coordinate system, the coordinate transformation relation set up, establishes the coupled motions relation of multiple rotor and blade and body.

3. rotor according to claim 1 and body Coupled Dynamics modal synthesis modeling method, is characterized in that: described employing finite element method sets up body and rotor blade and control structure system dynamics model by analyzing rigid motion and elastic movement.

4. rotor according to claim 1 and body Coupled Dynamics modal synthesis modeling method, it is characterized in that: the method involving inertia and the load such as pneumatic analyzing that body foundation motion produces comprises, first adopt the remaining nodal displacement of minor structure to derive, be transformed into body Modal Space afterwards.

5. rotor according to claim 1 and body Coupled Dynamics modal synthesis modeling method, is characterized in that: the ONERA aerodynamic model calculating profile lift, profile drag and the moment coefficient response characteristic in unsteady flow environment based on aerofoil profile aerodynamic characteristics tests data.

6. rotor according to claim 1 and body Coupled Dynamics modal synthesis modeling method, it is characterized in that: set up and to be coupled with body the model that stability controls and vibratory response or load control with the blade trailing edge winglet yaw motion lifting airscrew that is control variable, by design analysis trailing edge winglet the exhibition of every sheet blade to position, length and width change, to provide initiatively or Passive Control rotor is coupled with body stability and vibratory response or load.

7. rotor according to claim 1 and body Coupled Dynamics modal synthesis modeling method, it is characterized in that: set up and to derive complicated high-order nonlinear kinetics equation in matrix representation mode, according to the operation relation of derived relation formula, computer programming is utilized to realize formation and the calculating of the whole matrix of coefficients of kinetics equation.