Merge pull request #523 from WISDEM/develop

Prep for next release

.github/workflows/CI_WISDEM.yml

@@ -92,13 +92,11 @@ jobs:
       # Run coveralls
       - name: Run coveralls
         if: contains( matrix.os, 'ubuntu') && contains( matrix.python-version, '3.10')
-        uses: coverallsapp/github-action@v2
-        # This also works, https://github.com/AndreMiras/coveralls-python-action
-        #uses: AndreMiras/coveralls-python-action@develop
-        #env:
-        #  GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
-        #run: |
-        #  coveralls --service=github
+        #uses: coverallsapp/github-action@v2
+        env:
+          GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }}
+        run: |
+          coveralls --service=github
           
   build_pip:
     name: Pip Build-Test (${{ matrix.os }} Python ${{ matrix.python-version }})

docs/wisdem/drivetrainse/layout.rst

@@ -51,7 +51,7 @@ In addition to the user-defined dimensions, the other values are derived in the
    L_{nose}	&= L_{12} + L_{2n}\\
    L_{drive}    &= L_{h1} + L_{12} + L_{2n}\\
    L_{bedplate} &= L_{overhang} - L_{drive}\cos \gamma\\
-   H_{bedplate} &= H_{htt} - L_{drive}\sin \gamma
+   H_{bedplate} &= H_{htt} - (L_{drive} + D_{hub}/2)\sin \gamma
 
 Here the length from the hub flange to the generator rotor attachment, :math:`L_{grs}`, is assumed to be at the halfway point between the flange and the first main bearing, :math:`L_{h1}`.  Similarly, the distance between the second main bearing and the nose/turret interface with the bedplate, :math:`L_{2n}`, is twice the distance as that from the same interface to the generator stator attachment, :math:`L_{gsn}`.  After adding up the total length of the low speed shaft and nose/turret, the total drivetrain length from bedplate to hub can be determined.  Then, the bedplate dimensions are determined in order to meet the target overhang and hub-to-tower top height. To ensure that these layout dimensions are adequately satisfied during a design optimization, a constraint is enforced such that :math:`L_{bed} \geq 0.5 D_{top}`.
 
@@ -104,6 +104,6 @@ In addition to the user-defined dimensions, the other values are derived in the
    L_{lss}      &= L_{12} + L_{h1} + \delta\\
    L_{drive}    &= L_{lss} + L_{gearbox} + L_{hss} + L_{generator}\\
    L_{bedplate} &= L_{drive} \cos \gamma \\
-   H_{bedplate} &= H_{htt} - L_{drive} \sin \gamma
+   H_{bedplate} &= H_{htt} - (L_{drive} + D_{hub}/2) \sin \gamma
 
 The dimension, :math:`\delta` is the space between the second main bearing and the gearbox attachment where the shrink disk lies.  This is assumed to be 0.1 meters.  The bedplate height is sized to ensure that the desired height from tower top to hub is obtained.  To achieve the desired overhang distance, the tower is centered at the exact overhang distance from the hub and a constraint is enforced such that the drivetrain length is sufficient to extend past the tower, :math:`L_{drive} \cos \gamma - L_{overhang} \geq 0.5 D_{top}`.

environment.yml

@@ -11,16 +11,15 @@ dependencies:
   - matplotlib
   - meson
   - moorpy
-  - mpi4py
   - ninja
   - nlopt
   - numpy>=1.26
   - numpydoc
-  - openmdao<3.28
+  - openmdao>=3.31
   - openpyxl
   - pandas
   - pip
-  - pydoe2
+  - pydoe3
   - pyoptsparse
   - pytest
   - pytest-cov
@@ -41,3 +40,4 @@ dependencies:
 #  - libpython # [win]
 #  - compilers # [not win]
 #  - petsc4py # [not win]
+#  - mpi4py (only needed if local system has mpi built and is used)

examples/15_step_size_study/driver.py

@@ -22,7 +22,7 @@
 all_derivs = np.load("total_derivs.npy", allow_pickle=True)[()]
 
 # Select a key of interest; an (of, wrt) pair
-key = ("wt.towerse.post.constr_stress", "wt.wt_init.tower.diameter")
+key = ("towerse.post.constr_stress", "tower.diameter")
 
 # Collect data to plot later
 step_sizes = []
@@ -49,5 +49,5 @@
     plt.gca().invert_xaxis()
 
     plt.show()
-except:
+except Exception:
     pass

pyproject.toml

@@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta"
 
 [project]
 name = "wisdem"
-version = "3.14.0"
+version = "3.15.0"
 description = "Wind-Plant Integrated System Design & Engineering Model"
 readme = "README.md"
 requires-python = ">=3.9"
@@ -49,10 +49,10 @@ dependencies = [
   "moorpy",
   "nlopt",
   "numpy>=1.26",
-  "openmdao<3.28",
+  "openmdao>=3.31",
   "openpyxl",
   "pandas",
-  "pydoe2",
+  "pydoe3",
   "python-benedict",
   "pyyaml",
   "ruamel.yaml",

wisdem/ccblade/Polar.py

@@ -1684,8 +1684,11 @@ def _find_slope(x, y, xi=None, x0=None, window=None, method="max", opts=None, nI
             dx = x[im + 1] - x[im - 1]
             if np.abs(dx) > 1e-7:
                 a = (y[im + 1] - y[im - 1]) / dx
-                yi = np.interp(xi, x, y)
-                x0 = xi - yi / a
+                if a != 0.:
+                    yi = np.interp(xi, x, y)
+                    x0 = xi - yi / a
+                else:
+                    x0 = xi
             else:
                 a = np.inf
                 x0 = xi
@@ -2004,4 +2007,4 @@ def smooth_heaviside(x, k=1, rng=(-np.inf, np.inf), method="exp"):
 
 
 if __name__ == "__main__":
-    pass
+    pass

wisdem/ccblade/ccblade_component.py

@@ -62,8 +62,8 @@ def setup(self):
         self.declare_partials("presweepTip", "presweep_in", val=1.0, rows=[0], cols=[n_span - 1])
 
     def compute(self, inputs, outputs):
-        Rtip = inputs["Rtip"]
-        precone = inputs["precone"]
+        Rtip = inputs["Rtip"][0]
+        precone = inputs["precone"][0]
 
         outputs["precurveTip"] = inputs["precurve_in"][-1]
         outputs["presweepTip"] = inputs["presweep_in"][-1]
@@ -270,27 +270,27 @@ def compute(self, inputs, outputs, discrete_inputs, discrete_outputs):
         r = inputs["r"]
         chord = inputs["chord"]
         theta = inputs["theta"]
-        Rhub = inputs["Rhub"]
-        Rtip = inputs["Rtip"]
-        hub_height = inputs["hub_height"]
-        precone = inputs["precone"]
-        tilt = inputs["tilt"]
-        yaw = inputs["yaw"]
+        Rhub = inputs["Rhub"][0]
+        Rtip = inputs["Rtip"][0]
+        hub_height = inputs["hub_height"][0]
+        precone = inputs["precone"][0]
+        tilt = inputs["tilt"][0]
+        yaw = inputs["yaw"][0]
         precurve = inputs["precurve"]
-        precurveTip = inputs["precurveTip"]
+        precurveTip = inputs["precurveTip"][0]
         B = discrete_inputs["nBlades"]
-        rho = inputs["rho"]
-        mu = inputs["mu"]
-        shearExp = inputs["shearExp"]
+        rho = inputs["rho"][0]
+        mu = inputs["mu"][0]
+        shearExp = inputs["shearExp"][0]
         nSector = discrete_inputs["nSector"]
         tiploss = discrete_inputs["tiploss"]
         hubloss = discrete_inputs["hubloss"]
         wakerotation = discrete_inputs["wakerotation"]
         usecd = discrete_inputs["usecd"]
-        V_load = inputs["V_load"]
-        Omega_load = inputs["Omega_load"]
-        pitch_load = inputs["pitch_load"]
-        azimuth_load = inputs["azimuth_load"]
+        V_load = inputs["V_load"][0]
+        Omega_load = inputs["Omega_load"][0]
+        pitch_load = inputs["pitch_load"][0]
+        azimuth_load = inputs["azimuth_load"][0]
 
         if len(precurve) == 0:
             precurve = np.zeros_like(r)
@@ -559,28 +559,28 @@ def compute(self, inputs, outputs, discrete_inputs, discrete_outputs):
             inputs["chord"],
             np.zeros_like(inputs["chord"]),
             af,
-            inputs["Rhub"],
-            inputs["Rtip"],
+            inputs["Rhub"][0],
+            inputs["Rtip"][0],
             discrete_inputs["nBlades"],
-            inputs["rho"],
-            inputs["mu"],
-            inputs["precone"],
-            inputs["tilt"],
-            inputs["yaw"],
-            inputs["shearExp"],
-            inputs["hub_height"],
+            inputs["rho"][0],
+            inputs["mu"][0],
+            inputs["precone"][0],
+            inputs["tilt"][0],
+            inputs["yaw"][0],
+            inputs["shearExp"][0],
+            inputs["hub_height"][0],
             discrete_inputs["nSector"],
             inputs["precurve"],
-            inputs["precurveTip"],
+            inputs["precurveTip"][0],
             inputs["presweep"],
-            inputs["presweepTip"],
+            inputs["presweepTip"][0],
             discrete_inputs["tiploss"],
             discrete_inputs["hubloss"],
             discrete_inputs["wakerotation"],
             discrete_inputs["usecd"],
         )
 
-        Omega = inputs["tsr"] * inputs["Uhub"] / inputs["r"][-1] * 30.0 / np.pi
+        Omega = inputs["tsr"][0] * inputs["Uhub"][0] / inputs["r"][-1] * 30.0 / np.pi
 
         if self.options["opt_options"]["design_variables"]["blade"]["aero_shape"]["twist"]["inverse"]:
             if self.options["opt_options"]["design_variables"]["blade"]["aero_shape"]["twist"]["flag"]:
@@ -647,7 +647,7 @@ def compute(self, inputs, outputs, discrete_inputs, discrete_outputs):
         ccblade.inverse_analysis = False
 
         # Call ccblade at azimuth 0 deg
-        loads, _ = ccblade.distributedAeroLoads(inputs["Uhub"][0], Omega[0], inputs["pitch"][0], 0.0)
+        loads, _ = ccblade.distributedAeroLoads(inputs["Uhub"][0], Omega, inputs["pitch"][0], 0.0)
 
         # Call ccblade evaluate (averaging across azimuth)
         myout, _ = ccblade.evaluate([inputs["Uhub"]], [Omega], [inputs["pitch"]], coefficients=True)
@@ -768,28 +768,28 @@ def compute(self, inputs, outputs, discrete_inputs, discrete_outputs):
         r = inputs["r"]
         chord = inputs["chord"]
         theta = inputs["theta"]
-        Rhub = inputs["Rhub"]
-        Rtip = inputs["Rtip"]
-        hub_height = inputs["hub_height"]
-        precone = inputs["precone"]
-        tilt = inputs["tilt"]
-        yaw = inputs["yaw"]
+        Rhub = inputs["Rhub"][0]
+        Rtip = inputs["Rtip"][0]
+        hub_height = inputs["hub_height"][0]
+        precone = inputs["precone"][0]
+        tilt = inputs["tilt"][0]
+        yaw = inputs["yaw"][0]
         precurve = inputs["precurve"]
-        precurveTip = inputs["precurveTip"]
+        precurveTip = inputs["precurveTip"][0]
         presweep = inputs["presweep"]
-        presweepTip = inputs["presweepTip"]
+        presweepTip = inputs["presweepTip"][0]
         B = discrete_inputs["nBlades"]
-        rho = inputs["rho"]
-        mu = inputs["mu"]
-        shearExp = inputs["shearExp"]
+        rho = inputs["rho"][0]
+        mu = inputs["mu"][0]
+        shearExp = inputs["shearExp"][0]
         nSector = discrete_inputs["nSector"]
         tiploss = discrete_inputs["tiploss"]
         hubloss = discrete_inputs["hubloss"]
         wakerotation = discrete_inputs["wakerotation"]
         usecd = discrete_inputs["usecd"]
-        V_load = inputs["V_load"]
-        Omega_load = inputs["Omega_load"]
-        pitch_load = inputs["pitch_load"]
+        V_load = inputs["V_load"][0]
+        Omega_load = inputs["Omega_load"][0]
+        pitch_load = inputs["pitch_load"][0]
 
         if len(precurve) == 0:
             precurve = np.zeros_like(r)
@@ -846,28 +846,28 @@ def compute_partials(self, inputs, J, discrete_inputs):
         r = inputs["r"]
         chord = inputs["chord"]
         theta = inputs["theta"]
-        Rhub = inputs["Rhub"]
-        Rtip = inputs["Rtip"]
-        hub_height = inputs["hub_height"]
-        precone = inputs["precone"]
-        tilt = inputs["tilt"]
-        yaw = inputs["yaw"]
+        Rhub = inputs["Rhub"][0]
+        Rtip = inputs["Rtip"][0]
+        hub_height = inputs["hub_height"][0]
+        precone = inputs["precone"][0]
+        tilt = inputs["tilt"][0]
+        yaw = inputs["yaw"][0]
         precurve = inputs["precurve"]
-        precurveTip = inputs["precurveTip"]
+        precurveTip = inputs["precurveTip"][0]
         presweep = inputs["presweep"]
-        presweepTip = inputs["presweepTip"]
+        presweepTip = inputs["presweepTip"][0]
         B = discrete_inputs["nBlades"]
-        rho = inputs["rho"]
-        mu = inputs["mu"]
-        shearExp = inputs["shearExp"]
+        rho = inputs["rho"][0]
+        mu = inputs["mu"][0]
+        shearExp = inputs["shearExp"][0]
         nSector = discrete_inputs["nSector"]
         tiploss = discrete_inputs["tiploss"]
         hubloss = discrete_inputs["hubloss"]
         wakerotation = discrete_inputs["wakerotation"]
         usecd = discrete_inputs["usecd"]
-        V_load = inputs["V_load"]
-        Omega_load = inputs["Omega_load"]
-        pitch_load = inputs["pitch_load"]
+        V_load = inputs["V_load"][0]
+        Omega_load = inputs["Omega_load"][0]
+        pitch_load = inputs["pitch_load"][0]
 
         if len(precurve) == 0:
             precurve = np.zeros_like(r)

wisdem/commonse/cylinder_member.py

@@ -2146,6 +2146,17 @@ def setup(self):
         self.add_input("rho_full", np.zeros(n_full - 1), units="kg/m**3")
         self.add_input("sigma_y_full", np.zeros(n_full - 1), units="Pa")
 
+        self.add_input("section_A", np.zeros(n_full - 1), units="m**2")
+        self.add_input("section_Asx", np.zeros(n_full - 1), units="m**2")
+        self.add_input("section_Asy", np.zeros(n_full - 1), units="m**2")
+        self.add_input("section_Ixx", np.zeros(n_full - 1), units="kg*m**2")
+        self.add_input("section_Iyy", np.zeros(n_full - 1), units="kg*m**2")
+        self.add_input("section_J0", np.zeros(n_full - 1), units="kg*m**2")
+        self.add_input("section_rho", np.zeros(n_full - 1), units="kg/m**3")
+        self.add_input("section_E", np.zeros(n_full - 1), units="Pa")
+        self.add_input("section_G", np.zeros(n_full - 1), units="Pa")
+        self.add_input("section_L", np.zeros(n_full - 1), units="m")
+
         # Processed Frame3DD/OpenFAST outputs
         self.add_input("cylinder_Fz", val=np.zeros((n_full - 1, n_dlc)), units="N")
         self.add_input("cylinder_Vx", val=np.zeros((n_full - 1, n_dlc)), units="N")
@@ -2177,7 +2188,15 @@ def compute(self, inputs, outputs):
         h = np.diff(z, axis=0)
         d_sec, _ = util.nodal2sectional(d)
         r_sec = 0.5 * d_sec
-        itube = cs.Tube(d_sec, t)
+
+        # Geom properties
+        #itube = cs.Tube(d_sec, t)
+        Az = np.tile(inputs["section_A"], (n_dlc, 1)).T
+        Asx = np.tile(inputs["section_Asx"], (n_dlc, 1)).T
+        Asy = np.tile(inputs["section_Asy"], (n_dlc, 1)).T
+        Jz = np.tile(inputs["section_J0"], (n_dlc, 1)).T
+        Ixx = np.tile(inputs["section_Ixx"], (n_dlc, 1)).T
+        Iyy = np.tile(inputs["section_Iyy"], (n_dlc, 1)).T
 
         L_buckling = self.options["modeling_options"]["buckling_length"]
         gamma_f = self.options["modeling_options"]["gamma_f"]
@@ -2206,14 +2225,6 @@ def compute(self, inputs, outputs):
         M = np.sqrt(Mxx**2 + Myy**2)
         V = np.sqrt(Vx**2 + Vy**2)
 
-        # Geom properties
-        Az = itube.Area
-        Asx = itube.Asx
-        Asy = itube.Asy
-        Jz = itube.J0
-        Ixx = itube.Ixx
-        Iyy = itube.Iyy
-
         # See http://svn.code.sourceforge.net/p/frame3dd/code/trunk/doc/Frame3DD-manual.html#structuralmodeling
         outputs["axial_stress"] = axial_stress = Fz / Az + M * r_sec / Iyy
         outputs["shear_stress"] = shear_stress = np.abs(Mzz) / Jz * r_sec + V / Asx
@@ -2257,6 +2268,8 @@ def compute(self, inputs, outputs):
                 G=G[:, 0],
                 sigma_y=sigma_y[:, 0],
                 gamma=gamma_f * gamma_b,
+                A=Az[:,0],
+                I=Ixx[:,0],
             )
 
             for k in range(n_dlc):