From: Bill Kuker Date: Fri, 10 Apr 2009 15:08:12 +0000 (+0000) Subject: Changed to non-sampled x-section area X-Git-Url: https://git.gag.com/?a=commitdiff_plain;h=d24efc3712c345a96a1ab7439098bf45c9c827c6;p=sw%2Fmotorsim Changed to non-sampled x-section area --- diff --git a/src/com/billkuker/rocketry/motorsim/grain/ExtrudedGrain.java b/src/com/billkuker/rocketry/motorsim/grain/ExtrudedGrain.java index 70ae575..0924300 100644 --- a/src/com/billkuker/rocketry/motorsim/grain/ExtrudedGrain.java +++ b/src/com/billkuker/rocketry/motorsim/grain/ExtrudedGrain.java @@ -10,9 +10,11 @@ import java.awt.Rectangle; import java.awt.Shape; import java.awt.geom.AffineTransform; import java.awt.geom.Ellipse2D; +import java.awt.geom.GeneralPath; import java.awt.geom.PathIterator; import java.awt.geom.Rectangle2D; import java.util.HashSet; +import java.util.Iterator; import java.util.NoSuchElementException; import java.util.Set; import java.util.SortedMap; @@ -37,267 +39,292 @@ import com.billkuker.rocketry.motorsim.Grain; import com.billkuker.rocketry.motorsim.visual.Chart; public class ExtrudedGrain implements Grain, Grain.Graphical { - - Set plus = new HashSet(); - Set minus = new HashSet(); + + Set plus = new HashSet(); + + Set minus = new HashSet(); + Set inhibited = new HashSet(); - + Amount length = Amount.valueOf(25, SI.MILLIMETER); - + Amount rStep; - - private class RegEntry{ + + private class RegEntry { Amount surfaceArea; + Amount volume; } - + SortedMap, RegEntry> data = new TreeMap, RegEntry>(); - + Amount webThickness; - + { - /* Similar test grain*/ - Shape outside = new Ellipse2D.Double(50,50,30,30); - plus.add(outside); - minus.add(new Ellipse2D.Double(50,60,10,10)); - inhibited.add(outside); - length = Amount.valueOf(70, SI.MILLIMETER); - /*/ - - /*Big c-slot - Shape outside = new Ellipse2D.Double(0,0,30,30); + /* + * Similar test grain Shape outside = new Ellipse2D.Double(50,50,30,30); + * plus.add(outside); minus.add(new Ellipse2D.Double(50,60,10,10)); + * inhibited.add(outside); length = Amount.valueOf(70, SI.MILLIMETER); / + */ + + /* Big c-slot */ + Shape outside = new Ellipse2D.Double(0, 0, 30, 30); plus.add(outside); inhibited.add(outside); - minus.add(new Rectangle2D.Double(12,12,6,20)); + minus.add(new Rectangle2D.Double(12, 25, 5, 10)); + minus.add(new Ellipse2D.Double(10, 4, 5, 5)); length = Amount.valueOf(70, SI.MILLIMETER); - */ - - /*Plus sign - Shape outside = new Ellipse2D.Double(0,0,200,200); - plus.add(outside); - inhibited.add(outside); - minus.add(new Rectangle2D.Double(90,40,20,120)); - minus.add(new Rectangle2D.Double(40,90,120,20)); - */ - + /**/ + + /* + * Plus sign Shape outside = new Ellipse2D.Double(0,0,200,200); + * plus.add(outside); inhibited.add(outside); minus.add(new + * Rectangle2D.Double(90,40,20,120)); minus.add(new + * Rectangle2D.Double(40,90,120,20)); + */ + findWebThickness(); - + fillInData(); } - + @Override public Amount surfaceArea(Amount regression) { - if ( regression.isGreaterThan(webThickness) ) + if (regression.isGreaterThan(webThickness)) return Amount.valueOf(0, Area.UNIT); - Amount highKey = data.tailMap(regression).firstKey(); + Amount highKey = data.tailMap(regression).firstKey(); Amount lowKey; - try{ - lowKey = data.headMap(regression).lastKey(); - } catch ( NoSuchElementException e ){ + try { + lowKey = data.headMap(regression).lastKey(); + } catch (NoSuchElementException e) { return data.get(highKey).surfaceArea; } - - double lp = regression.minus(lowKey).divide(highKey.minus(lowKey)).to(Dimensionless.UNIT).doubleValue(Dimensionless.UNIT); - + + double lp = regression.minus(lowKey).divide(highKey.minus(lowKey)).to( + Dimensionless.UNIT).doubleValue(Dimensionless.UNIT); + Amount lowVal = data.get(lowKey).surfaceArea; Amount highVal = data.get(highKey).surfaceArea; - - return lowVal.times(1-lp).plus(highVal.times(lp)); + + return lowVal.times(1 - lp).plus(highVal.times(lp)); } - + @Override public Amount volume(Amount regression) { - if ( regression.isGreaterThan(webThickness) ) + if (regression.isGreaterThan(webThickness)) return Amount.valueOf(0, Volume.UNIT); - Amount highKey = data.tailMap(regression).firstKey(); + Amount highKey = data.tailMap(regression).firstKey(); Amount lowKey; - try{ - lowKey = data.headMap(regression).lastKey(); - } catch ( NoSuchElementException e ){ + try { + lowKey = data.headMap(regression).lastKey(); + } catch (NoSuchElementException e) { return data.get(highKey).volume; } - - double lp = regression.minus(lowKey).divide(highKey.minus(lowKey)).to(Dimensionless.UNIT).doubleValue(Dimensionless.UNIT); - + + double lp = regression.minus(lowKey).divide(highKey.minus(lowKey)).to( + Dimensionless.UNIT).doubleValue(Dimensionless.UNIT); + Amount lowVal = data.get(lowKey).volume; Amount highVal = data.get(highKey).volume; - - return lowVal.times(1-lp).plus(highVal.times(lp)); + + return lowVal.times(1 - lp).plus(highVal.times(lp)); } - - private void fillInData(){ + + private void fillInData() { double max = webThickness().doubleValue(SI.MILLIMETER); - double delta = max / 5; + double delta = max / 100; rStep = Amount.valueOf(delta, SI.MILLIMETER).divide(10); - for( double r = 0; r <= max+1; r += delta){ + for (double r = 0; r <= max + 1; r += delta) { RegEntry e = new RegEntry(); Amount regression = Amount.valueOf(r, SI.MILLIMETER); - Amount rLen = length.minus(regression.times(2)); - if ( rLen.isLessThan(Amount.valueOf(0, SI.MILLIMETER))) + if (rLen.isLessThan(Amount.valueOf(0, SI.MILLIMETER))) break; - + System.out.println("Calculating area for regression " + regression); - + java.awt.geom.Area burn = getArea(regression); - if ( burn.isEmpty() ) + if (burn.isEmpty()) break; - burn.subtract(getArea(regression.plus(Amount.valueOf(.001, SI.MILLIMETER)))); - - Amount xSection = crossSectionArea(getArea(regression)); + burn.subtract(getArea(regression.plus(Amount.valueOf(.001, + SI.MILLIMETER)))); + + //Amount xSection = crossSectionArea(getArea(regression)); + Amount xSection = crossSectionArea(regression); + e.volume = xSection.times(rLen).to(Volume.UNIT); - - e.surfaceArea = perimeter(burn).divide(2).times(rLen).plus(xSection.times(2)).to(Area.UNIT); - + + e.surfaceArea = perimeter(burn).divide(2).times(rLen).plus( + xSection.times(2)).to(Area.UNIT); + data.put(regression, e); } - + RegEntry e = new RegEntry(); e.surfaceArea = Amount.valueOf(0, Area.UNIT); e.volume = Amount.valueOf(0, Volume.UNIT); - data.put( webThickness(), e ); + data.put(webThickness(), e); } + - private Amount crossSectionArea(java.awt.geom.Area a ) { - Rectangle r = a.getBounds(); - int cnt = 0; - int total = 0; - double e = .1; - for( double x = r.getMinX(); x < r.getMaxX(); x += e){ - for( double y = r.getMinY(); y < r.getMaxY(); y += e){ - if ( a.contains(x, y) ) - cnt++; - total++; - } + private Amount crossSectionArea(Amount regression) { + //Get the PLUS shape and sum its area + java.awt.geom.Area plus = getPlus(regression); + Amount plusArea = Amount.valueOf(0, SI.SQUARE_METRE); + for (java.awt.geom.Area a : separate(plus)) { + plusArea = plusArea.plus(area(a)); } - System.out.println(cnt + " out of " + total); - return Amount.valueOf(r.getWidth() * r.getHeight() * cnt / total, SI.MILLIMETER.pow(2)).to(Area.UNIT); - + + //Get the MINUS shape, intersect it with PLUS to get just the parts + //that are removed, sum it's area + java.awt.geom.Area minus = getMinus(regression); + minus.intersect(plus); + Amount minusArea = Amount.valueOf(0, SI.SQUARE_METRE); + for (java.awt.geom.Area a : separate(minus)) { + minusArea = minusArea.plus(area(a)); + } + + //Subtract PLUS from MINUS and return + Amount area = plusArea.minus(minusArea); + System.out.println(area.to(SI.MILLIMETER.pow(2))); + return area; } - - - - @Override public Amount webThickness() { return webThickness; } - - private Amount perimeter(java.awt.geom.Area a){ + + private Amount perimeter(java.awt.geom.Area a) { + //TODO: I think I need to handle seg_close!! PathIterator i = a.getPathIterator(new AffineTransform(), .001); - double x=0, y=0; + double x = 0, y = 0; double len = 0; - while ( !i.isDone() ){ + while (!i.isDone()) { double coords[] = new double[6]; int type = i.currentSegment(coords); - if ( type == PathIterator.SEG_LINETO ){ - //System.out.println("Line"); + if (type == PathIterator.SEG_LINETO) { + // System.out.println("Line"); double nx = coords[0]; double ny = coords[1]; - //System.out.println(x+","+y+ " to " + nx+"," + ny); - len += Math.sqrt(Math.pow(x-nx, 2) + Math.pow(y-ny,2)); + // System.out.println(x+","+y+ " to " + nx+"," + ny); + len += Math.sqrt(Math.pow(x - nx, 2) + Math.pow(y - ny, 2)); x = nx; y = ny; - } else if ( type == PathIterator.SEG_MOVETO ){ - //System.out.println("Move"); + } else if (type == PathIterator.SEG_MOVETO) { + // System.out.println("Move"); x = coords[0]; y = coords[1]; } else { - //System.err.println("Got " + type); + // System.err.println("Got " + type); } i.next(); } return Amount.valueOf(len, SI.MILLIMETER); } - - private void findWebThickness(){ + + private void findWebThickness() { java.awt.geom.Area a = getArea(Amount.valueOf(0, SI.MILLIMETER)); Rectangle r = a.getBounds(); - double max = r.getWidth() regression){ + + private java.awt.geom.Area getArea(Amount regression) { + java.awt.geom.Area res = getPlus(regression); + res.subtract(getMinus(regression)); + return res; + } + + private java.awt.geom.Area getPlus(Amount regression) { + java.awt.geom.Area a = new java.awt.geom.Area(); + for (Shape s : plus) + a.add(new java.awt.geom.Area(regress(s, regression + .doubleValue(SI.MILLIMETER), true))); + return a; + } + + private java.awt.geom.Area getMinus(Amount regression) { java.awt.geom.Area a = new java.awt.geom.Area(); - for ( Shape s: plus ) - a.add(new java.awt.geom.Area(regress(s, regression.doubleValue(SI.MILLIMETER), true))); - for ( Shape s: minus ) - a.subtract(new java.awt.geom.Area(regress(s, regression.doubleValue(SI.MILLIMETER), false))); + for (Shape s : minus) + a.add(new java.awt.geom.Area(regress(s, regression + .doubleValue(SI.MILLIMETER), false))); return a; } - private Shape regress(Shape s, double mm, boolean plus){ - if ( inhibited.contains(s) ) + + private Shape regress(Shape s, double mm, boolean plus) { + if (inhibited.contains(s)) return s; - if ( s instanceof Ellipse2D ){ - Ellipse2D e = (Ellipse2D)s; - - double d = plus?-2*mm:2*mm; - + if (s instanceof Ellipse2D) { + Ellipse2D e = (Ellipse2D) s; + + double d = plus ? -2 * mm : 2 * mm; + double w = e.getWidth() + d; double h = e.getHeight() + d; - double x = e.getX() - d/2; - double y = e.getY() - d/2; - - return new Ellipse2D.Double(x,y,w,h); - } else if ( s instanceof Rectangle2D ){ - Rectangle2D r = (Rectangle2D)s; - - double d = plus?-2*mm:2*mm; - + double x = e.getX() - d / 2; + double y = e.getY() - d / 2; + + return new Ellipse2D.Double(x, y, w, h); + } else if (s instanceof Rectangle2D) { + Rectangle2D r = (Rectangle2D) s; + + double d = plus ? -2 * mm : 2 * mm; + double w = r.getWidth() + d; double h = r.getHeight() + d; - double x = r.getX() - d/2; - double y = r.getY() - d/2; - - return new Rectangle2D.Double(x,y,w,h); + double x = r.getX() - d / 2; + double y = r.getY() - d / 2; + + return new Rectangle2D.Double(x, y, w, h); } return null; } - public static void main(String args[]) throws Exception { ExtrudedGrain e = new ExtrudedGrain(); - new GrainPanel(e).show(); - + new GrainPanel(e).show(); } @Override public void draw(Graphics2D g2d, Amount regression) { - + java.awt.geom.Area reg = getArea(regression); java.awt.geom.Area burn = getArea(regression); - burn.subtract(getArea(regression.plus(Amount.valueOf(.001, SI.MILLIMETER)))); - java.awt.geom.Area noreg = getArea(Amount.valueOf(0,SI.MILLIMETER)); - + burn.subtract(getArea(regression.plus(Amount.valueOf(.001, + SI.MILLIMETER)))); + java.awt.geom.Area noreg = getArea(Amount.valueOf(0, SI.MILLIMETER)); + Rectangle bounds = noreg.getBounds(); - g2d.scale(200/bounds.getWidth(), 200/bounds.getHeight()); + g2d.scale(200 / bounds.getWidth(), 200 / bounds.getHeight()); g2d.translate(-bounds.getX(), -bounds.getY()); - g2d.setStroke(new BasicStroke(0.5f)); - + g2d.setColor(Color.GRAY); g2d.fill(reg); g2d.setColor(Color.RED); @@ -305,7 +332,116 @@ public class ExtrudedGrain implements Grain, Grain.Graphical { g2d.setColor(Color.BLACK); g2d.draw(noreg); + } + /* + * Separate an area into multiple distinct area. + * Area CAN NOT HAVE HOLES. HOLES WILL BE RETURNED AS AREAS, + * SO A DONUT WILL TURN INTO TWO CIRCLES. + */ + private Set separate(java.awt.geom.Area a) { + Set res = new HashSet(); + PathIterator i = a.getPathIterator(new AffineTransform()); + GeneralPath cur = null; + + while (!i.isDone()) { + double coords[] = new double[6]; + int type = i.currentSegment(coords); + switch (type) { + case PathIterator.SEG_CLOSE: + cur.closePath(); + if (cur != null ){ + java.awt.geom.Area area = new java.awt.geom.Area(cur); + if ( !a.isEmpty() ) + res.add(area); + } + cur = new GeneralPath(i.getWindingRule()); + break; + case PathIterator.SEG_MOVETO: + if (cur != null ){ + java.awt.geom.Area area = new java.awt.geom.Area(cur); + if ( !a.isEmpty() ) + res.add(area); + } + cur = new GeneralPath(i.getWindingRule()); + cur.moveTo(coords[0], coords[1]); + break; + case PathIterator.SEG_CUBICTO: + cur.curveTo(coords[0], coords[1], coords[2], coords[3], + coords[4], coords[5]); + break; + case PathIterator.SEG_LINETO: + cur.lineTo(coords[0], coords[1]); + break; + case PathIterator.SEG_QUADTO: + cur.quadTo(coords[0], coords[1], coords[2], coords[3]); + break; + + } + i.next(); + } + + return res; + } + + /* + * Return the Area of a singular polygon (NO HOLES OR DISJOINT PARTS). + * Coordinates assumed to be in MM. + * http://valis.cs.uiuc.edu/~sariel/research/CG/compgeom/msg00831.html + * http://stackoverflow.com/questions/451426/how-do-i-calculate-the-surface-area-of-a-2d-polygon + * http://www.wikihow.com/Calculate-the-Area-of-a-Polygon + */ + private Amount area(java.awt.geom.Area a) { + if ( !a.isSingular() ) + throw new IllegalArgumentException("Can not calculate area of non-singular shape!"); + PathIterator i = a.getPathIterator(new AffineTransform(), .001); + + + double x = 0, y = 0, sx = 0, sy = 0; + double nx, ny; + double area = 0; + while (!i.isDone()) { + double coords[] = new double[6]; + int type = i.currentSegment(coords); + switch( type ){ + case PathIterator.SEG_CLOSE: + //Go back to the start + nx = sx; + ny = sy; + area += x * ny; + area -= y * nx; + break; + case PathIterator.SEG_LINETO: + nx = coords[0]; + ny = coords[1]; + area += x * ny; + area -= y * nx; + + //Remember the last points + x = nx; + y = ny; + + break; + case PathIterator.SEG_MOVETO: + //Remember the starting point + x = sx = coords[0]; + y = sy = coords[1]; + break; + default: + System.err.println("Got " + type); + } + i.next(); + } + + area = area / 2.0; // Result so far is double the signed area + + if ( area < 0 ) //Depending on winding it could be negative + area = area * -1.0; + + System.out.println(area); + + return Amount.valueOf(area, SI.MILLIMETER.pow(2)).to(Area.UNIT); + } }