Skip to content
GitLab
Menu
Projects
Groups
Snippets
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in / Register
Toggle navigation
Menu
Open sidebar
Frederik Zahle
SEAMRotor
Commits
b6a73ba2
Commit
b6a73ba2
authored
Jun 26, 2015
by
Frederik Zahle
Browse files
a few minor adjustments
parent
b2465642
Changes
4
Hide whitespace changes
Inline
Side-by-side
rotor_unit.py
deleted
100644 → 0
View file @
b2465642
__author__
=
's127504'
from
decimal
import
*
import
matplotlib.pyplot
as
plt
from
math
import
pow
import
numpy
as
np
from
scipy.optimize
import
minimize
from
openmdao.main.api
import
Component
from
openmdao.lib.datatypes.api
import
Float
,
Array
,
Bool
class
SEAMRotor
(
Component
):
Nsections
=
Float
(
iotype
=
'in'
,
desc
=
'number of sections'
)
Neq
=
Float
(
1.e7
,
iotype
=
'in'
,
desc
=
''
)
WohlerExpFlap
=
Float
(
iotype
=
'in'
,
desc
=
''
)
PMtarget
=
Float
(
iotype
=
'in'
,
desc
=
''
)
Diameter
=
Float
(
iotype
=
'in'
,
units
=
'm'
,
desc
=
''
)
#[m]
# RootChord = Float(iotype='in', units='m', desc='') #[m]
# MaxChord = Float(iotype='in', units='m', desc='') #[m]
MaxChordrR
=
Float
(
iotype
=
'in'
,
units
=
'm'
,
desc
=
''
)
#[m]
OverallMaxFlap
=
Float
(
iotype
=
'in'
,
desc
=
''
)
OverallMaxEdge
=
Float
(
iotype
=
'in'
,
desc
=
''
)
TIF_FLext
=
Float
(
iotype
=
'in'
,
desc
=
''
)
# Tech Impr Factor _ flap extreme
TIF_EDext
=
Float
(
iotype
=
'in'
,
desc
=
''
)
FlapLEQ
=
Float
(
iotype
=
'in'
,
desc
=
''
)
EdgeLEQ
=
Float
(
iotype
=
'in'
,
desc
=
''
)
TIF_FLfat
=
Float
(
iotype
=
'in'
,
desc
=
''
)
sc_frac_flap
=
Float
(
iotype
=
'in'
,
desc
=
''
)
# sparcap fraction of chord flap
sc_frac_edge
=
Float
(
iotype
=
'in'
,
desc
=
''
)
# sparcap fraction of thickness edge
SF_blade
=
Float
(
iotype
=
'in'
,
desc
=
''
)
#[factor]
Slim_ext_blade
=
Float
(
iotype
=
'in'
,
units
=
'MPa'
,
desc
=
''
)
Slim_fat_blade
=
Float
(
iotype
=
'in'
,
units
=
'MPa'
,
desc
=
''
)
AddWeightFactorBlade
=
Float
(
iotype
=
'in'
,
desc
=
''
)
# Additional weight factor for blade shell
BladeDens
=
Float
(
iotype
=
'in'
,
units
=
'kg/m**3'
,
desc
=
'density of blades'
)
# [kg / m^3]
BladeCostPerMass
=
Float
(
iotype
=
'in'
,
desc
=
''
)
#[e/kg]
HubCostPerMass
=
Float
(
iotype
=
'in'
,
desc
=
''
)
#[e/kg]
SpinnerCostPerMass
=
Float
(
iotype
=
'in'
,
desc
=
''
)
#[e/kg]
BladeWeight
=
Float
(
iotype
=
'out'
,
units
=
'kg'
,
desc
=
'BladeMass'
)
def
execute
(
self
):
volumen
=
0.
# [m^3]
r
=
np
.
linspace
(
0
,
self
.
Diameter
/
2.
,
self
.
Nsections
)
C
=
np
.
zeros
(
r
.
shape
)
thick
=
np
.
zeros
(
r
.
shape
)
RootChord
=
(
self
.
Diameter
/
2.
)
/
25.
if
RootChord
<
1
:
RootChord
=
1
MaxChord
=
0.001
*
(
self
.
Diameter
/
2.
)
**
2
-
0.0354
*
(
self
.
Diameter
/
2.
)
+
1.6635
#Empirical
if
MaxChord
<
1.5
:
RootChord
=
1.5
for
i
in
range
(
0
,
int
(
self
.
Nsections
)):
# Calculating the chord and thickness
if
r
[
i
]
>
(
self
.
MaxChordrR
*
self
.
Diameter
/
2.
):
C
[
i
]
=
MaxChord
-
(
MaxChord
)
/
((
self
.
Diameter
/
2.
)
-
self
.
MaxChordrR
*
self
.
Diameter
/
2.
)
*
(
r
[
i
]
-
(
self
.
MaxChordrR
*
self
.
Diameter
/
2.
))
thick
[
i
]
=
0.4
*
MaxChord
-
(
0.4
*
MaxChord
-
0.05
*
MaxChord
)
/
(
self
.
Diameter
/
2.
-
self
.
MaxChordrR
*
self
.
Diameter
/
2.
)
*
(
r
[
i
]
-
self
.
MaxChordrR
*
self
.
Diameter
/
2.
)
if
thick
[
i
]
<
0.001
:
thick
[
i
]
=
0.001
if
C
[
i
]
<
0.001
:
C
[
i
]
=
0.001
else
:
C
[
i
]
=
MaxChord
-
(
MaxChord
)
/
((
self
.
Diameter
/
2.
)
-
self
.
MaxChordrR
*
(
self
.
Diameter
/
2.
))
*
(
self
.
MaxChordrR
*
(
self
.
Diameter
/
2.
)
-
r
[
i
])
thick
[
i
]
=
RootChord
if
thick
[
i
]
<
0.001
:
thick
[
i
]
=
0.001
if
C
[
i
]
<
0.001
:
C
[
i
]
=
0.001
# Printing radius, thickness and chord
print
'radius
\n
'
,
r
print
print
'thickness
\n
'
,
thick
print
print
'chord
\n
'
,
C
print
# Calculating load flapwise extreme
Mext_flap
=
(
self
.
OverallMaxFlap
-
1.75
*
self
.
OverallMaxFlap
*
r
/
(
self
.
Diameter
/
2.
)
+
0.75
*
self
.
OverallMaxFlap
*
r
/
(
self
.
Diameter
/
2.
)
*
r
/
(
self
.
Diameter
/
2.
))
*
self
.
TIF_FLext
# Calculating load edgewise extreme
Mext_edge
=
(
self
.
OverallMaxEdge
-
1.75
*
self
.
OverallMaxEdge
*
r
/
(
self
.
Diameter
/
2.
)
+
0.75
*
self
.
OverallMaxEdge
*
r
/
(
self
.
Diameter
/
2.
)
*
r
/
(
self
.
Diameter
/
2.
))
*
self
.
TIF_EDext
# Calculating load flapwise fatigue
Mfat_flap
=
(
self
.
FlapLEQ
-
1.75
*
self
.
FlapLEQ
*
r
/
(
self
.
Diameter
/
2.
)
+
0.75
*
self
.
FlapLEQ
*
r
/
(
self
.
Diameter
/
2.
)
*
r
/
(
self
.
Diameter
/
2.
))
*
self
.
TIF_FLfat
# Calculating load edgewise fatigue
Mfat_edge
=
(
self
.
EdgeLEQ
-
1.75
*
self
.
EdgeLEQ
*
r
/
(
self
.
Diameter
/
2.
)
+
0.75
*
self
.
EdgeLEQ
*
r
/
(
self
.
Diameter
/
2.
)
*
r
/
(
self
.
Diameter
/
2.
))
*
self
.
TIF_EDext
# Printing the 4 load cases
print
'Mext_flap
\n
'
,
Mext_flap
print
print
'Mext_edge
\n
'
,
Mext_edge
print
print
'Mfat_flap
\n
'
,
Mfat_flap
print
print
'Mfat_edge
\n
'
,
Mfat_edge
print
# Calculating thickness flapwise extreme
text_flap
=
np
.
zeros
(
self
.
Nsections
)
for
i
in
range
(
int
(
self
.
Nsections
)):
res
=
minimize
(
self
.
solve_text_flap
,
0.03
,
args
=
(
C
[
i
],
thick
[
i
],
Mext_flap
[
i
]),
bounds
=
None
,
method
=
'COBYLA'
,
tol
=
1.e-12
,
options
=
{
'maxiter'
:
10000
})
text_flap
[
i
]
=
res
[
'x'
]
#print 'done', text_flap[i]
# if not res['success']: print 'WARNING', res
# Calculating thickness edgewise extreme
text_edge
=
np
.
zeros
(
self
.
Nsections
)
for
i
in
range
(
int
(
self
.
Nsections
)):
res
=
minimize
(
self
.
solve_text_edge
,
0.03
,
args
=
(
C
[
i
],
thick
[
i
],
Mext_edge
[
i
]),
bounds
=
None
,
method
=
'COBYLA'
,
tol
=
1.e-8
)
text_edge
[
i
]
=
res
[
'x'
]
# if not res['success']: print 'WARNING', res
# Calculating thickness flapwise fatigue
tfat_flap
=
np
.
zeros
(
self
.
Nsections
)
for
i
in
range
(
int
(
self
.
Nsections
)):
res
=
minimize
(
self
.
solve_tfat_flap
,
0.03
,
args
=
(
C
[
i
],
thick
[
i
],
Mfat_flap
[
i
]),
bounds
=
None
,
method
=
'COBYLA'
,
tol
=
1.e-8
,
options
=
{
'maxiter'
:
10000
})
tfat_flap
[
i
]
=
res
[
'x'
]
#print 'done', tfat_flap[i]
# if not res['success']: print 'WARNING', res
# Calculating thickness edgewise fatigue
tfat_edge
=
np
.
zeros
(
self
.
Nsections
)
for
i
in
range
(
int
(
self
.
Nsections
)):
res
=
minimize
(
self
.
solve_tfat_edge
,
0.03
,
args
=
(
C
[
i
],
thick
[
i
],
Mfat_edge
[
i
]),
bounds
=
None
,
method
=
'COBYLA'
,
tol
=
1.e-8
)
tfat_edge
[
i
]
=
res
[
'x'
]
# if not res['success']: print 'WARNING', res
# Printing the 4 thicknesses
print
'text_flap
\n
'
,
text_flap
print
print
'text_edge
\n
'
,
text_edge
print
print
'tfat_flap
\n
'
,
tfat_flap
print
print
'tfat_edge
\n
'
,
tfat_edge
print
self
.
x
=
np
.
array
(
range
(
1
,
22
))
self
.
text_flap
=
text_flap
self
.
text_edge
=
text_edge
self
.
tfat_flap
=
tfat_flap
self
.
tfat_edge
=
tfat_edge
# Comparing found values for thickness with Kenneth's
datacompare_text_flap
=
(
abs
(
data
[:,
5
]
-
text_flap
))
datacompare_text_edge
=
(
abs
(
data
[:,
7
]
-
text_edge
))
datacompare_tfat_flap
=
(
abs
(
data
[:,
9
]
-
tfat_flap
))
datacompare_tfat_edge
=
(
abs
(
data
[:,
11
]
-
tfat_edge
))
# Calculating the final thickness for flap and edge direction
tfinal_flap
=
np
.
maximum
(
text_flap
,
tfat_flap
)
tfinal_edge
=
np
.
maximum
(
text_edge
,
tfat_edge
)
print
'tfinal_flap'
,
tfinal_flap
print
print
'tfinal_edge'
,
tfinal_edge
print
# Calculating different costs and mass
for
i
in
range
(
0
,
int
(
self
.
Nsections
)):
if
i
>
0
:
volumen
=
volumen
+
self
.
AddWeightFactorBlade
*
(
r
[
i
]
-
r
[
i
-
1
])
*
((
2
*
self
.
sc_frac_flap
*
C
[
i
]
*
tfinal_flap
[
i
]
+
2
*
self
.
sc_frac_edge
*
thick
[
i
]
*
tfinal_edge
[
i
])
+
(
2
*
self
.
sc_frac_flap
*
C
[
i
-
1
]
*
tfinal_flap
[
i
-
1
]
+
2
*
self
.
sc_frac_edge
*
thick
[
i
-
1
]
*
tfinal_edge
[
i
-
1
]))
/
2.
self
.
BladeWeight
=
self
.
BladeDens
*
volumen
print
'volumen'
,
volumen
print
'BladeWeigth'
,
self
.
BladeWeight
rotor
=
self
.
BladeWeight
*
self
.
BladeCostPerMass
/
1e6
# Meuro
print
'rotor'
,
rotor
# Scaling laws for hub, spinner, pitch
HubMass
=
0.954
*
self
.
BladeWeight
+
5680.3
print
'HubMass'
,
HubMass
HubCost
=
HubMass
*
self
.
HubCostPerMass
/
1e6
print
'HubCost'
,
HubCost
#PitchCost
SpinnerMass
=
18.5
*
self
.
Diameter
-
520.5
print
'SpinnerMass'
,
SpinnerMass
SpinnerCost
=
SpinnerMass
*
self
.
SpinnerCostPerMass
/
1e6
print
'SpinnerCost'
,
SpinnerCost
#Rotor = 3*(self.BladeWeight*self.BladeCostPerMass/1e6)+HubCost+PitchCost+SpinnerCost
# Plot illustrating the absolute differences for the thicknesses as function of the radius
plt
.
figure
(
1
)
plt
.
plot
(
r
,
datacompare_text_flap
,
'o-'
,
label
=
'text_flap'
)
plt
.
plot
(
r
,
datacompare_text_edge
,
'o-'
,
label
=
'text_edge'
)
plt
.
plot
(
r
,
datacompare_tfat_flap
,
'o-'
,
label
=
'tfat_flap'
)
plt
.
plot
(
r
,
datacompare_tfat_edge
,
'o-'
,
label
=
'tfat_edge'
)
plt
.
legend
()
plt
.
title
(
"abs differences between thickness's as function of radius"
)
plt
.
axis
([
0
,
90
,
0
,
0.0003
])
plt
.
xlabel
(
'radius [m]'
)
plt
.
ylabel
(
"abs difference between thickness's"
)
plt
.
grid
()
plt
.
show
()
# Plot illustrating the 4 thicknesess's as function of the radius
plt
.
figure
(
2
)
plt
.
plot
(
r
,
text_flap
,
'o-'
,
label
=
'text_flap'
)
plt
.
plot
(
r
,
text_edge
,
'o-'
,
label
=
'text_edge'
)
plt
.
plot
(
r
,
tfat_flap
,
'o-'
,
label
=
'tfat_flap'
)
plt
.
plot
(
r
,
tfat_edge
,
'o-'
,
label
=
'tfat_edge'
)
plt
.
legend
()
plt
.
title
(
'Thickness
\'
s as function of time'
)
plt
.
xlabel
(
'radius [m]'
)
plt
.
ylabel
(
'spar thickness [m]'
)
plt
.
grid
()
plt
.
show
()
# Defining functions in order to solve for the thickness
# Solving for t in flap direction, extremes
def
solve_text_flap
(
self
,
t
,
C
,
thick
,
Mext_flap
):
Ine
=
(
2.
/
3.
)
*
self
.
sc_frac_flap
*
C
*
t
**
3
-
self
.
sc_frac_flap
*
C
*
thick
*
t
**
2
+
self
.
sc_frac_flap
*
C
*
thick
**
2
/
2.
*
t
W
=
Ine
/
(
thick
/
2.
)
sext
=
self
.
SF_blade
*
1e3
*
Mext_flap
/
W
/
1e6
return
abs
(
sext
-
self
.
Slim_ext_blade
)
#Solving for t in edge direction, extremes
def
solve_text_edge
(
self
,
t
,
C
,
thick
,
Mext_edge
):
Ine
=
(
2
/
3.
)
*
self
.
sc_frac_edge
*
thick
*
t
**
3
-
self
.
sc_frac_edge
*
thick
*
C
*
t
**
2
+
self
.
sc_frac_edge
*
thick
*
C
**
2
/
2.
*
t
W
=
Ine
/
(
C
/
2.
);
sext
=
self
.
SF_blade
*
1.e3
*
Mext_edge
/
W
/
1.e6
return
abs
(
sext
-
self
.
Slim_ext_blade
)
# Solving for t in flap direction, fatigue
def
solve_tfat_flap
(
self
,
t
,
C
,
thick
,
Mfat_flap
):
Ine
=
(
2
/
3.
)
*
self
.
sc_frac_flap
*
C
*
t
**
3
-
self
.
sc_frac_flap
*
C
*
thick
*
t
**
2
+
self
.
sc_frac_flap
*
C
*
thick
**
2
/
2.
*
t
W
=
Ine
/
(
thick
/
2.
)
sfat
=
np
.
maximum
(
1.e-6
,
self
.
SF_blade
*
1.e3
*
Mfat_flap
/
W
/
1.e6
)
PM
=
self
.
Neq
/
(
pow
(
10
,
(
self
.
Slim_fat_blade
-
self
.
WohlerExpFlap
*
np
.
log10
(
sfat
))))
return
abs
(
PM
-
self
.
PMtarget
)
# Solving for t in edge direction, fatigue
def
solve_tfat_edge
(
self
,
t
,
C
,
thick
,
Mfat_edge
):
Ine
=
(
2
/
3.
)
*
self
.
sc_frac_edge
*
thick
*
t
**
3
-
self
.
sc_frac_edge
*
thick
*
C
*
t
**
2
+
self
.
sc_frac_edge
*
thick
*
C
**
2
/
2.
*
t
W
=
Ine
/
(
C
/
2.
)
sfat
=
np
.
maximum
(
1.e-6
,
self
.
SF_blade
*
1.e3
*
Mfat_edge
/
W
/
1.e6
)
PM
=
self
.
Neq
/
(
pow
(
10
,
(
self
.
Slim_fat_blade
-
self
.
WohlerExpFlap
*
np
.
log10
(
sfat
))))
return
abs
(
PM
-
self
.
PMtarget
)
#Just to check if the script works
if
__name__
==
'__main__'
:
top
=
SEAMRotor
()
top
.
Nsections
=
21.
top
.
Neq
=
1e7
top
.
WohlerExpFlap
=
10.0
top
.
PMtarget
=
1.0
top
.
Diameter
=
177.
#[m]
# top.RootChord = 2.0 #[m]
# top.MaxChord = 2.5 #[m]
top
.
MaxChordrR
=
0.2
#[m]
top
.
OverallMaxFlap
=
47225.
top
.
OverallMaxEdge
=
26712.
top
.
TIF_FLext
=
1.
# Tech Impr Factor _ flap extreme
top
.
TIF_EDext
=
1.
top
.
FlapLEQ
=
26975.
top
.
EdgeLEQ
=
24252.
top
.
TIF_FLfat
=
1.
top
.
sc_frac_flap
=
0.3
# sparcap fraction of chord flap
top
.
sc_frac_edge
=
0.8
# sparcap fraction of thickness edge
top
.
SF_blade
=
1.1
#[factor]
top
.
Slim_ext_blade
=
200.0
top
.
Slim_fat_blade
=
27
top
.
AddWeightFactorBlade
=
1.2
# Additional weight factor for blade shell
top
.
BladeDens
=
2100.
# [kg / m^3]
top
.
BladeCostPerMass
=
15.0
#[e/kg]
top
.
HubCostPerMass
=
3.5
#[e/kg]
top
.
SpinnerCostPerMass
=
4.5
#[e/kg]
data
=
np
.
loadtxt
(
'Rotor_outputfile.txt'
)
top
.
execute
()
src/SEAMRotor/SEAMRotor.py
View file @
b6a73ba2
__a
ll__
=
[
'SEAMRotor'
]
__a
uthor__
=
's127504'
from
openmdao.main.api
import
Component
from
openmdao.lib.datatypes.api
import
Float
from
decimal
import
*
import
matplotlib.pyplot
as
plt
from
math
import
pow
import
numpy
as
np
from
scipy.optimize
import
minimize
from
openmdao.main.api
import
Component
from
openmdao.lib.datatypes.api
import
Float
,
Array
,
Bool
# Make sure that your class has some kind of docstring. Otherwise
# the descriptions for your variables won't show up in the
# source ducumentation.
class
SEAMRotor
(
Component
):
"""
"""
# declare inputs and outputs here, for example:
#x = Float(0.0, iotype='in', desc='description for x')
#y = Float(0.0, iotype='out', desc='description for y')
Nsections
=
Float
(
iotype
=
'in'
,
desc
=
'number of sections'
)
Neq
=
Float
(
1.e7
,
iotype
=
'in'
,
desc
=
''
)
WohlerExpFlap
=
Float
(
iotype
=
'in'
,
desc
=
''
)
PMtarget
=
Float
(
iotype
=
'in'
,
desc
=
''
)
Diameter
=
Float
(
iotype
=
'in'
,
units
=
'm'
,
desc
=
''
)
#[m]
# RootChord = Float(iotype='in', units='m', desc='') #[m]
# MaxChord = Float(iotype='in', units='m', desc='') #[m]
MaxChordrR
=
Float
(
iotype
=
'in'
,
units
=
'm'
,
desc
=
''
)
#[m]
OverallMaxFlap
=
Float
(
iotype
=
'in'
,
desc
=
''
)
OverallMaxEdge
=
Float
(
iotype
=
'in'
,
desc
=
''
)
TIF_FLext
=
Float
(
iotype
=
'in'
,
desc
=
''
)
# Tech Impr Factor _ flap extreme
TIF_EDext
=
Float
(
iotype
=
'in'
,
desc
=
''
)
FlapLEQ
=
Float
(
iotype
=
'in'
,
desc
=
''
)
EdgeLEQ
=
Float
(
iotype
=
'in'
,
desc
=
''
)
TIF_FLfat
=
Float
(
iotype
=
'in'
,
desc
=
''
)
sc_frac_flap
=
Float
(
iotype
=
'in'
,
desc
=
''
)
# sparcap fraction of chord flap
sc_frac_edge
=
Float
(
iotype
=
'in'
,
desc
=
''
)
# sparcap fraction of thickness edge
SF_blade
=
Float
(
iotype
=
'in'
,
desc
=
''
)
#[factor]
Slim_ext_blade
=
Float
(
iotype
=
'in'
,
units
=
'MPa'
,
desc
=
''
)
Slim_fat_blade
=
Float
(
iotype
=
'in'
,
units
=
'MPa'
,
desc
=
''
)
AddWeightFactorBlade
=
Float
(
iotype
=
'in'
,
desc
=
''
)
# Additional weight factor for blade shell
BladeDens
=
Float
(
iotype
=
'in'
,
units
=
'kg/m**3'
,
desc
=
'density of blades'
)
# [kg / m^3]
BladeCostPerMass
=
Float
(
iotype
=
'in'
,
desc
=
''
)
#[e/kg]
HubCostPerMass
=
Float
(
iotype
=
'in'
,
desc
=
''
)
#[e/kg]
SpinnerCostPerMass
=
Float
(
iotype
=
'in'
,
desc
=
''
)
#[e/kg]
BladeWeight
=
Float
(
iotype
=
'out'
,
units
=
'kg'
,
desc
=
'BladeMass'
)
def
execute
(
self
):
""" do your calculations here """
pass
\ No newline at end of file
volumen
=
0.
# [m^3]
r
=
np
.
linspace
(
0
,
self
.
Diameter
/
2.
,
self
.
Nsections
)
C
=
np
.
zeros
(
r
.
shape
)
thick
=
np
.
zeros
(
r
.
shape
)
RootChord
=
(
self
.
Diameter
/
2.
)
/
25.
if
RootChord
<
1
:
RootChord
=
1
MaxChord
=
0.001
*
(
self
.
Diameter
/
2.
)
**
2
-
0.0354
*
(
self
.
Diameter
/
2.
)
+
1.6635
#Empirical
if
MaxChord
<
1.5
:
RootChord
=
1.5
for
i
in
range
(
0
,
int
(
self
.
Nsections
)):
# Calculating the chord and thickness
if
r
[
i
]
>
(
self
.
MaxChordrR
*
self
.
Diameter
/
2.
):
C
[
i
]
=
MaxChord
-
(
MaxChord
)
/
((
self
.
Diameter
/
2.
)
-
self
.
MaxChordrR
*
self
.
Diameter
/
2.
)
*
(
r
[
i
]
-
(
self
.
MaxChordrR
*
self
.
Diameter
/
2.
))
thick
[
i
]
=
0.4
*
MaxChord
-
(
0.4
*
MaxChord
-
0.05
*
MaxChord
)
/
(
self
.
Diameter
/
2.
-
self
.
MaxChordrR
*
self
.
Diameter
/
2.
)
*
(
r
[
i
]
-
self
.
MaxChordrR
*
self
.
Diameter
/
2.
)
if
thick
[
i
]
<
0.001
:
thick
[
i
]
=
0.001
if
C
[
i
]
<
0.001
:
C
[
i
]
=
0.001
else
:
C
[
i
]
=
MaxChord
-
(
MaxChord
)
/
((
self
.
Diameter
/
2.
)
-
self
.
MaxChordrR
*
(
self
.
Diameter
/
2.
))
*
(
self
.
MaxChordrR
*
(
self
.
Diameter
/
2.
)
-
r
[
i
])
thick
[
i
]
=
RootChord
if
thick
[
i
]
<
0.001
:
thick
[
i
]
=
0.001
if
C
[
i
]
<
0.001
:
C
[
i
]
=
0.001
# Printing radius, thickness and chord
print
'radius
\n
'
,
r
print
print
'thickness
\n
'
,
thick
print
print
'chord
\n
'
,
C
print
self
.
thick
=
thick
self
.
chord
=
C
# Calculating load flapwise extreme
Mext_flap
=
(
self
.
OverallMaxFlap
-
1.75
*
self
.
OverallMaxFlap
*
r
/
(
self
.
Diameter
/
2.
)
+
0.75
*
self
.
OverallMaxFlap
*
r
/
(
self
.
Diameter
/
2.
)
*
r
/
(
self
.
Diameter
/
2.
))
*
self
.
TIF_FLext
# Calculating load edgewise extreme
Mext_edge
=
(
self
.
OverallMaxEdge
-
1.75
*
self
.
OverallMaxEdge
*
r
/
(
self
.
Diameter
/
2.
)
+
0.75
*
self
.
OverallMaxEdge
*
r
/
(
self
.
Diameter
/
2.
)
*
r
/
(
self
.
Diameter
/
2.
))
*
self
.
TIF_EDext
# Calculating load flapwise fatigue
Mfat_flap
=
(
self
.
FlapLEQ
-
1.75
*
self
.
FlapLEQ
*
r
/
(
self
.
Diameter
/
2.
)
+
0.75
*
self
.
FlapLEQ
*
r
/
(
self
.
Diameter
/
2.
)
*
r
/
(
self
.
Diameter
/
2.
))
*
self
.
TIF_FLfat
# Calculating load edgewise fatigue
Mfat_edge
=
(
self
.
EdgeLEQ
-
1.75
*
self
.
EdgeLEQ
*
r
/
(
self
.
Diameter
/
2.
)
+
0.75
*
self
.
EdgeLEQ
*
r
/
(
self
.
Diameter
/
2.
)
*
r
/
(
self
.
Diameter
/
2.
))
*
self
.
TIF_EDext
# Printing the 4 load cases
print
'Mext_flap
\n
'
,
Mext_flap
print
print
'Mext_edge
\n
'
,
Mext_edge
print
print
'Mfat_flap
\n
'
,
Mfat_flap
print
print
'Mfat_edge
\n
'
,
Mfat_edge
print
# Calculating thickness flapwise extreme
text_flap
=
np
.
zeros
(
self
.
Nsections
)
for
i
in
range
(
int
(
self
.
Nsections
)):
res
=
minimize
(
self
.
solve_text_flap
,
0.03
,
args
=
(
C
[
i
],
thick
[
i
],
Mext_flap
[
i
]),
bounds
=
None
,
method
=
'COBYLA'
,
tol
=
1.e-12
,
options
=
{
'maxiter'
:
10000
})
text_flap
[
i
]
=
res
[
'x'
]
#print 'done', text_flap[i]
# if not res['success']: print 'WARNING', res
# Calculating thickness edgewise extreme
text_edge
=
np
.
zeros
(
self
.
Nsections
)
for
i
in
range
(
int
(
self
.
Nsections
)):
res
=
minimize
(
self
.
solve_text_edge
,
0.03
,
args
=
(
C
[
i
],
thick
[
i
],
Mext_edge
[
i
]),
bounds
=
None
,
method
=
'COBYLA'
,
tol
=
1.e-8
)
text_edge
[
i
]
=
res
[
'x'
]
# if not res['success']: print 'WARNING', res
# Calculating thickness flapwise fatigue
tfat_flap
=
np
.
zeros
(
self
.
Nsections
)
for
i
in
range
(
int
(
self
.
Nsections
)):
res
=
minimize
(
self
.
solve_tfat_flap
,
0.03
,
args
=
(
C
[
i
],
thick
[
i
],
Mfat_flap
[
i
]),
bounds
=
None
,
method
=
'COBYLA'
,
tol
=
1.e-8
,
options
=
{
'maxiter'
:
10000
})
tfat_flap
[
i
]
=
res
[
'x'
]
#print 'done', tfat_flap[i]
# if not res['success']: print 'WARNING', res
# Calculating thickness edgewise fatigue
tfat_edge
=
np
.
zeros
(
self
.
Nsections
)
for
i
in
range
(
int
(
self
.
Nsections
)):
res
=
minimize
(
self
.
solve_tfat_edge
,
0.03
,
args
=
(
C
[
i
],
thick
[
i
],
Mfat_edge
[
i
]),
bounds
=
None
,
method
=
'COBYLA'
,
tol
=
1.e-8
)
tfat_edge
[
i
]
=
res
[
'x'
]
# if not res['success']: print 'WARNING', res
# Printing the 4 thicknesses
print
'text_flap
\n
'
,
text_flap
print
print
'text_edge
\n
'
,
text_edge
print
print
'tfat_flap
\n
'
,
tfat_flap
print
print
'tfat_edge
\n
'
,
tfat_edge
print
self
.
x
=
np
.
array
(
range
(
1
,
22
))
self
.
text_flap
=
text_flap
self
.
text_edge
=
text_edge
self
.
tfat_flap
=
tfat_flap
self
.
tfat_edge
=
tfat_edge
self
.
r
=
r
# Calculating the final thickness for flap and edge direction
tfinal_flap
=
np
.
maximum
(
text_flap
,
tfat_flap
)
tfinal_edge
=
np
.
maximum
(
text_edge
,
tfat_edge
)
print
'tfinal_flap'
,
tfinal_flap
print
print
'tfinal_edge'
,
tfinal_edge
print
# Calculating different costs and mass
for
i
in
range
(
0
,
int
(
self
.
Nsections
)):
if
i
>
0
:
volumen
=
volumen
+
self
.
AddWeightFactorBlade
*
(
r
[
i
]
-
r
[
i
-
1
])
*
((
2
*
self
.
sc_frac_flap
*
C
[
i
]
*
tfinal_flap
[
i
]
+
2
*
self
.
sc_frac_edge
*
thick
[
i
]
*
tfinal_edge
[
i
])
+
(
2
*
self
.
sc_frac_flap
*
C
[
i
-
1
]
*
tfinal_flap
[
i
-
1
]
+
2
*
self
.
sc_frac_edge
*
thick
[
i
-
1
]
*
tfinal_edge
[
i
-
1
]))
/
2.
self
.
BladeWeight
=
self
.
BladeDens
*
volumen
print
'volumen'
,
volumen
print
'BladeWeigth'
,
self
.
BladeWeight
rotor
=
self
.
BladeWeight
*
self
.
BladeCostPerMass
/
1e6
# Meuro
print
'rotor'
,
rotor
# Scaling laws for hub, spinner, pitch
HubMass
=
0.954
*
self
.
BladeWeight
+
5680.3
print
'HubMass'
,
HubMass
HubCost
=
HubMass
*
self
.
HubCostPerMass
/
1e6
print
'HubCost'
,
HubCost
#PitchCost
SpinnerMass
=
18.5
*
self
.
Diameter
-
520.5
print
'SpinnerMass'
,
SpinnerMass
SpinnerCost
=
SpinnerMass
*
self
.
SpinnerCostPerMass
/
1e6
print
'SpinnerCost'
,
SpinnerCost
#Rotor = 3*(self.BladeWeight*self.BladeCostPerMass/1e6)+HubCost+PitchCost+SpinnerCost
# Defining functions in order to solve for the thickness
# Solving for t in flap direction, extremes
def
solve_text_flap
(
self
,
t
,
C
,
thick
,
Mext_flap
):
Ine
=
(
2.
/
3.
)
*
self
.
sc_frac_flap
*
C
*
t
**
3
-
self
.
sc_frac_flap
*
C
*
thick
*
t
**
2
+
self
.
sc_frac_flap
*
C
*
thick
**
2
/
2.
*
t
W
=
Ine
/
(
thick
/
2.
)
sext
=
self
.
SF_blade
*
1e3
*
Mext_flap
/
W
/
1e6
return
abs
(
sext
-
self
.
Slim_ext_blade
)
#Solving for t in edge direction, extremes
def
solve_text_edge
(
self
,
t
,
C
,
thick
,
Mext_edge
):
Ine
=
(
2
/
3.
)
*
self
.
sc_frac_edge
*
thick
*
t
**
3
-
self
.
sc_frac_edge
*
thick
*
C
*
t
**
2
+
self
.
sc_frac_edge
*
thick
*
C
**
2
/
2.
*
t
W
=
Ine
/
(
C
/
2.
);
sext
=
self
.
SF_blade
*
1.e3
*
Mext_edge
/
W
/
1.e6
return
abs
(
sext
-
self
.
Slim_ext_blade
)
# Solving for t in flap direction, fatigue
def
solve_tfat_flap
(
self
,
t
,
C
,
thick
,
Mfat_flap
):
Ine
=
(
2
/
3.
)
*
self
.
sc_frac_flap
*
C
*
t
**
3
-
self
.
sc_frac_flap
*
C
*
thick
*
t
**
2
+
self
.
sc_frac_flap
*
C
*
thick
**
2
/
2.
*
t
W
=
Ine
/
(
thick
/
2.
)
sfat
=
np
.
maximum
(
1.e-6
,
self
.
SF_blade
*
1.e3
*
Mfat_flap
/
W
/
1.e6
)
PM
=
self
.
Neq
/
(
pow
(
10
,
(
self
.
Slim_fat_blade
-
self
.
WohlerExpFlap
*
np
.
log10
(
sfat
))))
return
abs
(
PM
-
self
.
PMtarget
)
# Solving for t in edge direction, fatigue
def
solve_tfat_edge
(
self
,
t
,
C
,
thick
,
Mfat_edge
):
Ine
=
(
2
/
3.
)
*
self
.
sc_frac_edge
*
thick
*
t
**
3
-
self
.
sc_frac_edge
*
thick
*
C
*
t
**
2
+
self
.
sc_frac_edge
*
thick
*
C
**
2
/
2.
*
t
W
=
Ine
/
(
C
/
2.
)
sfat
=
np
.
maximum
(
1.e-6
,
self
.
SF_blade
*
1.e3
*
Mfat_edge
/
W
/
1.e6
)
PM
=
self
.
Neq
/
(
pow
(
10
,
(
self
.
Slim_fat_blade
-
self
.
WohlerExpFlap
*
np
.
log10
(
sfat
))))
return
abs
(
PM
-
self
.
PMtarget
)
Rotor_outputfile.txt
→
src/SEAMRotor/test/
Rotor_outputfile.txt
View file @
b6a73ba2
File moved
src/SEAMRotor/test/test_SEAMRotor.py
View file @
b6a73ba2
import
numpy
as
np
import
unittest
from
SEAMRotor.SEAMRotor
import
SEAMRotor
class
SEAMRotorTestCase
(
unittest
.
TestCase
):
...
...
@@ -16,5 +18,82 @@ class SEAMRotorTestCase(unittest.TestCase):
#pass
if
__name__
==
"__main__"
:
unittest
.
main
()
\ No newline at end of file
# unittest.main()
top
=
SEAMRotor
()
top
.
Nsections
=
21.
top
.
Neq
=
1e7
top
.
WohlerExpFlap
=
10.0
top
.
PMtarget
=
1.0