windSimSuite.blade

1 #!/usr/bin/env python 2 # -*- coding: utf-8 -*- 3 # 4 # This file is part of MBDyn sim suite. 5 # Copyright (C) 2007 André ESPAZE, as part of a Master thesis supervised by 6 # Martin O.L.Hansen (DTU) and Nicolas Chauvat (Logilab) 7 8 # MBDyn sim suite is free software; you can redistribute it and/or modify 9 # it under the terms of the GNU General Public License as published by 10 # the Free Software Foundation; either version 2 of the License, or 11 # (at your option) any later version. 12 # 13 # This program is distributed in the hope that it will be useful, 14 # but WITHOUT ANY WARRANTY; without even the implied warranty of 15 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 16 # GNU General Public License for more details. 17 # 18 # You should have received a copy of the GNU General Public License 19 # along with this program; if not, write to the Free Software 20 # Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. 21 # 22 23 """The blade definition. The blade is first described in an 24 aerodynamics context for the BEM calculations. 25 The final L{Blade} class will also have MBDyn nodes and 26 elements as attributes. The aerodynamics calculations 27 are in fact done by C{Sections} described in 28 L{windSimSuite.sections} while the coupling with MBDyn is assured 29 by the C{BemNode} class from L{windSimSuite.bem_node}. 30 """ 31 import numpy as N 32 from windSimSuite.common import ObjectWithNodesAndElems, Angle 33 from windSimSuite.references import REFERENCE_FRAME, ReferenceFrames 34 from windSimSuite.sections import Section 35 from windSimSuite.bem_node import BemNode 36 37

38 -class AerodynamicsBlade:

39 """The Blade Element Momentum part concerning the blade 40 """ 41

42 - def __init__(self):

43 self.angle = {"wing" : Angle(0., "deg"), 44 "cone" : Angle(0., "deg"), 45 "pitch" : Angle(0., "deg")} 46 self.offset_wing_angle = Angle(0., "deg") 47 self.references = ReferenceFrames() 48 self._add_reference_frames() 49 50 self.radius = None 51 self.correction_radius = None 52 self.skew_section = None 53 54 # Redefined by Blade that manages the user interface 55 self.sections = [] 56 self.nb_sections = 0. 57 58 self.tangential_force = 0. 59 self.normal_force = 0. 60 self.torque = 0. 61 self.power = 0. 62 63 self.force = {} 64 self.normal_vector_wake = N.array([ [1.], 65 [0.], 66 [0.] ])

67

68 - def _add_reference_frames(self):

69 """Add the reference frames for the wing, 70 cone and pitch angles""" 71 for key in ["wing", "cone", "pitch"]: 72 ref = REFERENCE_FRAME[key]() 73 ref.set_angle(self.angle[key]) 74 ref.calculate() 75 self.references.add(ref)

76

77 - def set_radius(self, value):

78 """Set the blade radius""" 79 self.radius = value

80

81 - def set_wing_angle(self, rotor_wing, unit_key):

82 """Set the wind angle for the blade and recalculate 83 the corresponding reference frame. 84 85 @type rotor_wing: a L{Angle} instance 86 @param rotor_wing: the rotor wing angle 87 88 @type unit_key: a string 89 @param unit_key: the 'rad' or 'deg' value""" 90 self.angle["wing"][unit_key] = rotor_wing[unit_key] + \ 91 self.offset_wing_angle[unit_key] 92 self.references["wing"].calculate()

93

94 - def set_pitch_angle(self, value, unit_key):

95 """Set the pitch angle for the blade and recalculate 96 the corresponding reference frame. 97 98 @type value: a float 99 @param value: the pitch value 100 101 @type unit_key: a string 102 @param unit_key: the 'rad' or 'deg' value""" 103 self.angle["pitch"][unit_key] = value 104 self.references["pitch"].calculate()

105

106 - def set_offset_angle(self, value, unit_key="deg"):

107 """Set the blade offset angle relative to the rotor. 108 109 @type value: a float 110 @param value: the offset value 111 112 @type unit_key: a string 113 @param unit_key: the 'rad' or 'deg' value""" 114 self.offset_wing_angle[unit_key] = value

115

116 - def set_correction_radius(self, value):

117 """The correction radius when a problem occurs with 118 the Prandlt's correction""" 119 self.correction_radius = value

120

121 - def set_cone_angle(self, value, unit_key="deg"):

122 """Set the cone angle for the blade and recalculate 123 the corresponding reference frame. 124 125 @type value: a float 126 @param value: the cone value 127 128 @type unit_key: a string 129 @param unit_key: the 'rad' or 'deg' value""" 130 self.angle["cone"][unit_key] = value 131 self.references["cone"].calculate()

132

133 - def set_skew_section(self, skew_section_id):

134 """Set the skew section, used for applying the yaw model 135 on the rotor.""" 136 self.skew_section = self.sections[skew_section_id]

137

138 - def compute_forces_and_power(self, rotational_speed):

139 """Get the normal and tangential forces as well 140 as the torque and power applied on 141 the blade. The calculations are performed from the sections""" 142 for fkey in ["normal", "tangential"]: 143 self.force[fkey] = 0. 144 for section_id in range(self.nb_sections - 1): 145 sec2 = self.sections[section_id +1] 146 sec1 = self.sections[section_id] 147 delta_r = sec2.radial_position - sec1.radial_position 148 acoeff = (sec2.press[fkey] - sec1.press[fkey]) / delta_r 149 bcoeff = (sec1.press[fkey] * sec2.radial_position - \ 150 sec2.press[fkey] * sec1.radial_position) / delta_r 151 square_delta_r = sec2.radial_position**2 - \ 152 sec1.radial_position**2 153 self.force[fkey] += 0.5 * acoeff * square_delta_r \ 154 + bcoeff * delta_r 155 self.normal_force = self.force["normal"] 156 self.tangential_force = self.force["tangential"] 157 158 self.torque = 0. 159 for section_id in range(self.nb_sections - 1): 160 sec2 = self.sections[section_id +1] 161 sec1 = self.sections[section_id] 162 delta_r = sec2.radial_position - sec1.radial_position 163 acoeff = (sec2.press[fkey] - sec1.press[fkey]) / delta_r 164 bcoeff = (sec1.press[fkey] * sec2.radial_position - \ 165 sec2.press[fkey] * sec1.radial_position) / delta_r 166 square_delta_r = sec2.radial_position**2 - \ 167 sec1.radial_position**2 168 threep_delta_r = sec2.radial_position**3 - \ 169 sec1.radial_position**3 170 self.torque += 1./3. * acoeff * threep_delta_r + \ 171 0.5 * bcoeff * square_delta_r 172 self.power = self.torque * rotational_speed["rad/s"]

173

174 - def get_yaw_model_coeff(self, rotor):

175 """Return the coefficient of the yaw model. Used 176 to get the restoring moment when the turbine is yawed""" 177 scalar = N.dot(self.normal_vector_wake.transpose(), 178 self.skew_section.velocity_prime) 179 skew_angle = N.arccos(float(scalar) / \ 180 self.skew_section.velocity_prime_norm) 181 return rotor.skew_radius_ratio * N.tan(skew_angle/2.) * \ 182 N.cos(self.angle["wing"]["rad"] - rotor.deepest_pt_angle) + 1.

183 184

185 -class Blade(AerodynamicsBlade, ObjectWithNodesAndElems):

186 """The blade interface manipulated by the user. 187 """ 188

189 - def __init__(self, idx=0):

190 AerodynamicsBlade.__init__(self) 191 self.idx = idx 192 ObjectWithNodesAndElems.__init__(self, "Blade %i" % idx) 193 194 self.sections = [] 195 self.section_paras = [] 196 self.section_radius = [] 197 self.nb_sections = 0 198 199 self.bem_nodes = [] 200 self.bem_node_paras = [] 201 self.nb_bem_nodes = 0 202 203 self.save_pressures = False 204 205 self.own_para_names += ["idx", "section_paras", 206 "radius", 207 "correction_radius", 208 "angle", 209 "save_pressures", 210 "bem_node_paras"] 211 212 self.coupling_torque = 0. 213 self.will_save("tangential_force", 214 "normal_force", 215 "torque", 216 "power", 217 "coupling_torque")

218

219 - def add_section(self, section):

220 """Add an aerodynamic section to the blade""" 221 self.sections.append(section) 222 self.section_radius.append(section.radial_position) 223 self.nb_sections += 1

224

225 - def add_bem_node(self, bem_node):

226 """Add a bem node to the blade""" 227 self.bem_nodes.append(bem_node) 228 self.nb_bem_nodes += 1

229

230 - def will_save_pressures(self, boolean):

231 """Will save the pressures on the blades""" 232 self.save_pressures = boolean

233

234 - def _set_force_ranges(self):

235 """Set the force range of actions, according 236 to the bem nodes present""" 237 238 if self.nb_bem_nodes == 0: 239 return 240 241 elif self.nb_bem_nodes == 1: 242 self.bem_nodes[0].create_force() 243 self.bem_nodes[0].create_force_areas_from(self.sections) 244 245 elif self.nb_bem_nodes == 2: 246 one = self.bem_nodes[0] 247 two = self.bem_nodes[1] 248 rsep = (two.radial_position - one.radial_position) / 2. \ 249 + one.radial_position 250 251 one.set_lower_radius(0.) 252 one.set_upper_radius(rsep) 253 254 two.set_lower_radius(rsep) 255 two.set_upper_radius(self.radius) 256 257 else: 258 for idx in range(1, self.nb_bem_nodes - 1): 259 lower = self.bem_nodes[idx - 1] 260 current = self.bem_nodes[idx] 261 upper = self.bem_nodes[idx + 1] 262 263 current.set_boundaries(lower.radial_position, 264 upper.radial_position) 265 266 # Only the first and the last have not been current 267 first = self.bem_nodes[0] 268 first.set_boundaries(0., 269 self.bem_nodes[1].lower_radius) 270 271 last = self.bem_nodes[-1] 272 last.set_boundaries(self.bem_nodes[-2].upper_radius, 273 self.radius) 274 275 section_radius = N.array(self.section_radius) 276 for bem_node in self.bem_nodes: 277 bem_node.find_sections(section_radius, 278 self.sections)

279

280 - def init_coupling(self):

281 """Init the coupling between MBDyn and the unsteady BEM 282 code, as well as the result objects.""" 283 if self.nb_sections > 0: 284 self._set_force_ranges()

285

286 - def init_results(self):

287 """Initialize the results objects for the 288 parameters that will be saved.""" 289 self.common_init_results() 290 for section in self.sections: 291 if self.save_pressures: 292 section.will_save("tangential_press", 293 "normal_press") 294 section.init_results()

295

296 - def save(self):

297 """Save the blade status""" 298 self.coupling_torque = 0. 299 #for bem_node in self.bem_nodes: 300 # if bem_node.force: 301 # self.coupling_torque += bem_node.force.tangential_force * \ 302 # bem_node.radial_position 303 self.save_results() 304 for section in self.sections: 305 section.save() 306 for bem_node in self.bem_nodes: 307 bem_node.save()

308

309 - def collect_parameters(self):

310 """Collect the parameters into a dictionary, for saving 311 the blade results in a file""" 312 self.collect_own_parameters() 313 self.collect_node_parameters() 314 self.collect_elt_parameters() 315 for section in self.sections: 316 section.collect_parameters() 317 self.section_paras.append(section.para) 318 for bem_node in self.bem_nodes: 319 bem_node.collect_parameters() 320 self.bem_node_paras.append(bem_node.para)

321

322 - def set_parameters(self, para):

323 """Set the results from C{para} on the blade. 324 This method will also create the sections and 325 bem nodes.""" 326 self.set_own_parameters(para) 327 self.set_node_parameters() 328 self.set_elt_parameters() 329 330 for section_para in self.section_paras: 331 section = Section() 332 section.set_parameters(section_para) 333 self.add_section(section) 334 335 for bem_node_para in self.bem_node_paras: 336 bem_node = BemNode() 337 bem_node.set_parameters(bem_node_para) 338 self.add_bem_node(bem_node)

339

Source Code for Module windSimSuite.blade