Skip to content
Snippets Groups Projects
Commit 4d431cc6 authored by Mads M. Pedersen's avatar Mads M. Pedersen
Browse files

some stuff

parent 1772f3e3
No related branches found
No related tags found
No related merge requests found
...@@ -14,6 +14,7 @@ classifiers = Development Status :: 4 - Beta, ...@@ -14,6 +14,7 @@ classifiers = Development Status :: 4 - Beta,
Programming Language :: Python :: 3, Programming Language :: Python :: 3,
Programming Language :: Python :: 3.3, Programming Language :: Python :: 3.3,
Programming Language :: Python :: 3.4, Programming Language :: Python :: 3.4,
Programming Language :: Python :: 3.5,
Environment :: Console, Environment :: Console,
Intended Audience :: Education, Intended Audience :: Education,
Intended Audience :: Science/Research, Intended Audience :: Science/Research,
......
...@@ -4,9 +4,9 @@ from __future__ import division ...@@ -4,9 +4,9 @@ from __future__ import division
from __future__ import absolute_import from __future__ import absolute_import
from future import standard_library from future import standard_library
standard_library.install_aliases() standard_library.install_aliases()
import pkg_resources test = "TEST"
test = "TEST"
try: try:
__version__ = pkg_resources.get_distribution(__name__).version import pkg_resources
__version__ = pkg_resources.safe_version(pkg_resources.get_distribution(__name__).version)
except: except:
__version__ = 'unknown' __version__ = 'unknown'
...@@ -85,7 +85,7 @@ class Simulation(object): ...@@ -85,7 +85,7 @@ class Simulation(object):
if not os.path.isabs(htcfilename): if not os.path.isabs(htcfilename):
htcfilename = os.path.join(modelpath, htcfilename) htcfilename = os.path.join(modelpath, htcfilename)
self.filename = os.path.basename(htcfilename) self.filename = os.path.basename(htcfilename)
self.htcFile = H2aeroHTCFile(htcfilename) self.htcFile = HTCFile(htcfilename)
self.time_stop = self.htcFile.simulation.time_stop[0] self.time_stop = self.htcFile.simulation.time_stop[0]
self.hawc2exe = hawc2exe self.hawc2exe = hawc2exe
self.copy_turbulence = copy_turbulence self.copy_turbulence = copy_turbulence
......
...@@ -22,8 +22,8 @@ class StFile(object): ...@@ -22,8 +22,8 @@ class StFile(object):
- r : curved length distance from main_body node 1 [m] - r : curved length distance from main_body node 1 [m]
- m : mass per unit length [kg/m] - m : mass per unit length [kg/m]
- xm : xc2-coordinate from C1/2 to mass center [m] - x_cg : xc2-coordinate from C1/2 to mass center [m]
- ym : yc2-coordinate from C1/2 to mass center [m] - y_cg : yc2-coordinate from C1/2 to mass center [m]
- rix : radius of gyration related to elastic center. Corresponds to rotation about principal bending xe axis [m] - rix : radius of gyration related to elastic center. Corresponds to rotation about principal bending xe axis [m]
- riy : radius of gyration related to elastic center. Corresponds to rotation about principal bending ye axis [m] - riy : radius of gyration related to elastic center. Corresponds to rotation about principal bending ye axis [m]
- xs : xc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion. - xs : xc2-coordinate from C1/2 to shear center [m]. The shear center is the point where external forces only contributes to pure bending and no torsion.
......
...@@ -8,6 +8,7 @@ from __future__ import print_function ...@@ -8,6 +8,7 @@ from __future__ import print_function
from __future__ import division from __future__ import division
from __future__ import absolute_import from __future__ import absolute_import
from future import standard_library from future import standard_library
import sys
standard_library.install_aliases() standard_library.install_aliases()
import inspect import inspect
...@@ -80,7 +81,11 @@ def cache_function(f): ...@@ -80,7 +81,11 @@ def cache_function(f):
self.cache_attr_lst.add(name) self.cache_attr_lst.add(name)
return getattr(self, name) return getattr(self, name)
if 'reload' in inspect.getargspec(f)[0]: version = sys.version_info
if version >= (3,3):
if 'reload' in inspect.signature(f).parameters.values():
raise AttributeError("Functions decorated with cache_function are not allowed to take a parameter called 'reload'")
elif 'reload' in inspect.getargspec(f)[0]:
raise AttributeError("Functions decorated with cache_function are not allowed to take a parameter called 'reload'") raise AttributeError("Functions decorated with cache_function are not allowed to take a parameter called 'reload'")
return wrap return wrap
...@@ -151,7 +151,7 @@ def cython_import(import_module_name, compiler=None): ...@@ -151,7 +151,7 @@ def cython_import(import_module_name, compiler=None):
fid.close() fid.close()
# compile, import compiled module and delete temporary files # compile, import compiled module and delete temporary files
module_relname = os.path.relpath(eval(module_name).__file__, str(os.getcwd())).replace(os.path.sep, ".")[:-3] module_relname = eval(module_name).__package__ + "." + module_name #"os.path.relpath(eval(module_name).__file__, str(os.getcwd())).replace(os.path.sep, ".")[:-3]
return compile_and_cleanup(import_module_name, pyx_filename, module_relname, compiler) return compile_and_cleanup(import_module_name, pyx_filename, module_relname, compiler)
return eval(module_name) return eval(module_name)
......
...@@ -162,7 +162,7 @@ def rpm2rads(rpm): ...@@ -162,7 +162,7 @@ def rpm2rads(rpm):
return rpm * 2 * np.pi / 60 return rpm * 2 * np.pi / 60
def abvrel2xyz(alpha, beta, vrel): def abvrel2xyz_old(alpha, beta, vrel):
"""Convert pitot tube alpha, beta and relative velocity to local Cartesian wind speed velocities """Convert pitot tube alpha, beta and relative velocity to local Cartesian wind speed velocities
Parameters Parameters
...@@ -203,3 +203,57 @@ def abvrel2xyz(alpha, beta, vrel): ...@@ -203,3 +203,57 @@ def abvrel2xyz(alpha, beta, vrel):
z[m] = sign_vsz[m] * np.sqrt(vrel[m] ** 2 / ((1 / np.tan(beta[m])) ** 2 + 1 + (np.tan(alpha[m]) / np.tan(beta[m])) ** 2)) z[m] = sign_vsz[m] * np.sqrt(vrel[m] ** 2 / ((1 / np.tan(beta[m])) ** 2 + 1 + (np.tan(alpha[m]) / np.tan(beta[m])) ** 2))
return x, y, z return x, y, z
def abvrel2xyz(alpha, beta, vrel):
"""Convert pitot tube alpha, beta and relative velocity to local Cartesian wind speed velocities
x : parallel to pitot tube, direction pitot tube root to tip, i.e. normal flow gives negative x\n
y : component in alpha plane
z : component in beta plane
For typical usage where pitot tube is mounted on leading edge:\n
x: Opposite rotational direction\n
y: Direction of mean wind\n
z: From blade root to tip\n
Parameters
----------
alpha : array_like
Pitot tube angle of attack [rad]. Zero for flow towards pitot tube. Positive around z-axis. I.e.
negative alpha (normal flow) gives positive y component
beta : array_like
Pitot tube side slip angle [rad]. Zero for flow towards pitot tube. Positive around y-axis. I.e.
Positive beta (normal flow due to expansion and position in front of blade) gives positive z
vrel : array_like
Pitot tube relative velocity. Positive: flow towards pitot tube
Returns
-------
x : array_like
Wind component away from pitot tube (positive for postive vrel and -90<beta<90)
y : array_like
Wind component in alpha plane (positive for positive alpha)
z : array_like
Wind component in beta plane (positive for negative beta)
"""
alpha = np.array(alpha, dtype=np.float)
beta = np.array(beta, dtype=np.float)
vrel = np.array(vrel, dtype=np.float)
sign_vsx = ((np.abs(beta) > np.pi / 2) * 2 - 1) # -1 for |beta| < 90, +1 for |beta|>90
sign_vsy = -np.sign(alpha) #- for alpha > 0
sign_vsz = np.sign(beta) # for beta>0
x = sign_vsx * np.sqrt(vrel ** 2 / (1 + np.tan(alpha) ** 2 + np.tan(beta) ** 2))
m = alpha != 0
y = np.zeros_like(alpha)
y[m] = sign_vsy[m] * np.sqrt(vrel[m] ** 2 / ((1 / np.tan(alpha[m])) ** 2 + 1 + (np.tan(beta[m]) / np.tan(alpha[m])) ** 2))
m = beta != 0
z = np.zeros_like(alpha)
z[m] = sign_vsz[m] * np.sqrt(vrel[m] ** 2 / ((1 / np.tan(beta[m])) ** 2 + 1 + (np.tan(alpha[m]) / np.tan(beta[m])) ** 2))
return np.array([x, y, z]).T
...@@ -14,7 +14,7 @@ import unittest ...@@ -14,7 +14,7 @@ import unittest
import wetb.gtsdf import wetb.gtsdf
import numpy as np import numpy as np
from wetb.utils.geometry import rad, deg, mean_deg, sind, cosd, std_deg, xyz2uvw, \ from wetb.utils.geometry import rad, deg, mean_deg, sind, cosd, std_deg, xyz2uvw, \
wsp_dir2uv, wsp_dir_tilt2uvw, tand wsp_dir2uv, wsp_dir_tilt2uvw, tand, abvrel2xyz
import os import os
...@@ -117,6 +117,15 @@ class TestGeometry(unittest.TestCase): ...@@ -117,6 +117,15 @@ class TestGeometry(unittest.TestCase):
np.testing.assert_array_almost_equal(wsp_dir_tilt2uvw(wsp, dir, tilt, wsp_horizontal=True, dir_ref=180), np.array([x, -y, z]), 5) np.testing.assert_array_almost_equal(wsp_dir_tilt2uvw(wsp, dir, tilt, wsp_horizontal=True, dir_ref=180), np.array([x, -y, z]), 5)
np.testing.assert_array_almost_equal(xyz2uvw(*wsp_dir_tilt2uvw(wsp, dir, tilt, wsp_horizontal=True), left_handed=False), xyz2uvw(x, y, z)) np.testing.assert_array_almost_equal(xyz2uvw(*wsp_dir_tilt2uvw(wsp, dir, tilt, wsp_horizontal=True), left_handed=False), xyz2uvw(x, y, z))
def test_abvrel2xyz(self):
abvrel = np.array([(0., 0, 1), (30, 0, 1), (-30, 0, 1), (0, 30, 1), (0, -30, 1), (30, 30, 1), (30, 30, 2)])
abvrel[:, :2] = rad(abvrel[:, :2])
for (x, y, z), (a, b, vrel) in zip(abvrel2xyz(*abvrel.T), abvrel):
#print (deg(a), deg(b), vrel, x, y, z)
self.assertAlmostEqual(a, np.arctan(y / x))
self.assertAlmostEqual(b, np.arctan(-z / x))
self.assertAlmostEqual(vrel, np.sqrt(x ** 2 + y ** 2 + z ** 2))
if __name__ == "__main__": if __name__ == "__main__":
#import sys;sys.argv = ['', 'Test.test_rad'] #import sys;sys.argv = ['', 'Test.test_rad']
......
...@@ -19,7 +19,7 @@ def _z_u(z_u_lst): ...@@ -19,7 +19,7 @@ def _z_u(z_u_lst):
u = np.array([np.mean(np.array([u])[:]) for _, u in z_u_lst]) u = np.array([np.mean(np.array([u])[:]) for _, u in z_u_lst])
return z, u return z, u
def power_shear(alpha, z_ref, u_ref, z): def power_shear(alpha, z_ref, u_ref):
"""Power shear """Power shear
Parameters Parameters
...@@ -28,22 +28,20 @@ def power_shear(alpha, z_ref, u_ref, z): ...@@ -28,22 +28,20 @@ def power_shear(alpha, z_ref, u_ref, z):
The alpha shear parameter The alpha shear parameter
z_ref : int or float z_ref : int or float
The reference height The reference height
z : int, float or array_like
The heights of interest
u_ref : int or float u_ref : int or float
The wind speed of the reference height The wind speed of the reference height
Returns Returns
------- -------
power_shear : array-like power_shear : function
Wind speeds at height z Function returning for wsp at input heights: f(height) -> wsp
Example Example
-------- --------
>>> power_shear(.5, 70, [20,50,70], 9) >>> power_shear(.5, 70, 9)([20,50,70])
[ 4.81070235 7.60638829 9. ] [ 4.81070235 7.60638829 9. ]
""" """
return u_ref * (np.array(z) / z_ref) ** alpha return lambda z : u_ref * (np.array(z) / z_ref) ** alpha
def fit_power_shear(z_u_lst): def fit_power_shear(z_u_lst):
...@@ -72,7 +70,7 @@ def fit_power_shear(z_u_lst): ...@@ -72,7 +70,7 @@ def fit_power_shear(z_u_lst):
alpha, _ = np.polyfit(np.log(z / z_hub), np.log(u / u_hub), 1) alpha, _ = np.polyfit(np.log(z / z_hub), np.log(u / u_hub), 1)
return alpha return alpha
def fit_power_shear_ref(z_u_lst, z_ref): def fit_power_shear_ref(z_u_lst, z_ref, plt=None):
"""Estimate power shear parameter, alpha, from two or more specific reference heights using polynomial fit. """Estimate power shear parameter, alpha, from two or more specific reference heights using polynomial fit.
Parameters Parameters
...@@ -84,6 +82,8 @@ def fit_power_shear_ref(z_u_lst, z_ref): ...@@ -84,6 +82,8 @@ def fit_power_shear_ref(z_u_lst, z_ref):
- u_z2: Wind speeds or mean wind speeds at z2 - u_z2: Wind speeds or mean wind speeds at z2
z_ref : float or int z_ref : float or int
Reference height (hub height) Reference height (hub height)
plt : matplotlib.pyplot (or similar) or None
Used to plot result if not None
Returns Returns
------- -------
...@@ -101,7 +101,13 @@ def fit_power_shear_ref(z_u_lst, z_ref): ...@@ -101,7 +101,13 @@ def fit_power_shear_ref(z_u_lst, z_ref):
alpha, u_ref = x alpha, u_ref = x
return np.sum([(u - u_ref * (z / z_ref) ** alpha) ** 2 for z, u in z_u_lst]) return np.sum([(u - u_ref * (z / z_ref) ** alpha) ** 2 for z, u in z_u_lst])
z_u_lst = [(z, np.mean(u)) for z, u in z_u_lst] z_u_lst = [(z, np.mean(u)) for z, u in z_u_lst]
return fmin(shear_error, (.1, 10), (z_u_lst, z_ref), disp=False) alpha, u_ref = fmin(shear_error, (.1, 10), (z_u_lst, z_ref), disp=False)
if plt:
z, u = list(zip(*z_u_lst))
plt.plot(u, z, '.')
z = np.linspace(min(z), max(z), 100)
plt.plot(power_shear(alpha, z_ref, u_ref)(z), z)
return alpha, u_ref
...@@ -188,7 +194,7 @@ if __name__ == '__main__': ...@@ -188,7 +194,7 @@ if __name__ == '__main__':
from matplotlib.pyplot import plot, show from matplotlib.pyplot import plot, show
z = np.arange(0, 211) z = np.arange(0, 211)
for alpha, c in zip([0.00001, 1, 2], ['r', 'b', 'g']): for alpha, c in zip([0.00001, 1, 2], ['r', 'b', 'g']):
u = power_shear(alpha, 120, 10, z) u = power_shear(alpha, 120, 10)(z)
plot(u, z, c) plot(u, z, c)
plot(u.mean(), 120, c + '.') plot(u.mean(), 120, c + '.')
......
...@@ -98,8 +98,10 @@ class TestShear(unittest.TestCase): ...@@ -98,8 +98,10 @@ class TestShear(unittest.TestCase):
wsp53 = data[:, 2] wsp53 = data[:, 2]
wsp21 = data[:, 4] wsp21 = data[:, 4]
import matplotlib.pyplot as plt
alpha, u_ref = fit_power_shear_ref([(85, wsp85), (21, wsp21)], 87.13333) alpha, u_ref = fit_power_shear_ref([(85, wsp85), (21, wsp21)], 87.13333, plt)
if 0:
plt.show()
self.assertAlmostEqual(alpha, .5, delta=.001) self.assertAlmostEqual(alpha, .5, delta=.001)
self.assertAlmostEqual(u_ref, 9, delta=.01) self.assertAlmostEqual(u_ref, 9, delta=.01)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment