Skip to content
Snippets Groups Projects
Commit 025f62ca authored by Malo Rosemeier's avatar Malo Rosemeier
Browse files

wrapper works, verification of cross section props required.

parent 71decd4a
No related branches found
No related tags found
No related merge requests found
import numpy as np
from finstrip_wrapper.core.main import FINSTRIPWrapper, status
from PGL.components.airfoil import AirfoilShape
class FINSTRIPWrapperBlade(FINSTRIPWrapper):
......@@ -12,18 +13,10 @@ class FINSTRIPWrapperBlade(FINSTRIPWrapper):
# data base output name
self.name = 'blade'
# loads dummy
self.loads = np.array([])
#self.aligns = np.array([])
self.core_mats = [] # list of materials to be identified as core
self.spanpos = 0 # position along span of analysed cross-section
# flags
#self.post_rst = True
# standard buckling analysis variables
self.buckle_modes = 10 # number of bucling modes
self.buckle_load_start = 0.0 # start boundary load
self.buckle_load_end = 10.0 # end boundary load
# self.buckle_method = 'LANB' # 'SUBSP'
self.use_becas_stiffness = False
for k, w in kwargs.iteritems():
try:
......@@ -31,56 +24,248 @@ class FINSTRIPWrapperBlade(FINSTRIPWrapper):
except:
pass
self._init_half_wavelength_steps()
self._init_finstrip_input_files_names()
def init_dicts(self, pf):
def init_dicts(self, cs2d):
# dictionaries written by FINSTRIPBladeStructure according to FUSED-Wind format
# NOTE: absolut dimensions are required!
self.pf = pf
def init_full_scale_load_case(self, loads, aligns):
# TODO: here we need the conversion from load distro to point loads
# and the calculation of the load introduction positions
# By now it takes only external point loads
# TODO: remove hardcoded loads
# position of master nodes and force vector,and actuator positions
# x,y,z,rotx,roty,rotz,Fx,Fy,Fz,Mx,My,Mz,ax,ay,az
# actuator positions in the test_tig locinate system
self.loads = loads
self.nloads = self.loads.shape[0]
# tilt, pitch, x_orig, y_orig, z_orig
# tilt_angle rotation about root x-axis in deg
# pitch_angle rotation about root z-axis in deg
# test_rig_origin use distance in y-direction for vertical pull
self.aligns = aligns
self.cs2d = cs2d
def init_becas_stiffness(self, csprops, csprops_ref):
# stiffness properties as described in BECASWrapper
self.csprops = csprops
self.csprops_ref = csprops_ref
def init_cs_load_cases(self, csloads_rcs, twist):
''' csload is assumed to be provided in root cs
:param csloads: cross section loads (Fx,Fy,Fz,Mx,My,Mz) size ((ncase,6))
:param twist: geometric twist angle of cross section w.r.t root coordinate system size ((1))
note: so far only bending moments mx and my are used
'''
Mx_rcs = csloads_rcs[:, 3]
My_rcs = csloads_rcs[:, 4]
# transformation matrix from root cs to cross section cs
sin_twist = np.sin(np.deg2rad(twist))
cos_twist = np.cos(np.deg2rad(twist))
self.T_twist = np.matrix(
[[cos_twist, sin_twist], [-sin_twist, cos_twist]])
# transform from root cs to cross section cs
v_Mrcs = np.array([[Mx_rcs], [My_rcs]])
v_Mcs = self.T_twist * v_Mrcs
M_x = np.squeeze(np.asarray(v_Mcs[:, 0]))
M_y = np.squeeze(np.asarray(v_Mcs[:, 1]))
self.csloads = np.zeros_like(csloads_rcs)
self.csloads[:, 3] = M_x
self.csloads[:, 4] = M_y
self.csloads[:, 3:] = self.csloads[:, 3:] * 1E3 # convert to Nmm
self.ncase = len(self.csloads[:, 0])
def setup(self):
''' Sets up the FINSTRIP execution.
'''
super(FINSTRIPWrapperBlade, self).setup()
self._det_maxlength_element()
self._write_finstrip_input_files()
self._write_finstrip_analysis_files()
def _det_maxlength_element(self):
# calculate arc length of airfoil surface
af = AirfoilShape()
af.initialize(points=self.cs2d['coords'])
# maximum size of generated element in m
self.maxlength_element = af.smax * self.elsize_factor
def _write_finstrip_input_files(self):
super(FINSTRIPWrapperBlade, self)._write_finstrip_input_files()
self._write_file_buc()
def _write_file_buc(self):
pass
out_str = []
n = '\n'
# write header
out_str.append(
'// Finstrip input file, created by FINSTRIPWrapper' + n)
# write vars converted to mm
for var_name, var_value in zip(['maxlength_element',
'lambda_steps',
'lambda_min',
'lambda_max'],
[self.maxlength_element * 1E3,
self.lambda_steps,
self.lambda_min * 1E3,
self.lambda_max * 1E3,
]):
out_str.append('%s, %0.2f' % (var_name, var_value) + n)
# write materials
mat_names = self.cs2d['materials'].keys()
matprops = self.cs2d['matprops']
def _write_iso(out_str, n, label, E1, nu12, rho):
out_str.append('iso, %s, %g, %g, %g' %
(label, E1, nu12, rho) + n)
def _write_ortho(out_str, n, label, E1, E2, nu12, G12, rho):
out_str.append('ortho, %s, %g, %g, %g, %g, %g' % (
label, E1, E2, nu12, G12, rho) + n)
def _write_core(out_str, n, label, E1, E2, nu12, G12, G13, G23, rho):
out_str.append('core, %s, %g, %g, %g, %g, %g,%g, %g' % (
label, E1, E2, nu12, G12, G13, G23, rho) + n)
for i, label in enumerate(mat_names):
E1 = matprops[i, 0] * 1E-6 # to N/mm*2
E2 = matprops[i, 1] * 1E-6 # to N/mm*2
E3 = matprops[i, 2] * 1E-6 # to N/mm*2
nu12 = matprops[i, 3]
nu13 = matprops[i, 4]
nu23 = matprops[i, 5]
G12 = matprops[i, 6] * 1E-6 # to N/mm*2
G13 = matprops[i, 7] * 1E-6 # to N/mm*2
G23 = matprops[i, 8] * 1E-6 # to N/mm*2
rho = matprops[i, 9] * 1E-9 # to kg/mm*3
# deterimne material type
if label in self.core_mats:
# mat_type = 'core'
_write_core(
out_str, n, label, E1, E2, nu12, G12, G13, G23, rho)
elif E1 == E2 and nu12 == nu13:
# mat_type = 'iso'
_write_iso(out_str, n, label, E1, nu12, rho)
else:
# mat_type = 'ortho'
_write_ortho(out_str, n, label, E1, E2, nu12, G12, rho)
# write blade name and span pos
self.spanpos = self.cs2d['s'] * 1E3 # to mm
out_str.append('blade, %s, %g' % (self.name, self.spanpos) + n)
# write first (and only) section
out_str.append('Section, %g' % (self.spanpos) + n)
# write coords
coords = self.cs2d['coords'] * 1E3 # to mm
for i, c in enumerate(coords):
out_str.append('k, %g, %g' % (c[0], c[1]) + n)
# write region layup and dps
for ireg, _ in enumerate(self.cs2d['regions']):
layers = self.cs2d['regions'][ireg]['layers']
thicknesses = self.cs2d['regions'][
ireg]['thicknesses'] * 1E3 # to mm
angles = self.cs2d['regions'][ireg]['angles']
dp0 = self.cs2d['DPs'][ireg]
dp1 = self.cs2d['DPs'][ireg + 1]
for l, t, a in zip(layers, thicknesses, angles):
out_str.append('Layer, %s, %g, %g, %g, %g' %
(l[:-2], t, a, dp0, dp1) + n)
# write web layup and dps
for ireg, _ in enumerate(self.cs2d['webs']):
layers = self.cs2d['webs'][ireg]['layers']
thicknesses = self.cs2d['webs'][
ireg]['thicknesses'] * 1E3 # to mm
angles = self.cs2d['regions'][ireg]['angles']
# check if web exists at this cross section
if not np.sum(thicknesses) == 0.0:
dp0 = self.cs2d['DPs'][self.cs2d['web_def'][ireg][0]]
dp1 = self.cs2d['DPs'][self.cs2d['web_def'][ireg][1]]
out_str.append('Shearweb, %g, %g' % (dp0, dp1) + n)
wdp0 = 0.0 # web layup does not change along width
wdp1 = 1.0
for l, t, a in zip(layers, thicknesses, angles):
out_str.append('Layer, %s, %g, %g, %g, %g' %
(l[:-2], t, a, wdp0, wdp1) + n)
if self.use_becas_stiffness:
x_e = self.csprops_ref[2] * 1E3 # to mm
y_e = self.csprops_ref[3] * 1E3 # to mm
I_xx = self.csprops_ref[7]
I_yy = self.csprops_ref[8]
I_xy = self.csprops_ref[9]
nu = self.csprops_ref[17] + np.pi / 2
# TODO: check if sign is correct in I_xy!
# A = self.csprops_ref[15]
# Gere 12-30
theta = 0.5 * np.arctan(-2 * I_xy / (I_xx - I_yy))
# Gere 12-32
R = np.sqrt((0.5 * (I_xx - I_yy))**2 + I_xy**2)
# Gere 12-31a
theta1 = 0.5 * np.arccos((I_xx - I_yy) / (2 * R))
# Gere 12-31b
theta2 = 0.5 * np.arcsin(-I_xy / R)
# if not np.sign(theta) == np.sign(nu):
# raise RuntimeError('BECAS input not consistent')
config = {}
config['BECASWrapper'] = {}
config['BECASWrapper']['hawc2_FPM'] = True
from becas_wrapper.becas_bladestructure import set_cs_size
cs_size, _ = set_cs_size(config)
if cs_size == len(self.csprops):
#hawc2_FPM = True
AE = self.csprops[20] * 1E-3 # to N/mm
else:
#hawc2_FPM = False
E = self.csprops[8]
An = self.csprops[15]
AE = E * An * 1E-3 # to N/mm
# AE = A * E # to N/mm
EI_xx = E * I_xx * 1E3 # to Nmm**2/mm
EI_yy = E * I_yy * 1E3 # to Nmm**2/mm
EI_xy = E * I_xy * 1E3 # to Nmm**2/mm
out_str.append('farob, %g, %g, %g, %g, %g, %g' %
(x_e, y_e, AE, EI_xx, EI_yy, EI_xy) + n)
fid = open(self.file_input, 'w')
for i, line in enumerate(out_str):
fid.write(line)
fid.close()
def _write_finstrip_analysis_files(self):
self._write_file_loads()
self._write_file_lcd()
def _write_file_lcd(self):
'''
np.savetxt(self.file_loads, self.loads, fmt='%.10e', delimiter=',')
np.savetxt(
self.file_aligns, self.aligns, fmt='%.10e', delimiter=',')
'''
out_str = []
for case in range(self.ncase):
fx = self.csloads[case, 0]
fy = self.csloads[case, 1]
fz = self.csloads[case, 2]
mx = self.csloads[case, 3]
my = self.csloads[case, 4]
mz = self.csloads[case, 5]
n = '\n'
# write header
out_str.append(
'loads, %s, %i, %g, %g, %g, %g, %g' % (self.spanpos, fx, fy, fz, mx, my, mz) + n)
fid = open(self.file_load, 'w')
for i, line in enumerate(out_str):
fid.write(line)
fid.close()
def execute(self):
super(FINSTRIPWrapperBlade, self).execute()
......@@ -91,16 +276,19 @@ class FINSTRIPWrapperBlade(FINSTRIPWrapper):
self._read_rst_mode_shape()
def _read_rst_eigenvalue(self):
# proces_nr, section_nr, load_nr, wavelength, eigenvalue
self.f_eig = np.loadtxt(self.file_rst_eigen, skiprows=1)
# convert wavelength to m
self.f_eig[:, -2] = self.f_eig[:, -2] * 1E-3
self.f_eig = np.zeros((self.ncase, self.nmode, 2))
for case in range(self.ncase):
# proces_nr, section_nr, load_nr, wavelength, eigenvalue
self.f_eig[case, :, :] = np.loadtxt(
self.file_rst_eigen, skiprows=1)[case * self.nmode:(case + 1) * self.nmode, -2:]
# convert wavelength to m
self.f_eig[:, :, 0] = self.f_eig[:, :, 0] * 1E-3
def _read_rst_mode_shape(self):
# determine number of elements
self.nelm = len(np.loadtxt('R' +
'%i' % (self.spanpos) +
'%05i' % (self.spanpos) +
'_LC' +
'%03d' % 0 +
'_mode'
......@@ -108,41 +296,43 @@ class FINSTRIPWrapperBlade(FINSTRIPWrapper):
self.txt_sfx,
skiprows=2)[:, 0])
# get mode shapes
self.mode_shape = np.zeros((self.nmode, self.nelm, 8))
self.mode_shape = np.zeros((self.ncase, self.nmode, self.nelm, 8))
# x-loc p1, y-loc p1, x-loc p2, y-loc p2,
# x-u p1, y-u p1, x-u p2, y-u p2
for case in range(self.ncase):
for mode in range(self.nmode):
# transform to m
self.mode_shape[mode, :, :] = 1E-3 * np.loadtxt('R' +
'%i' % (self.spanpos) +
'_LC' +
'%03d' % (case) +
'_mode'
'%03d' % (mode) +
self.txt_sfx,
skiprows=2)
self.mode_shape[case, mode, :, :] = np.loadtxt('R' +
'%05i' % (self.spanpos) +
'_LC' +
'%03d' % (case) +
'_mode'
'%03d' % (mode) +
self.txt_sfx,
skiprows=2)
idxs = []
# identify surface and web nodes
for elm in range(self.nelm - 1):
elm1_node2 = self.mode_shape[0, elm, 2:4]
elm2_node1 = self.mode_shape[0, elm + 1, :2]
elm1_node2 = self.mode_shape[0, 0, elm, 2:4]
elm2_node1 = self.mode_shape[0, 0, elm + 1, :2]
if elm1_node2[0] != elm2_node1[0] and elm1_node2[1] != elm2_node1[1]:
idxs.append(elm)
# append last elm index
idxs.append(len(self.mode_shape[0, :, 0]) - 1)
idxs.append(len(self.mode_shape[0, 0, :, 0]) - 1)
# obtain surface nodes
self.nnode_loc_surf = idxs[0] + 1 + 1
self.node_loc_surf = np.zeros((self.nnode_loc_surf, 2))
self.nnode_surf = idxs[0] + 1 + 1
self.node_loc_surf = np.zeros((self.nnode_surf, 2))
self.node_loc_surf[
:-1, :] = self.mode_shape[0, :idxs[0] + 1, :2]
:-1, :] = self.mode_shape[0, 0, :idxs[0] + 1, :2]
# add second node of last element
self.node_loc_surf[-1:,
:] = self.mode_shape[0, idxs[0], 2:4]
:] = self.mode_shape[0, 0, idxs[0], 2:4]
# obtain web nodes
# actual amount of webs in post processing files can differ from input
# NOTE: if the flatback web is oriented such that it closes the cross section,
# then it is identified as if the web was part of the surface
self.nweb = len(idxs[1:])
self.nnode_web = []
......@@ -154,29 +344,29 @@ class FINSTRIPWrapperBlade(FINSTRIPWrapper):
(self.nweb, self.nnode_web_max, 2))
for w in range(self.nweb):
self.node_loc_web[
w, :self.nnode_web[w] - 1, :] = self.mode_shape[0, idxs[w] + 1:idxs[w + 1] + 1, :2]
w, :self.nnode_web[w] - 1, :] = self.mode_shape[0, 0, idxs[w] + 1:idxs[w + 1] + 1, :2]
# add second node of last element
self.node_loc_web[w, self.nnode_web[w] - 1, :] = self.mode_shape[
0, idxs[w + 1], 2:4]
0, 0, idxs[w + 1], 2:4]
# obtain surface displacements
self.node_u_surf = np.zeros(
(self.nmode, self.nnode_loc_surf, 2))
(self.ncase, self.nmode, self.nnode_surf, 2))
self.node_u_surf[
:, :-1, :] = self.mode_shape[:, :idxs[0] + 1, 4:6]
:, :, :-1, :] = self.mode_shape[:, :, :idxs[0] + 1, 4:6]
# add second node of last element
self.node_u_surf[
:, -1, :] = self.mode_shape[:, idxs[0], 6:]
:, :, -1, :] = self.mode_shape[:, :, idxs[0], 6:]
# obtain web displacements
self.node_u_web = np.zeros(
(self.nweb, self.nmode, self.nnode_web_max, 2))
(self.ncase, self.nweb, self.nmode, self.nnode_web_max, 2))
for w in range(self.nweb):
self.node_u_web[
w, :, :self.nnode_web[w] - 1, :] = self.mode_shape[:, idxs[w] + 1:idxs[w + 1] + 1, 4:6]
:, w, :, :self.nnode_web[w] - 1, :] = self.mode_shape[:, :, idxs[w] + 1:idxs[w + 1] + 1, 4:6]
# add second node of last element
self.node_u_web[w, :, self.nnode_web[w] - 1, :] = self.mode_shape[
:, idxs[w + 1], 6:]
self.node_u_web[:, w, :, self.nnode_web[w] - 1, :] = self.mode_shape[
:, :, idxs[w + 1], 6:]
plot = False
if plot:
......@@ -199,5 +389,3 @@ class FINSTRIPWrapperBlade(FINSTRIPWrapper):
self.node_loc_web[
w, :self.nnode_web[w], 0] + sf * self.node_u_web[w, mode, :self.nnode_web[w], 0],
self.node_loc_web[w, :self.nnode_web[w], 1] + sf * self.node_u_web[w, mode, :self.nnode_web[w], 1], '--')
plt.show()
......@@ -34,23 +34,24 @@ class FINSTRIPWrapper(object):
def __init__(self, **kwargs):
super(FINSTRIPWrapper, self).__init__()
# finstrip program parameters
# maximum size of generated element in mm
self.maxlength_element = 0.252E+03
self.lambda_steps = 20 # number of half wavelength steps
self.lambda_min = 0.01E+03 # minimum half wavelength in mm
self.lambda_max = 15.126E+03 # maximum half wavelength in mm
# element length as fraction of arc length of airfoil surface
self.elsize_factor = 0.015
self.spanpos = 0 # position along span of analysed cross-section
# finstrip program parameters
self.lambda_min = 0.01 # minimum half wavelength in m
self.lambda_max = 15.126 # maximum half wavelength in m
self.nmode = self.lambda_steps + 1
self.ncoord = 79
self.ncase = 1
self.nmode = 21
#self.ncoord = 79
#self.ncase = 1
# flags
self.debug_mode = True
self.dry_run = False
def _init_half_wavelength_steps(self):
self.lambda_steps = self.nmode - 1 # number of half wavelength steps
def _init_finstrip_input_files_names(self):
''' Initializes standard input file names or extension
'''
......@@ -59,10 +60,11 @@ class FINSTRIPWrapper(object):
# input file names
self.buc_sfx = '.buc'
self.lcd_sfx = '.lcd'
self.lcd_sfx = '.ldc'
# self.fdb_prfx = 'db' # why no "." ??
self.txt_sfx = '.txt'
self.file_input = self.name + self.buc_sfx
self.file_load = self.name + self.lcd_sfx
#self.file_output = self.name + self.fcdb_prfx
#self.file_mat = 'MAT' + self.fd_sfx
#self.file_failmat = 'FAILMAT' + self.fd_sfx
......
This diff is collapsed.
......@@ -23,6 +23,8 @@ from becas_wrapper.becas_bladestructure import BECASBeamStructureKM,\
from becas_wrapper.becas_stressrecovery import BECASStressRecovery
from distutils.spawn import find_executable
from finstrip_wrapper.omdao.finstrip_blade import FINSTRIPBeamStructure
from robe.lib.load import convert_point_forces_to_moment_distro
_matlab_installed = find_executable('matlab')
......@@ -43,12 +45,12 @@ else:
DRYRUN = False
def configure_BECASBeamStructure(nsec, exec_mode, path_data, dry_run=False, FPM=False, with_sr=False, param2=False):
def configure_BECASBeamStructure(nsec, exec_mode, path_data, dry_run=False, FPM=False, with_sr=False, param2=False, with_becas=False):
p = Problem(impl=impl, root=Group())
p.root.add('blade_length_c', IndepVarComp(
'blade_length', 86.366), promotes=['*'])
'blade_length', 86.366, units='m'), promotes=['*'])
pf = read_blade_planform(
os.path.join(path_data, 'DTU_10MW_RWT_blade_axis_prebend.dat'))
......@@ -116,30 +118,79 @@ def configure_BECASBeamStructure(nsec, exec_mode, path_data, dry_run=False, FPM=
cfg['debug_mode'] = False
config['BECASWrapper'] = cfg
p.root.add('stiffness', BECASBeamStructure(p.root, config, st3dn,
(200, nsec_st, 3)), promotes=['*'])
if with_sr:
p.root.add('stress_recovery', BECASStressRecovery(
config, s_st, 2), promotes=['*'])
if with_becas:
p.root.add('stiffness', BECASBeamStructure(p.root, config, st3dn,
(200, nsec_st, 3)), promotes=['*'])
# if with_sr:
# p.root.add('stress_recovery', BECASStressRecovery(
# config, s_st, 2), promotes=['*'])
cfg = {}
cfg['nmode'] = 20
cfg['core_mats'] = ['balsa']
cfg['elsize_factor'] = 0.015
cfg['lambda_min'] = 0.01 # minimum half wavelength in m
cfg['lambda_max'] = 15.126 # maximum half wavelength in m
cfg['use_becas_stiffness'] = with_becas
config['FINSTRIPWrapper'] = cfg
ncase = 2
p.root.add('stability', FINSTRIPBeamStructure(p.root, config, st3dn,
(200, nsec_st, 3), ncase), promotes=['*'])
p.setup()
p['hub_radius'] = 2.8
for k, v in pf.iteritems():
if k in p.root.pf_splines.params.keys():
p.root.pf_splines.params[k] = v
if with_sr:
# set some arbitrary values in the load vectors used to compute strains
for i, x in enumerate(s_st):
try:
p['load_cases_sec%03d' % i] = np.ones((2, 6))
except:
pass
# assuming we have an external shear force loading along the blade, which
# need to be transformed to internal section load case
nload = 9
cfg = {}
cfg['load_cases'] = np.zeros((ncase, nload, 15))
hr = 2.8
cfg['load_cases'][0, :, :] = np.array(
[[0.26658496, -0.064059925, 20.073 - hr, 0.0, 0.0, 0.0,
0.0, 230000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.49336469, -0.04793707, 30.437 - hr, 0.0, 0.0,
0.0, 0.0, 270000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.47505332, -0.053173575, 33.028 - hr, 0.0,
0.0, 0.0, 290000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.37001989, -0.01194097, 47.71 - hr, 0.0, 0.0, 0.0,
180000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.33257294, 0.00242491, 52.892 - hr, 0.0, 0.0, 0.0,
0.0, 250000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.26347361, 0.01553866, 62.393 - hr, 0.0, 0.0, 0.0,
130000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.23826286, 0.01686549, 65.847 - hr, 0.0, 0.0, 0.0,
0.0, 220000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.16251535, 0.017204025, 76.211 - hr, 0.0, 0.0, 0.0,
18000, 190000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
[0.09939912, 0.014232575, 84.848 - hr, 0.0, 0.0, 0.0,
25000, 165000, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]])
csmoments = np.zeros((len(s_st), 3))
moment_positions = s_st * p['blade_length']
for cs, _ in enumerate(csmoments):
csmoments[cs, :] = convert_point_forces_to_moment_distro(
moment_position=moment_positions[cs],
force_positions=cfg['load_cases'][0, :, 2],
forces=cfg['load_cases'][0, :, 6:9])
for i, _ in enumerate(s_st):
p['load_cases_sec%03d' % i] = np.zeros((ncase, 6))
# Fx,Fy,Fz,Mx,My,Mz
p['load_cases_sec%03d' % i][0, 3] = csmoments[i, 0]
p['load_cases_sec%03d' % i][0, 4] = csmoments[i, 1]
# make second dummy load case
p['load_cases_sec%03d' % i][1, 3] = -2 * csmoments[i, 0]
p['load_cases_sec%03d' % i][1, 4] = -2 * csmoments[i, 1]
if with_becas:
p['hub_radius'] = hr
return p
class BECASBeamStructureTestCase(unittest.TestCase):
class FINSTRIPBeamStructureTestCase(unittest.TestCase):
def setUp(self):
pass
......@@ -147,6 +198,7 @@ class BECASBeamStructureTestCase(unittest.TestCase):
def tearDown(self):
pass
#@unittest.skip("skipping")
@unittest.skipIf(not _matlab_installed,
"Matlab not available on this system")
def test_standard_matlab(self):
......@@ -155,8 +207,50 @@ class BECASBeamStructureTestCase(unittest.TestCase):
dry_run=DRYRUN,
FPM=False,
with_sr=True,
param2=False)
param2=False,
with_becas=True)
p.run()
ncase = 2
nsec = len(p['s_st'])
f_eig_min = np.zeros((ncase, nsec))
for case in range(ncase):
for sec in range(nsec):
f_eig_min[case, sec] = np.min(
p['blade_f_eig'][case, sec, :, 1])
f_eig_min_exp = np.array([[0.727, 0.833, 1.385, 1.],
[4.734, 3.066, 1.694, 1.]])
for case in range(ncase):
for i, f in enumerate(f_eig_min[case, :]):
self.assertAlmostEqual(
f, f_eig_min_exp[case, i], places=3)
#@unittest.skip("skipping")
def test_finstrip(self):
p = configure_BECASBeamStructure(4, 'matlab', 'data',
dry_run=DRYRUN,
FPM=False,
with_sr=True,
param2=False,
with_becas=False)
p.run()
ncase = 2
nsec = len(p['s_st'])
f_eig_min = np.zeros((ncase, nsec))
for case in range(ncase):
for sec in range(nsec):
f_eig_min[case, sec] = np.min(
p['blade_f_eig'][case, sec, :, 1])
f_eig_min_exp = np.array([[0.996, 1.877, 2.09, 1.],
[6.484, 5.243, 2.011, 1.]])
for case in range(ncase):
for i, f in enumerate(f_eig_min[case, :]):
self.assertAlmostEqual(
f, f_eig_min_exp[case, i], places=3)
if __name__ == "__main__":
......
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