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                         this.setFilled(preset.get(ComponentPreset.FILLED));
 156                 }
 157
 158                 super.loadFromPreset(preset);
 159         }
 160
 161
 162
 163
 164         /**
 165          * Calculate volume of the component by integrating over the length of the component.
 166          * The method caches the result, so subsequent calls are instant.  Subclasses may
 167          * override this method for simple shapes and use this method as necessary.
 168          *
 169          * @return  The volume of the component.
 170          */
 171         @Override
 172         public double getComponentVolume() {
 173                 if (volume < 0)
 174                         integrate();
 175                 return volume;
 176         }
 177
 178
 179         /**
 180          * Calculate full (filled) volume of the component by integrating over the length
 181          * of the component.  The method caches the result, so subsequent calls are instant.
 182          * Subclasses may override this method for simple shapes and use this method as
 183          * necessary.
 184          *
 185          * @return  The filled volume of the component.
 186          */
 187         public double getFullVolume() {
 188                 if (fullVolume < 0)
 189                         integrate();
 190                 return fullVolume;
 191         }
 192
 193
 194         /**
 195          * Calculate the wetted area of the component by integrating over the length
 196          * of the component.  The method caches the result, so subsequent calls are instant.
 197          * Subclasses may override this method for simple shapes and use this method as
 198          * necessary.
 199          *
 200          * @return  The wetted area of the component.
 201          */
 202         public double getComponentWetArea() {
 203                 if (wetArea < 0)
 204                         integrate();
 205                 return wetArea;
 206         }
 207
 208
 209         /**
 210          * Calculate the planform area of the component by integrating over the length of
 211          * the component.  The method caches the result, so subsequent calls are instant.
 212          * Subclasses may override this method for simple shapes and use this method as
 213          * necessary.
 214          *
 215          * @return  The planform area of the component.
 216          */
 217         public double getComponentPlanformArea() {
 218                 if (planArea < 0)
 219                         integrate();
 220                 return planArea;
 221         }
 222
 223
 224         /**
 225          * Calculate the planform center X-coordinate of the component by integrating over
 226          * the length of the component.  The planform center is defined as
 227          * <pre>   integrate(x*2*r(x)) / planform area  </pre>
 228          * The method caches the result, so subsequent calls are instant.  Subclasses may
 229          * override this method for simple shapes and use this method as necessary.
 230          *
 231          * @return  The planform center of the component.
 232          */
 233         public double getComponentPlanformCenter() {
 234                 if (planCenter < 0)
 235                         integrate();
 236                 return planCenter;
 237         }
 238
 239
 240         /**
 241          * Calculate CG of the component by integrating over the length of the component.
 242          * The method caches the result, so subsequent calls are instant.  Subclasses may
 243          * override this method for simple shapes and use this method as necessary.
 244          *
 245          * @return  The CG+mass of the component.
 246          */
 247         @Override
 248         public Coordinate getComponentCG() {
 249                 if (cg == null)
 250                         integrate();
 251                 return cg;
 252         }
 253
 254
 255         @Override
 256         public double getLongitudinalUnitInertia() {
 257                 if (longitudinalInertia < 0) {
 258                         if (getComponentVolume() > 0.0000001) // == 0.1cm^3
 259                                 integrateInertiaVolume();
 260                         else
 261                                 integrateInertiaSurface();
 262                 }
 263                 return longitudinalInertia;
 264         }
 265
 266
 267         @Override
 268         public double getRotationalUnitInertia() {
 269                 if (rotationalInertia < 0) {
 270                         if (getComponentVolume() > 0.0000001)
 271                                 integrateInertiaVolume();
 272                         else
 273                                 integrateInertiaSurface();
 274                 }
 275                 return rotationalInertia;
 276         }
 277
 278
 279
 280         /**
 281          * Performs integration over the length of the component and updates the cached variables.
 282          */
 283         private void integrate() {
 284                 double x, r1, r2;
 285                 double cgx;
 286
 287                 // Check length > 0
 288                 if (length <= 0) {
 289                         wetArea = 0;
 290                         planArea = 0;
 291                         planCenter = 0;
 292                         volume = 0;
 293                         cg = Coordinate.NUL;
 294                         return;
 295                 }
 296
 297
 298                 // Integrate for volume, CG, wetted area and planform area
 299
 300                 final double l = length / DIVISIONS;
 301                 final double pil = Math.PI * l; // PI * l
 302                 final double pil3 = Math.PI * l / 3; // PI * l/3
 303                 r1 = getRadius(0);
 304                 x = 0;
 305                 wetArea = 0;
 306                 planArea = 0;
 307                 planCenter = 0;
 308                 fullVolume = 0;
 309                 volume = 0;
 310                 cgx = 0;
 311
 312                 for (int n = 1; n <= DIVISIONS; n++) {
 313                         /*
 314                          * r1 and r2 are the two radii
 315                          * x is the position of r1
 316                          * hyp is the length of the hypotenuse from r1 to r2
 317                          * height if the y-axis height of the component if not filled
 318                          */
 319
 320                         r2 = getRadius(x + l);
 321                         final double hyp = MathUtil.hypot(r2 - r1, l);
 322
 323
 324                         // Volume differential elements
 325                         final double dV;
 326                         final double dFullV;
 327
 328                         dFullV = pil3 * (r1 * r1 + r1 * r2 + r2 * r2);
 329                         if (filled || r1 < thickness || r2 < thickness) {
 330                                 // Filled piece
 331                                 dV = dFullV;
 332                         } else {
 333                                 // Hollow piece
 334                                 final double height = thickness * hyp / l;
 335                                 dV = MathUtil.max(pil * height * (r1 + r2 - height), 0);
 336                         }
 337
 338                         // Add to the volume-related components
 339                         volume += dV;
 340                         fullVolume += dFullV;
 341                         cgx += (x + l / 2) * dV;
 342
 343                         // Wetted area ( * PI at the end)
 344                         wetArea += hyp * (r1 + r2);
 345
 346                         // Planform area & center
 347                         final double p = l * (r1 + r2);
 348                         planArea += p;
 349                         planCenter += (x + l / 2) * p;
 350
 351                         // Update for next iteration
 352                         r1 = r2;
 353                         x += l;
 354                 }
 355
 356                 wetArea *= Math.PI;
 357
 358                 if (planArea > 0)
 359                         planCenter /= planArea;
 360
 361                 if (volume < 0.0000000001) { // 0.1 mm^3
 362                         volume = 0;
 363                         cg = new Coordinate(length / 2, 0, 0, 0);
 364                 } else {
 365                         // getComponentMass is safe now
 366                         // Use super.getComponentMass() to ensure only the transition shape mass
 367                         // is used, not the shoulders
 368                         cg = new Coordinate(cgx / volume, 0, 0, super.getComponentMass());
 369                 }
 370         }
 371
 372
 373         /**
 374          * Integrate the longitudinal and rotational inertia based on component volume.
 375          * This method may be used only if the total volume is zero.
 376          */
 377         private void integrateInertiaVolume() {
 378                 double x, r1, r2;
 379
 380                 final double l = length / DIVISIONS;
 381                 final double pil = Math.PI * l; // PI * l
 382                 final double pil3 = Math.PI * l / 3; // PI * l/3
 383
 384                 r1 = getRadius(0);
 385                 x = 0;
 386                 longitudinalInertia = 0;
 387                 rotationalInertia = 0;
 388
 389                 double vol = 0;
 390
 391                 for (int n = 1; n <= DIVISIONS; n++) {
 392                         /*
 393                          * r1 and r2 are the two radii, outer is their average
 394                          * x is the position of r1
 395                          * hyp is the length of the hypotenuse from r1 to r2
 396                          * height if the y-axis height of the component if not filled
 397                          */
 398                         r2 = getRadius(x + l);
 399                         final double outer = (r1 + r2) / 2;
 400
 401
 402                         // Volume differential elements
 403                         final double inner;
 404                         final double dV;
 405
 406                         if (filled || r1 < thickness || r2 < thickness) {
 407                                 inner = 0;
 408                                 dV = pil3 * (r1 * r1 + r1 * r2 + r2 * r2);
 409                         } else {
 410                                 final double hyp = MathUtil.hypot(r2 - r1, l);
 411                                 final double height = thickness * hyp / l;
 412                                 dV = pil * height * (r1 + r2 - height);
 413                                 inner = Math.max(outer - height, 0);
 414                         }
 415
 416                         rotationalInertia += dV * (pow2(outer) + pow2(inner)) / 2;
 417                         longitudinalInertia += dV * ((3 * (pow2(outer) + pow2(inner)) + pow2(l)) / 12
 418                                         + pow2(x + l / 2));
 419
 420                         vol += dV;
 421
 422                         // Update for next iteration
 423                         r1 = r2;
 424                         x += l;
 425                 }
 426
 427                 if (MathUtil.equals(vol, 0)) {
 428                         integrateInertiaSurface();
 429                         return;
 430                 }
 431
 432                 rotationalInertia /= vol;
 433                 longitudinalInertia /= vol;
 434
 435                 // Shift longitudinal inertia to CG
 436                 longitudinalInertia = Math.max(longitudinalInertia - pow2(getComponentCG().x), 0);
 437         }
 438
 439
 440
 441         /**
 442          * Integrate the longitudinal and rotational inertia based on component surface area.
 443          * This method may be used only if the total volume is zero.
 444          */
 445         private void integrateInertiaSurface() {
 446                 double x, r1, r2;
 447
 448                 final double l = length / DIVISIONS;
 449
 450                 r1 = getRadius(0);
 451                 System.out.println(r1);
 452                 x = 0;
 453
 454                 longitudinalInertia = 0;
 455                 rotationalInertia = 0;
 456
 457                 double surface = 0;
 458
 459                 for (int n = 1; n <= DIVISIONS; n++) {
 460                         /*
 461                          * r1 and r2 are the two radii, outer is their average
 462                          * x is the position of r1
 463                          * hyp is the length of the hypotenuse from r1 to r2
 464                          * height if the y-axis height of the component if not filled
 465                          */
 466                         r2 = getRadius(x + l);
 467                         final double hyp = MathUtil.hypot(r2 - r1, l);
 468                         final double outer = (r1 + r2) / 2;
 469
 470                         final double dS = hyp * (r1 + r2) * Math.PI;
 471
 472                         rotationalInertia += dS * pow2(outer);
 473                         longitudinalInertia += dS * ((6 * pow2(outer) + pow2(l)) / 12 + pow2(x + l / 2));
 474
 475                         surface += dS;
 476
 477                         // Update for next iteration
 478                         r1 = r2;
 479                         x += l;
 480                 }
 481
 482                 if (MathUtil.equals(surface, 0)) {
 483                         longitudinalInertia = 0;
 484                         rotationalInertia = 0;
 485                         return;
 486                 }
 487
 488                 longitudinalInertia /= surface;
 489                 rotationalInertia /= surface;
 490
 491                 // Shift longitudinal inertia to CG
 492                 longitudinalInertia = Math.max(longitudinalInertia - pow2(getComponentCG().x), 0);
 493         }
 494
 495
 496
 497
 498         /**
 499          * Invalidates the cached volume and CG information.
 500          */
 501         @Override
 502         protected void componentChanged(ComponentChangeEvent e) {
 503                 super.componentChanged(e);
 504                 if (!e.isOtherChange()) {
 505                         wetArea = -1;
 506                         planArea = -1;
 507                         planCenter = -1;
 508                         volume = -1;
 509                         fullVolume = -1;
 510                         longitudinalInertia = -1;
 511                         rotationalInertia = -1;
 512                         cg = null;
 513                 }
 514         }
 515
 516
 517
 518         ///////////   Auto radius helper methods
 519
 520
 521         /**
 522          * Returns the automatic radius for this component towards the
 523          * front of the rocket.  The automatics will not search towards the
 524          * rear of the rocket for a suitable radius.  A positive return value
 525          * indicates a preferred radius, a negative value indicates that a
 526          * match was not found.
 527          */
 528         protected abstract double getFrontAutoRadius();
 529
 530         /**
 531          * Returns the automatic radius for this component towards the
 532          * end of the rocket.  The automatics will not search towards the
 533          * front of the rocket for a suitable radius.  A positive return value
 534          * indicates a preferred radius, a negative value indicates that a
 535          * match was not found.
 536          */
 537         protected abstract double getRearAutoRadius();
 538
 539
 540
 541         /**
 542          * Return the previous symmetric component, or null if none exists.
 543          * NOTE: This method currently assumes that there are no external
 544          * "pods".
 545          *
 546          * @return      the previous SymmetricComponent, or null.
 547          */
 548         protected final SymmetricComponent getPreviousSymmetricComponent() {
 549                 RocketComponent c;
 550                 for (c = this.getPreviousComponent(); c != null; c = c.getPreviousComponent()) {
 551                         if (c instanceof SymmetricComponent) {
 552                                 return (SymmetricComponent) c;
 553                         }
 554                         if (!(c instanceof Stage) &&
 555                                         (c.relativePosition == RocketComponent.Position.AFTER))
 556                                 return null; // Bad component type as "parent"
 557                 }
 558                 return null;
 559         }
 560
 561         /**
 562          * Return the next symmetric component, or null if none exists.
 563          * NOTE: This method currently assumes that there are no external
 564          * "pods".
 565          *
 566          * @return      the next SymmetricComponent, or null.
 567          */
 568         protected final SymmetricComponent getNextSymmetricComponent() {
 569                 RocketComponent c;
 570                 for (c = this.getNextComponent(); c != null; c = c.getNextComponent()) {
 571                         if (c instanceof SymmetricComponent) {
 572                                 return (SymmetricComponent) c;
 573                         }
 574                         if (!(c instanceof Stage) &&
 575                                         (c.relativePosition == RocketComponent.Position.AFTER))
 576                                 return null; // Bad component type as "parent"
 577                 }
 578                 return null;
 579         }
 580
 581 }