FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-1

POLE-ZERO MAP OF A DYNAMIC SYSTEM MODEL WITH PLOT

CUSTOMIZATION

AIM:
To draw a Pole-Zero Map of a Dynamic System Model

THEORY

A system is often defined in terms of the poles and zeros of its transfer function. The poles of
a transfer function are the values of the Laplace transform variable, s, that cause the transfer
function to become infinite. The zeros of a transfer function are the values of the Laplace
transform variable, s, that cause the transfer function to become zero

A system can be defined by its transfer function, which is a ratio of polynomials in the Laplace
variable "s.

The above transfer function H(s) can be represented in the form of Pole-Zero representation
as

The constant k=b0/a0. The zi terms are the zeros of the transfer function; as s→zi the numerator
polynomial goes to zero, so the transfer function also goes to zero. The pi terms are the poles
of the transfer function; as s→pi the denominator polynomial is zero, so the transfer function
goes to infinity.

A pole–zero plot is a graphical representation of a rational transfer function in the complex

plane which helps to convey certain properties of the system such as stability, Causality etc.,.

In Mat lab, pzplot (sys) computes the poles and (transmission) zeros of the dynamic system
model sys and plots them in the complex plane. The poles are plotted as X's and the zeros are
plotted as O's. It also returns the plot handle h. You can use this handle to customize the plot
with the get options and set options commands

1 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

Example :Plot the poles and zeros of the continuous-time system represented by the
following transfer function:

MATLAB CODE:

clear all;
clc;
sys = tf([2 5 1],[1 3 5]);
h = pzplot(sys);
grid on
p = getoptions(h);
p.Title.Color = [1,0,0];
setoptions(h,p);

Problem 1: Plot the poles and zeros of the continuous-time system represented by the
following transfer function with plot customization on Xaxis and Yaxis. Sketch the results

Problem 2: Plot the poles and zeros of the continuous-time system represented by the
following transfer function with plot customization. Sketch the results
2s2 + s + 3
G(s) =
3s3 + 4s2 + 1

VIVA QUESTIONS

1. What is MatLab?
2. Name some of the basic commands in MatLab?
3. What are all the functions of the commands in the above program?
4. How MatLab is superior to other software simulation?
5. What is a Pole and Zero?
6. What is the prediction obtained from a Pole zero plot?

2 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-3

ROOT LOCUS FOR A DYNAMIC SYSTEM.

AIM:
a. To Plot the root locus with variables in transfer function
b. To Plot the root locus of a dynamic system

THEORY

Root locus analysis is a graphical method for examining how the roots of a system change with
variation of a certain system parameter, commonly a gain within a feedback system. This is a
technique used as a stability criterion in the field of classical control theory developed by
Walter R. Evans which can determine stability of the system. The root locus plots the poles of
the closed loop transfer function in the complex s-plane as a function of a gain parameter.

Locus is defined as a set of all points satisfying a set of conditions, the term root refers to the
roots of the characteristic equation which are the poles of the closed loop transfer function.
These poles define the time response of the system and hence the performance and stability of
the system.

The root locus of a feedback system is the graphical representation in the complex s-plane of
the possible locations of its closed-loop poles for varying values of a certain system parameter.
The points that are part of the root locus satisfy the angle condition. The value of the parameter
for a certain point of the root locus can be obtained using the magnitude condition.

Suppose there is a feedback system with input signal X(s) and output signal Y(s). The forward
path transfer function is G(s), the feedback path transfer function is H(s)

For this system, the closed-loop transfer function is given by

3 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

Thus, the closed-loop poles of the closed-loop transfer function are the roots of the
characteristic equation 1+G(s)H(s)=0. The roots of this equation may be found wherever
G(s)H(s)=-1

In systems without pure delay, the product G(s)H(s) is a rational polynomial function and may
be expressed as

Where num is the numerator of the polynomial and den is the denominator polynomial and K
is the gain (K>0). The vector K contains all the gain values for which the closed loop poles are
to be computed

The root loci is plotted by using the MATLAB Command

rlocus(num,den)

The gain vector K is supplied by the user

Example: A plant to be controlled is described by a transfer function

+5
()= 2 + 7 + 25

Obtain the root locus plot using MATLAB

MATLAB CODE

clear all;
clc;
num=[1 5];
den=[1 7 25];
G=tf(num,den);
rlocus(G);

VIVA QUESTIONS

1. What is root locus?

2. What can be interpreted from the root locus analysis?
3. What is stability?
4. What are the steps involved in root locus method?
5. Which command is used to obtain the root locus of any dynamic system in MATLAB

4 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-5

SIMULATE A SPRING MASS – DAMPER SYSTEM WITH AND

WITHOUT A FORCING FUNCTION THROUGH SIMULINK.
AIM:
To generate and analyze the response of a spring mass motion to a step and ramp input

THEORY

The Mathematical model of the system is described by:

1
x=M (− Bx − Kx + f (t ))
Let Mass M=2kg; Stiffness Constant K=2 N/m; B = 2 Ns/m. Then,

1
x= 2 (−2 x − 2x + f (t ))

PROCEDURE FOR SIMULATION:

Set f(t)=0 for free oscillation response or set it to a constant value. For f(t)=0 the equation
becomes:
̈+ +
̇ =0

STEP 1: Open MATLAB, Type “simulink” in command window.

STEP 2: Go to File, Open a New Model.
STEP 3: Drag and place the following blocks from the corresponding library

5 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

Select BLOCK set Location in Simulink Library

Step Sources

Sum Math Operation

Gain Math Operation

Integrator Continuous

Scope Sinks

STEP 4: Make Connections

STEP 5: Set Parameters

6 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

STEP 6: Running Simulations

1 Open “simulation parameters”

2
RUN
Set “Stop time” =

Simulatio

STEP 7: Viewing Results

7 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

VIVA QUESTIONS

1. What is Simulink?
2. What are types of math operators are used in MatLab?
3. What is a mathematical model and how it is helpful?
4. What are the different types of inputs given to the system?
5. What is damping ratio and what does it signifies?
6. What is stiffness?
7. What is viscous damping?

8 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-6

SIMULATION OF SIMPLE SERVOMECHANISM FEEDBACK SYSTEM

AIM:
To simulate a simple servo-mechanism feedback system in “S” domain and plot its response
to a step input.

THEORY:

Consider a servomechanism system as shown in the following figure.

r(t) e(t) R θ(t)

T
A L SH
+ c(t)

Where,
A = Gain of Amplifier
R = Resistance of Motor Winding
L = Inductance of Motor Winding
SH = Shaft

It consists of the following elements,

(a) An electrical input which acts as a “reference” of shaft angle. It is referred as
r(t). In the „S‟ domain this can be represented by R(s).
(b) A differential amplifier which compares the output c(t) i.e the controlled shaft angle
withy input reference proportional to shaft angle.
Mathematical equations for this may be described
as, e(t) = r(t)-c(t)
where r(t) is the reference input, c(t) is the controlled shaft angle and e(t) is the
output of the differential amplifier which acts as an input to the amplifier.
In the „S‟ domain this can be represented by,
E(s) = R(s)-C(s)
(c) An output amplifier amplifies the output of this differential amplifier. Its input is e(t)
and output is v(t). Its mathematical equation may be described as,
v(t) = A e(t)

9 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

In the „S‟ domain this can be represented by,

V(s) = AE(s)
(d) An electrical d.c shunt motor which controls the controlled shaft. The torque T,
generated by the motor proportional to the field current and represented by,
T/I = K
Where K is a constant determined by the characteristics of the motor. The ratio of the
field current I and amplifier output is,
V(t) = R.I + L (dI/dt)
where L and R are the inductance and resistance respectively of the motor field
winding. In the „S‟ domain this can be represented by,
V(s) = R. I(s) + L S I(s)
(e) The equation of the shaft angle θ(t) can be written by the application of Newton‟s
Second Law as,
T = L (d2 θ(t) / dt2) + F (d θ(t)/dt)
In the „S‟ domain this can be represented by,
T(s) = JS2θ(s) + FS θ(s)
= (JS2 + FS) θ(s)
The output voltage is c(t) = B θ(t), and in „S‟ domain this can be written
as C(s) = B θ(s)
From the above equations we can obtain,
C(s) / T(s) = (1/S) . (K1 / (1+STL))
T(s) / V(s) = K2 / (1+STM)
Therefore the above diagram can be reduced to its block diagram representation in “s”
domain
as,

10 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

PROCEDURE FOR SIMULATION:

STEP 1: Open MATLAB, Type “simulink” in command window.

STEP 2: Go to File, Open a New Model.

STEP 3: Drag and place the following blocks from the corresponding library
Select BLOCK set Location in Simulink Library
Step Sources
Sum Math Operation
Transfer Function Continuous
Scope Sinks

STEP 4: Make Connections as shown in the reduced block diagram above

STEP 6: Running Simulations

1 Open “simulation parameters”

2
RUN
Set “Stop time” =

Simulation

STEP 7: View the response plot by double clicking on the scope block

11 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

VIVA QUESTIONS

1. What is a frequency Response?

2. What is the gain of an amplifier?
3. What is the inductance of motor winding?
4. What are the different types of controllers used?
5. What is the advantage of taking Laplace Transform of parameter in Control system?
6. What is the requirement of ideal control system?
7. What is an open loop and closed loop control system? Differentiate.

12 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-7

SIMULATE A BOMB DROP FROM AN AIRCRAFT ON A MOVING

TANK FOR PURE-PURSUIT MOTION
AIM:
To simulate a bomb drop from an aircraft on a moving tank for pure pursuit motion and plot
the trajectory.

THEORY:

Assumptions:
1. Target and the pursuer are moving in the same horizontal plane when the pursuer first
sights the tank.
2. The Pursuer speed is constant.
3. Forward speed of the bomb remains constant.

Let us assume a pursuer (fighter) at a height yo from the ground moving in the same plane as
that of the target and also in the same direction. Let the pure pursuit velocity be vo. Let the
angle of elevation with respect to the horizon of the fighter be θ. If t is the total time of flight
for the trajectory, then equation of motion can be represented by,
−g
y = t 2 + (V sin )t
O  2  O

The above equation is of the form s = ut + (1/2)at2. Thus solving the above equation we get
the total time of trajectory for the bomb drop. The velocity vector along the x axis is given by
• •

x = VO cos( ) and that along the y axis is given by y = VO sin( ) . Therefore the new position
•
vector at any instant along the x axis is given by x = x t and that along the y axis is given by,
 −g2 •
y= t + y t + yO
 2 

13 Department of AE, MITE.

 FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

FLOW CHART:

Start

Initialize Fighter altitude, velocity, angle of

elevation

Obtain Trajectory time from the standard equation of motion. Obtain velocity vector along the x
axis and y axis Over different instants of trajectory time, obtain the position of bomb at respective
x and y coordinates

Plot the coordinates

End

MATLAB CODE:
clc;
clear all;
y0=5; % initial height in m
v0=10; % initial velocity in m/s
theta=0; % trajectory elevation w.r.t horizon in deg
g=9.81; % acceleration due to gravity in m/s^2
b=v0*sin(pi*(theta/180)); a=-g/2;

c=y0;
t_flight=(-b-sqrt(b^2-4*a*c))/(2*a);
t=linspace(0,t_flight,30);
xdot0=v0*cos(pi*(theta/180));
ydot0=v0*sin(pi*(theta/180));
x=xdot0*t+x0;
y=-(g/2)*t.^2+ydot0*t+y0;
plot(x,y);
xlabel('Distance (m)');
ylabel('Height (m)');
title('Bomb Drop Trajectory for pure pursuit Motion');

14 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

Result:
Case1: v0=10;
Case2: v0=0;
VIVA QUESTIONS
1. What is the initial velocity of the bomb when it is dropped?
2. What are the forces acting on a falling bomb?
3. What are the directions of the above forces?
4. What are the uncertain forces acting on the falling body?
5. What are the effects of these forces?
6. What is the shape traced by the bomb falling from aircraft under ideal condition?
7. Explain pure pursuit motion?
8. What is the function of linspace in MatLab?

15 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-8

DEVELOP A STRAIGHT AND LEVEL FLIGHT SIMULATION

PROGRAM USING MATLAB

AIM:
To Develop straight and level flight simulation by comparing forces and moment acting on
the aircraft.

THEORY:
For steady flight along a straight line, whether level or not, it is not only necessary to balance the
four forces so that they produce no resultant force. Their lines of action must also be arranged so
that they produce no resultant moment, otherwise the aircraft will rotate either nose-up or nose-
down. When there is no resultant moment, the aircraft is said to be trimmed. Mathematically, if the
sum of the moments is M, then the condition for trim is that M=0. One way to achieve this would
be to arrange that all of the forces act through a single point, however this is not generally practical,
as there are many factors that tend to alter the line of action, apart from those already stated.
For example, lowering the undercarriage tends to shift the line of the resultant drag down, and on
the floatplane variant of a light aircraft the position of the drag resultant would be very much lower
than for the original float less design.
It is possible to balance the moment produced by the drag and thrust being out of line, by arranging
the lift and weight forces to be out of line by an amount that causes them to exactly produce the
necessary balancing moment, as in Fig. However, because of the way that the lines of action of the
forces tend to change according to aircraft attitude, fuel weight etc., there is no simple design
solution to ensure that the resultant moment will always be zero. Some active involvement of the
pilot is required in order to keep the aircraft trimmed.

1. Tail-Less Airplane:

This type has had followers from the very early days of flying – and among birds from prehistoric
times – and although the reasons for its adoption have changed somewhat, a common feature has
been a large degree of sweepback, or even delta-shaped wings, so that although this type may
appear to have no tail, the exact equivalent is found at its wing tips, the wings being, in fact, swept
back so that the tip portion can fulfil the functions of the tail plane in the orthodox airplane.

16 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

2. Tail Planes:

Where the four main forces can be satisfactorily balanced in themselves, the duty of the tail plane
is merely to act as a “stand-by”. Therefore, it will usually be set at such an angle, that at cruise
speed it will be at zero angle of attack, thereby producing no lift. At flight speeds higher than the
cruise speed, the lift coefficient must be lowered to compensate for the higher dynamic pressure
otherwise the lift would be greater than the weight. This means that the aircraft must be trimmed a
little more nose-down. In doing so, the centre of lift of the wing will move back, giving a nose-
down pitching moment.

However, if the angle of attack of the tail plane was zero for cruise, then it will now become
negative, and the tail will generate a down-force producing a nose-up pitching moment which will
tend to more than counteract the nose-down moment produced by the wing. In fact, as the flight
speed increases, the pilot normally has to make a small nose-down trim adjustment.
Correspondingly, at low speeds, the nose of the aircraft must be raised in order to increase the angle
of attack. This means that the tail lift will now become positive (in Fig). As the tail plane is equally
likely to carry an upward or a downward force, it is usually of symmetrical camber, and therefore
produces no lift at zero angle of attack. On a tail-first or canard aircraft, the foreplane is set at a
slightly higher angle of attack than the wings for reasons of stability, and both wings produce lift
in normal flight.

Conditions of balance considered numerically:

1. Lift =Weight, i.e. L =W.

2. Thrust=Drag, i.e. T =D.

3. The ‗nose-down„ pitching moment of L and W must balance the ‗taildown„ pitching moment of
T and D.

The two forces, L and W, are two equal and opposite parallel forces, i.e. a couple; their
moment is measured by ‗one of the forces multiplied by the perpendicular distance
between them‘. So, if the distance between L and W is x metres, the moment is Lx or Wx
newton-metres. Similarly, T and D form a couple and, if the distance between them is y
metres, their moment is Ty or Dy newton-metres. Therefore the third condition is that –

Lx(or Wx)= Ty(or Dy)

Case1: Without Tail Effect

17 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

Case 2: With Tail Effect

MATLAB Code:
clc;
fprintf('Without Tail Effect')
lift=20000;
thrust=2000;
z=thrust/lift
fprintf('With tail effect')
lift2=10000;
thrust2=1250;
x=0.025;
y1=0.15;
t=6;
drag=1250;
syms y P
e=y+P==lift2;
e2=(x*y)+(t*P)==(y1*drag);
sol=solve([e e2],[y P]);
wing_lift= sol.y
tail_lift= sol.P
fprintf('negative sign indicates the assumed direction
is opposite\n')

Viva Questions

1. Explain Straight and Level Flight?

2. Explain the Condition for level flight?
3. Explain the parameter affecting the level flight?
4. Explain Static Margin?
5. Explain the Procedure to balance Moment in aircraft?
6. Give the example of Tail-less planes?
7. Explain the procedure to develop straight level flight in tail and tail-less aircrafts?

18 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-9

SIMULATE AIRCRAFT TAKE-OFF AND LANDING WITH

TRAJECTORY TRACING

AIM:
To calculate Lift – off and Landing distance of Aircraft and also study the effect due to
other parameters.

THEORY:
1. Takeoff Performance:
The total takeoff distance, as defined in the Federal Aviation Requirements
(FAR), is the sum of the ground roll distance SLO and the distance (measured
along the ground) to clear a 35-ft height (for jet-powered civilian transports) or a
50-ft height (for all other airplanes).

Lift off distance

1. Liftoff distance is very sensitive to the weight of the airplane, varying directly as W 2.
If the weight is doubled, the ground roll of the airplane is quadrupled.
2. Liftoff distance is dependent on the ambient density If we assume that thrust is directly
proportional to

3. The liftoff distance can be decreased by increasing the wing area, increasing
CL(max), and increasing the thrust, all of which simply make common sense.

Velocity Liftoff:

Lift:

Drag:

19 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

2. Landing Performance:
The total landing distance, as defined in FAR, is the sum of the ground roll distance
and the istance (measured along the ground) to achieve touchdown in a glide from a 50-ft
height. During the landing ground roll, the pilot is applying brakes. Modern jet transports
utilize thrust reversal during the landing ground roll. Thrust reversal is created by ducting air
from the jet engines and blowing it in the upstream direction. Opposite to the usual downstream
direction when normal thrust is produced.

Landing Distance:

TR is the reserved thrust

MATLAB Code:
clc;
w=10000:20:15000;
k=10^-5;
Slo=k*(w.^2);
plot(Slo,w);
title('takeoff distance v/s weight:')
xlabel('Liftoff Distance')
ylabel('Weight')
figure(2)
s=90:1:100;
k=10^5;
Slo=k./s;
plot(Slo,s);
title('takeoff distance v/s wing aera:')
xlabel('Liftoff Distance')
ylabel('Wing Area')
figure(3)
w=10000:20:15000;
k=10^-5;
Slo=k*(w.^2);
plot(Slo,w);
title('Landing distance v/s weight:')
xlabel('Landing Distance')
ylabel('Weight')
figure(4)
s=90:1:100;
k=10^5;
Slo=k./s;
plot(Slo,s);
title('Landing distance v/s wing aera:')
xlabel('Landing Distance')

20 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

ylabel('Wing Area')
Viva Questions:

1. Explain the modes of flight?

2. Explain Take-off performance?
3. Explain the factors affecting takeoff performance?
4. Explain the minimum power required for Takeoff?
5. Explain landing performance?

6. Explain the factors affecting Landing distance?

7. Explain the clearance height required for landing and takeoff?

21 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-10

SIMULATE STALL OF AIRCRAFT AND SHOWS THE EFFECT OF

VARIATION IN STATIC MARGIN ON STALL CHARACTERISTICS

AIM:

To Simulate stall Characteristics of aircraft and the effect of variation in Static Margin

THEORY:

1.Stall:

The factors that affect the stalling characteristics of the airplane are: balance (load distribution),
bank (wing loading), pitch attitude (critical angle of attack), coordination (control movement), drag
(gear or flaps), and power. The pilot should learn the effect of each on the stall characteristics of
the airplane being flown and what should be done to affect the proper correction. It should be
reemphasized here that a stall can occur at any airspeed, in any attitude, or at any power setting,
depending on the total number of factors affecting the particular airplane.

A number of factors may be induced as the result of other factors. For example, when the airplane
is in a nose high turning attitude, the angle of bank has a tendency to increase. This occurs because
with the airspeed decreasing, the airplane begins flying in a smaller and smaller arc. Since the outer
wing is moving in a larger radius and thus traveling faster than the inner wing, it has more lift and
causes an overbanking tendency. At the same time, because of the decreasing airspeed and
decreasing lift on both wings, the pitch attitude tends to lower. In addition, since the airspeed is
decreasing while the power setting remains constant, the effect of torque becomes more prominent,
causing the airplane to yaw.

During the practice of nose high turning stalls, to compensate for these factors and to maintain a
constant flight attitude until the stall occurs, aileron pressure must be continually adjusted to keep
the bank attitude constant. At the same time, back elevator pressure must be continually increased
to maintain the pitch attitude, as well as right rudder pressure increased to prevent adverse yaw
from changing the turn rate. If the bank is allowed to become too steep, the vertical component of
lift decreases and makes it even more difficult to maintain a constant pitch attitude.

Whenever practicing turning stalls, a constant pitch and bank attitude should be maintained until
the stall occurs. Whatever control pressures are necessary should be applied even though the
controls appear to be crossed (aileron pressure in one direction, rudder pressure in the opposite
direction). During the entry to a power on turning stall to the right in particular, the controls will
be crossed to some extent. This is due to right rudder pressure being used to overcome torque and
left aileron pressure being used to prevent the bank from increasing.

2.Static Margin:
If an aircraft in flight suffers a disturbance in pitch that causes an increase (or decrease) in angle
of attack, it is desirable that the aerodynamic forces on the aircraft cause a decrease (or increase)

22 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

in angle of attack so that the disturbance does not cause a continuous increase (or decrease) in
angle of attack. This is longitudinal static stability.

Static margin is a concept used to characterize the static longitudinal stability and controllability
of aircraft and missiles.

In aircraft analysis, static margin is defined as the distance between the center of gravity and the
neutral point of the aircraft, expressed as a percentage of the mean aerodynamic chord of the wing.
The greater this distance and the narrower the wing, the more stable the aircraft.

Conventionally, the neutral point is aft of the c.g., although in rare cases (computer-controlled
fighter aircraft) it may be forward of the c.g., i.e. slightly unstable, to obtain quickness of response
in combat. Too great longitudinal stability makes the aircraft "stiff" in pitch, resulting in such
undesirable features as difficulty in obtaining the necessary stalled nose-up pitch when landing.

The position of the neutral point is found by taking the algebraic net moment of all horizontal
surfaces, measured from the nose of the aircraft, in the same manner as the c.g. is determined, i.e.
the sum of all such moments divided by their total area. The stabilizer and elevator dominate this
result, but it is necessary to account for all surfaces such as fuselage, landing gear, prop-normal,
etc. It is also necessary to take account of the center of pressure of the wing, which can move a
good deal fore and aft as angle of attack of a flat-bottom wing section (Clark Y) changes, or not at
all in the case of self-stabilizing sections such as the M6.

The neutral point in conventional aircraft is a short distance behind the c.g. ("The feathers of the
arrow must be at the back"); but in unconventional aircraft such as canards and those with dual-
wings, such as the Quickie, this will not be so. The overall rule stated above must hold, i.e. the
neutral point must be aft of the c.g., wherever that may be.

Static margin (SM) reflects the distance between the

• center of mass (CM) and the net center of pressure (CP)
• Body Axes
• Normalized by mean aerodynamic chord
• Does not reflect Z position of center of pressure
• Positive SM if CP is behind CM.

Effect of Static Margin on Pitching Coefficient:

• Zero crossing determines trim angle of attack, i.e., sum of moments = 0
• Negative slope required for static stability
• Slope, varies with static margin.

23 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

MATLAB Code:
clc;
close all;
[v,T,vT]=xlsread('C:\Users\rohit\Documents\MATLAB\covsalpha.xlsx'
);
% 'xlsx' for
exell %v: Double
%T and vT : cell
%use v containing
numbers
alpha=v(:,1);cl=v(:,2);
cd=v(:,3); cm=v(:,4);
%if u have to plot second colone depending on
first: plot(alpha,cl)
grid on
title('CL v/s Alpha
(NACA0012)') xlabel('Alpha')
ylabel('Cl')
figure(2)
plot(alpha,cd
) grid on
title('Cd v/s Alpha
(NACA0012)') xlabel('Alpha')
ylabel('Cd')
figure(3)
plot(alpha,cm
) grid on
title('Cm v/s Alpha
(NACA0012)') xlabel('Alpha')
ylabel('Cm')

Viva Questions:
Q1. Explain Stall of an aircraft?
Q2. Explain the factor affecting Stalling in the aircraft?
Q3. Explain the precautions to be taken to avoid stall?
Q4. Explain Static Margin?
Q5. What is Neutral Point?
Q6. Explain longitudinal stability?
Q7. Explain how static margin affect the stability of an aircraft?
Q8. Explain the effect of static margin on stall characteristics?

24 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-11

SIMULATE AIRCRAFT LONGITUDINAL MOTION AND

DEMONSTRATE THE EFFECT OF STATIC MARGIN VARIATION FOR
A PULSE INPUT IN PITCH THAT IS INTENDED TO BLEED THE
AIRSPEED

AIM:
To simulate aircraft longitudinal motion and demonstrate the effect of static margin variation for
a pulse input in pitch that is intended to bleed the airspeed.

THEORY:
The simulation of the aircraft longitudinal motion would be obtained from the longitudinal
control transfer functions for elevator deflections a control inputs. Hence the first step towards
that is to obtain the said transfer functions.
The linearized small-disturbance longitudinal rigid body equation of motion is given by:

The expressions of stability parameters estimated for, are as follows:

25 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

The transfer function gives the relationship between the output of and input to a system. In the case
of aircraft dynamics it specifies the relationship between the motion variables and the control input.
The transfer function is defined as the ratio of the Laplace transform of the output to the Laplace
transform of the input with all the initial conditions set to 0. (i.e, the system is assumed to be
initially in equilibrium).

In this case the linearized equations of motion with the help of the stability parameters, will be
used to generate the longitudinal control transfer functions.

The longitudinal control transfer function will be having a general form of

Where,
= Numerator polynomial of the transfer function, forward velocity vs elevator deflection.
= Denominator polynomial of the longitudinal transfer function, which remains same for all
longitudinal control transfer functions.

Where in which the numerator and denominator polynomials can be written generally as:

The coefficients A, B, C, D and E belongs to the denominator polynomial and coefficients Au,
Bu, Cu and Du are of the numerator polynomial. Similar representation is available for the transfer
() ()

function such as ()
and ()

The table below shows the simplified expressions to estimate the coefficients of the transfer
functions.

Upon using the provided relationships with the basic flight parameters, stability derivatives and
aerodynamics coefficients the required transfer function can be estimated.

26 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

Impulse Response:
The aircraft time response can be summarized as follows:

Here, ( )is the input function, upon multiplying the input function with the system transfer function the Laplace transform of the output, U(s), will be
obtained. The inverse Laplace transform of the output would give the time function of the output u(t).

In this experiment, the system transfer function connects the input variable ( ), elevator deflection, and the controlled variable is the U(s) , the forward
velocity. Here, the input is an input function δ(t). Defined as:

The unit impulse is defined as a function of time that is zero everywhere, except for an
infinitesimally small neighbourhood around the origin, in which the function attains unbounded
values. However, its time integral from −∞ to +∞ is exactly +1. Thus, if we denote by 0− and 0+
the instants just before and just after 0, respectively, the unit-impulse function. The height of the
arrow, then, denotes the time integral of the associated impulse function on the whole real axis.
The Laplace transform of the function is equivalent to 1.

Case 1:
Static Margin = 0.2211
clc
A =1; % Coefficients in the characteristic equation A - E
B = 1.2632;
C = 0.9618;
D = 0.09786;
E = 0.0127;
Au = 0.9710; % Coefficients in the Numerator polynomial Au,
Bu, Cu, & Du
Bu = 0.0903;
Cu = 1.6882;
Du = 10.5012;
DTF = [A B C D E]; % Denominator polynomial of the TF
Nu = [Au Bu Cu Du]; % Numerator polynomial of the TF
TFu = tf([Nu], [DTF]) % TF: forward velocity vs elevator
deflection [U(s)/de(s)]

27 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

figure
impulseplot(TFu) % Impulse response of the TF
[U(s)/de(s)] grid on

Case 2:
Static Margin = 0.6633.
clc
A =1; % Coefficients in the characteristic equation A - E
B = 1.2632;
C = 2.0199;
D = 0.2121;
E = 0.0310;
Au = 0.9710; % Coefficients in the Numerator polynomial Au,
Bu, Cu, & Du
Bu = 0.0903;
Cu = 2.7156;
Du = 9.3144;
DTF = [A B C D E]; % Denominator polynomial of the TF
Nu = [Au Bu Cu Du]; % Numerator polynomial of the TF
TFu = tf([Nu], [DTF]) % TF: forward velocity vs elevator
deflection [U(s)/de(s)]
figure
impulseplot(TFu) % Impulse response of the TF
[U(s)/de(s)] grid on

Note:
The reduction in the forward velocity, u(t), for the impulse response upon varying, increasing, the
static margin is significant and easily observable.

Viva Questions:

1. How the Linearized Equations of Motion for an aircraft are obtained.

2. What are stability derivatives?
3. What are stability parameters
4. What do you understand by Time response of a system?
5. State few standard test signals.
6. What do you understand by Short period and Phugoid motion?
7. How is the time response of a first order system differ from a second or higher
order system, for a unit impulse input?
8. What effect the change in static margin bring in the forward velocity of an aircraft
for certain elevator control deflection.
9. Mention the physical meaning of implementation of an Impulse function as an
elevation command.
10. State few time response specifications applicable for higher order system and mention
its physical significance as well.

28 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-12
SIMULATE AIRCRAFT LONGITUDINAL MOTION AND
DEMONSTRATE THE EFFECT OF STATIC MARGIN VARIATION
FOR A DOUBLET INPUT IN PITCH

AIM:
To simulate aircraft longitudinal motion and demonstrate the effect of static margin variation for
a pulse input in pitch that is intended to bleed the airspeed.

THEORY:
The simulation of the aircraft longitudinal motion would be obtained from the longitudinal
control transfer functions for elevator deflections a control inputs. Hence the first step towards
that is to obtain the said transfer functions.
The linearized small-disturbance longitudinal rigid body equation of motion is given by:

The expressions of stability parameters estimated for, are as follows:

29 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

The transfer function gives the relationship between the output of and input to a system. In the case
of aircraft dynamics it specifies the relationship between the motion variables and the control input.
The transfer function is defined as the ratio of the Laplace transform of the output to the Laplace
transform of the input with all the initial conditions set to 0. (i.e, the system is assumed to be
initially in equilibrium).

In this case the linearized equations of motion with the help of the stability parameters, will be
used to generate the longitudinal control transfer functions.

The longitudinal control transfer function will be having a general form of

Where,
= Numerator polynomial of the transfer function, forward velocity vs elevator deflection.
= Denominator polynomial of the longitudinal transfer function, which remains same for all
longitudinal control transfer functions.

Where in which the numerator and denominator polynomials can be written generally as:

The coefficients A, B, C, D and E belongs to the denominator polynomial and coefficients Au,
Bu, Cu and Du are of the numerator polynomial. Similar representation is available for the transfer
() ()

function such as ()
and ()

The table below shows the simplified expressions to estimate the coefficients of the transfer
functions.

Upon using the provided relationships with the basic flight parameters, stability derivatives and
aerodynamics coefficients the required transfer function can be estimated.

30 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

Impulse Response:
The aircraft time response can be summarized as follows:

Here, ( )is the input function, upon multiplying the input function with the system transfer function the Laplace transform of the output, U(s), will be
obtained. The inverse Laplace transform of the output would give the time function of the output u(t).

In this experiment, the system transfer function connects the input variable ( ), elevator deflection, and the controlled variable is the U(s) , the
forward velocity. Here, the input is an input function δ(t). The time derivative of a unit impulse, δ (t), is called the doublet function. It is
formally defined as:

Using the definition of the time-derivative:

With ε being defined as: 2 ε ≡ 0+ − 0−

Therefore, the physical interpretation of a doublet is the limiting case of two impulses of ∞
amplitude, one applied at t = 0− and the other at t = 0+, the latter being the negative of the
former. As shown below:

The Laplace transform of doublet function, δ„(t), will be equals to ‗s„. Hence, upon applying the
doublet function as the input, the output definition would be simply the system transfer function with
the order of the numerator polynomial increases by 1, with the coefficient of s0 will be 0.

31 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

Note:
The method adopted in this experiment to obtain the doublet response of the system is by
increasing the order of the numerator polynomial, in the same order of defined coefficients,
and the coefficient of s0 becoming 0, and taking the impulse response of this modified transfer
function. The Procedure in depicted in the following section.

MATLAB Code
Case 1:
Static Margin = 0.2211.
%-Boeing B747-8, Sea Level Condition, M = 0.25, Sm = 0.2211-
% clc
A =1; % Coefficients in the characteristic equation A - E
B = 1.2632;
C = 0.9618;
D = 0.0979;
E = 0.0127;
Atheta = -0.5602; % Coefficients in the Numerator polynomial Au,
Bu, Cu, & Du
Btheta = -0.3867;
Ctheta = -0.0440;
DTF = [A B C D E]; % Denominator polynomial of the TF
Ntheta = [Atheta Btheta Ctheta 0]; % TF Numerator
polynomial TFtheta = tf([Ntheta],[DTF]) % TF:
[theta(s)/de(s)] figure
impulseplot(TFtheta) % Impulse response of the
TF grid on

Case 2:
Static Margin = 0.6633.
%----- Boeing B747-8, Sea Level Condition, M = 0.25, Sm = 0.6633-
------%
clc
A =1; % Coefficients in the characteristic equation A - E
B = 1.2632;
C = 2.0199;
D = 0.2121;
E = 0.0310;
Atheta = -0.5602; % Coefficients in the Numerator polynomial Au,
Bu, Cu, & Du
Btheta = -0.3498;
Ctheta = -0.0395;
DTF = [A B C D E]; % Denominator polynomial of the TF
Ntheta = [Atheta Btheta Ctheta 0]; % Numerator polynomial of the
TF
TFtheta = tf([Ntheta],[DTF]) % TF: forward velocity vs elevator
deflection [theta(s)/de(s)]
figure
impulseplot(TFtheta) % Impulse response of the TF
[theta(s)/de(s)]

32 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

grid on

Viva Questions:

1. Which all are the force equations?

2. Which all are the kinematic equations?
3. What do you understand by Euler angles?
4. What do understand by time response of a system?
5. How do you differentiate Shor period and Phugoid motion?
6. What do you understand by ‗Small-Disturbance Theory„?
7. What are stability parameters?
8. What do understand by order and type of a system?
9. How is the time response of a first order and higher order system differ for a doublet input?
10. What do you understand by Static Margin and how it influences the longitudinal stability of
an aircraft?

33 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-13

GIVEN A QUADRATIC EQUATION, DETERMINE TWO

QUADRATICS THAT SHALL RESULT IN SHORT PERIOD
OSCILLATION AND POLES OF PHUGOID. VARY THE
COEFFICIENTS OF POLYNOMIAL TO STUDY THE MOVEMENT OF
POLES.

AIM:
To investigate the poles of short period oscillation and phugoid mode.
THEORY:
The fourth order characteristic equation for longitudinal motion can be written as the product of
two second order polynomials as:

The fourth order characteristic equation can be convert into Bi-quadratic equations where one will
define as short period and other as phugoid period mode. The coefficients of each characteristic
equation change with flight condition, airplane mass, mass distribution airplane geometry and
aerodynamics characteristics. These changes translate to zeta and Wn
1. Short period Mode:

Characterized by complex conjugate roots with moderate to relatively high damping ratio, high
natural and damping frequency. The control action can be easily demonstrated by first trimming
the aircraft then disturbing it from trim with a forward aft natural stick input (Commonly called
doublet).
The resulting response back to trim maybe first order (exponential delay) or second order
(oscillatory) significant vibrations to angle of attack ( α) and pitch attitude ( ) longitudinal motion
occurs while the airspeed remains fairly constant.
Trim is generally regained in few seconds, thus the descriptive name short period and the small
vibration in airspeed.

2. Phugoid Mode:
Characterized by complex conjugate roots, with a relatively low damping ratio, natural and
damping frequency.
Control action can be demonstrated by trimming the aircraft in level flight, then input aftstick for
2-3 seconds bleed off some airspped and then returning the stick to the natural (trimmed) position.
The resulting response is usually oscillatory with significant vibrations in pitch attitude and
airspeed, while angle of attack remains relatively constant. As the oscillation starts the airspeed
decreases while the airplane gain altitude (Pitch angle is positive). The aircraft then begins to lose
altitude and airspeed increases while pitch angle decreases. The period of phugoid is typically quite
long (30-120s).

Characteristics equation

34 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

MATLAB CODE

CASE 1:
clc;
x=[1 4.5898 21.653584 0.251526 0.18792];
sys_g=tf([-0.08924 -29.9945 -0.31579805],x);
figure(1)
h=pzplot(sys_g,'r');
p=getoptions(h);
p.Title.Color=[0,0,1];
setoptions(h,p);
y=roots(x);
z1=[1 -(y(1)+y(2)) y(1)*y(2)];
z2=[1 -(y(3)+y(4)) y(3)*y(4)];
wn1=(y(1)*y(2))^0.5;
wn2=(y(3)*y(4))^0.5;
if (wn1>wn2)
sp = poly2sym(sym(z2))
wn_sp=wn2
eta_sp=(z2(2)/(2*wn_sp))
ph = poly2sym(sym(z1))
wn_ph=wn1
eta_ph=(z1(2)/(2*wn_ph))
figure(2)
sys1=tf([1],z2);
sys2=tf([1],z1);
h1=pzplot(sys1,'r',sys2,'b');
q=getoptions(h1);
q.Title.Color=[1,0,0];
setoptions(h1,p);
else
sp = poly2sym(sym(z1))
wn_sp=wn1
eta_sp=(z1(2)/(2*wn_sp))
ph = poly2sym(sym(z2))
wn_ph=wn2
eta_ph=(z2(2)/(2*wn_ph))
figure(2)
sys1=tf([1],z2);
sys2=tf([1],z1);
h1=pzplot(sys1,'b',sys2,'r');
q=getoptions(h1);
q.Title.Color=[1,0,0];
setoptions(h1,p);
end

CASE 2:
Substitute the value of x in case 1:
x=[1 8.2 21.653584 0.251526 0.18792];

35 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

VIVA QUESTIONS:

1. Explain Short Period Mode?

2. Explain Phugoid period?
3. What is meant by damping ratio?
4. Explain the difference between natural and damped frequency?
5. Explain Longitudinal stability?
6. Explain the Parameters affect the Short Period Mode?
7. Explain the Parameters affect the Phugoid Mode?

36 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

EXPERIMENT NO-14

GIVEN A QUARTIC CHARACTERISTICS EQUATION, DETERMINE

POLES AND TIME CONSTANT FOR ROLL MODE, SPIRAL MODE
AND DUTCH ROLL. VARY THE COEFFICIENTS OF POLYNOMIAL
TO STUDY THE MOVEMENTS OF POLES

AIM:

To determine poles and time constant for Roll mode, Spiral Mode and Dutch roll from a given
quartic equation and study the movement of poles.

THEORY:
Three-degree-of-freedom analysis of the lateral-directional modes of motion

The preceding development of transfer functions for the three lateral-directional motion variables
, β, ϕ and ψ leads to a 3-DOF solution for lateral-directional motion. Normally, the fourth-order
characteristic equation for lateral-directional motion is written as the product of one second-order
(oscillatory) and two first-order (non-oscillatory) polynomials.

Dutch Roll, Spiral and Roll three lateral-directional dynamic modes. Each of these polynominals
can be thought of as a separate characteristic equation that defines the dynamic characteristics of
its respective mode.

As with the longitudinal case, the coefficients (and roots) of each characteristic equation change
with flight condition, airplane mass, mass distribution, airplane geometry, and aerodynamic
characteristics. These changes translate to changes in but the fundamental presence of the Dutch
roll, roll, and spiral modes is maintained.

The dutch roll mode is a second-order response (complex conjugate roots) usually characterized
by concurrent oscillations in the three lateral-directional motion variables β, ϕ and ψ.

The roll mode has a real root and a first-order (non-oscillatory) response that involves almost a
pure rolling motion about the x stability axis. It is usually stable at low and moderate angles of
attack but may be unstable at high angles of attack. The roll mode can be excited by a disturbance
or an aileron input.

The spiral mode is a first-order response (real root) that involves a relatively slow roll and yawing
motion of the aircraft. It may be stable or unstable. The spiral is usually initiated by a displacement
in roll angle and appears as a descending turn with increasing roll angle if unstable. If the spiral is
stable, the aircraft simply returns to wings level after a roll angle displacement. The primary motion
variables during the spiral are ϕ and ψ, while β remains close to zero. A high
degree of lateral stability will tend to make the spiral stable.

The fourth order characteristic equation can be defined as:

37 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

DR- Dutch Roll

r- Roll Mode
s- Spiral Mode

MATLAB Code:

clc;
close all;
x=[1 102 202 201 100]; %characteristic
equation sys_g=tf([1],x);% Transfer Function
figure(1)
h=pzplot(sys_g,'r');% Pole-Zero plot for transfer function
p=getoptions(h);
p.Title.Color=[0,0,1];
setoptions(h,p);
y=roots(x);
z1=abs(y(1));
z2=abs(y(2));
if (z1>z2)
roll_mode= poly2sym(sym([1 z1]))% Equation for Roll mode
ploes_rollmode=y(1)
timeconstant_roll=(1/y(1))
spiral_mode= poly2sym(sym([1 z2]))% Equation for Spiral mode
ploes_spiralmode=y(2)
timeconstant_spiral=(1/y(2))
else
roll_mode= poly2sym(sym([1 z2]))% Equation for Roll mode
ploes_rollmode=y(2)
timeconstant_roll=(1/y(2))
spiral_mode= poly2sym(sym([1 z1]))% Equation for Spiral mode
ploes_spiralmode=y(1)
timeconstant_spiral=(1/y(1))
end
z3=[1 -(y(3)+y(4)) y(3)*y(4)];
dutchroll_mode= poly2sym(sym(z3))% Equation for Dutchroll mode
poles_dutchroll=y(3)
poles_dutchroll2=y(4)
y5=abs((y(3)+y(4))/2);
timeconstant_dutchroll=(1/y5)

VIVA QUESTIONS:
1. Explain Roll Mode?
2. Explain Spiral Mode?
3. Explain Dutch Roll Mode?
4. Derive the damping ratio and natural frequency for the Dutch roll?

38 Department of AE, MITE.

FLIGHT SIMULATION LAB – BAEL606 AERONAUTICAL ENGINEERING 2024-25

5. Explain the time constant for first and second order system?
6. Explain Lateral-Directional Motion?
7. Explain the factors affecting lateral-directional Motion?

39 Department of AE, MITE.