00001
00002
00003 #ifndef TURBINE_MECHANISM_H
00004 #define TURBINE_MECHANISM_H
00005
00006 #define MAXBLADES 10
00007
00008
00009 #define TURBINE_LOG_ON
00010
00011 #include "Point.h"
00012 #include "Assembly.h"
00013 #include "DH_Chain.h"
00014 #include "Scene.h"
00015
00021 class TurbineControl : public controlable {
00022 public:
00023 TurbineControl() : yaw_amp(0.), yaw_initial(180.), yaw_rate(0.), ontop(1) {}
00024
00025 ~TurbineControl() {}
00026
00027 virtual void register_at(ControlDevice& Ctrl) {}
00028
00029 virtual void register_at(ControlDevice& Ctrl, const std::string& prefix, int num) {
00030 char VariableName[32];
00031
00032 idNum = num;
00033 std::sprintf(VariableName,"Turbine(%02d)",num+1);
00034
00035 #ifdef TURBINE_DEBUG_ON
00036 std::cout << "Registering " << VariableName << std::endl;
00037 #endif
00038
00039 TLocCtrl = Ctrl.getSubDevice(std::string(VariableName));
00040 RegisterAt(TLocCtrl,"Turbine_name",turbine_name);
00041 RegisterAt(TLocCtrl,"Tower_name",tower_name);
00042 RegisterAt(TLocCtrl,"Nacelle_name",nacelle_name);
00043 RegisterAt(TLocCtrl,"Hub_name",hub_name);
00044 RegisterAt(TLocCtrl,"Tower_file",tower_file);
00045 RegisterAt(TLocCtrl,"Nacelle_file",nacelle_file);
00046 RegisterAt(TLocCtrl,"Hub_file",hub_file);
00047 RegisterAt(TLocCtrl,"Tower_filetype",tower_filetype);
00048 RegisterAt(TLocCtrl,"Nacelle_filetype",nacelle_filetype);
00049 RegisterAt(TLocCtrl,"Hub_filetype",hub_filetype);
00050
00051 RegisterAt(TLocCtrl,"offset",part_offset);
00052
00053 RegisterAt(TLocCtrl,"num_blades",num_blades);
00054 for (register int i = 0; i < MAXBLADES; i++) {
00055 std::sprintf(VariableName,"Blade(%02d)_name",i+1);
00056 RegisterAt(TLocCtrl,std::string(VariableName),blade_name[i]);
00057 std::sprintf(VariableName,"Blade(%02d)_file",i+1);
00058 RegisterAt(TLocCtrl,std::string(VariableName),blade_file[i]);
00059 std::sprintf(VariableName,"Blade(%02d)_filetype",i+1);
00060 RegisterAt(TLocCtrl,std::string(VariableName),blade_filetype[i]);
00061 }
00062
00063 towerBase.register_at(TLocCtrl,"towerBase");
00064 tower_nacelle.register_at(TLocCtrl,"tower-nacelle");
00065 nacelleBase.register_at(TLocCtrl,"nacelleBase");
00066 nacelle_tower.register_at(TLocCtrl,"nacelle_tower");
00067 nacelle_rotor.register_at(TLocCtrl,"nacelle_rotor");
00068 for (register int i=0; i < MAXBLADES; i++) {
00069 std::sprintf(VariableName,"BladeJoint(%02d)",i+1);
00070 blades[i].register_at(TLocCtrl,std::string(VariableName));
00071 }
00072
00073 scene_tower.register_at(TLocCtrl,"scene_tower");
00074 RegisterAt(TLocCtrl,"ontop",ontop);
00075 RegisterAt(TLocCtrl,"xpos",xpos);
00076 RegisterAt(TLocCtrl,"ypos",ypos);
00077
00078 RegisterAt(TLocCtrl,"yaw_initial",yaw_initial);
00079 RegisterAt(TLocCtrl,"yaw_rate",yaw_rate);
00080 RegisterAt(TLocCtrl,"yaw_amplitude",yaw_amp);
00081 RegisterAt(TLocCtrl,"rotation_initial",rotation_initial);
00082 RegisterAt(TLocCtrl,"rotation_rate",rotation_rate);
00083 RegisterAt(TLocCtrl,"pitch_control",pitch_control);
00084 RegisterAt(TLocCtrl,"v_minPitch",v_minPitch);
00085 RegisterAt(TLocCtrl,"alpha",alpha);
00086
00087 for (register int i=0; i < MAXBLADES; i++) {
00088 std::sprintf(VariableName,"theta(%02d)",i+1);
00089 RegisterAt(TLocCtrl,std::string(VariableName),theta[i]);
00090 std::sprintf(VariableName,"pitch_initial(%02d)",i+1);
00091 RegisterAt(TLocCtrl,std::string(VariableName),pitch_initial[i]);
00092 std::sprintf(VariableName,"pitch_rate(%02d)",i+1);
00093 RegisterAt(TLocCtrl,std::string(VariableName),pitch_rate[i]);
00094 std::sprintf(VariableName,"pitch_I(%02d)",i+1);
00095 RegisterAt(TLocCtrl,std::string(VariableName),pitch_I[i]);
00096 }
00097
00098
00099 RegisterAt(TLocCtrl,"I_rotor",I_rotor);
00100 RegisterAt(TLocCtrl,"brake_set",brake_set);
00101 RegisterAt(TLocCtrl,"brake_release",brake_release);
00102
00103 RegisterAt(TLocCtrl,"Wind_avg",wind_avg);
00104 RegisterAt(TLocCtrl,"blade_length",blade_length);
00105
00106 }
00107
00108 void print() {
00109 ( mklog() << "TURBINE CONTROL DEVICE\n"
00110 << "\tnum_blades = " << num_blades << " tower_name = " << tower_name << std::endl).flush();
00111 }
00112
00113 public:
00114 double LastTime;
00115 int Outputs, OutputEvery, CheckPointEvery;
00116 int idNum;
00117
00118 std::string turbine_name, tower_name, tower_file, nacelle_name, nacelle_file, hub_name, hub_file;
00119 std::string blade_name[MAXBLADES], blade_file[MAXBLADES];
00120 int num_blades, pitch_control;
00121 int tower_filetype, nacelle_filetype, hub_filetype, blade_filetype[MAXBLADES];
00122
00123 int ontop;
00124 DataType xpos, ypos, part_offset;
00125
00126 DataType yaw_initial, yaw_rate, yaw_amp, rotation_initial, rotation_rate, theta[MAXBLADES];
00127 DataType pitch_initial[MAXBLADES], pitch_rate[MAXBLADES], pitch_I[MAXBLADES], v_minPitch, alpha;
00128
00129 JointControl towerBase, tower_nacelle, nacelleBase, nacelle_tower, nacelle_rotor, rotorBase, blades[MAXBLADES];
00130
00131 LinkControl scene_tower;
00132
00133 ControlDevice TLocCtrl;
00134
00135 DataType I_rotor, I_driveTrain, I_nacelleTower, I_rotorTower, I_DriveTower;
00136 DataType brake_set, brake_release, wind_avg, blade_length;
00137
00138 };
00139
00141 class SolidSolverSpecific;
00142
00143 class Turbine {
00144 typedef Assembly<DataType,3> AssemblyType;
00145 typedef Scene<DataType,3> SceneType;
00146 typedef Loft<DataType> LoftType;
00147 typedef Joint<DataType,3> JointType;
00148 typedef DH_Link<DataType,3> LinkType;
00149 typedef DH_Chain<DataType,3> ChainType;
00150
00151 public:
00152 friend class Terrain;
00153
00154 enum pitchType { free, collective, linear, relWind };
00155
00156 Turbine() {
00157 p_rotor = 0.;
00158 p_nacelle = 0.;
00159 v_rotor = 0.;
00160 v_nacelle = 0;
00161
00162 for (register int i = 0; i < MAXBLADES; i++) {
00163 a_blades[i] = 0.;
00164 v_blades[i] = 0.;
00165 p_blades[i] = 0.;
00166 }
00167
00168 }
00169
00170 ~Turbine() {}
00171
00172 bool build(int i, std::string name, std::string file, int type) {
00173 Turbine0pt->GetNthPart(i)->read(type,file);
00174 Turbine0pt->GetNthPart(i)->setName(name);
00175 return 1;
00176 }
00177
00178 bool specificSetup(SceneType& ThisScene, TurbineControl * TurbineCtrl ) {
00179 std::string nameTmp;
00180 ( mklog() << "\n\n*********\t Turbine specificSetup \tBegin\n*********" << std::endl ).flush();
00181
00183 nameTmp = "Turbine_Assembly0";
00184 Turbine0pt = new AssemblyType;
00185 Turbine0pt->setName(TurbineCtrl->turbine_name);
00186
00187 ThisScene.AddAssembly(*Turbine0pt);
00188 ThisScene.GetNthSceneAssembly(TurbineCtrl->idNum).print(0);
00189
00190 Turbine0pt = ThisScene.GetNthAssembly(TurbineCtrl->idNum);
00191
00192 Tower0 = new LoftType;
00193 Turbine0pt->AddPart(Tower0);
00194
00195 Turbine0pt->GetNthPart(0)->read(TurbineCtrl->tower_filetype,TurbineCtrl->tower_file);
00196 Turbine0pt->GetNthPart(0)->setName(TurbineCtrl->tower_name);
00197
00198 Tower0 = (LoftType *) Turbine0pt->GetNthPart(0);
00199 Tower0->measure();
00200
00201 Nacel0 = new LoftType;
00202 Turbine0pt->AddPart(Nacel0);
00203 Turbine0pt->GetNthPart(1)->read(TurbineCtrl->nacelle_filetype,TurbineCtrl->nacelle_file);
00204 nameTmp = "Nacel_Loft0";
00205 Turbine0pt->GetNthPart(1)->setName(TurbineCtrl->nacelle_name);
00206 Nacel0 = (LoftType *) Turbine0pt->GetNthPart(1);
00207
00208 Rotor0 = new AssemblyType;
00209 nameTmp = "Rotor_Assembly0";
00210 char VariableName[16];
00211 std::sprintf(VariableName,"Rotor_%02d",TurbineCtrl->idNum+1);
00212 Rotor0->setName(std::string(VariableName));
00213
00214 Turbine0pt->AddSubAssembly(*Rotor0);
00215 Rotor0 = Turbine0pt->GetNthSubAssembly(0);
00216
00217 for (register int ib = 0; ib < TurbineCtrl->num_blades; ib++) {
00218 Blades[ib] = new LoftType;
00219 Rotor0->AddPart(Blades[ib]);
00220 Rotor0->GetNthPart(ib)->read(TurbineCtrl->blade_filetype[ib],TurbineCtrl->blade_file[ib]);
00221 Rotor0->GetNthPart(ib)->setName(TurbineCtrl->blade_name[ib]);
00222 Blades[ib] = (LoftType *) Rotor0->GetNthPart(ib);
00223 }
00224
00225 Hub0 = new LoftType;
00226 Rotor0->AddPart(Hub0);
00227 Rotor0->GetNthPart(Rotor0->GetNumParts()-1 )->read(TurbineCtrl->hub_filetype,TurbineCtrl->hub_file);
00228 nameTmp = "Hub_Loft0";
00229 Rotor0->GetNthPart(Rotor0->GetNumParts()-1 )->setName(TurbineCtrl->hub_name);
00230 Hub0 = (LoftType *) Rotor0->GetNthPart(Rotor0->GetNumParts()-1);
00231
00232 ( mklog() << "\nparts added\n").flush();
00233
00234
00235 PType o(0.), x(0.), y(0.), z(0.);
00236
00237 o(0) = 0.; o(1) = 0.; o(2) = 0.;
00238 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00239 y(0) = 0.; y(1) = 1.; y(2) = 0.;
00240 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00241
00242 o(0) = 0.; o(1) = 0.; o(2) = 0.;
00243 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00244 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00245 JointType * zero;
00246 zero = new JointType(TurbineCtrl->towerBase.origin,TurbineCtrl->towerBase.zaxis,TurbineCtrl->towerBase.xaxis,
00247 JointType::revolute,Tower0->IdTag());
00248
00249 L0 = new LinkType(*ThisScene.GetCoordFrame(),*zero, TurbineCtrl->scene_tower.r, TurbineCtrl->scene_tower.a, TurbineCtrl->scene_tower.d, TurbineCtrl->scene_tower.t);
00250 L0->setId(-1);
00251
00252 o(0) = 0.; o(1) = 0.; o(2) = 31.5;
00253 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00254 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00255 JointType * one;
00256 one = new JointType(o,z,x,JointType::sphericalWrist,Tower0->IdTag());
00257 LinkType * L1;
00258 L1 = new LinkType(*zero,*one);
00259
00260 o(0) = -0.05; o(1) = 0.; o(2) = 0.;
00261 z(0) = 1.; z(1) = 0.; z(2) = 0.;
00262 x(0) = 0.; x(1) = 0.; x(2) = 1.;
00263 JointType * two;
00264 two = new JointType(o,z,x,JointType::revolute,Nacel0->IdTag());
00265 LinkType * L2;
00266 L2 = new LinkType(*one,*two, 0., 0., 0.2, 0.);
00267
00268 o(0) = 0.; o(1) = 0.; o(2) = 2.;
00269 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00270 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00271 JointType * three;
00272 three = new JointType(o,z,x,JointType::revolute,Nacel0->IdTag());
00273 LinkType * L3;
00274 L3 = new LinkType(*two,*three);
00275
00276 o(0) = 0.; o(1) = 0.; o(2) = 1.6;
00277 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00278 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00279 JointType * four;
00280 four = new JointType(o,z,x,JointType::prismatic,Nacel0->IdTag());
00281 LinkType * L4;
00282 L4 = new LinkType(*three,*four);
00283
00284 o(0) = 0.; o(1) = 0.; o(2) = 0.;
00285 z(0) = 1.; z(1) = 0.; z(2) = 0.;
00286 x(0) = 0.; x(1) = 0.; x(2) = 1.;
00287 JointType * five;
00288 five = new JointType(o,z,x,JointType::revolute,Rotor0->IdTag());
00289 LinkType * L5;
00290 L5 = new LinkType(*three,*five, 0., 0., 2.5, 0.);
00292
00293 o(0) = 0.0; o(1) = 0.; o(2) = 0.2;
00294 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00295 x(0) = 0.; x(1) = 1.; x(2) = 0.;
00296 JointType * six;
00297 six = new JointType(o,z,x,JointType::revolute,Rotor0->GetNthPart(0)->IdTag());
00298 LinkType * L6;
00299 L6 = new LinkType(*five,*six);
00300
00301 o(0) = 0.; o(1) = 0.; o(2) = 0.;
00302 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00303 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00304 JointType * six_n;
00305 six_n = new JointType(o,z,x,JointType::revolute,Rotor0->GetNthPart(0)->IdTag());
00306 LinkType * L6_n;
00307 L6_n = new LinkType(*six,*six_n);
00309
00310 o(0) = 0.0; o(1) = 0.; o(2) = 0.2;
00311 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00312 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00313 JointType * seven;
00314 seven = new JointType(o,z,x,JointType::revolute,Rotor0->GetNthPart(1)->IdTag());
00315 LinkType * L7;
00316 L7 = new LinkType(*five,*seven);
00317
00318 o(0) = 0.; o(1) = 0.; o(2) = 0.;
00319 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00320 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00321 JointType * seven_n;
00322 seven_n = new JointType(o,z,x,JointType::revolute,Rotor0->GetNthPart(1)->IdTag());
00323 LinkType * L7_n;
00324 L7_n = new LinkType(*seven,*seven_n);
00326
00327 o(0) = 0.0; o(1) = 0.; o(2) = 0.2;
00328 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00329 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00330 JointType * eight;
00331 eight = new JointType(o,z,x,JointType::revolute,Rotor0->GetNthPart(2)->IdTag());
00332 LinkType * L8;
00333 L8 = new LinkType(*five,*eight);
00334
00335 o(0) = 0.; o(1) = 0.; o(2) = 0.;
00336 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00337 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00338 JointType * eight_n;
00339 eight_n = new JointType(o,z,x,JointType::revolute,Rotor0->GetNthPart(2)->IdTag());
00340 LinkType * L8_n;
00341 L8_n = new LinkType(*eight,*eight_n);
00343
00344 o(0) = 0.0; o(1) = 0.; o(2) = 0.;
00345 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00346 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00347 JointType * nine;
00348 nine = new JointType(o,z,x,JointType::revolute,Rotor0->GetNthPart(3)->IdTag());
00349 LinkType * L9;
00350 L9 = new LinkType(*five,*nine);
00351
00352 o(0) = 0.; o(1) = 0.; o(2) = -0.25;
00353 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00354 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00355 JointType * nine_n;
00356 nine_n = new JointType(o,z,x,JointType::revolute,Rotor0->GetNthPart(3)->IdTag());
00357 LinkType * L9_n;
00358 L9_n = new LinkType(*nine,*nine_n);
00360
00361 ThisScene.addLink(L0);
00362 L0_index = ThisScene.GetNumLinks();
00363
00364 C1 = new ChainType;
00365 C1->AddJoint(zero);
00366 C1->AddJoint(one);
00367 C1->AddJoint(two);
00368 C1->AddJoint(three);
00369 C1->AddJoint(four);
00370 C1->AddJoint(five);
00371 C1->addLink(L1);
00372 C1->addLink(L2);
00373 C1->addLink(L3);
00374 C1->addLink(L5);
00375 C1->addNativeLink(L4);
00376 C1->addTernLink(2,0);
00377 C1->setParentChain(-L0_index);
00378 C1->GetNthLink(0)->setQ_pre(0.0);
00379
00380 ThisScene.AddChain(C1);
00381 C1_index = ThisScene.GetNumChains()-1;
00382
00383 C2 = new ChainType;
00384 C2->AddJoint(five);
00385 C2->AddJoint(six);
00386 C2->AddJoint(six_n);
00387 C2->addLink(L6);
00388 C2->addNativeLink(L6_n);
00389 C2->addTernLink(0,0);
00390 C2->setParentChain(C1_index);
00391
00392 C3 = new ChainType;
00393 C3->AddJoint(five);
00394 C3->AddJoint(seven);
00395 C3->AddJoint(seven_n);
00396 C3->addLink(L7);
00397 C3->addNativeLink(L7_n);
00398 C3->addTernLink(0,0);
00399 C3->setParentChain(C1_index);
00400
00401 C4 = new ChainType;
00402 C4->AddJoint(five);
00403 C4->AddJoint(eight);
00404 C4->AddJoint(eight_n);
00405 C4->addLink(L8);
00406 C4->addNativeLink(L8_n);
00407 C4->addTernLink(0,0);
00408 C4->setParentChain(C1_index);
00409
00410 C5 = new ChainType;
00411 C5->AddJoint(five);
00412 C5->AddJoint(nine);
00413 C5->AddJoint(nine_n);
00414 C5->addLink(L9);
00415 C5->addNativeLink(L9_n);
00416 C5->addTernLink(0,0);
00417 C5->setParentChain(C1_index);
00418
00419 ThisScene.AddChain(C2);
00420 ThisScene.AddChain(C3);
00421 ThisScene.AddChain(C4);
00422 ThisScene.AddChain(C5);
00423
00425 yaw(TurbineCtrl->yaw_initial);
00426
00428 for ( register int i = 0; i < TurbineCtrl->num_blades; i++ ) {
00429 positionBlade(i,TurbineCtrl->theta[i]);
00430 pitch(i,TurbineCtrl->pitch_initial[i]);
00431 }
00433 rotate(TurbineCtrl->rotation_initial );
00434
00436 ( mklog() << "Turbine "<<TurbineCtrl->turbine_name<<"specific setup end\n" ).flush();
00437 return true;
00438 }
00439
00440
00441 void yaw(DataType val) {
00442 C1->GetNthLink(1)->actuate(val,0.,0.,1.);
00443 }
00444
00445 void ontop(SceneType& ThisScene, TurbineControl * TurbineCtrl, DataType offset ) {
00446 PType o(0.), x(0.), y(0.), z(0.), tmp(0.);
00447
00448 o(0) = 0.; o(1) = 0.; o(2) = 0.;
00449 x(0) = 1.; x(1) = 0.; x(2) = 0.;
00450 z(0) = 0.; z(1) = 0.; z(2) = 1.;
00451 ( mklog() << "=================== Joint at base of" << TurbineCtrl->tower_name << " ===================\n").flush();
00452 PType targ;
00453 targ(0) = TurbineCtrl->xpos; targ(1) = TurbineCtrl->ypos;
00454 targ(2) = 100.;
00455 DataType zmin = 0.;
00456 int flag = 0;
00457
00458 o = ThisScene.GetNthSceneAssembly( ThisScene.GetNumAssemblies()-1 ).surfacePoint(targ,zmin,flag,offset);
00459 L0 = ThisScene.GetNthLink(L0_index-1);
00460 ( mklog() << " PLACING ONTOP " << o << std::endl).flush();
00461 L0->getChild()->setOrigin(o);
00462 y = o - ThisScene.GetCoordFrame()->getOrigin();
00463 y = normalize(y);
00464 L0->update();
00465
00466 }
00467
00468 void positionBlade(int n, DataType val) {
00469 PType norm;
00470 norm(0) = 1.; norm(1) = 0.; norm(2) = 0.;
00471
00472 reflect( Blades[n]->GetCPoints(), Blades[n]->GetCPoints_org(), norm );
00473 std::copy (Blades[n]->GetCPoints_org().begin (), Blades[n]->GetCPoints_org().end (), std::back_inserter(Blades[n]->GetCPoints()) );
00474
00475 switch (n) {
00476 case 0 : {
00477 C2->GetNthLink(0)->actuate(val);
00478 break;
00479 }
00480 case 1 : {
00481 C3->GetNthLink(0)->actuate(val);
00482 break;
00483 }
00484 case 2 : {
00485 C4->GetNthLink(0)->actuate(val);
00486 break;
00487 }
00488 }
00489
00490 }
00491
00492 void rotate(DataType time, DataType period) {
00493 C1->GetNthLink(3)->actuate(-(time*DataType(360.0/period)));
00494 }
00495
00496 void rotate(DataType val) {
00497 C1->GetNthLink(3)->actuate(val);
00498 }
00499
00500
00501 void pitch(int n, DataType val) {
00502 switch (n) {
00503 case 0 : {
00504 C2->GetNthNative(0)->actuate(val);
00505 break;
00506 }
00507 case 1 : {
00508 C3->GetNthNative(0)->actuate(val);
00509 break;
00510 }
00511 case 2 : {
00512 C4->GetNthNative(0)->actuate(val);
00513 break;
00514 }
00515 }
00516
00517 }
00518
00519 AssemblyType * GetRotor() { return Rotor0; }
00520
00521 void response(Turbine * ThisTurbine, TurbineControl * TurbineCtrl, double& t, double& dt, int nsteps, int logevery) {
00522
00523 PType rotorCenter, arm, ftmp, ratmp;
00524 DataType rtmp;
00525 MType mtmp;
00526
00527 if ( (t-TurbineCtrl->brake_release) < dt) {
00528 p_rotor = TurbineCtrl->rotation_initial;
00529 for ( register int i = 0; i < MAXBLADES; i++) {
00530 p_blades[i] = TurbineCtrl->pitch_initial[i];
00531 }
00532 }
00533
00534 if ( (TurbineCtrl->brake_set > 0 && t >= TurbineCtrl->brake_release) || (TurbineCtrl->brake_set <= 0 ) ) {
00535 #ifdef TURBINE_DEBUG_ON
00536 if (nsteps % logevery == 0 )
00537 ( mklog() << "\t" << TurbineCtrl->turbine_name << " FSI response at time = " << t << " , dt = " << dt << std::endl ).flush();
00538 #endif
00539 rotorCenter(0) = 0.; rotorCenter(1) = 0.; rotorCenter(2) = 0.;
00540 for (register int b = 0; b < TurbineCtrl->num_blades; b++) {
00541 rotorCenter += ThisTurbine->Blades[b]->Centroid();
00542 }
00543 rotorCenter *= DataType(1.)/DataType(TurbineCtrl->num_blades);
00544 T_rotor = 0.;
00545 for (register int i = 0; i < TurbineCtrl->num_blades; i++) {
00546 arm = ThisTurbine->Blades[i]->Centroid() - rotorCenter;
00547 ftmp = ThisTurbine->Blades[i]->GetNetP();
00548 T_rotor += cross3D(arm, ftmp);
00549 }
00550 ratmp = ThisTurbine->Blades[0]->Centroid() - rotorCenter;
00551 ratmp = cross3D(arm,ratmp);
00552 T_rotor = project(T_rotor,ratmp);
00553 rtmp = projectScalar(T_rotor,ratmp);
00554 a_rotor = ( rtmp/TurbineCtrl->I_rotor );
00555
00556 if ( std::fabs(a_rotor/d2r) > 2.e3) {
00557 #ifdef TURBINE_LOG_ON
00558 ( mklog() << "\t _ACCELERATION LIMIT " << a_rotor/d2r << " reset to " << sgn(a_rotor)*2.e3 << std::endl ).flush();
00559 #endif
00560 a_rotor = sgn(a_rotor)*2.e3*d2r;
00561 }
00562
00563 v_rotor += a_rotor*dt;
00564 p_rotor += (DataType(0.5)*a_rotor*dt*dt + v_rotor*dt)/d2r;
00565
00567 ThisTurbine->rotate(-p_rotor);
00568
00569 #ifdef MOTION_ONLY
00570 DataType t_off = t - TurbineCtrl->brake_release;
00571 ThisTurbine->rotate(TurbineCtrl->rotation_initial+TurbineCtrl->rotation_rate*t_off);
00572 ThisTurbine->yaw(TurbineCtrl->yaw_initial+TurbineCtrl->yaw_amp*std::sin(TurbineCtrl->yaw_rate*d2r*t_off) );
00573 #endif
00574
00576 DataType pBlade = 0., vBlade = 0.;
00577
00579 if ( TurbineCtrl->pitch_control == linear ) {
00580 #ifdef TURBINE_DEBUG_ON
00581 if (nsteps % logevery == 0 )
00582 ( mklog() << "\n\t=-linear Pitch Control v_minPitch " << TurbineCtrl->v_minPitch;
00583 #endif
00584 DataType vRtmp = std::fabs(v_rotor/d2r);
00585 for (register int i = 0; i < TurbineCtrl->num_blades; i++) {
00586 if ( vRtmp > TurbineCtrl->v_minPitch ) {
00587 pBlade = TurbineCtrl->pitch_initial[i] + TurbineCtrl->pitch_rate[i]*(vRtmp-TurbineCtrl->v_minPitch);
00588 #ifdef TURBINE_DEBUG_ON
00589 if (nsteps % logevery == 0 )
00590 ( mklog() << "\t\t" << i << " pitching : " << TurbineCtrl->pitch_initial[i] << " + "
00591 << TurbineCtrl->pitch_rate[i] << " * " << (vRtmp-TurbineCtrl->v_minPitch) << " = " << pBlade << std::endl ).flush();
00592 #endif
00593 ThisTurbine->pitch(i,pBlade);
00594 p_blades[i] = pBlade;
00595 }
00596 }
00597 }
00598 if ( TurbineCtrl->pitch_control == relWind ) {
00599 #ifdef TURBINE_DEBUG_ON
00600 if (nsteps % logevery == 0 )
00601 ( mklog() << "\t| Relative Wind Pitch Control\n";
00602 #endif
00603
00604 DataType wind = TurbineCtrl->wind_avg;
00605 DataType blade_length = TurbineCtrl->blade_length;
00606 DataType blade_speed = blade_length*std::fabs(v_rotor);
00607 pBlade = 180. - std::atan(wind/blade_speed)/d2r;
00608
00609 for (register int i = 0; i < TurbineCtrl->num_blades; i++) {
00610 p_blades[i] = DataType(0.5)*(pBlade + p_blades[i]);
00611 ThisTurbine->pitch(i,p_blades[i]);
00612
00613 }
00614
00615 #ifdef TURBINE_DEBUG_ON
00616 if (nsteps % logevery == 0 )
00617 ( mklog() << "\t blade tip speed [m/s] " << blade_speed << "\t tip speed ratio " << blade_speed/wind
00618 << "\t pitch [deg] : " << p_blades[0] << std::endl ).flush();
00619 #endif
00620
00621 }
00622
00624 else if ( TurbineCtrl->pitch_control == free ) {
00625 for (register int i = 0; i < TurbineCtrl->num_blades; i++) {
00626
00627 arm = ThisTurbine->Blades[i]->Centroid() - rotorCenter;
00628 ftmp = ThisTurbine->Blades[i]->GetMomP();
00629 rtmp = projectScalar(ftmp, arm);
00630 a_blades[i] = rtmp / TurbineCtrl->pitch_I[i];
00631 v_blades[i] += a_blades[i]*dt;
00632 p_blades[i] += DataType(0.5)*a_blades[i]*dt*dt + v_blades[i]*dt;
00633 ThisTurbine->pitch(i,p_blades[i]);
00634 #ifdef TURBINE_DEBUG_ON
00635 if (nsteps % logevery == 0 )
00636 ( mklog() << "\t " << i << " ThisTurbine->Blades[i]->GetMomP() " << ThisTurbine->Blades[i]->GetMomP()
00637 << " \t arm " << arm << " rtmp " << rtmp << std::endl
00638 << "\t\t" << i << " pitching : " << p_blades[i] << std::endl ).flush();
00639 #endif
00640 }
00641 }
00643 else if ( TurbineCtrl->pitch_control == collective ) {
00644 pBlade = 0.; vBlade = 0.;
00645 for (register int i = 0; i < TurbineCtrl->num_blades; i++) {
00646
00647
00648 arm = ThisTurbine->Blades[i]->Centroid() - rotorCenter;
00649 ftmp = ThisTurbine->Blades[i]->GetMomP();
00650 rtmp = projectScalar(ftmp, arm);
00651 a_blades[i] = rtmp / TurbineCtrl->pitch_I[i];
00652 v_blades[i] += a_blades[i]*dt;
00653 p_blades[i] += (DataType(0.5)*a_blades[i]*dt*dt + v_blades[i]*dt)/d2r;
00654 vBlade += v_blades[i];
00655 pBlade += p_blades[i];
00656 #ifdef TURBINE_DEBUG_ON
00657 ( mklog() << "\t p_blades["<<i<<"] " << p_blades[i] << " Blades["<<i<<"]->GetMomP() " << ThisTurbine->Blades[i]->GetMomP()
00658 << " \t arm " << arm << " rtmp " << rtmp << std::endl ).flush();
00659 #endif
00660 }
00661 pBlade /= TurbineCtrl->num_blades;
00662 vBlade /= TurbineCtrl->num_blades;
00663
00664 for (register int i = 0; i < TurbineCtrl->num_blades; i++) {
00665 v_blades[i] = vBlade;
00666 p_blades[i] = pBlade;
00667 ThisTurbine->pitch(i,pBlade);
00668 }
00669 #ifdef TURBINE_DEBUG_ON
00670 ( mklog() << "\n\t\t pitching collectively : " << pBlade << std::endl ).flush();
00671 #endif
00672 }
00673
00674 if ( nsteps % logevery == 0 ) {
00675 char fname[256];
00676 sprintf(fname,"motion_Kinetics.dat");
00677 std::ofstream ofs;
00678
00679 if ( ! fexists(fname) ) {
00680 ofs.open(fname, std::ios::out);
00681 ofs << "# 1 NSteps\t2 time\t3 TurbineName\t4 a_rotor\t5 v_rotor\6 tp_rotor";
00682 for (int b=0; b < TurbineCtrl->num_blades; b++) {
00683 ofs << "\t" << 6+b << " pitch[" << b << "]";
00684 }
00685 ofs << std::endl;
00686 }
00687 else {
00688 ofs.open(fname, std::ios::app);
00689 }
00690 ofs << nsteps << "\t" << t << "\t" << TurbineCtrl->turbine_name << "\t"
00691 << a_rotor/d2r << "\t" << v_rotor/d2r << "\t" << p_rotor << "\t";
00692 #ifdef TURBINE_DEBUG_ON
00693 ( mklog() << "\t " << nsteps << "\t time " << t << "\t" << TurbineCtrl->turbine_name << "\t accel. [deg s^-2] "
00694 << a_rotor/d2r << "\t vel. [deg/s] " << v_rotor/d2r << "\t pos. " << p_rotor << "\n";
00695 #endif
00696 for (int b=0; b < TurbineCtrl->num_blades; b++) {
00697 ofs << "\t" << p_blades[b];
00698 }
00699 ofs << std::endl;
00700 ofs.close();
00701 }
00702
00703 if ( nsteps % logevery == 0 ) {
00704
00705 #ifdef TURBINE_DEBUG_ON
00706 ( mklog() << "\t " << nsteps << "\t time " << t << "\t" << TurbineCtrl->turbine_name << "\t accel. [deg s^-2] "
00707 << a_rotor/d2r << "\t vel. [deg/s] " << v_rotor/d2r << "\t pos. " << p_rotor << "\n";
00708 #endif
00709 }
00710 }
00711 }
00712
00713 virtual void logKinetics(Turbine * ThisTurbine, TurbineControl * TurbineCtrl, std::ofstream& ofs, DataType time_, int steps ) {
00714
00715 PType cen, ftmp, axis;
00716
00717 ofs << "\t\t\"" << TurbineCtrl->turbine_name << "\": {\n";
00718 ofs << "\t\t\t\"rotor\": {\n\t\t\t\t\"acc\": " << a_rotor/d2r;
00719 ofs << ", \"vel\": " << v_rotor/d2r;
00720 ofs << ", \"pos\": " << p_rotor/d2r;
00721 ofs << ", \"colPitch\": " << p_blades[0] << ",";
00722 for (int b=0; b < TurbineCtrl->num_blades; b++) {
00723 cen = ThisTurbine->Blades[b]->Centroid();
00724 ftmp = ThisTurbine->Blades[b]->GetNetP();
00725 axis = ThisTurbine->Blades[b]->GetXaxis();
00726 ofs << "\n\t\t\t\t\"blade" << b << "\": {\"lift\": " << projectScalar(ftmp,axis);
00727 axis = ThisTurbine->Blades[b]->GetYaxis();
00728 ofs << ", \"drag\": " << projectScalar(ftmp,axis);
00729 axis = ThisTurbine->Blades[b]->GetZaxis();
00730 ofs << ", \"radial\": " << projectScalar(ftmp,axis);
00731 ofs << ", \"cenX\": " << cen(0);
00732 ofs << ", \"cenY\": " << cen(1);
00733 ofs << ", \"cenZ\": " << cen(2);
00734 ofs << ", \"ZaxiX\": " << axis(0);
00735 ofs << ", \"ZaxiY\": " << axis(1);
00736 ofs << ", \"ZaxiZ\": " << axis(2) << " }";
00737 if (b<TurbineCtrl->num_blades-1) ofs << ",";
00738 }
00739 ofs << "\n\t\t\t},";
00740 cen = ThisTurbine->Hub0->Centroid();
00741 ftmp = ThisTurbine->Hub0->GetNetP();
00742 axis = ThisTurbine->Hub0->GetXaxis();
00743 ofs << "\n\t\t\t\"hub\": { \"lateral\": " << projectScalar(ftmp,axis);
00744 axis = ThisTurbine->Hub0->GetYaxis();
00745 ofs << ", \"drag\": " << projectScalar(ftmp,axis);
00746 axis = ThisTurbine->Hub0->GetZaxis();
00747 ofs << ", \"lift\": " << projectScalar(ftmp,axis);
00748 ofs << ", \"cenX\": " << cen(0);
00749 ofs << ", \"cenY\": " << cen(1);
00750 ofs << ", \"cenZ\": " << cen(2);
00751 ofs << ", \"axiX\": " << axis(0);
00752 ofs << ", \"axiY\": " << axis(1);
00753 ofs << ", \"axiZ\": " << axis(2) << " },";
00754 cen = ThisTurbine->Nacel0->Centroid();
00755 ftmp = ThisTurbine->Nacel0->GetNetP();
00756 axis = ThisTurbine->Nacel0->GetXaxis();
00757 ofs << "\n\t\t\t\"nacelle\": { \"lateral\": " << projectScalar(ftmp,axis);
00758 axis = ThisTurbine->Nacel0->GetYaxis();
00759 ofs << ", \"drag\": " << projectScalar(ftmp,axis);
00760 axis = ThisTurbine->Nacel0->GetZaxis();
00761 ofs << ", \"lift\": " << projectScalar(ftmp,axis);
00762 ofs << ", \"cenX\": " << cen(0);
00763 ofs << ", \"cenY\": " << cen(1);
00764 ofs << ", \"cenZ\": " << cen(2);
00765 ofs << ", \"ZaxiX\": " << axis(0);
00766 ofs << ", \"ZaxiY\": " << axis(1);
00767 ofs << ", \"ZaxiZ\": " << axis(2) << " },";
00768 cen = ThisTurbine->Tower0->Centroid();
00769 ftmp = ThisTurbine->Tower0->GetNetP();
00770 axis = ThisTurbine->Tower0->GetXaxis();
00771 ofs << "\n\t\t\t\"tower\": { \"lateral\": " << projectScalar(ftmp,axis);
00772 axis = ThisTurbine->Tower0->GetYaxis();
00773 ofs << ", \"drag\": " << projectScalar(ftmp,axis);
00774 axis = ThisTurbine->Tower0->GetZaxis();
00775 ofs << ", \"lift\": " << projectScalar(ftmp,axis);
00776 ofs << ", \"cenX\": " << cen(0);
00777 ofs << ", \"cenY\": " << cen(1);
00778 ofs << ", \"cenZ\": " << cen(2);
00779 ofs << ", \"ZaxiX\": " << axis(0);
00780 ofs << ", \"ZaxiY\": " << axis(1);
00781 ofs << ", \"ZaxiZ\": " << axis(2) << " }";
00782 ofs << "\n\t\t}";
00783
00784 }
00785
00786
00787
00788 void movement(Turbine * ThisTurbine, TurbineControl * TurbineCtrl, double& t, double& dt) {
00789
00790 #ifdef TURBINE_DEBUG_ON
00791 ( mklog() << "\n\n*********\tSolidSolverSpecific movement aDVANCE bEGIN\n*********" << std::endl ).flush();
00792 ( mklog() << " time = " << t << " time*DataType(360.0/6.) = " << t*DataType(360.0/6.)
00793 << " time*DataType(360.0/6.)*d2r = " << t*DataType(360.0/6.)*d2r << std::endl
00794 << " dt = " << dt << std::endl ).flush();
00795 #endif
00796
00797 if ( (TurbineCtrl->brake_set > 0 && t >= TurbineCtrl->brake_release) || (TurbineCtrl->brake_set <= 0 ) ) {
00798 DataType t_off = t - TurbineCtrl->brake_release;
00800 ThisTurbine->yaw(TurbineCtrl->yaw_initial+TurbineCtrl->yaw_amp*std::sin(TurbineCtrl->yaw_rate*d2r*t_off) );
00801
00803 ThisTurbine->rotate(TurbineCtrl->rotation_initial+TurbineCtrl->rotation_rate*t_off);
00804
00805 }
00806
00807 #ifdef TURBINE_DEBUG_ON
00808 ( mklog() << "\n\n*********\tSolidSolverSpecific movement aDVANCE fINISHED\n*********" << std::endl ).flush();
00809 #endif
00810
00811 }
00812
00813
00814 virtual void Restart(std::ifstream& ifs, int& pos, double& t, double& dt) {
00815 #ifdef CHECKRESTART_DEBUG_ON
00816 ( mklog() << "\nRestarting TURBINEMECHANISM SOLVER t= " << t << " dt = " << dt << " pos = " << pos << std::endl ).flush();
00817 #endif
00818
00819 L0->Restart(ifs,pos,t,dt);
00820 pos = ifs.tellg();
00821 C1->Restart(ifs,pos,t,dt);
00822 pos = ifs.tellg();
00823 C2->Restart(ifs,pos,t,dt);
00824 pos = ifs.tellg();
00825 C3->Restart(ifs,pos,t,dt);
00826 pos = ifs.tellg();
00827 C4->Restart(ifs,pos,t,dt);
00828 pos = ifs.tellg();
00829
00830 ifs.read((char*)&L0_index,sizeof(int));
00831 ifs.read((char*)&C1_index,sizeof(int));
00832
00833 ifs.read((char*)&T_rotor(0),sizeof(DataType));
00834 ifs.read((char*)&T_rotor(1),sizeof(DataType));
00835 ifs.read((char*)&T_rotor(2),sizeof(DataType));
00836
00837 ifs.read((char*)&T_nacelle(0),sizeof(DataType));
00838 ifs.read((char*)&T_nacelle(1),sizeof(DataType));
00839 ifs.read((char*)&T_nacelle(2),sizeof(DataType));
00840
00841 ifs.read((char*)&a_rotor,sizeof(DataType));
00842 ifs.read((char*)&v_rotor,sizeof(DataType));
00843 ifs.read((char*)&p_rotor,sizeof(DataType));
00844
00845 ifs.read((char*)&a_nacelle,sizeof(DataType));
00846 ifs.read((char*)&v_nacelle,sizeof(DataType));
00847 ifs.read((char*)&p_nacelle,sizeof(DataType));
00848
00849 pos = ifs.tellg();
00850 #ifdef CHECKRESTART_DEBUG_ON
00851 ( mklog() << "\nDone Restarting TURBINEMECHANISM t = " << t << " dt = " << dt << " pos = " << pos << std::endl ).flush();
00852 #endif
00853 }
00854
00855 virtual void Checkpointing(std::ofstream& ofs) {
00856 #ifdef CHECKRESTART_DEBUG_ON
00857 ( mklog() << "\nCHECKPOINTING TURBINEMECHANISM SOLVER pos = " << ofs.tellp() << std::endl ).flush();
00858 #endif
00859
00860 L0->Checkpointing(ofs);
00861 C1->Checkpointing(ofs);
00862 C2->Checkpointing(ofs);
00863 C3->Checkpointing(ofs);
00864 C4->Checkpointing(ofs);
00865
00866 ofs.write((char*)&L0_index,sizeof(int));
00867 ofs.write((char*)&C1_index,sizeof(int));
00868
00869 ofs.write((char*)&T_rotor(0),sizeof(DataType));
00870 ofs.write((char*)&T_rotor(1),sizeof(DataType));
00871 ofs.write((char*)&T_rotor(2),sizeof(DataType));
00872
00873 ofs.write((char*)&T_nacelle(0),sizeof(DataType));
00874 ofs.write((char*)&T_nacelle(1),sizeof(DataType));
00875 ofs.write((char*)&T_nacelle(2),sizeof(DataType));
00876
00877 ofs.write((char*)&a_rotor,sizeof(DataType));
00878 ofs.write((char*)&v_rotor,sizeof(DataType));
00879 ofs.write((char*)&p_rotor,sizeof(DataType));
00880
00881 ofs.write((char*)&a_nacelle,sizeof(DataType));
00882 ofs.write((char*)&v_nacelle,sizeof(DataType));
00883 ofs.write((char*)&p_nacelle,sizeof(DataType));
00884
00885 #ifdef CHECKRESTART_DEBUG_ON
00886 ( mklog() << "\n+++ CHECKPOINTING TURBINEMECHANISM COMPLETE pos = " << ofs.tellp() <<std::endl ).flush();
00887 #endif
00888 }
00889
00890 protected:
00891
00892 AssemblyType * Turbine0pt;
00893 AssemblyType * Rotor0;
00894 LoftType * Tower0;
00895 LoftType * Nacel0;
00896 LoftType * Blades[MAXBLADES];
00897 LoftType * Hub0;
00898 LinkType * L0;
00899
00900 ChainType * C1;
00901 ChainType * C2;
00902 ChainType * C3;
00903 ChainType * C4;
00904 ChainType * C5;
00905
00906 int L0_index, C1_index;
00907
00908 PType T_rotor, T_nacelle;
00909
00910 DataType a_rotor, a_nacelle;
00911 DataType v_rotor, v_nacelle;
00912 DataType p_rotor, p_nacelle;
00913
00914 DataType a_blades[MAXBLADES];
00915 DataType v_blades[MAXBLADES];
00916 DataType p_blades[MAXBLADES];
00917
00918 };
00919
00920 #endif