Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
def test():
"""
Test each different engine
"""
#dict of initial guesses
x0 = {
'W_{engine}': 1e4*units('N'),
'P_{t_0}': 1e1*units('kPa'),
'T_{t_0}': 1e3*units('K'),
'h_{t_0}': 1e6*units('J/kg'),
'P_{t_{1.8}}': 1e1*units('kPa'),
'T_{t_{1.8}}': 1e3*units('K'),
'h_{t_{1.8}}': 1e6*units('J/kg'),
'P_{t_2}': 1e1*units('kPa'),
'T_{t_2}': 1e3*units('K'),
'h_{t_2}': 1e6*units('J/kg'),
'P_{t_2.1}': 1e3*units('K'),
'T_{t_2.1}': 1e3*units('K'),
'h_{t_2.1}': 1e6*units('J/kg'),
'P_{t_{2.5}}': 1e3*units('kPa'),
'T_{t_{2.5}}': 1e3*units('K'),
'h_{t_{2.5}}': 1e6*units('J/kg'),
'P_{t_3}': 1e4*units('kPa'),
def test():
m = Mission(SimPleAC(),4)
m.substitutions.update({
'h_{cruise_m}' :5000*units('m'),
'Range_m' :3000*units('km'),
'W_{p_m}' :3000*units('N'),
'\\rho_{p_m}' :1500*units('kg/m^3'),
'C_m' :120*units('1/hr'),
'V_{min_m}' :35*units('m/s'),
'T/O factor_m' :2,
})
m.cost = m['W_{f_m}']*units('1/N') + m['C_m']*m['t_m']
sol = m.localsolve(verbosity=0)
def __init__(self):
CL = Variable('C_L', '-', 'Lift coefficient')
CLmax = Variable('C_{L_{max}}', '-', 'Max lift coefficient')
CD = Variable('C_D', '-', 'Drag coefficient')
D = Variable('D', 'N', 'Total aircraft drag (cruise)')
Dfuse = Variable('D_{fuse}', 'N', 'Fuselage drag')
Dht = Variable('D_{ht}', 'N', 'Horizontal tail drag')
Dvt = Variable('D_{vt}', 'N', 'Vertical tail drag')
Dwing = Variable('D_{wing}', 'N', 'Wing drag')
LD = Variable('\\frac{L}{D}', '-', 'Lift/drag ratio')
Lh = Variable('L_h', 'N', 'Horizontal tail downforce')
Lw = Variable('L_w', 'N', 'Wing lift')
M = Variable('M', '-', 'Cruise Mach number')
R = Variable('Range', 'nautical_miles', 'Range')
Sw = Variable('S_w', 'm**2', 'Wing reference area')
Te = Variable('T_e', 'N', 'Engine thrust at takeoff')
TSFC = Variable('TSFC', 'lb/lbf/hr',
'Thrust specific fuel consumption')
V = Variable('V_{\\infty}', 'm/s', 'Cruise velocity')
VTO = Variable('V_{TO}', 'm/s', 'Takeoff speed')
W = Variable('W', 'N', 'Total aircraft weight')
Weng = Variable('W_{eng}', 'N', 'Engine weight')
Wfuel = Variable('W_{fuel}', 'N', 'Fuel weight')
Wfuse = Variable('W_{fuse}', 'N', 'Fuselage weight')
Wht = Variable('W_{ht}', 'N', 'Horizontal tail weight')
def SimPleAC():
"Creates SimpleAC model"
# Env. constants
g = Variable("g", 9.81, "m/s^2", "gravitational acceleration")
mu = Variable("\\mu", 1.775e-5, "kg/m/s", "viscosity of air")
rho = Variable("\\rho", 1.23, "kg/m^3", "density of air")
rho_f = Variable("\\rho_f", 817, "kg/m^3", "density of fuel")
# Non-dimensional constants
C_Lmax = Variable("C_{L,max}", 1.6, "-", "max CL with flaps down")
e = Variable("e", 0.92, "-", "Oswald efficiency factor")
k = Variable("k", 1.17, "-", "form factor")
N_ult = Variable("N_{ult}", 3.3, "-", "ultimate load factor")
S_wetratio = Variable("(\\frac{S}{S_{wet}})", 2.075, "-",
"wetted area ratio")
tau = Variable("\\tau", 0.12, "-", "airfoil thickness to chord ratio")
W_W_coeff1 = Variable("W_{W_{coeff1}}", 2e-5, "1/m",
"wing weight coefficent 1") # 12e-5 originally
W_W_coeff2 = Variable("W_{W_{coeff2}}", 60, "Pa",
"wing weight coefficent 2")
# Dimensional constants
Range = Variable("Range", 3000, "km", "aircraft range")
TSFC = Variable("TSFC", 0.6, "1/hr", "thrust specific fuel consumption")
V_min = Variable("V_{min}", 25, "m/s", "takeoff speed")
W_0 = Variable("W_0", 6250, "N", "aircraft weight excluding wing")
# Free Variables
LoD = Variable("L/D", "-", "lift-to-drag ratio")
D = Variable("D", "N", "total drag force")
V = Variable("V", "m/s", "cruising speed")
W = Variable("W", "N", "total aircraft weight")
aircraft = Aircraft(Nclimb, Ncruise, enginestate)
# vectorize
with Vectorize(Nclimb):
climb = ClimbSegment(aircraft, Nclimb)
with Vectorize(Ncruise):
cruise = CruiseSegment(aircraft, Nclimb)
statelinking = StateLinking(climb.state, cruise.state, enginestate, Nclimb, Ncruise)
# declare new variables
W_ftotal = Variable('W_{f_{total}}', 'lbf', 'Total Fuel Weight')
W_fclimb = Variable('W_{f_{climb}}', 'lbf',
'Fuel Weight Burned in Climb')
W_fcruise = Variable('W_{f_{cruise}}', 'lbf',
'Fuel Weight Burned in Cruise')
W_total = Variable('W_{total}', 'lbf', 'Total Aircraft Weight')
## W_dry = Variable('W_{dry}', 'lbf', 'Dry Aircraft Weight')
CruiseAlt = Variable('CruiseAlt', 'ft', 'Cruise Altitude [feet]')
ReqRng = Variable('ReqRng', 'nautical_miles', 'Required Cruise Range')
h = climb.state['h']
hftClimb = climb.state['hft']
dhft = climb.climbP['dhft']
hftCruise = cruise.state['hft']
# make overall constraints
constraints = []
constraints.extend([
# weight constraints
#pressure ratios
pihpt = Variable('\pi_{HPT}', '-', 'HPT Pressure Ratio')
pilpt = Variable('\pi_{LPT}', '-', 'LPT Pressure Ratio')
#turbine efficiences
etaht = Variable('\eta_{ht}', 0.9, '-', 'Polytropic Efficiency of HPT')
etalt = Variable('\eta_{lt}', 0.9, '-', 'Polytropic Efficiency of LPT')
#flow faction f
f = Variable('f', '-', 'Fuel Air Mass Flow Fraction')
#BPR
alpha = Variable('alpha', '-', 'By Pass Ratio')
#relavent pressure ratio
pitn = Variable('\pi_{tn}', '-', 'Turbine Nozzle Pressure Ratio')
with SignomialsEnabled():
constraints = [
#HPT shafter power balance
#SIGNOMIAL
SignomialEquality((1+f)*(ht41-ht45),ht3 - ht25), #B.161
#LPT shaft power balance
#SIGNOMIAL
SignomialEquality((1+f)*(ht49 - ht45),-((ht25-ht18)+alpha*(ht21 - ht2))), #B.165
#HPT Exit states (station 4.5)
Pt45 == pihpt * Pt41,
pihpt == (Tt45/Tt41)**(4.60517), #turbine efficiency is 0.9
ht45 == Cpt1 * Tt45,
def __init__(self, ht, fuse, wing, fuseP, wingP, state, **kwargs):
self.ht = ht
self.fuse = fuse
self.wing = wing
self.fuseP = fuseP
self.wingP = wingP
#variables
alpha = Variable('\\alpha', '-', 'Horizontal tail angle of attack')
D = Variable('D_{ht}', 'N', 'Horizontal tail drag')
Lh = Variable('L_h', 'N', 'Horizontal tail downforce')
SM = VectorVariable(2, 'S.M.', '-', 'Stability margin')
Rec = Variable('Re_{c_h}', '-',
'Cruise Reynolds number (Horizontal tail)')
CLah = Variable('C_{L_{ah}}', '-', 'Lift curve slope (htail)')
CLah0 = Variable('C_{L_{ah_0}}', '-',
'Isolated lift curve slope (htail)')
CLaw = Variable('C_{L_{aw}}', '-', 'Lift curve slope (wing)')
CLh = Variable('C_{L_h}', '-', 'Lift coefficient (htail)')
dxlead = Variable('\\Delta x_{{lead}_h}', 'm',
'Distance from CG to horizontal tail leading edge')
dxtrail = Variable('\\Delta x_{{trail}_h}', 'm',
'Distance from CG to horizontal tail trailing edge')
dxw = VectorVariable(2, '\\Delta x_w', 'm',
'Distance from aerodynamic centre to CG')
def setup(self, aircraft, state, **kwargs):
#submodels
self.aircraft = aircraft
self.aircraftP = AircraftP(aircraft, state)
self.wingP = self.aircraftP.wingP
self.fuseP = self.aircraftP.fuseP
#variable definitions
theta = Variable('\\theta', '-', 'Aircraft Climb Angle')
excessP = Variable('excessP', 'W', 'Excess Power During Climb')
RC = Variable('RC', 'feet/min', 'Rate of Climb/Decent')
dhft = Variable('dhft', 'feet', 'Change in Altitude Per Climb Segment [feet]')
RngClimb = Variable('RngClimb', 'nautical_miles', 'Down Range Covered in Each Climb Segment')
#constraints
constraints = []
constraints.extend([
RC == excessP/self.aircraftP['W_{avg}'],
RC >= 500*units('ft/min'),
#make the small angle approximation and compute theta
theta * state['V'] == RC,
dhft == self.aircraftP['tmin'] * RC,
k = Variable("k", 1.2, "-", "form factor")
e = Variable("e", 0.95, "-", "Oswald efficiency factor")
mu = Variable("\\mu", 1.78e-5, "kg/m/s", "viscosity of air")
pi = Variable("\\pi", np.pi, "-", "half of the circle constant")
rho = Variable("\\rho", 1.23, "kg/m^3", "density of air")
tau = Variable("\\tau", 0.12, "-", "airfoil thickness to chord ratio")
N_ult = Variable("N_{ult}", 3.8, "-", "ultimate load factor")
V_min = Variable("V_{min}", 22, "m/s", "takeoff speed")
C_Lmax = Variable("C_{L,max}", 1.5, "-", "max CL with flaps down")
S_wetratio = Variable("(\\frac{S}{S_{wet}})", 2.05, "-",
"wetted area ratio")
W_W_coeff1 = Variable("W_{W_{coeff1}}", 8.71e-5, "1/m",
"Wing Weight Coefficent 1")
W_W_coeff2 = Variable("W_{W_{coeff2}}", 45.24, "Pa",
"Wing Weight Coefficent 2")
CDA0 = Variable("(CDA0)", 0.031, "m^2", "fuselage drag area")
W_0 = Variable("W_0", 4940.0, "N", "aircraft weight excluding wing")
# Free Variables
D = Variable("D", "N", "total drag force")
A = Variable("A", "-", "aspect ratio")
S = Variable("S", "m^2", "total wing area")
V = Variable("V", "m/s", "cruising speed")
W = Variable("W", "N", "total aircraft weight")
Re = Variable("Re", "-", "Reynold's number")
C_D = Variable("C_D", "-", "Drag coefficient of wing")
C_L = Variable("C_L", "-", "Lift coefficent of wing")
C_f = Variable("C_f", "-", "skin friction coefficient")
W_w = Variable("W_w", "N", "wing weight")
constraints = []
taucone = Variable('\\tau_{cone}', 'N/m^2', 'Shear stress in cone')
taufloor = Variable('\\tau_{floor}', 'N/m^2',
'Max allowable shear web stress')
tcone = Variable('t_{cone}', 'm', 'Cone thickness')
tshell = Variable('t_{shell}', 'm', 'Shell thickness')
tskin = Variable('t_{skin}', 'm', 'Skin thickness')
waisle = Variable('w_{aisle}', 'm', 'Aisle width')
wfloor = Variable('w_{floor}', 'm', 'Floor width')
wfuse = Variable('w_{fuse}', 'm', 'Fuselage width')
wseat = Variable('w_{seat}', 'm', 'Seat width')
wsys = Variable('w_{sys}', 'm',
'Width between cabin and skin for systems')
xCGfu = Variable('x_{CG_{fu}}', 'm', 'x-location of fuselage CG')
xVbulk = Variable('xVbulk', 'm^4', 'Volume moment of bulkhead')
xVcyl = Variable('xVcyl', 'm^4', 'Volume moment of cylinder')
xVnose = Variable('xVnose', 'm^4', 'Volume moment of nose')
xWapu = Variable('xWapu', 'N*m', 'Moment of APU')
xWcone = Variable('xWcone', 'N*m', 'Moment of cone')
xWfix = Variable('xWfix', 'N*m', 'Moment of fixed weights')
xWfloor = Variable('xWfloor', 'N*m', 'Moment of floor weight')
xWfuse = Variable('xWfuse', 'N*m', 'Fuselage moment')
xWinsul = Variable('xWinsul', 'N*m', 'Moment of insulation material')
xWpadd = Variable('xWpadd', 'N*m', 'Moment of misc weights')
xWseat = Variable('xWseat', 'N*m', 'Moment of seats')
xWshell = Variable('xWshell', 'N*m', 'Mass moment of shell')
xWskin = Variable('xWskin', 'N*m', 'Mass moment of skin')
xWwindow = Variable('xWwindow', 'N*m', 'Mass moment of windows')
x_upswp = Variable('x_{up}', 'm', 'Fuselage upsweep point')
xapu = Variable('xapu', 'ft', 'x-location of APU')
xconend = Variable('xconend', 'm', 'x-location of cone end')
xfix = Variable('xfix', 'm', 'x-location of fixed weight')
xshell1 = Variable('x_{shell1}', 'm', 'Start of cylinder section')