001/* 002* NACA4ModUncambered -- An arbitrary uncambered NACA 4 digit modified airfoil. 003* 004* Copyright (C) 2000-2012, by Joseph A. Huwaldt 005* All rights reserved. 006* 007* This library is free software; you can redistribute it and/or 008* modify it under the terms of the GNU Lesser General Public 009* License as published by the Free Software Foundation; either 010* version 2.1 of the License, or (at your option) any later version. 011* 012* This library is distributed in the hope that it will be useful, 013* but WITHOUT ANY WARRANTY; without even the implied warranty of 014* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 015* Lesser General Public License for more details. 016* 017* You should have received a copy of the GNU Lesser General Public License 018* along with this program; if not, write to the Free Software 019* Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. 020* Or visit: http://www.gnu.org/licenses/lgpl.html 021**/ 022package jahuwaldt.aero.airfoils; 023 024import java.util.List; 025import java.awt.geom.Point2D; 026import java.text.DecimalFormat; 027import java.text.NumberFormat; 028 029 030/** 031* <p> This class represents an arbitrary uncambered NACA 032* modified 4 digit airfoil section such as a NACA 0012-34 033* airfoil. The 1st two digits indicate a symmetric airfoil, 034* the 2nd two, the thickness-chord ratio. The 1st digit 035* after the hyphen is the leading edge radius index and 036* the last digit after the hyphen is the position of 037* maximum thickness in tenths of chord. 038* </p> 039* 040* <p> Ported from FORTRAN "NACA4.FOR" to Java by: 041* Joseph A. Huwaldt, October 19, 2000 </p> 042* 043* <p> Original FORTRAN "NACA4" code had the following note: </p> 044* 045* <pre> 046* AUTHORS - Charles L.Ladson and Cuyler W. Brooks, NASA Langley 047* Liam Hardy, NASA Ames 048* Ralph Carmichael, Public Domain Aeronautical Software 049* Last FORTRAN version: 8Aug95 1.7 RLC 050* 051* NOTES - This program has been known by the names ANALIN, FOURDIGIT and NACA4. 052* REFERENCES- NASA Technical Memorandum TM X-3284 (November, 1975), 053* "Development of a Computer Program to Obtain Ordinates for 054* NACA 4-Digit, 4-Digit Modified, 5-Digit and 16-Digit Airfoils", 055* by Charles L. Ladson and Cuyler W. Brooks, Jr., 056* NASA Langley Research Center. 057* 058* NASA Technical Memorandum TM 4741 (December 1996), 059* "Computer Program to Obtain Ordinates for NACA Airfoils", 060* by Charles L. Ladson, Cuyler W. Brooks, Jr., and Acquilla S. Hill, 061* NASA Langley Research Center and 062* Darrell W. Sproles, Computer Sciences Corporation, Hampton, VA. 063* 064* "Theory of Wing Sections", by Ira Abbott and Albert Von Doenhoff. 065* </pre> 066* 067* <p> Modified by: Joseph A. Huwaldt </p> 068* 069* @author Joseph A. Huwaldt Date: October 19, 2000 070* @version September 15, 2012 071**/ 072public class NACA4ModUncambered extends NACA4Uncambered { 073 074 /** 075 * Table relating t/c ratio and LE index to LE radius. 076 **/ 077 private static final double[][] tableRle = { 078 { 0.000077, 0.000306, 0.000689, 0.001224, 0.001913, 0.002755, 0.003749, 0.004897 }, 079 { 0.000110, 0.000441, 0.000992, 0.001763, 0.002755, 0.003967, 0.005399, 0.007052 }, 080 { 0.000150, 0.000600, 0.001350, 0.002400, 0.003749, 0.005399, 0.007349, 0.009599 }, 081 { 0.000196, 0.000784, 0.001763, 0.003134, 0.004897, 0.007052, 0.009599, 0.012537 }, 082 { 0.000248, 0.000992, 0.002231, 0.003967, 0.006198, 0.008925, 0.012148, 0.015867 }, 083 { 0.000306, 0.001224, 0.002755, 0.004897, 0.007652, 0.011019, 0.014998, 0.019589 }, 084 { 0.000370, 0.001481, 0.003333, 0.005925, 0.009259, 0.013333, 0.018147, 0.023703 }, 085 { 0.000441, 0.001763, 0.003967, 0.007052, 0.011019, 0.015867, 0.021597, 0.028208 }, 086 { 0.000517, 0.002069, 0.004655, 0.008276, 0.012932, 0.018622, 0.025346, 0.033105 }, 087 { 0.000600, 0.002400, 0.005399, 0.009599, 0.014998, 0.021597, 0.029395, 0.038394 }, 088 { 0.000689, 0.002755, 0.006198, 0.011019, 0.017217, 0.024792, 0.033745, 0.044075 }, 089 { 0.000784, 0.003134, 0.007052, 0.012537, 0.019589, 0.028208, 0.038394, 0.050147 }, 090 { 0.000885, 0.003438, 0.007961, 0.014153, 0.022114, 0.031844, 0.043343, 0.056612 }, 091 { 0.000992, 0.003967, 0.008925, 0.015867, 0.024792, 0.035701, 0.048592, 0.063468 }, 092 { 0.001105, 0.004420, 0.009944, 0.017679, 0.027623, 0.039778, 0.054142, 0.070716 }, 093 { 0.001224, 0.004897, 0.011019, 0.019589, 0.030608, 0.044075, 0.059991, 0.078355 }, 094 { 0.001350, 0.005399, 0.012148, 0.021597, 0.033745, 0.048592, 0.066140, 0.086387 } }; 095 096 private static final double[] toc = {0.05, 0.06, 0.07, 0.08, 0.09, 097 0.10, 0.11, 0.12, 0.13, 0.14, 0.15, 0.16, 0.17, 0.18, 0.19, 098 0.20, 0.21 }; 099 100 /** 101 * The leading edge radius in percent chord length. 102 **/ 103 protected double RLE; 104 105 /** 106 * The leading edge radius index. 107 **/ 108 protected int leIndex = 0; 109 110 /** 111 * The coefficients of the profile equation. 112 **/ 113 private double a0, a1, a2, a3; 114 private double d0, d1, d2, d3; // d1 = trailing edge slope. 115 116 117 /** 118 * Create an uncambered NACA modified 4 digit airfoil with the 119 * specified parameters. For example: a NACA 0012-34 airfoil 120 * translates to: thickness = 0.12, Rle = 0.003967, xThickness = 0.40. 121 * 122 * @param thickness The thickness to chord ratio (e.g.: 0.09 ==> 9% t/c). 123 * @param Rle The radius of the leading edge of the airfoil in 124 * percent chord. 125 * @param xThickness The position of maximum thickness in percent chord 126 * (e.g.: 0.40 => 40% t/c). 127 * @param length The chord length. 128 **/ 129 public NACA4ModUncambered(double thickness, double Rle, double xThickness, double length) { 130 super(thickness, length); 131 132 xMaxT = xThickness; 133 RLE = Rle; 134 135 calcCoefficients(); 136 137 } 138 139 /** 140 * Create an uncambered NACA modified 4 digit airfoil with the 141 * specified parameters. For example: a NACA 0012-34 airfoil 142 * translates to: thickness = 0.12, LEindex = 3, xThickness = 0.40. 143 * 144 * @param thickness The thickness to chord ratio (e.g.: 0.09 ==> 9% t/c). 145 * @param LEindex The leading edge radius index for this airfoil. Should 146 * be a number from 1 to 8. 147 * @param xThickness The position of maximum thickness in percent chord 148 * (e.g.: 0.40 => 40% t/c). 149 * @param length The chord length. 150 **/ 151 public NACA4ModUncambered(double thickness, int LEindex, double xThickness, double length) { 152 super(thickness, length); 153 154 xMaxT = xThickness; 155 RLE = lookupRle(thickness, LEindex); 156 leIndex = LEindex; 157 158 calcCoefficients(); 159 160 } 161 162 163 /** 164 * Returns a string representation of this airfoil 165 * (the NACA designation of this airfoil). 166 **/ 167 @Override 168 public String toString() { 169 StringBuilder buffer = new StringBuilder("NACA 00"); 170 171 if (TOC < 0.1) 172 buffer.append("0"); 173 buffer.append((int)(TOC*100)); 174 175 buffer.append("-"); 176 177 if (leIndex == 0) 178 buffer.append("?"); 179 else 180 buffer.append(leIndex); 181 182 buffer.append((int)(xMaxT*10)); 183 184 if (chord != 1.0) { 185 buffer.append(" c="); 186 buffer.append(chord); 187 } 188 return buffer.toString(); 189 } 190 191 /** 192 * Method to calculate the ordinates of the uncambered airfoil 193 * forward of the maximum thickness point. This implementation 194 * computes a set of ordinates that are appropriate for a modified 195 * NACA 4 digit airfoil. 196 * 197 * @param x The x/c location currently being calculated. 198 * @param o The ordinate data structure. The following are set: y, yp, ypp. 199 **/ 200 @Override 201 protected void calcOrdinateForward(double x, Ordinate o) { 202 // Uncambered ordinate. 203 double sqrtx = Math.sqrt(x); 204 double x2 = x*x; 205 o.y = a0*sqrtx + a1*x + a2*x2 + a3*x2*x; 206 o.yp = 0.5*a0/sqrtx + a1 + 2*a2*x + 3*a3*x2; 207 o.ypp = -0.5*0.5*a0/Math.sqrt(x2*x) + 2*a2 + 3*2*a3*x; 208 } 209 210 /** 211 * Method to calculate the ordinates of the uncambered airfoil 212 * aft of the maximum thickness point. This implementation 213 * computes a set of ordinates that are appropriate for a modified 214 * NACA 4 digit airfoil. 215 * 216 * @param x The x/c location currently being calculated. 217 * @param o The ordinate data structure. The following are set: y, yp, ypp. 218 **/ 219 @Override 220 protected void calcOrdinateAft(double x, Ordinate o) { 221 double oneMx = 1 - x; 222 double oneMx2 = oneMx*oneMx; 223 o.y = d0 + d1*oneMx + d2*oneMx2 + d3*oneMx2*oneMx; 224 o.yp = -d1 - 2*d2*oneMx - 3*d3*oneMx2; 225 o.ypp = 2*d2 + 6*d3*oneMx; 226 } 227 228 //----------------------------------------------------------------------------------- 229 230 /** 231 * Calculate the ordinate equation coefficients required 232 * for a modified, uncambered, NACA 4 digit airfoil. 233 **/ 234 private void calcCoefficients() { 235 a0 = Math.sqrt(2*RLE)*0.2/TOC; 236 d0 = 0.002; 237 double xM2 = xMaxT*xMaxT; 238 d1 = 0.1*(2.24 - 5.42*xMaxT + 12.3*xM2)/(1 - 0.878*xMaxT); 239 double oneMxM = 1 - xMaxT; 240 double oneMxM2 = oneMxM*oneMxM; 241 d3 = (3*d1 - 0.588/oneMxM)/(3*oneMxM2); 242 d2 = -1.5*oneMxM*d3 - 0.5*d1/oneMxM; 243 a3 = 0.1/xM2/xMaxT + (2*d1*oneMxM - 0.588)/(2*xMaxT*oneMxM2) - 244 3*a0/(8*Math.pow(xMaxT, 2.5)); 245 a2 = -0.10/xM2 + 0.5*a0/Math.pow(xMaxT, 1.5) - 2*xMaxT*a3; 246 a1 = -0.5*a0/Math.sqrt(xMaxT) - 2*xMaxT*a2 - 3*xM2*a3; 247 } 248 249 /** 250 * Lookup the appropriate RLE value from the LE index and t/c. 251 **/ 252 private static double lookupRle(double thickness, int LEindex) { 253 --LEindex; 254 255 // Do a range check. 256 if (LEindex < 0) 257 LEindex = 0; 258 else if (LEindex > 7) 259 LEindex = 7; 260 261 // Find the indexes of the toc values that bound the thickness. 262 int length = toc.length; 263 int i=0; 264 for (; i < length; ++i) 265 if (toc[i] > thickness) 266 break; 267 int indexL = 0; 268 if (i == length) 269 indexL = i - 2; 270 else if (i != 0) 271 indexL = i - 1; 272 int indexH = indexL + 1; 273 274 // Do the interpolation into the table. 275 double value = lineInterp(toc[indexH], tableRle[indexH][LEindex], 276 toc[indexL], tableRle[indexL][LEindex], thickness); 277 278 return value; 279 } 280 281 /** 282 * Simple linear 1D interpolation between two points. 283 * 284 * @param x1,y1 Coordinates of the 1st point (the high point). 285 * @param x2,y2 Coordinates of the 2nd point (the low point). 286 * @param x The X coordinate of the point for which we want to 287 * interpolate to determine a Y coordinate. Will 288 * extrapolate if X is outside of the bounds of the 289 * point arguments. 290 * @return The interpolated Y value corresponding to the input X 291 * value is returned. 292 **/ 293 private static double lineInterp( double x1, double y1, double x2, double y2, double x ) { 294 return ((y2 - y1)/(x2 - x1)*(x - x1) + y1); 295 } 296 297 298 /** 299 * Simple method to test this class. 300 **/ 301 public static void main(String[] args) { 302 303 DecimalFormat nf = (DecimalFormat)NumberFormat.getInstance(); 304 nf.setMaximumFractionDigits(5); 305 nf.setMinimumFractionDigits(5); 306 307 System.out.println("Start NACA4ModUncambered..."); 308 309 System.out.println("Creating a NACA 0010-64 airfoil..."); 310 Airfoil af = new NACA4ModUncambered(0.10, 6, 0.4, 1); 311 312 System.out.println("Airfoil = " + af.toString()); 313 314 // Output the upper surface of the airfoil. 315 List<Point2D> upper = af.getUpper(); 316 List<Double> ypArr = af.getUpperYp(); 317 318 System.out.println(" X \t Y \t dy/dx"); 319 int length = upper.size(); 320 for (int i=0; i < length; ++i) { 321 Point2D o = upper.get(i); 322 System.out.println(" " + nf.format(o.getX()) + "\t" + nf.format(o.getY()) + 323 "\t" + nf.format(ypArr.get(i))); 324 } 325 326 System.out.println("# ordinates = " + length); 327 System.out.println("Done!"); 328 } 329 330} 331 332