Commit 3d0805a7 authored by Frederik Zahle's avatar Frederik Zahle
Browse files

Merge branch 'linear_dist' into 'master'

Linear dist and more robust interp_x

See merge request !6
parents c11de8d4 490abb77
......@@ -2,7 +2,7 @@
import numpy as np
import copy
from scipy.optimize import minimize
from scipy.interpolate import pchip, Akima1DInterpolator
from scipy.interpolate import pchip, Akima1DInterpolator, interp1d
from PGL.main.geom_tools import curvature
from PGL.main.curve import Curve
......@@ -76,7 +76,8 @@ class AirfoilShape(Curve):
return ds1
def redistribute(self, ni, even=False, dist=None, dLE=False, dTE=-1., close_te=0):
def redistribute(self, ni, even=False, dist=None, dLE=False, dTE=-1.,
close_te=0, linear=False):
"""
redistribute the points on the airfoil using fusedwind.lib.distfunc
......@@ -129,7 +130,7 @@ class AirfoilShape(Curve):
dist = [[0., minTE, 1],
[self.sLE, dist_LE, main_ni / 2],
[1., minTE, main_ni]]
main_seg.redistribute(dist=dist)
main_seg.redistribute(dist=dist, linear=linear)
points = lte_seg.points.copy()
points = np.append(points, main_seg.points[1:], axis=0)
points = np.append(points, ute_seg.points[1:], axis=0)
......@@ -145,7 +146,7 @@ class AirfoilShape(Curve):
else:
dist = [[0., dTE, 1], [self.sLE, self.leading_edge_dist(ni), ni / 2], [1., dTE, ni]]
super(AirfoilShape, self).redistribute(dist)
super(AirfoilShape, self).redistribute(dist, linear=linear)
return self
......@@ -262,26 +263,41 @@ class AirfoilShape(Curve):
"""
interpolate s(x) for lower or upper side
"""
if self.LE[0] < self.TE[0]:
iLE = np.argmin(self.points[:, 0])
iTEl = np.argmax(self.points[:iLE, 0])
iTEu = np.argmax(self.points[iLE:, 0]) + iLE
if x < self.points[iLE, 0]:
return self.points[iLE, 1]
elif x > self.TE[0]:
return 0.
else:
iLE = np.argmax(self.points[:, 0])
iTEl = np.argmin(self.points[:iLE, 0])
iTEu = np.argmin(self.points[iLE:, 0]) + iLE
if x > self.points[iLE, 0]:
return self.sLE
elif x < self.TE[0]:
return 0.
if side == 'lower':
# interpolate pressure side coordinates
s = self.points[20:self.iLE-20, 0]
s = self.points[iTEl:iLE, 0]
if s[0] > s[-1]:
s = s[::-1]
y = self.s[20:self.iLE-20][::-1]
y = self.s[iTEl:iLE][::-1]
else:
y = self.s[20:self.iLE-20]
y = self.s[iTEl:iLE]
spl = NaturalCubicSpline(s, y)
return spl(x)
if side == 'upper':
# interpolate pressure side coordinates
s = self.points[self.iLE+20:, 0]
s = self.points[iLE:iTEu, 0]
if s[0] > s[-1]:
s = s[::-1]
y = self.s[self.iLE+20:][::-1]
y = self.s[iLE:iTEu][::-1]
else:
y = self.s[self.iLE+20:]
y = self.s[iLE:iTEu]
spl = NaturalCubicSpline(s, y)
return spl(x)
......@@ -333,6 +349,8 @@ class BlendAirfoilShapes(object):
def initialize(self):
if self.spline == 'linear':
self._spline = interp1d
if self.spline == 'pchip':
self._spline = pchip
elif self.spline == 'cubic':
......
......@@ -105,7 +105,7 @@ class Curve(object):
for j in range(self.points.shape[1]):
self._splines.append(interp1d(self.s, self.points[:, j], kind='linear'))
def redistribute(self, dist=None, s=None, ni=100):
def redistribute(self, dist=None, s=None, ni=100, linear=False):
"""
redistribute the points on the curve using distfunc
or a user-supplied distribution
......@@ -129,8 +129,13 @@ class Curve(object):
redistributed evenly using ni points.
"""
if dist is not None:
self.s = distfunc(np.real(dist))
# self.s = distfunc(dist, dist[-1][-1])
ni = dist[-1][-1]
if linear:
self.s = np.zeros(ni)
for d0, d1 in zip(dist[:-1], dist[1:]):
self.s[d0[2]-1:d1[2]] = np.linspace(d0[0], d1[0], d1[2]-d0[2]+1)
else:
self.s = distfunc(np.real(dist))
elif s is not None:
self.s = s
else:
......@@ -521,9 +526,9 @@ class SegmentedCurve(Curve):
_last = seg
self.points = np.append(points,np.array([_last.points[-1,:]]),0)
def redistribute(self, dist=None, s=None, ni=None):
def redistribute(self, dist=None, s=None, ni=None, linear=False):
super(SegmentedCurve,self).redistribute(dist, s, ni)
super(SegmentedCurve,self).redistribute(dist, s, ni, linear)
class Line(Curve):
......
......@@ -8,29 +8,30 @@ class AirfoilShapeTests(unittest.TestCase):
def setUp(self):
unittest.TestCase.setUp(self)
self.af = AirfoilShape()
self.af.initialize(points=np.loadtxt('data/ffaw3480.dat'))
self.af.initialize(points=np.loadtxt('data/cylinder.dat'))
# self.af.initialize(points=np.loadtxt('data/ffaw3480.dat'))
self.dps_s01 = np.array([0.0, 0.25, 0.5, 0.75, 1.0])
def test_s_to_11_10(self):
dps_s11 = np.array([self.af.s_to_11(s) for s in self.dps_s01])
dps_s01n = np.array([self.af.s_to_01(s) for s in dps_s11])
self.assertEqual(np.testing.assert_allclose(dps_s01n, self.dps_s01, 1E-06), None)
def test_s_to_11_10_rotate(self):
afi = self.af.copy()
afi.rotate_z(-45.0)
dps_s11 = np.array([afi.s_to_11(s) for s in self.dps_s01])
afii = self.af.copy()
afii.rotate_z(+45.0)
dps_s01n = np.array([afii.s_to_01(s) for s in dps_s11])
self.assertEqual(np.testing.assert_allclose(dps_s01n, self.dps_s01, 1E-06), None)
def test_s_to_11_10_scale_equal(self):
afi = self.af.copy()
afi.scale(1.5)
points = afi.points
......@@ -39,12 +40,12 @@ class AirfoilShapeTests(unittest.TestCase):
afii.scale(1.5)
pointsn = afii.points
dps_s01n = np.array([afii.s_to_01(s) for s in dps_s11])
self.assertEqual(np.testing.assert_allclose(dps_s01n, self.dps_s01, 1E-06), None)
self.assertEqual(np.testing.assert_allclose(pointsn, points, 1E-06), None)
def test_s_to_11_10_scale_not_equal(self):
afi = self.af.copy()
afi.scale(1.1)
#points = afi.points
......@@ -53,10 +54,45 @@ class AirfoilShapeTests(unittest.TestCase):
afii.scale(1.5)
#pointsn = afii.points
dps_s01n = np.array([afii.s_to_01(s) for s in dps_s11])
self.assertEqual(np.testing.assert_allclose(dps_s01n, self.dps_s01, 1E-06), None)
#self.assertEqual(np.testing.assert_allclose(pointsn, points, 1E-06), None)
class AirfoilShapeTests2(unittest.TestCase):
''' Tests for components.airfoil.py
'''
def test_redistribute_linear(self):
af = AirfoilShape(points=np.loadtxt('data/ffaw3241.dat'))
ds = np.array([ 0.11712651, 0.117116 , 0.11720435, 0.11725923, 0.11719773,
0.11713939, 0.11709662, 0.11696592, 0.11212361, 0.02098937,
0.02108704, 0.02109582, 0.02110256, 0.0211067 , 0.02110718,
0.02110731, 0.02110809, 0.02110783, 0.02110793, 0.08436049,
0.08439066, 0.08441201, 0.0844209 , 0.08442478, 0.08442662,
0.08442745, 0.08442748, 0.08442729, 0.08442751])
af.redistribute(ni=30, dist=[[0, 0.01, 1],
[0.5, 0.01, 10],
[0.6, 0.01, 20],
[1, 0.01, 30]], linear=True)
self.assertEqual(np.testing.assert_allclose(af.ds, ds, 1E-06), None)
def test_redistribute_distfunc(self):
af = AirfoilShape(points=np.loadtxt('data/ffaw3241.dat'))
ds = np.array([ 0.03187173, 0.0662534 , 0.1239282 , 0.19230617, 0.22543051,
0.19175545, 0.1237126 , 0.06599569, 0.03104575, 0.02098937,
0.02108704, 0.02109582, 0.02110256, 0.0211067 , 0.02110718,
0.02110731, 0.02110809, 0.02110783, 0.02110793, 0.02879817,
0.05032088, 0.0816139 , 0.1176486 , 0.1436471 , 0.1436658 ,
0.11769158, 0.0816484 , 0.05033041, 0.02879825])
af.redistribute(ni=30, dist=[[0, 0.01, 1],
[0.5, 0.01, 10],
[0.6, 0.01, 20],
[1, 0.01, 30]])
self.assertEqual(np.testing.assert_allclose(af.ds, ds, 1E-06), None)
if __name__ == '__main__':
unittest.main()
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment