# -*- coding: utf-8 -*-
"""
pyBADA
Generic BADA4 aircraft performance module
Developed @EUROCONTROL (EIH)
2024
"""
from math import sqrt, pow, pi, isnan, sin, asin
import numpy as np
import os
from datetime import date
import xml.etree.ElementTree as ET
from scipy.optimize import fmin, fminbound
import pandas as pd
from pyBADA import constants as const
from pyBADA import conversions as conv
from pyBADA import atmosphere as atm
from pyBADA import configuration as configuration
from pyBADA.aircraft import Airplane, BadaFamily, Bada
[docs]
def proper_round(num, dec=0):
num = str(num)[: str(num).index(".") + dec + 2]
if num[-1] >= "5":
return float(num[: -2 - (not dec)] + str(int(num[-2 - (not dec)]) + 1))
return float(num[:-1])
[docs]
def checkArgument(argument, **kwargs):
if kwargs.get(argument) is not None:
return kwargs.get(argument)
else:
raise TypeError("Missing " + argument + " argument")
[docs]
class Parser:
"""This class implements the BADA4 parsing mechanism to parse xml and GPF(xml) BADA4 files."""
def __init__(self):
pass
[docs]
@staticmethod
def readMappingFile(filePath):
"""
Parses the BADA4 mapping XML file and stores a dictionary of aircraft code names and their corresponding XML file paths.
This function processes the BADA4 aircraft model mapping XML file to create a dictionary that maps the aircraft code names
to the XML file paths for their corresponding BADA models. The mapping file contains information about available models
for the specified BADA version.
:param filePath: The path to the directory containing the BADA4 mapping XML file.
:type filePath: str
:raises IOError: If the XML file cannot be found or parsed.
:return: A dictionary with aircraft code names as keys and corresponding file names as values.
:rtype: dict
"""
filename = os.path.join(filePath, "aircraft_model_default.xml")
code_fileName = {}
if os.path.isfile(filename):
try:
tree = ET.parse(filename)
root = tree.getroot()
except:
raise IOError(filename + " not found or in correct format")
for child in root.iter("MAP"):
code = child.find("code").text
file = child.find("file").text
code_fileName[code] = file
return code_fileName
[docs]
@staticmethod
def parseMappingFile(filePath, acName):
"""
Retrieves the file name for a given aircraft name from the parsed BADA4 mapping file.
This function uses the readMappingFile method to parse the BADA4 XML mapping file and returns the file name
associated with the given aircraft name (acName).
:param filePath: The path to the directory containing the BADA4 mapping XML file.
:param acName: The aircraft code name for which the corresponding file is being requested.
:type filePath: str
:type acName: str
:return: The file name corresponding to the aircraft code, or None if the aircraft code is not found.
:rtype: str or None
"""
code_fileName = Parser.readMappingFile(filePath)
if acName in code_fileName:
fileName = code_fileName[acName]
return fileName
else:
return None
[docs]
@staticmethod
def parseXML(filePath, acName):
"""
Parses the BADA4 XML file for a specific aircraft model and extracts various parameters.
This function parses the BADA4 aircraft XML file for a given aircraft model (acName). It retrieves
general information about the aircraft, engine type, aerodynamic configurations, performance parameters, and more.
:param filePath: The path to the folder containing the BADA4 XML file.
:param acName: The aircraft code name for which the XML file is being parsed.
:type filePath: str
:type acName: str
:raises IOError: If the XML file cannot be found or parsed.
:return: A pandas DataFrame containing the parsed data for the specified aircraft.
:rtype: pd.DataFrame
"""
acXmlFile = os.path.join(filePath, acName, acName) + ".xml"
try:
tree = ET.parse(acXmlFile)
root = tree.getroot()
except:
raise IOError(acXmlFile + " not found or in correct format")
# Parse general aircraft data
model = root.find("model").text # aircraft model
engineType = root.find("type").text # aircraft type
engines = root.find("engine").text # engine type
ICAO_desig = {} # ICAO designator and WTC
ICAO = root.find("ICAO").find("designator").text
WTC = root.find("ICAO").find("WTC").text
# Parse engine data
PFM = root.find("PFM") # get PFM
MREF = float(PFM.find("MREF").text) # reference mass
WREF = MREF * const.g
LHV = float(PFM.find("LHV").text)
n_eng = int(PFM.find("n_eng").text) # number of engines
rho = float(PFM.find("rho").text)
TFA = None
if PFM.find("TFA") is not None:
TFA = float(PFM.find("TFA").text)
# parameters introduced with BADA 4.3
p_delta = None
if PFM.find("p_delta") is not None:
p_delta = float(PFM.find("p_delta").text)
p_theta = None
if PFM.find("p_theta") is not None:
p_theta = float(PFM.find("p_theta").text)
TFM = PFM.find("TFM") # get TFM
# set all the parameters to NONE as a default
a = None
f = None
b = None
c = None
ti = None
fi = None
throttle = None
prop_dia = None
max_eff = None
p = None
Hd_turbo = None
CPSFC = None
P = None
kink = {}
max_power = {}
if engineType == "JET":
CT = TFM.find("CT") # Thrust coefficients
a = []
for i in CT.findall("a"):
a.append(float(i.text)) # C_T polynomial coefficients
CF = TFM.find("CF") # Fuel flow coefficients
f = []
for i in CF.findall("f"):
f.append(float(i.text)) # FF polynomial coefficients
b = {}
c = {}
for rating in ["MCMB", "MCRZ", "MTKF"]:
ENG = TFM.find(rating)
if ENG is not None:
flat_rating = ENG.find("flat_rating")
temp_rating = ENG.find("temp_rating")
kink[rating] = float(
ENG.find("kink").text
) # kink point for Max Climb
b[rating] = []
for i in flat_rating.findall("b"):
b[rating].append(float(i.text))
c[rating] = []
for i in temp_rating.findall("c"):
c[rating].append(float(i.text))
# Idle data
LIDL = TFM.find("LIDL")
CT = LIDL.find("CT")
ti = []
for i in CT.findall("ti"):
ti.append(float(i.text)) # idle thrust coefficients
CF = LIDL.find("CF")
fi = []
for i in CF.findall("fi"):
fi.append(float(i.text)) # idle fuel flow coefficients
throttle = {}
throttle["low"] = float(TFM.find("throttle").find("low").text)
throttle["high"] = float(TFM.find("throttle").find("high").text)
elif engineType == "TURBOPROP":
TPM = PFM.find("TPM") # get TFM
prop_dia = float(TPM.find("prop_dia").text)
max_eff = float(TPM.find("max_eff").text)
p = {}
CP = TPM.find("CP") # Thrust coefficients
a = []
for i in CP.findall("a"):
a.append(float(i.text)) # C_P polynomial coefficients
CF = TPM.find("CF") # Fuel flow coefficients
f = []
for i in CF.findall("f"):
f.append(float(i.text)) # FF polynomial coefficients
for rating in ["MCMB", "MCRZ"]:
ENG = TPM.find(rating)
if ENG is not None:
r = ENG.find("rating")
max_power[rating] = float(ENG.find("max_power").text)
p[rating] = []
for i in r.findall("p"):
p[rating].append(float(i.text))
# Idle data
LIDL = TPM.find("LIDL")
CT = LIDL.find("CT")
ti = []
for i in CT.findall("ti"):
ti.append(float(i.text)) # idle thrust coefficients
CF = LIDL.find("CF")
fi = []
for i in CF.findall("fi"):
fi.append(float(i.text)) # idle fuel flow coefficients
throttle = {}
throttle["low"] = float(TPM.find("throttle").find("low").text)
throttle["high"] = float(TPM.find("throttle").find("high").text)
elif engineType == "PISTON":
PEM = PFM.find("PEM")
prop_dia = float(PEM.find("prop_dia").text)
max_eff = float(PEM.find("max_eff").text)
Hd_turbo = float(PEM.find("Hd_turbo").text)
CPSFC = float(PEM.find("CPSFC").text)
P = float(PEM.find("P").text)
# Parse aerodynamic data
AFCM = root.find("AFCM") # get AFCM
S = float(AFCM.find("S").text)
HLPosition = {}
configName = {}
VFE = {}
d = {}
CL_max = {}
bf = []
HLids = []
Mmin = None
Mmax = None
CL_Mach0 = None
CL_clean = None
for conf in AFCM.findall("Configuration"):
HLid = float(conf.get("HLid"))
HLids.append(str(HLid))
HLPosition[HLid] = float(conf.find("HLPosition").text)
configName[HLid] = conf.find("name").text
VFE[HLid] = float(conf.find("vfe").text)
d[HLid] = {}
LGUP = conf.find("LGUP")
LGDN = conf.find("LGDN")
if LGUP is not None:
DPM = LGUP.find("DPM_nonclean")
if DPM is not None: # DPM is not clean
CD = DPM.find("CD_nonclean")
if CD is not None:
if CD.find("d") is not None:
d[HLid]["LGUP"] = []
for i in CD.findall("d"):
d[HLid]["LGUP"].append(float(i.text))
else: # DPM is clean
DPM = LGUP.find("DPM_clean")
if DPM.find("M_max") is not None:
M_max = float(DPM.find("M_max").text)
if DPM.find("scalar") is not None:
scalar = float(DPM.find("scalar").text)
CD = DPM.find("CD_clean")
if CD is not None:
if CD.find("d") is not None:
d[HLid]["LGUP"] = []
for i in CD.findall("d"):
d[HLid]["LGUP"].append(float(i.text))
BLM = LGUP.find("BLM")
if BLM is not None: # BLM is not clean
if BLM.find("CL_max") is not None:
if HLid not in CL_max:
CL_max[HLid] = {}
CL_max[HLid]["LGUP"] = float(BLM.find("CL_max").text)
else: # BLM is clean
BLM = LGUP.find("BLM_clean")
if BLM.find("Mmin") is not None:
Mmin = float(BLM.find("Mmin").text)
if BLM.find("Mmax") is not None:
Mmax = float(BLM.find("Mmax").text)
if BLM.find("CL_Mach0") is not None:
CL_Mach0 = float(BLM.find("CL_Mach0").text)
CL_clean = BLM.find("CL_clean")
CL_clean = BLM.find("CL_clean")
if CL_clean is not None:
for i in CL_clean.findall("bf"):
bf.append(float(i.text))
if (
LGDN is not None
): # Landing gear NOT allowed in clean configuration
DPM = LGDN.find("DPM_nonclean")
if DPM is not None: # DPM is not clean
CD = DPM.find("CD_nonclean")
if CD.find("d") is not None:
d[HLid]["LGDN"] = []
for i in CD.findall("d"):
d[HLid]["LGDN"].append(float(i.text))
BLM = LGDN.find("BLM")
if BLM is not None: # BLM is not clean
if HLid not in CL_max:
CL_max[HLid] = {}
CL_max[HLid]["LGDN"] = float(BLM.find("CL_max").text)
ALM = root.find("ALM") # get ALM
DLM = ALM.find("DLM") # get DLM
MTOW = float(DLM.find("MTOW").text)
OEW = float(DLM.find("OEW").text)
MFL = float(DLM.find("MFL").text)
MTW = None
MZFW = None
MPL = None
MLW = None
if DLM.find("MTW") is not None:
MTW = float(DLM.find("MTW").text)
if DLM.find("MZFW") is not None:
MZFW = float(DLM.find("MZFW").text)
if DLM.find("MPL") is not None:
MPL = float(DLM.find("MPL").text)
if DLM.find("MLW") is not None:
MLW = float(DLM.find("MLW").text)
GLM = ALM.find("GLM") # get GLM
hmo = float(GLM.find("hmo").text)
if GLM.find("mfa") is not None:
mfa = float(GLM.find("mfa").text)
else:
mfa = None
KLM = ALM.find("KLM") # get GLM
MMO = None
MLE = None
VLE = None
VMO = None
if KLM.find("mmo") is not None:
MMO = float(KLM.find("mmo").text)
if KLM.find("mle") is not None:
MLE = float(KLM.find("mle").text)
if KLM.find("vmo") is not None:
VMO = float(KLM.find("vmo").text)
if KLM.find("vle") is not None:
VLE = float(KLM.find("vle").text)
ground = root.find("Ground")
dimensions = ground.find("Dimensions")
span = float(dimensions.find("span").text)
length = float(dimensions.find("length").text)
# ARPM model
ARPM = root.find("ARPM")
aeroConfSchedule = ARPM.find("AeroConfSchedule")
# all aerodynamic configurations
aeroConfig = {}
for conf in aeroConfSchedule.findall("AeroPhase"):
name = conf.find("name").text
HLid = float(conf.find("HLid").text)
LG = "LG" + conf.find("LG").text
aeroConfig[name] = {"HLid": HLid, "LG": LG}
speedScheduleList = ARPM.find("SpeedScheduleList")
SpeedSchedule = speedScheduleList.find("SpeedSchedule")
# all phases of flight
speedSchedule = {}
for phaseOfFlight in SpeedSchedule.findall("SpeedPhase"):
name = phaseOfFlight.find("name").text
CAS1 = float(phaseOfFlight.find("CAS1").text)
CAS2 = float(phaseOfFlight.find("CAS2").text)
M = float(phaseOfFlight.find("M").text)
speedSchedule[name] = {"CAS1": CAS1, "CAS2": CAS2, "M": M}
# Single row dataframe
data = {
"acName": [acName],
"model": [model],
"engineType": [engineType],
"engines": [engines],
"ICAO": [ICAO],
"WTC": [WTC],
"MREF": [MREF],
"WREF": [WREF],
"LHV": [LHV],
"n_eng": [n_eng],
"rho": [rho],
"TFA": [TFA],
"p_delta": [p_delta],
"p_theta": [p_theta],
"a": [a],
"f": [f],
"b": [b],
"c": [c],
"ti": [ti],
"fi": [fi],
"throttle": [throttle],
"prop_dia": [prop_dia],
"max_eff": [max_eff],
"p": [p],
"Hd_turbo": [Hd_turbo],
"CPSFC": [CPSFC],
"P": [P],
"kink": [kink],
"max_power": [max_power],
"S": [S],
"HLPosition": [HLPosition],
"configName": [configName],
"VFE": [VFE],
"d": [d],
"CL_max": [CL_max],
"bf": [bf],
"HLids": [HLids],
"Mmin": [Mmin],
"Mmax": [Mmax],
"CL_Mach0": [CL_Mach0],
"CL_clean": [CL_clean],
"M_max": [M_max],
"scalar": [scalar],
"MTOW": [MTOW],
"OEW": [OEW],
"MFL": [MFL],
"MTW": [MTW],
"MZFW": [MZFW],
"MPL": [MPL],
"MLW": [MLW],
"hmo": [hmo],
"mfa": [mfa],
"MMO": [MMO],
"MLE": [MLE],
"VLE": [VLE],
"VMO": [VMO],
"span": [span],
"length": [length],
"aeroConfig": [aeroConfig],
"speedSchedule": [speedSchedule],
}
df_single = pd.DataFrame(data)
return df_single
[docs]
@staticmethod
def parseGPF(filePath):
"""
Parses the BADA4 GPF XML file and extracts key performance factors.
This function processes the BADA4 GPF XML file to extract data related to various flight performance factors, such as minimum
climb/descent speeds, maximum altitude limits for different phases of flight, and speed schedules for climb and descent.
:param filePath: The path to the directory containing the BADA4 GPF XML file.
:type filePath: str
:raises IOError: If the GPF XML file cannot be found or parsed.
:return: A pandas DataFrame containing the parsed GPF data.
:rtype: pd.DataFrame
"""
filename = os.path.join(filePath, "GPF.xml")
tree = ET.parse(filename)
try:
tree = ET.parse(filename)
root = tree.getroot()
except:
raise IOError(filename + " not found or in correct format")
CVminTO = float(root.find("CVminTO").text) # CVminTO
CVmin = float(root.find("CVmin").text) # CVmin
HmaxList = root.find("HmaxList")
HmaxPhase = {} # phase of flight
for Hmax_Phase in HmaxList.findall("HmaxPhase"):
phaseName = Hmax_Phase.find("Phase").text
Hmax = float(Hmax_Phase.find("Hmax").text)
HmaxPhase[phaseName] = Hmax
V_des = {} # V_des
V_cl = {} # V_cl
VdList = root.find("VdList")
for VdPhase in VdList.findall("VdPhase"):
Phase = VdPhase.find("Phase")
name = Phase.find("name").text
index = int(Phase.find("index").text)
Vd = float(VdPhase.find("Vd").text)
if name == "CL":
V_cl[index] = Vd
elif name == "DES":
V_des[index] = Vd
# Single row dataframe
data = {
"CVminTO": [CVminTO],
"CVmin": [CVmin],
"HmaxPhase": [HmaxPhase],
"V_des": [V_des],
"V_cl": [V_cl],
}
df_single = pd.DataFrame(data)
return df_single
[docs]
@staticmethod
def combineXML_GPF(XMLDataFrame, GPFDataframe):
"""
Combines the parsed aircraft XML DataFrame with the parsed GPF DataFrame.
This function merges two DataFrames, one containing the parsed aircraft-specific data (from the XML file)
and the other containing the parsed GPF data. This combination provides a unified set of aircraft performance
data along with general performance factors.
:param XMLDataFrame: A DataFrame containing the parsed aircraft XML data.
:param GPFDataframe: A DataFrame containing the parsed GPF data.
:type XMLDataFrame: pd.DataFrame
:type GPFDataframe: pd.DataFrame
:return: A combined DataFrame with both aircraft and GPF data.
:rtype: pd.DataFrame
"""
# Combine data with GPF data (temporary solution)
combined_df = pd.concat(
[
XMLDataFrame.reset_index(drop=True),
GPFDataframe.reset_index(drop=True),
],
axis=1,
)
return combined_df
[docs]
@staticmethod
def parseAll(badaVersion, filePath=None):
"""
Parses all BADA4 XML files and combines the data into a single DataFrame.
This function parses both the BADA4 aircraft XML files and the GPF (General Performance Factors) file.
It combines the data from both sources and creates a unified DataFrame that contains all aircraft and
performance-related data for the specified BADA version.
:param badaVersion: The version of BADA (e.g., "4.3") to be parsed.
:param filePath: The path to the folder containing the BADA4 XML files and GPF file.
If None, the default path from the configuration will be used.
:type badaVersion: str
:type filePath: str, optional
:raises IOError: If any of the necessary XML files cannot be found or parsed.
:return: A DataFrame containing the combined data from all parsed aircraft and GPF files.
:rtype: pd.DataFrame
"""
if filePath == None:
filePath = configuration.getBadaVersionPath(
badaFamily="BADA4", badaVersion=badaVersion
)
else:
filePath = filePath
# parsing GPF file
GPFparsedDataframe = Parser.parseGPF(filePath)
# retrieving mapping data
code_fileName = Parser.readMappingFile(filePath)
# get names of all the folders in the main BADA model folder to search for XML files
subfolders = configuration.list_subfolders(filePath)
# Initialize an empty list to collect DataFrames
mapping_dfs = []
if code_fileName:
for code in code_fileName:
file = code_fileName[code]
if file in subfolders:
# Parse the original XML of a model
df = Parser.parseXML(filePath, file)
# Rename 'acName' in the DataFrame to match the code model name
df.at[0, "acName"] = code
# Combine data with GPF data (temporary solution)
combined_df = Parser.combineXML_GPF(df, GPFparsedDataframe)
# Drop columns that are all NaN
combined_df = combined_df.dropna(axis=1, how="all")
# Check if combined_df is not empty
if not combined_df.empty:
mapping_dfs.append(combined_df)
# Concatenate all collected DataFrames
if mapping_dfs:
merged_mapping_df = pd.concat(mapping_dfs, ignore_index=True)
else:
merged_mapping_df = pd.DataFrame()
# Initialize an empty list to collect DataFrames
original_dfs = []
for file in subfolders:
# Parse the original XML of a model
df = Parser.parseXML(filePath, file)
# Combine data with GPF data (temporary solution)
combined_df = pd.concat(
[
df.reset_index(drop=True),
GPFparsedDataframe.reset_index(drop=True),
],
axis=1,
)
# Drop columns that are all NaN
combined_df = combined_df.dropna(axis=1, how="all")
# Check if combined_df is not empty
if not combined_df.empty:
original_dfs.append(combined_df)
# Concatenate all collected DataFrames
if original_dfs:
merged_original_df = pd.concat(original_dfs, ignore_index=True)
else:
merged_original_df = pd.DataFrame()
# Merge mapping and original aircraft models
merged_final_df = pd.concat(
[merged_original_df, merged_mapping_df], ignore_index=True
)
return merged_final_df
[docs]
class BADA4(Airplane, Bada):
"""This class implements the part of BADA4 performance model that will be used in other classes following the BADA4 manual.
:param AC: Aircraft object {BADA4}.
:type AC: bada4Aircraft.
"""
def __init__(self, AC):
"""
Initializes the BADA4 class by inheriting from the parent Airplane class
and assigns the parsed aircraft data to the class instance.
:param AC: Aircraft object {BADA4}.
:type AC: bada4Aircraft.
"""
super().__init__()
self.AC = AC
[docs]
def CL(self, delta, mass, M, nz=1.0):
"""
Computes the aircraft's lift coefficient based on the current Mach number (M),
normalized air pressure (delta), aircraft mass, and load factor (nz).
:param M: Mach number [-].
:param delta: Normalized air pressure [-].
:param mass: Aircraft mass [kg].
:param nz: Load factor [-]. Default is 1.0 (straight and level flight).
:type M: float
:type delta: float
:type mass: float
:type nz: float
:return: Lift coefficient (CL) [-].
:rtype: float
"""
return (
2
* mass
* const.g
* nz
/ (delta * const.p_0 * const.Agamma * pow(M, 2) * self.AC.S)
)
[docs]
def CLPoly(self, M):
"""
Computes the lift coefficient polynomial for the given Mach number (M).
This method uses a 5th-degree polynomial defined by the coefficients in the parsed aircraft data (self.AC.bf).
:param M: Mach number [-].
:type M: float
:return: Lift coefficient (polynomial approximation) [-].
:rtype: float
"""
CLpoly = 0.0
for i in range(5):
CLpoly += self.AC.bf[i] * pow(M, i)
return CLpoly
[docs]
def CLmax(self, M, HLid, LG):
"""
Computes the maximum lift coefficient (CLmax) for the given Mach number (M),
high-lift device (HLid) position, and landing gear (LG) configuration.
If the aircraft is in a clean configuration (HLid == 0 and LG == "LGUP"),
the method interpolates or extrapolates the lift coefficient based on Mach number.
:param M: Mach number [-].
:param HLid: High-lift device position [-].
:param LG: Landing gear position [LGUP/LGDN] [-].
:type M: float
:type HLid: float
:type LG: str
:return: Maximum lift coefficient (CLmax) [-].
:rtype: float
"""
CLmax = 0.0
# CLmax available - non clean configuration
if HLid in self.AC.CL_max and LG in self.AC.CL_max[HLid]:
CLmax = self.AC.CL_max[HLid][LG]
# CLmax unavailable - clean configuration
elif HLid == 0 and LG == "LGUP":
if self.AC.CL_clean is not None:
if M < self.AC.Mmin:
CLmax = self.AC.CL_Mach0 + (M / self.AC.Mmin) * (
self.AC.CLPoly(self.AC.Mmin) - self.AC.CL_Mach0
)
elif M > self.AC.Mmax:
CLder = (
self.AC.bf[1]
+ 2 * self.AC.bf[2] * self.AC.Mmax
+ 3 * self.AC.bf[3] * self.AC.Mmax * self.AC.Mmax
+ 4
* self.AC.bf[4]
* self.AC.Mmax
* self.AC.Mmax
* self.AC.Mmax
)
CLmax = (
self.CLPoly(self.AC.Mmax)
+ self.CLPoly(M - self.AC.Mmax) * CLder
)
else:
CLmax = self.CLPoly(M)
else:
CLmax = self.AC.CL_max[0]["LGUP"]
return CLmax
[docs]
def CF_idle(self, delta, theta, M):
"""
Computes the fuel flow coefficient at idle throttle for JET and TURBOPROP engines.
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param M: Mach speed [-].
:type delta: float
:type theta: float
:type M: float
:return: Idle fuel flow coefficient [-].
:rtype: float
"""
if self.AC.engineType == "JET":
CF_idle = 0.0
for i in range(0, 3):
for j in range(0, 3):
CF_idle += self.AC.fi[i * 3 + j] * (delta**j) * (M**i)
if self.AC.BADAVersion == "4.2":
CF_idle = CF_idle * pow(delta, -1) * pow(theta, -0.5)
elif (
self.AC.BADAVersion == "4.3" or self.AC.BADAVersion == "DUMMY"
):
CF_idle = CF_idle * pow(delta, -1)
elif self.AC.engineType == "TURBOPROP":
CF_idle = 0.0
for i in range(0, 3):
for j in range(0, 3):
CF_idle += self.AC.fi[i * 3 + j] * (delta**j) * (M**i)
CF_idle += (
self.AC.fi[9] * theta
+ self.AC.fi[10] * (theta**2)
+ self.AC.fi[11] * M * theta
+ self.AC.fi[12] * M * delta * sqrt(theta)
+ self.AC.fi[13] * M * delta * theta
)
if self.AC.BADAVersion == "4.2":
CF_idle = CF_idle * pow(delta, -1) * pow(theta, -0.5)
elif (
self.AC.BADAVersion == "4.3" or self.AC.BADAVersion == "DUMMY"
):
CF_idle = CF_idle * pow(delta, -1)
return CF_idle
[docs]
def CF(self, delta, theta, DeltaTemp, **kwargs):
"""
Computes the fuel flow coefficient (CF) for JET, TURBOPROP, and PISTON engines.
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param DeltaTemp: Temperature deviation from ISA [K].
:param kwargs: Optional parameters including 'rating', 'deltaT', 'CT', or 'M'.
:type delta: float
:type theta: float
:type DeltaTemp: float
:return: Fuel flow coefficient [-].
:rtype: float
"""
if self.AC.engineType == "JET":
M = checkArgument("M", **kwargs)
# when idle rating is used
CF_idle = self.CF_idle(delta=delta, theta=theta, M=M)
# for adaptive thrust calcualation if CT is an input
if "CT" in kwargs:
CT = kwargs.get("CT")
CF_gen_rating = 0.0
for i in range(0, 5):
for j in range(0, 5):
CF_gen_rating += (
self.AC.f[i * 5 + j] * (CT**j) * (M**i)
)
CF = max(CF_gen_rating, CF_idle)
# rating as input parameter
elif "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
# in case MCRZ rating is not defined, we switch to MCMB
if rating == "MCRZ" and rating not in self.AC.kink.keys():
rating = "MCMB"
if rating not in self.AC.kink.keys() and rating != "LIDL":
raise ValueError("Unknown engine rating " + rating)
if rating == "LIDL":
CF = CF_idle
elif rating in self.AC.kink.keys():
# when non-idle rating is used
CF_gen_rating = 0.0
CT_rating = self.CT(
rating=rating,
theta=theta,
delta=delta,
M=M,
DeltaTemp=DeltaTemp,
)
for i in range(0, 5):
for j in range(0, 5):
CF_gen_rating += (
self.AC.f[i * 5 + j] * (CT_rating**j) * (M**i)
)
CF = max(CF_gen_rating, CF_idle)
# deltaT - direct throttle as input parameter
elif "deltaT" in kwargs:
# when no rating is used
CF_gen_deltaT = 0.0
CT_deltaT = self.CT(
deltaT=deltaT,
theta=theta,
delta=delta,
M=M,
DeltaTemp=DeltaTemp,
)
for i in range(0, 5):
for j in range(0, 5):
CF_gen_deltaT += (
self.AC.f[i * 5 + j] * (CT_deltaT**j) * (M**i)
)
CF = max(CF_gen_deltaT, CF_idle)
elif self.AC.engineType == "TURBOPROP":
M = checkArgument("M", **kwargs)
# when idle rating is used
CF_idle = self.CF_idle(delta=delta, theta=theta, M=M)
# for adaptive thrust calcualation if CT is an input
if "CT" in kwargs:
CT = kwargs.get("CT")
CP = self.CP(CT=CT, M=M)
CF_gen_rating = 0.0
for i in range(0, 5):
for j in range(0, 5):
CF_gen_rating += (
self.AC.f[i * 5 + j] * (CP**j) * (M**i)
)
CF = max(CF_gen_rating, CF_idle)
# rating as input parameter
elif "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
# in case MCRZ rating is not defined, we switch to MCMB
if rating == "MCRZ" and rating not in self.AC.max_power.keys():
rating = "MCMB"
if rating not in self.AC.max_power.keys() and rating != "LIDL":
raise ValueError("Unknown engine rating " + rating)
if rating == "LIDL":
CF = CF_idle
elif rating in self.AC.max_power.keys():
# when non-idle rating is used
CF_gen_rating = 0.0
CP_rating = self.CP(
rating=rating, theta=theta, delta=delta, M=M
)
for i in range(0, 5):
for j in range(0, 5):
CF_gen_rating += (
self.AC.f[i * 5 + j] * (CP_rating**j) * (M**i)
)
CF = max(CF_gen_rating, CF_idle)
# deltaT - direct throttle as input parameter
elif "deltaT" in kwargs:
# when non-idle rating is used
CF_gen_deltaT = 0.0
CP_deltaT = self.CP(deltaT=deltaT, M=M)
CT_deltaT = self.CT(
deltaT=deltaT,
theta=theta,
delta=delta,
M=M,
DeltaTemp=DeltaTemp,
)
for i in range(0, 5):
for j in range(0, 5):
CF_gen_deltaT += (
self.AC.f[i * 5 + j] * (CP_deltaT**j) * (M**i)
)
CF = max(CF_gen_deltaT, CF_idle)
elif self.AC.engineType == "PISTON":
sigma = atm.sigma(theta=theta, delta=delta)
# for adaptive thrust calcualation if CT is an input
if "CT" in kwargs:
CT = kwargs.get("CT")
M = checkArgument("M", **kwargs)
deltaT_vec = np.arange(0.01, 1.01, 0.01)
CT_diff = []
for k in range(len(deltaT_vec)):
CT_k = self.CT(
theta=theta,
delta=delta,
deltaT=deltaT_vec[k],
M=M,
DeltaTemp=DeltaTemp,
)
CT_diff.append(abs(CT_k - CT))
CT_min_idx = CT_diff.index(min(CT_diff))
deltaT = deltaT_vec[CT_min_idx]
CP = self.CP(
theta=theta,
delta=delta,
sigma=sigma,
deltaT=deltaT,
DeltaTemp=DeltaTemp,
)
# rating as input parameter
elif "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
CP = self.CP(
rating=rating,
theta=theta,
delta=delta,
sigma=sigma,
DeltaTemp=DeltaTemp,
)
CF = self.AC.CPSFC * CP / (delta * sqrt(theta))
return CF
[docs]
def CT(self, delta, **kwargs):
"""
Computes the thrust coefficient (CT) based on engine type, throttle setting,
or engine rating. The thrust coefficient is calculated differently for JET, TURBOPROP,
and PISTON engines based on normalized pressure, temperature, Mach number, and other inputs.
:param delta: Normalized pressure [-].
:type delta: float
:param kwargs: Optional parameters:
- 'rating': Engine rating {MCMB, MCRZ, MTKF, LIDL}.
- 'deltaT': Direct throttle parameter [-].
- 'Thrust': Direct thrust value [N].
- 'M': Mach number [-].
- 'theta': Normalized temperature [-].
- 'DeltaTemp': Deviation from the standard ISA temperature [K].
:return: Thrust coefficient (CT) [-].
:rtype: float
:raises: ValueError: If an invalid rating is provided.
"""
if "Thrust" in kwargs:
Thrust = kwargs.get("Thrust")
CT = Thrust / (delta * self.AC.WREF)
return CT
else:
theta = checkArgument("theta", **kwargs)
DeltaTemp = checkArgument("DeltaTemp", **kwargs)
M = checkArgument("M", **kwargs)
if self.AC.engineType == "JET":
# rating as input parameter
if "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
# in case MCRZ rating is not defined, we switch to MCMB
if rating == "MCRZ" and rating not in self.AC.kink.keys():
rating = "MCMB"
if rating == "MTKF" and rating not in self.AC.kink.keys():
rating = "MCMB"
if rating not in self.AC.kink.keys() and rating != "LIDL":
raise ValueError("Unknown engine rating " + rating)
if rating == "LIDL":
CT = self.CT_LIDL(delta=delta, M=M)
elif rating in self.AC.kink.keys():
CT = self.CT_nonLIDL(
rating=rating,
theta=theta,
delta=delta,
M=M,
DeltaTemp=DeltaTemp,
)
# deltaT - direct throttle as input parameter
elif "deltaT" in kwargs:
deltaT = checkArgument("deltaT", **kwargs)
CT = 0.0
for i in range(0, 6):
for j in range(0, 6):
CT += self.AC.a[i * 6 + j] * (M**j) * (deltaT**i)
# limit CT with CT_LIDL and CT_MCMB
if CT > self.CT_nonLIDL(
rating="MCMB",
theta=theta,
delta=delta,
M=M,
DeltaTemp=DeltaTemp,
):
raise ValueError(
"Throttle parameter value result in CT > CT_MCMB"
+ deltaT
)
elif CT < self.CT_LIDL(delta=delta, M=M):
raise ValueError(
"Throttle parameter value result in CT < CT_LIDL"
+ deltaT
)
elif self.AC.engineType == "TURBOPROP":
# rating as input parameter
if "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
# in case MCRZ rating is not defined, we switch to MCMB
if rating == "MCRZ" and rating not in self.AC.max_power.keys():
rating = "MCMB"
if rating == "MTKF" and rating not in self.AC.kink.keys():
rating = "MCMB"
if rating not in self.AC.max_power.keys() and rating != "LIDL":
raise ValueError("Unknown engine rating " + rating)
if rating == "LIDL":
CT = self.CT_LIDL(theta=theta, delta=delta, M=M)
elif rating in self.AC.max_power.keys():
CT = self.CT_nonLIDL(
rating=rating, theta=theta, delta=delta, M=M
)
# deltaT - direct throttle as input parameter
elif "deltaT" in kwargs:
deltaT = checkArgument("deltaT", **kwargs)
CP = self.CP(deltaT=deltaT, M=M)
CT = CP / M
if CT > self.CT_nonLIDL(
rating="MCMB", theta=theta, delta=delta, M=M
):
raise ValueError(
"Throttle parameter value result in CT > CT_MCMB"
+ deltaT
)
elif CT < self.CT_LIDL(theta=theta, delta=delta, M=M):
raise ValueError(
"Throttle parameter value result in CT < CT_LIDL"
+ deltaT
)
elif self.AC.engineType == "PISTON":
sigma = atm.sigma(theta=theta, delta=delta)
# rating as input parameter
if "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
if rating not in ["MCMB", "MCRZ"] and rating != "LIDL":
raise ValueError("Unknown engine rating " + rating)
CP = self.CP(
rating=rating,
theta=theta,
delta=delta,
sigma=sigma,
DeltaTemp=DeltaTemp,
)
# deltaT - direct throttle as input parameter
elif "deltaT" in kwargs:
deltaT = checkArgument("deltaT", **kwargs)
CP = self.CP(
theta=theta,
delta=delta,
sigma=sigma,
deltaT=deltaT,
DeltaTemp=DeltaTemp,
)
CT = self.CT_nonLIDL(theta=theta, delta=delta, M=M, CP=CP)
return CT
[docs]
def CT_LIDL(self, **kwargs):
"""
Computes the thrust coefficient (CT) for the LIDL (idle) rating.
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param M: Mach number [-].
:param CP: Power coefficient (for PISTON engines) [-].
:type delta: float.
:type theta: float.
:type M: float.
:type CP: float.
:return: Idle thrust coefficient (CT) [-].
:rtype: float
"""
if self.AC.engineType == "JET":
delta = checkArgument("delta", **kwargs)
M = checkArgument("M", **kwargs)
CT = 0.0
for i in range(0, 3):
for j in range(0, 4):
CT += self.AC.ti[i * 4 + j] * pow(delta, j - 1) * (M**i)
elif self.AC.engineType == "TURBOPROP":
theta = checkArgument("theta", **kwargs)
delta = checkArgument("delta", **kwargs)
M = checkArgument("M", **kwargs)
CT = 0.0
for i in range(0, 3):
for j in range(0, 4):
CT += self.AC.ti[i * 4 + j] * pow(delta, j - 1) * (M**i)
CT += (
self.AC.ti[12] * sqrt(theta)
+ self.AC.ti[13] * theta
+ self.AC.ti[14] / sqrt(theta)
+ self.AC.ti[15] * theta**2
)
CT += (
self.AC.ti[16] / delta
+ self.AC.ti[17] * delta
+ self.AC.ti[18] * delta**2
+ self.AC.ti[19] * M
+ self.AC.ti[20] * M**2
) / sqrt(theta)
CT += (
self.AC.ti[21] / M
+ self.AC.ti[22] * delta / M
+ self.AC.ti[23] * M**3
)
CT += (
self.AC.ti[24] * M
+ self.AC.ti[25] * M**2
+ self.AC.ti[26]
+ self.AC.ti[27] * M / delta
) / theta
CT += (
self.AC.ti[28] * M / (delta * theta**2)
+ self.AC.ti[29] * M**2 / (delta * theta**2)
+ self.AC.ti[30] * M**2 / (delta * sqrt(theta))
+ self.AC.ti[31] * delta / theta
)
elif self.AC.engineType == "PISTON":
theta = checkArgument("theta", **kwargs)
delta = checkArgument("delta", **kwargs)
CP = checkArgument("CP", **kwargs)
M = checkArgument("M", **kwargs)
CT = self.CT_nonLIDL(theta=theta, delta=delta, M=M, CP=CP)
return CT
[docs]
def CT_nonLIDL(self, theta, delta, M, **kwargs):
"""
Computes the thrust coefficient (CT) for non-LIDL ratings {MCMB, MCRZ}.
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param M: Mach number [-].
:param CP: Power coefficient (for PISTON engines) [-].
:type delta: float.
:type theta: float.
:type M: float.
:type CP: float.
:return: Thrust coefficient (CT) [-].
:rtype: float
"""
if self.AC.engineType == "JET":
rating = checkArgument("rating", **kwargs)
DeltaTemp = checkArgument("DeltaTemp", **kwargs)
# deltaT below kink point -> flat-rated area
deltaTFlat = 0.0
for i in range(0, 6):
for j in range(0, 6):
deltaTFlat += (
self.AC.b[rating][i * 6 + j] * (M**j) * (delta**i)
)
# deltaT above kink point -> temperature-rated area
thetaT = theta * (1 + (M**2) * (const.Agamma - 1.0) / 2.0)
deltaTTemp = 0.0
for i in range(0, 5):
for j in range(0, 5):
deltaTTemp += (
self.AC.c[rating][i * 5 + j] * (M**j) * (thetaT**i)
)
for i in range(5, 9):
for j in range(0, 5):
deltaTTemp += (
self.AC.c[rating][i * 5 + j]
* (M**j)
* (delta ** (i - 4))
)
# compute deltaT according to DeltaTemp with respect to kink point
if DeltaTemp <= self.AC.kink[rating]:
deltaT = deltaTFlat
else:
deltaT = deltaTTemp
CT = 0.0
for i in range(0, 6):
for j in range(0, 6):
CT += self.AC.a[i * 6 + j] * (M**j) * (deltaT**i)
elif self.AC.engineType == "TURBOPROP":
rating = checkArgument("rating", **kwargs)
CP = self.CP(rating=rating, theta=theta, delta=delta, M=M)
CT = CP / M
elif self.AC.engineType == "PISTON":
CP = checkArgument("CP", **kwargs)
sigma = atm.sigma(theta=theta, delta=delta)
Wp = self.AC.WREF * const.a_0 * CP
TAS = atm.mach2Tas(Mach=M, theta=theta)
propEff = self.propEfficiency(
Wp=Wp * self.AC.n_eng, sigma=sigma, tas=TAS
)
CT = propEff * const.a_0 * CP / (delta * TAS)
return CT
[docs]
def CPmax(self, rating, delta, theta, M):
"""
Computes the maximum engine power coefficient (CPmax) for TURBOPROP engines.
:param rating: Throttle setting {MCMB, MCRZ}.
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param M: Mach number [-].
:type rating: str.
:type delta: float.
:type theta: float.
:type M: float.
:return: Maximum engine power coefficient (CPmax) [-].
:rtype: float.
:raises: ValueError: If the rating is unknown.
"""
if self.AC.engineType == "TURBOPROP":
if rating not in self.AC.max_power.keys():
raise ValueError("Unknown engine rating " + rating)
Wpmax = self.AC.max_power[rating]
aSound = atm.aSound(theta=theta)
tas = atm.mach2Tas(Mach=M, theta=theta)
sigma = atm.sigma(theta=theta, delta=delta)
propEff = self.propEfficiency(Wp=Wpmax, sigma=sigma, tas=tas)
if propEff is None:
return None
CPmax = Wpmax * propEff / (delta * self.AC.WREF * aSound)
else:
raise ValueError("CPmax implemented only for turboprop")
return CPmax
[docs]
def CP(self, **kwargs):
"""
Computes the power coefficient (CP) for TURBOPROP and PISTON engines.
:param rating: Throttle setting {MCMB, MCRZ, LIDL}.
:param deltaT: Direct throttle parameter [-].
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param sigma: Normalized density [-] (for piston engines).
:param M: Mach number [-].
:type rating: str.
:type deltaT: float.
:type delta: float.
:type theta: float.
:type sigma: float.
:type M: float.
:return: Power coefficient (CP) [-].
:rtype: float.
:raises: ValueError if an unknown rating is provided.
"""
if self.AC.engineType == "TURBOPROP":
M = checkArgument("M", **kwargs)
# CT as input parameter
# computes the power coefficient from thrust coefficient assuming efficiency of 1
if "CT" in kwargs:
CT = kwargs.get("CT")
CP = CT * M
return CP
# rating as input parameter
elif "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
# in case MCRZ rating is not defined, we switch to MCMB
if rating == "MCRZ" and rating not in self.AC.max_power.keys():
rating = "MCMB"
delta = checkArgument("delta", **kwargs)
theta = checkArgument("theta", **kwargs)
deltaT = 0.0
for i in range(0, 6):
for j in range(0, 6):
deltaT += (
self.AC.p[rating][i * 6 + j] * (M**j) * (theta**i)
)
CP = 0.0
for i in range(0, 6):
for j in range(0, 6):
CP += self.AC.a[i * 6 + j] * (M**j) * (deltaT**i)
CPmax = self.CPmax(
rating=rating, theta=theta, delta=delta, M=M
)
CP = min(CP, CPmax)
# deltaT - direct throttle as input parameter
elif "deltaT" in kwargs:
deltaT = checkArgument("deltaT", **kwargs)
CP = 0.0
for i in range(0, 6):
for j in range(0, 6):
CP += self.AC.a[i * 6 + j] * (M**j) * (deltaT**i)
elif self.AC.engineType == "PISTON":
delta = checkArgument("delta", **kwargs)
theta = checkArgument("theta", **kwargs)
sigma = checkArgument("sigma", **kwargs)
DeltaTemp = checkArgument("DeltaTemp", **kwargs)
if self.AC.Hd_turbo <= 0:
theta_turbo = atm.theta(
conv.ft2m(self.AC.Hd_turbo), DeltaTemp=DeltaTemp
)
delta_turbo = atm.delta(
conv.ft2m(self.AC.Hd_turbo), DeltaTemp=DeltaTemp
)
# ensure that all real sigmas are smaller than sigma_turbo
sigma_turbo = float("Inf")
else:
sigma_turbo = atm.sigma(theta=theta_turbo, delta=delta_turbo)
CPmaxStdMSL = (
conv.hp2W(self.AC.P)
* self.AC.n_eng
/ (self.AC.WREF * const.a_0)
)
# deltaT - direct throttle as input parameter
if "deltaT" in kwargs:
deltaT = checkArgument("deltaT", **kwargs)
# rating as input parameter
elif "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
if rating == "LIDL":
if self.AC.BADAVersion == "4.2":
deltaT = 0.0
elif (
self.AC.BADAVersion == "4.3"
or self.AC.BADAVersion == "DUMMY"
):
deltaT = 0.1
elif rating == "MCMB" or rating == "MCRZ":
deltaT = 1.0
CPstdMSL = CPmaxStdMSL * deltaT
if sigma >= sigma_turbo:
CP = CPstdMSL
else:
CP = min(
CPstdMSL,
CPmaxStdMSL
* delta
* sqrt(theta_turbo)
/ (delta_turbo * sqrt(theta)),
)
return CP
[docs]
def CDClean(self, CL, M):
"""
Computes the drag coefficient (CD) in a clean configuration based on the Mach number (M) and the lift coefficient (CL).
:param M: Mach number [-].
:param CL: Lift coefficient [-].
:type M: float.
:type CL: float.
:return: Drag coefficient (CD) in clean configuration [-].
:rtype: float
"""
param = 1 - M * M
d = self.AC.d[0]["LGUP"]
C0 = (
d[0]
+ d[1] / sqrt(param)
+ d[2] / (param)
+ d[3] / pow(param, 3.0 / 2.0)
+ d[4] / pow(param, 2.0)
)
C2 = (
d[5]
+ d[6] / pow(param, 3.0 / 2.0)
+ d[7] / pow(param, 3.0)
+ d[8] / pow(param, 9.0 / 2.0)
+ d[9] / pow(param, 6.0)
)
C6 = (
d[10]
+ d[11] / pow(param, 7.0)
+ d[12] / pow(param, 15.0 / 2.0)
+ d[13] / pow(param, 8.0)
+ d[14] / pow(param, 17.0 / 2.0)
)
CD_clean = self.AC.scalar * (C0 + C2 * CL * CL + C6 * pow(CL, 6))
return CD_clean
[docs]
def CD(
self,
HLid,
LG,
CL,
M,
speedBrakes={"deployed": False, "value": 0.03},
**kwargs,
):
"""
Computes the drag coefficient (CD) based on the Mach number (M), lift coefficient (CL), high lift devices (HLid),
landing gear position (LG), and speed brakes status. The drag coefficient is calculated for both clean and non-clean configurations.
:param M: Mach number [-].
:param CL: Lift coefficient [-].
:param HLid: High lift devices position [-].
:param LG: Landing gear position, [LGUP/LGDN] [-].
:param speedBrakes: Dictionary indicating if speed brakes are deployed and their value. Default: {"deployed": False, "value": 0.03}.
:type M: float.
:type CL: float.
:type HLid: float.
:type LG: string.
:type speedBrakes: dict.
:returns: Drag coefficient [-].
:rtype: float.
"""
# clean configuration
if HLid == 0 and LG == "LGUP":
# below Mmax
if M <= self.AC.M_max:
CD = self.AC.CDClean(M=M, CL=CL)
# above Mmax (accounting for air compresibility)
else:
CD = self.CDClean(M=self.AC.M_max - 0.01, CL=CL) + pow(
(M - (self.AC.M_max - 0.01)) / 0.01, 3 / 2
) * (
self.CDClean(M=self.AC.M_max, CL=CL)
- self.CDClean(M=self.AC.M_max - 0.01, CL=CL)
)
# non-clean configuration
else:
CD = (
self.AC.d[HLid][LG][0]
+ self.AC.d[HLid][LG][1] * CL
+ self.AC.d[HLid][LG][2] * CL * CL
)
# implementation of a simple speed brakes model
if speedBrakes["deployed"]:
if speedBrakes["value"] is not None:
CD = CD + speedBrakes["value"]
else:
CD = CD + 0.03
return CD
# calculation of drag coefficient in transition for HLid assuming LG is not changing
if "HLid_init" in kwargs and "HLid_final" in kwargs:
HLid_init = checkArgument("HLid_init", **kwargs)
HLid_final = checkArgument("HLid_final", **kwargs)
LG_init = LG
LG_final = LG
if HLid_init == 0 and LG_init == "LGUP":
CD_init = self.AC.CDClean(M=M, CL=CL)
else:
CD_init = (
self.AC.d[HLid_init][LG_init][0]
+ self.AC.d[HLid_init][LG_init][1] * CL
+ self.AC.d[HLid_init][LG_init][2] * CL
)
if HLid_final == 0 and LG_final == "LGUP":
CD_final = self.CDClean(M=M, CL=CL)
else:
CD_final = (
self.AC.d[HLid_final][LG_final][0]
+ self.AC.d[HLid_final][LG_final][1] * CL
+ self.AC.d[HLid_final][LG_final][2] * CL
)
# linear interpolation
xp = [HLid_init, HLid_final]
fp = [CD_init, CD_final]
CD = np.interp(HLid, xp, fp)
# calculation of drag coefficient in transition for LG assuming HLid is not changing
if "LG_init" in kwargs and "LG_final" in kwargs:
LG_init = checkArgument("LG_init", **kwargs)
LG_final = checkArgument("LG_final", **kwargs)
HLid_init = HLid
HLid_final = HLid
if HLid_init == 0 and LG_init == "LGUP":
CD_init = self.CDClean(M=M, CL=CL)
else:
CD_init = (
self.AC.d[HLid_init][LG_init][0]
+ self.AC.d[HLid_init][LG_init][1] * CL
+ self.AC.d[HLid_init][LG_init][2] * CL
)
if HLid_final == 0 and LG_final == "LGUP":
CD_final = self.CDClean(M=M, CL=CL)
else:
CD_final = (
self.AC.d[HLid_final][LG_final][0]
+ self.AC.d[HLid_final][LG_final][1] * CL
+ self.AC.d[HLid_final][LG_final][2] * CL
)
# linear interpolation
xp = [HLid_init, HLid_final]
fp = [CD_init, CD_final]
CD = np.interp(HLid, xp, fp)
return CD
[docs]
def L(self, delta, M, CL):
"""
Computes the aerodynamic lift based on the Mach number (M), normalized air pressure (delta), and lift coefficient (CL).
:param M: Mach number [-].
:param delta: Normalized air pressure [-].
:param CL: Lift coefficient [-].
:type M: float.
:type delta: float.
:type CL: float.
:returns: Aerodynamic lift [N].
:rtype: float.
"""
return 0.5 * delta * const.p_0 * const.Agamma * M * M * self.AC.S * CL
[docs]
def D(self, delta, M, CD):
"""
Computes the thrust based on throttle settings, normalized air pressure (delta), and other flight parameters.
:param rating: Throttle setting {MCMB, MCRZ, LIDL}.
:param deltaT: Direct throttle parameter [-].
:param delta: Normalized air pressure [-].
:param theta: Normalized temperature [-].
:param M: Mach number [-].
:param DeltaTemp: Temperature deviation with respect to ISA [K].
:type rating: string.
:type deltaT: float.
:type delta: float.
:type theta: float.
:type M: float.
:type DeltaTemp: float.
:returns: Thrust [N].
:rtype: float.
"""
return 0.5 * delta * const.p_0 * const.Agamma * M * M * self.AC.S * CD
[docs]
def Thrust(self, delta, **kwargs):
"""
Computes the maximum thrust produced by the aircraft based on throttle settings,
normalized air pressure (delta), and other flight parameters like Mach number and temperature deviation.
:param rating: Throttle setting {MCMB (Max Climb), MCRZ (Max Cruise), LIDL (Idle)}.
:param deltaT: Direct throttle parameter for intermediate throttle levels [-].
:param delta: Normalized air pressure relative to sea-level pressure (ISA) [-].
:param theta: Normalized temperature relative to sea-level temperature (ISA) [-].
:param M: Mach number, the ratio of the aircraft's speed to the speed of sound [-].
:param DeltaTemp: Temperature deviation from the International Standard Atmosphere (ISA) [K].
:type rating: string (optional).
:type deltaT: float (optional).
:type delta: float.
:type theta: float (optional).
:type M: float (optional).
:type DeltaTemp: float (optional).
:returns: Thrust force produced by the engine [N].
:rtype: float.
"""
CT = self.CT(delta=delta, **kwargs)
return delta * self.AC.WREF * CT
[docs]
def ff(self, delta, theta, DeltaTemp, **kwargs):
"""
Computes the fuel flow (FF) of the aircraft based on engine throttle settings,
normalized air pressure (delta), normalized temperature (theta), and other relevant parameters.
:param rating: Throttle setting {MCMB (Max Climb), MCRZ (Max Cruise), LIDL (Idle), TAXI}.
If 'TAXI' is selected, the taxi fuel allowance is used for ground operations.
:param deltaT: Direct throttle parameter for intermediate throttle levels [-].
:param delta: Normalized air pressure relative to sea-level pressure (ISA) [-].
:param theta: Normalized temperature relative to sea-level temperature (ISA) [-].
:param M: Mach number, the ratio of the aircraft's speed to the speed of sound [-].
:param DeltaTemp: Temperature deviation from the International Standard Atmosphere (ISA) [K].
:type rating: string (optional).
:type deltaT: float (optional).
:type delta: float.
:type theta: float.
:type M: float.
:type DeltaTemp: float.
:returns: Fuel flow rate in kilograms per second [kg s^-1].
:rtype: float.
"""
if "rating" in kwargs:
rating = checkArgument("rating", **kwargs)
if rating == "TAXI":
if self.AC.TFA is not None:
return self.AC.TFA / 60
else:
return None
CF = self.CF(delta=delta, theta=theta, DeltaTemp=DeltaTemp, **kwargs)
if self.AC.BADAVersion == "4.2":
return (
delta
* pow(theta, 0.5)
* self.AC.WREF
* const.a_0
* CF
/ self.AC.LHV
)
elif self.AC.BADAVersion == "4.3" or self.AC.BADAVersion == "DUMMY":
return (
pow(delta, self.AC.p_delta)
* pow(theta, self.AC.p_theta)
* self.AC.WREF
* const.a_0
* CF
/ self.AC.LHV
)
[docs]
def ROCD(self, T, D, v, mass, ESF, h, DeltaTemp):
"""
Computes the Rate of Climb or Descent (ROCD) of the aircraft based on
the provided thrust, drag, true airspeed, mass, and Energy Share Factor (ESF).
:param T: Thrust produced by the engines [N].
:param D: Drag acting on the aircraft [N].
:param v: True airspeed (TAS) of the aircraft [m/s].
:param mass: Current aircraft mass [kg].
:param ESF: Energy Share Factor, which controls the allocation of excess thrust
between climb/descent and acceleration [-].
:param h: Altitude of the aircraft above sea level [m].
:param DeltaTemp: Temperature deviation from the International Standard Atmosphere (ISA) [K].
:type T: float.
:type D: float.
:type v: float.
:type mass: float.
:type ESF: float.
:type h: float.
:type DeltaTemp: float.
:returns: Rate of climb or descent (ROCD) in meters per second [m/s].
:rtype: float.
"""
theta = atm.theta(h=h, DeltaTemp=DeltaTemp)
temp = theta * const.temp_0
ROCD = (
((temp - DeltaTemp) / temp) * (T - D) * v * ESF / (mass * const.g)
)
return ROCD
[docs]
def controlLawThrust(self, ROCD, D, v, mass, ESF, h, DeltaTemp):
"""
Computes the required thrust based on the TEM (Thrust-Equilibrium Method) control law.
This calculation takes into account the aircraft's rate of climb/descent (ROCD),
drag (D), true airspeed (v), and energy share factor (ESF), along with altitude and
temperature deviations.
:param h: Altitude of the aircraft above sea level [m].
:param ROCD: Rate of climb or descent of the aircraft [m/s].
Positive value for climb, negative for descent.
:param D: Drag force acting on the aircraft [N].
:param v: True airspeed (TAS) of the aircraft [m/s].
:param mass: Current aircraft mass [kg].
:param ESF: Energy Share Factor, determining the fraction of excess thrust used for acceleration or climb/descent [-].
:param DeltaTemp: Temperature deviation from the International Standard Atmosphere (ISA) [K].
:type h: float.
:type ROCD: float.
:type D: float.
:type v: float.
:type mass: float.
:type ESF: float.
:type DeltaTemp: float.
:returns: Thrust required to maintain the specified rate of climb or descent [N].
:rtype: float.
"""
theta = atm.theta(h=h, DeltaTemp=DeltaTemp)
temp = theta * const.temp_0
if ROCD == 0.0 or ESF == 0.0:
thrust = (temp / (temp - DeltaTemp)) * (
ROCD * mass * const.g
) / v + D
else:
thrust = (temp / (temp - DeltaTemp)) * (ROCD * mass * const.g) / (
ESF * v
) + D
return thrust
[docs]
def propEfficiency(self, Wp, sigma, tas):
"""
Computes the propeller efficiency of a piston or turboprop engine using momentum theory.
The calculation estimates the efficiency based on power, air density, and true airspeed.
:param Wp: Total engine power output for all engines [W].
:param sigma: Normalized air density relative to sea-level density [-].
:param tas: True airspeed (TAS) of the aircraft [m/s].
:type Wp: float.
:type sigma: float.
:type tas: float.
:returns: Propeller efficiency as a dimensionless ratio [-].
:rtype: float.
"""
if self.AC.engineType == "TURBOPROP" or self.AC.engineType == "PISTON":
a1 = (
2
* (Wp / self.AC.n_eng)
/ (
sigma
* const.rho_0
* self.AC.prop_dia
* self.AC.prop_dia
* pi
* tas
* tas
* tas
* self.AC.max_eff
)
)
a2 = 0.0
a3 = 1.0
a4 = -1 * (self.AC.max_eff)
coef = np.array([a1, a2, a3, a4])
roots = np.roots(coef)
for root in roots:
if not np.iscomplex(root):
eff = float(np.real(root))
return eff
[docs]
class FlightEnvelope(BADA4):
"""This class is a BADA4 aircraft subclass and implements the flight envelope caclulations
following the BADA4 manual.
:param AC: Aircraft object {BADA4}.
:type AC: bada4Aircraft.
"""
def __init__(self, AC):
super().__init__(AC)
[docs]
def maxM(self, LG):
"""
Computes the maximum allowable Mach speed (Mmax) based on the kinematic limitations
of the aircraft and the position of the landing gear (LG).
:param LG: Landing gear position, either "LGUP" (gear up) or "LGDN" (gear down).
:type LG: str
:return: The maximum allowable Mach number (Mmax) based on the aircraft's limitations.
:rtype: float
"""
if LG == "LGUP":
Mmax = self.AC.MMO
else:
if self.AC.MLE is not None:
Mmax = self.AC.MLE
else:
Mmax = self.AC.MMO
return Mmax
[docs]
def maxCAS(self, HLid, LG):
"""
Computes the maximum allowable Calibrated Airspeed (CASmax) based on the position of
the high-lift devices (HLid) and landing gear (LG), following the kinematic limitations
of the aircraft.
:param HLid: Position of high-lift devices (0 for clean configuration, >0 for extended).
:type HLid: float
:param LG: Landing gear position, either "LGUP" (gear up) or "LGDN" (gear down).
:type LG: str
:return: The maximum allowable Calibrated Airspeed (CASmax) in meters per second [m/s].
:rtype: float
"""
if HLid == 0 and LG == "LGUP":
CASmax = conv.kt2ms(self.AC.VMO)
elif HLid > 0 and LG == "LGUP":
if self.AC.VFE[HLid] is not None:
CASmax = conv.kt2ms(self.AC.VFE[HLid])
else:
CASmax = conv.kt2ms(self.AC.VMO)
elif HLid == 0 and LG == "LGDN":
if self.AC.VLE is not None:
CASmax = conv.kt2ms(self.AC.VLE)
else:
CASmax = self.AC.VMO
elif HLid > 0 and LG == "LGDN":
if self.AC.VLE is not None:
CASmax = conv.kt2ms(min(self.AC.VFE[HLid], self.AC.VLE))
else:
CASmax = self.AC.VMO
return CASmax
[docs]
def VMax(self, h, HLid, LG, delta, theta, mass, nz=1.0):
"""
Computes the maximum speed (CAS)
:param h: Altitude in meters [m].
:param HLid: High-lift devices position [-].
:param LG: Landing gear position, either "LGUP" (gear up) or "LGDN" (gear down) [-].
:param theta: Normalized temperature [-].
:param delta: Normalized pressure [-].
:param mass: Aircraft mass in kilograms [kg].
:param nz: Load factor [-], default is 1.0.
:type h: float.
:type HLid: float.
:type LG: str.
:type theta: float.
:type delta: float.
:type mass: float.
:type nz: float.
:return: Maximum allowable calibrated airspeed (CAS) in meters per second [m/s].
:rtype: float.
"""
if self.AC.MMO is not None:
crossoverAlt = atm.crossOver(
cas=self.maxCAS(HLid=HLid, LG=LG), Mach=self.maxM(LG=LG)
)
if h >= crossoverAlt:
M = self.maxM(LG=LG)
M_buffet = self.maxMbuffet(
HLid=HLid, LG=LG, delta=delta, mass=mass, nz=nz
)
if M_buffet is None:
return None
sigma = atm.sigma(theta=theta, delta=delta)
VMax = atm.mach2Cas(
Mach=min(M, M_buffet),
theta=theta,
delta=delta,
sigma=sigma,
)
# if M_buffet == float('-inf'):
# VMax = float('-inf')
else:
VMax = self.maxCAS(HLid=HLid, LG=LG)
else:
VMax = self.maxCAS(HLid=HLid, LG=LG)
return VMax
[docs]
def Vmax_thrustLimited(self, h, mass, DeltaTemp, rating, config):
"""
Computes the maximum speed (CAS) within the certified flight envelope while accounting
for thrust limitations at a given altitude, temperature deviation, and configuration.
:param h: Altitude in meters [m].
:param mass: Aircraft mass in kilograms [kg].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:param rating: Engine rating, options include "MTKF", "MCMB", "MCRZ" [-].
:param config: Aircraft aerodynamic configuration, such as "TO", "IC", "CR" [-].
:type h: float.
:type mass: float.
:type DeltaTemp: float.
:type rating: str.
:type config: str.
:return: Maximum thrust-limited calibrated airspeed (CAS) in meters per second [m/s].
:rtype: float.
"""
[HLid, LG] = self.getAeroConfig(config=config)
[theta, delta, sigma] = atm.atmosphereProperties(
h=h, DeltaTemp=DeltaTemp
)
VmaxCertified = self.VMax(
h=h, HLid=HLid, LG=LG, delta=delta, theta=theta, mass=mass, nz=1.0
)
VminCertified = self.VStall(
theta=theta, delta=delta, mass=mass, HLid=HLid, LG=LG, nz=1.0
)
maxCASList = []
for CAS in np.linspace(
VminCertified, VmaxCertified, num=200, endpoint=True
):
[M, CAS, TAS] = atm.convertSpeed(
v=conv.ms2kt(CAS),
speedType="CAS",
theta=theta,
delta=delta,
sigma=sigma,
)
maxThrust = self.Thrust(
delta=delta,
theta=theta,
M=M,
rating=rating,
DeltaTemp=DeltaTemp,
)
CL = self.CL(delta=delta, mass=mass, M=M, nz=1.0)
CD = self.CD(HLid=HLid, LG=LG, CL=CL, M=M)
Drag = self.D(delta=delta, M=M, CD=CD)
if maxThrust >= Drag:
maxCASList.append(CAS)
if not maxCASList:
return None
else:
return max(maxCASList)
[docs]
def maxMbuffet(self, HLid, LG, delta, mass, nz=1.0):
"""
Computes the maximum allowable Mach number (M) under buffet limitations, where the lift
coefficient (CL) cannot exceed the maximum lift coefficient (CL_max) for the given Mach
number and configuration.
:param HLid: High-lift devices position [-].
:param LG: Landing gear position, either "LGUP" (gear up) or "LGDN" (gear down) [-].
:param delta: Normalized pressure [-].
:param mass: Aircraft mass in kilograms [kg].
:param nz: Load factor [-], default is 1.0.
:type HLid: float.
:type LG: str.
:type delta: float.
:type mass: float.
:type nz: float.
:return: Maximum allowable Mach number (M) limited by buffet conditions.
:rtype: float.
"""
# if CLMax model exist for aircraft, additional limitation apply
if self.AC.CL_clean is not None:
M_list = np.arange(0.01, self.AC.MMO + 0.001, 0.001)
# start from maximum value, since we are looking for max M
idx = -1
M = M_list[idx]
while True:
CL = self.CL(delta=delta, mass=mass, M=M, nz=nz)
CL_max = self.CLmax(M=M, HLid=HLid, LG=LG)
if CL_max - CL > 0:
return M
else:
# didn't find any M satisfying the CL < CLmax condition
if abs(idx) == len(M_list):
return None
else:
idx -= 1
M = M_list[idx]
else:
return self.AC.MMO
[docs]
def minMbuffet(self, HLid, LG, theta, delta, mass, nz=1.0):
"""
Computes the minimum Mach number (M) applying buffet limitations, where the lift coefficient (CL)
must not exceed the maximum lift coefficient (CL_max) for the given Mach number and aerodynamic configuration.
:param HLid: High-lift devices position [-].
:param LG: Landing gear position, either "LGUP" (gear up) or "LGDN" (gear down) [-].
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param mass: Aircraft mass in kilograms [kg].
:param nz: Load factor [-], default is 1.0.
:type HLid: float.
:type LG: str.
:type delta: float.
:type theta: float.
:type mass: float.
:type nz: float.
:returns: Minimum Mach number (M) limited by buffet conditions.
:rtype: float.
"""
if HLid in self.AC.CL_max and LG in self.AC.CL_max[HLid]:
CLmax = self.AC.CL_max[HLid][LG]
# estimation of min M where CLmax = CL
Mmin = sqrt(
2
* mass
* const.g
/ (delta * const.p_0 * const.Agamma * CLmax * self.AC.S)
)
return Mmin
else:
if self.AC.MMO is not None:
MMO_max = self.AC.MMO
else:
sigma = atm.sigma(theta=theta, delta=delta)
MMO_max = atm.cas2Mach(
cas=conv.kt2ms(self.AC.VMO),
theta=theta,
delta=delta,
sigma=sigma,
)
M_list = np.arange(0.1, MMO_max + 0.001, 0.001)
idx = 0
M = M_list[idx]
while True:
CL = self.CL(delta=delta, mass=mass, M=M, nz=nz)
CL_max = self.CLmax(M=M, HLid=HLid, LG=LG)
if CL_max - CL > 0:
return M
else:
# didn't find any M satisfying the CL < CLmax condition
if idx == len(M_list) - 1:
return None
else:
idx += 1
M = M_list[idx]
[docs]
def VMin(self, config, theta, delta, mass):
"""
Computes the minimum speed (CAS) for the given aerodynamic configuration, accounting for
stall speed and other configuration-based limitations.
:param config: Aircraft configuration, options include "CR", "IC", "TO", "AP", "LD" [-].
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param mass: Aircraft mass in kilograms [kg].
:type config: str.
:type delta: float.
:type theta: float.
:type mass: float.
:returns: Minimum calibrated airspeed (CAS) in meters per second [m/s].
:rtype: float.
"""
aeroConf = self.getAeroConfig(config=config)
HLid = aeroConf[0]
LG = aeroConf[1]
if (HLid == 0 and LG == "LGUP") and self.AC.CL_clean is not None:
Vmin = self.VStall(
theta=theta, delta=delta, mass=mass, HLid=HLid, LG=LG, nz=1.2
)
else:
if config == "TO":
Vmin = self.AC.CVminTO * self.VStall(
theta=theta,
delta=delta,
mass=mass,
HLid=HLid,
LG=LG,
nz=1.0,
)
else:
Vmin = self.AC.CVmin * self.VStall(
theta=theta,
delta=delta,
mass=mass,
HLid=HLid,
LG=LG,
nz=1.0,
)
return Vmin
[docs]
def VStall(self, mass, HLid, LG, nz=1.0, **kwargs):
"""
Calculates the stall speed (CAS) for the given aerodynamic configuration and load factor,
taking into account altitude and temperature deviations from ISA.
:param HLid: High-lift devices position [-].
:param LG: Landing gear position, either "LGUP" (gear up) or "LGDN" (gear down) [-].
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param mass: Aircraft mass in kilograms [kg].
:param nz: Load factor [-], default is 1.0.
:param h: Altitude above mean sea level (AMSL) in meters [m].
:param DeltaTemp: Temperature deviation from ISA in Kelvin [K].
:type HLid: float.
:type LG: str.
:type delta: float.
:type theta: float.
:type mass: float.
:type nz: float.
:type h: float, optional.
:type DeltaTemp: float, optional.
:returns: Stall speed in calibrated airspeed (CAS) [m/s].
:rtype: float.
"""
if "h" in kwargs:
h = kwargs.get("h")
DeltaTemp = checkArgument("DeltaTemp", **kwargs)
delta = atm.delta(h, DeltaTemp)
theta = atm.theta(h, DeltaTemp)
else:
theta = checkArgument("theta", **kwargs)
delta = checkArgument("delta", **kwargs)
sigma = atm.sigma(theta=theta, delta=delta)
minM = self.minMbuffet(
theta=theta, delta=delta, mass=mass, HLid=HLid, LG=LG, nz=nz
)
if minM is None:
return None
minCAS = atm.mach2Cas(Mach=minM, theta=theta, delta=delta, sigma=sigma)
return minCAS
[docs]
def maxAltitude(self, HLid, LG, M, DeltaTemp, mass, nz=1.0):
"""
Computes the maximum altitude taking into account buffet limitations. The altitude is calculated
based on the aerodynamic configuration and the available buffet boundary conditions.
:param HLid: High-lift devices position [-].
:param LG: Landing gear position, either "LGUP" (gear up) or "LGDN" (gear down) [-].
:param M: Mach airspeed [-].
:param mass: Aircraft mass [kg].
:param nz: Load factor [-], default is 1.0.
:param DeltaTemp: Temperature deviation from ISA [K].
:type HLid: float.
:type LG: str.
:type M: float.
:type mass: float.
:type nz: float.
:type DeltaTemp: float.
:returns: Maximum altitude [m].
:rtype: float.
"""
if HLid > 0:
if self.AC.mfa is not None:
hMax = self.AC.mfa
else:
hMax = self.AC.hmo
else:
hMax = self.AC.hmo
if self.AC.CL_clean is not None:
def f(H):
delta = atm.delta(h=H[0], DeltaTemp=DeltaTemp)
CL = self.CL(delta=delta, mass=mass, M=M, nz=nz)
CL_max = self.CLmax(M=M, HLid=HLid, LG=LG)
return -CL - CL_max
hMax = float(
fminbound(
f,
x1=np.array([0]),
x2=np.array([conv.ft2m(hMax)]),
disp=False,
)
)
return hMax
[docs]
def getConfig(self, phase, h, mass, v, DeltaTemp=0.0, hRWY=0.0):
"""
Returns the aircraft's aerodynamic configuration based on altitude, speed, and phase of flight.
:param phase: Phase of flight [Climb, Cruise, Descent] [-].
:param h: Altitude above mean sea level (AMSL) [m].
:param mass: Aircraft mass [kg].
:param v: Calibrated airspeed (CAS) [m/s].
:param DeltaTemp: Deviation from ISA temperature [K], default is 0.0.
:param hRWY: Runway elevation above mean sea level (AMSL) [m], default is 0.0.
:type phase: str.
:type h: float.
:type mass: float.
:type v: float.
:type DeltaTemp: float.
:type hRWY: float.
:returns: Aircraft aerodynamic configuration [TO/IC/CR/AP/LD] [-].
:rtype: str.
:raises: TypeError if unable to determine the configuration.
"""
config = None
[theta, delta, sigma] = atm.atmosphereProperties(
h=h, DeltaTemp=DeltaTemp
)
# aircraft AGL altitude assuming being close to the RWY [m]
h_AGL = h - hRWY
HmaxTO_AGL = conv.ft2m(self.AC.HmaxPhase["TO"]) - hRWY
HmaxIC_AGL = conv.ft2m(self.AC.HmaxPhase["IC"]) - hRWY
HmaxAPP_AGL = conv.ft2m(self.AC.HmaxPhase["AP"]) - hRWY
HmaxLD_AGL = conv.ft2m(self.AC.HmaxPhase["LD"]) - hRWY
if phase == "Climb" and h_AGL <= HmaxTO_AGL:
config = "TO"
return config
elif phase == "Climb" and (h_AGL > HmaxTO_AGL and h_AGL <= HmaxIC_AGL):
config = "IC"
return config
elif (
phase == "Cruise"
or (phase == "Climb" and h_AGL >= HmaxIC_AGL)
or (phase == "Descent" and h_AGL >= HmaxAPP_AGL)
):
config = "CR"
return config
else:
vMinCR = self.VMin(
config="CR", mass=mass, theta=theta, delta=delta
)
vMinAPP = self.VMin(
config="AP", mass=mass, theta=theta, delta=delta
)
ep = 1e-6
if (
phase == "Descent"
and (h_AGL + ep) < HmaxLD_AGL
and (v + ep) < (vMinAPP + conv.kt2ms(10))
):
config = "LD"
elif (
phase == "Descent"
and h_AGL >= HmaxLD_AGL
and (h_AGL + ep) < HmaxAPP_AGL
and v < (vMinCR + conv.kt2ms(10))
) or (
phase == "Descent"
and (h_AGL + ep) < HmaxLD_AGL
and (
(v + ep) < (vMinCR + conv.kt2ms(10))
and v >= (vMinAPP + conv.kt2ms(10))
)
):
config = "AP"
elif (
phase == "Descent"
and (h_AGL + ep) < HmaxAPP_AGL
and v >= (vMinCR + conv.kt2ms(10))
):
config = "CR"
if config is None:
raise TypeError("Unable to determine aircraft configuration")
return config
[docs]
def getAeroConfig(self, config):
"""
Returns the aircraft aerodynamic configuration based on the provided configuration ID.
This includes the high-lift device (HLID) position and landing gear (LG) position.
If the configuration is not found, it returns None.
:param config: Aircraft configuration (TO/IC/CR/AP/LD).
:type config: str
:returns: Aerodynamic configuration as a combination of HLID and LG.
:rtype: [float, str]
"""
configDict = self.AC.aeroConfig.get(config)
return [configDict["HLid"], configDict["LG"]]
[docs]
def getSpeedSchedule(self, phase):
"""
Returns the speed schedule based on the phase of flight (Climb, Cruise, Descent).
The schedule includes two CAS values (CAS1, CAS2) and a Mach number (M).
:param phase: Aircraft phase of flight (Climb, Cruise, Descent).
:type phase: str
:returns: Speed schedule as a combination of CAS1, CAS2 (in m/s) and Mach number (M).
:rtype: [float, float, float]
"""
speedScheduleDict = self.AC.speedSchedule[phase]
return [
conv.kt2ms(speedScheduleDict["CAS1"]),
conv.kt2ms(speedScheduleDict["CAS2"]),
speedScheduleDict["M"],
]
[docs]
def checkConfigurationContinuity(
self, phase, previousConfig, currentConfig
):
"""
Ensures the continuity of the aerodynamic configuration changes based on the phase of flight.
It prevents sudden or improper configuration transitions, ensuring the aerodynamic configuration
does not change in the wrong direction during Climb, Cruise, or Descent.
:param phase: Aircraft phase of flight (Climb, Cruise, Descent).
:param previousConfig: Previous aerodynamic configuration.
:param currentConfig: Current aerodynamic configuration.
:type phase: str
:type previousConfig: str
:type currentConfig: str
:returns: The appropriate configuration for the current flight phase.
:rtype: str
"""
newConfig = ""
# previous configuration is NOT empty/unknown
if previousConfig is not None:
if phase == "Descent":
if currentConfig == "CR" and (
previousConfig == "AP" or previousConfig == "LD"
):
newConfig = previousConfig
elif currentConfig == "AP" and previousConfig == "LD":
newConfig = previousConfig
else:
newConfig = currentConfig
elif phase == "Climb":
if currentConfig == "TO" and (
previousConfig == "IC" or previousConfig == "CR"
):
newConfig = previousConfig
elif currentConfig == "IC" and previousConfig == "CR":
newConfig = previousConfig
else:
newConfig = currentConfig
elif phase == "Cruise":
newConfig = currentConfig
# previous configuration is empty/unknown
else:
newConfig = currentConfig
return newConfig
[docs]
class ARPM(BADA4):
"""This class is a BADA4 aircraft subclass and implements the Airline Procedure Model (ARPM)
following the BADA4 user manual.
:param AC: Aircraft object {BADA4}.
:type AC: bada4Aircraft.
"""
def __init__(self, AC):
super().__init__(AC)
self.flightEnvelope = FlightEnvelope(AC)
[docs]
def climbSpeed(
self,
theta,
delta,
mass,
h,
hRWY=0.0,
speedSchedule_default=None,
procedure="BADA",
config=None,
NADP1_ALT=3000,
NADP2_ALT=[1000, 3000],
DeltaTemp=0.0,
):
"""
Computes the climb speed schedule (CAS) for the given altitude based on various procedures and aircraft parameters.
:param theta: Normalized air temperature [-].
:param delta: Normalized air pressure [-].
:param mass: Aircraft mass in kilograms [kg].
:param h: Altitude in meters [m].
:param hRWY: Runway elevation AMSL in meters [m].
:param speedSchedule_default: Optional, a default speed schedule that overrides the BADA schedule. It should be in the form [Vcl1, Vcl2, Mcl].
:param procedure: Climb procedure to be followed, e.g., 'BADA', 'NADP1', 'NADP2'. Default is 'BADA'.
:param config: Optional, current aircraft aerodynamic configuration (TO/IC/CR/AP/LD).
:param NADP1_ALT: Altitude in feet for NADP1 procedure. Default is 3000 feet.
:param NADP2_ALT: Altitude range in feet for NADP2 procedure. Default is [1000, 3000].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type theta: float
:type delta: float
:type mass: float
:type h: float
:type hRWY: float, optional
:type speedSchedule_default: list[float], optional
:type procedure: str
:type config: str, optional
:type NADP1_ALT: float, optional
:type NADP2_ALT: list[float], optional
:type DeltaTemp: float, optional
:returns: A tuple containing the climb calibrated airspeed (CAS) in meters per second [m/s] and a status flag indicating whether the calculated CAS is constrained ('C'), unconstrained ('V' or 'v'), or not altered ('').
:rtype: tuple[float, str]
This function computes the climb speed schedule for different phases of flight and aircraft types.
It supports BADA, NADP1, and NADP2 procedures for both jet and turboprop/piston/electric aircraft.
The climb schedule uses specific speed profiles depending on altitude and aircraft model. For jet engines, the speed is constrained
below 250 knots below 10,000 feet, and then it follows a defined speed schedule, either from BADA or NADP procedures.
Additionally, the function applies speed limits based on the aircraft's flight envelope, adjusting the calculated climb speed if necessary.
- For procedure='BADA', it uses the BADA climb speed schedule.
- For procedure='NADP1', it implements the Noise Abatement Departure Procedure 1.
- For procedure='NADP2', it implements the Noise Abatement Departure Procedure 2.
The function also ensures that the calculated CAS remains within the bounds of the aircraft's minimum and maximum speeds.
"""
# aircraft AGL altitude assuming being close to the RWY [m]
h_AGL = h - hRWY
phase = "Climb"
acModel = self.AC.engineType
[HLidTO, LGTO] = self.flightEnvelope.getAeroConfig(config="TO")
VstallTO = self.flightEnvelope.VStall(
h=h_AGL,
mass=mass,
HLid=HLidTO,
LG=LGTO,
nz=1.0,
DeltaTemp=DeltaTemp,
)
[HLidCR, LGCR] = self.flightEnvelope.getAeroConfig(config="CR")
VstallCR = self.flightEnvelope.VStall(
h=h_AGL,
mass=mass,
HLid=HLidCR,
LG=LGCR,
nz=1.0,
DeltaTemp=DeltaTemp,
)
[Vcl1, Vcl2, Mcl] = self.flightEnvelope.getSpeedSchedule(phase=phase)
if speedSchedule_default is not None:
Vcl1 = speedSchedule_default[0]
Vcl2 = speedSchedule_default[1]
Mcl = speedSchedule_default[2]
crossOverAlt = atm.crossOver(cas=Vcl2, Mach=Mcl)
if procedure == "BADA":
if acModel == "JET":
speed = list()
speed.append(min(Vcl1, conv.kt2ms(250)))
speed.append(
self.AC.CVmin * VstallTO + conv.kt2ms(self.AC.V_cl[5])
)
speed.append(
self.AC.CVmin * VstallTO + conv.kt2ms(self.AC.V_cl[4])
)
speed.append(
self.AC.CVmin * VstallTO + conv.kt2ms(self.AC.V_cl[3])
)
speed.append(
self.AC.CVmin * VstallTO + conv.kt2ms(self.AC.V_cl[2])
)
speed.append(
self.AC.CVmin * VstallTO + conv.kt2ms(self.AC.V_cl[1])
)
n = 1
while n < len(speed):
if speed[n] > speed[n - 1]:
speed[n] = speed[n - 1]
n = n + 1
if h < conv.ft2m(1500):
cas = speed[5]
elif h >= conv.ft2m(1500) and h < conv.ft2m(3000):
cas = speed[4]
elif h >= conv.ft2m(3000) and h < conv.ft2m(4000):
cas = speed[3]
elif h >= conv.ft2m(4000) and h < conv.ft2m(5000):
cas = speed[2]
elif h >= conv.ft2m(5000) and h < conv.ft2m(6000):
cas = speed[1]
elif h >= conv.ft2m(6000) and h < conv.ft2m(10000):
cas = speed[0]
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vcl2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mcl, theta=theta, delta=delta, sigma=sigma
)
elif acModel == "TURBOPROP" or acModel == "PISTON":
speed = list()
speed.append(min(Vcl1, conv.kt2ms(250)))
speed.append(
self.AC.CVmin * VstallTO + conv.kt2ms(self.AC.V_cl[8])
)
speed.append(
self.AC.CVmin * VstallTO + conv.kt2ms(self.AC.V_cl[7])
)
speed.append(
self.AC.CVmin * VstallTO + conv.kt2ms(self.AC.V_cl[6])
)
n = 1
while n < len(speed):
if speed[n] > speed[n - 1]:
speed[n] = speed[n - 1]
n = n + 1
if h < conv.ft2m(500):
cas = speed[3]
elif h >= conv.ft2m(500) and h < conv.ft2m(1000):
cas = speed[2]
elif h >= conv.ft2m(1000) and h < conv.ft2m(1500):
cas = speed[1]
elif h >= conv.ft2m(1500) and h < conv.ft2m(10000):
cas = speed[0]
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vcl2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mcl, theta=theta, delta=delta, sigma=sigma
)
elif procedure == "NADP1":
if acModel == "JET":
speed = list()
speed.append(min(Vcl1, conv.kt2ms(250)))
speed.append(
self.AC.CVminTO * VstallTO + conv.kt2ms(self.AC.V_cl[2])
)
n = 1
while n < len(speed):
if speed[n] > speed[n - 1]:
speed[n] = speed[n - 1]
n = n + 1
if h < conv.ft2m(NADP1_ALT):
cas = speed[1]
elif h >= conv.ft2m(NADP1_ALT) and h < conv.ft2m(10000):
cas = speed[0]
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vcl2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mcl, theta=theta, delta=delta, sigma=sigma
)
elif acModel == "TURBOPROP" or acModel == "PISTON":
speed = list()
speed.append(min(Vcl1, conv.kt2ms(250)))
speed.append(
self.AC.CVminTO * VstallTO + conv.kt2ms(self.AC.V_cl[1])
)
n = 1
while n < len(speed):
if speed[n] > speed[n - 1]:
speed[n] = speed[n - 1]
n = n + 1
if h < conv.ft2m(NADP1_ALT):
cas = speed[1]
elif h >= conv.ft2m(NADP1_ALT) and h < conv.ft2m(10000):
cas = speed[0]
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vcl2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mcl, theta=theta, delta=delta, sigma=sigma
)
elif procedure == "NADP2":
if acModel == "JET":
speed = list()
speed.append(min(Vcl1, conv.kt2ms(250)))
speed.append(
self.AC.CVmin * VstallCR + conv.kt2ms(self.AC.V_cl[2])
)
speed.append(
self.AC.CVminTO * VstallTO + conv.kt2ms(self.AC.V_cl[2])
)
n = 1
while n < len(speed):
if speed[n] > speed[n - 1]:
speed[n] = speed[n - 1]
n = n + 1
if h < conv.ft2m(NADP2_ALT[0]):
cas = speed[2]
elif h >= conv.ft2m(NADP2_ALT[0]) and h < conv.ft2m(
NADP2_ALT[1]
):
cas = speed[1]
elif h >= conv.ft2m(NADP2_ALT[1]) and h < conv.ft2m(10000):
cas = speed[0]
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vcl2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mcl, theta=theta, delta=delta, sigma=sigma
)
elif acModel == "TURBOPROP" or acModel == "PISTON":
speed = list()
speed.append(min(Vcl1, conv.kt2ms(250)))
speed.append(
self.AC.CVmin * VstallCR + conv.kt2ms(self.AC.V_cl[2])
)
speed.append(
self.AC.CVminTO * VstallTO + conv.kt2ms(self.AC.V_cl[1])
)
n = 1
while n < len(speed):
if speed[n] > speed[n - 1]:
speed[n] = speed[n - 1]
n = n + 1
if h < conv.ft2m(NADP2_ALT[0]):
cas = speed[2]
elif h >= conv.ft2m(NADP2_ALT[0]) and h < conv.ft2m(
NADP2_ALT[1]
):
cas = speed[1]
elif h >= conv.ft2m(NADP2_ALT[1]) and h < conv.ft2m(10000):
cas = speed[0]
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vcl2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mcl, theta=theta, delta=delta, sigma=sigma
)
# check if the speed is within the limits of minimum and maximum speed from the flight envelope, if not, overwrite calculated speed with flight envelope min/max speed
if config is None:
config = self.flightEnvelope.getConfig(
h=h,
phase=phase,
v=cas,
mass=mass,
DeltaTemp=DeltaTemp,
hRWY=hRWY,
)
minSpeed = self.flightEnvelope.VMin(
config=config, mass=mass, theta=theta, delta=delta
)
[HLid, LG] = self.flightEnvelope.getAeroConfig(config=config)
maxSpeed = self.flightEnvelope.VMax(
h=h, HLid=HLid, LG=LG, theta=theta, delta=delta, mass=mass, nz=1.2
)
eps = 1e-6 # float calculation precision
# empty envelope - keep the original calculated CAS speed
if minSpeed is None or maxSpeed is None:
return [cas, "vV"]
if maxSpeed < minSpeed:
if (cas - eps) > maxSpeed and (cas - eps) > minSpeed:
return [cas, "V"]
elif (cas + eps) < minSpeed and (cas + eps) < maxSpeed:
return [cas, "v"]
else:
return [cas, "vV"]
if minSpeed > (cas + eps):
return [minSpeed, "C"]
if maxSpeed < (cas - eps):
return [maxSpeed, "C"]
return [cas, ""]
[docs]
def cruiseSpeed(
self,
theta,
delta,
mass,
h,
hRWY=0.0,
speedSchedule_default=None,
config=None,
DeltaTemp=0.0,
):
"""
Computes the cruise speed schedule (CAS) for the given altitude based on various aircraft parameters.
:param theta: Normalized air temperature [-].
:param delta: Normalized air pressure [-].
:param mass: Aircraft mass in kilograms [kg].
:param h: Altitude in meters [m].
:param hRWY: Runway elevation AMSL in meters [m].
:param speedSchedule_default: Optional, a default speed schedule that overrides the BADA schedule. It should be in the form [Vcr1, Vcr2, Mcr].
:param config: Optional, current aircraft aerodynamic configuration (TO/IC/CR/AP/LD).
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type theta: float
:type delta: float
:type mass: float
:type h: float
:type hRWY: float, optional
:type speedSchedule_default: list[float], optional
:type config: str, optional
:type DeltaTemp: float, optional
:returns: A tuple containing the cruise calibrated airspeed (CAS) in meters per second [m/s] and a status flag indicating whether the calculated CAS is constrained ('C'), unconstrained ('V' or 'v'), or not altered ('').
:rtype: tuple[float, str]
This function computes the cruise speed schedule for different phases of flight and aircraft types. It uses either the default
speed schedule or the BADA speed schedule based on the aircraft model and altitude.
The cruise speed schedule varies depending on the altitude and type of engine (JET, TURBOPROP, or PISTON).
The function also applies speed limits based on the aircraft's flight envelope, ensuring the calculated cruise speed remains within
the aircraft's minimum and maximum allowable speeds.
The function ensures the calculated CAS remains within the aircraft's operational speed limits, adjusting the speed if necessary.
"""
phase = "Cruise"
acModel = self.AC.engineType
[Vcr1, Vcr2, Mcr] = self.flightEnvelope.getSpeedSchedule(phase=phase)
if speedSchedule_default is not None:
Vcr1 = speedSchedule_default[0]
Vcr2 = speedSchedule_default[1]
Mcr = speedSchedule_default[2]
crossOverAlt = atm.crossOver(cas=Vcr2, Mach=Mcr)
if acModel == "JET":
if h < conv.ft2m(3000):
cas = min(Vcr1, conv.kt2ms(170))
elif h >= conv.ft2m(3000) and h < conv.ft2m(6000):
cas = min(Vcr1, conv.kt2ms(220))
elif h >= conv.ft2m(6000) and h < conv.ft2m(14000):
cas = min(Vcr1, conv.kt2ms(250))
elif h >= conv.ft2m(14000) and h < crossOverAlt:
cas = Vcr2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mcr, theta=theta, delta=delta, sigma=sigma
)
elif acModel == "TURBOPROP" or acModel == "PISTON":
if h < conv.ft2m(3000):
cas = min(Vcr1, conv.kt2ms(150))
elif h >= conv.ft2m(3000) and h < conv.ft2m(6000):
cas = min(Vcr1, conv.kt2ms(180))
elif h >= conv.ft2m(6000) and h < conv.ft2m(10000):
cas = min(Vcr1, conv.kt2ms(250))
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vcr2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mcr, theta=theta, delta=delta, sigma=sigma
)
# check if the speed is within the limits of minimum and maximum speed from the flight envelope, if not, overwrite calculated speed with flight envelope min/max speed
if config is None:
config = self.flightEnvelope.getConfig(
h=h,
phase=phase,
v=cas,
mass=mass,
DeltaTemp=DeltaTemp,
hRWY=hRWY,
)
minSpeed = self.flightEnvelope.VMin(
config=config, mass=mass, theta=theta, delta=delta
)
[HLid, LG] = self.flightEnvelope.getAeroConfig(config=config)
maxSpeed = self.flightEnvelope.VMax(
h=h, HLid=HLid, LG=LG, theta=theta, delta=delta, mass=mass, nz=1.2
)
eps = 1e-6 # float calculation precision
# empty envelope - keep the original calculated CAS speed
if minSpeed is None or maxSpeed is None:
return [cas, "vV"]
if maxSpeed < minSpeed:
if (cas - eps) > maxSpeed and (cas - eps) > minSpeed:
return [cas, "V"]
elif (cas + eps) < minSpeed and (cas + eps) < maxSpeed:
return [cas, "v"]
else:
return [cas, "vV"]
if minSpeed > (cas + eps):
return [minSpeed, "C"]
if maxSpeed < (cas - eps):
return [maxSpeed, "C"]
return [cas, ""]
[docs]
def descentSpeed(
self,
theta,
delta,
mass,
h,
hRWY=0.0,
speedSchedule_default=None,
config=None,
DeltaTemp=0.0,
):
"""
Computes the descent speed schedule (CAS) for the given altitude based on various aircraft parameters.
:param theta: Normalized air temperature [-].
:param delta: Normalized air pressure [-].
:param mass: Aircraft mass in kilograms [kg].
:param h: Altitude in meters [m].
:param hRWY: Runway elevation AMSL in meters [m].
:param speedSchedule_default: Optional, a default speed schedule that overrides the BADA schedule. It should be in the form [Vdes1, Vdes2, Mdes].
:param config: Optional, current aircraft aerodynamic configuration (TO/IC/CR/AP/LD).
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type theta: float
:type delta: float
:type mass: float
:type h: float
:type hRWY: float, optional
:type speedSchedule_default: list[float], optional
:type config: str, optional
:type DeltaTemp: float, optional
:returns: A tuple containing the descent calibrated airspeed (CAS) in meters per second [m/s] and a status flag indicating whether the calculated CAS is constrained ('C'), unconstrained ('V' or 'v'), or not altered ('').
:rtype: tuple[float, str]
This function computes the descent speed schedule for different phases of flight and aircraft types.
It uses either the default speed schedule or the BADA speed schedule based on the aircraft model and altitude.
The descent speed schedule varies depending on the altitude and type of engine (JET, TURBOPROP, or PISTON).
The function ensures the calculated CAS remains within the aircraft's operational speed limits, adjusting the speed if necessary.
"""
# aircraft AGL altitude assuming being close to the RWY [m]
h_AGL = h - hRWY
phase = "Descent"
acModel = self.AC.engineType
[HLid, LG] = self.flightEnvelope.getAeroConfig(config="LD")
VstallDES = self.flightEnvelope.VStall(
h=h_AGL, mass=mass, HLid=HLid, LG=LG, nz=1.0, DeltaTemp=DeltaTemp
)
[Vdes1, Vdes2, Mdes] = self.flightEnvelope.getSpeedSchedule(
phase=phase
)
if speedSchedule_default is not None:
Vdes1 = speedSchedule_default[0]
Vdes2 = speedSchedule_default[1]
Mdes = speedSchedule_default[2]
crossOverAlt = atm.crossOver(cas=Vdes2, Mach=Mdes)
if acModel == "JET" or acModel == "TURBOPROP":
speed = []
speed.append(min(Vdes1, conv.kt2ms(220)))
speed.append(
self.AC.CVmin * VstallDES + conv.kt2ms(self.AC.V_des[4])
)
speed.append(
self.AC.CVmin * VstallDES + conv.kt2ms(self.AC.V_des[3])
)
speed.append(
self.AC.CVmin * VstallDES + conv.kt2ms(self.AC.V_des[2])
)
speed.append(
self.AC.CVmin * VstallDES + conv.kt2ms(self.AC.V_des[1])
)
n = 1
while n < len(speed):
if speed[n] > speed[n - 1]:
speed[n] = speed[n - 1]
n = n + 1
epsilon = 1e-6
if h < conv.ft2m(1000):
cas = speed[4]
elif h >= conv.ft2m(1000) and h < conv.ft2m(1500):
cas = speed[3]
elif h >= conv.ft2m(1500) and h < conv.ft2m(2000):
cas = speed[2]
elif h >= conv.ft2m(2000) and h < conv.ft2m(3000):
cas = speed[1]
elif h >= conv.ft2m(3000) and h < conv.ft2m(6000):
cas = min(Vdes1, conv.kt2ms(220))
elif h >= conv.ft2m(6000) and h < conv.ft2m(10000):
cas = min(Vdes1, conv.kt2ms(250))
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vdes2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mdes, theta=theta, delta=delta, sigma=sigma
)
elif acModel == "PISTON":
speed = list()
speed.append(Vdes1)
speed.append(
self.AC.CVmin * VstallDES + conv.kt2ms(self.AC.V_des[7])
)
speed.append(
self.AC.CVmin * VstallDES + conv.kt2ms(self.AC.V_des[6])
)
speed.append(
self.AC.CVmin * VstallDES + conv.kt2ms(self.AC.V_des[5])
)
n = 1
while n < len(speed):
if speed[n] > speed[n - 1]:
speed[n] = speed[n - 1]
n = n + 1
if h < conv.ft2m(500):
cas = speed[3]
elif h >= conv.ft2m(500) and h < conv.ft2m(1000):
cas = speed[2]
elif h >= conv.ft2m(1000) and h < conv.ft2m(1500):
cas = speed[1]
elif h >= conv.ft2m(1500) and h < conv.ft2m(10000):
cas = speed[0]
elif h >= conv.ft2m(10000) and h < crossOverAlt:
cas = Vdes2
elif h >= crossOverAlt:
sigma = atm.sigma(theta=theta, delta=delta)
cas = atm.mach2Cas(
Mach=Mdes, theta=theta, delta=delta, sigma=sigma
)
# check if the speed is within the limits of minimum and maximum speed from the flight envelope, if not, overwrite calculated speed with flight envelope min/max speed
if config is None:
config = self.flightEnvelope.getConfig(
h=h,
phase=phase,
v=cas,
mass=mass,
DeltaTemp=DeltaTemp,
hRWY=hRWY,
)
minSpeed = self.flightEnvelope.VMin(
config=config, mass=mass, theta=theta, delta=delta
)
[HLid, LG] = self.flightEnvelope.getAeroConfig(config=config)
maxSpeed = self.flightEnvelope.VMax(
h=h, HLid=HLid, LG=LG, theta=theta, delta=delta, mass=mass, nz=1.2
)
eps = 1e-6 # float calculation precision
# empty envelope - keep the original calculated CAS speed
if minSpeed is None or maxSpeed is None:
return [cas, "vV"]
if maxSpeed < minSpeed:
if (cas - eps) > maxSpeed and (cas - eps) > minSpeed:
return [cas, "V"]
elif (cas + eps) < minSpeed and (cas + eps) < maxSpeed:
return [cas, "v"]
else:
return [cas, "vV"]
if minSpeed > (cas + eps):
return [minSpeed, "C"]
if maxSpeed < (cas - eps):
return [maxSpeed, "C"]
return [cas, ""]
[docs]
class Optimization(BADA4):
"""This class implements the BADA4 optimization following the BADA4 manual.
:param AC: Aircraft object {BADA4}.
:type AC: bada4Aircraft.
"""
def __init__(self, AC):
super().__init__(AC)
self.flightEnvelope = FlightEnvelope(AC)
[docs]
def CCI(self, theta, delta, cI):
"""
Computes the cost index coefficient for given flight conditions.
:param cI: Cost index in kilograms per minute [kg min^-1].
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:type cI: float
:type delta: float
:type theta: float
:returns: Cost index coefficient [-].
:rtype: float
"""
if self.AC.BADAVersion == "4.2":
return (
(cI / 60.0)
* self.AC.LHV
/ (self.AC.MTOW * delta * const.g * const.a_0 * sqrt(theta))
)
elif self.AC.BADAVersion == "4.3" or self.AC.BADAVersion == "DUMMY":
return (
(cI / 60.0)
* self.AC.LHV
/ (
self.AC.MTOW
* pow(delta, self.AC.p_delta)
* const.g
* const.a_0
* pow(theta, self.AC.p_theta)
)
)
[docs]
def CW(self, mass, delta):
"""
Computes the weight coefficient at a given mass and pressure.
:param mass: Aircraft mass in kilograms [kg].
:param delta: Normalized pressure [-].
:type mass: float
:type delta: float
:returns: Weight coefficient [-].
:rtype: float
The weight coefficient is used to represent the aircraft's weight relative to its maximum takeoff weight (MTOW)
under given atmospheric conditions.
"""
return mass * const.g / (self.AC.MTOW * delta * const.g)
[docs]
def SR(self, M, CF):
"""
Computes the specific range (SR) for given flight conditions.
:param M: Mach ground speed [-].
:param CF: Fuel coefficient [-].
:type M: float
:type CF: float
:returns: Specific range in nautical miles per kilogram [NM kg^-1].
:rtype: float
Specific range is a measure of the distance that can be flown per unit of fuel mass.
It is calculated as the ratio of Mach number to fuel flow coefficient.
"""
return M / CF
[docs]
def econMach(self, theta, delta, mass, DeltaTemp, cI, wS):
"""
Computes the economic Mach number for a given flight condition and cost index.
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param mass: Aircraft mass in kilograms [kg].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:param cI: Cost index in kilograms per minute [kg min^-1].
:param wS: Longitudinal wind speed in meters per second [m/s].
:type delta: float
:type theta: float
:type mass: float
:type DeltaTemp: float
:type cI: float
:type wS: float
:returns: Maximum Range Cruise (MRC) Mach number [-].
:rtype: float
This function calculates the economic Mach number by balancing the fuel consumption and time costs,
taking into account the aircraft’s flight conditions and the provided cost index (cI).
It iterates through possible Mach numbers within the aircraft's flight envelope (bounded by Mmin and Mmax),
computing the drag, thrust, and fuel flow, then selects the Mach that maximizes the cost efficiency.
- **Mmin**: Lower bound of Mach speed, limited by buffet constraints.
- **Mmax**: Upper bound of Mach speed, limited by buffet constraints.
- **SR**: Specific Range, used to calculate the most efficient cruise speed in terms of fuel consumption.
"""
# clean configuration during CR
HLid = 0
LG = "LGUP"
ccI = self.CCI(cI=cI, delta=delta, theta=theta)
Mws = atm.tas2Mach(v=wS, theta=theta)
# min/max M speed limitation
Mmin = self.flightEnvelope.minMbuffet(
theta=theta, delta=delta, mass=mass, HLid=HLid, LG=LG
)
Mmax = self.flightEnvelope.maxMbuffet(
delta=delta, mass=mass, HLid=HLid, LG=LG
)
epsilon = 0.001
M_list = np.arange(Mmin, Mmax + epsilon, epsilon)
M_econ = []
cost_econ = []
for M in M_list:
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
Thrust = Drag
ThrustMax = self.Thrust(
rating="MCRZ",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
# max Thrust limitation
if Thrust > ThrustMax:
continue
CT = self.CT(Thrust=Thrust, delta=delta)
CF = self.CF(
CT=CT, delta=delta, theta=theta, M=M, DeltaTemp=DeltaTemp
)
# maximize the cost function
cost = self.SR(M=M + Mws, CF=ccI + CF)
M_econ.append(M)
cost_econ.append(cost)
if not cost_econ:
return float("Nan")
econM = M_econ[cost_econ.index(max(cost_econ))]
return proper_round(econM, 10)
# def f(M):
# CL = self.CL(M=M[0], delta=delta, mass=mass)
# CD = self.CD(M=M[0], CL=CL, HLid=HLid, LG=LG)
# Drag = self.D(M=M[0], delta=delta, CD=CD)
# CT = self.CT(Thrust=Drag, delta=delta)
# CF = self.CF(CT=CT, delta=delta, theta=theta, M=M[0], DeltaTemp=DeltaTemp)
# maximize the cost function -> to minimize, change the sign to -1 what was originally a maximization
# cost = - (self.SR(M=M[0]+Mws, CF=ccI+CF))
# return cost
# bnds = Bounds([Mmin],[Mmax])
# ThrustMax - Thrust >= 0
# cons = ({'type': 'ineq','fun': lambda M: self.Thrust(rating='MCRZ', delta=delta, theta=theta, M=M[0], DeltaTemp=DeltaTemp) - self.D(M=M[0], delta=delta, CD=self.CD(M=M[0], CL=self.CL(M=M[0], delta=delta, mass=mass), HLid=HLid, LG=LG))})
# econ = minimize(f, np.array([Mmin]), method='SLSQP', bounds=bnds, constraints=cons)
# return float(econ.x)
[docs]
def MRC(self, theta, delta, mass, DeltaTemp, wS):
"""
Computes the Mach number representing Maximum Range Cruise (MRC) for the given flight conditions.
:param theta: Normalized air temperature [-].
:param delta: Normalized air pressure [-].
:param mass: Aircraft mass in kilograms [kg].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:param wS: Longitudinal wind speed in meters per second [m/s].
:type theta: float
:type delta: float
:type mass: float
:type DeltaTemp: float
:type wS: float
:returns: Maximum Range Cruise (MRC) Mach number [-].
:rtype: float
This function calculates the Mach number that corresponds to the Maximum Range Cruise (MRC),
which maximizes the specific range. The calculation assumes clean configuration during the cruise phase.
It uses the `econMach` function with a cost index (cI) of zero to find the MRC.
If the calculated MRC is invalid, it returns NaN.
"""
mrcM = self.econMach(
cI=0.0,
theta=theta,
delta=delta,
mass=mass,
DeltaTemp=DeltaTemp,
wS=wS,
)
if isnan(mrcM):
return float("Nan")
return mrcM
[docs]
def LRC(self, theta, delta, mass, DeltaTemp, wS):
"""
Computes the Mach number representing Long Range Cruise (LRC) for the given flight conditions.
:param theta: Normalized air temperature [-].
:param delta: Normalized air pressure [-].
:param mass: Aircraft mass in kilograms [kg].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:param wS: Longitudinal wind speed in meters per second [m/s].
:type theta: float
:type delta: float
:type mass: float
:type DeltaTemp: float
:type wS: float
:returns: Long Range Cruise (LRC) Mach number [-].
:rtype: float
The Long Range Cruise (LRC) is defined as the speed where fuel efficiency is 99% of the Maximum Range Cruise (MRC).
This function calculates the LRC based on the MRC and iterates through possible Mach numbers to find the one that
minimizes the difference between the specific range (SR) at LRC and 99% of the SR at MRC.
If no valid LRC is found, it returns NaN.
"""
Mws = atm.tas2Mach(v=wS, theta=theta)
MRC = self.MRC(
theta=theta, delta=delta, mass=mass, DeltaTemp=DeltaTemp, wS=wS
)
if isnan(MRC):
return float("Nan")
# clean configuration during CR
HLid = 0
LG = "LGUP"
CL = self.CL(M=MRC, delta=delta, mass=mass)
CD = self.CD(M=MRC, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=MRC, delta=delta, CD=CD)
CT = self.CT(Thrust=Drag, delta=delta)
CF = self.CF(
CT=CT, delta=delta, theta=theta, M=MRC, DeltaTemp=DeltaTemp
)
SR_LRC = 0.99 * self.SR(M=MRC + Mws, CF=CF)
# min/max M speed limitation
Mmax = self.flightEnvelope.maxMbuffet(
delta=delta, mass=mass, HLid=HLid, LG=LG
)
# LRC > MRC
epsilon = 0.001
M_list = np.arange(MRC, Mmax + epsilon, epsilon)
M_LRC = []
cost_LRC = []
for M in M_list:
CL = self.CL(M=M, delta=delta, mass=mass)
CL_max = self.CLmax(M=M, HLid=HLid, LG=LG)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
Thrust = Drag
ThrustMax = self.Thrust(
rating="MCRZ",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
# max Thrust limitation
if Thrust > ThrustMax:
continue
if CL > CL_max:
continue
CT = self.CT(Thrust=Thrust, delta=delta)
CF = self.CF(
CT=CT, delta=delta, theta=theta, M=M, DeltaTemp=DeltaTemp
)
# specific range for LRC (definition)
SR = self.SR(M=M + Mws, CF=CF)
# minimize the cost function
cost_LRC.append(sqrt((SR - SR_LRC) ** 2))
M_LRC.append(M)
lrcM = M_LRC[cost_LRC.index(min(cost_LRC))]
return lrcM
# def f(M):
# CL = self.CL(delta=delta, mass=mass, M=M[0])
# CD = self.CD(M=M[0], CL=CL, HLid=HLid, LG=LG)
# Drag = self.D(M=M[0], delta=delta, CD=CD)
# CT = self.CT(Thrust=Drag, delta=delta)
# CF = self.CF(CT=CT, delta=delta, theta=theta, M=M[0], DeltaTemp=DeltaTemp)
# SR = self.SR(M=M[0]+Mws, CF=CF)
# return sqrt((SR - SR_LRC)**2)
# bnds = Bounds([MRC],[Mmax])
# ThrustMax - Thrust >= 0
# cons = ({'type': 'ineq','fun': lambda M: self.Thrust(rating='MCRZ', delta=delta, theta=theta, M=M[0], DeltaTemp=DeltaTemp) - self.D(M=M[0], delta=delta, CD=self.CD(M=M[0], CL=self.CL(M=M[0], delta=delta, mass=mass), HLid=HLid, LG=LG))})
# lrc = minimize(f, np.array([0.1]), method='SLSQP', bounds=bnds, constraints=cons)
# return float(lrc.x)
# minimum = float(fmin(f, x0=np.array([MRC]), disp=False))
# return minimum
[docs]
def MEC(self, theta, delta, mass, DeltaTemp, wS):
"""
Computes the Mach number representing Maximum Endurance Cruise (MEC) for given flight conditions.
:param delta: Normalized pressure [-].
:param theta: Normalized temperature [-].
:param mass: Aircraft mass in kilograms [kg].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:param wS: Longitudinal wind speed in meters per second [m/s].
:type delta: float
:type theta: float
:type mass: float
:type DeltaTemp: float
:type wS: float
:returns: Maximum Endurance Cruise (MEC) in Mach [-].
:rtype: float
The Maximum Endurance Cruise (MEC) Mach is the Mach number that minimizes the fuel consumption rate,
maximizing the endurance of the flight. This function iterates over a range of possible Mach numbers
within the flight envelope and returns the Mach number with the lowest fuel coefficient (CF).
The calculation is subject to thrust limitations, and the function ensures that the resulting Mach
does not exceed the maximum thrust available.
"""
# clean configuration during CR
HLid = 0
LG = "LGUP"
Mws = atm.tas2Mach(v=wS, theta=theta)
# min/max M speed limitation
Mmin = self.flightEnvelope.minMbuffet(
theta=theta, delta=delta, mass=mass, HLid=HLid, LG=LG
)
Mmax = self.flightEnvelope.maxMbuffet(
delta=delta, mass=mass, HLid=HLid, LG=LG
)
epsilon = 0.001
M_list = np.arange(Mmin, Mmax + epsilon, epsilon)
M_mec = []
CF_mec = []
for M in M_list:
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
Thrust = Drag
ThrustMax = self.Thrust(
rating="MCRZ",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
# max Thrust limitation
if Thrust > ThrustMax:
continue
CT = self.CT(Thrust=Thrust, delta=delta)
CF = self.CF(
CT=CT, delta=delta, theta=theta, M=M, DeltaTemp=DeltaTemp
)
# minimize the cost function
CF_mec.append(CF)
M_mec.append(M)
if not CF_mec:
return float("Nan")
mecM = M_mec[CF_mec.index(min(CF_mec))]
return proper_round(mecM, 10)
# def f(M):
# CL = self.CL(M=M[0], delta=delta, mass=mass)
# CD = self.CD(M=M[0], CL=CL, HLid=HLid, LG=LG)
# Drag = self.D(M=M[0], delta=delta, CD=CD)
# CT = self.CT(Thrust=Drag, delta=delta)
# CF = self.CF(CT=CT, delta=delta, theta=theta, M=M[0], DeltaTemp=DeltaTemp)
# return CF
# bnds = Bounds([Mmin],[Mmax + 1e-8])
# ThrustMax - Thrust >= 0
# cons = ({'type': 'ineq','fun': lambda M: self.Thrust(rating='MCRZ', delta=delta, theta=theta, M=M[0], DeltaTemp=DeltaTemp) - self.D(M=M[0], delta=delta, CD=self.CD(M=M[0], CL=self.CL(M=M[0], delta=delta, mass=mass), HLid=HLid, LG=LG))})
# mecM = minimize(f, np.array([Mmin]), method='SLSQP', bounds=bnds, constraints=cons)
# return float(mecM.x)
[docs]
def optAltitude(self, M, mass, DeltaTemp):
"""
Computes the optimum altitude for a given flight condition and Mach number.
:param M: Mach number [-].
:param mass: Aircraft mass in kilograms [kg].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type M: float
:type mass: float
:type DeltaTemp: float
:returns: Optimum altitude in feet [ft].
:rtype: float
The optimum altitude is the altitude where the aircraft achieves the maximum efficiency for a given
Mach number and mass, subject to thrust and buffet limitations. The function iterates over a range
of altitudes and returns the one with the best fuel efficiency.
The function also ensures that the optimum altitude is bounded by a minimum of 2000 feet.
"""
# clean configuration during CR
HLid = 0
LG = "LGUP"
Hmax = self.flightEnvelope.maxAltitude(
HLid=HLid, LG=LG, M=M, DeltaTemp=DeltaTemp, mass=mass, nz=1.0
)
epsilon = 100
H_list = np.arange(0, Hmax, epsilon)
H_opt = []
cost_opt = []
for H in H_list:
theta = atm.theta(H, DeltaTemp=DeltaTemp)
delta = atm.delta(H, DeltaTemp=DeltaTemp)
# min/max M speed limitation
Mmin = self.flightEnvelope.minMbuffet(
theta=theta, delta=delta, mass=mass, HLid=HLid, LG=LG
)
Mmax = self.flightEnvelope.maxMbuffet(
delta=delta, mass=mass, HLid=HLid, LG=LG
)
if M < Mmin or M > Mmax:
continue
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
Thrust = Drag
ThrustMax = self.Thrust(
rating="MCRZ",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
# max Thrust limitation
if Thrust > ThrustMax:
continue
CT = self.CT(Thrust=Thrust, delta=delta)
CF = self.CF(
CT=CT, delta=delta, theta=theta, M=M, DeltaTemp=DeltaTemp
)
ff = self.ff(
CT=CT, delta=delta, theta=theta, M=M, DeltaTemp=DeltaTemp
)
a = atm.aSound(theta=theta)
# maximize the cost function
cost = (M * a) / ff
# cost = CL/CD
H_opt.append(H)
cost_opt.append(cost)
if not cost_opt:
return float("Nan")
optH = conv.m2ft(H_opt[cost_opt.index(max(cost_opt))])
# bound the optimum altitude at 2000ft
if optH < 2000.0:
return 2000.0
return proper_round(optH, 10)
# def f(H):
# theta = atm.theta(h=H[0],DeltaTemp=DeltaTemp)
# delta = atm.delta(h=H[0],DeltaTemp=DeltaTemp)
# min/max M speed limitation
# Mmin = self.flightEnvelope.minMbuffet(theta=theta, delta=delta, mass=mass, HLid=HLid, LG=LG)
# Mmax = self.flightEnvelope.maxMbuffet(delta=delta, mass=mass, HLid=HLid, LG=LG)
# if M < Mmin or M > Mmax:
# return float('Inf')
# CL = self.CL(M=M, delta=delta, mass=mass)
# CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
# Drag = self.D(M=M, delta=delta, CD=CD)
# ThrustMax = self.Thrust(rating='MCRZ', delta=delta, theta=theta, M=M, DeltaTemp=DeltaTemp)
# max Thrust limitation
# if Thrust > ThrustMax:
# return float('Inf')
# CT = self.CT(Thrust=Drag, delta=delta)
# ff = self.ff(CT=CT, delta=delta, theta=theta, M=M, DeltaTemp=DeltaTemp)
# a = atm.aSound(theta=theta)
# maximize the cost function
# return -((M*a) / ff)
# optAlt = conv.m2ft(float(fminbound(f, x1=np.array([0]), x2=np.array([Hmax]), disp=False)))
# bound the optimum altitude at 2000ft
# print(optAlt)
# if optAlt < 1e-3:
# return float('Nan')
# if optAlt < 2000:
# return 2000
# return optAlt
[docs]
def getOPTParam(self, optParam, var_1, var_2=None):
"""
Returns the value of an optimization (OPT) parameter from a BADA4 OPT file for various flight conditions.
The OPT file contains values for parameters like Long Range Cruise (LRC), Maximum Endurance Cruise (MEC),
Maximum Range Cruise (MRC), ECON speed, or OPTALT altitude, and this function interpolates the appropriate
value based on input variables.
Note:
The array used in this function is expected to be sorted, as per the design of OPT files.
:param optParam: Name of the optimization parameter file {LRC, MEC, MRC, ECON, OPTALT}.
:param var_1: First optimizing factor, typically a flight condition like mass or altitude.
:param var_2: (Optional) Second optimizing factor if the parameter depends on more than one factor.
:type optParam: str
:type var_1: float
:type var_2: float, optional
:returns: Interpolated optimization parameter value or NaN if not found.
:rtype: float
"""
filename = os.path.join(
self.AC.filePath,
self.AC.acName,
optParam + ".OPT",
)
def findNearest(value, array):
"""
Returns the indices of the nearest value(s) in the array.
If the value is lower/higher than the lowest/highest value in the array, only one index is returned.
If the value is between two values, two closest indices (left and right) are returned.
Note:
The array is expected to be sorted.
:param value: Value to be compared to the array elements.
:param array: Sorted array of values to search.
:type value: float
:type array: list[float]
:returns: A list containing one or two nearest indices.
:rtype: list[int]
"""
nearestIdx = list()
idx = np.searchsorted(array, value, side="left")
if idx == len(array):
nearestIdx = idx - 1
elif idx == 0 or value == array[idx]:
nearestIdx = idx
elif value < array[idx] or value > array[idx]:
nearestIdx = [idx - 1, idx]
return nearestIdx
def parseOPT(filename):
"""
Parses a BADA4 OPT file and populates variables for later interpolation.
:param filename: Path to the OPT file.
:type filename: str
"""
file = open(filename, "r")
lines = file.readlines()
self.tableTypes = lines[8].split(":")[1].strip()
self.tableDimension = lines[10].split(":")[1].strip()
self.var_1 = list()
self.var_2 = list()
self.var_3 = list()
if self.tableTypes == "2D":
self.tableDimensionColumns = int(
self.tableDimension.split("x")[0]
)
self.tableDimensionRows = int(
self.tableDimension.split("x")[1]
)
self.var_2 = [
float(i)
for i in list(
filter(
None, lines[13].split("|")[1].strip().split(" ")
)
)
]
for j in range(16, 16 + self.tableDimensionRows, 1):
self.var_1.append(float(lines[j].split("|")[0].strip()))
self.var_3.extend(
[
float(i) if i != "---" else float("nan")
for i in list(
filter(
None,
lines[j].split("|")[1].strip().split(" "),
)
)
]
)
if self.tableTypes == "1D":
self.tableDimensionColumns = int(
self.tableDimension.split("x")[1]
)
self.tableDimensionRows = int(
self.tableDimension.split("x")[0]
)
for j in range(15, 15 + self.tableDimensionRows, 1):
self.var_1.append(float(lines[j].split("|")[0].strip()))
self.var_2.append(float(lines[j].split("|")[1].strip()))
parseOPT(filename=filename)
if self.tableTypes == "2D":
if var_2 is None:
return float("NaN")
nearestIdx_1 = np.array(findNearest(var_1, self.var_1))
nearestIdx_2 = np.array(findNearest(var_2, self.var_2))
# if nearestIdx_1 & nearestIdx_2 [1] [1]
if (nearestIdx_1.size == 1) & (nearestIdx_2.size == 1):
return self.var_3[
nearestIdx_1 * (self.tableDimensionColumns) + nearestIdx_2
]
# if nearestIdx_1 & nearestIdx_2 [1] [1,2]
if (nearestIdx_1.size == 1) & (nearestIdx_2.size == 2):
varTemp_1 = self.var_3[
nearestIdx_1 * (self.tableDimensionColumns)
+ nearestIdx_2[0]
]
varTemp_2 = self.var_3[
nearestIdx_1 * (self.tableDimensionColumns)
+ nearestIdx_2[1]
]
# interpolation between the 2 found points
interpVar = np.interp(
var_2,
[self.var_2[nearestIdx_2[0]], self.var_2[nearestIdx_2[1]]],
[varTemp_1, varTemp_2],
)
return interpVar
# if nearestIdx_1 & nearestIdx_2 [1,2] [1]
if (nearestIdx_1.size == 2) & (nearestIdx_2.size == 1):
varTemp_1 = self.var_3[
nearestIdx_1[0] * (self.tableDimensionColumns)
+ nearestIdx_2
]
varTemp_2 = self.var_3[
nearestIdx_1[1] * (self.tableDimensionColumns)
+ nearestIdx_2
]
# interpolation between the 2 found points
interpVar = np.interp(
var_1,
[self.var_1[nearestIdx_1[0]], self.var_1[nearestIdx_1[1]]],
[varTemp_1, varTemp_2],
)
return interpVar
# if nearestIdx_1 & nearestIdx_2 [1,2] [1,2]
if (nearestIdx_1.size == 2) & (nearestIdx_2.size == 2):
varTemp_1 = self.var_3[
nearestIdx_1[0] * (self.tableDimensionColumns)
+ nearestIdx_2[0]
]
varTemp_2 = self.var_3[
nearestIdx_1[0] * (self.tableDimensionColumns)
+ nearestIdx_2[1]
]
varTemp_3 = self.var_3[
nearestIdx_1[1] * (self.tableDimensionColumns)
+ nearestIdx_2[0]
]
varTemp_4 = self.var_3[
nearestIdx_1[1] * (self.tableDimensionColumns)
+ nearestIdx_2[1]
]
# interpolation between the 4 found points
interpVar_1 = np.interp(
var_2,
[self.var_2[nearestIdx_2[0]], self.var_2[nearestIdx_2[1]]],
[varTemp_1, varTemp_2],
)
interpVar_2 = np.interp(
var_2,
[self.var_2[nearestIdx_2[0]], self.var_2[nearestIdx_2[1]]],
[varTemp_3, varTemp_4],
)
interpVar_3 = np.interp(
var_1,
[self.var_1[nearestIdx_1[0]], self.var_1[nearestIdx_1[1]]],
[interpVar_1, interpVar_2],
)
return interpVar_3
if self.tableTypes == "1D":
nearestIdx_1 = np.array(findNearest(var_1, self.var_1))
# if nearestIdx_1 & nearestIdx_2 [1] [1]
if nearestIdx_1.size == 1:
return self.var_2[nearestIdx_1]
if nearestIdx_1.size == 2:
varTemp_1 = self.var_2[nearestIdx_1[0]]
varTemp_2 = self.var_2[nearestIdx_1[1]]
interpVar = np.interp(
var_1,
[self.var_1[nearestIdx_1[0]], self.var_1[nearestIdx_1[1]]],
[varTemp_1, varTemp_2],
)
return interpVar
[docs]
class PTD(BADA4):
"""This class implements the PTD file creator for BADA4 aircraft following BADA4 manual.
:param AC: Aircraft object {BADA4}.
:type AC: bada4Aircraft.
"""
def __init__(self, AC):
super().__init__(AC)
self.flightEnvelope = FlightEnvelope(AC)
self.ARPM = ARPM(AC)
[docs]
def create(self, DeltaTemp, saveToPath):
"""
Creates the BADA4 Performance Table Data (PTD) file, calculating climb, cruise, and descent profiles
for different mass levels and altitudes, and saves the output to the specified path.
:param saveToPath: Directory path where the PTD file should be stored.
:param DeltaTemp: Deviation from the ISA (International Standard Atmosphere) temperature in Kelvin [K].
:type saveToPath: str
:type DeltaTemp: float
:returns: None
:rtype: None
The function generates PTD files by computing the performance for three mass levels:
- 120% of the Operating Empty Weight (OEW),
- OEW plus 70% of the difference between Maximum Take-Off Weight (MTOW) and OEW,
- and Maximum Take-Off Weight (MTOW).
For each mass level, climb, cruise, and descent performance is calculated at various altitudes
(up to the maximum operating altitude or a limit of 51,000 feet). The data is then saved to a PTD file.
"""
# 3 different mass levels [kg]
massList = [
1.2 * self.AC.OEW,
self.AC.OEW + 0.7 * (self.AC.MTOW - self.AC.OEW),
self.AC.MTOW,
]
max_alt_ft = self.AC.hmo
if max_alt_ft > 51000:
max_alt_ft = 51000
# original PTF altitude list
altitudeList = list(range(0, 2000, 500))
altitudeList.extend(range(2000, 4000, 1000))
if int(max_alt_ft) < 30000:
altitudeList.extend(range(4000, int(max_alt_ft), 2000))
altitudeList.append(max_alt_ft)
else:
altitudeList.extend(range(4000, 30000, 2000))
altitudeList.extend(range(29000, int(max_alt_ft), 2000))
altitudeList.append(max_alt_ft)
CRList = []
CLList = []
DESList = []
for mass in massList:
CLList.append(
self.PTD_climb(
mass=mass, altitudeList=altitudeList, DeltaTemp=DeltaTemp
)
)
CRList.append(
self.PTD_cruise(
mass=mass, altitudeList=altitudeList, DeltaTemp=DeltaTemp
)
)
DESList.append(
self.PTD_descent(
mass=mass, altitudeList=altitudeList, DeltaTemp=DeltaTemp
)
)
self.save2PTD(
saveToPath=saveToPath,
CLList=CLList,
CRList=CRList,
DESList=DESList,
DeltaTemp=DeltaTemp,
)
[docs]
def save2PTD(self, saveToPath, CLList, CRList, DESList, DeltaTemp):
"""
Saves climb, cruise, and descent performance data into a PTD (Performance Table Data) file.
:param saveToPath: Directory path where the PTD file will be stored.
:param CLList: List of climb performance data.
:param CRList: List of cruise performance data.
:param DESList: List of descent performance data.
:param DeltaTemp: Deviation from ISA (International Standard Atmosphere) temperature in Kelvin [K].
:type saveToPath: str
:type CLList: list
:type CRList: list
:type DESList: list
:type DeltaTemp: float
:returns: None
:rtype: None
"""
def Nan2Zero(list):
# replace NAN values by 0 for printing purposes
for n in range(len(list)):
for k in range(len(list[n])):
for m in range(len(list[n][k])):
if isinstance(list[n][k][m], float):
if isnan(list[n][k][m]):
list[n][k][m] = 0
return list
newpath = saveToPath
if not os.path.exists(newpath):
os.makedirs(newpath)
if DeltaTemp == 0.0:
ISA = ""
elif DeltaTemp > 0.0:
ISA = "+" + str(int(DeltaTemp))
elif DeltaTemp < 0.0:
ISA = str(int(DeltaTemp))
filename = saveToPath + self.AC.acName + "_ISA" + ISA + ".PTD"
file = open(filename, "w")
file.write("BADA PERFORMANCE FILE RESULTS\n")
file = open(filename, "a")
file.write(
"=============================\n=============================\n\n"
)
file.write("Low mass CLIMB\n")
file.write("==============\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# replace NAN values by 0 for printing purposes
CLList = Nan2Zero(CLList)
CRList = Nan2Zero(CRList)
DESList = Nan2Zero(DESList)
# low mass
list_mass = CLList[0]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
file.write("\n\nMedium mass CLIMB\n")
file.write("=================\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# medium mass
list_mass = CLList[1]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
file.write("\n\nHigh mass CLIMB\n")
file.write("===============\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# high mass
list_mass = CLList[2]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
file.write("\n\nLow mass DESCENT\n")
file.write("================\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# low mass
list_mass = DESList[0]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
file.write("\n\nMedium mass DESCENT\n")
file.write("===================\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# medium mass
list_mass = DESList[1]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
file.write("\n\nHigh mass DESCENT\n")
file.write("=================\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# high mass
list_mass = DESList[2]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
file.write("\n\nLow mass CRUISE\n")
file.write("===============\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# low mass
list_mass = CRList[0]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
file.write("\n\nMedium mass CRUISE\n")
file.write("==================\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# medium mass
list_mass = CRList[1]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
file.write("\n\nHigh mass CRUISE\n")
file.write("================\n\n")
file.write(
" FL T p rho a TAS CAS M mass Thrust Drag Fuel ESF ROCD gamma Conf Lim\n"
)
file.write(
"[-] [K] [Pa] [kg/m3] [m/s] [kt] [kt] [-] [kg] [N] [N] [kgm] [-] [fpm] [deg] [-] \n"
)
# high mass
list_mass = CRList[2]
for k in range(0, len(list_mass[0])):
file.write(
"%3d %7.2f %7.0f %6.3f %6.1f %7.2f %7.2f %6.3f %7.0f %7.0f %7.0f %9.2f %6.3f %6.0f %7.2f %s %s\n"
% (
list_mass[0][k],
list_mass[1][k],
list_mass[2][k],
list_mass[3][k],
list_mass[4][k],
list_mass[5][k],
list_mass[6][k],
list_mass[7][k],
list_mass[8][k],
list_mass[9][k],
list_mass[10][k],
list_mass[11][k],
list_mass[12][k],
list_mass[13][k],
list_mass[14][k],
list_mass[15][k],
list_mass[16][k],
)
)
[docs]
def PTD_climb(self, mass, altitudeList, DeltaTemp):
"""
Calculates the BADA4 PTD (Performance Table Data) for the CLIMB phase of flight.
:param mass: Aircraft mass in kilograms [kg].
:param altitudeList: List of aircraft altitudes in feet [ft].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type mass: float
:type altitudeList: list of int
:type DeltaTemp: float
:returns: List of PTD CLIMB data.
:rtype: list
"""
FL_complet = []
T_complet = []
p_complet = []
rho_complet = []
a_complet = []
TAS_complet = []
CAS_complet = []
M_complet = []
mass_complet = []
Thrust_complet = []
Drag_complet = []
ff_comlet = []
ESF_complet = []
ROCD_complet = []
gamma_complet = []
conf_complet = []
Lim_complet = []
phase = "Climb"
[Vcl1, Vcl2, Mcl] = self.flightEnvelope.getSpeedSchedule(phase=phase)
Vcl1 = min(Vcl1, conv.kt2ms(250))
crossAlt = atm.crossOver(cas=Vcl2, Mach=Mcl)
for h in altitudeList:
H_m = conv.ft2m(h) # altitude [m]
[theta, delta, sigma] = atm.atmosphereProperties(
h=H_m, DeltaTemp=DeltaTemp
)
[cas, speedUpdated] = self.ARPM.climbSpeed(
h=H_m,
mass=mass,
theta=theta,
delta=delta,
DeltaTemp=DeltaTemp,
speedSchedule_default=[Vcl1, Vcl2, Mcl],
)
tas = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
M = atm.tas2Mach(v=tas, theta=theta)
a = atm.aSound(theta=theta)
FL = h / 100
# add limitation that has been applied (if some has been applied)
limitation = speedUpdated
ff = (
self.ff(
rating="MCMB",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
* 60
)
Thrust = self.Thrust(
rating="MCMB",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
config = self.flightEnvelope.getConfig(
h=H_m,
phase=phase,
v=cas,
mass=mass,
DeltaTemp=DeltaTemp,
)
# ensure the continuity of aerodyamic configuration during Climb phase of flight
if conf_complet:
prevConf = conf_complet[-1]
if config == "TO" and (prevConf == "IC" or prevConf == "CR"):
config = prevConf
elif config == "IC" and prevConf == "CR":
config = prevConf
[HLid, LG] = self.flightEnvelope.getAeroConfig(config=config)
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
if H_m < crossAlt:
ESF = self.esf(
h=H_m, flightEvolution="constCAS", M=M, DeltaTemp=DeltaTemp
)
else:
ESF = self.esf(
h=H_m, flightEvolution="constM", M=M, DeltaTemp=DeltaTemp
)
ROCD = (
conv.m2ft(
self.ROCD(
h=H_m,
T=Thrust,
D=Drag,
v=tas,
mass=mass,
ESF=ESF,
DeltaTemp=DeltaTemp,
)
)
* 60
)
temp_const = (theta * const.temp_0) / (
theta * const.temp_0 - DeltaTemp
)
dhdt = (conv.ft2m(ROCD / 60)) * temp_const
gamma = conv.rad2deg(asin(dhdt / tas))
FL_complet.append(proper_round(FL))
T_complet.append(theta * const.temp_0)
p_complet.append(delta * const.p_0)
rho_complet.append(sigma * const.rho_0)
a_complet.append(a)
TAS_complet.append(conv.ms2kt(tas))
CAS_complet.append(conv.ms2kt(cas))
M_complet.append(M)
mass_complet.append(mass)
Thrust_complet.append(Thrust)
Drag_complet.append(Drag)
ff_comlet.append(ff)
ESF_complet.append(ESF)
ROCD_complet.append(ROCD)
gamma_complet.append(gamma)
conf_complet.append(config)
Lim_complet.append(limitation)
CLList = [
FL_complet,
T_complet,
p_complet,
rho_complet,
a_complet,
TAS_complet,
CAS_complet,
M_complet,
mass_complet,
Thrust_complet,
Drag_complet,
ff_comlet,
ESF_complet,
ROCD_complet,
gamma_complet,
conf_complet,
Lim_complet,
]
return CLList
[docs]
def PTD_descent(self, mass, altitudeList, DeltaTemp):
"""
Calculates the BADA4 PTD (Performance Table Data) for the DESCENT phase of flight.
:param mass: Aircraft mass in kilograms [kg].
:param altitudeList: List of aircraft altitudes in feet [ft].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type mass: float
:type altitudeList: list of int
:type DeltaTemp: float
:returns: List of PTD DESCENT data.
:rtype: list
"""
FL_complet = []
T_complet = []
p_complet = []
rho_complet = []
a_complet = []
TAS_complet = []
CAS_complet = []
M_complet = []
mass_complet = []
Thrust_complet = []
Drag_complet = []
ff_complet = []
ESF_complet = []
ROCD_complet = []
gamma_complet = []
conf_complet = []
Lim_complet = []
phase = "Descent"
[Vdes1, Vdes2, Mdes] = self.flightEnvelope.getSpeedSchedule(
phase=phase
)
Vdes1 = min(Vdes1, conv.kt2ms(250))
crossAlt = atm.crossOver(cas=Vdes2, Mach=Mdes)
for h in reversed(altitudeList):
H_m = conv.ft2m(h) # altitude [m]
[theta, delta, sigma] = atm.atmosphereProperties(
h=H_m, DeltaTemp=DeltaTemp
)
[cas, speedUpdated] = self.ARPM.descentSpeed(
h=H_m,
mass=mass,
theta=theta,
delta=delta,
DeltaTemp=DeltaTemp,
speedSchedule_default=[Vdes1, Vdes2, Mdes],
)
tas = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
M = atm.tas2Mach(v=tas, theta=theta)
a = atm.aSound(theta=theta)
FL = h / 100
# add limitation that has been applied (if some has been applied)
limitation = speedUpdated
Thrust = self.Thrust(
rating="LIDL",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
config = self.flightEnvelope.getConfig(
h=H_m, phase="Descent", v=cas, mass=mass, DeltaTemp=DeltaTemp
)
# ensure the continuity of aerodyamic configuration during Descent phase of flight
if conf_complet:
prevConf = conf_complet[0]
if config == "CR" and (prevConf == "AP" or prevConf == "LD"):
config = prevConf
elif config == "AP" and prevConf == "LD":
config = prevConf
[HLid, LG] = self.flightEnvelope.getAeroConfig(config=config)
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
ff = (
self.ff(
rating="LIDL",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
* 60
)
if H_m < crossAlt:
ESF = self.esf(
h=H_m, flightEvolution="constCAS", M=M, DeltaTemp=DeltaTemp
)
else:
ESF = self.esf(
h=H_m, flightEvolution="constM", M=M, DeltaTemp=DeltaTemp
)
ROCD = (
conv.m2ft(
self.ROCD(
h=H_m,
T=Thrust,
D=Drag,
v=tas,
mass=mass,
ESF=ESF,
DeltaTemp=DeltaTemp,
)
)
* 60
)
temp_const = (theta * const.temp_0) / (
theta * const.temp_0 - DeltaTemp
)
dhdt = (conv.ft2m(ROCD / 60)) * temp_const
gamma = conv.rad2deg(asin(dhdt / tas))
minSpeed = self.flightEnvelope.VMin(
config=config, mass=mass, theta=theta, delta=delta
)
[HLid, LG] = self.flightEnvelope.getAeroConfig(config=config)
maxSpeed = self.flightEnvelope.VMax(
h=h, HLid=HLid, LG=LG, theta=theta, delta=delta, mass=mass
)
# in case of AP & LD thrust is computed to fly a 3deg slope
if config == "AP" or config == "LD":
gamma = -3.0
temp_const = (theta * const.temp_0) / (
theta * const.temp_0 - DeltaTemp
)
ROCD_gamma = sin(conv.deg2rad(gamma)) * tas * (1 / temp_const)
ROCD = conv.m2ft(ROCD_gamma) * 60 # [ft/min]
n = 1.0 # aircraft.loadFactor(gamma) - use this in case of L = W * (cos(gamma))
CL = self.CL(M=M, delta=delta, mass=mass, nz=n)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
Thrust = (ROCD_gamma * mass * const.g) * temp_const / (
ESF * tas
) + Drag
CT = self.CT(Thrust=Thrust, delta=delta)
ff = (
self.ff(
CT=CT,
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
* 60
)
FL_complet = [proper_round(FL)] + FL_complet
T_complet = [theta * const.temp_0] + T_complet
p_complet = [delta * const.p_0] + p_complet
rho_complet = [sigma * const.rho_0] + rho_complet
a_complet = [a] + a_complet
TAS_complet = [conv.ms2kt(tas)] + TAS_complet
CAS_complet = [conv.ms2kt(cas)] + CAS_complet
M_complet = [M] + M_complet
mass_complet = [mass] + mass_complet
Thrust_complet = [Thrust] + Thrust_complet
Drag_complet = [Drag] + Drag_complet
ff_complet = [ff] + ff_complet
ESF_complet = [ESF] + ESF_complet
ROCD_complet = [-1 * ROCD] + ROCD_complet
gamma_complet = [gamma] + gamma_complet
conf_complet = [config] + conf_complet
Lim_complet = [limitation] + Lim_complet
DESList = [
FL_complet,
T_complet,
p_complet,
rho_complet,
a_complet,
TAS_complet,
CAS_complet,
M_complet,
mass_complet,
Thrust_complet,
Drag_complet,
ff_complet,
ESF_complet,
ROCD_complet,
gamma_complet,
conf_complet,
Lim_complet,
]
return DESList
[docs]
def PTD_cruise(self, mass, altitudeList, DeltaTemp):
"""
Calculates the BADA4 PTD (Performance Table Data) for the CRUISE phase of flight.
:param mass: Aircraft mass in kilograms [kg].
:param altitudeList: List of aircraft altitudes in feet [ft].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type mass: float
:type altitudeList: list of int
:type DeltaTemp: float
:returns: List of PTD CRUISE data.
:rtype: list
"""
FL_complet = []
T_complet = []
p_complet = []
rho_complet = []
a_complet = []
TAS_complet = []
CAS_complet = []
M_complet = []
mass_complet = []
Thrust_complet = []
Drag_complet = []
ff_comlet = []
ESF_complet = []
ROCD_complet = []
gamma_complet = []
conf_complet = []
Lim_complet = []
phase = "Cruise"
[Vcr1, Vcr2, Mcr] = self.flightEnvelope.getSpeedSchedule(phase=phase)
Vcr1 = min(Vcr1, conv.kt2ms(250))
for h in altitudeList:
H_m = conv.ft2m(h) # altitude [m]
[theta, delta, sigma] = atm.atmosphereProperties(
h=H_m, DeltaTemp=DeltaTemp
)
[cas, speedUpdated] = self.ARPM.cruiseSpeed(
h=H_m,
mass=mass,
theta=theta,
delta=delta,
speedSchedule_default=[Vcr1, Vcr2, Mcr],
)
tas = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
M = atm.tas2Mach(v=tas, theta=theta)
a = atm.aSound(theta=theta)
FL = h / 100
# add limitation that has been applied (if some has been applied)
limitation = ""
config = "CR"
HLid = 0
LG = "LGUP"
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
Thrust = Drag
ThrustMax = self.Thrust(
rating="MCRZ",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
CT = self.CT(Thrust=Thrust, delta=delta)
ff = (
self.ff(
CT=CT, delta=delta, theta=theta, M=M, DeltaTemp=DeltaTemp
)
* 60
)
if Thrust > ThrustMax:
# "(T)" - as thrust limited
limitation += "T"
limitation += speedUpdated
ESF = 0.0
ROCD = 0.0
gamma = 0.0
FL_complet.append(proper_round(FL))
T_complet.append(theta * const.temp_0)
p_complet.append(delta * const.p_0)
rho_complet.append(sigma * const.rho_0)
a_complet.append(a)
TAS_complet.append(conv.ms2kt(tas))
CAS_complet.append(conv.ms2kt(cas))
M_complet.append(M)
mass_complet.append(mass)
Thrust_complet.append(Thrust)
Drag_complet.append(Drag)
ff_comlet.append(ff)
ESF_complet.append(ESF)
ROCD_complet.append(ROCD)
gamma_complet.append(gamma)
conf_complet.append(config)
Lim_complet.append(limitation)
CRList = [
FL_complet,
T_complet,
p_complet,
rho_complet,
a_complet,
TAS_complet,
CAS_complet,
M_complet,
mass_complet,
Thrust_complet,
Drag_complet,
ff_comlet,
ESF_complet,
ROCD_complet,
gamma_complet,
conf_complet,
Lim_complet,
]
return CRList
[docs]
class PTF(BADA4):
"""This class implements the PTF file creator for BADA4 aircraft following BADA4 manual.
:param AC: Aircraft object {BADA4}.
:type AC: bada4Aircraft.
"""
def __init__(self, AC):
super().__init__(AC)
self.flightEnvelope = FlightEnvelope(AC)
self.ARPM = ARPM(AC)
[docs]
def create(self, DeltaTemp, saveToPath):
"""
Creates the BADA4 PTF and saves it to the specified directory.
:param saveToPath: Path to the directory where the PTF should be stored.
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type saveToPath: str
:type DeltaTemp: float
:returns: None
:rtype: None
"""
# 3 different mass levels [kg]
massList = [
1.2 * self.AC.OEW,
self.AC.OEW + 0.7 * (self.AC.MTOW - self.AC.OEW),
self.AC.MTOW,
]
max_alt_ft = self.AC.hmo
if max_alt_ft > 51000:
max_alt_ft = 51000
# original PTF altitude list
altitudeList = list(range(0, 2000, 500))
altitudeList.extend(range(2000, 4000, 1000))
if int(max_alt_ft) < 30000:
altitudeList.extend(range(4000, int(max_alt_ft), 2000))
altitudeList.append(max_alt_ft)
else:
altitudeList.extend(range(4000, 30000, 2000))
altitudeList.extend(range(29000, int(max_alt_ft), 2000))
altitudeList.append(max_alt_ft)
CRList = self.PTF_cruise(
massList=massList, altitudeList=altitudeList, DeltaTemp=DeltaTemp
)
CLList = self.PTF_climb(
massList=massList, altitudeList=altitudeList, DeltaTemp=DeltaTemp
)
DESList = self.PTF_descent(
massList=massList, altitudeList=altitudeList, DeltaTemp=DeltaTemp
)
self.save2PTF(
saveToPath=saveToPath,
altitudeList=altitudeList,
massList=massList,
CRList=CRList,
CLList=CLList,
DESList=DESList,
DeltaTemp=DeltaTemp,
)
[docs]
def save2PTF(
self,
saveToPath,
CLList,
CRList,
DESList,
DeltaTemp,
massList,
altitudeList,
):
"""
Saves the BADA4 performance data to a PTF format.
:param saveToPath: Path to the directory where the PTF should be stored.
:param CLList: List of PTD data for CLIMB.
:param CRList: List of PTD data for CRUISE.
:param DESList: List of PTD data for DESCENT.
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:param massList: List of aircraft mass levels [kg].
:param altitudeList: List of altitudes [ft].
:type saveToPath: str
:type CLList: list
:type CRList: list
:type DESList: list
:type DeltaTemp: float
:type massList: list(float)
:returns: None
:rtype: None
This function formats and writes the climb, cruise, and descent data for different mass levels
and altitudes into a .PTF file, adhering to the BADA4 performance file format.
"""
def Nan2Zero(list):
# replace NAN values by 0 for printing purposes
for k in range(len(list)):
for m in range(len(list[k])):
if isinstance(list[k][m], float):
if isnan(list[k][m]):
list[k][m] = 0
return list
newpath = saveToPath
if not os.path.exists(newpath):
os.makedirs(newpath)
if DeltaTemp == 0.0:
ISA = ""
elif DeltaTemp > 0.0:
ISA = "+" + str(int(DeltaTemp))
elif DeltaTemp < 0.0:
ISA = str(int(DeltaTemp))
filename = saveToPath + self.AC.acName + "_ISA" + ISA + ".PTF"
[Vcl1, Vcl2, Mcl] = self.flightEnvelope.getSpeedSchedule(phase="Climb")
[Vcr1, Vcr2, Mcr] = self.flightEnvelope.getSpeedSchedule(
phase="Cruise"
)
[Vdes1, Vdes2, Mdes] = self.flightEnvelope.getSpeedSchedule(
phase="Descent"
)
V1cl = min(250, conv.ms2kt(Vcl1))
V2cl = conv.ms2kt(Vcl2)
V1des = min(250, conv.ms2kt(Vdes1))
V2des = conv.ms2kt(Vdes2)
V1cr = min(250, conv.ms2kt(Vcr1))
V2cr = conv.ms2kt(Vcr2)
today = date.today()
d3 = today.strftime("%b %d %Y")
acModel = self.AC.model
file = open(filename, "w")
file.write(
"BADA PERFORMANCE FILE %s\n\n"
% (d3)
)
file = open(filename, "a")
file.write("AC/Type: %s\n\n" % (acModel))
file.write(
" Speeds: CAS(LO/HI) Mach Mass Levels [kg] Temperature: ISA%s\n"
% (ISA)
)
file.write(
" climb - %3d/%3d %4.3f low - %.0f\n"
% (V1cl, V2cl, Mcl, massList[0])
)
file.write(
" cruise - %3d/%3d %4.3f nominal - %-5.0f Max Alt. [ft]:%7d\n"
% (V1cr, V2cr, Mcr, massList[1], altitudeList[-1])
)
file.write(
" descent - %3d/%3d %4.3f high - %0.f\n"
% (V1des, V2des, Mdes, massList[2])
)
file.write(
"======================================================================================================\n"
)
file.write(
" FL | CRUISE | CLIMB | DESCENT \n"
)
file.write(
" | TAS fuel | TAS ROCD fuel | TAS ROCD fuel \n"
)
file.write(
" | [kts] [kg/min] | [kts] [fpm] [kg/min] | [kts] [fpm] [kg/min]\n"
)
file.write(
" | lo nom hi | lo nom hi nom | lo nom hi nom \n"
)
file.write(
"======================================================================================================\n"
)
# replace NAN values by 0 for printing purposes
CLList = Nan2Zero(CLList)
DESList = Nan2Zero(DESList)
for k in range(0, len(altitudeList)):
FL = proper_round(altitudeList[k] / 100)
file.write(
"%3.0f | %s %s %s %s | %3.0f %5.0f %5.0f %5.0f %5.1f | %3.0f %5.0f %5.0f %5.0f %5.1f\n"
% (
FL,
CRList[0][k],
CRList[1][k],
CRList[2][k],
CRList[3][k],
CLList[0][k],
CLList[1][k],
CLList[2][k],
CLList[3][k],
CLList[4][k],
DESList[0][k],
DESList[1][k],
DESList[2][k],
DESList[3][k],
DESList[4][k],
)
)
file.write(
" | | | \n"
)
file.write(
"======================================================================================================\n"
)
[docs]
def PTF_cruise(self, massList, altitudeList, DeltaTemp):
"""
Calculates the BADA4 PTF for the CRUISE phase of flight.
:param massList: List of aircraft mass levels in kilograms [kg].
:param altitudeList: List of aircraft altitudes in feet [ft].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type massList: list
:type altitudeList: list of int
:type DeltaTemp: float
:returns: List of PTF CRUISE data.
:rtype: list
"""
TAS_CR_complet = []
FF_CR_LO_complet = []
FF_CR_NOM_complet = []
FF_CR_HI_complet = []
massNominal = massList[1]
[Vcr1, Vcr2, Mcr] = self.flightEnvelope.getSpeedSchedule(
phase="Cruise"
)
Vcr1 = min(Vcr1, conv.kt2ms(250))
for h in altitudeList:
H_m = conv.ft2m(h) # altitude [m]
[theta, delta, sigma] = atm.atmosphereProperties(
h=H_m, DeltaTemp=DeltaTemp
)
[cas, speedUpdated] = self.ARPM.cruiseSpeed(
h=H_m,
mass=massNominal,
theta=theta,
delta=delta,
speedSchedule_default=[Vcr1, Vcr2, Mcr],
DeltaTemp=DeltaTemp,
)
tas_nominal = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
FL = h / 100
ff = []
for mass in massList:
[cas, speedUpdated] = self.ARPM.cruiseSpeed(
h=H_m,
mass=mass,
theta=theta,
delta=delta,
speedSchedule_default=[Vcr1, Vcr2, Mcr],
DeltaTemp=DeltaTemp,
)
tas = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
M = atm.tas2Mach(v=tas, theta=theta)
config = "CR"
HLid = 0
LG = "LGUP"
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
Thrust = Drag
ThrustMax = self.Thrust(
rating="MCRZ",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
CL = self.flightEnvelope.CL(
delta=delta, mass=mass, M=M, nz=1.2
)
CL_max = self.flightEnvelope.CLmax(M=M, HLid=HLid, LG=LG)
epsilon = 0.01
if Thrust > ThrustMax:
# "(T)" - as thrust limited
ff.append("(T)")
elif CL > (CL_max + epsilon):
# "(B)" - as buffet limited
ff.append("(B)")
else:
CT = self.CT(Thrust=Thrust, delta=delta)
ff.append(
self.ff(
CT=CT,
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
* 60
)
TAS_CR_complet.append(f"{conv.ms2kt(tas_nominal):3.0f}")
if isinstance(ff[0], str):
FF_CR_LO_complet.append(" " + ff[0] + " ")
else:
FF_CR_LO_complet.append(f"{ff[0]:5.1f}")
if isinstance(ff[1], str):
FF_CR_NOM_complet.append(" " + ff[1] + " ")
else:
FF_CR_NOM_complet.append(f"{ff[1]:5.1f}")
if isinstance(ff[2], str):
FF_CR_HI_complet.append(" " + ff[2] + " ")
else:
FF_CR_HI_complet.append(f"{ff[2]:5.1f}")
CRList = [
TAS_CR_complet,
FF_CR_LO_complet,
FF_CR_NOM_complet,
FF_CR_HI_complet,
]
return CRList
[docs]
def PTF_climb(self, massList, altitudeList, DeltaTemp):
"""
Calculates the BADA4 PTF for the CLIMB phase of flight.
:param massList: List of aircraft mass levels in kilograms [kg].
:param altitudeList: List of aircraft altitudes in feet [ft].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type massList: list
:type altitudeList: list of int
:type DeltaTemp: float
:returns: List of PTF CLIMB data.
:rtype: list
"""
TAS_CL_complet = []
ROCD_CL_LO_complet = []
ROCD_CL_NOM_complet = []
ROCD_CL_HI_complet = []
FF_CL_NOM_complet = []
conf_LO_complet = []
conf_NOM_complet = []
conf_HI_complet = []
conf_complet = {}
for mass in massList:
conf_complet[str(mass)] = []
massNominal = massList[1]
[Vcl1, Vcl2, Mcl] = self.flightEnvelope.getSpeedSchedule(phase="Climb")
Vcl1 = min(Vcl1, conv.kt2ms(250))
crossAlt = atm.crossOver(cas=Vcl2, Mach=Mcl)
for h in altitudeList:
H_m = conv.ft2m(h) # altitude [m]
[theta, delta, sigma] = atm.atmosphereProperties(
h=H_m, DeltaTemp=DeltaTemp
)
[cas, speedUpdated] = self.ARPM.climbSpeed(
h=H_m,
mass=massNominal,
theta=theta,
delta=delta,
DeltaTemp=DeltaTemp,
speedSchedule_default=[Vcl1, Vcl2, Mcl],
)
tas_nominal = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
FL = h / 100
M_nominal = atm.tas2Mach(v=tas_nominal, theta=theta)
ff_nominal = (
self.ff(
rating="MCMB",
delta=delta,
theta=theta,
M=M_nominal,
DeltaTemp=DeltaTemp,
)
* 60
)
ROC = []
for mass in massList:
[cas, speedUpdated] = self.ARPM.climbSpeed(
h=H_m,
mass=mass,
theta=theta,
delta=delta,
DeltaTemp=DeltaTemp,
speedSchedule_default=[Vcl1, Vcl2, Mcl],
)
tas = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
M = atm.tas2Mach(v=tas, theta=theta)
Thrust = self.Thrust(
rating="MCMB",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
config = self.flightEnvelope.getConfig(
h=H_m,
phase="Climb",
v=cas,
mass=mass,
DeltaTemp=DeltaTemp,
)
# ensure the continuity of aerodyamic configuration during CLimb phase of flight
if conf_complet[str(mass)]:
prevConf = conf_complet[str(mass)][-1]
if config == "TO" and (
prevConf == "IC" or prevConf == "CR"
):
config = prevConf
elif config == "IC" and prevConf == "CR":
config = prevConf
[HLid, LG] = self.flightEnvelope.getAeroConfig(config=config)
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
if H_m < crossAlt:
ESF = self.esf(
h=H_m,
flightEvolution="constCAS",
M=M,
DeltaTemp=DeltaTemp,
)
else:
ESF = self.esf(
h=H_m,
flightEvolution="constM",
M=M,
DeltaTemp=DeltaTemp,
)
ROC_val = (
conv.m2ft(
self.ROCD(
h=H_m,
T=Thrust,
D=Drag,
v=tas,
mass=mass,
ESF=ESF,
DeltaTemp=DeltaTemp,
)
)
* 60
)
if ROC_val < 0:
ROC_val = float("Nan")
ROC.append(ROC_val)
conf_complet[str(mass)].append(config)
TAS_CL_complet.append(conv.ms2kt(tas_nominal))
ROCD_CL_LO_complet.append(ROC[0])
ROCD_CL_NOM_complet.append(ROC[1])
ROCD_CL_HI_complet.append(ROC[2])
FF_CL_NOM_complet.append(ff_nominal)
CLList = [
TAS_CL_complet,
ROCD_CL_LO_complet,
ROCD_CL_NOM_complet,
ROCD_CL_HI_complet,
FF_CL_NOM_complet,
]
return CLList
[docs]
def PTF_descent(self, massList, altitudeList, DeltaTemp):
"""
Calculates the BADA4 PTF for the DESCENT phase of flight.
:param massList: List of aircraft mass levels in kilograms [kg].
:param altitudeList: List of aircraft altitudes in feet [ft].
:param DeltaTemp: Deviation from ISA temperature in Kelvin [K].
:type massList: list
:type altitudeList: list of int
:type DeltaTemp: float
:returns: List of PTF DESCENT data.
:rtype: list
"""
TAS_DES_complet = []
ROCD_DES_LO_complet = []
ROCD_DES_NOM_complet = []
ROCD_DES_HI_complet = []
FF_DES_NOM_complet = []
conf_complet = {}
for mass in massList:
conf_complet[str(mass)] = []
massNominal = massList[1]
[Vdes1, Vdes2, Mdes] = self.flightEnvelope.getSpeedSchedule(
phase="Descent"
)
Vdes1 = min(Vdes1, conv.kt2ms(250))
crossAlt = atm.crossOver(cas=Vdes2, Mach=Mdes)
# for h in altitudeList:
for h in reversed(altitudeList):
H_m = conv.ft2m(h) # altitude [m]
[theta, delta, sigma] = atm.atmosphereProperties(
h=H_m, DeltaTemp=DeltaTemp
)
[cas, speedUpdated] = self.ARPM.descentSpeed(
h=H_m,
mass=massNominal,
theta=theta,
delta=delta,
DeltaTemp=DeltaTemp,
speedSchedule_default=[Vdes1, Vdes2, Mdes],
)
tas_nominal = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
FL = h / 100
M_nominal = atm.tas2Mach(v=tas_nominal, theta=theta)
ff_nominal = (
self.ff(
rating="LIDL",
delta=delta,
theta=theta,
M=M_nominal,
DeltaTemp=DeltaTemp,
)
* 60
)
ROD = []
ff_gamma_list = []
for mass in massList:
[cas, speedUpdated] = self.ARPM.descentSpeed(
h=H_m,
mass=mass,
theta=theta,
delta=delta,
DeltaTemp=DeltaTemp,
speedSchedule_default=[Vdes1, Vdes2, Mdes],
)
tas = atm.cas2Tas(cas=cas, delta=delta, sigma=sigma)
M = atm.tas2Mach(v=tas, theta=theta)
Thrust = self.Thrust(
rating="LIDL",
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
config = self.flightEnvelope.getConfig(
h=H_m,
phase="Descent",
v=cas,
mass=mass,
DeltaTemp=DeltaTemp,
)
# ensure the continuity of aerodyamic configuration during Descent phase of flight
if conf_complet[str(mass)]:
prevConf = conf_complet[str(mass)][0]
if config == "CR" and (
prevConf == "AP" or prevConf == "LD"
):
config = prevConf
elif config == "AP" and prevConf == "LD":
config = prevConf
[HLid, LG] = self.flightEnvelope.getAeroConfig(config=config)
CL = self.CL(M=M, delta=delta, mass=mass)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
if H_m < crossAlt:
ESF = self.esf(
h=H_m,
flightEvolution="constCAS",
M=M,
DeltaTemp=DeltaTemp,
)
else:
ESF = self.esf(
h=H_m,
flightEvolution="constM",
M=M,
DeltaTemp=DeltaTemp,
)
ROCD = (
conv.m2ft(
self.ROCD(
h=H_m,
T=Thrust,
D=Drag,
v=tas,
mass=mass,
ESF=ESF,
DeltaTemp=DeltaTemp,
)
)
* 60
)
# in case of AP & LD thrust is computed to fly a 3? slope
if config == "AP" or config == "LD":
gamma = -3.0
temp_const = (theta * const.temp_0) / (
theta * const.temp_0 - DeltaTemp
)
ROCD_gamma = (
sin(conv.deg2rad(gamma)) * tas * (1 / temp_const)
)
ROCD = conv.m2ft(ROCD_gamma) * 60 # [ft/min]
n = 1.0 # aircraft.loadFactor(gamma) - use this in case of L = W * (cos(gamma))
CL = self.CL(M=M, delta=delta, mass=mass, nz=n)
CD = self.CD(M=M, CL=CL, HLid=HLid, LG=LG)
Drag = self.D(M=M, delta=delta, CD=CD)
Thrust = (ROCD_gamma * mass * const.g) * temp_const / (
ESF * tas
) + Drag
CT = self.CT(Thrust=Thrust, delta=delta)
ff_gamma = (
self.ff(
CT=CT,
delta=delta,
theta=theta,
M=M,
DeltaTemp=DeltaTemp,
)
* 60
)
ff_gamma_list.append(ff_gamma)
else:
ff_gamma_list.append(float("Nan"))
ROD.append(ROCD)
conf_complet[str(mass)] = [config] + conf_complet[str(mass)]
if not isnan(ff_gamma_list[1]):
ff_nominal = ff_gamma_list[1]
TAS_DES_complet = [conv.ms2kt(tas_nominal)] + TAS_DES_complet
ROCD_DES_LO_complet = [-1 * ROD[0]] + ROCD_DES_LO_complet
ROCD_DES_NOM_complet = [-1 * ROD[1]] + ROCD_DES_NOM_complet
ROCD_DES_HI_complet = [-1 * ROD[2]] + ROCD_DES_HI_complet
FF_DES_NOM_complet = [ff_nominal] + FF_DES_NOM_complet
DESList = [
TAS_DES_complet,
ROCD_DES_LO_complet,
ROCD_DES_NOM_complet,
ROCD_DES_HI_complet,
FF_DES_NOM_complet,
]
return DESList
[docs]
class Bada4Aircraft(BADA4):
"""
This class encapsulates the BADA4 performance model for an aircraft, extending the BADA4 base class.
:param badaVersion: The version of the BADA4 model being used.
:param acName: The ICAO designation or name of the aircraft.
:param filePath: (Optional) Path to the BADA4 XML file. If not provided, a default path is used.
:param allData: (Optional) Dataframe containing pre-loaded aircraft data, typically used to
initialize the aircraft parameters without needing to parse XML files.
:type badaVersion: str
:type acName: str
:type filePath: str, optional
:type allData: pd.DataFrame, optional
This class initializes the aircraft's performance model using data from a dataframe or by
reading from XML files in the BADA4 format.
"""
def __init__(self, badaVersion, acName, filePath=None, allData=None):
"""
Initializes the BADA4Aircraft class by loading aircraft-specific data.
- If `allData` is provided and contains the aircraft's information, it will be used to
initialize various parameters such as engine type, mass, thrust settings, and performance
data.
- If the aircraft is not found in `allData`, the class will search for the corresponding
BADA4 XML file or synonym file (if applicable) in the specified or default file path.
- Once the aircraft data is found, the class initializes various performance modules such
as the flight envelope, aerodynamic model, and performance optimizations.
:param badaVersion: Version of the BADA4 model (e.g., "4.2", "4.3").
:param acName: ICAO aircraft designation or model name.
:param filePath: (Optional) Custom file path to load the aircraft data. If not provided,
a default directory is used.
:param allData: (Optional) Dataframe containing pre-loaded aircraft data for initialization.
"""
super().__init__(self)
self.BADAFamily = BadaFamily(BADA4=True)
self.BADAFamilyName = "BADA4"
self.BADAVersion = badaVersion
self.acName = acName
if filePath == None:
self.filePath = configuration.getBadaVersionPath(
badaFamily="BADA4", badaVersion=badaVersion
)
else:
self.filePath = filePath
# check if the aircraft is in the allData dataframe data
if allData is not None and acName in allData["acName"].values:
filtered_df = allData[allData["acName"] == acName]
self.model = configuration.safe_get(filtered_df, "model", None)
self.engineType = configuration.safe_get(
filtered_df, "engineType", None
)
self.engines = configuration.safe_get(filtered_df, "engines", None)
self.WTC = configuration.safe_get(filtered_df, "WTC", None)
self.ICAO = configuration.safe_get(filtered_df, "ICAO", None)
self.MREF = configuration.safe_get(filtered_df, "MREF", None)
self.WREF = configuration.safe_get(filtered_df, "WREF", None)
self.LHV = configuration.safe_get(filtered_df, "LHV", None)
self.n_eng = configuration.safe_get(filtered_df, "n_eng", None)
self.rho = configuration.safe_get(filtered_df, "rho", None)
self.TFA = configuration.safe_get(filtered_df, "TFA", None)
self.p_delta = configuration.safe_get(filtered_df, "p_delta", None)
self.p_theta = configuration.safe_get(filtered_df, "p_theta", None)
self.kink = configuration.safe_get(filtered_df, "kink", None)
self.b = configuration.safe_get(filtered_df, "b", None)
self.c = configuration.safe_get(filtered_df, "c", None)
self.max_power = configuration.safe_get(
filtered_df, "max_power", None
)
self.p = configuration.safe_get(filtered_df, "p", None)
self.a = configuration.safe_get(filtered_df, "a", None)
self.f = configuration.safe_get(filtered_df, "f", None)
self.ti = configuration.safe_get(filtered_df, "ti", None)
self.fi = configuration.safe_get(filtered_df, "fi", None)
self.throttle = configuration.safe_get(
filtered_df, "throttle", None
)
self.prop_dia = configuration.safe_get(
filtered_df, "prop_dia", None
)
self.max_eff = configuration.safe_get(filtered_df, "max_eff", None)
self.Hd_turbo = configuration.safe_get(
filtered_df, "Hd_turbo", None
)
self.CPSFC = configuration.safe_get(filtered_df, "CPSFC", None)
self.P = configuration.safe_get(filtered_df, "P", None)
self.S = configuration.safe_get(filtered_df, "S", None)
self.HLPosition = configuration.safe_get(
filtered_df, "HLPosition", None
)
self.configName = configuration.safe_get(
filtered_df, "configName", None
)
self.VFE = configuration.safe_get(filtered_df, "VFE", None)
self.d = configuration.safe_get(filtered_df, "d", None)
self.CL_max = configuration.safe_get(filtered_df, "CL_max", None)
self.bf = configuration.safe_get(filtered_df, "bf", None)
self.HLids = configuration.safe_get(filtered_df, "HLids", None)
self.CL_clean = configuration.safe_get(
filtered_df, "CL_clean", None
)
self.M_max = configuration.safe_get(filtered_df, "M_max", None)
self.scalar = configuration.safe_get(filtered_df, "scalar", None)
self.Mmin = configuration.safe_get(filtered_df, "Mmin", None)
self.Mmax = configuration.safe_get(filtered_df, "Mmax", None)
self.CL_Mach0 = configuration.safe_get(
filtered_df, "CL_Mach0", None
)
self.MTOW = configuration.safe_get(filtered_df, "MTOW", None)
self.OEW = configuration.safe_get(filtered_df, "OEW", None)
self.MFL = configuration.safe_get(filtered_df, "MFL", None)
self.MTW = configuration.safe_get(filtered_df, "MTW", None)
self.MZFW = configuration.safe_get(filtered_df, "MZFW", None)
self.MPL = configuration.safe_get(filtered_df, "MPL", None)
self.MLW = configuration.safe_get(filtered_df, "MLW", None)
self.hmo = configuration.safe_get(filtered_df, "hmo", None)
self.mfa = configuration.safe_get(filtered_df, "mfa", None)
self.MMO = configuration.safe_get(filtered_df, "MMO", None)
self.MLE = configuration.safe_get(filtered_df, "MLE", None)
self.VLE = configuration.safe_get(filtered_df, "VLE", None)
self.VMO = configuration.safe_get(filtered_df, "VMO", None)
self.span = configuration.safe_get(filtered_df, "span", None)
self.length = configuration.safe_get(filtered_df, "length", None)
self.aeroConfig = configuration.safe_get(
filtered_df, "aeroConfig", None
)
self.speedSchedule = configuration.safe_get(
filtered_df, "speedSchedule", None
)
# GPF data (temporary)
self.CVminTO = configuration.safe_get(filtered_df, "CVminTO", None)
self.CVmin = configuration.safe_get(filtered_df, "CVmin", None)
self.HmaxPhase = configuration.safe_get(
filtered_df, "HmaxPhase", None
)
self.V_des = configuration.safe_get(filtered_df, "V_des", None)
self.V_cl = configuration.safe_get(filtered_df, "V_cl", None)
self.flightEnvelope = FlightEnvelope(self)
self.ARPM = ARPM(self)
self.OPT = Optimization(self)
self.PTD = PTD(self)
self.PTF = PTF(self)
else:
self.ACModelAvailable = False
self.synonymFileAvailable = False
self.ACinSynonymFile = False
# check if SYNONYM file exist - since for BADA4 this is not a standard procedure (yet)
synonymFile = os.path.join(
self.filePath,
"aircraft_model_default.xml",
)
if os.path.isfile(synonymFile):
self.synonymFileAvailable = True
# if SYNONYM exist - look for synonym based on defined acName
self.SearchedACName = Parser.parseMappingFile(
filePath=self.filePath,
acName=acName,
)
# if cannot find - look for full name (in sub folder names) based on acName (may not be ICAO designator)
if self.SearchedACName == None:
self.SearchedACName = acName
else:
self.ACinSynonymFile = True
else:
# if it doesn't exist - look for full name (in sub folder names) based on acName (may not be ICAO designator)
self.SearchedACName = acName
acXmlFile = (
os.path.join(
self.filePath,
self.SearchedACName,
self.SearchedACName,
)
+ ".xml"
)
# look for either found synonym or original full BADA4 model name designator
if self.SearchedACName is not None:
if os.path.isfile(acXmlFile):
self.ACModelAvailable = True
XMLDataFrame = Parser.parseXML(
filePath=self.filePath,
acName=self.SearchedACName,
)
GPFDataframe = Parser.parseGPF(filePath=self.filePath)
combined_df = Parser.combineXML_GPF(
XMLDataFrame, GPFDataframe
)
self.acName = self.SearchedACName
self.model = configuration.safe_get(
combined_df, "model", None
)
self.engineType = configuration.safe_get(
combined_df, "engineType", None
)
self.engines = configuration.safe_get(
combined_df, "engines", None
)
self.WTC = configuration.safe_get(combined_df, "WTC", None)
self.ICAO = configuration.safe_get(
combined_df, "ICAO", None
)
self.MREF = configuration.safe_get(
combined_df, "MREF", None
)
self.WREF = configuration.safe_get(
combined_df, "WREF", None
)
self.LHV = configuration.safe_get(combined_df, "LHV", None)
self.n_eng = configuration.safe_get(
combined_df, "n_eng", None
)
self.rho = configuration.safe_get(combined_df, "rho", None)
self.TFA = configuration.safe_get(combined_df, "TFA", None)
self.p_delta = configuration.safe_get(
combined_df, "p_delta", None
)
self.p_theta = configuration.safe_get(
combined_df, "p_theta", None
)
self.kink = configuration.safe_get(
combined_df, "kink", None
)
self.b = configuration.safe_get(combined_df, "b", None)
self.c = configuration.safe_get(combined_df, "c", None)
self.max_power = configuration.safe_get(
combined_df, "max_power", None
)
self.p = configuration.safe_get(combined_df, "p", None)
self.a = configuration.safe_get(combined_df, "a", None)
self.f = configuration.safe_get(combined_df, "f", None)
self.ti = configuration.safe_get(combined_df, "ti", None)
self.fi = configuration.safe_get(combined_df, "fi", None)
self.throttle = configuration.safe_get(
combined_df, "throttle", None
)
self.prop_dia = configuration.safe_get(
combined_df, "prop_dia", None
)
self.max_eff = configuration.safe_get(
combined_df, "max_eff", None
)
self.Hd_turbo = configuration.safe_get(
combined_df, "Hd_turbo", None
)
self.CPSFC = configuration.safe_get(
combined_df, "CPSFC", None
)
self.P = configuration.safe_get(combined_df, "P", None)
self.S = configuration.safe_get(combined_df, "S", None)
self.HLPosition = configuration.safe_get(
combined_df, "HLPosition", None
)
self.configName = configuration.safe_get(
combined_df, "configName", None
)
self.VFE = configuration.safe_get(combined_df, "VFE", None)
self.d = configuration.safe_get(combined_df, "d", None)
self.CL_max = configuration.safe_get(
combined_df, "CL_max", None
)
self.bf = configuration.safe_get(combined_df, "bf", None)
self.HLids = configuration.safe_get(
combined_df, "HLids", None
)
self.CL_clean = configuration.safe_get(
combined_df, "CL_clean", None
)
self.M_max = configuration.safe_get(
combined_df, "M_max", None
)
self.scalar = configuration.safe_get(
combined_df, "scalar", None
)
self.Mmin = configuration.safe_get(
combined_df, "Mmin", None
)
self.Mmax = configuration.safe_get(
combined_df, "Mmax", None
)
self.CL_Mach0 = configuration.safe_get(
combined_df, "CL_Mach0", None
)
self.MTOW = configuration.safe_get(
combined_df, "MTOW", None
)
self.OEW = configuration.safe_get(combined_df, "OEW", None)
self.MFL = configuration.safe_get(combined_df, "MFL", None)
self.MTW = configuration.safe_get(combined_df, "MTW", None)
self.MZFW = configuration.safe_get(
combined_df, "MZFW", None
)
self.MPL = configuration.safe_get(combined_df, "MPL", None)
self.MLW = configuration.safe_get(combined_df, "MLW", None)
self.hmo = configuration.safe_get(combined_df, "hmo", None)
self.mfa = configuration.safe_get(combined_df, "mfa", None)
self.MMO = configuration.safe_get(combined_df, "MMO", None)
self.MLE = configuration.safe_get(combined_df, "MLE", None)
self.VLE = configuration.safe_get(combined_df, "VLE", None)
self.VMO = configuration.safe_get(combined_df, "VMO", None)
self.span = configuration.safe_get(
combined_df, "span", None
)
self.length = configuration.safe_get(
combined_df, "length", None
)
self.aeroConfig = configuration.safe_get(
combined_df, "aeroConfig", None
)
self.speedSchedule = configuration.safe_get(
combined_df, "speedSchedule", None
)
# GPF data (temporary)
self.CVminTO = configuration.safe_get(
combined_df, "CVminTO", None
)
self.CVmin = configuration.safe_get(
combined_df, "CVmin", None
)
self.HmaxPhase = configuration.safe_get(
combined_df, "HmaxPhase", None
)
self.V_des = configuration.safe_get(
combined_df, "V_des", None
)
self.V_cl = configuration.safe_get(
combined_df, "V_cl", None
)
# BADA4.__init__(self, AC_parsed)
self.flightEnvelope = FlightEnvelope(self)
self.ARPM = ARPM(self)
self.OPT = Optimization(self)
self.PTD = PTD(self)
self.PTF = PTF(self)
else:
# AC name cannot be found
raise ValueError(
acName + " Cannot be found at path " + self.filePath
)
def __str__(self):
return f"(BADA4, AC_name: {self.acName}, searched_AC_name: {self.SearchedACName}, model_ICAO: {self.ICAO}, ID: {id(self.AC)})"