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

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