Coverage for docs/source/examples/sp_to_gp_sweep.py: 100%
Shortcuts on this page
r m x toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
Shortcuts on this page
r m x toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1"example of an SP that turns into a GP"
2import numpy as np
3from gpkit import Model, Variable, SignomialsEnabled, units
4from gpkit.constraints.tight import Tight
7def SimPleAC():
8 "Creates SimpleAC model"
9 # Env. constants
10 g = Variable("g", 9.81, "m/s^2", "gravitational acceleration")
11 mu = Variable("\\mu", 1.775e-5, "kg/m/s", "viscosity of air")
12 rho = Variable("\\rho", 1.23, "kg/m^3", "density of air")
13 rho_f = Variable("\\rho_f", 817, "kg/m^3", "density of fuel")
15 # Non-dimensional constants
16 C_Lmax = Variable("C_{L,max}", 1.6, "-", "max CL with flaps down")
17 e = Variable("e", 0.92, "-", "Oswald efficiency factor")
18 k = Variable("k", 1.17, "-", "form factor")
19 N_ult = Variable("N_{ult}", 3.3, "-", "ultimate load factor")
20 S_wetratio = Variable("(\\frac{S}{S_{wet}})", 2.075, "-",
21 "wetted area ratio")
22 tau = Variable("\\tau", 0.12, "-", "airfoil thickness to chord ratio")
23 W_W_coeff1 = Variable("W_{W_{coeff1}}", 2e-5, "1/m",
24 "wing weight coefficent 1") # 12e-5 originally
25 W_W_coeff2 = Variable("W_{W_{coeff2}}", 60, "Pa",
26 "wing weight coefficent 2")
28 # Dimensional constants
29 Range = Variable("Range", 3000, "km", "aircraft range")
30 TSFC = Variable("TSFC", 0.6, "1/hr", "thrust specific fuel consumption")
31 V_min = Variable("V_{min}", 25, "m/s", "takeoff speed")
32 W_0 = Variable("W_0", 6250, "N", "aircraft weight excluding wing")
34 # Free Variables
35 LoD = Variable("L/D", "-", "lift-to-drag ratio")
36 D = Variable("D", "N", "total drag force")
37 V = Variable("V", "m/s", "cruising speed")
38 W = Variable("W", "N", "total aircraft weight")
39 Re = Variable("Re", "-", "Reynold's number")
40 CDA0 = Variable("(CDA0)", "m^2", "fuselage drag area") # 0.035 originally
41 C_D = Variable("C_D", "-", "drag coefficient")
42 C_L = Variable("C_L", "-", "lift coefficient of wing")
43 C_f = Variable("C_f", "-", "skin friction coefficient")
44 W_f = Variable("W_f", "N", "fuel weight")
45 V_f = Variable("V_f", "m^3", "fuel volume")
46 V_f_avail = Variable("V_{f_{avail}}", "m^3", "fuel volume available")
47 T_flight = Variable("T_{flight}", "hr", "flight time")
49 # Free variables (fixed for performance eval.)
50 A = Variable("A", "-", "aspect ratio")
51 S = Variable("S", "m^2", "total wing area")
52 W_w = Variable("W_w", "N", "wing weight")
53 W_w_strc = Variable("W_w_strc", "N", "wing structural weight")
54 W_w_surf = Variable("W_w_surf", "N", "wing skin weight")
55 V_f_wing = Variable("V_f_wing", "m^3", "fuel volume in the wing")
56 V_f_fuse = Variable("V_f_fuse", "m^3", "fuel volume in the fuselage")
58 objective = W_f
60 constraints = []
62 # Weight and lift model
63 constraints += [
64 W >= W_0 + W_w + W_f,
65 W_0 + W_w + 0.5 * W_f <= 0.5 * rho * S * C_L * V ** 2,
66 W <= 0.5 * rho * S * C_Lmax * V_min ** 2,
67 T_flight >= Range / V,
68 LoD == C_L/C_D]
70 # Thrust and drag model
71 C_D_fuse = CDA0 / S
72 C_D_wpar = k * C_f * S_wetratio
73 C_D_ind = C_L ** 2 / (np.pi * A * e)
74 constraints += [
75 W_f >= TSFC * T_flight * D,
76 D >= 0.5 * rho * S * C_D * V ** 2,
77 C_D >= C_D_fuse + C_D_wpar + C_D_ind,
78 V_f_fuse <= 10*units("m")*CDA0,
79 Re <= (rho / mu) * V * (S / A) ** 0.5,
80 C_f >= 0.074 / Re ** 0.2]
82 # Fuel volume model
83 with SignomialsEnabled():
84 constraints += [
85 V_f == W_f / g / rho_f,
86 # linear with b and tau, quadratic with chord
87 V_f_wing**2 <= 0.0009*S**3/A*tau**2,
88 V_f_avail <= V_f_wing + V_f_fuse, # [SP]
89 Tight([V_f_avail >= V_f])]
91 # Wing weight model
92 constraints += [
93 W_w_surf >= W_W_coeff2 * S,
94 W_w_strc**2 >= W_W_coeff1**2/tau**2 * N_ult**2*A**3*(V_f_fuse*g*rho_f
95 + W_0)*W*S,
96 W_w >= W_w_surf + W_w_strc]
98 m = Model(objective, constraints)
100 return m
103sa = SimPleAC()
104sa.substitutions.update({"V_f_wing": ("sweep", np.linspace(0.1, 0.5, 3)),
105 "V_f_fuse": 0.5})
106sol = sa.solve(verbosity=0)
107print(sol.summary())