Understanding the control characteristics of electric vertical take-off and landing (eVTOL) aircraft for urban air mobility

The aim of the present paper is to present an understanding of the control characteristics for representative eVTOL (Electric Vertical Takeoff and Landing) quadcopters having the same properties of interest as the eVTOL conﬁgurations developed presently in Urban Air Mobility. In this context, the paper concentrates on eVTOL quadcopters with pure rpm control, pure propeller pitch control and both rpm and propeller pitch control. Dynamic effects of electric motors can potentially have signiﬁcant effects on the ﬂight characteristics of these vehicles. The paper will demonstrate for a pitch manoeuvre that while varying purely the rpm of the rotors results in a system that behaves as an acceleration-control system, varying purely the propeller pitch corresponds more to a velocity-control system. The rotor pitch control is more sensitive to the coupling rotor-motor and its transients dynamics. Varying both the rpm and the rotor pitch for controlling the pitch manoeuvre is the best option as this brings the system more towards velocity-control behaviour and dampens the transients due to the motor dynamics. © 2022 The Author. Published by Elsevier Masson SAS. This is an open access article


Introduction
History has shown that each century that human society has passed is highly dependent on the technology that people have developed.A strategic scan of the aerospace environment at the beginning of this 21st century strongly suggests that the world might be approaching a new age of airpower-the era of electrified/hybrid aircraft propulsion.Undeniably, starting from the Wright Brothers piston engine flight in 1903, the jet engine of the 1960s, to the space age of today, one can say that leaps in propulsion technology have marked the different ages of human flight.The technological advancements, brought at the beginning of 21st century by the revolution in data exchange, computational power, sensors, wireless communication, internet, and autonomy, contributed to the development of a new vision towards a new era of aviation.
The World eVTOL (Electric Vertical Takeoff and Landing) Aircraft Directory (https://evtol .news/aircraft) started in 2016 when there were only a half-dozen known eVTOL designs.Today it counts more than four hundred fifty companies which are in the process of developing prototypes in a fierce competition between startups, including Kitty Hawk (US), Lilium (Germany), Joby Aviation (US), E-Hang (China), Volocopter (Germany), as well as large firms like Airbus (with a special A3 by Airbus located in Silicon Val-E-mail address: m.d.pavel@tudelft.nl.ley), Boeing (US), Bell (US), Embraer (Brazil), and Uber (US).Additionally, the "Big Four" technology companies from US -Amazon, Google, Apple, and Facebook -known for the disruption of wellestablished industries through technological innovation are moving into the eVTOL aircraft endeavours, placing vast venture capital and highly-talented human capital into these efforts.These new eVTOLs belong to a new mobility -Urban Air Mobility (UAM)which will bring our mobility in the cities from 2 to 3 dimensions.A preliminary investigation of individuals' perceptions and expectations towards the adoption of flying cars as a 3 rd dimension revealed that various individual-specific socio-demographic, behavioural and driving attributes, as well as individuals' attitudinal perspectives towards the cost, safety, security and environmental implications of the flying cars, affect their willingness to adopt this new emerging transportation technology [20].
Thinking about the advancements in aircraft technology for UAM, eVTOL aircraft have been prototyped and are to be flown in the coming years with shapes that had previously been relegated to science fiction.These configurations require new flight control schemes such as adaptive control [22] or fault-tolerant control [23] as it corresponds to aircraft with large operating envelopes.This freedom to design new configurations emerges directly from the availability of reliable electric motors.These motors are controlled by fast, closed-loop electronic feedback controllers and sensors.Those electronic control systems can quickly command speed changes to the motors, which are very responsive.In-deed, Electric motors offer thrust capabilities for alternative modes of aircraft manoeuvrability through distributed electric propulsion.However, the envisioned UAM concept brings several technological challenges w.r.t.flight control systems, namely: unconventional command and control modules applied to new systems with various degrees of automation and autonomy in unconventional displays and inceptors.This is why, these new controllers need to be designed and programmed properly in order to make eVTOL aircraft fly with stability and precision.
Previous assessment of multirotor control systems for eVTOL aircraft, especially variable speed-controlled models has been performed primarily on small scale UAVs in order to be primarily used for surveillance, package delivery, etc.In such configurations, the electric motor is decoupled from the rotor analyses for performance, loads, or vibratory investigations [1,6,10,15,16,21].However, when upscaling such configurations to multirotor configurations sized for passengers in UAM mobility, the transient dynamics effects between the rotor and the motor cannot be neglected.Recent literature shows that dynamic effects of electric motors can potentially have significant effects on the rotor and power system design [3,5,7,11,12,17,18].It should be noticed that controlling a multirotor VTOL aircraft can be achieved considering either 1) fixed-pitch propellers where the control of the aircraft is achieved through rpm variation or 2) variable-pitch propellers where one can change the blades pitch.The driving force of an electric motor is torque-not horsepower.The power produced by the electric motor depends on the speed and torque requirements at each operational point.Motor's performance curves defined as profiles giving the change in torque versus change in speed characteristics affect the energy storage of the battery, the power distribution and finally the eVTOL flight performance.The objective of this paper is to investigate the rotor-electric motor coupling for a VTOL in a pitch manoeuvre understanding the control characteristics of eVTOL aircraft in the pitch axis.The paper will investigate quadcopters with either variable speed control or pitch control demonstrating the need to understand such transient dynamics in order to build reliable controllers.For the beginning the paper will analyse the case of transient dynamics in the pitch axis for a helicopter case.Then, analysis of this transient dynamics and rotor-motor couplings for the quadcopter case will be performed.Finally the paper will conclude on the rotor-motor characteristics for the UAM mobility as compared to the helicopter case.

A simple manoeuvre, analysed in the time domain
In order to get some physical feeling for the problem of rotormotor coupling for eVTOL aircraft, a very simple manoeuvre is used as an example, i.e. the first few instants during the transition from hover to forward flight.This manoeuvre will be first investigated for the helicopter case.Studying the problem first for a helicopter and then for a quadcopter will help understanding how VTOL response types depend on rotor characteristics.This is important as eVTOL handling qualities are under development and a lot of knowledge can be gained from understanding how the overall mission performance was approached at helicopters.Fig. 1 presents a simple model representing the pitching motion for a helicopter.One can see that for a helicopter the pitching motion is controlled by controlling the main rotor, tilting it forward in order to build the forward speed.As the pilot applies longitudinal cyclic pitch θ 1s tilting the helicopter swashplate in order to build the vehicle pitch rate q, the rotor will flap back with a disc tilt angle a 1 .
One may assume that simply a pitching motion of the helicopter occurs at the very beginning of this manoeuvre, before forward speed builds up and begins to have an influence.In classical treatments of the subject, the rotor disctilt is often assumed to respond instantaneously to control inputs, as well as to pitching motion and helicopter velocity.This in fact is equivalent to neglecting the transient flapping motion, which indeed damps out very quickly after a disturbance [4].In this classical approach, only the quasi steady response of the rotor disc is taken into account.In the case considered, backward tilt of the rotor disc with respect to the control plane (no-feathering plane) is given by [14]: indicating the classical damping effect that the flap motion has on the building of the fuselage pitch rate q.The notations in eq.(1) correspond to γ the Lock number, the rotor rpm and q the pitch rate (rad/sec, deg/sec).This equation may be combined with the equation describing the pitching of the body: where M y is the fuselage pitch moment (Nm) and I y is the helicopter moment of inertia on y axis (Nm 2 ), q is the pitch rate (rad/sec, deg/sec) and θ f is the fuselage pitch tilt (rad, deg), positive up.
From Fig. 1 it follows that the fuselage pitch moment with respect to the disc plane (DP) 1 can be expressed as: where h is the distance between body CG and rotor hub (m), T is the rotor thrust force (N), θ 1s is the pilot longitudinal cyclic input (rad, deg) and a 1 is the longitudinal flapping angle (more correct the longitudinal disc-tilt angle).Substituting (3) into (2), using the small-angle approximation and with the notation K h = T h/I y gives the final pitch equation as: 1 Looking at the rotor system, there are three distinctive rotor planes which can be used to characterize the helicopter dynamics: 1) the shaft Plane (SP) defined as the plane perpendicular to the rotor shaft (hub), 2) the control plane (CP) defined as the plane containing pilot pitch control input and 3) the disc plane (DP) defined as the plane containing the final stabilized rotor position.One may conclude that the helicopter pitch motion in this simple case for helicopter will depend on the pilot input and the rotor flapping.Substituting now (1) into (4) yields: which is a first order differential equation in time that needs to be solved in order to determine the helicopter motion.Considering the initial solution as q(0) = 0 results in the solution of the pitch motion for the helicopter manoeuvre as: The control type characteristics of the helicopter can be changed through rotor design.Namely, when the rotor is mounted flexible on the shaft, the equivalent "rigid-blade concept" can be used for rotor modelling.In the rigid blade concept, the flexible blades are assumed to be rigid in bending and torsion, the blade flexibility being concentrated in virtual hinges by means of springs (flexible lumps).In the case of a flapping hinge this means that the blade flapping motion can be characterized by a specific spring constant K β .A fuselage moment is therefore introduced by the spring that needs to be accounted for, see Fig. 2. The restrained moment acting on the blade flapping hinge is: The total elastic moment acting on the body pitch axis is obtained by projecting each blade elastic moment on the nonrotating fuselage pitch axis, multiplying it with the number of the blades N, and averaging over the azimuth.Observing that β s = (θ 1s − a 1 ) is the longitudinal disc tilt yields the elastic moment on the pitch axis as: where N is the number of blades and K β is the rotor spring constant (Nm).The fuselage pitch moment w.r.t.DP, eq. ( 3) becomes now: Using the small-angle approximation, the pitch equation of motion will still have the form of eq. ( 4) but now the term K h will have a different expression: showing that the fuselage pitch motion depends on the pilot input, rotor flapping and on the hinge spring properties.Consider as numerical example the helicopter data as given in the Appendix and assume a 1 deg pilot longitudinal cyclic pitch input (θ 1s = 1 deg).Fig. 3, left-hand-side, represents fuselage pitch rate responses for values of K β typical for a teeterrotor (K β = 0 Nm) and a semi-rigid rotor configuration respectively (K β = 46000 Nm).
In handling qualities assessments, the aircraft response to an input is characterized by its response-type.There are two different response-types which affect the pilots point of view on the aircraft: acceleration-control response and velocity-control response.In an acceleration-control system the pilot command affects system's acceleration and therefore it takes more time to correct the aircraft position.An acceleration-response type would be described by the pilot as "sluggish to control input" or "tends to wallow" (low bandwidth).A velocity-response type corresponds to the pilot command being affected by the system's velocity.Typical pilot commentary for a velocity-type response would be more as "crisp" or "rapid and well damped (high bandwidth).Looking at Fig. 3 one can read important motion characteristics through which the system response characteristics are defined that will define the future performance requirements of the controller.
First, the rise time can be defined as the amount of time the system takes to go from 10% to 90% of the steady-state, or final, value.From Fig. 3 it follows that while for the teeter rotor the fuselage time response is, t r = 4 sec for the semi-rigid rotor is around t r = 1 sec.This means that the teeter rotor response is slower showing a response characteristic typical for acceleration-control system.An acceleration-control system requires pilot anticipation in order to predict fuselage behaviour.In an acceleration-control system, the command to the wished position requires two integration levels for pilot anticipation and is characterized by higher pilot workload.The semi-rigid rotor is faster in response showing a response characteristic typical for velocity-control system.In a velocity-control system (rate-control) the command to the wished position requires one integration level from the pilot and therefore less pilot anticipation is required.Second, the time constant of a transient motion can be defined as the time needed for the system's step response to reach 1 − 1/e ≈ 63.2% of its final (asymptotic) value (or equivalently the time needed for the system to reduce its transient to 1/e ≈ 36.2% of its initial value).The time constant for the helicopter case is according to eq. ( 5) This demonstrates that, for the example considered in this paper, the teetering rotor has a time constant of approximately 2.1 sec while the semi-rigid rotor time constant is 0.4 sec.Again, this shows that the teetering rotor is much slower than the semi-rigid rotor (or equivalently, the semi-rigid rotor is much faster than the teetering rotor) as the time needed for the system to reduce its transient motion is much longer in the case of the teetering rotor when compared to the semi-rigid rotor.
Both, the rise time and the time constant of the system can be connected to the helicopter bandwidth.More precisely, the rise time/time constant of the system are inverse proportional to the system bandwidth.In the example of this paper this means that the teetering rotor has a smaller bandwidth w.r.t.semi-rigid rotor.In practice, the bandwidth is defined as that frequency beyond which system closed-loop stability is threatened, or in other words, the bandwidth frequency is defined as the highest frequency at which the pilot can double his gain or allow a 135 degree phase lag between control input and aircraft attitude response without causing instability [13].The higher the bandwidth, the larger will be the aircraft's safety margin in high gain tracking tasks.For the case of small amplitude, high-frequency pitch motions, bandwidth criteria have been developed for both fixed-and rotary-wing aircraft as part of aircraft flying and handling qualities.Therefore, rise time and time constant of the system connect aircraft's response to the flying qualities characteristics of the aircraft.
The art of a control system to be developed is to find the optimum trade-off between cost, rise-time, bandwidth, overshoot and stability.Therefore, the goal of the control system to be designed will be to adjust the process variable to a desired state: either increasing the rise time (slow down the closed loop control system as a whole) if the designer wants the system inputs to flow into it gradually and not too fast so that the controller saturates; either decreasing the rise time (speed up the closed loop control system) if the designer wants to be able to change the system variables as fast as possible.
Next, consider that some refinement of eq. ( 1) is introduced in the model accounting for the tilting dynamics of the rotor disc: where τ β in the time constant of the flapping motion and ȧ1 = da 1 dt is the first-order disctilt dynamics.According to the theory, this kind of extension of the equation for the disctilt corresponds to taking into account the low-frequency flapping mode dynamics (the so-called regressing flapping mode), on top of the steady solution.Fig. 3, right-hand side, represents the helicopter pitch response when considering first order disctilt dynamics (eq.( 11) added to the body motion eq. ( 10)) and varying the flapping time constant τ β between 0 and 0.3 seconds.Looking at this figure one can see that the first-order disctilt dynamics ȧ1 does influence the response of the semi-rigid system rather profoundly, creating an overshoot in response, in such a way that it will probably be noticeable to the pilot.Percent Overshoot as system characteristic is the amount that the process variable overshoots the final value, expressed as a percentage of the final value.On the other hand, in the case of the teetering rotor, the additional dynamics due to ȧ1 creates no overshoot in the response and is hardly noticeable for the pilot.Overshoot and number of oscillations till the settling time are parameters which should be carefully investigated in their relation to model uncertainty.They may lead to poorer control performance and potentially instability in response [8].The above-discussion on control characteristics of a helicopter in the pitch axis will be next used for understanding the control characteristics of the aircraft developed in the UAM future.

Pitch control characteristics of eVTOL quadcopter in urban air mobility
Consider next the case of an electric eVTOL quadcopter.Fig. 4 presents the pitching motion for the quadcopter case.One can see that in the case of a quadcopter the swashplate is eliminated and the pilot can use either 1) rpm variation as in small scale UAV case by increasing the rpm on the back rotor and decrease the rpm in the front rotor or 2) vary the blade pitch of the propellers increasing the blade pitch of the back propeller and decreasing the blade pitch on the front rotor.However, the difference between controlling a small scale UAV and a larger scale quadcopter is that quadcopter's rotors are mounted flexible on the shaft and therefore involve also flapping motion.The fuselage pitch moment is given by varying the rpm of rotor 3 and 1, namely, the rpm of rotor 3 increases and the rpm of rotor 1 decreases: Using the notation of thrust coefficient C T = T ρ( R) 2 π R 2 results in the following expression for the pitch moment: with l the horizontal distance between CG and rotor hub and Consider the pitch moment with respect to the Oxyz system of the quadcopter.Consider a pitch θ of the fuselage and the flapping angles on the two rotors a 13 on the rotor 3 and a 11 on the rotor 1.The pitch moment w.r.t. the Oxyz system can be determined as (after making small angle approximation for the flapping angle): The contribution of the pitch rate to the flapping motion will be for the rotor 3 and a 11 ∼ = − 16 γ q 1 therefore eq. ( 14) becomes: The pitch equation of motion of the quadcopter becomes then: ) sin θ cos θ (16) Considering now an equilibrium point where 3 = 1 = 0 and a variation δ of the angular velocity in order to pitch the quadcopter, positive on rotor 3 and negative on rotor 1, i.e. 3 = 0 + δ and 1 = 0 − δ , the pitch motion becomes: Consider next that the quadcopter rotors have a hinge spring which introduces a hinge spring moment on the fuselage as represented in eq. ( 8).The longitudinal disc tilt is in the case of the quadcopter β s = −a 1 so that the total elastic pitch moment is a combination of rotor elastic moments of rotor 3 and 1, i.e.
where and adding eq. ( 18) to the pitching moment ( 14) results in the pitching motion as: A general rule is that the motors should be able to provide twice as much torque at maximum, therefore δ ∈ [0, 2 0 ].Assume that C T 3 = C T 1 = C T for the quadcopter with constant pitch.Equation (19) for the pitch response of a quadcopter with constant pitch and variable rpm becomes: Assume that the quadcopter is a variable pitch configuration.Equation (19) for the pitch response of a quadcopter with variable pitch is then:

Coupled rotor -motor dynamics for a pitch manoeuvre
The new freedom of design in UAM emerges directly from the availability of reliable brushless direct current (DC) electric motors − known as BLDC 2 motors or BL motors − controlled by fast, closed-loop electronic feedback controllers and sensors.
Those electronic control systems can quickly command speed changes to the motors, which are very responsive.The BLDC motors, which are used to drive rotor or propeller in eVTOL, pursue abilities of rapid speed response and disturbance rejection.The comprehensive model integrated by bare BLDC dynamics and electronic speed controller (ESC) model is depicted in Fig. 5 from reference [2].The Brushless DC motor equivalent model consists of two parts: 1) one is an electrical part which calculates electromagnetic torque and current of the motor 2) the other one is a mechanical part composed of the BLDC motor (a fan blade attached to a permanent magnet rotor that surrounds the electromagnetic coils of the stator and associated control electronics) and an inverter for speed control of BLDC motor.Assuming the variations of the stator self-inductance with rotor position and the mutual inductance between the stator windings are negligible, the electrical dynamics of the BLDC motor can be described as [2,5]: where i a is the motor armature current, L a is the equivalent circuit armature inductance, R a is the equivalent resistance, e b is the motor-back electromotive force (or just the motor-back EMF 3 ) K e is the back-EMF constant in units of V/(rad/s), r is the drive system gear ratio and V a is the voltage applied at the armature.The coupling of permanent magnet synchronous motors with flexible rotors has been of concern for the eVTOL quadcopter [9,12].The coupled rotor-motor mechanical equation of motion can be written according to [12] as: where is the rotor rpm, I r is the rotor rotational moment of inertia, J is the inertia of the high-speed drive components (motor and coupled transmission components), r drive system gear ratio, K m is the motor torque constant, i a is the motor armature current, Q A is the rotor aerodynamic torque with drive system and B is 2 Brushless DC motors are similar to alternating current (AC) motors powered by a 3-phase waveform. 3The term back electromotive force, or just the back-EMF, is most commonly used to refer to the voltage that occurs in electric motors where there is relative motion between the armature of the motor and the magnetic field from the motor's field magnets, or windings.
a linear representation of mechanical friction or viscous losses in the drive system.The motor torque constant K m can be related to K e back-EMF constant through the relationship K m = cK e , where proportionality constant c is the torque SI unit conversion constant (assumed 0.7374 lb-ft/Nm as in ref. [5]).
For current state-of-the-art BLDC motors, expected inductance values are on the order of a few micro Henry and result in extremely low time constants of the current dynamics [5].Therefore, one can make the assumption that the electric current response is instantaneous.This is equivalent with assuming that, for direct drive response, the inductance contribution to the tension of the circuit is small and therefore negligible, i.e.L a ≈ 0. This yields the expression of R a i a = K e r − V a in eq. ( 22) which then substituted in (23) provides the fundamental dynamic response of the coupled rotor-motor direct drive system as: Equation ( 24) can be rewritten as: wherein the time constant of the motor is: Looking at the value of this time constant one can see that this scales proportionally with the rotor and motor inertia.Increasing the values of the ratio K 2 e /R a will reduce this time constant and can therefore help reduce the rotor speed change in time.The mechanical advantage of the gearbox between the motor and rotor needs to be selected to maximize the motor efficiency, given the nominal rotor design operating conditions, while minimizing the component weight.One can see that increased gearbox ratios can further help reducing the time constant.The aerodynamic torque can be expressed as a function as its torque coefficient 2 .In hover the torque coefficient is given by where σ is the blade solidity, C Dm is blade medium drag coefficient and λ i is the induced velocity which in hover can be expressed in actuator disc theory as λ i = W 2ρπ R 2 .Equation ( 25), (26) combined with eq. ( 19), (20) or (21) can be used to represent the response in a pitch manoeuvre of a quadcopter with both speed and pitch control, quadcopter with rpm control and respectively a quadcopter with pitch control.Note that rotor longitudinal disctilt a 1 and lateral disctilt b 1 will introduce one extra component of the rotor torque in the pitch and roll di- 2 a 1 in eq. ( 19), ( 20) or ( 21) which was neglected in this paper as the lateral dynamics was not analysed in the present investigation.Furthermore, higherorder yaw unbalance may arise, considering the nonlinearity of the torque-angular velocity and torque-blade pitch relationships, this was not included in the paper.
Consider as numerical example the quadcopters data as given in the Appendix (data correspond to the 4-PAX quadcopters used in ref. [12]).Fig. 6 presents some calculations for eVTOL response to a step input in the pitch when only rpm control was used, when only blade pitch control was used and when both rpm and blade pitch control was used.Linearizing equation ( 21 19) and combining this with eq.(25) one can obtain the time constant of the transient pitch motion when both rpm and pitch control are applied, τ q−rpm = γ 0 16 I • 2 ).Looking at Fig. 6, for a control step input given in pitch, rpm or both pitch and rpm, one can see that the rpm control corresponds more to an acceleration-control vehicle and the variable rotor pitch control corresponds to a velocity-control system.The time constant of the motor plays an important role in the case of rpm control and pitch+rpm control and yields an overshoot in the pitch response.Varying rotor pitch results in faster changes in pitch response as compared to the pure rpm control obtained through varying motor voltage.This can be seen also in the dead time (i.e. the elapsed time before a response is discernible in the axis) of the initial response of the pure rpm control system.In designing V/STOL vehicles it is known that controlling a system with high time constants requires pilot full attention, the pilot tending to overcontrol at time constants above 1.2 seconds [19].Also, it is known that an allowance for a 0.3 extra second delay with regard to vertical thrust response can result in a degradation of approximately one unit in handling qualities pilot rating [19].It follows that system's time constants are the simplest and the most powerful basis for revealing general patterns in system behaviour.Varying both motor speed and rotor pitch reduces the effects of motor dynamics and eliminates the overshoot in response but it increases the time constant so it slows down the vehicle response.

Conclusions
The goal of the present paper was to present an understanding of the control characteristics for the eVTOL quadcopter concentrating on the first seconds of a pitch manoeuvre executed from hover.Comparing the case of an eVTOL quadcopter to the case of a helicopter, the paper demonstrated that while varying purely the rpm of the rotors results in a system that behaves as an accelerationcontrol system, varying purely the propeller pitch corresponds more to a velocity-control system.An acceleration-control system requires pilot anticipation in order to predict fuselage behaviour.In an acceleration-control system, the command to the wished position requires two integration levels for pilot anticipation and is characterized by higher pilot workload.In a velocity-control system (rate-control) the command to the wished position requires one integration level from the pilot and therefore less pilot anticipation is required.However, the rotor pitch control is more sensitive to the coupling rotor-motor and its transients dynamics.Varying both the motor speed and the rotor pitch for controlling the pitch manoeuvre is the best option as this brings the system more towards velocity-control behaviour and dampens the transients due to the motor dynamics, however it slows down the vehicle response.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Fig. 1 .
Fig. 1.Modelling the pitch motion of a helicopter after a longitudinal cyclic step input.
0 )I bl was written as a function of flap frequency assuming an rpm constant for the rotor 0 .Assuming a 13 = a 11 ∼ = − 16 γ q 0

Fig. 6 .
Fig. 6.Quadcopter eVTOL response in a pitch manoeuvre for the case of a pure rpm control, pure blade pitch control and mixed rpm and pitch control.

2 N 2 0
) in pitch gives the time constant of the transient pitch motion as τ q = γ 0 16 I y • 2 N•I bl .Linearizing equation (20) in pitch and combining this with eq.(25) for pure rpm control results in the time constant of the transient pitch motion as τ rpm = γ 0 16 I y • • l • K ).Finally, linearizing eq. (