create gust loading from chordsparl

gpkitmodels/GP/aircraft/wing/chord_spar_loading.py

@@ -10,16 +10,8 @@ class ChordSparL(Model):
 
     def new_qbarFun(self, c):
         " define qbar model for chord loading "
-        return [f(c) for f in self.static.substitutions["\\bar{c}"]]
-
-    def new_SbarFun(self, c):
-        " define Sbar model for chord loading "
-        Sb = [1e-10]*self.static.N
-        for i in flip(range(self.static.N-1), 0):
-            Sb[i] = (Sb[i+1] + self.static.substitutions["d\\eta"][i](c)
-                     * (self.new_qbarFun(c)[i]
-                        + self.new_qbarFun(c)[i+1])/2)
-        return Sb
+        barc = self.static["\\bar{c}"]
+        return [f(c) for f in self.static.substitutions[barc]]
 
     def setup(self, static, Wcent):
 
@@ -31,25 +23,25 @@ def setup(self, static, Wcent):
         kappa = Variable("\\kappa", 0.2, "-", "max tip deflection ratio")
         W = Variable("W", "lbf", "loading weight")
 
-        with Vectorize(static.N-1):
+        with Vectorize(self.static.N-1):
             Mr = Variable("M_r", "N*m", "wing section root moment")
 
 
         Beam.qbarFun = self.new_qbarFun
-        Beam.SbarFun = self.new_SbarFun
-        beam = Beam(static.N)
+        self.beam = Beam(self.static.N)
 
         constraints = [
             # dimensionalize moment of inertia and young's modulus
-            beam["dx"] == static["d\\eta"],
+            self.beam["dx"] == self.static["d\\eta"],
             W == Wcent,
-            beam["\\bar{EI}"] <= (8*static["E"]*static["I"]/Nmax
-                                  / W/static["b"]**2),
-            Mr == (beam["\\bar{M}"][:-1]*W*Nmax*static["b"]/4),
-            sigmacfrp >= Mr/static["S_y"],
-            beam["\\bar{\\delta}"][-1] <= kappa,
-            taucfrp >= (beam["\\bar{S}"][-1]*W*Nmax/4/static["t_{shear}"]
-                        / static["c_{ave}"]/static["\\tau"])
+            self.beam["\\bar{EI}"] <= (8*self.static["E"]*self.static["I"]/Nmax
+                                       / W/self.static["b"]**2),
+            Mr == (self.beam["\\bar{M}"][:-1]*W*Nmax*self.static["b"]/4),
+            sigmacfrp >= Mr/self.static["S_y"],
+            self.beam["\\bar{\\delta}"][-1] <= kappa,
+            taucfrp >= (self.beam["\\bar{S}"][-1]*W*Nmax/4
+                        / self.static["t_{shear}"]/self.static["c_{ave}"]
+                        / self.static["\\tau"])
             ]
 
-        return beam, constraints
+        return self.beam, constraints

gpkitmodels/GP/aircraft/wing/gustloading.py

@@ -2,25 +2,21 @@
 import os
 from numpy import pi, hstack, cos
 import pandas as pd
-from gpkit import Model, Variable, Vectorize
-from gpkitmodels.GP.beam.beam import Beam
+from gpkit import Variable, Vectorize
 from gpfit.fit_constraintset import FitCS
+from .chord_spar_loading import ChordSparL
 
-#pylint: disable=invalid-name, no-member, too-many-locals
+#pylint: disable=invalid-name, no-member, arguments-differ
+#pylint: disable=attribute-defined-outside-init
 
-class GustL(Model):
+class GustL(ChordSparL):
     "spar loading model"
-    def setup(self, static, Wcent, Wwing, V, CL):
-
-        Nmax = Variable("N_{max}", 2, "-", "load safety factor")
+    new_qbarFun = None
+    new_SbarFun = None
 
-        sigmacfrp = Variable("\\sigma_{CFRP}", 1700e6, "Pa", "CFRP max stress")
-        taucfrp = Variable("\\tau_{CFRP}", 450e6, "Pa", "CFRP fabric stress")
-        kappa = Variable("\\kappa", 0.2, "-", "max tip deflection ratio")
-
-        with Vectorize(static.N-1):
-            Mr = Variable("M_r", "N*m", "wing section root moment")
+    def setup(self, static, Wcent, Wwing, V, CL):
 
+        self.load = ChordSparL.setup(self, static, Wcent)
         vgust = Variable("V_{gust}", 10, "m/s", "gust velocity")
 
         with Vectorize(static.N):
@@ -30,27 +26,15 @@ def setup(self, static, Wcent, Wwing, V, CL):
             cosminus1 = Variable("1-cos(\\eta)", return_cosm1,
                                  "-", "1 minus cosine factor")
 
-        Beam.qbarFun = None
-        Beam.SbarFun = None
-        beam = Beam(static.N)
-        path = os.path.abspath(__file__).replace(os.path.basename(__file__), "")
+        path = os.path.dirname(os.path.abspath(__file__))
         df = pd.read_csv(path + os.sep + "arctan_fit.csv").to_dict(
             orient="records")[0]
 
         constraints = [
             # fit for arctan from 0 to 1, RMS = 0.044
             FitCS(df, agust, [cosminus1*vgust/V]),
-            beam["\\bar{q}"] >= static["\\bar{c}"]*(
+            self.beam["\\bar{q}"] >= self.static["\\bar{c}"]*(
                 1 + 2*pi*agust/CL*(1+Wwing/Wcent)),
-            beam["dx"] == static["d\\eta"],
-            # dimensionalize moment of inertia and young's modulus
-            beam["\\bar{EI}"] <= (8*static["E"]*static["I"]/Nmax
-                                  / Wcent/static["b"]**2),
-            Mr == (beam["\\bar{M}"][:-1]*Wcent*Nmax*static["b"]/4),
-            sigmacfrp >= Mr/static["S_y"],
-            beam["\\bar{\\delta}"][-1] <= kappa,
-            taucfrp >= (beam["\\bar{S}"][-1]*Wcent*Nmax/4/static["t_{shear}"]
-                        / static["c_{ave}"]/static["\\tau"])
             ]
 
-        return beam, constraints
+        return self.load, constraints

gpkitmodels/GP/beam/beam.py

@@ -17,26 +17,28 @@ def setup(self, N):
             dx = Variable("dx", "-", "normalized length of element")
 
         with Vectorize(N):
-            qbar = Variable("\\bar{q}", self.qbarFun, "-", "normalized loading")
             Sbar = Variable("\\bar{S}", self.SbarFun, "-", "normalized shear")
             Mbar = Variable("\\bar{M}", self.MbarFun, "-", "normalized moment")
             th = Variable("\\theta", "-", "deflection slope")
             dbar = Variable("\\bar{\\delta}", "-", "normalized displacement")
 
 
-        Sbartip = Variable("\\bar{S}_{tip}", 1e-10, "-", "Tip loading")
-        Mbartip = Variable("\\bar{M}_{tip}", 1e-10, "-", "Tip moment")
         throot = Variable("\\theta_{root}", 1e-10, "-", "Base angle")
         dbarroot = Variable("\\bar{\\delta}_{root}", 1e-10, "-",
                             "Base deflection")
 
         constraints = []
         if self.SbarFun is None:
+            with Vectorize(N):
+                qbar = Variable("\\bar{q}", self.qbarFun, "-",
+                                "normalized loading")
+            Sbartip = Variable("\\bar{S}_{tip}", 1e-10, "-", "Tip loading")
             constraints.extend([
                 Sbar[:-1] >= Sbar[1:] + 0.5*dx*(qbar[:-1] + qbar[1:]),
                 Sbar[-1] >= Sbartip])
 
         if self.MbarFun is None:
+            Mbartip = Variable("\\bar{M}_{tip}", 1e-10, "-", "Tip moment")
             constraints.extend([
                 Mbar[:-1] >= Mbar[1:] + 0.5*dx*(Sbar[:-1] + Sbar[1:]),
                 Mbar[-1] >= Mbartip])