Quadcoptor Platform Equations of Motion: Dynamic Model15 min read

Quadcoptor Platform Equations of Motion: Dynamic Model15 min read

image_pdfimage_print

It’s been over one month since I last posted. Back then I took a journey into fluid mechanics to offer a glimpse into how propellers work. I was motivated to do this as I was previously working-out the electric motor and propeller dynamics when I started to recall my mechanics studies.

It’s a good time to revisit my method for going at this, “system” problem of quadcoptor design.

Quadcoptor Project Review:

Basic Design Objectives

I added an objective to try a design that doesn’t assume the quadcoptor will always be level.

  • Have you noticed that camera drones don’t bank? There are some simplifications to the design if you want to maintain a horizontal platform at all times.
  • I’m curious to model a system that can operate away from the horizontal plane, at least in theory. I’m designing for a larger, “flight envelope”.

Sub-System Breakdown

  • When we get to a physical build there will be electrical, mechanical, and software sub-systems but at this point I’m looking at the system on paper:
  • How does the problem break-down into clear input-output relationships? Think of these as layers:
    • One level operates with certain data and makes a decision to command another, lower layer.
    • That lower layer operates with the above command and other information (sensors) to command what is it responsible for.
    • …and so on producing…
    • Physical action of the quadcoptor that is sensed again by each layer in this process which repeats itself.

This repeating is what makes these layers, “loops”. They process inputs, command outputs, sense the environment as inputs…and go again in a, “loop”. The sampling-rates are the time between successive trips through each loop at each layer.

Layers or, “Loops” I am working within

  1. “Flight Planning”: I’m not even into this yet so ignore it for now:
    1. Whether we use a joystick or upload some coordinates our craft will follow based on GPS.
  2. “Platform Control”
    1. The Quadcopter flight dynamics: the complete system in the sky that can move up, down, sideways, and rotate.
    2. Perhaps bank and perform some complicated actions predictably. This will be a, “flight envelope” or range-of-motion that likely exceeds a commercial, “camera drone” and for which we can come to understand the mechanics and the limitations.
    3. Sensor inputs to the platform. This will involve sensors for position and orientation: we’re not there yet but it will involve measuring and calculating positions and angles we’re getting into now.
    4. Outputs from the, “platform controller”: Propulsion
  3. “Propulsion”: I committed to a, “quad-copter” so we’ll have four propulsion units
    1.  Kinematics: how the placement and geometry of the 4 propulsion units affect the platform with their output thrust.
    2. Kinetics: producing the thrust from each output unit.

Earlier Posts Covered Propulsion System Modeling

If you look back you’ll see we covered motor-propeller electro-mechanics including a foray into fluid mechanics of lift from propellers. We now know the thrust producer is going to be a function of propeller speed. Our, “platform controller” is going to calculate how it wants to move, determine what combinations of the four propeller speed changes it needs to make to attempt this movement. Then it will command each motor-propeller propulsion unit accordingly. This, “loop” will likely run at near 10Hz or-so (this will be explained in a future post).

These four propulsion units each have a speed control loop to maintain the commanded propeller speed. These are the innermost, “loops” of control. This loop will run faster than the platform control loop (some tens-of-Hz. Again, we’ll cover this later). The specifics will depend on the motor speed controller we buy. We will not likely make our own motor speed controller.

No motors, propellers, sensors, or controllers for each propulsion unit have been selected. We understand the electro-mechanics though. We’re not going to shop for any parts until we get all the way through the system design. We’re headed there now!

Quadcoptor Platform Control

In my reading and on paper I’ve been studying in this area from the beginning even though I’ve only been writing about the propulsion unit. System design is an iterative process. In the mind it is all mixed-up and everything is interrelated.

Quadcoptor Equations of Motion

Degrees-of-Freedom

Our quad operates with 6 degrees-of-freedom (DOF). Its center-of-mass (historical term is center-of-gravity CG so this will be our abbreviation) moves up, down, and sideways relative to an “inertial reference frame”. It also tilts-and-twists as described by three angles below. That term, “Inertial reference frame” just means fixed relative to the environment our object will be moving: not accelerating.

A typical choice of reference frame for a hobby-grade quadcoptor is the, “take-off” point as the zero-point reference, an X-axis pointing north, Y-axis pointing east, and Z-axis pointing down into the earth via the right-hand-rule relative to X and Y axis. Point your right arm north with palm facing east. Curl your fingers towards the east. That is the right-hand-rule from X-to-Y axis. Extend your thumb. It points down. That’s your local, “Z” axis. Let’s call this our North-East-Down (NED) reference frame.

Reference Frame Assumptions

Our craft will operate within sight of this point. We will assume a flat earth. As our quadcoptor moves in the sky the CG position will be determined as an (X,Y,Z) coordinate relative to the inertial reference frame we chose to pin to the ground: this NED frame. Technically, this point on the earth is accelerating as the earth rotates but we’re assuming the entire local flight space is rotating through the day as we all do here on earth, and we accept this as a fair, “inertial frame”.

The quadcoptor is not a simple point mass. If you hold a ball-bearing and move it around you can think of it as a point-mass, and even though you know a bearing twists and rolls ignore that and assume it doesn’t rotate at all. The quadcoptor ball-bearing equivalent places all the mass of the object at the CG for the purposes of the above (X, Y, Z) position relative the the inertial NED reference frame. However, the tilts and rotations of the quadcoptor define the, “orientation” of the entire assembly in the sky.

For aircraft there are three angular rotations commonly used to describe orientation:

\phi = roll

\theta = pitch

\psi = yaw

These gives us our remaining 3 degrees-of-freedom, for a total of 6 relative to the earth-fixed NED reference frame:

(X,\: Y,\: Z,\: \phi,\: \theta,\: \psi)

The, “Body Frame”

Figure 1-2 near the end of the Blakelock Chapter 1 extract below offers a helpful illustration and legend describing the body frame, fixed NED frame, and Euler angles. Here’s the figure. See the PDF below for a detailed description on pages 15 and 16.

The sketch below illustrates a coordinate system fixed to the body of the quadcoptor. The X-Axis is defined to point always from the CG to the location of the first propeller. The Y-axis points towards what I label as the fourth propeller axis, and again by right-hand-rule Z is orthogonal to body X-Y plane by the right-hand-rule.

I’m making a design assumption that our quadcoptor will be a square cross configuration and not an, “X-wing” form-factor. I don’t have a particular reason other than wanting to sketch it this way and keep the diagrams and coordinate reference frames as simple as possible.

ref_frames

Angular Equations of Motion

The Bouabdallah paper refers to the Lagrangian method before summarizing the equations of motion for the 3 rotational axes as described above. I had some difficulty following the author’s use of Euler Angles and I’m not sure the derivations use these angles consistently within the paper.

I’m going to try a Newton-Euler method (Morin p.376) because although the Lagrangian method is elegantly simple I find equating the rate-of-change of body momentum to applied forces and torques straightforward. In this case I get less confused between using body angles and Euler angles (body angles relative to the fixed-frame. I hope to repeat the derivations using the Lagrangian at some point.

Recommended Reading

I was referred to the following text as a, “classic”. The first edition was published in 1965, second edition in 1991. I acquired a copy. I see why it is revered. The first 10 pages of the first chapter lay-out the equations of motion for aircraft more concisely than any other text, paper, or online reference I have found. The clarity I need is Blakelock’s labeling and naming of the various body displacements, angles, rates-of-change of position and rotation, linear and rotational momentum, and his description of the Euler angles. In my opinion many other sources are not as clear on these points.

Retired USAF Colonel Blakelock died in 2015. Click his picture to read his obituary. I am thankful he left behind this great book for us.

It is confusing to think back-and-forth from the, “inertial” frame and the “body frame”. Specifically, my mind jumps to the rotations on the body that an eventual gyro package will measure as (p,q,r) where…

p = Rotational rate of the body rolling about the body X-axis.

q = Rotational rate of the body pitching about the body Y-axis.

r = Rotational rate of the body yawing (twisting) about the body Z-axis

Then there are, “Euler Angles” (\Phi , \Theta , \Psi )  representing roll, pitch, and yaw, respectively. These are relative to the fixed-frame: the local NED coordinate system in our case.

In some papers and texts capital and lowercase Greek letters, introduction of other symbols, presence or absence of the (p,q,r) definitions and other information can create confusion and doubt that one is indeed treating the accelerating coordinate frame(s) correctly at the same time accounting for vehicle dynamics.

Blakelock’s first 10 pages offer a clean and clear description that we will later see gets us to the equations in our main reference for the quadcoptor: the Bouabdallah paper.

Blakelock’s straightforward equations-of-motion derivation is included here:

blakelock_eom2

Timeout: Angular Momentum Clarification

In a linear system we know the simple statement, “momentum equals mass times velocity”. We know that linear motion in three dimensions means that velocity is a vector with x-, y-, and z-components. Mass for linear momentum is lumped to a point (our ball-bearing example above). It is a scalar. So, in a linear system we’re left with a momentum vector: the mass times the 3 components of it’s velocity.

A rotational system in three-dimensions will have a rotational momentum vector as well, but in this system, “mass” is represented as, “moment of inertia (I)” and velocity is angular velocity of the entire body represented as a three-coordinate vector \omega or \Omega  (capital or lowercase may vary but it is almost always this Greek letter, “omega”that is used).

At any instant you can think of this \omega-vector as the axle on which the body spins.

To complicate matters further, what was a scalar mass for the linear system becomes not simply a vector of inertia about various axes but a 3×3 matrix!

we’re going to ultimately arrive at a 3-element vector of angular momentum when we multiply our, ‘I’ matrix by our body angular velocity vector \vec\omega .

\vec L\:=\:\b I\vec\omega

First let’s get more comfortable with the inertia matrix, ‘I’.

Platform Inertia Matrix

The Blakelock text above skips over the details of the, “Inertia Tensor”: a matrix that represents to moment of inertia of the body. These are the, “I” terms. Good treatment of this topic is provided by David Morin in Ch. 9 of, “Classical Mechanics with Problems and solutions“. This Chapter is extracted here for your reference. The matrix describing a body’s moment of inertia is covered beginning on P. 376.

MorinCh9

Angular Momentum Changes

Page 393 in Morin above states, “Angular momentum, ‘L’ in the body frame changes relative to the fixed-frame (our NED frame) by two effects”..

  1. Body coordinates change relative to the body frame.
  2. The body frame itself rotates with respect to the fixed-frame.

Morin sets \vec L_0 = \vec L  as the momentum vector of the body at any given instant and imagines painting it on the body at this instant. The rate-of-change of momentum with respect to the fixed-frame is then…

\frac{d\vec L}{dt} = \frac{d(\vec L- \vec L_0)}{dt} + \frac{d\vec L_0}{dt}

That first term is the rate-of-change of momentum with-respect-to the body frame. This is what someone sitting on the rigid body would measure:

\frac{\delta \vec L}{\delta t} = \frac{d(\vec L- \vec L_0)}{dt}

The remaining \frac{dL_0}{dt}  term is the rate-of-change of the body-fixed momentum vector which is the body momentum at this instant as defined above. The momentum vector in this instant is assumed, “painted” on the body like our body axes. This temporarily, “constant” momentum is rotating with respect to the fixed NED frame on the ground at this instant. This is the \frac{dL_0}{dt}  term.

Morin reminds us that the fancy sum above that gives us total change in momentum \frac{d\vec L}{dt} is saying just this:

The total rate-of-change of angular momentum is the rate-of-change of momentum relative to the moving frame PLUS the rate-of-change of the momentum of the moving frame with respect to the fixed frame.

When we have any number of, “reference frames” for vectors this is how we map any vector from one frame to another. There is nothing magic or new here, even though the nomenclature looks fancy.

Body-Fixed Momentum Vector Relative to NED Frame

Just as Blakelock on page 13, Morin p. 93 gives us,

\frac{d\vec L_0}{dt} = \vec\omega \times \vec L

In Morin chapter 9 you can see derivation of the, “Inertia tensor” but it simplifies to a diagonal matrix when we select the body-frame axis to be principal axes of the body.

I =\begin{bmatrix} I_x & 0 & 0\\ 0 & I_y & 0\\ 0 & 0 & I_z \end{bmatrix}

Blakelock doesn’t zero the lower-left and upper-right terms of the inertial matrix. For a fixed-wing aircraft he maintains a cross-coupling of the inertia terms associated with roll and yaw. If you study his figure 1-2 in the book section above you can make sense of this for an aircraft that does not appear symmetrical front-to-back. However, our quadcoptor is symmetrical here too. Our inertia matrix simplifies to the diagonal matrix above.

This leaves us with L_i = I_iw_i for each body axis x, y, and z. The cross-product  \vec\omega \times \vec L  can be computed via the determinant of the matrix below. The determinant of a 3-by-3 matrix can be had by augmenting it with the first two columns as shown and subtracting left-diagonal products from right-diagonal products.

\begin{vmatrix} \hat x & \hat y & \hat z\\ w_x & w_y & w_z\\ I_xw_x & I_yw_y & I_zw_z \end{vmatrix} \begin{matrix} \hat x &\hat y \\ w_x &w_y \\ I_xw_x & I_yw_y \end{matrix}

Below are the hand calculations for the  \vec\omega \times \vec L  cross-product. The black lines represent the, “right diagonal” products. the red lines represent the left-diagonal products.

These hand calculations agree with Blakelock equations 1-32 (the central term with the principal axes moments of inertia in parenthesis). Our rotational rates match Blakelock’s as follows:

\omega_x = P \\ \omega_y = Q \\ \omega_z = R

These are the body rates-of-rotation about the principal axes of the body: the frame, “painted” on the body. We will soon get into details for an Inertial measurement unit (IMU) with gyros that will give us (P,Q,R) data for our platform controller.

Rate-of-change of body momentum relative to the body frame: \frac{\delta \vec L}{\delta t}

This is what we would measure if we were sitting still on the body.

As in Blakelock’s equation 1-32 we will assume the primary-axis inertia terms are constant. For aircraft, missiles, and the like fuel burn changes mass and distribution of this mass within the body. Our battery-powered quad will not experience these changes, but it worthwhile to consider this effect should we desire to deliver a payload with our quadrotor. When we drop the load the mass and moments would change.

Blakelock indicates that even when this is relevant for fuel-thirsty aircraft the above assumption works over time intervals much longer than our tens-of-Hz platform control loop. We would implement a, “sliding” controller or, “gain scheduling” system that effectively remodels our design periodically as fuel is consumed or payload is released. between those slower updates we would take our moment-of-intertia, “I” terms as constants.

Thus, we have

\frac{\delta \vec L}{\delta t} = (\dot PI_x, \: \dot QI_y,\:\dot RI_z)

Combining this with the result above our equations of rotational motion become…

\sum_{}^{} T_x = I_x\dot\omega_x\:+\:(I_z-Iy)\omega_y\omega_z \\\sum_{}^{} T_y = I_x\dot\omega_y\:+\:(I_x-Iz)\omega_x\omega_z \\\sum_{}^{} T_z = I_x\dot\omega_z\:+\:(I_y-Ix)\omega_x\omega_y

If you assume no applied Torques T (natural rotational response of the system) you can rearrange the above equations to match the Bouabdallah paper equation (6). Note that the, “I” subtractions are reversed because assuming zero applied torque we move terms to the other side of the equal sign and solve for the angular acceleration.

These match when we relate Bouabdallah’s terms to our terms above as follows…

\dot \phi = \omega_x \\ \dot \theta = \omega_y \\ \dot \psi = \omega_z \\ \ddot \phi = \dot \omega_x \\ \ddot \theta = \dot \omega_y \\ \ddot \psi = \dot \omega_z \\

Equations of Linear Motion

The Bouabdallah paper describes a test-stand for which various various control strategies are studied. The test stand pins the CG to a point, constraining it to angular motion only. The paper does not cover details for linear motion.

Blakelock chapter 1 above offers good descriptions for modelling linear motion. Below is his conclusion. If you look back at his Equation 1.15 you can see how each primary axis velocity couples into acceleration of the other axes via another cross-product. I’m going to leave this here for now because depending on how we measure position and velocity relative to the body and/or fixed NED frame these terms may or may not be relevant in our design.

 

Concluding Remarks…

It looks like we have a good handle on the equations of motion now. I am more confident in nomenclature and with the mental gymnastics involved with thinking within the two reference frames:

  • The North-East-Down (NED) frame pinned to the earth at our take-off point.
  • The Body frame we can think of as permanently painted along primary (x, y, z) axes of our quadcopter.

I’m going to pause here, and pick-up the applied torque and force elements of the equations-of-motion in the next post. These are the following terms in the equations above.

\sum_{}^{} T_{(x,y,z)}\\ \sum_{}^{} F_{(x,y,z)}\:

They are the torques and forces that will act on the quadcoptor. They will include the desired outputs from the propulsion units, and we can also model external factors such as wind. We will need a controller design that is robust to some applied, “disturbances” from wind.

I’m looking forward to getting into sensor and control system designs!