This repository has been archived by the owner on May 9, 2019. It is now read-only.
-
Notifications
You must be signed in to change notification settings - Fork 16
/
Copy pathrotor_components.py
120 lines (90 loc) · 6.5 KB
/
rotor_components.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
from openmdao.main.api import Component, VariableTree
from openmdao.lib.datatypes.api import Array, Float, Bool, Int, List, Str, VarTree
import numpy as np
from scipy import interp
## Components that interpolate a predefined CP/CT curve and apply a yaw correction
class windSpeedToCPCT(VariableTree):
# variable tree that defines known wind speed to CP/CT curve
def __init__(self, datasize=0):
super(windSpeedToCPCT, self).__init__()
self.add('wind_speed', Array(np.zeros(datasize), iotype='in', units='m/s', desc='range of wind speeds'))
self.add('CP', Array(np.zeros(datasize), iotype='out', desc='power coefficients'))
self.add('CT', Array(np.zeros(datasize), iotype='out', desc='thrust coefficients'))
class CPCT_Interpolate(Component):
pP = Float(3.0, iotype='in')
def __init__(self, nTurbines, datasize=0):
#print 'in CPCT_Interpolate init'
super(CPCT_Interpolate, self).__init__()
self.nTurbines = nTurbines
self.datasize = datasize
self.add('windSpeedToCPCT', VarTree(windSpeedToCPCT(datasize), iotype='in', desc='pre-calculated CPCT'))
self.add('yaw', Array(np.zeros(nTurbines), iotype='in', desc='yaw error', units='deg'))
self.add('wind_speed_hub', Array(np.zeros(nTurbines), iotype='in', units='m/s', desc='hub height wind speed')) # Uhub
self.add('CP', Array(np.zeros(nTurbines), iotype='out'))
self.add('CT', Array(np.zeros(nTurbines), iotype='out'))
def execute(self):
#print 'in CPCT_Interpolate'
print "first speed, first CT", self.windSpeedToCPCT.wind_speed[0], self.windSpeedToCPCT.CT[0], self.windSpeedToCPCT.CT[1]
print "last speed, last CT", self.windSpeedToCPCT.wind_speed[-1], self.windSpeedToCPCT.CP[-1]
print "pP: ", self.pP
wind_speed_ax = np.cos(self.yaw*np.pi/180.0)**(self.pP/3.0)*self.wind_speed_hub
# use interpolation on precalculated CP-CT curve
wind_speed_ax = np.maximum(wind_speed_ax, self.windSpeedToCPCT.wind_speed[0])
wind_speed_ax = np.minimum(wind_speed_ax, self.windSpeedToCPCT.wind_speed[-1])
self.CP = interp(wind_speed_ax, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CP)
self.CT = interp(wind_speed_ax, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CT)
# normalize on incoming wind speed to correct coefficients for yaw
self.CP = self.CP * np.cos(self.yaw*np.pi/180.0)**self.pP
self.CT = self.CT * np.cos(self.yaw*np.pi/180.0)**2
print "in rotor, Cp = ", self.CP
print "in rotor, Ct = ", self.CT
print "in rotor, pP = ", self.pP
# print 'in CPCT interp, wind_speed_hub = ', self.wind_speed_hub
def list_deriv_vars(self):
return ('yaw', 'wind_speed_hub'), ('CP', 'CT')
def provideJ(self):
#print 'in CPCT_Interpolate - provideJ'
# standard central differencing
# set step size for finite differencing
h = 1e-6
# calculate upper and lower function values
wind_speed_ax_high_yaw = np.cos((self.yaw+h)*np.pi/180.0)**(self.pP/3.0)*self.wind_speed_hub
wind_speed_ax_low_yaw = np.cos((self.yaw-h)*np.pi/180.0)**(self.pP/3.0)*self.wind_speed_hub
wind_speed_ax_high_wind = np.cos(self.yaw*np.pi/180.0)**(self.pP/3.0)*(self.wind_speed_hub+h)
wind_speed_ax_low_wind = np.cos(self.yaw*np.pi/180.0)**(self.pP/3.0)*(self.wind_speed_hub-h)
# use interpolation on precalculated CP-CT curve
wind_speed_ax_high_yaw = np.maximum(wind_speed_ax_high_yaw, self.windSpeedToCPCT.wind_speed[0])
wind_speed_ax_low_yaw = np.maximum(wind_speed_ax_low_yaw, self.windSpeedToCPCT.wind_speed[0])
wind_speed_ax_high_wind = np.maximum(wind_speed_ax_high_wind, self.windSpeedToCPCT.wind_speed[0])
wind_speed_ax_low_wind = np.maximum(wind_speed_ax_low_wind, self.windSpeedToCPCT.wind_speed[0])
wind_speed_ax_high_yaw = np.minimum(wind_speed_ax_high_yaw, self.windSpeedToCPCT.wind_speed[-1])
wind_speed_ax_low_yaw = np.minimum(wind_speed_ax_low_yaw, self.windSpeedToCPCT.wind_speed[-1])
wind_speed_ax_high_wind = np.minimum(wind_speed_ax_high_wind, self.windSpeedToCPCT.wind_speed[-1])
wind_speed_ax_low_wind = np.minimum(wind_speed_ax_low_wind, self.windSpeedToCPCT.wind_speed[-1])
CP_high_yaw = interp(wind_speed_ax_high_yaw, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CP)
CP_low_yaw = interp(wind_speed_ax_low_yaw, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CP)
CP_high_wind = interp(wind_speed_ax_high_wind, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CP)
CP_low_wind = interp(wind_speed_ax_low_wind, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CP)
CT_high_yaw = interp(wind_speed_ax_high_yaw, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CT)
CT_low_yaw = interp(wind_speed_ax_low_yaw, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CT)
CT_high_wind = interp(wind_speed_ax_high_wind, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CT)
CT_low_wind = interp(wind_speed_ax_low_wind, self.windSpeedToCPCT.wind_speed, self.windSpeedToCPCT.CT)
# normalize on incoming wind speed to correct coefficients for yaw
CP_high_yaw = CP_high_yaw * np.cos((self.yaw+h)*np.pi/180.0)**self.pP
CP_low_yaw = CP_low_yaw * np.cos((self.yaw-h)*np.pi/180.0)**self.pP
CP_high_wind = CP_high_wind * np.cos((self.yaw)*np.pi/180.0)**self.pP
CP_low_wind = CP_low_wind * np.cos((self.yaw)*np.pi/180.0)**self.pP
CT_high_yaw = CT_high_yaw * np.cos((self.yaw+h)*np.pi/180.0)**2
CT_low_yaw = CT_low_yaw * np.cos((self.yaw-h)*np.pi/180.0)**2
CT_high_wind = CT_high_wind * np.cos((self.yaw)*np.pi/180.0)**2
CT_low_wind = CT_low_wind * np.cos((self.yaw)*np.pi/180.0)**2
# compute derivative via central differencing and arrange in sub-matrices of the Jacobian
dCP_dyaw = np.eye(self.nTurbines)*(CP_high_yaw-CP_low_yaw)/(2.0*h)
dCP_dwind = np.eye(self.nTurbines)*(CP_high_wind-CP_low_wind)/(2.0*h)
dCT_dyaw = np.eye(self.nTurbines)*(CT_high_yaw-CT_low_yaw)/(2.0*h)
dCT_dwind = np.eye(self.nTurbines)*(CT_high_wind-CT_low_wind)/(2.0*h)
# compile full Jacobian from sub-matrices
dCP = np.hstack((dCP_dyaw, dCP_dwind))
dCT = np.hstack((dCT_dyaw, dCT_dwind))
J = np.vstack((dCP, dCT))
return J