git.gag.com Git - debian/openrocket/blob - src/net/sf/openrocket/rocketcomponent/SymmetricComponent.java

   1 package net.sf.openrocket.rocketcomponent;
   2
   3 import static net.sf.openrocket.util.MathUtil.pow2;
   4
   5 import java.util.ArrayList;
   6 import java.util.Collection;
   7 import java.util.List;
   8
   9 import net.sf.openrocket.util.Coordinate;
  10 import net.sf.openrocket.util.MathUtil;
  11
  12
  13 /**
  14  * Class for an axially symmetric rocket component generated by rotating
  15  * a function y=f(x) >= 0 around the x-axis (eg. tube, cone, etc.)
  16  *
  17  * @author Sampo Niskanen <sampo.niskanen@iki.fi>
  18  */
  19
  20 public abstract class SymmetricComponent extends BodyComponent implements RadialParent {
  21         public static final double DEFAULT_RADIUS = 0.025;
  22         public static final double DEFAULT_THICKNESS = 0.002;
  23
  24         private static final int DIVISIONS = 100; // No. of divisions when integrating
  25
  26         protected boolean filled = false;
  27         protected double thickness = DEFAULT_THICKNESS;
  28
  29
  30         // Cached data, default values signify not calculated
  31         private double wetArea = -1;
  32         private double planArea = -1;
  33         private double planCenter = -1;
  34         private double volume = -1;
  35         private double fullVolume = -1;
  36         private double longitudinalInertia = -1;
  37         private double rotationalInertia = -1;
  38         private Coordinate cg = null;
  39
  40
  41
  42         public SymmetricComponent() {
  43                 super();
  44         }
  45
  46
  47         /**
  48          * Return the component radius at position x.
  49          * @param x Position on x-axis.
  50          * @return  Radius of the component at the given position, or 0 if outside
  51          *          the component.
  52          */
  53         public abstract double getRadius(double x);
  54
  55         @Override
  56         public abstract double getInnerRadius(double x);
  57
  58         public abstract double getForeRadius();
  59
  60         public abstract boolean isForeRadiusAutomatic();
  61
  62         public abstract double getAftRadius();
  63
  64         public abstract boolean isAftRadiusAutomatic();
  65
  66
  67         // Implement the Radial interface:
  68         @Override
  69         public final double getOuterRadius(double x) {
  70                 return getRadius(x);
  71         }
  72
  73
  74         @Override
  75         public final double getRadius(double x, double theta) {
  76                 return getRadius(x);
  77         }
  78
  79         @Override
  80         public final double getInnerRadius(double x, double theta) {
  81                 return getInnerRadius(x);
  82         }
  83
  84
  85
  86         /**
  87          * Return the component wall thickness.
  88          */
  89         public double getThickness() {
  90                 if (filled)
  91                         return Math.max(getForeRadius(), getAftRadius());
  92                 return Math.min(thickness, Math.max(getForeRadius(), getAftRadius()));
  93         }
  94
  95
  96         /**
  97          * Set the component wall thickness.  Values greater than the maximum radius are not
  98          * allowed, and will result in setting the thickness to the maximum radius.
  99          */
 100         public void setThickness(double thickness) {
 101                 if ((this.thickness == thickness) && !filled)
 102                         return;
 103                 this.thickness = MathUtil.clamp(thickness, 0, Math.max(getForeRadius(), getAftRadius()));
 104                 filled = false;
 105                 fireComponentChangeEvent(ComponentChangeEvent.MASS_CHANGE);
 106         }
 107
 108
 109         /**
 110          * Returns whether the component is set as filled.  If it is set filled, then the
 111          * wall thickness will have no effect.
 112          */
 113         public boolean isFilled() {
 114                 return filled;
 115         }
 116
 117
 118         /**
 119          * Sets whether the component is set as filled.  If the component is filled, then
 120          * the wall thickness will have no effect.
 121          */
 122         public void setFilled(boolean filled) {
 123                 if (this.filled == filled)
 124                         return;
 125                 this.filled = filled;
 126                 fireComponentChangeEvent(ComponentChangeEvent.MASS_CHANGE);
 127         }
 128
 129
 130         /**
 131          * Adds component bounds at a number of points between 0...length.
 132          */
 133         @Override
 134         public Collection<Coordinate> getComponentBounds() {
 135                 List<Coordinate> list = new ArrayList<Coordinate>(20);
 136                 for (int n = 0; n <= 5; n++) {
 137                         double x = n * length / 5;
 138                         double r = getRadius(x);
 139                         addBound(list, x, r);
 140                 }
 141                 return list;
 142         }
 143
 144
 145
 146         /**
 147          * Calculate volume of the component by integrating over the length of the component.
 148          * The method caches the result, so subsequent calls are instant.  Subclasses may
 149          * override this method for simple shapes and use this method as necessary.
 150          *
 151          * @return  The volume of the component.
 152          */
 153         @Override
 154         public double getComponentVolume() {
 155                 if (volume < 0)
 156                         integrate();
 157                 return volume;
 158         }
 159
 160
 161         /**
 162          * Calculate full (filled) volume of the component by integrating over the length
 163          * of the component.  The method caches the result, so subsequent calls are instant.
 164          * Subclasses may override this method for simple shapes and use this method as
 165          * necessary.
 166          *
 167          * @return  The filled volume of the component.
 168          */
 169         public double getFullVolume() {
 170                 if (fullVolume < 0)
 171                         integrate();
 172                 return fullVolume;
 173         }
 174
 175
 176         /**
 177          * Calculate the wetted area of the component by integrating over the length
 178          * of the component.  The method caches the result, so subsequent calls are instant.
 179          * Subclasses may override this method for simple shapes and use this method as
 180          * necessary.
 181          *
 182          * @return  The wetted area of the component.
 183          */
 184         public double getComponentWetArea() {
 185                 if (wetArea < 0)
 186                         integrate();
 187                 return wetArea;
 188         }
 189
 190
 191         /**
 192          * Calculate the planform area of the component by integrating over the length of
 193          * the component.  The method caches the result, so subsequent calls are instant.
 194          * Subclasses may override this method for simple shapes and use this method as
 195          * necessary.
 196          *
 197          * @return  The planform area of the component.
 198          */
 199         public double getComponentPlanformArea() {
 200                 if (planArea < 0)
 201                         integrate();
 202                 return planArea;
 203         }
 204
 205
 206         /**
 207          * Calculate the planform center X-coordinate of the component by integrating over
 208          * the length of the component.  The planform center is defined as
 209          * <pre>   integrate(x*2*r(x)) / planform area  </pre>
 210          * The method caches the result, so subsequent calls are instant.  Subclasses may
 211          * override this method for simple shapes and use this method as necessary.
 212          *
 213          * @return  The planform center of the component.
 214          */
 215         public double getComponentPlanformCenter() {
 216                 if (planCenter < 0)
 217                         integrate();
 218                 return planCenter;
 219         }
 220
 221
 222         /**
 223          * Calculate CG of the component by integrating over the length of the component.
 224          * The method caches the result, so subsequent calls are instant.  Subclasses may
 225          * override this method for simple shapes and use this method as necessary.
 226          *
 227          * @return  The CG+mass of the component.
 228          */
 229         @Override
 230         public Coordinate getComponentCG() {
 231                 if (cg == null)
 232                         integrate();
 233                 return cg;
 234         }
 235
 236
 237         @Override
 238         public double getLongitudinalUnitInertia() {
 239                 if (longitudinalInertia < 0) {
 240                         if (getComponentVolume() > 0.0000001) // == 0.1cm^3
 241                                 integrateInertiaVolume();
 242                         else
 243                                 integrateInertiaSurface();
 244                 }
 245                 return longitudinalInertia;
 246         }
 247
 248
 249         @Override
 250         public double getRotationalUnitInertia() {
 251                 if (rotationalInertia < 0) {
 252                         if (getComponentVolume() > 0.0000001)
 253                                 integrateInertiaVolume();
 254                         else
 255                                 integrateInertiaSurface();
 256                 }
 257                 return rotationalInertia;
 258         }
 259
 260
 261
 262         /**
 263          * Performs integration over the length of the component and updates the cached variables.
 264          */
 265         private void integrate() {
 266                 double x, r1, r2;
 267                 double cgx;
 268
 269                 // Check length > 0
 270                 if (length <= 0) {
 271                         wetArea = 0;
 272                         planArea = 0;
 273                         planCenter = 0;
 274                         volume = 0;
 275                         cg = Coordinate.NUL;
 276                         return;
 277                 }
 278
 279
 280                 // Integrate for volume, CG, wetted area and planform area
 281
 282                 final double l = length / DIVISIONS;
 283                 final double pil = Math.PI * l; // PI * l
 284                 final double pil3 = Math.PI * l / 3; // PI * l/3
 285                 r1 = getRadius(0);
 286                 x = 0;
 287                 wetArea = 0;
 288                 planArea = 0;
 289                 planCenter = 0;
 290                 fullVolume = 0;
 291                 volume = 0;
 292                 cgx = 0;
 293
 294                 for (int n = 1; n <= DIVISIONS; n++) {
 295                         /*
 296                          * r1 and r2 are the two radii
 297                          * x is the position of r1
 298                          * hyp is the length of the hypotenuse from r1 to r2
 299                          * height if the y-axis height of the component if not filled
 300                          */
 301
 302                         r2 = getRadius(x + l);
 303                         final double hyp = MathUtil.hypot(r2 - r1, l);
 304
 305
 306                         // Volume differential elements
 307                         final double dV;
 308                         final double dFullV;
 309
 310                         dFullV = pil3 * (r1 * r1 + r1 * r2 + r2 * r2);
 311                         if (filled || r1 < thickness || r2 < thickness) {
 312                                 // Filled piece
 313                                 dV = dFullV;
 314                         } else {
 315                                 // Hollow piece
 316                                 final double height = thickness * hyp / l;
 317                                 dV = MathUtil.max(pil * height * (r1 + r2 - height), 0);
 318                         }
 319
 320                         // Add to the volume-related components
 321                         volume += dV;
 322                         fullVolume += dFullV;
 323                         cgx += (x + l / 2) * dV;
 324
 325                         // Wetted area ( * PI at the end)
 326                         wetArea += hyp * (r1 + r2);
 327
 328                         // Planform area & center
 329                         final double p = l * (r1 + r2);
 330                         planArea += p;
 331                         planCenter += (x + l / 2) * p;
 332
 333                         // Update for next iteration
 334                         r1 = r2;
 335                         x += l;
 336                 }
 337
 338                 wetArea *= Math.PI;
 339
 340                 if (planArea > 0)
 341                         planCenter /= planArea;
 342
 343                 if (volume < 0.0000000001) { // 0.1 mm^3
 344                         volume = 0;
 345                         cg = new Coordinate(length / 2, 0, 0, 0);
 346                 } else {
 347                         // getComponentMass is safe now
 348                         // Use super.getComponentMass() to ensure only the transition shape mass
 349                         // is used, not the shoulders
 350                         cg = new Coordinate(cgx / volume, 0, 0, super.getComponentMass());
 351                 }
 352         }
 353
 354
 355         /**
 356          * Integrate the longitudinal and rotational inertia based on component volume.
 357          * This method may be used only if the total volume is zero.
 358          */
 359         private void integrateInertiaVolume() {
 360                 double x, r1, r2;
 361
 362                 final double l = length / DIVISIONS;
 363                 final double pil = Math.PI * l; // PI * l
 364                 final double pil3 = Math.PI * l / 3; // PI * l/3
 365
 366                 r1 = getRadius(0);
 367                 x = 0;
 368                 longitudinalInertia = 0;
 369                 rotationalInertia = 0;
 370
 371                 double volume = 0;
 372
 373                 for (int n = 1; n <= DIVISIONS; n++) {
 374                         /*
 375                          * r1 and r2 are the two radii, outer is their average
 376                          * x is the position of r1
 377                          * hyp is the length of the hypotenuse from r1 to r2
 378                          * height if the y-axis height of the component if not filled
 379                          */
 380                         r2 = getRadius(x + l);
 381                         final double outer = (r1 + r2) / 2;
 382
 383
 384                         // Volume differential elements
 385                         final double inner;
 386                         final double dV;
 387
 388                         if (filled || r1 < thickness || r2 < thickness) {
 389                                 inner = 0;
 390                                 dV = pil3 * (r1 * r1 + r1 * r2 + r2 * r2);
 391                         } else {
 392                                 final double hyp = MathUtil.hypot(r2 - r1, l);
 393                                 final double height = thickness * hyp / l;
 394                                 dV = pil * height * (r1 + r2 - height);
 395                                 inner = Math.max(outer - height, 0);
 396                         }
 397
 398                         rotationalInertia += dV * (pow2(outer) + pow2(inner)) / 2;
 399                         longitudinalInertia += dV * ((3 * (pow2(outer) + pow2(inner)) + pow2(l)) / 12
 400                                         + pow2(x + l / 2));
 401
 402                         volume += dV;
 403
 404                         // Update for next iteration
 405                         r1 = r2;
 406                         x += l;
 407                 }
 408
 409                 if (MathUtil.equals(volume, 0)) {
 410                         integrateInertiaSurface();
 411                         return;
 412                 }
 413
 414                 rotationalInertia /= volume;
 415                 longitudinalInertia /= volume;
 416
 417                 // Shift longitudinal inertia to CG
 418                 longitudinalInertia = Math.max(longitudinalInertia - pow2(getComponentCG().x), 0);
 419         }
 420
 421
 422
 423         /**
 424          * Integrate the longitudinal and rotational inertia based on component surface area.
 425          * This method may be used only if the total volume is zero.
 426          */
 427         private void integrateInertiaSurface() {
 428                 double x, r1, r2;
 429
 430                 final double l = length / DIVISIONS;
 431
 432                 r1 = getRadius(0);
 433                 System.out.println(r1);
 434                 x = 0;
 435
 436                 longitudinalInertia = 0;
 437                 rotationalInertia = 0;
 438
 439                 double surface = 0;
 440
 441                 for (int n = 1; n <= DIVISIONS; n++) {
 442                         /*
 443                          * r1 and r2 are the two radii, outer is their average
 444                          * x is the position of r1
 445                          * hyp is the length of the hypotenuse from r1 to r2
 446                          * height if the y-axis height of the component if not filled
 447                          */
 448                         r2 = getRadius(x + l);
 449                         final double hyp = MathUtil.hypot(r2 - r1, l);
 450                         final double outer = (r1 + r2) / 2;
 451
 452                         final double dS = hyp * (r1 + r2) * Math.PI;
 453
 454                         rotationalInertia += dS * pow2(outer);
 455                         longitudinalInertia += dS * ((6 * pow2(outer) + pow2(l)) / 12 + pow2(x + l / 2));
 456
 457                         surface += dS;
 458
 459                         // Update for next iteration
 460                         r1 = r2;
 461                         x += l;
 462                 }
 463
 464                 if (MathUtil.equals(surface, 0)) {
 465                         longitudinalInertia = 0;
 466                         rotationalInertia = 0;
 467                         return;
 468                 }
 469
 470                 longitudinalInertia /= surface;
 471                 rotationalInertia /= surface;
 472
 473                 // Shift longitudinal inertia to CG
 474                 longitudinalInertia = Math.max(longitudinalInertia - pow2(getComponentCG().x), 0);
 475         }
 476
 477
 478
 479
 480         /**
 481          * Invalidates the cached volume and CG information.
 482          */
 483         @Override
 484         protected void componentChanged(ComponentChangeEvent e) {
 485                 super.componentChanged(e);
 486                 if (!e.isOtherChange()) {
 487                         wetArea = -1;
 488                         planArea = -1;
 489                         planCenter = -1;
 490                         volume = -1;
 491                         fullVolume = -1;
 492                         longitudinalInertia = -1;
 493                         rotationalInertia = -1;
 494                         cg = null;
 495                 }
 496         }
 497
 498
 499
 500         ///////////   Auto radius helper methods
 501
 502
 503         /**
 504          * Returns the automatic radius for this component towards the
 505          * front of the rocket.  The automatics will not search towards the
 506          * rear of the rocket for a suitable radius.  A positive return value
 507          * indicates a preferred radius, a negative value indicates that a
 508          * match was not found.
 509          */
 510         protected abstract double getFrontAutoRadius();
 511
 512         /**
 513          * Returns the automatic radius for this component towards the
 514          * end of the rocket.  The automatics will not search towards the
 515          * front of the rocket for a suitable radius.  A positive return value
 516          * indicates a preferred radius, a negative value indicates that a
 517          * match was not found.
 518          */
 519         protected abstract double getRearAutoRadius();
 520
 521
 522
 523         /**
 524          * Return the previous symmetric component, or null if none exists.
 525          * NOTE: This method currently assumes that there are no external
 526          * "pods".
 527          *
 528          * @return      the previous SymmetricComponent, or null.
 529          */
 530         protected final SymmetricComponent getPreviousSymmetricComponent() {
 531                 RocketComponent c;
 532                 for (c = this.getPreviousComponent(); c != null; c = c.getPreviousComponent()) {
 533                         if (c instanceof SymmetricComponent) {
 534                                 return (SymmetricComponent) c;
 535                         }
 536                         if (!(c instanceof Stage) &&
 537                                         (c.relativePosition == RocketComponent.Position.AFTER))
 538                                 return null; // Bad component type as "parent"
 539                 }
 540                 return null;
 541         }
 542
 543         /**
 544          * Return the next symmetric component, or null if none exists.
 545          * NOTE: This method currently assumes that there are no external
 546          * "pods".
 547          *
 548          * @return      the next SymmetricComponent, or null.
 549          */
 550         protected final SymmetricComponent getNextSymmetricComponent() {
 551                 RocketComponent c;
 552                 for (c = this.getNextComponent(); c != null; c = c.getNextComponent()) {
 553                         if (c instanceof SymmetricComponent) {
 554                                 return (SymmetricComponent) c;
 555                         }
 556                         if (!(c instanceof Stage) &&
 557                                         (c.relativePosition == RocketComponent.Position.AFTER))
 558                                 return null; // Bad component type as "parent"
 559                 }
 560                 return null;
 561         }
 562
 563 }