2021 International Conference on Computer Communication and Informatics (ICCCI -2021), Jan.

27 – 29, 2021, Coimbatore, INDIA

Modelling, Simulation, and Implementation of

PIDController on Quadrotors

Abhishek Iyer Dr. Hari Om Bansal

Department of Electrical and Electronics Engineering Department of Electrical and Electronics Engineering
BITS Pilani, Pilani Campus BITS Pilani, Pilani Campus
Jhunjhunu, India Jhunjhunu, India
[email protected] [email protected]
2021 International Conference on Computer Communication and Informatics (ICCCI) | 978-1-7281-5875-4/21/$31.00 ©2021 IEEE | DOI: 10.1109/ICCCI50826.2021.9402301

Abstract— A quadrotor is a type of Unmanned Aerial Vehicle is a flow chart diagram explaining the approach taken to
(UAV) that is popular for its vast applications due to its small size implement the proportional integral and derivative (PID)
and maneuverability. UAVs fall under the category of drones. controller to successfully to achieve attitude and altitude control
Drones are unmanned aircrafts or ships guided by remote control on the quadrotor.
or on-board computers. They are high in demand because of their
deployment for surveillance, product delivery, aerial photography
and defence applications. Implementing a control system on a
quadrotor can be very challenging due to its nonlinear functioning.
Thus the motivation of this paper is to successfully implement the
PID algorithm on the plant to achieve altitude and attitude control.
This paper also compares a linearized plant system of the
quadrotor with that of a nonlinear plant system.

A mathematical model for the quadrotor, progressing from a

uni dimensional to a three dimensional configuration has been
derived. Newton-Euler equations have been used to develop the
dynamic model. The simulation are carried out in
MATLAB/SIMULINK environment and different PID controller
tuning methods have been looked into for smoother control
Fig. 1. Control System overview flowchart.
Keywords- Modelling, Matlab Simulink, PID Controller, Relevant research papers have been referred to, for primary
Quadrotors, Simulation information along with other articles, videos and courses which
I. INTRODUCTION have been referenced [2] [3] [4].

Quadrotors are a class of Unmanned Aerial Vehicles II. SCOPE OF THIS PAPER
(UAVs) which are small and offer great manoeuvrability. They The goal of this paper is to utilize existing work on the
are also classified as Vertical Take-Off and Landing (VTOL) mechanics and dynamics of quadrotors [5] [6] to successfully
vehicles. Quadrotors are built for several different purposes model and simulate a quadrotor in Simulink Matlab [7]. We will
encompassing carrying payloads and deliveries, surveillance, be ignoring some of the more complicated effects while
recreational purposes, photography, videography, and for modelling the quadrotor. These can be found in [8] Furthermore,
defence applications to name a few. While the structure of the this paper will begin with modelling and simulating a quadrotor
quadrotor changes according to its purpose, the general in 1-dimension and build to a 3-dimensional simulation to allow
structure remains the same. A quadrotor, aptly suggested by its for proper understanding on which the PID [9] algorithm will be
name, has 4 motors which are equidistant from its centre, which implemented to achieve attitude and altitude control in each
is where the centre of mass lies. They make use of the 4 motors case. Gain tuning will be discussed briefly for PD and PID
for take-off, landing and other kinds of manoeuvring in the air. models in each case. Furthermore, few graphs will be shown to
As estimated in a study by Goldman Sachs, the drone highlight the difference between a nonlinear plant and a
market is estimated to be worth $100 billion by the end of 2020, linearized one. Finally, desired outputs versus actual outputs will
with the defence sector being responsible for $70 billion, be looked at and the degree of success of the modelling and
consumer market being worth $17 billion and the commercial simulation will be determined.
and civil market being worth $13 billion. All major companies
are trying to automate their processes with robotics for
enhanced efficiency and accuracy. Quadrotors, being reliable,
are being used in numerous sectors. Further details on the
application and scope of the industry can be found in [1]. Below

978-1-7281-5875-4/21/$31.00 ©2021 IEEE

Authorized licensed use limited to: BOURNEMOUTH UNIVERSITY. Downloaded on June 26,2021 at 18:51:30 UTC from IEEE Xplore. Restrictions apply.
 2021 International Conference on Computer Communication and Informatics (ICCCI -2021), Jan. 27 – 29, 2021, Coimbatore, INDIA

III. 1-D MODELLING AND SIMULATION OF QUADCOPTER We will use the provided block and after tweaking the gain
coefficients a bit, the following results were the ones that lead
A. Mechanics and Quadrotor Dynamics
to an acceptable response.
For one dimensional analysis of the drone, we shall
investigate the altitude control such that the quadrotor is only D. System Response using PD and PID controllers
allowed to move in the Z direction. Furthermore, the drone will
be considered symmetrical around the centre of mass and the X
and Y axes. The equation of motion along the Z axis for a
quadrotor can be derived from the figure below.
TABLE I. PID GAINS OF SYSTEM ABOVE

Constants Gains

Kp 50
Ki 25
Kd 10

Fig. 2. 1 Dimensional quadrotor force body diagram

′′ = 1 − (1)

Here U1 is the thrust generated by the drone motors. Mg is

the force applied by gravity and Ma is the net acceleration
caused. z’’ is the net acceleration in z direction.
B. Simulink Implementation
To successfully simulate a system, there are two parts one
must take care of, the plant system, and its control. The plant
describes the behaviour of the quadrotor and will be replicated
by using equation (1), which will give the Z value for Fig. 5. Step response with PD implementation
calculating the error and feeding it to the PID controller. The
controller will be tuned to minimize error between the desired The system response using PD and PID controllers is shown
and actual Z value. in figure 6 and 7 respectively. While the PD controller leaves a
steady state error, a PID controller causes the peak overshoot to
be higher and a decrease in the rise time. To compensate for
this, often the PID controller has a lower P gain than a PD
controller.
Compared to the PD graph above, the PID graph barely
Fig. 3. Implementation of (1) has any steady state error and a slightly higher peak overshoot
With the PID controller implemented, the overall system looks which can be combated by using a lower P gain or a reduced
like this. D gain.

Fig. 4. 1D PID implementation on Quadrotor plant

C. Gains and PID

Equation (1) will be used as the plant equation and a control
system will be developed around it to get the desirable output.
PID algorithm or the PID controller is a controller that uses Fig. 6. Step response with PID implementation
the sum of the error multiplied by the coefficients Kp, Ki and Kd
to elicit a response well within the design parameters of the IV. 2-D MODELLING AND SIMULATION OF QUADROTOR
system. MATLAB Simulink provides a readymade PID block, A. Mechanics and Quadrotor Dynamics
but one could easily make a block as per the applications.
In this section, the mechanics of a quadrotor fixed in the ZY

Authorized licensed use limited to: BOURNEMOUTH UNIVERSITY. Downloaded on June 26,2021 at 18:51:30 UTC from IEEE Xplore. Restrictions apply.
 2021 International Conference on Computer Communication and Informatics (ICCCI -2021), Jan. 27 – 29, 2021, Coimbatore, INDIA

plane will be explored. All our previous discussions about the

symmetry and behaviour still hold. Although now, the
quadrotor can move freely in the ZY plane which leaves us with
three variables z, y and ɸ (Roll Angle).

Fig. 7. Planar Quadrotor diagram Fig. 8. Implementation of nonlinear plant system

From the free body diagram in fig 7, one can derive the Finally, after implementing the PID controller on the plant
equations for y’’, z’’ and ɸ’’ in which U1 is the thrust input system we can observe the step responses.
and U2 is the moments input.
D. Gains and PID
′′= − 1 sin(ɸ) (3)
In the 2D system, there are 3 PID blocks to control the
error of all the 3 variables that are being worked with – z, y
′′ = − + 1 cos(ɸ) (4) and ɸ. The aim will be to control z and y and let ɸ vary
independently to provide the desired y output since y and ɸ
ɸ′′ = 2
(5) are dependent on each other. Another important thing to take
care of in nested PID loops is that the internal loop should
−sin(ɸ)
0 always have higher gains than the external ones to provide
′′ 0 effective responsiveness. Below are the gains for the PID
cos(ɸ) 1
′′ = − + 0 (6) loops.
2
ɸ′′ 0 1
0
TABLE II. PID GAINS OF SYSTEM ABOVE
Extending (3), (4) and (5) : ɸ y z

2 Kp 100 45 40
= ′′ = − 1 sin ɸ (7) Ki 0 20 30
2 Kd 100 10 15
2
= ′′ = − + 1 cos ɸ (8)
2 E. System Response of Linear and Nonlinear Systems
2ɸ
For the simulation described above, the system responses
= ɸ′′ = 2 (9) are shown. Below are the system responses for the altitude
2
control of the linear system.
B. Linearization
In the situation of the quadrotor hovering, the angle ɸ
changes nominally. Hence one can linearize the above
equations by applying the small angle approximations:
2
= ′′ = − ɸ (10)
22
= ′′ = − + 1
(11)
2

Fig. 9. Step response of z in linear system

Equation (9) is already in the linearized form.
C. Simulink Implementation We can compare this approximation to the actual nonlinear
Now we shall compare with the linear model and will implementation. Below is the system response for altitude
provide an insight into how well the approximation holds. For control for the nonlinear implementation of the plant. As we
the linear plant system, (9), (10), and (11) will be implemented can see the settling time is higher and the responsiveness of
and for the nonlinear plant system, equations (7), (8), and (9) the system seems a bit more delayed.
will be implemented. Below is the nonlinear implementation.

Authorized licensed use limited to: BOURNEMOUTH UNIVERSITY. Downloaded on June 26,2021 at 18:51:30 UTC from IEEE Xplore. Restrictions apply.
 2021 International Conference on Computer Communication and Informatics (ICCCI -2021), Jan. 27 – 29, 2021, Coimbatore, INDIA

Now one can look at the moments acting on the quadrotor (Mb)

i = 1 2Ω 2
i (15)
2

i =1 2Ω 2
i (16)
2

Where:

ρ = Air Density

A = Area swept by the blades

Ct, Cd = Coefficients

Fig. 10. Step response of z in nonlinear system R = Radius of the blades

The step response of y is almost identical in both the Ωi = Angular Velocity of motor
nonlinear and linear implementation and for the sake of
brevity is not included here. Equations (15) and (16) can be written as

V. 3-D MODELLING AND SIMULATION OF QUADROTOR i = ƒ Ωi 2 (17)

A. Mechanics and Quadrotor Dynamics = Ω i2 (18)
i
In this section, the mechanics of how a quadrotor functions
in the 3-dimensional space is explored. The position of the
quadrotor can be defined by [x y z] and the orientation of the
quadrotor can be described by [ɸ θ ψ] which constitute the roll,
pitch and yaw and represent rotations about the X, Y and Z axes,
respectively. The distance from the body frame to the inertial
frame will be represented by r = [x y z] T. If one describes a
movement with roll first, followed by pitch and finally yaw then
the multiplication of each of the individual matrices will lead to
the following matrix.

ƒ ɸ ƒ ɸ ƒ+ ɸ ƒ
= ƒ ɸ ƒ+ ɸ ƒ ɸ ƒ− ƒ (12)
− ɸ ɸ Fig. 11. Forces and Moments acting on the Quadrotor

c is the shorthand for cosine and s is the shorthand for sine here. Now calculating the moment about the X,Y and Z axes from
The rotational equation of motion can be written as follows Fig.13

= ′ + + (13) = ( 4 − 2) = ƒ (Ω4
2 − Ω2 2 ) (19)
= ( − )= (Ω 2 −Ω 2) (20)
1 3 ƒ 1 3
Where: = 1- 2+ 3 - 4 =
J = Inertial Matrix which because of the symmetry around the (Ω12 − Ω22 + Ω32 − Ω42 ) (21)
centre is Combining (19), (20), and (21)
Ixx 0 0
ƒ (Ω4 − Ω2 2 )
2
0 Iyy 0
Mb = ƒ (Ω1 − Ω3 )
2 2 (22)
0 0 Izz
(Ω12 − Ω22 + Ω32 − Ω42 )
ω = Angular Rates of Velocity
The translational motion of equation take the general form:
Mg = [0 0 JrΩr]T = Gyroscopic moments (14)
0
Ωr = Relative motor speed. ′′ = 0 + (23)
Mb = Moments acting on the body frame

Authorized licensed use limited to: BOURNEMOUTH UNIVERSITY. Downloaded on June 26,2021 at 18:51:30 UTC from IEEE Xplore. Restrictions apply.
 2021 International Conference on Computer Communication and Informatics (ICCCI -2021), Jan. 27 – 29, 2021, Coimbatore, INDIA

Fb represents all the non-gravitational forces acting in the body These provide the rotational equations. The translational
frame and is given by equations can be derived from expanding (23) and using (12),

0 0 ƒ ɸ ƒ ɸ ƒ+ ɸ ƒ 0
= 0 + ƒ ɸ ƒ+ ɸ ƒ ɸ ƒ− ƒ 0
Fb = 0 (24) − ɸ ɸ − 1
- ƒ (Ω1
2 + Ω2 2 + Ω3 2 + Ω 4 2 ) Expanding and rewriting,

Fbz = F1 + F2 + F3 + F4 (25) ′′ = − 1 ( o ɸ i o ƒ + i ɸ i ƒ) (36)

which converts to the equation above. Positive z direction is
′′ = − 1 ( o ɸ i i ƒ+ i i ƒ) (37)
assumed in direction of gravity.

B. State Space Representation ′′ = − 1( o ɸ o ) (38)

Now that the forces and moments acting on the body of the
drone have been discussed, one can convert these equations into
state space form which will be implemented to simulate the
behaviour. State vector X can be represented as [x1 x2 x3 x4 x5 This can be rewritten in state space form as,
x6 x7 x8 x9 x10 x11 x12]T which can be mapped to X = [ɸ ɸ’ θ θ’
ψ ψ’ x x’ y y’ z z’]T ′′ = − 1 ( o 1 i 3 o 5 + i 1 i 5) (39)
The control input U = [U1 U2 U3 U4] can be written as
shown below ′′ = − 1 ( o 1 i 3 i 5 + i 3 i 5) (40)
U1 = ƒ (Ω1
2 + Ω2 2 + Ω3 2 + Ω4 2 ) (26) ′′ = − 1 ( o o ) (41)
1 3

U2 = ƒ (Ω4
2 − Ω2 2 ) (27) C. Simulink Implementation
U3 = ƒ (Ω1
2 − Ω3 2 ) (28) The parameters used for the simulation are listed below.
Parameter Value Unit
U4 = (Ω12 − Ω22 + Ω32 − Ω42 ) (29)
L 0.24 m
From equations (26) to (29), a coefficient matrix can be formed M 1 kg
g 9.81 m/s2
ƒ ƒ ƒ ƒ Ω12 Ixx 0.01 kg m2
1
0 − ƒ 0 Ω22 Iyy 0.01 kg m2
2 ƒ
= (30) Izz 0.02 kg m2
3 ƒ 0 − ƒ 0 Ω32 -6
Kf 1.87 x 10 -
4 − − Ω42 Km 1.87 x 10-6 -
Equation (30) can be rearranged so that the coefficient matrix Equations (32), (33), and (34) will be implemented to
achieve the roll pitch and yaw accelerations which in turn will
for motor speeds can be found from control vector U. provide the roll, pitch and yaw angles.
Ω 12 1/4Kf 0 1/2Kf 1/4Km U1

Ω22 = 1/4Kf -1/2Kf 0 -1/4Km U2 (31)

Ω 32 1/4Kf 0 -1/2Kf 1/4Km U3
Ω42 1/4Kf 1/2Kf 0 -1/4Km U4

Expanding and solving (13),(14), and (22)

ɸ′′ = 1 ( 2 − 4Ω +( − ) 4 6) (32)

′′ = 1
( 3 − 2Ω +( − ) 2 6) (33) Fig. 12. Implementation to get roll acceleration, ɸ’’.

Similarly one can implement the other equations for angular

ƒ′′ = 1 ( 4 +( − ) 2 4) (34) and linear accelerations. Equations (36), (37), and (38) will be
Izz implemented to give us the outputs x, y and z accelerations
which can be integrated twice to provide the actual positions.

Authorized licensed use limited to: BOURNEMOUTH UNIVERSITY. Downloaded on June 26,2021 at 18:51:30 UTC from IEEE Xplore. Restrictions apply.
 2021 International Conference on Computer Communication and Informatics (ICCCI -2021), Jan. 27 – 29, 2021, Coimbatore, INDIA

Note that the following are the nonlinear implementations and

to implement a linearized linear accelerations subsystem, one
can apply the small angle approximations described previously.

Fig. 16. Step response of pitch angle, θ

Fig. 14. Implementation to get relative motor speed, Ωr from control inputs,
Ui.
Now one can connect these subsystems with each other so
that the behaviour of the drone i.e., the plant system, is now
entirely modelled. Once the plant system is ready, the PID
control system will be implemented on the plant to achieve
altitude and attitude control. The gains used are mentioned in
the next section.

D. Gains and PID

To achieve altitude control a PID controller has been Fig. 17 Step response of yaw angle, ψ
applied on U1 and to achieve attitude control 3 PID controllers
have been applied on U2, U3 and U4 to get the desired ɸ, θ, and
ψ values. The gains for all the PID controller loops have been The responses for altitude control are given in the figures 18
shown below. These gains can be different based on what and 19 for linear and nonlinear implementation, respectively.
parameters are acceptable to the user. These controllers can be
just PD controllers as well, which was demonstrated earlier in
2D simulation.
TABLE III. PID GAINS OF SYSTEM ABOVE
z ɸ θ ψ
Kp 800 110 16 8
Ki 300 50 10 10
Kd 25 3 10 2

E. System Response for Attitude and Altitude Control

This section contains the system responses from the
implementation of the 3D plant and one can observe the
efficiency of the implemented control system used to achieve Fig. 18. Step response of Z in linear case
altitude and attitude control of the quadrotor.

Fig. 15. Step response of roll angle, ɸ Fig. 19. Step response of Z in non-linear case

Authorized licensed use limited to: BOURNEMOUTH UNIVERSITY. Downloaded on June 26,2021 at 18:51:30 UTC from IEEE Xplore. Restrictions apply.
 2021 International Conference on Computer Communication and Informatics (ICCCI -2021), Jan. 27 – 29, 2021, Coimbatore, INDIA

It is observed that attitude control implementations are Lastly, motion planning and environment tracking remain areas
identical in both the linear and nonlinear models. However, that need to be explored further for nuanced control of the drone
altitude control has a huge difference and so the small angle in its surrounding environment.
approximation is not one that is appropriate here. Further
methods or techniques must be looked into to provide a REFERENCES
smoother and more robust altitude controller for the 3D
nonlinear plant system. The results from the simulation are
tabulated below. [1] L. Camilli, “Emerging technologies, applications, regulations, and market
challenges in the consumer aerial drone industry,” Technical Report, San
TABLE IV. RESULTS OF SIMULATIONS Francisco State University College of Business, 2015.
Attribute Settling Rise Time Peak [2] O. M. Hamdy, and A. O. El-Nady, “Modeling, Simulation and Control of
Time (ms) (ms) Overshoot (%) Quadcopter using PID controller,” 4th International under-Graduate
Research Conference, IUGRC 2019, Cairo, 2019.
Nonlinear 413.0 179.7 14.2 [3] A. T. Bayisa, “Controlling Quadcopter Altitude using PID-Control
System,” in International Journal of Engineering Research & Technology
Altitude (IJERT), Vol 8, pp. 195-199, 2019.
Linear Altitude 214.9 35.8 13.0 [4] Q. Huang, “Mathematical Modeling of Quadcopter Dynamics,” Report,
Rose Human Institute of Technology, 30 May 2016.
Roll 195.7 52.5 12.1
[5] “MathWorks,” [Online]. Available: www.mathworks.com [Accessed
Pitch 837.2 9.3 30.2 May 2020]
[6] University of Pennsylvania, “Robotics: Aerial Robotics,” Coursera,
Yaw 499.5 18.8 21.7 [Online]. Available: https://www.coursera.org/learn/robotics-
flight?specialization=robotics. [Accessed June 2020].
[7] V. Desh, “MATLAB & Simulink Tutorial: Quadrotor UAV Trajectory
VI. CONCLUSION AND FUTURE RESEARCH and Control Design (PID + Cascaded),” YouTube, [Online]. Available:
https://www.youtube.com/channel/UCRaiKJHsh_Ii4H-86uqXT_g.
As seen, one can quite accurately simulate the behavior of a [Accessed July 2020]
quadrotor discounting a few effects that would unnecessarily [8] T. Lukkonen, “Modelling and Control of quadcopter,” Thesis, Aalto
complicate the physics. While the simulations are fairly University - School of Science, 22 August 2011.
accurate and the PID algorithm is good enough to elicit the [9] A. Zulu, and S. John “A Review of Control Algorithms for Autonomous
desired response from the drone, this only provides us with a Quadrotors,” Open Journal of Applied Sciences, pp. 547-556, 2014.
starting point while building and testing a drone. Furthermore, [10] J. Anand, "Automatic Traffic Control Technologies for Remote
the nonlinear altitude control needs to be improved. It would Monitoring of Unmanned Railway Gates," Journal of Electronics, pp. 30-
37, 2020.
require a lot of further tweaking of the gains as the physics that
[11] E. Kamel, and K. Kamel, "PLC Batch Process Control Design and
has been used for the simulation involve approximations. Implementation Fundamentals," Journal of Electronics, pp. 155-161,
Upcoming research will be based on the different optimization 2020.
tools to precisely predict the gain values such as reinforced
learning, genetic algorithms, and swarm optimization etc.

Authorized licensed use limited to: BOURNEMOUTH UNIVERSITY. Downloaded on June 26,2021 at 18:51:30 UTC from IEEE Xplore. Restrictions apply.