1 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
6 Purpose: Encapsulates the rotor object
8 ------------- Copyright (C) 2000 Jon S. Berndt (jon@jsbsim.org) -------------
10 This program is free software; you can redistribute it and/or modify it under
11 the terms of the GNU Lesser General Public License as published by the Free Software
12 Foundation; either version 2 of the License, or (at your option) any later
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
20 You should have received a copy of the GNU Lesser General Public License along with
21 this program; if not, write to the Free Software Foundation, Inc., 59 Temple
22 Place - Suite 330, Boston, MA 02111-1307, USA.
24 Further information about the GNU Lesser General Public License can also be found on
25 the world wide web at http://www.gnu.org.
27 FUNCTIONAL DESCRIPTION
28 --------------------------------------------------------------------------------
31 --------------------------------------------------------------------------------
33 01/01/10 T.Kreitler test implementation
34 11/15/10 T.Kreitler treated flow solver bug, flow and torque calculations
35 simplified, tiploss influence removed from flapping angles
36 01/10/11 T.Kreitler changed to single rotor model
38 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
40 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
47 #include "models/FGPropagate.h"
48 #include "models/FGAtmosphere.h"
49 #include "models/FGAuxiliary.h"
50 #include "models/FGMassBalance.h"
52 #include "input_output/FGXMLElement.h"
59 static const char *IdSrc = "$Id: FGRotor.cpp,v 1.11 2011/01/17 22:09:59 jberndt Exp $";
60 static const char *IdHdr = ID_ROTOR;
62 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
64 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
66 static inline double sqr(double x) { return x*x; }
68 /*%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
70 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%*/
73 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
77 FGRotor::FGRotor(FGFDMExec *exec, Element* rotor_element, int num)
78 : FGThruster(exec, rotor_element, num),
82 dt(0.0), rho(0.002356),
84 // configuration parameters
85 Radius(0.0), BladeNum(0),
87 Sense(1.0), NominalRPM(0.0), ExternalRPM(0), RPMdefinition(0), ExtRPMsource(NULL),
89 BladeChord(0.0), LiftCurveSlope(0.0), BladeTwist(0.0), HingeOffset(0.0),
90 BladeFlappingMoment(0.0), BladeMassMoment(0.0), PolarMoment(0.0),
94 GroundEffectExp(0.0), GroundEffectShift(0.0),
97 LockNumberByRho(0.0), Solidity(0.0),
100 RPM(0.0), Omega(0.0),
103 a0(0.0), a_1(0.0), b_1(0.0), a_dw(0.0), a1s(0.0), b1s(0.0),
105 H_drag(0.0), J_side(0.0), Torque(0.0), C_T(0.0),
107 lambda(-0.001), mu(0.0), nu(0.001), v_induced(0.0),
108 theta_downwash(0.0), phi_downwash(0.0),
111 ControlMap(eMainCtrl),
112 CollectiveCtrl(0.0), LateralCtrl(0.0), LongitudinalCtrl(0.0)
115 FGColumnVector3 location(0.0, 0.0, 0.0), orientation(0.0, 0.0, 0.0);
116 Element *thruster_element;
118 // initialise/set remaining variables
119 SetTransformType(FGForce::tCustom);
120 PropertyManager = exec->GetPropertyManager();
124 dt = exec->GetDeltaT();
125 for (int i=0; i<5; i++) R[i] = 0.0;
126 for (int i=0; i<5; i++) B[i] = 0.0;
129 thruster_element = rotor_element->GetParent()->FindElement("sense");
130 if (thruster_element) {
131 double s = thruster_element->GetDataAsNumber();
133 Sense = -1.0; // 'CW' as seen from above
134 } else if (s < 0.1) {
135 Sense = 0.0; // 'coaxial'
137 Sense = 1.0; // 'CCW' as seen from above
141 thruster_element = rotor_element->GetParent()->FindElement("location");
142 if (thruster_element) {
143 location = thruster_element->FindElementTripletConvertTo("IN");
145 cerr << "No thruster location found." << endl;
148 thruster_element = rotor_element->GetParent()->FindElement("orient");
149 if (thruster_element) {
150 orientation = thruster_element->FindElementTripletConvertTo("RAD");
152 cerr << "No thruster orientation found." << endl;
155 SetLocation(location);
156 SetAnglesToBody(orientation);
157 InvTransform = Transform().Transposed();
160 ControlMap = eMainCtrl;
161 if (rotor_element->FindElement("controlmap")) {
162 string cm = rotor_element->FindElementValue("controlmap");
165 ControlMap = eTailCtrl;
166 } else if (cm == "TANDEM") {
167 ControlMap = eTandemCtrl;
169 cerr << "# found unknown controlmap: '" << cm << "' using main rotor config." << endl;
173 // ExternalRPM -- is the RPM dictated ?
174 if (rotor_element->FindElement("ExternalRPM")) {
176 RPMdefinition = (int) rotor_element->FindElementValueAsNumber("ExternalRPM");
179 // configure the rotor parameters
180 Configure(rotor_element);
182 // shaft representation - a rather simple transform,
183 // but using a matrix is safer.
184 TboToHsr.InitMatrix( 0.0, 0.0, 1.0,
187 HsrToTbo = TboToHsr.Transposed();
189 // smooth out jumps in hagl reported, otherwise the ground effect
190 // calculation would cause jumps too. 1Hz seems sufficient.
191 damp_hagl = Filter(1.0,dt);
193 // enable import-export
200 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
207 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
209 // 5in1: value-fetch-convert-default-return function
211 double FGRotor::ConfigValueConv( Element* el, const string& ename, double default_val,
212 const string& unit, bool tell)
216 double val = default_val;
218 string pname = "*No parent element*";
221 e = el->FindElement(ename);
222 pname = el->GetName();
227 val = e->GetDataAsNumber();
229 val = el->FindElementValueAsNumberConvertTo(ename,unit);
233 cerr << pname << ": missing element '" << ename <<
234 "' using estimated value: " << default_val << endl;
241 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
243 double FGRotor::ConfigValue(Element* el, const string& ename, double default_val, bool tell)
245 return ConfigValueConv(el, ename, default_val, "", tell);
248 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
250 // 1. read configuration and try to fill holes, ymmv
251 // 2. calculate derived parameters and transforms
252 void FGRotor::Configure(Element* rotor_element)
256 const bool yell = true;
257 const bool silent = false;
260 Radius = 0.5 * ConfigValueConv(rotor_element, "diameter", 42.0, "FT", yell);
261 Radius = Constrain(1e-3, Radius, 1e9);
263 BladeNum = (int) ConfigValue(rotor_element, "numblades", 3 , yell);
265 GearRatio = ConfigValue(rotor_element, "gearratio", 1.0, yell);
267 // make sure that v_tip (omega*r) is below 0.7mach ~ 750ft/s
268 estimate = (750.0/Radius)/(2.0*M_PI) * 60.0; // 7160/Radius
269 NominalRPM = ConfigValue(rotor_element, "nominalrpm", estimate, yell);
271 estimate = Constrain(0.07, 2.0/Radius , 0.14); // guess solidity
272 estimate = estimate * M_PI*Radius/BladeNum;
273 BladeChord = ConfigValueConv(rotor_element, "chord", estimate, "FT", yell);
275 LiftCurveSlope = ConfigValue(rotor_element, "liftcurveslope", 6.0); // "1/RAD"
276 BladeTwist = ConfigValueConv(rotor_element, "twist", -0.17, "RAD");
278 HingeOffset = ConfigValueConv(rotor_element, "hingeoffset", 0.05 * Radius, "FT" );
280 estimate = sqr(BladeChord) * sqr(Radius - HingeOffset) * 0.57;
281 BladeFlappingMoment = ConfigValueConv(rotor_element, "flappingmoment", estimate, "SLUG*FT2");
282 BladeFlappingMoment = Constrain(0.001, BladeFlappingMoment, 1e9);
284 // guess mass from moment of a thin stick, and multiply by the blades cg distance
285 estimate = ( 3.0 * BladeFlappingMoment / sqr(Radius) ) * (0.45 * Radius) ;
286 BladeMassMoment = ConfigValue(rotor_element, "massmoment", estimate); // unit is slug-ft
287 BladeMassMoment = Constrain(0.001, BladeMassMoment, 1e9);
289 TipLossB = ConfigValue(rotor_element, "tiplossfactor", 1.0, silent);
291 estimate = 1.1 * BladeFlappingMoment * BladeNum;
292 PolarMoment = ConfigValueConv(rotor_element, "polarmoment", estimate, "SLUG*FT2");
293 PolarMoment = Constrain(0.001, PolarMoment, 1e9);
295 InflowLag = ConfigValue(rotor_element, "inflowlag", 0.2, yell); // fixme, depends on size
296 InflowLag = Constrain(1e-6, InflowLag, 2.0);
300 if (rotor_element->FindElement("cgroundeffect")) {
302 cge = rotor_element->FindElementValueAsNumber("cgroundeffect");
303 cge = Constrain(1e-9, cge, 1.0);
304 gee = 1.0 / ( 2.0*Radius * cge );
305 cerr << "# *** 'cgroundeffect' is defunct." << endl;
306 cerr << "# *** use 'groundeffectexp' with: " << gee << endl;
309 GroundEffectExp = ConfigValue(rotor_element, "groundeffectexp", 0.0);
310 GroundEffectShift = ConfigValueConv(rotor_element, "groundeffectshift", 0.0, "FT");
312 // precalc often used powers
313 R[0]=1.0; R[1]=Radius; R[2]=R[1]*R[1]; R[3]=R[2]*R[1]; R[4]=R[3]*R[1];
314 B[0]=1.0; B[1]=TipLossB; B[2]=B[1]*B[1]; B[3]=B[2]*B[1]; B[4]=B[3]*B[1];
316 // derived parameters
317 LockNumberByRho = LiftCurveSlope * BladeChord * R[4] / BladeFlappingMoment;
318 Solidity = BladeNum * BladeChord / (M_PI * Radius);
323 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
325 // calculate control-axes components of total airspeed at the hub.
326 // sets rotor orientation angle (beta) as side effect. /SH79/ eqn(19-22)
328 FGColumnVector3 FGRotor::hub_vel_body2ca( const FGColumnVector3 &uvw,
329 const FGColumnVector3 &pqr,
330 double a_ic, double b_ic)
332 FGColumnVector3 v_r, v_shaft, v_w;
335 pos = fdmex->GetMassBalance()->StructuralToBody(GetActingLocation());
338 v_shaft = TboToHsr * InvTransform * v_r;
340 beta_orient = atan2(v_shaft(eV),v_shaft(eU));
342 v_w(eU) = v_shaft(eU)*cos(beta_orient) + v_shaft(eV)*sin(beta_orient);
344 v_w(eW) = v_shaft(eW) - b_ic*v_shaft(eU) - a_ic*v_shaft(eV);
349 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
351 // express fuselage angular velocity in control axes /SH79/ eqn(30,31)
353 FGColumnVector3 FGRotor::fus_angvel_body2ca( const FGColumnVector3 &pqr)
355 FGColumnVector3 av_s_fus, av_w_fus;
358 // av_s_fus = BodyToShaft * pqr; /SH79/
359 // BodyToShaft = TboToHsr * InvTransform
360 av_s_fus = TboToHsr * InvTransform * pqr;
362 av_w_fus(eP)= av_s_fus(eP)*cos(beta_orient) + av_s_fus(eQ)*sin(beta_orient);
363 av_w_fus(eQ)= - av_s_fus(eP)*sin(beta_orient) + av_s_fus(eQ)*cos(beta_orient);
364 av_w_fus(eR)= av_s_fus(eR);
370 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
372 // The calculation is a bit tricky because thrust depends on induced velocity,
375 // The flow_scale parameter (ranging from 0.5-1.0) is used to approximate a
376 // reduction of inflow if the helicopter is close to the ground, yielding to
377 // higher thrust, see /TA77/ eqn(10a).
379 void FGRotor::calc_flow_and_thrust( double theta_0, double Uw, double Ww,
383 double ct_over_sigma = 0.0;
384 double c0, ct_l, ct_t0, ct_t1;
387 mu = Uw/(Omega*Radius); // /SH79/ eqn(24)
390 ct_t0 = (1.0/3.0*B[3] + 1.0/2.0 * TipLossB*mu2 - 4.0/(9.0*M_PI) * mu*mu2 ) * theta_0;
391 ct_t1 = (1.0/4.0*B[4] + 1.0/4.0 * B[2]*mu2) * BladeTwist;
393 ct_l = (1.0/2.0*B[2] + 1.0/4.0 * mu2) * lambda; // first time
395 c0 = (LiftCurveSlope/2.0)*(ct_l + ct_t0 + ct_t1) * Solidity;
396 c0 = c0 / ( 2.0 * sqrt( sqr(mu) + sqr(lambda) ) + 1e-15);
398 // replacement for /SH79/ eqn(26).
399 // ref: dnu/dt = 1/tau ( Ct / (2*sqrt(mu^2+lambda^2)) - nu )
400 // taking mu and lambda constant, this integrates to
402 nu = flow_scale * ((nu - c0) * exp(-dt/InflowLag) + c0);
404 // now from nu to lambda, C_T, and Thrust
406 lambda = Ww/(Omega*Radius) - nu; // /SH79/ eqn(25)
408 ct_l = (1.0/2.0*B[2] + 1.0/4.0 * mu2) * lambda;
410 ct_over_sigma = (LiftCurveSlope/2.0)*(ct_l + ct_t0 + ct_t1); // /SH79/ eqn(27)
412 Thrust = BladeNum*BladeChord*Radius*rho*sqr(Omega*Radius) * ct_over_sigma;
414 C_T = ct_over_sigma * Solidity;
415 v_induced = nu * (Omega*Radius);
420 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
422 // The coning angle doesn't apply for teetering rotors, but calculating
423 // doesn't hurt. /SH79/ eqn(29)
425 void FGRotor::calc_coning_angle(double theta_0)
427 double lock_gamma = LockNumberByRho * rho;
429 double a0_l = (1.0/6.0 + 0.04 * mu*mu*mu) * lambda;
430 double a0_t0 = (1.0/8.0 + 1.0/8.0 * mu*mu) * theta_0;
431 double a0_t1 = (1.0/10.0 + 1.0/12.0 * mu*mu) * BladeTwist;
432 a0 = lock_gamma * ( a0_l + a0_t0 + a0_t1);
436 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
438 // Flapping angles relative to control axes /SH79/ eqn(32)
440 void FGRotor::calc_flapping_angles(double theta_0, const FGColumnVector3 &pqr_fus_w)
442 double lock_gamma = LockNumberByRho * rho;
445 double mu2_2 = sqr(mu)/2.0;
446 double t075 = theta_0 + 0.75 * BladeTwist; // common approximation for rectangular blades
448 a_1 = 1.0/(1.0 - mu2_2) * (
449 (2.0*lambda + (8.0/3.0)*t075)*mu
450 + pqr_fus_w(eP)/Omega
451 - 16.0 * pqr_fus_w(eQ)/(lock_gamma*Omega)
454 b_1 = 1.0/(1.0 + mu2_2) * (
456 - pqr_fus_w(eQ)/Omega
457 - 16.0 * pqr_fus_w(eP)/(lock_gamma*Omega)
460 // used in force calc
461 a_dw = 1.0/(1.0 - mu2_2) * (
462 (2.0*lambda + (8.0/3.0)*t075)*mu
463 - 24.0 * pqr_fus_w(eQ)/(lock_gamma*Omega)
464 * ( 1.0 - ( 0.29 * t075 / (C_T/Solidity) ) )
470 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
474 void FGRotor::calc_drag_and_side_forces(double theta_0)
476 double cy_over_sigma ;
477 double t075 = theta_0 + 0.75 * BladeTwist;
479 H_drag = Thrust * a_dw;
482 0.75*b_1*lambda - 1.5*a0*mu*lambda + 0.25*a_1*b_1*mu
483 - a0*a_1*sqr(mu) + (1.0/6.0)*a0*a_1
484 - (0.75*mu*a0 - (1.0/3.0)*b_1 - 0.5*sqr(mu)*b_1)*t075
486 cy_over_sigma *= LiftCurveSlope/2.0;
488 J_side = BladeNum * BladeChord * Radius * rho * sqr(Omega*Radius) * cy_over_sigma;
493 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
495 // Simplified version of /SH79/ eqn(36). Uses an estimate for blade drag
496 // (a new config parameter to come...).
497 // From "Bramwell's Helicopter Dynamics" Â second edition, eqn(3.43) and (3.44)
499 void FGRotor::calc_torque(double theta_0)
501 // estimate blade drag
502 double delta_dr = 0.009 + 0.3*sqr(6.0*C_T/(LiftCurveSlope*Solidity));
504 Torque = rho * BladeNum * BladeChord * delta_dr * sqr(Omega*Radius) * R[2] *
505 (1.0+4.5*sqr(mu))/8.0
506 - (Thrust*lambda + H_drag*mu)*Radius;
511 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
513 // transform rotor forces from control axes to shaft axes, and express
514 // in body axes /SH79/ eqn(40,41)
516 FGColumnVector3 FGRotor::body_forces(double a_ic, double b_ic)
519 - H_drag*cos(beta_orient) - J_side*sin(beta_orient) + Thrust*b_ic,
520 - H_drag*sin(beta_orient) + J_side*cos(beta_orient) + Thrust*a_ic,
523 return HsrToTbo * F_s;
526 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
528 // calculates the additional moments due to hinge offset and handles
531 FGColumnVector3 FGRotor::body_moments(double a_ic, double b_ic)
533 FGColumnVector3 M_s, M_hub, M_h;
536 // cyclic flapping relative to shaft axes /SH79/ eqn(43)
537 a1s = a_1*cos(beta_orient) + b_1*sin(beta_orient) - b_ic;
538 b1s = b_1*cos(beta_orient) - a_1*sin(beta_orient) + a_ic;
540 mf = 0.5 * HingeOffset * BladeNum * Omega*Omega * BladeMassMoment;
544 M_s(eN) = Torque * Sense ;
546 return HsrToTbo * M_s;
549 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
551 void FGRotor::CalcStatePart1(void)
553 double A_IC; // lateral (roll) control in radians
554 double B_IC; // longitudinal (pitch) control in radians
555 double theta_col; // rotor collective pitch in radians
559 FGColumnVector3 UVW_h, PQR_h;
560 FGColumnVector3 vHub_ca, avFus_ca;
562 double h_agl_ft, filtered_hagl = 0.0;
563 double ge_factor = 1.0;
565 // fetch needed values from environment
566 Vt = fdmex->GetAuxiliary()->GetVt(); // total vehicle velocity including wind
567 dt = fdmex->GetDeltaT();
568 rho = fdmex->GetAtmosphere()->GetDensity(); // slugs/ft^3.
569 UVW_h = fdmex->GetAuxiliary()->GetAeroUVW();
570 PQR_h = fdmex->GetAuxiliary()->GetAeroPQR();
571 h_agl_ft = fdmex->GetPropagate()->GetDistanceAGL();
572 // update InvTransform, the rotor orientation could have been altered
573 InvTransform = Transform().Transposed();
575 // handle RPM requirements, calc omega.
576 if (ExternalRPM && ExtRPMsource) {
577 RPM = ExtRPMsource->getDoubleValue() / GearRatio;
580 if (RPM < 1.0) { // kludge, otherwise calculations go bananas
584 Omega = (RPM/60.0)*2.0*M_PI;
586 // set control inputs
588 B_IC = LongitudinalCtrl;
589 theta_col = CollectiveCtrl;
592 if (GroundEffectExp > 1e-5) {
593 if (h_agl_ft<0.0) h_agl_ft = 0.0; // clamp
594 filtered_hagl = damp_hagl.execute(h_agl_ft) + GroundEffectShift;
595 // actual/nominal factor avoids absurd scales at startup
596 ge_factor -= exp(-filtered_hagl*GroundEffectExp) * (RPM / NominalRPM);
597 if (ge_factor<0.5) ge_factor=0.5; // clamp
600 // all set, start calculations
602 vHub_ca = hub_vel_body2ca(UVW_h, PQR_h, A_IC, B_IC);
604 avFus_ca = fus_angvel_body2ca(PQR_h);
606 calc_flow_and_thrust(theta_col, vHub_ca(eU), vHub_ca(eW), ge_factor);
608 calc_coning_angle(theta_col);
610 calc_flapping_angles(theta_col, avFus_ca);
612 calc_drag_and_side_forces(theta_col);
614 calc_torque(theta_col);
616 // Fixme: only valid for a 'decent' rotor
617 theta_downwash = atan2( - UVW_h(eU), v_induced - UVW_h(eW));
618 phi_downwash = atan2( UVW_h(eV), v_induced - UVW_h(eW));
620 vFn = body_forces(A_IC, B_IC);
621 vMn = Transform() * body_moments(A_IC, B_IC);
625 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
627 void FGRotor::CalcStatePart2(double PowerAvailable)
631 double ExcessTorque = PowerAvailable / Omega;
632 double deltaOmega = ExcessTorque / PolarMoment * dt;
633 RPM += deltaOmega/(2.0*M_PI) * 60.0;
634 if (RPM < 0.0) RPM = 0.0; // Engine won't turn backwards
638 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
640 double FGRotor::GetPowerRequired(void)
643 PowerRequired = Torque * Omega;
644 return PowerRequired;
647 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
649 double FGRotor::Calculate(double PowerAvailable)
651 CalcStatePart2(PowerAvailable);
656 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
659 bool FGRotor::BindModel(void)
661 string property_name, base_property_name;
662 base_property_name = CreateIndexedPropertyName("propulsion/engine", EngineNum);
664 property_name = base_property_name + "/rotor-rpm";
665 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetRPM );
667 property_name = base_property_name + "/x-engine-rpm"; // used for RPM eXchange
668 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetEngineRPM );
670 property_name = base_property_name + "/rotor-thrust-lbs"; // might be redundant - check!
671 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetThrust );
673 property_name = base_property_name + "/a0-rad";
674 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetA0 );
676 property_name = base_property_name + "/a1-rad";
677 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetA1 );
679 property_name = base_property_name + "/b1-rad";
680 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetB1 );
682 property_name = base_property_name + "/inflow-ratio";
683 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetLambda );
685 property_name = base_property_name + "/advance-ratio";
686 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetMu );
688 property_name = base_property_name + "/induced-inflow-ratio";
689 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetNu );
691 property_name = base_property_name + "/vi-fps";
692 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetVi );
694 property_name = base_property_name + "/thrust-coefficient";
695 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetCT );
697 property_name = base_property_name + "/torque-lbsft";
698 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetTorque );
700 property_name = base_property_name + "/theta-downwash-rad";
701 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetThetaDW );
703 property_name = base_property_name + "/phi-downwash-rad";
704 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetPhiDW );
706 switch (ControlMap) {
708 property_name = base_property_name + "/antitorque-ctrl-rad";
709 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetCollectiveCtrl, &FGRotor::SetCollectiveCtrl);
712 property_name = base_property_name + "/tail-collective-ctrl-rad";
713 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetCollectiveCtrl, &FGRotor::SetCollectiveCtrl);
714 property_name = base_property_name + "/lateral-ctrl-rad";
715 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetLateralCtrl, &FGRotor::SetLateralCtrl);
716 property_name = base_property_name + "/longitudinal-ctrl-rad";
717 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetLongitudinalCtrl, &FGRotor::SetLongitudinalCtrl);
719 default: // eMainCtrl
720 property_name = base_property_name + "/collective-ctrl-rad";
721 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetCollectiveCtrl, &FGRotor::SetCollectiveCtrl);
722 property_name = base_property_name + "/lateral-ctrl-rad";
723 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetLateralCtrl, &FGRotor::SetLateralCtrl);
724 property_name = base_property_name + "/longitudinal-ctrl-rad";
725 PropertyManager->Tie( property_name.c_str(), this, &FGRotor::GetLongitudinalCtrl, &FGRotor::SetLongitudinalCtrl);
729 if (RPMdefinition == -1) {
730 property_name = base_property_name + "/x-rpm-dict";
731 ExtRPMsource = PropertyManager->GetNode(property_name, true);
732 } else if (RPMdefinition >= 0 && RPMdefinition != EngineNum) {
733 string ipn = CreateIndexedPropertyName("propulsion/engine", RPMdefinition);
734 property_name = ipn + "/x-engine-rpm";
735 ExtRPMsource = PropertyManager->GetNode(property_name, false);
736 if (! ExtRPMsource) {
737 cerr << "# Warning: Engine number " << EngineNum << "." << endl;
738 cerr << "# No 'x-engine-rpm' property found for engine " << RPMdefinition << "." << endl;
739 cerr << "# Please check order of engine definitons." << endl;
742 cerr << "# Engine number " << EngineNum;
743 cerr << ", given ExternalRPM value '" << RPMdefinition << "' unhandled." << endl;
750 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
752 string FGRotor::GetThrusterLabels(int id, string delimeter)
757 buf << Name << " RPM (engine " << id << ")";
763 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
765 string FGRotor::GetThrusterValues(int id, string delimeter)
776 //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
777 // The bitmasked value choices are as follows:
778 // unset: In this case (the default) JSBSim would only print
779 // out the normally expected messages, essentially echoing
780 // the config files as they are read. If the environment
781 // variable is not set, debug_lvl is set to 1 internally
782 // 0: This requests JSBSim not to output any messages
784 // 1: This value explicity requests the normal JSBSim
786 // 2: This value asks for a message to be printed out when
787 // a class is instantiated
788 // 4: When this value is set, a message is displayed when a
789 // FGModel object executes its Run() method
790 // 8: When this value is set, various runtime state variables
791 // are printed out periodically
792 // 16: When set various parameters are sanity checked and
793 // a message is printed out when they go out of bounds
795 void FGRotor::Debug(int from)
797 string ControlMapName;
799 if (debug_lvl <= 0) return;
801 if (debug_lvl & 1) { // Standard console startup message output
802 if (from == 0) { // Constructor
803 cout << "\n Rotor Name: " << Name << endl;
804 cout << " Diameter = " << 2.0 * Radius << " ft." << endl;
805 cout << " Number of Blades = " << BladeNum << endl;
806 cout << " Gear Ratio = " << GearRatio << endl;
807 cout << " Sense = " << Sense << endl;
808 cout << " Nominal RPM = " << NominalRPM << endl;
811 if (RPMdefinition == -1) {
812 cout << " RPM is controlled externally" << endl;
814 cout << " RPM source set to engine " << RPMdefinition << endl;
818 cout << " Blade Chord = " << BladeChord << endl;
819 cout << " Lift Curve Slope = " << LiftCurveSlope << endl;
820 cout << " Blade Twist = " << BladeTwist << endl;
821 cout << " Hinge Offset = " << HingeOffset << endl;
822 cout << " Blade Flapping Moment = " << BladeFlappingMoment << endl;
823 cout << " Blade Mass Moment = " << BladeMassMoment << endl;
824 cout << " Polar Moment = " << PolarMoment << endl;
825 cout << " Inflow Lag = " << InflowLag << endl;
826 cout << " Tip Loss = " << TipLossB << endl;
827 cout << " Lock Number = " << LockNumberByRho * 0.002356 << " (SL)" << endl;
828 cout << " Solidity = " << Solidity << endl;
830 switch (ControlMap) {
831 case eTailCtrl: ControlMapName = "Tail Rotor"; break;
832 case eTandemCtrl: ControlMapName = "Tandem Rotor"; break;
833 default: ControlMapName = "Main Rotor";
835 cout << " Control Mapping = " << ControlMapName << endl;
839 if (debug_lvl & 2 ) { // Instantiation/Destruction notification
840 if (from == 0) cout << "Instantiated: FGRotor" << endl;
841 if (from == 1) cout << "Destroyed: FGRotor" << endl;
843 if (debug_lvl & 4 ) { // Run() method entry print for FGModel-derived objects
845 if (debug_lvl & 8 ) { // Runtime state variables
847 if (debug_lvl & 16) { // Sanity checking
849 if (debug_lvl & 64) {
850 if (from == 0) { // Constructor
851 cout << IdSrc << endl;
852 cout << IdHdr << endl;
859 } // namespace JSBSim