git.gag.com Git - debian/openrocket/blob - src/net/sf/openrocket/util/Analysis.java

   1 package net.sf.openrocket.util;
   2
   3 import static net.sf.openrocket.aerodynamics.AtmosphericConditions.GAMMA;
   4 import static net.sf.openrocket.aerodynamics.AtmosphericConditions.R;
   5
   6 import java.io.File;
   7 import java.io.PrintStream;
   8 import java.util.Arrays;
   9
  10 import net.sf.openrocket.aerodynamics.AerodynamicCalculator;
  11 import net.sf.openrocket.aerodynamics.AerodynamicForces;
  12 import net.sf.openrocket.aerodynamics.AtmosphericConditions;
  13 import net.sf.openrocket.aerodynamics.BarrowmanCalculator;
  14 import net.sf.openrocket.aerodynamics.ExactAtmosphericConditions;
  15 import net.sf.openrocket.aerodynamics.FlightConditions;
  16 import net.sf.openrocket.document.OpenRocketDocument;
  17 import net.sf.openrocket.file.GeneralRocketLoader;
  18 import net.sf.openrocket.file.RocketLoadException;
  19 import net.sf.openrocket.file.RocketLoader;
  20 import net.sf.openrocket.rocketcomponent.Configuration;
  21
  22 public class Analysis {
  23
  24         private static final double MACH_MIN = 0.01;
  25         private static final double MACH_MAX = 5.00001;
  26         private static final double MACH_STEP = 0.02;
  27
  28         private static final double AOA_MACH = 0.6;
  29         private static final double AOA_MIN = 0;
  30         private static final double AOA_MAX = 15.00001*Math.PI/180;
  31         private static final double AOA_STEP = 0.5*Math.PI/180;
  32
  33         private static final double REYNOLDS = 9.8e6;
  34         private static final double STAG_TEMP = 330;
  35
  36
  37         private final RocketLoader loader = new GeneralRocketLoader();
  38         private final AerodynamicCalculator calculator = new BarrowmanCalculator();
  39
  40         private final FlightConditions conditions;
  41         private final double length;
  42
  43         private final Configuration config;
  44
  45         private final AtmosphericConditions atmosphere;
  46
  47
  48
  49         private Analysis(String filename) throws RocketLoadException {
  50
  51                 OpenRocketDocument doc = loader.load(new File(filename));
  52                 config = doc.getRocket().getDefaultConfiguration();
  53
  54                 calculator.setConfiguration(config);
  55
  56                 conditions = new FlightConditions(config);
  57                 System.out.println("Children: " + Arrays.toString(config.getRocket().getChildren()));
  58                 System.out.println("Children: " + Arrays.toString(config.getRocket().getChild(0).getChildren()));
  59                 length = config.getLength();
  60                 System.out.println("Rocket length: " + (length*1000)+"mm");
  61
  62                 atmosphere = new ExactAtmosphericConditions();
  63
  64         }
  65
  66
  67         private double computeVelocityAndAtmosphere(double mach, double reynolds, double stagTemp) {
  68                 final double temperature;
  69                 final double pressure;
  70
  71
  72                 temperature = stagTemp / (1 + (GAMMA-1)/2 * MathUtil.pow2(mach));
  73
  74                 // Speed of sound
  75                 double c = 331.3 * Math.sqrt(1 + (temperature - 273.15)/273.15);
  76
  77                 // Free-stream velocity
  78                 double v0 = c * mach;
  79
  80 //              kin.visc. = (3.7291e-06 + 4.9944e-08 * temperature) / density
  81                 pressure = reynolds * (3.7291e-06 + 4.9944e-08 * temperature) * R * temperature /
  82                                         (v0 * length);
  83
  84                 atmosphere.pressure = pressure;
  85                 atmosphere.temperature = temperature;
  86                 conditions.setAtmosphericConditions(atmosphere);
  87                 conditions.setVelocity(v0);
  88
  89                 if (Math.abs(conditions.getMach() - mach) > 0.001) {
  90                         System.err.println("Computed mach: "+conditions.getMach() + " requested "+mach);
  91 //                      System.exit(1);
  92                 }
  93
  94                 return v0;
  95         }
  96
  97
  98
  99         private void computeVsMach(PrintStream stream) {
 100
 101                 conditions.setAOA(0);
 102                 conditions.setTheta(45*Math.PI/180);
 103                 stream.println("% Mach, Caxial, CP, , CNa, Croll");
 104
 105                 for (double mach = MACH_MIN; mach <= MACH_MAX; mach += MACH_STEP) {
 106
 107                         computeVelocityAndAtmosphere(mach, REYNOLDS, STAG_TEMP);
 108 //                      conditions.setMach(mach);
 109
 110
 111                         AerodynamicForces forces = calculator.getAerodynamicForces(0, conditions, null);
 112
 113
 114                         double Re = conditions.getVelocity() *
 115                                         calculator.getConfiguration().getLength() /
 116                                         conditions.getAtmosphericConditions().getKinematicViscosity();
 117                         if (Math.abs(Re - REYNOLDS) > 1) {
 118                                 throw new RuntimeException("Re="+Re);
 119                         }
 120                         stream.printf("%f, %f, %f, %f, %f\n", mach, forces.Caxial, forces.cp.x, forces.CNa,
 121                                         forces.Croll);
 122                 }
 123
 124         }
 125
 126
 127
 128         private void computeVsAOA(PrintStream stream, double thetaDeg) {
 129
 130                 computeVelocityAndAtmosphere(AOA_MACH, REYNOLDS, STAG_TEMP);
 131                 conditions.setTheta(thetaDeg * Math.PI/180);
 132                 stream.println("% AOA, CP, CN, Cm   at theta = "+thetaDeg);
 133
 134                 for (double aoa = AOA_MIN; aoa <= AOA_MAX; aoa += AOA_STEP) {
 135
 136                         conditions.setAOA(aoa);
 137                         AerodynamicForces forces = calculator.getAerodynamicForces(0, conditions, null);
 138
 139
 140                         double Re = conditions.getVelocity() *
 141                                         calculator.getConfiguration().getLength() /
 142                                         conditions.getAtmosphericConditions().getKinematicViscosity();
 143                         if (Math.abs(Re - REYNOLDS) > 1) {
 144                                 throw new RuntimeException("Re="+Re);
 145                         }
 146                         stream.printf("%f, %f, %f, %f\n", aoa*180/Math.PI, forces.cp.x, forces.CN, forces.Cm);
 147                 }
 148
 149         }
 150
 151
 152
 153
 154         public static void main(String arg[]) throws Exception {
 155
 156                 if (arg.length != 2) {
 157                         System.err.println("Arguments:  <rocket file> <output prefix>");
 158                         System.exit(1);
 159                 }
 160
 161                 Analysis a = new Analysis(arg[0]);
 162                 final String prefix = arg[1];
 163
 164
 165                 String name;
 166                 double v0 = a.computeVelocityAndAtmosphere(0.6, 9.8e6, 322);
 167                 System.out.printf("Sanity test: mach = %.1f v=%.1f temp=%.1f pres=%.0f c=%.1f " +
 168                                 "ref.length=%.1fmm\n",
 169                                 a.conditions.getMach(), v0, a.atmosphere.temperature, a.atmosphere.pressure,
 170                                 a.atmosphere.getMachSpeed(), a.conditions.getRefLength()*1000);
 171                 System.out.println();
 172
 173
 174                 // CA, CP, Croll vs. Mach  at AOA=0
 175                 name = prefix + "-CA-CP-CNa-Croll-vs-Mach.csv";
 176                 System.out.println("Computing CA, CP, CNa, Croll vs. Mach to file "+name);
 177                 a.computeVsMach(new PrintStream(name));
 178
 179
 180                 // CN & Cm vs. AOA  at M=0.6
 181                 name = prefix + "-CP-CN-Cm-vs-AOA-0.csv";
 182                 System.out.println("Computing CP, CN, Cm vs. AOA at theta=0 to file "+name);
 183                 a.computeVsAOA(new PrintStream(name), 0);
 184
 185                 // CN & Cm vs. AOA  at M=0.6
 186                 name = prefix + "-CP-CN-Cm-vs-AOA-22.5.csv";
 187                 System.out.println("Computing CP, CN, Cm vs. AOA at theta=22.5 to file "+name);
 188                 a.computeVsAOA(new PrintStream(name), 0);
 189
 190                 // CN & Cm vs. AOA  at M=0.6
 191                 name = prefix + "-CP-CN-Cm-vs-AOA-45.csv";
 192                 System.out.println("Computing CP, CN, Cm vs. AOA at theta=45 to file "+name);
 193                 a.computeVsAOA(new PrintStream(name), 0);
 194
 195
 196                 System.out.println("Done.");
 197         }
 198
 199
 200
 201
 202
 203 }