This is the first example o PyFME. The main purpose of this example is to check if the aircraft trimmed in a given state maintains the trimmed flight condition.
The aircraft used is a Cessna 310, ISA1976 integrated with Flat Earth (euler angles).
Example with trimmed aircraft: stationary, horizontal, symmetric, wings level flight.
In [100]:
# -*- coding: utf-8 -*-
Import python libreries needed.
In [101]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
Import PyFME classes
In [102]:
from pyfme.aircrafts import Cessna310
from pyfme.environment.environment import Environment
from pyfme.environment.atmosphere import ISA1976
from pyfme.environment.gravity import VerticalConstant
from pyfme.environment.wind import NoWind
from pyfme.models.systems import EulerFlatEarth
from pyfme.simulator import BatchSimulation
from pyfme.utils.trimmer import steady_state_flight_trimmer
Initialize variables
In [103]:
aircraft = Cessna310()
atmosphere = ISA1976()
gravity = VerticalConstant()
wind = NoWind()
environment = Environment(atmosphere, gravity, wind)
Initial conditions
In [104]:
TAS = 312.5 * 0.3048 # m/s
h0 = 8000 * 0.3048 # m
psi0 = 1.0 # rad
x0, y0 = 0, 0 # m
turn_rate = 0.0 # rad/s
gamma0 = 0.0 # rad
In [105]:
system = EulerFlatEarth(lat=0, lon=0, h=h0, psi=psi0, x_earth=x0, y_earth=y0)
not_trimmed_controls = {'delta_elevator': 0.05,
'hor_tail_incidence': 0.00,
'delta_aileron': 0.01 * np.sign(turn_rate),
'delta_rudder': 0.01 * np.sign(turn_rate),
'delta_t': 0.5}
controls2trim = ['delta_elevator', 'delta_aileron', 'delta_rudder', 'delta_t']
In [106]:
trimmed_ac, trimmed_sys, trimmed_env, results = steady_state_flight_trimmer(
aircraft, system, environment, TAS=TAS, controls_0=not_trimmed_controls,
controls2trim=controls2trim, gamma=gamma0, turn_rate=turn_rate, verbose=2)
In [107]:
print(results)
In [108]:
my_simulation = BatchSimulation(trimmed_ac, trimmed_sys, trimmed_env)
In [109]:
tfin = 150 # seconds
N = tfin * 100 + 1
time = np.linspace(0, tfin, N)
initial_controls = trimmed_ac.controls
In [110]:
controls = {}
for control_name, control_value in initial_controls.items():
controls[control_name] = np.ones_like(time) * control_value
my_simulation.set_controls(time, controls)
In [111]:
par_list = ['x_earth', 'y_earth', 'height',
'psi', 'theta', 'phi',
'u', 'v', 'w',
'v_north', 'v_east', 'v_down',
'p', 'q', 'r',
'alpha', 'beta', 'TAS',
'F_xb', 'F_yb', 'F_zb',
'M_xb', 'M_yb', 'M_zb']
In [112]:
my_simulation.set_par_dict(par_list)
my_simulation.run_simulation()
In [113]:
plt.style.use('ggplot')
In [114]:
for ii in range(len(par_list) // 3):
three_params = par_list[3*ii:3*ii+3]
fig, ax = plt.subplots(3, 1, sharex=True)
for jj, par in enumerate(three_params):
ax[jj].plot(time, my_simulation.par_dict[par])
ax[jj].set_ylabel(par)
ax[jj].set_xlabel('time (s)')
fig.tight_layout()
fig = plt.figure()
ax = Axes3D(fig)
ax.plot(my_simulation.par_dict['x_earth'],
my_simulation.par_dict['y_earth'],
my_simulation.par_dict['height'])
ax.plot(my_simulation.par_dict['x_earth'],
my_simulation.par_dict['y_earth'],
my_simulation.par_dict['height'] * 0)
ax.set_xlabel('x_earth')
ax.set_ylabel('y_earth')
ax.set_zlabel('z_earth')
plt.show()