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