src/FDM/YASim/Wing.cpp

   1 #include "Math.hpp"
   2 #include "Surface.hpp"
   3 #include "Wing.hpp"
   4 namespace yasim {
   5
   6 Wing::Wing()
   7 {
   8     _mirror = false;
   9     _base[0] = _base[1] = _base[2] = 0;
  10     _length = 0;
  11     _chord = 0;
  12     _taper = 0;
  13     _sweep = 0;
  14     _dihedral = 0;
  15     _stall = 0;
  16     _stallWidth = 0;
  17     _stallPeak = 0;
  18     _twist = 0;
  19     _camber = 0;
  20     _incidence = 0;
  21     _inducedDrag = 1;
  22     _dragScale = 1;
  23     _liftRatio = 1;
  24     _flap0Start = 0;
  25     _flap0End = 0;
  26     _flap0Lift = 0;
  27     _flap0Drag = 0;
  28     _flap1Start = 0;
  29     _flap1End = 0;
  30     _flap1Lift = 0;
  31     _flap1Drag = 0;
  32     _spoilerStart = 0;
  33     _spoilerEnd = 0;
  34     _spoilerLift = 0;
  35     _spoilerDrag = 0;
  36     _slatStart = 0;
  37     _slatEnd = 0;
  38     _slatAoA = 0;
  39     _slatDrag = 0;
  40 }
  41
  42 Wing::~Wing()
  43 {
  44     int i;
  45     for(i=0; i<_surfs.size(); i++) {
  46         SurfRec* s = (SurfRec*)_surfs.get(i);
  47         delete s->surface;
  48         delete s;
  49     }
  50 }
  51
  52 int Wing::numSurfaces()
  53 {
  54     return _surfs.size();
  55 }
  56
  57 Surface* Wing::getSurface(int n)
  58 {
  59     return ((SurfRec*)_surfs.get(n))->surface;
  60 }
  61
  62 float Wing::getSurfaceWeight(int n)
  63 {
  64     return ((SurfRec*)_surfs.get(n))->weight;
  65 }
  66
  67 void Wing::setMirror(bool mirror)
  68 {
  69     _mirror = mirror;
  70 }
  71
  72 void Wing::setBase(float* base)
  73 {
  74     int i;
  75     for(i=0; i<3; i++) _base[i] = base[i];
  76 }
  77
  78 void Wing::setLength(float length)
  79 {
  80     _length = length;
  81 }
  82
  83 void Wing::setChord(float chord)
  84 {
  85     _chord = chord;
  86 }
  87
  88 void Wing::setTaper(float taper)
  89 {
  90     _taper = taper;
  91 }
  92
  93 void Wing::setSweep(float sweep)
  94 {
  95     _sweep = sweep;
  96 }
  97
  98 void Wing::setDihedral(float dihedral)
  99 {
 100     _dihedral = dihedral;
 101 }
 102
 103 void Wing::setStall(float aoa)
 104 {
 105     _stall = aoa;
 106 }
 107
 108 void Wing::setStallWidth(float angle)
 109 {
 110     _stallWidth = angle;
 111 }
 112
 113 void Wing::setStallPeak(float fraction)
 114 {
 115     _stallPeak = fraction;
 116 }
 117
 118 void Wing::setTwist(float angle)
 119 {
 120     _twist = angle;
 121 }
 122
 123 void Wing::setCamber(float camber)
 124 {
 125     _camber = camber;
 126 }
 127
 128 void Wing::setIncidence(float incidence)
 129 {
 130     _incidence = incidence;
 131     int i;
 132     for(i=0; i<_surfs.size(); i++)
 133         ((SurfRec*)_surfs.get(i))->surface->setIncidence(incidence);
 134 }
 135
 136 void Wing::setFlap0(float start, float end, float lift, float drag)
 137 {
 138     _flap0Start = start;
 139     _flap0End = end;
 140     _flap0Lift = lift;
 141     _flap0Drag = drag;
 142 }
 143
 144 void Wing::setFlap1(float start, float end, float lift, float drag)
 145 {
 146     _flap1Start = start;
 147     _flap1End = end;
 148     _flap1Lift = lift;
 149     _flap1Drag = drag;
 150 }
 151
 152 void Wing::setSlat(float start, float end, float aoa, float drag)
 153 {
 154     _slatStart = start;
 155     _slatEnd = end;
 156     _slatAoA = aoa;
 157     _slatDrag = drag;
 158 }
 159
 160 void Wing::setSpoiler(float start, float end, float lift, float drag)
 161 {
 162     _spoilerStart = start;
 163     _spoilerEnd = end;
 164     _spoilerLift = lift;
 165     _spoilerDrag = drag;
 166 }
 167
 168 void Wing::setFlap0(float lval, float rval)
 169 {
 170     lval = Math::clamp(lval, -1, 1);
 171     rval = Math::clamp(rval, -1, 1);
 172     int i;
 173     for(i=0; i<_flap0Surfs.size(); i++) {
 174         ((Surface*)_flap0Surfs.get(i))->setFlap(lval);
 175         if(_mirror) ((Surface*)_flap0Surfs.get(++i))->setFlap(rval);
 176     }
 177 }
 178
 179 void Wing::setFlap1(float lval, float rval)
 180 {
 181     lval = Math::clamp(lval, -1, 1);
 182     rval = Math::clamp(rval, -1, 1);
 183     int i;
 184     for(i=0; i<_flap1Surfs.size(); i++) {
 185         ((Surface*)_flap1Surfs.get(i))->setFlap(lval);
 186         if(_mirror) ((Surface*)_flap1Surfs.get(++i))->setFlap(rval);
 187     }
 188 }
 189
 190 void Wing::setSpoiler(float lval, float rval)
 191 {
 192     lval = Math::clamp(lval, 0, 1);
 193     rval = Math::clamp(rval, 0, 1);
 194     int i;
 195     for(i=0; i<_spoilerSurfs.size(); i++) {
 196         ((Surface*)_spoilerSurfs.get(i))->setSpoiler(lval);
 197         if(_mirror) ((Surface*)_spoilerSurfs.get(++i))->setSpoiler(rval);
 198     }
 199 }
 200
 201 void Wing::setSlat(float val)
 202 {
 203     val = Math::clamp(val, 0, 1);
 204     int i;
 205     for(i=0; i<_slatSurfs.size(); i++)
 206         ((Surface*)_slatSurfs.get(i))->setSlat(val);
 207 }
 208
 209 float Wing::getGroundEffect(float* posOut)
 210 {
 211     int i;
 212     for(i=0; i<3; i++) posOut[i] = _base[i];
 213     float span = _length * Math::cos(_sweep) * Math::cos(_dihedral);
 214     span = 2*(span + Math::abs(_base[2]));
 215     return span;
 216 }
 217
 218 void Wing::getTip(float* tip)
 219 {
 220     tip[0] = -Math::tan(_sweep);
 221     tip[1] = Math::cos(_dihedral);
 222     tip[2] = Math::sin(_dihedral);
 223     Math::unit3(tip, tip);
 224     Math::mul3(_length, tip, tip);
 225     Math::add3(_base, tip, tip);
 226 }
 227
 228 bool Wing::isMirrored()
 229 {
 230     return _mirror;
 231 }
 232
 233 void Wing::compile()
 234 {
 235     // Have we already been compiled?
 236     if(_surfs.size() != 0) return;
 237
 238     // Assemble the start/end coordinates of all control surfaces
 239     // and the wing itself into an array, sort them,
 240     // and remove duplicates.  This gives us the boundaries of our
 241     // segments.
 242     float bounds[10];
 243     bounds[0] = _flap0Start;   bounds[1] = _flap0End;
 244     bounds[2] = _flap1Start;   bounds[3] = _flap1End;
 245     bounds[4] = _spoilerStart; bounds[5] = _spoilerEnd;
 246     bounds[6] = _slatStart;    bounds[7] = _slatEnd;
 247     //and don't forget the root and the tip of the wing itself
 248     bounds[8] = 0;             bounds[9] = 1;
 249
 250     // Sort in increasing order
 251     int i;
 252     for(i=0; i<10; i++) {
 253         int minIdx = i;
 254         float minVal = bounds[i];
 255         int j;
 256         for(j=i+1; j<10; j++) {
 257             if(bounds[j] < minVal) {
 258                 minIdx = j;
 259                 minVal = bounds[j];
 260             }
 261         }
 262         float tmp = bounds[i];
 263         bounds[i] = minVal; bounds[minIdx] = tmp;
 264     }
 265
 266     // Uniqify
 267     float last = bounds[0];
 268     int nbounds = 1;
 269     for(i=1; i<10; i++) {
 270         if(bounds[i] != last)
 271             bounds[nbounds++] = bounds[i];
 272         last = bounds[i];
 273     }
 274
 275     // Calculate a "nominal" segment length equal to an average chord,
 276     // normalized to lie within 0-1 over the length of the wing.
 277     float segLen = _chord * (0.5f*(_taper+1)) / _length;
 278
 279     // Generating a unit vector pointing out the left wing.
 280     float left[3];
 281     left[0] = -Math::tan(_sweep);
 282     left[1] = Math::cos(_dihedral);
 283     left[2] = Math::sin(_dihedral);
 284     Math::unit3(left, left);
 285
 286     // Calculate coordinates for the root and tip of the wing
 287     float root[3], tip[3];
 288     Math::set3(_base, root);
 289     Math::set3(left, tip);
 290     Math::mul3(_length, tip, tip);
 291     Math::add3(root, tip, tip);
 292
 293     // The wing's Y axis will be the "left" vector.  The Z axis will
 294     // be perpendicular to this and the local (!) X axis, because we
 295     // want motion along the local X axis to be zero AoA (i.e. in the
 296     // wing's XY plane) by definition.  Then the local X coordinate is
 297     // just Y cross Z.
 298     float orient[9], rightOrient[9];
 299     float *x = orient, *y = orient+3, *z = orient+6;
 300     x[0] = 1; x[1] = 0; x[2] = 0;
 301     Math::set3(left, y);
 302     Math::cross3(x, y, z);
 303     Math::unit3(z, z);
 304     Math::cross3(y, z, x);
 305
 306     if(_mirror) {
 307         // Derive the right side orientation matrix from this one.
 308         int i;
 309         for(i=0; i<9; i++)  rightOrient[i] = orient[i];
 310
 311         // Negate all Y coordinates, this gets us a valid basis, but
 312         // it's left handed!  So...
 313         for(i=1; i<9; i+=3) rightOrient[i] = -rightOrient[i];
 314
 315         // Change the direction of the Y axis to get back to a
 316         // right-handed system.
 317         for(i=3; i<6; i++)  rightOrient[i] = -rightOrient[i];
 318     }
 319
 320     // Now go through each boundary and make segments
 321     for(i=0; i<(nbounds-1); i++) {
 322         float start = bounds[i];
 323         float end = bounds[i+1];
 324         float mid = (start+end)/2;
 325
 326         bool flap0=0, flap1=0, slat=0, spoiler=0;
 327         if(_flap0Start   < mid && mid < _flap0End)   flap0 = 1;
 328         if(_flap1Start   < mid && mid < _flap1End)   flap1 = 1;
 329         if(_slatStart    < mid && mid < _slatEnd)    slat = 1;
 330         if(_spoilerStart < mid && mid < _spoilerEnd) spoiler = 1;
 331
 332         // FIXME: Should probably detect an error here if both flap0
 333         // and flap1 are set.  Right now flap1 overrides.
 334
 335         int nSegs = (int)Math::ceil((end-start)/segLen);
 336         if (_twist != 0 && nSegs < 8) // more segments if twisted
 337             nSegs = 8;
 338         float segWid = _length * (end - start)/nSegs;
 339
 340         int j;
 341         for(j=0; j<nSegs; j++) {
 342             float frac = start + (j+0.5f) * (end-start)/nSegs;
 343             float pos[3];
 344             interp(root, tip, frac, pos);
 345
 346             float chord = _chord * (1 - (1-_taper)*frac);
 347
 348             Surface *s = newSurface(pos, orient, chord,
 349                                     flap0, flap1, slat, spoiler);
 350
 351             SurfRec *sr = new SurfRec();
 352             sr->surface = s;
 353             sr->weight = chord * segWid;
 354             s->setTotalDrag(sr->weight);
 355             s->setTwist(_twist * frac);
 356             _surfs.add(sr);
 357
 358             if(_mirror) {
 359                 pos[1] = -pos[1];
 360                 s = newSurface(pos, rightOrient, chord,
 361                                flap0, flap1, slat, spoiler);
 362                 sr = new SurfRec();
 363                 sr->surface = s;
 364                 sr->weight = chord * segWid;
 365                 s->setTotalDrag(sr->weight);
 366                 s->setTwist(_twist * frac);
 367                 _surfs.add(sr);
 368             }
 369         }
 370     }
 371
 372     // Last of all, re-set the incidence in case setIncidence() was
 373     // called before we were compiled.
 374     setIncidence(_incidence);
 375 }
 376
 377 float Wing::getDragScale()
 378 {
 379     return _dragScale;
 380 }
 381
 382 void Wing::setDragScale(float scale)
 383 {
 384     _dragScale = scale;
 385     int i;
 386     for(i=0; i<_surfs.size(); i++) {
 387         SurfRec* s = (SurfRec*)_surfs.get(i);
 388         s->surface->setTotalDrag(scale * s->weight);
 389     }
 390 }
 391
 392 void Wing::setLiftRatio(float ratio)
 393 {
 394     _liftRatio = ratio;
 395     int i;
 396     for(i=0; i<_surfs.size(); i++)
 397         ((SurfRec*)_surfs.get(i))->surface->setZDrag(ratio);
 398 }
 399
 400 float Wing::getLiftRatio()
 401 {
 402     return _liftRatio;
 403 }
 404
 405 Surface* Wing::newSurface(float* pos, float* orient, float chord,
 406                           bool flap0, bool flap1, bool slat, bool spoiler)
 407 {
 408     Surface* s = new Surface();
 409
 410     s->setPosition(pos);
 411     s->setOrientation(orient);
 412     s->setChord(chord);
 413
 414     // Camber is expressed as a fraction of stall peak, so convert.
 415     s->setBaseZDrag(_camber*_stallPeak);
 416
 417     // The "main" (i.e. normal) stall angle
 418     float stallAoA = _stall - _stallWidth/4;
 419     s->setStall(0, stallAoA);
 420     s->setStallWidth(0, _stallWidth);
 421     s->setStallPeak(0, _stallPeak);
 422
 423     // The negative AoA stall is the same if we're using an uncambered
 424     // airfoil, otherwise a "little badder".
 425     if(_camber > 0) {
 426         s->setStall(1, stallAoA * 0.8f);
 427         s->setStallWidth(1, _stallWidth * 0.5f);
 428     } else {
 429         s->setStall(1, stallAoA);
 430         s->setStall(1, _stallWidth);
 431     }
 432
 433     // The "reverse" stalls are unmeasurable junk.  Just use 13deg and
 434     // "sharp".
 435     s->setStallPeak(1, 1);
 436     int i;
 437     for(i=2; i<4; i++) {
 438         s->setStall(i, 0.2267f);
 439         s->setStallWidth(i, 0.01);
 440     }
 441
 442     if(flap0)   s->setFlapParams(_flap0Lift, _flap0Drag);
 443     if(flap1)   s->setFlapParams(_flap1Lift, _flap1Drag);
 444     if(slat)    s->setSlatParams(_slatAoA, _slatDrag);
 445     if(spoiler) s->setSpoilerParams(_spoilerLift, _spoilerDrag);
 446
 447     if(flap0)   _flap0Surfs.add(s);
 448     if(flap1)   _flap1Surfs.add(s);
 449     if(slat)    _slatSurfs.add(s);
 450     if(spoiler) _spoilerSurfs.add(s);
 451
 452     s->setInducedDrag(_inducedDrag);
 453
 454     return s;
 455 }
 456
 457 void Wing::interp(float* v1, float* v2, float frac, float* out)
 458 {
 459     out[0] = v1[0] + frac*(v2[0]-v1[0]);
 460     out[1] = v1[1] + frac*(v2[1]-v1[1]);
 461     out[2] = v1[2] + frac*(v2[2]-v1[2]);
 462 }
 463
 464 }; // namespace yasim