diff --git a/src/fusedwind/examples/turbine/cs_area_example.py b/src/fusedwind/examples/turbine/cs_area_example.py new file mode 100644 index 0000000..7226ec2 --- /dev/null +++ b/src/fusedwind/examples/turbine/cs_area_example.py @@ -0,0 +1,116 @@ + +# --- 1 ----- + +from fusedwind.turbine.structure_vt import CrossSectionAreasVT + +#top = Assembly() + + + +cs_areas = CrossSectionAreasVT() + +uniax = cs_areas.add_material('uniax') +uniax.E1 = 41.63e9 +uniax.E2 = 14.93e9 +uniax.E3 = 14.93e9 +uniax.nu12 = 0.241 +uniax.nu13 = 0.241 +uniax.nu23 = 0.241 +uniax.G12 = 5.047e9 +uniax.G13 = 5.047e9 +uniax.G23 = 5.047e9 +uniax.rho = 1915.5 + +biax = cs_areas.add_material('biax') +biax.E1 = 13.92e9 +biax.E2 = 13.92e9 +biax.E3 = 13.92e9 +biax.nu12 = 0.533 +biax.nu13 = 0.533 +biax.nu23 = 0.533 +biax.G12 = 11.5e9 +biax.G13 = 4.539e9 +biax.G23 = 4.539e9 +biax.rho = 1845 + + +# add keypoints +KP000 = cs_areas.add_kp('KP000') +KP000.x = -0.1 +KP000.y = 0.1 +KP001 = cs_areas.add_kp('KP001') +KP001.x = 0.1 +KP001.y = 0.1 +KP002 = cs_areas.add_kp('KP002') +KP002.x = -0.1 +KP002.y = 0.05 +KP003 = cs_areas.add_kp('KP003') +KP003.x = -0.1 +KP003.y = 0.05 +KP004 = cs_areas.add_kp('KP004') +KP004.x = -0.05 +KP004.y = 0.05 +KP005 = cs_areas.add_kp('KP005') +KP005.x = +0.05 +KP005.y = 0.05 +KP006 = cs_areas.add_kp('KP006') +KP006.x = -0.1 +KP006.y = -0.05 +KP007 = cs_areas.add_kp('KP007') +KP007.x = +0.1 +KP007.y = -0.05 +KP008 = cs_areas.add_kp('KP008') +KP008.x = +0.1 +KP008.y = -0.1 +KP009 = cs_areas.add_kp('KP009') +KP009.x = -0.1 +KP009.y = -0.1 +KP010 = cs_areas.add_kp('KP010') +KP010.x = -0.05 +KP010.y = -0.05 +KP011 = cs_areas.add_kp('KP011') +KP011.x = +0.05 +KP011.y = -0.05 + +# add areas +cap_up = cs_areas.add_area('cap_up') +cap_up.KPs = ['KP003', + 'KP002', + 'KP001', + 'KP000' + ] +cap_up.mat = 'uniax' +cap_up.fiber_plane_angle = .0 +cap_up.fiber_dir_angle = .0 + +cap_low = cs_areas.add_area('cap_low') +cap_low.KPs = ['KP009', + 'KP009', + 'KP007', + 'KP006' + ] +cap_low.mat = 'uniax' +cap_low.fiber_plane_angle = .0 +cap_low.fiber_dir_angle = .0 + +web = cs_areas.add_area('web') +web.KPs = ['KP010', + 'KP011', + 'KP005', + 'KP004' + ] +web.mat = 'biax' +web.fiber_plane_angle = -90.0 +web.fiber_dir_angle = -90.0 + +# HOWTO access single list object +kp = getattr(cs_areas, 'KP005') +print kp.x +print kp.y + +# HOWTO access list obejct in for loop +for kpname in cs_areas.KPs: + print kpname + kp = getattr(cs_areas, kpname) + print kp.x + print kp.y \ No newline at end of file diff --git a/src/fusedwind/turbine/blade_structure.py b/src/fusedwind/turbine/blade_structure.py index 8f702b5..0f6dc5f 100644 --- a/src/fusedwind/turbine/blade_structure.py +++ b/src/fusedwind/turbine/blade_structure.py @@ -9,7 +9,9 @@ from fusedwind.turbine.geometry_vt import BladeSurfaceVT, BladePlanformVT, Curve, AirfoilShape from fusedwind.turbine.geometry import RedistributedBladePlanform, SplineComponentBase, FFDSplineComponentBase -from fusedwind.turbine.structure_vt import BladeStructureVT3D, CrossSectionStructureVT, BeamStructureVT +from fusedwind.turbine.structure_vt import BladeStructureVT3D, CrossSectionStructureVT, BeamStructureVT,\ + CrossSectionMeshVT, ResultVectorArray, CrossSectionElementStressRecoveryVT,\ + CrossSectionAreasVT, KeyPointsVT, CrossSectionAreasVT3D, MeshProps from fusedwind.turbine.rotoraero_vt import LoadVectorCaseList from fusedwind.interface import base, implement_base from fusedwind.lib.geom_tools import curvature @@ -37,6 +39,14 @@ class ModifyBladeStructureBase(Component): st3dOut = VarTree(BladeStructureVT3D(), iotype='out', desc='Vartree containing modified discrete definition of blade structure') +@base +class ModifyCrossSectionAreasBase(Component): + + areas3dIn = VarTree(CrossSectionAreasVT3D(), iotype='in', + desc='Vartree containing initial discrete definition of blade cross section areas') + areas3dOut = VarTree(CrossSectionAreasVT3D(), iotype='out', + desc='Vartree containing modified discrete definition of blade cross section areas') + @implement_base(BladeStructureReaderBase) class BladeStructureReader(Component): @@ -651,7 +661,6 @@ def _post_execute(self): for layer in region.layers: region.thickness += np.maximum(0., getattr(region, layer).thickness) - class BladeStructureProperties(Component): surface = VarTree(BladeSurfaceVT(), iotype='in', desc='Stacked blade surface object') @@ -855,6 +864,99 @@ def execute(self): self.cs2d.append(st2d) +@base +class CrossSectionAreasCSBuilder(Component): + """ + Class that generates a series of 2D cross-sectional area + vartrees (CrossSectionAreasVT3D) used by a mesher like FEPROC + """ + areas3d = VarTree(CrossSectionAreasVT3D(), iotype='in', desc = 'Vartree containing initial discrete definition of blade cross section areas') + areas2d_names = List(iotype='out', desc = 'List of Area names') + + areas2d = List(iotype='out', desc='List of cross-sectional area' + 'vartrees') + + def execute(self): + """ + generate cross sections at every spanwise node of the areas3d vartree + """ + # clear list of outputs! + self.areas2d = [] + self.areas2d_names = self.areas3d.areas + + ni = self.areas3d.s.shape[0] + for i in range(ni): + + area2d = CrossSectionAreasVT() + area2d.s = self.areas3d.s[i] + for kp_name in self.areas3d.KPs: + kp3d = getattr(self.areas3d, kp_name) + kp2dIn = kp3d.extract_kp2d(i) # extract kp object + kp2dOut = area2d.add_kp(kp_name) # add kp object + kp2dOut.x = kp2dIn.x # copy kp coordinates + kp2dOut.y = kp2dIn.y # copy kp coordinates + area2d.materials = self.areas3d.materials.copy() + area2d.areas = self.areas3d.areas # copy area name list + for area_name in self.areas3d.areas: + area_obj3d = getattr(self.areas3d, area_name) + + area_obj2d = area2d.add_area(area_name) + area_obj2d.add_kps(area_obj3d.KPs) + area_obj2d.mat = area_obj3d.mat + area_obj2d.fiber_plane_angle = area_obj3d.fiber_plane_angle + area_obj2d.fiber_dir_angle = area_obj3d.fiber_dir_angle + + self.areas2d.append(area2d) + + + +@base +class BeamMeshCSCode(Component): + """ + Base class for a code such as FEPROC computing beam cross-sectional meshes + suitable for 2D FE codes such as PreComp, BECAS or VABS. + """ + mesh_props = VarTree(MeshProps(), iotype='in', desc = 'mesh properties') + + areas2d = List(iotype='in', desc='List of cross-sectional area' + 'vartrees (CrossSectionAreasVT)') + + meshes2d = List(iotype='out',desc='List of cross-sectional mesh' + 'vartrees (CrossSectionMeshVT)') + + + +@base +class BeamStructureCSCodeWOMesher(Component): + """ + Base class for computing beam structural properties using a cross-sectional + code such as PreComp, BECAS or VABS. + + Complementary to BeamStructureCSCode + """ + + cs_mesh_3D = List(CrossSectionMeshVT, desc='List of CrossSectionMeshVTs', iotype='in') + + beam_structure = List(BeamStructureVT, iotype='out', desc='List Structural beam properties') + + + +@base +class StressRecoveryCSStressStrain(Component): + """ + Base class for performing cross sectional stress/strain analysis + using codes like BECAS and VABS. + + This analysis will typically be in a workflow preceeded by + a call to a BeamStructureCSCode or BeamStructureCSCodeWOMesher. It is assumed that the list of + LoadVectorCaseList vartrees are interpolated onto the structural grid. + """ + + load_cases = List(LoadVectorCaseList, iotype='in', + desc='List of lists of section load vectors for each radial section' + 'used to perform stress/strain analysis') + + cs_res_3D = List(CrossSectionElementStressRecoveryVT, desc='List of element result dictionary VTs', iotype='out') @base class BeamStructureCSCode(Component): diff --git a/src/fusedwind/turbine/geometry_vt.py b/src/fusedwind/turbine/geometry_vt.py index 1b6ffce..dce786c 100644 --- a/src/fusedwind/turbine/geometry_vt.py +++ b/src/fusedwind/turbine/geometry_vt.py @@ -72,7 +72,7 @@ def redistribute(self, dist=None, s=None): points[:, i] = self._splines[i](self.s) self.initialize(points) - + def interp_s(self, s): """ interpolate (x,y) at some curve fraction s diff --git a/src/fusedwind/turbine/rotoraero_vt.py b/src/fusedwind/turbine/rotoraero_vt.py index 0b83fb1..11d9736 100644 --- a/src/fusedwind/turbine/rotoraero_vt.py +++ b/src/fusedwind/turbine/rotoraero_vt.py @@ -302,6 +302,7 @@ class LoadVectorArrayCaseList(VariableTree): List of load vector cases as function of span """ cases = List(LoadVectorArray, desc='List of load cases') + # maro: why is this not a list of load case names according to the convention i.e.in structure_vt/BladeStructureVT3D def _interp_s(self, s): """ diff --git a/src/fusedwind/turbine/structure_vt.py b/src/fusedwind/turbine/structure_vt.py index 594f782..335b0d8 100644 --- a/src/fusedwind/turbine/structure_vt.py +++ b/src/fusedwind/turbine/structure_vt.py @@ -1,4 +1,4 @@ - +import numpy as np from openmdao.main.api import VariableTree from openmdao.lib.datatypes.api import Int, Float, Array, List, Str, Enum, Bool, VarTree, Slot, Dict @@ -147,7 +147,7 @@ def compute_thickness(self): @base -class CrossSectionStructureVT(VariableTree): +class CrossSectionStructureVT(VariableTree): # more appropriate name: CrossSectionLayupVT """ Container for a cross-sectional definition of the internal structure of a blade. @@ -218,7 +218,7 @@ def add_layer(self, name): @base -class BladeStructureVT3D(VariableTree): +class BladeStructureVT3D(VariableTree): # more appropriate name: BladeLayupVT3D """ Variable tree for the structural definition of a blade. """ @@ -240,7 +240,7 @@ def configure_regions(self, nr, names=[]): for i in range(nr): try: - name = rnames[i] + name = names[i] except: name = 'region%02d' % i self.add_region(name) @@ -251,7 +251,7 @@ def configure_webs(self, nw, iwebs, names=[]): for i in range(nw): try: - name = wnames[i] + name = names[i] except: name = 'web%02d' % i self.add_web(name) @@ -290,3 +290,209 @@ def add_material(self, name): self.add(name, VarTree(mat)) self.materials[name] = getattr(self, name) return getattr(self, name) + + +@base +class ResultVectorArray(VariableTree): + ''' + Container for element results + ''' + id = Str(desc = 'Result identifier') + comp_11 = Array(desc = 'Result component 11') + comp_22 = Array(desc = 'Result component 22') + comp_33 = Array(desc = 'Result component 33') + comp_12 = Array(desc = 'Result component 12') + comp_13 = Array(desc = 'Result component 13') + comp_23 = Array(desc = 'Result component 23') + + def _fromarray(self, d): + + self.comp_11 = d[:,0] + self.comp_22 = d[:,1] + self.comp_33 = d[:,2] + self.comp_12 = d[:,3] + self.comp_13 = d[:,4] + self.comp_23 = d[:,5] + + def _toarray(self): + + return np.array([self.comp_11, self.comp_22, self.comp_33, self.comp_12, self.comp_13, + self.comp_23]).T + +@base +class KeyPointsVT(VariableTree): + kp_name = Str('cs_xx',desc='Pointer for key point array') + kp_coords = Array([[0.,0.],[1.,2.]],desc='Array containing arrays of key point coordinates') + +@base +class KeyPoint(VariableTree): + x = Float(desc = 'x coordinates of cs keypoints as function of span', units = 'm') + y = Float(desc = 'y coordinates of cs keypoints as function of span', units = 'm') + + +@base +class Area(VariableTree): + ''' + Area of a cross section. + ''' + #no_KPs,ET,MAT,ASYS,REAL,CUT_flag + KPs = List(desc = 'List of keypoint names') + kp_ni = Int(desc='Number of area keypoints') + mat = Str(desc='Material name') + fiber_plane_angle = Float(units='deg', desc='Material (fiber) plane orientation angle (BECAS alpha)') + fiber_dir_angle = Float(units='deg', desc='Material (fiber) direction angle (BECAS beta)') + #csys_name = Str(desc='Coordinate system object name') + + def add_kps(self, kps): + """ + add a list of keypoints to area + + parameters + ----------- + kps: list + list of kp names + """ + self.KPs = kps + self.kp_ni = len(kps) + +@base +class CrossSectionAreasVT(VariableTree): + """ + Container for a cross section discretized by areas. An area is represented by + at least 3 keypoints. + """ + s = Float(desc = 'position along span') + KPs = List(desc = 'List of keypoint names') + areas = List(desc = 'List of Area names') + materials = Dict(desc='Dictionary of MaterialProps vartrees') + + def add_kp(self, name): + self.add(name, VarTree(KeyPoint())) + self.KPs.append(name) + return getattr(self, name) + + def add_material(self, name): + mat = MaterialProps() + mat.materialname = name + self.add(name, VarTree(mat)) + self.materials[name] = getattr(self, name) + return getattr(self, name) + + def add_area(self, name): + self.add(name, VarTree(Area())) + self.areas.append(name) + return getattr(self, name) + + def kp_coords(self): + ''' + Returns the keypoint coordinates as np array [[x0, y0], [x1, y1], ...] + ''' + kp_coords = [] + for kp_name in self.KPs: + kp_obj = getattr(self, kp_name) + kp_coords.append([kp_obj.x, kp_obj.y]) + return np.r_[kp_coords] + +@base +class KeyPoint3D(VariableTree): + x = Array(desc = 'x coordinates of cs keypoints as function of span', units = 'm') + y = Array(desc = 'y coordinates of cs keypoints as function of span', units = 'm') + + def append_kp2d(self, kp): + ''' + appends a 2D keypoint to KeyPoint3D + + parameters + ----------- + kp: KeyPoint + keypoint + ''' + # most efficient method to append items to a numpy array + x = self.x.tolist() + y = self.y.tolist() + x.append(kp.x) + y.append(kp.y) + self.x = np.r_[x] + self.y = np.r_[y] + + def extract_kp2d(self, kp_index): + ''' + extracts a 2D keypoint from KeyPoint3D + + parameters + ----------- + kp_index: int + index of the key point to be extracted + ''' + kp2d = VarTree(KeyPoint()) + kp2d.x = self.x[kp_index] + kp2d.y = self.y[kp_index] + return kp2d + +@base +class CrossSectionAreasVT3D(VariableTree): + """ + Container for a cross section discretized by areas as function of span. + """ + + s = Array(desc = 'spanwise discretization of blade') + KPs = List(desc = 'List of keypoint names') + materials = Dict(desc = 'Dictionary of MaterialProps vartrees') + areas = List(desc = 'List of Area names') + + def add_kp(self, name): + self.add(name, VarTree(KeyPoint3D())) + self.KPs.append(name) + return getattr(self, name) + + def add_material(self, name): + mat = MaterialProps() + mat.materialname = name + self.add(name, VarTree(mat)) + self.materials[name] = getattr(self, name) + return getattr(self, name) + + def add_area(self, name): + self.add(name, VarTree(Area())) + self.areas.append(name) + return getattr(self, name) + + + +@base +class MeshProps(VariableTree): + ''' + Mesh properties for a mesher + ''' + etype = Str(desc='Element type') + + +@base +class Elset(VariableTree): + ''' + Element set container + ''' + el_numbers = Array(desc = 'Element numbers belonging to the set' ) + +@base +class CrossSectionMeshVT(VariableTree): + ''' + Container for a 2D cross sectional mesh. + Note: Is this too BECAS specific? + ''' + s = Float(desc= 'psoition of mesh along span' ) + nl_2d = Array(desc='Nodal points (node nr, x, y, z)') + defs = '(Element number, node 1, n2, n3, ..., n8)' + el_2d = Array(desc='Elements %s' % defs) + defs = '(Element nr, material nr, fiber angle, fiberplane angle)' + emat = Array(desc='Material per element %s' % defs) + matprops = Array(desc='Material properties (see docs)') + elsets = Dict(desc='Dictionary of Elset vartrees') + +@base +class CrossSectionElementStressRecoveryVT(VariableTree): + ''' + Container for cross sectional element results returned from any FE code. + ''' + el_stresses = Dict(desc='element stresses') # List of ResultVectors + el_strains = Dict(desc='element strains') # List of ResultVectors \ No newline at end of file