gradients now returned as dictionary and as full Jacobian. Makes them…

… easier to access and more consistent. Fixed a transpose issue for dNp_dprecurve, etc. updated unit tests. updated version number

CHANGELOG.md

@@ -1,5 +1,24 @@
 # CCBlade Changelog
 
+## 1.1.0 (Jan 10, 2013)
+
+Andrew Ning <andrew.ning@nrel.gov>
+
+[NEW]:
+
+- added the remaining analytic gradients I didn't think I would need (yaw, azimuth, Uinf, Omega, pitch)
+- added analytic gradients for parked cases (Omega==0)
+- additional unit tests for the new gradients
+
+[CHANGE]:
+
+- gradients are now returned as a dictionary, and full Jacobians are returned.  This makes it much easier to access the gradients, because you don't have to worry about what the ordering is.  Also returning 2D arrays for every value maintains consistency so you don't have to remember which gradients were really just the nonzero diagonal entries of the corresponding Jacobian.  However, this makes the API not backwards compatible.  I apologize if it causes a lot of changes on your end.
+- AirfoilPrep.py was moved to a separate repository (https://github.com/NREL-WISDEM/AirfoilPreppy) to more easily accommodate changes and contributions.  The setup.py script should still install from the repository automatically.  For some reason the setup script throws an error saying it could not find AirfoilPrep however it does find it and installs properly.
+
+[FIX]:
+
+- As compared to what I wrote in the docstring, I was returning the transpose of dNp_dprecurve, dTp_dprecurve, dNp_dpresweep, and dTp_dpresweep.  This is now fixed.
+
 ## 1.0.2 (Nov 21, 2013)
 
 Andrew Ning <andrew.ning@nrel.gov>

docs/conf.py

@@ -22,8 +22,8 @@
 codeName = u'CCBlade'
 authorName = u'S. Andrew Ning'
 copyright = u'2013, NREL'
-version = '1.0'
-release = '1.0.5'
+version = '1.1'
+release = '1.1.0'
 
 # -- General configuration -----------------------------------------------------
 

docs/examples/gradients.py

@@ -94,53 +94,64 @@
 
 # 5 ----------
 
-Np, Tp, dNp_dX, dTp_dX, dNp_dprecurve, dTp_dprecurve = \
-    rotor.distributedAeroLoads(Uinf, Omega, pitch, azimuth)
-
-# X = [r, chord, theta, Rhub, Rtip, presweep, precone,
-#      tilt, hubHt, yaw, azimuth, Uinf, Omega, pitch]
-
-dNp_dr = dNp_dX[0, :]
-dNp_dchord = dNp_dX[1, :]
-dNp_dtheta = dNp_dX[2, :]
-dNp_dRhub = dNp_dX[3, :]
-dNp_dRtip = dNp_dX[4, :]
-dNp_dpresweep = dNp_dX[5, :]
-dNp_dprecone = dNp_dX[6, :]
-dNp_dtilt = dNp_dX[7, :]
-dNp_dhubHt = dNp_dX[8, :]
-dNp_dyaw = dNp_dX[9, :]
-dNp_dazimuth = dNp_dX[10, :]
-dNp_dUinf = dNp_dX[11, :]
-dNp_dOmega = dNp_dX[12, :]
-dNp_dpitch = dNp_dX[13, :]
+Np, Tp, dNp, dTp \
+    = rotor.distributedAeroLoads(Uinf, Omega, pitch, azimuth)
+
+n = len(r)
+
+# n x n (diagonal)
+dNp_dr = dNp['dr']
+dNp_dchord = dNp['dchord']
+dNp_dtheta = dNp['dtheta']
+dNp_dpresweep = dNp['dpresweep']
+
+# n x n (tridiagonal)
+dNp_dprecurve = dNp['dprecurve']
+
+# n x 1
+dNp_dRhub = dNp['dRhub']
+dNp_dRtip = dNp['dRtip']
+dNp_dprecone = dNp['dprecone']
+dNp_dtilt = dNp['dtilt']
+dNp_dhubHt = dNp['dhubHt']
+dNp_dyaw = dNp['dyaw']
+dNp_dazimuth = dNp['dazimuth']
+dNp_dUinf = dNp['dUinf']
+dNp_dOmega = dNp['dOmega']
+dNp_dpitch = dNp['dpitch']
 
 # 5 ----------
 
 
 # 6 ----------
 
 
-P, T, Q, dP_ds, dT_ds, dQ_ds, dP_dv, dT_dv, dQ_dv = \
-    rotor.evaluate([Uinf], [Omega], [pitch])
+P, T, Q, dP, dT, dQ \
+    = rotor.evaluate([Uinf], [Omega], [pitch])
 
-# scalars = [precone, tilt, hubHt, Rhub, Rtip,
-#            precurvetip, presweeptip, yaw,
-#            Uinf, Omega, pitch]
-# vectors = [r, chord, theta, precurve, presweep]
+npts = len(P)
 
-# unpack scalar derivatives
-dP_dprecone, dP_dtilt, dP_dhubHt, dP_dRhub, dP_dRtip, \
-    dP_dprecurvetip, dP_dpresweeptip, dP_dyaw, \
-    dP_dUinf, dP_dOmega, dP_dpitch = dP_ds[0, :]
+# npts x 1
+dP_dprecone = dP['dprecone']
+dP_dtilt = dP['dtilt']
+dP_dhubHt = dP['dhubHt']
+dP_dRhub = dP['dRhub']
+dP_dRtip = dP['dRtip']
+dP_dprecurveTip = dP['dprecurveTip']
+dP_dpresweepTip = dP['dpresweepTip']
+dP_dyaw = dP['dyaw']
 
-# unpack vector derivatives
-dP_dr = dP_dv[0, 0, :]
-dP_dchord = dP_dv[0, 1, :]
-dP_dtheta = dP_dv[0, 2, :]
-dP_dprecurve = dP_dv[0, 3, :]
-dP_dpresweep = dP_dv[0, 4, :]
+# npts x npts
+dP_dUinf = dP['dUinf']
+dP_dOmega = dP['dOmega']
+dP_dpitch = dP['dpitch']
 
+# npts x n
+dP_dr = dP['dr']
+dP_dchord = dP['dchord']
+dP_dtheta = dP['dtheta']
+dP_dprecurve = dP['dprecurve']
+dP_dpresweep = dP['dpresweep']
 
 
 # 6 ----------
@@ -163,7 +174,7 @@
 dNp_dr_fd = (Npd - Np) / delta
 dTp_dr_fd = (Tpd - Tp) / delta
 
-print '(analytic) dNp_i/dr_i =', dNp_dr[idx]
+print '(analytic) dNp_i/dr_i =', dNp_dr[idx, idx]
 print '(fin diff) dNp_i/dr_i =', dNp_dr_fd[idx]
 print
 # 7 ----------
@@ -185,7 +196,7 @@
 dQ_dprecone_fd = (Qd - Q) / delta
 dP_dprecone_fd = (Pd - P) / delta
 
-print '(analytic) dP/dprecone =', dP_dprecone
+print '(analytic) dP/dprecone =', dP_dprecone[0, 0]
 print '(fin diff) dP/dprecone =', dP_dprecone_fd[0]
 print
 # 8 ----------
@@ -208,7 +219,7 @@
 dQ_dr_fd = (Qd - Q) / delta
 dP_dr_fd = (Pd - P) / delta
 
-print '(analytic) dP/dr_i =', dP_dr[idx]
+print '(analytic) dP/dr_i =', dP_dr[0, idx]
 print '(fin diff) dP/dr_i =', dP_dr_fd[0]
 print
 # 9 ----------

docs/tutorial.rst

@@ -131,15 +131,13 @@ CCBlade optinally provides analytic gradients of every output with respect to al
     :start-after: # 3 ---
     :end-before: # 3 ---
 
-Now when we ask for the distributed loads, we also get the gradients.  The gradients are returned as a 2D array and unpacked as follows.
+Now when we ask for the distributed loads, we also get the gradients.  The gradients are returned as a dictionary containing 2D arrays.  These can be accessed as follows:
 
 .. literalinclude:: examples/gradients.py
     :start-after: # 5 ---
     :end-before: # 5 ---
 
-For each vector variable (with the exception of precurve) a change in the local quantity affects only the local normal load.  In other words: :math:`dNp_i/dr_j = 0 \text{ for all } i \neq j`.  Thus only the diagonal terms are returned in a 1D vector (``dNp_dr[3] =`` :math:`dNp_3 / dr_3`).  Precurve affects the sections adjacent to it, so the returned matrix is tridiagonal.
-
-We can compare against finite differencing as follows (with a randomly chosen station along the blade):
+Even though many of the matrices are diagonal, the full Jacobian is returned for consistency.  We can compare against finite differencing as follows (with a randomly chosen station along the blade):
 
 .. literalinclude:: examples/gradients.py
     :start-after: # 7 ---
@@ -150,7 +148,7 @@ The output is:
 >>> (analytic) dNp_i/dr_i = 107.680395098
 >>> (fin diff) dNp_i/dr_i = 107.680370762
 
-Similarly, when we compute thrust, torque, and power we also get the gradients (for either the nondimensional or dimensional form).  The gradients are returned as a multidimensional arrays.  The first index corresponds to the different input cases.  In the following example, we use only one input case (Uinf is a 1D array) so the first index is 0.  The derivatives with respect to scalar quantities and with respect to vector quantities are returned separately.
+Similarly, when we compute thrust, torque, and power we also get the gradients (for either the nondimensional or dimensional form).  The gradients are also returned as a dictionary containing 2D Jacobians.
 
 .. literalinclude:: examples/gradients.py
     :start-after: # 6 ---
@@ -174,3 +172,5 @@ Finally, we compare the derivative of power against finite differencing for one
 
 >>> (analytic) dP/dr_i = 848.368037506
 >>> (fin diff) dP/dr_i = 848.355994992
+
+For more comprehensive comparison to finite differencing, see the unit tests contained in test/test_gradients.py.

setup.py

@@ -8,7 +8,7 @@
 
 setup(
     name='CCBlade',
-    version='1.0.2',
+    version='1.1.0',
     description='Blade element momentum aerodynamics for wind turbines',
     author='S. Andrew Ning',
     author_email='andrew.ning@nrel.gov',