Airfoils¶
µ∩δeR Ɔonsτrµc†ioԥ
In [2]:
import numpy as np
import matplotlib.pyplot as plt
data = np.genfromtxt('NACA0015.txt')
fig,ax=plt.subplots(figsize=(10,10))
ax.set_aspect('equal', adjustable='box')
ax.plot(data[:,0],data[:,1])
ax.grid()
plt.show()
data = np.genfromtxt('NACA0012.txt')
fig,ax=plt.subplots(figsize=(10,10))
ax.set_aspect('equal', adjustable='box')
ax.plot(data[:,0],data[:,1])
ax.grid()
plt.show()
In [ ]:
from abaqus import *
from abaqusConstants import *
import numpy as np
def profile(file, scale, dx, dy, z):
data = np.loadtxt(file)
xs = data[:,0]*scale + dx
ys = data[:,1]*scale + dy
coordinates = []
for x,y in zip(xs,ys):
coordinates.append((x,y,z))
return coordinates
def airfoilB():
import part
modelname = 'M1'
profiles = [{'file':'NACA0015.txt', 'sf':200.0, 'dx': 0.0, 'dy':0.0, 'z': 0.0},
{'file':'NACA0015.txt', 'sf':100.0, 'dx':100.0, 'dy':0.0, 'z':500.0} ]
pxs = (0.1, 0.2, 0.4, 0.6, 0.8)
pzs = (0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9)
mod = mdb.Model(name=modelname)
prt = mod.Part(name='P1', dimensionality=THREE_D, type=DEFORMABLE_BODY)
featureNames = []
for pro in profiles:
crds = profile(pro['file'], pro['sf'], pro['dx'], pro['dy'], pro['z'])
fn1 = prt.WireSpline(points=tuple(crds), meshable=ON, smoothClosedSpline=OFF).name
fn2 = prt.WirePolyLine(points=((crds[0], crds[-1]), ), mergeType=IMPRINT, meshable=ON).name
featureNames.append((fn1,fn2))
loftsections=[]
for fns in featureNames:
pedges = [e for e in prt.edges if e.featureName in fns]
loftsections.append(pedges)
prt.SolidLoft(loftsections=loftsections, startCondition=NONE, endCondition=NONE)
ext = prt.edges.getBoundingBox()
minx, miny, minz = ext['low']
maxx, maxy, maxz = ext['high']
sf1, dx1, z1 = profiles[0]['sf'], profiles[0]['dx'], profiles[0]['z']
sf2, dx2, z2 = profiles[-1]['sf'], profiles[-1]['dx'], profiles[-1]['z']
for px in pxs:
id1 = prt.DatumPointByCoordinate(coords=(sf1*px+dx1, miny-100, z1)).id
id2 = prt.DatumPointByCoordinate(coords=(sf1*px+dx1, maxy+100, z1)).id
id3 = prt.DatumPointByCoordinate(coords=(sf2*px+dx2, 0.0, z2)).id
d=prt.datums
prt.PartitionCellByPlaneThreePoints(point1=d[id1], point2=d[id2], point3=d[id3], cells=prt.cells)
for pz in pzs:
id1 = prt.DatumPointByCoordinate(coords=(minx-100, miny-100, (z2-z1)*pz )).id
id2 = prt.DatumPointByCoordinate(coords=(minx-100, maxy+100, (z2-z1)*pz )).id
id3 = prt.DatumPointByCoordinate(coords=(maxx+100, 0.0, (z2-z1)*pz )).id
d=prt.datums
prt.PartitionCellByPlaneThreePoints(point1=d[id1], point2=d[id2], point3=d[id3], cells=prt.cells)
prt.RemoveCells(cellList = prt.cells)
session.viewports['Viewport: 1'].setValues(displayedObject=prt)
airfoilB()
