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