Merge commit '1dcda881a577e669c7abc4a9db740a278f51e781' into no_torsi…

…on_modes

WISDEM/docs/wisdem/nrelcsm/theory.rst

@@ -42,12 +42,12 @@ To obtain the blade mass in kilograms and cost in USD from the rotor diameter in
    b   &= (see below)\\
    k_c &= 14.6
 
-Where :math:`D_{rotor}` is the rotor diameter and :math:`k_c` is determined by:
+Where :math:`D_{rotor}` is the rotor diameter and :math:`b` is determined by:
 
-* If turbine class I and blade DOES have carbon fiber spar caps, :math:`k_c=2.47`
-* If turbine class I and blade DOES NOT have carbon fiber spar caps, :math:`k_c=2.54`
-* If turbine class II+ and blade DOES have carbon fiber spar caps, :math:`k_c=2.44`
-* If turbine class II+ and blade DOES NOT have carbon fiber spar caps, :math:`k_c=2.50`
+* If turbine class I and blade DOES have carbon fiber spar caps, :math:`b=2.47`
+* If turbine class I and blade DOES NOT have carbon fiber spar caps, :math:`b=2.54`
+* If turbine class II+ and blade DOES have carbon fiber spar caps, :math:`b=2.44`
+* If turbine class II+ and blade DOES NOT have carbon fiber spar caps, :math:`b=2.50`
 * User override of exponent value
 
 For variable names access to override the default values see the :ref:`csmsource`.

WISDEM/wisdem/commonse/utilities.py

@@ -31,7 +31,7 @@ def get_modal_coefficients(x, y, deg=6):
 
     return p6
 
-def get_mode_shapes(r, freqs, xdsp, ydsp, zdsp, xmpf, ympf, zmpf):
+def get_xy_mode_shapes(r, freqs, xdsp, ydsp, zdsp, xmpf, ympf, zmpf):
     # Number of frequencies and modes
     nfreq = len(freqs)
 
@@ -61,8 +61,6 @@ def get_mode_shapes(r, freqs, xdsp, ydsp, zdsp, xmpf, ympf, zmpf):
             mshapes_y[iy,:] = ypolys[m,:]
             freq_y[   iy  ] = freqs[m]
             iy += 1
-        else:
-            print('Warning: Unknown mode shape')
 
     return freq_x, freq_y, mshapes_x, mshapes_y
 

WISDEM/wisdem/commonse/vertical_cylinder.py

@@ -545,7 +545,8 @@ def compute(self, inputs, outputs):
         Mmethod = 1
         lump = 0
         shift = 0.0
-        cylinder.enableDynamics(NFREQ, Mmethod, lump, frame3dd_opt['tol'], shift)
+        # Run twice the number of modes to ensure that we can ignore the torsional modes and still get the desired number of fore-aft, side-side modes
+        cylinder.enableDynamics(2*NFREQ, Mmethod, lump, frame3dd_opt['tol'], shift)
         # ----------------------------
 
         # ------ static load case 1 ------------
@@ -592,15 +593,16 @@ def compute(self, inputs, outputs):
         # natural frequncies
         outputs['f1']    = modal.freq[0]
         outputs['f2']    = modal.freq[1]
-        outputs['freqs'] = modal.freq
+        outputs['freqs'] = modal.freq[:NFREQ]
 
         # Get all mode shapes in batch
-        freq_x, freq_y, mshapes_x, mshapes_y = util.get_mode_shapes(z, modal.freq, modal.xdsp, modal.ydsp, modal.zdsp, modal.xmpf, modal.ympf, modal.zmpf)
-        outputs['x_mode_freqs']  = freq_x
-        outputs['y_mode_freqs']  = freq_y
-        outputs['x_mode_shapes'] = mshapes_x
-        outputs['y_mode_shapes'] = mshapes_y
-
+        NFREQ2 = int(NFREQ/2)
+        freq_x, freq_y, mshapes_x, mshapes_y = util.get_xy_mode_shapes(z, modal.freq, modal.xdsp, modal.ydsp, modal.zdsp, modal.xmpf, modal.ympf, modal.zmpf)
+        outputs['x_mode_freqs']  = freq_x[:NFREQ2]
+        outputs['y_mode_freqs']  = freq_y[:NFREQ2]
+        outputs['x_mode_shapes'] = mshapes_x[:NFREQ2,:]
+        outputs['y_mode_shapes'] = mshapes_y[:NFREQ2,:]
+
         # deflections due to loading (from cylinder top and wind/wave loads)
         outputs['top_deflection'] = displacements.dx[iCase, n-1]  # in yaw-aligned direction
 

WISDEM/wisdem/glue_code/glue_code.py

@@ -295,7 +295,6 @@ def setup(self):
                 self.connect('nacelle.lss_diameter',              'drivese.bear1.D_shaft', src_indices=[0])
                 self.connect('nacelle.lss_diameter',              'drivese.bear2.D_shaft', src_indices=[-1])
             self.connect('nacelle.uptower',                   'drivese.uptower')
-            #self.connect('nacelle.gearbox_efficiency',        'drivese.gearbox_efficiency')
             self.connect('nacelle.brake_mass_user',           'drivese.brake_mass_user')
             self.connect('nacelle.hvac_mass_coeff',           'drivese.hvac_mass_coeff')
             self.connect('nacelle.converter_mass_user',       'drivese.converter_mass_user')

WISDEM/wisdem/rotorse/rotor_loads_defl_strains.py

@@ -316,7 +316,8 @@ def rotate(x,y):
         lump    = 0 # 0= consistent mass matrix, 1= lumped mass matrix
         tol     = 1e-9 # frequency convergence tolerance
         shift   = 0.0 # frequency shift-factor for rigid body modes, make 0 for pos.def. [K]
-        blade.enableDynamics(self.n_freq, Mmethod, lump, tol, shift)
+        # Run twice the number of modes to ensure that we can ignore the torsional modes and still get the desired number of fore-aft, side-side modes
+        blade.enableDynamics(2*self.n_freq, Mmethod, lump, tol, shift)
         # ----------------------------
 
         # ------ load case 1, blade 1 ------------
@@ -358,8 +359,13 @@ def rotate(x,y):
         dz = displacements.dz[iCase,:]
 
         # Mode shapes and frequencies
-        freq_x, freq_y, mshapes_x, mshapes_y = util.get_mode_shapes(r, modal.freq, modal.xdsp, modal.ydsp, modal.zdsp, modal.xmpf, modal.ympf, modal.zmpf)
-
+        n_freq2 = int(self.n_freq/2)
+        freq_x, freq_y, mshapes_x, mshapes_y = util.get_xy_mode_shapes(r, modal.freq, modal.xdsp, modal.ydsp, modal.zdsp, modal.xmpf, modal.ympf, modal.zmpf)
+        freq_x = freq_x[:n_freq2]
+        freq_y = freq_y[:n_freq2]
+        mshapes_x = mshapes_x[:n_freq2,:]
+        mshapes_y = mshapes_y[:n_freq2,:]
+
         # shear and bending, one per element (convert from local to global c.s.)
         Fz = np.r_[-forces.Nx[iCase,0],  forces.Nx[iCase, 1::2]]
         Vy = np.r_[-forces.Vy[iCase,0],  forces.Vy[iCase, 1::2]]
@@ -396,7 +402,7 @@ def strain(xu, yu, xl, yl):
         # Store outputs
         outputs['root_F'] = -1.0 * np.array([reactions.Fx.sum(), reactions.Fy.sum(), reactions.Fz.sum()])
         outputs['root_M'] = -1.0 * np.array([reactions.Mxx.sum(), reactions.Myy.sum(), reactions.Mzz.sum()])
-        outputs['freqs'] = modal.freq
+        outputs['freqs'] = modal.freq[:self.n_freq]
         outputs['edge_mode_shapes'] = mshapes_y
         outputs['flap_mode_shapes'] = mshapes_x
         # Dense numpy command that interleaves and alternates flap and edge modes

WISDEM/wisdem/test/test_commonse/test_utilities.py

@@ -33,6 +33,23 @@ def testModalCoefficients(self):
         pp  = util.get_modal_coefficients(x, y)
         npt.assert_almost_equal(p[2:]/p[2:].sum(), pp)
 
+    def testGetXYModes(self):
+        r = np.linspace(0,1,20)
+        n = 10
+        n2 = int(n/2)
+        dx = dy = dz = np.tile(np.r_[r**2 + 10.], (n,1))
+        freqs = np.arange(n)
+        xm = np.zeros(n)
+        ym = np.zeros(n)
+        zm = np.zeros(n)
+        xm[0] = xm[3] = xm[6] = xm[9] = 1
+        ym[1] = ym[4] = ym[7] = 1
+        zm[2] = zm[5] = zm[8] = 1
+
+        freq_x, freq_y, _, _ = util.get_xy_mode_shapes(r, freqs, dx, dy, dz, xm, ym, zm)
+        npt.assert_array_equal(freq_x, np.r_[0, 3, 6, 9, np.zeros(n2-4)])
+        npt.assert_array_equal(freq_y, np.r_[1, 4, 7, np.zeros(n2-3)])
+
     def testRotateI(self):
         I  = np.arange(6)+1
         th = np.deg2rad(45)