Skip to content

Commit

Permalink
Updated tests; added MIL constraint to multi-element propeller
Browse files Browse the repository at this point in the history
  • Loading branch information
courtin committed Apr 17, 2018
1 parent a043a68 commit da89baf
Show file tree
Hide file tree
Showing 4 changed files with 97 additions and 58 deletions.
13 changes: 6 additions & 7 deletions gpkitmodels/GP/aircraft/motor/motor.py
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
"Electric motor model "
from gpkit import Model, parse_variables
from gpkit.constraints.tight import Tight as TCS
from gpkitmodels.GP.aircraft.prop.propeller import Propeller
from gpkitmodels.GP.aircraft.prop.propeller import Propeller, Actuator_Propeller
from gpkitmodels import g

class MotorPerf(Model):
Expand Down Expand Up @@ -37,7 +37,7 @@ class Motor(Model):
Variables
---------
Qstar 1 [kg/(N*m)] motor specific torque
Qstar .5 [kg/(N*m)] motor specific torque
W [lbf] motor weight
Qmax [N*m] motor max. torque
V_max 300 [V] motor max voltage
Expand All @@ -64,8 +64,8 @@ class PropulsorPerf(Model):

def setup(self, static, state):
exec parse_variables(PropulsorPerf.__doc__)
self.prop = parent.prop.flight_model(static.prop, state)
self.motor = parent.motor.flight_model(static.motor, state)
self.prop = static.prop.flight_model(static.prop, state)
self.motor = static.motor.flight_model(static.motor, state)

self.components = [self.prop, self.motor]

Expand Down Expand Up @@ -103,14 +103,13 @@ class Actuator_Propulsor(Model):
W [lbf] propulsor weight
"""
flight_model = Propulsor_Performance
flight_model = PropulsorPerf

def setup(self):
exec parse_variables(Propulsor.__doc__)

#self.prop = Actuator_Propeller()
self.prop = Actuator_Propeller()
self.motor = ElecMotor()
self.motor = Motor()

components = [self.prop, self.motor]

Expand Down
6 changes: 3 additions & 3 deletions gpkitmodels/GP/aircraft/motor/motor_test.py
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
from gpkit import Model, parse_variables, SignomialsEnabled, SignomialEquality, units
from motor import Actuator_Propulsor,Propulsor, ElecMotor, ElecMotor_Performance
from gpkitmodels.GP.aircraft.prop.propeller import Actuator_Propeller
from motor import Actuator_Propulsor,Propulsor, Motor, MotorPerf
from gpkitmodels.GP.aircraft.prop.propeller import Propeller
from gpkitmodels.GP.aircraft.wing.wing_test import FlightState

class Propulsor_Test(Model):
Expand All @@ -27,7 +27,7 @@ def setup(self):
pp = p.flight_model(p,fs)
pp.substitutions[pp.prop.T] = 100
#pp.substitutions[pp.prop.AR_b] = 15
self.cost = pp.motor.Pelec/(1000*units('W')) + p.W/(1000*units('lbf')) + 1./pp.prop.eta
self.cost = pp.motor.Pelec/(1000*units('W')) + p.W/(1000*units('lbf'))

return fs,p,pp
def actuator_propulsor_test():
Expand Down
25 changes: 13 additions & 12 deletions gpkitmodels/GP/aircraft/prop/prop_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,10 +11,11 @@ def simpleprop_test():
Propeller.flight_model = ActuatorProp
p = Propeller()
pp = p.flight_model(p, fs)
m = Model(1/pp.eta + pp.Q/(100.*units("N*m")) + p.W/(100.*units("lbf")),
m = Model(1/pp.eta + p.W/(100.*units("lbf"))+ pp.Q/(100.*units("N*m")),
[fs, p, pp])
m.substitutions.update({"rho": 1.225, "V": 50, "T": 100, "T_m": 40})
m.solve()
m.substitutions.update({"rho": 1.225, "V": 50, "T": 100, "omega":1000})
sol = m.solve()
#print sol.table()

def OE_eta_test():
" simple qprop test "
Expand All @@ -23,27 +24,29 @@ def OE_eta_test():
Propeller.flight_model = SimpleQProp
p = Propeller()
pp = p.flight_model(p, fs)
pp.substitutions[pp.T] = 1000
pp.substitutions[pp.T] = 100
pp.substitutions[pp.AR_b] = 12
m = Model(1/pp.eta + pp.Q/(10.*units("N*m")) + p.W/(100.*units("lbf")),
m = Model(1/pp.eta + pp.Q/(10.*units("N*m"))+ p.W/(100.*units("lbf")),
[pp, p])
m.localsolve()
#m.debug()
sol = m.localsolve()
#print sol.table()

def ME_eta_test():

fs = FlightState()
p = Multi_Element_Propeller()
p.substitutions[p.T_m] = 100
p.substitutions[p.R] = 2
#p.substitutions[p.R] = 2
pp = p.flight_model(p,fs)
pp.substitutions[pp.T] = 100

#pp.substitutions[pp.omega] = 500

#pp.substitutions[pp.omega] = 1000
pp.cost = 1./pp.eta + pp.Q/(100.*units("N*m"))
sol = pp.localsolve(iteration_limit = 400)
#sol = pp.debug()
sol = pp.localsolve(iteration_limit = 400)
print sol.table()

#def qprop_test():
Expand All @@ -57,11 +60,9 @@ def ME_eta_test():

def test():
"tests"
#eta_test()
#OE_eta_test()
OE_eta_test()
simpleprop_test()
ME_eta_test()
#simpleprop_test()

if __name__ == "__main__":
test()

111 changes: 75 additions & 36 deletions gpkitmodels/GP/aircraft/prop/propeller.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,8 +13,8 @@ class ActuatorProp(Model):
---------
T [lbf] thrust
Tc [-] coefficient of thrust
etaadd 0.7 [-] swirl and nonuniformity losses
etav 0.85 [-] viscous losses
etaadd .7 [-] swirl and nonuniformity losses
etav .85 [-] viscous losses
etai [-] inviscid losses
eta [-] overall efficiency
z1 self.helper [-] efficiency helper 1
Expand All @@ -33,7 +33,7 @@ class ActuatorProp(Model):
def helper(self, c):
return 2. - 1./c[self.etaadd]

def setup(self, static, state):
def setup(self, static, state = None, rho = None, V = None):
exec parse_variables(ActuatorProp.__doc__)

V = state.V
Expand All @@ -50,7 +50,8 @@ def setup(self, static, state):
eta >= CT*lam/CP,
omega <= omega_max,
P_shaft == Q*omega,
(M_tip*a)**2 >= (omega*R)**2 + V**2
(M_tip*a)**2 >= (omega*R)**2 + V**2,
static.T_m >= T
]
return constraints

Expand All @@ -60,7 +61,7 @@ class Blade_Element_Performance(Model):
Variables
---------
dT [lbf] thrust
eta [-] local efficiency
eta_b [-] local efficiency
dQ [N*m] torque
omega [rpm] propeller rotation rate
Wa [m/s] Axial total relative velocity
Expand All @@ -69,21 +70,28 @@ class Blade_Element_Performance(Model):
va [m/s] Axial induced velocity
vt [m/s] Tangential induced velocity
G [m^2/s] Circulation
c .1 [m] local chord
cl [-] local lift coefficient
c [m] local chord
cl [-] local lift coefficient
cd [-] local drag coefficient
B 2 [-] number of blades
cdp [-] local drag coefficient
B 3 [-] number of blades
r [m] local radius
lam_w [-] advance ratio
eps [-] blade efficiency
AR_b [-] blade aspect ratio
AR_b_max 18 [-] max blade aspect ratio
AR_b [-] blade aspect ratio
AR_b_max 50 [-] max blade aspect ratio
Ut [m/s] tangential freestreem
Re [-] blade reynolds number
f [-] intermediate tip loss variable
F .05 [-] Prandtl tip loss factor
cl_max 1.3 [-] max airfoil cl
F [-] Prandtl tip loss factor
cl_max .9 [-] max airfoil cl
dr [m] length of blade element
FF [-] form factor
t_c .1 [-] airfoil t/c
C_f [-] skin friction coefficient
M [-] Mach number
a 295 [m/s] Speed of sound at altitude
"""

def setup(self,parent, state):
Expand All @@ -101,39 +109,47 @@ def setup(self,parent, state):
#Ut == omega*r,
TCS([G == (1./2.)*Wr*c*cl]),
#TCS([f+(r/R)*B/(2*lam_w) <= (B/2.)*(1./lam_w)]), #.2 = 1-.8 = 1-r/R
#F == (2./pi)*(1.02232*f**.0458729)**(1./.1), #This is the GP fit of arccos(exp(-f))
#TCS([vt == (B*G/(4.*pi*r*F))]),
F == (2./pi)*(1.02232*f**.0458729)**(1./.1), #This is the GP fit of arccos(exp(-f))
M == Wr/a,
lam_w == (r/R)*(Wa/Wt),
va == vt*(Wt/Wa),
eps == cd/cl,
#TCS([1 >= etai*(1+va/V)+vt/Ut]),
#TCS([1 >= etap*(1+eps*Wt/Wa)+eps*Wt/Wa]),
#TCS([dQ >= rho*B*G*(Wa+eps*Wt)*r*dr]),
Wa == V,
TCS([dQ >= rho*B*G*(Wa+eps*Wt)*r*dr]),
#Wa == V,
#Wt == omega*r,
#TCS([(M_tip*a)**2 >= Wa**2 + Wt**2]),
AR_b == R/c,
AR_b <= AR_b_max,
Re == Wr*c*rho/mu,
#eta == T*V/(Q*omega),
eta_b == (V/(omega*r))*(Wt/Wa),

#TCS([cd >= cl**2/(pi*AR_b*.85) + cdp]),
#TCS([FF >= 1. + .3*t_c + (100*t_c)**.4]),
#C_f == .027/(Re**(1./7.)),
#cdp >= FF*C_f*(2.+.4*t_c),

XfoilFit(fd, cd, [cl,Re], name="polar"),

cl <= cl_max
]
with SignomialsEnabled():
constraints += [#TCS([Wa<=V + va]),
constraints += [TCS([Wa>=V + va]),
#SignomialEquality(Wa, V+va),
SignomialEquality(Wr**2,(Wa**2+Wt**2)),
#TCS([Wr**2 >= (Wa**2+Wt**2)]),
SignomialEquality(Wa,V+va),
SignomialEquality(Wt,omega*r-vt),
#TCS([Wt>=omega*r-vt]),
#TCS([dT <= rho*B*G*(Wt-eps*Wa)*dr]),
#TCS([Wt<=omega*r-vt]),
TCS([dT <= rho*B*G*(Wt-eps*Wa)*dr]),
SignomialEquality(vt**2*F**2*(1.+(4.*lam_w*R/(pi*B*r))**2),(B*G/(4.*pi*r))**2),
#TCS([vt**2*F**2*(1.+(4.*lam_w*R/(pi*B*r))**2)<=(B*G/(4.*pi*r))**2]),

#TCS([f+(r/R)*B/(2*lam_w) <= (B/2.)*(1./lam_w)]),
SignomialEquality(f+(r/R)*B/(2*lam_w),(B/2.)*(1./lam_w))
#SignomialEquality(1,etap*(1+eps*Wt/Wa)+eps*Wa/Wt),
#SignomialEquality(1,etai*(1+va/V)+vt/Ut),
SignomialEquality(dQ,rho*B*G*(Wa+eps*Wt)*r*dr),
SignomialEquality(dT,rho*B*G*(Wt-eps*Wa)*dr),
#SignomialEquality(dQ,rho*B*G*(Wa+eps*Wt)*r*dr),
#SignomialEquality(dT,rho*B*G*(Wt-eps*Wa)*dr),

]
return constraints, state
Expand All @@ -145,36 +161,35 @@ class Multi_Element_Propeller_Performance(Model):
Variables
---------
T [lbf] overall thrust
eta [-] overall efficiency
Q [N*m] overall torque
omega [rpm] propeller rotation rate
Mtip .5 [-] tip mach number
"""
def setup(self,parent, state, N = 5):
def setup(self,parent, state, N = 4):
#exec parse_variables(Multi_Element_Propeller_Performance.__doc__)
T = self.T = Variable('T', 'lbf', 'Overall thrust')
Q = self.Q = Variable('Q', 'N*m', 'Overall torque')
eta = self.eta = Variable('eta', '-', 'Overall efficiency')
omega = self.omega = Variable('omega', 'rpm', 'rotation rate')
omega_max = Variable('omega_max',10000, 'rpm', 'rotation rate')

Mtip = Variable('Mtip', .5, '-', 'Tip Mach')
with Vectorize(N):
blade = Blade_Element_Performance(parent, state)


constraints = [blade.dr == parent.R/(N),
blade.omega == omega,
blade["r"][0] == parent.R/(2*N)]
blade["r"][0] == parent.R/(2.*N)]


with SignomialsEnabled():
for n in range(1,N):
constraints += [#TCS([blade["r"][n] >= blade["r"][n-1] + parent.R/N]),
SignomialEquality(blade["r"][n],blade["r"][n-1] + parent.R/N)
constraints += [TCS([blade["r"][n] >= blade["r"][n-1] + parent.R/N]),
#SignomialEquality(blade["r"][n],blade["r"][n-1] + parent.R/N),
blade["eta_b"][n] == blade["eta_b"][n-1]
]
constraints += [TCS([T <= sum(blade["dT"][n] for n in range(N))]),
TCS([Q >= sum(blade["dQ"][n] for n in range(N))]),
eta == state.V*T/(omega*Q),
blade["M"][-1] <= Mtip,
#parent.T_m>=T
omega <= omega_max
]
Expand Down Expand Up @@ -219,6 +234,7 @@ class SimpleQProp(Model):
f [-] intermediate tip loss variable
F [-] Prandtl tip loss factor
cl_max 1.3 [-] max airfoil cl
eta_check [-] efficiency check
"""

def helper(self, c):
Expand Down Expand Up @@ -256,9 +272,9 @@ def setup(self, static, state):
AR_b == R/c,
AR_b <= AR_b_max,
Re == Wr*c*rho/mu,
#eta == T*V/(Q*omega),
eta_check == T*V/(Q*omega),
XfoilFit(fd, cd, [cl,Re], name="polar"),
parent.T_m >= T,
static.T_m >= T,
cl <= cl_max
]
with SignomialsEnabled():
Expand All @@ -274,7 +290,7 @@ def setup(self, static, state):

#SignomialEquality(1,etap*(1+eps*Wt/Wa)+eps*Wa/Wt),
#SignomialEquality(1,etai*(1+va/V)+vt/Ut),
SignomialEquality(Q,rho*B*G*(Wa+eps*Wt)*R**2/2.),
SignomialEquality(Q,rho*B*G*(Wa+eps*Wt)*R**2),
SignomialEquality(T,rho*B*G*(Wt-eps*Wa)*R),
]
return constraints, state
Expand All @@ -288,9 +304,32 @@ class Propeller(Model):
W [lbf] prop weight
K 4e-4 [1/ft^2] prop weight scaling factor
T_m [lbf] prop max static thrust
Variables of length N
---------------------
c [ft] prop chord
"""

flight_model = SimpleQProp

def setup(self, N = 1):
exec parse_variables(Propeller.__doc__)

return [W >= K*T_m*R**2]

class Actuator_Propeller(Model):
""" Propeller Model
Variables
---------
R [ft] prop radius
W [lbf] prop weight
K 4e-4 [1/ft^2] prop weight scaling factor
T_m [lbf] prop max static thrust
"""

flight_model = ActuatorProp

def setup(self):
exec parse_variables(Propeller.__doc__)

Expand Down

0 comments on commit da89baf

Please sign in to comment.