Code-Eli  0.3.6
four_digit_test_suite.hpp
Go to the documentation of this file.
1 /*********************************************************************************
2 * Copyright (c) 2014 David D. Marshall <ddmarsha@calpoly.edu>
3 *
4 * All rights reserved. This program and the accompanying materials
5 * are made available under the terms of the Eclipse Public License v1.0
6 * which accompanies this distribution, and is available at
7 * http://www.eclipse.org/legal/epl-v10.html
8 *
9 * Contributors:
10 * David D. Marshall - initial code and implementation
11 ********************************************************************************/
12 
13 #ifndef four_digit_test_suite_hpp
14 #define four_digit_test_suite_hpp
15 
16 #include <cmath> // std::pow, std::exp
17 
18 #include <typeinfo> // typeid
19 #include <string> // std::string
20 #include <sstream> // std::stringstream
21 #include <iomanip> // std::setw
22 #include <limits> // std::numeric_limits
23 
24 #include "eli/constants/math.hpp"
25 #include "eli/mutil/fd/d1o2.hpp"
26 #include "eli/mutil/fd/d2o2.hpp"
28 
29 template<typename data__>
30 class four_digit_test_suite : public Test::Suite
31 {
32  private:
33  typedef data__ data_type;
37 
38  protected:
39  void AddTests(const float &)
40  {
41  // add the tests
50  }
51  void AddTests(const double &)
52  {
53  // add the tests
62  }
63  void AddTests(const long double &)
64  {
65  // add the tests
74  }
75 
76  public:
78  {
79  AddTests(data__());
80  }
82  {
83  }
84 
85  private:
87  {
88  airfoil_type af;
89  airfoil_thickness_coefficient_type a, a_ref;
90 
91  af.set_sharp_trailing_edge(false);
92  a_ref << static_cast<data_type>(0.2969), static_cast<data_type>(-0.1260), static_cast<data_type>(-0.3516), static_cast<data_type>(0.2843), static_cast<data_type>(-0.1015);
94  TEST_ASSERT((a-a_ref).norm()<1e-4);
95 
96  af.set_sharp_trailing_edge(true);
97  a_ref << static_cast<data_type>(0.2983), static_cast<data_type>(-0.1325), static_cast<data_type>(-0.3286), static_cast<data_type>(0.2442), static_cast<data_type>(-0.0815);
99  TEST_ASSERT((a-a_ref).norm()<1e-4);
100  }
101 
103  {
104  airfoil_type af;
105  data_type th, cam, cam_loc;
106  bool rtn;
107 
108  // set airfoil thickness
109  th=24;
110  rtn=af.set_thickness(th);
111  TEST_ASSERT(rtn);
112  TEST_ASSERT(af.get_thickness()==th);
113 
114  // set airfoil camber
115  cam=2;
116  cam_loc=3;
117  rtn=af.set_camber(cam, cam_loc);
118  TEST_ASSERT(rtn);
119  TEST_ASSERT(af.get_maximum_camber()==cam);
120  TEST_ASSERT(af.get_maximum_camber_location()==cam_loc);
121 
122  // test the name
123  std::string name, name_ref;
124 
125  name_ref="NACA "+std::to_string(static_cast<int>(std::round(cam)))
126  +std::to_string(static_cast<int>(std::round(cam_loc)))
127  +std::to_string(static_cast<int>(std::round(th)));
128  name=af.get_name();
129  TEST_ASSERT(name==name_ref);
130  }
131 
133  {
134  // non-sharp trailing edge
135  {
136  airfoil_type af;
137  std::vector<data_type> x(18), y(18);
138  airfoil_point_type xout;
139  const data_type tol(static_cast<data_type>(1e-5));
140 
141  // create 0021 airfoil to test against Abbott & von Doenhoff data NACA Report 824 p. 330
142  af.set_thickness(21);
143  x[0] = static_cast<data_type>(0); y[0] = static_cast<data_type>(0);
144  x[1] = static_cast<data_type>(1.25e-2); y[1] = static_cast<data_type>(3.315e-2);
145  x[2] = static_cast<data_type>(2.5e-2); y[2] = static_cast<data_type>(4.576e-2);
146  x[3] = static_cast<data_type>(5e-2); y[3] = static_cast<data_type>(6.221e-2);
147  x[4] = static_cast<data_type>(7.5e-2); y[4] = static_cast<data_type>(7.350e-2);
148  x[5] = static_cast<data_type>(10e-2); y[5] = static_cast<data_type>(8.195e-2);
149  x[6] = static_cast<data_type>(15e-2); y[6] = static_cast<data_type>(9.354e-2);
150  x[7] = static_cast<data_type>(20e-2); y[7] = static_cast<data_type>(10.040e-2);
151  x[8] = static_cast<data_type>(25e-2); y[8] = static_cast<data_type>(10.397e-2);
152  x[9] = static_cast<data_type>(30e-2); y[9] = static_cast<data_type>(10.504e-2);
153  x[10] = static_cast<data_type>(40e-2); y[10] = static_cast<data_type>(10.156e-2);
154  x[11] = static_cast<data_type>(50e-2); y[11] = static_cast<data_type>(9.265e-2);
155  x[12] = static_cast<data_type>(60e-2); y[12] = static_cast<data_type>(7.986e-2);
156  x[13] = static_cast<data_type>(70e-2); y[13] = static_cast<data_type>(6.412e-2);
157  x[14] = static_cast<data_type>(80e-2); y[14] = static_cast<data_type>(4.591e-2);
158  x[15] = static_cast<data_type>(90e-2); y[15] = static_cast<data_type>(2.534e-2);
159  x[16] = static_cast<data_type>(95e-2); y[16] = static_cast<data_type>(1.412e-2);
160  x[17] = static_cast<data_type>(1); y[17] = static_cast<data_type>(0.221e-2);
161 
162  // test the airfoil coordinates
163  for (size_t i=0; i<x.size(); ++i)
164  {
165  xout = af.f(x[i]);
166  TEST_ASSERT_DELTA(x[i], xout(0), tol);
167  TEST_ASSERT_DELTA(y[i], xout(1), tol);
168  }
169  }
170 
171  // sharp trailing edge
172  {
173  airfoil_type af;
174  std::vector<data_type> x(18), y(18);
175  airfoil_point_type xout;
176  const data_type tol(static_cast<data_type>(1e-5));
177 
178  // create 0021 airfoil to test against Abbott & von Doenhoff data NACA Report 824 p. 330
179  af.set_thickness(21);
180  af.set_sharp_trailing_edge(true);
181  x[0] = static_cast<data_type>(0); y[0] = static_cast<data_type>(0);
182  x[1] = static_cast<data_type>(1.25e-2); y[1] = static_cast<data_type>(3.323e-2);
183  x[2] = static_cast<data_type>(2.5e-2); y[2] = static_cast<data_type>(4.584e-2);
184  x[3] = static_cast<data_type>(5e-2); y[3] = static_cast<data_type>(6.226e-2);
185  x[4] = static_cast<data_type>(7.5e-2); y[4] = static_cast<data_type>(7.352e-2);
186  x[5] = static_cast<data_type>(10e-2); y[5] = static_cast<data_type>(8.195e-2);
187  x[6] = static_cast<data_type>(15e-2); y[6] = static_cast<data_type>(9.352e-2);
188  x[7] = static_cast<data_type>(20e-2); y[7] = static_cast<data_type>(10.039e-2);
189  x[8] = static_cast<data_type>(25e-2); y[8] = static_cast<data_type>(10.397e-2);
190  x[9] = static_cast<data_type>(30e-2); y[9] = static_cast<data_type>(10.503e-2);
191  x[10] = static_cast<data_type>(40e-2); y[10] = static_cast<data_type>(10.150e-2);
192  x[11] = static_cast<data_type>(50e-2); y[11] = static_cast<data_type>(9.241e-2);
193  x[12] = static_cast<data_type>(60e-2); y[12] = static_cast<data_type>(7.929e-2);
194  x[13] = static_cast<data_type>(70e-2); y[13] = static_cast<data_type>(6.308e-2);
195  x[14] = static_cast<data_type>(80e-2); y[14] = static_cast<data_type>(4.434e-2);
196  x[15] = static_cast<data_type>(90e-2); y[15] = static_cast<data_type>(2.333e-2);
197  x[16] = static_cast<data_type>(95e-2); y[16] = static_cast<data_type>(1.196e-2);
198  x[17] = static_cast<data_type>(1); y[17] = static_cast<data_type>(0);
199 
200  // test the airfoil coordinates
201  for (size_t i=0; i<x.size(); ++i)
202  {
203  xout = af.f(x[i]);
204  TEST_ASSERT_DELTA(x[i], xout(0), tol);
205  TEST_ASSERT_DELTA(y[i], xout(1), tol);
206 // std::cout << "i=" << i << "\tycal=" << xout(1) << "\tyref=" << y[i] << "\tdy=" << xout(1)-y[i] << std::endl;
207  }
208  }
209  }
210 
212  {
213  // test the known derivatives of thick trailing edge geometry
214  {
215  airfoil_type af;
216  airfoil_point_type x, xp, x_ref, xp_ref;
217  data_type xi;
218 
219  af.set_thickness(20);
220 
221  // lower surface trailing edge
222  xi=static_cast<data_type>(-1);
223  x_ref << -xi, static_cast<data_type>(-0.002);
224  xp_ref << static_cast<data_type>(-1), static_cast<data_type>(-0.234);
225  x=af.f(xi);
226  xp=af.fp(xi);
227  TEST_ASSERT((x-x_ref).norm()<2e-4);
228  TEST_ASSERT((xp-xp_ref).norm()<2e-4);
229 
230  // lower surface max thickness
231  xi=static_cast<data_type>(-0.3);
232  x_ref << -xi, static_cast<data_type>(-0.1);
233  xp_ref << static_cast<data_type>(-1), static_cast<data_type>(0);
234  x=af.f(xi);
235  xp=af.fp(xi);
236  TEST_ASSERT((x-x_ref).norm()<3e-5);
237  TEST_ASSERT((xp-xp_ref).norm()<2e-4);
238 
239  // leading edge
240  xi=static_cast<data_type>(0);
241  x_ref << xi, static_cast<data_type>(0);
242  xp_ref << static_cast<data_type>(0), static_cast<data_type>(1)/std::numeric_limits<data_type>::epsilon();
243  x=af.f(xi);
244  xp=af.fp(xi);
245  TEST_ASSERT(x==x_ref);
246  TEST_ASSERT(xp==xp_ref);
247 
248  // upper surface max thickness
249  xi=static_cast<data_type>(0.3);
250  x_ref << xi, static_cast<data_type>(0.1);
251  xp_ref << static_cast<data_type>(1), static_cast<data_type>(0);
252  x=af.f(xi);
253  xp=af.fp(xi);
254  TEST_ASSERT((x-x_ref).norm()<3e-5);
255  TEST_ASSERT((xp-xp_ref).norm()<2e-4);
256 
257  // upper surface trailing edge
258  xi=static_cast<data_type>(1);
259  x_ref << xi, static_cast<data_type>(0.002);
260  xp_ref << static_cast<data_type>(1), static_cast<data_type>(-0.234);
261  x=af.f(xi);
262  xp=af.fp(xi);
263  TEST_ASSERT((x-x_ref).norm()<2e-4);
264  TEST_ASSERT((xp-xp_ref).norm()<2e-4);
265 // std::cout << "xi=" << xi;
266 // std::cout << "\tx=" << x << "\tx_ref=" << x_ref;
267 // std::cout << "\tdx=" << x-x_ref << "\tdx.norm()=" << (x-x_ref).norm();
268 // std::cout << "\txp=" << xp << "\txp_ref=" << xp_ref;
269 // std::cout << "\tdxp=" << xp-xp_ref << "\tdxp.norm()=" << (xp-xp_ref).norm();
270 // std::cout << std::endl;
271  }
272 
273  // use finite differences to calculate the derivatives
274  {
275  airfoil_type af;
276  airfoil_point_type xr, xp, xpp, xp_ref, xpp_ref;
277  data_type xi, x[3], y[3];
278  const data_type dxi(static_cast<data_type>(1e2)*std::sqrt(std::numeric_limits<data_type>::epsilon()));
281 
282  af.set_thickness(12);
283 
284  // lower surface
285  xi=static_cast<data_type>(-0.8);
286  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
287  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
288  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
289  d1_calc.evaluate(xp_ref(0), x, dxi);
290  d1_calc.evaluate(xp_ref(1), y, dxi);
291  d2_calc.evaluate(xpp_ref(0), x, dxi);
292  d2_calc.evaluate(xpp_ref(1), y, dxi);
293  xp=af.fp(xi);
294  xpp=af.fpp(xi);
295  if (typeid(data_type)==typeid(float))
296  {
297  TEST_ASSERT((xp-xp_ref).norm()<1e-5);
298  TEST_ASSERT((xpp-xpp_ref).norm()<3e-3);
299  }
300  else
301  {
302  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
303  TEST_ASSERT((xpp-xpp_ref).norm()<5e-5);
304  }
305 
306  xi=static_cast<data_type>(-0.25);
307  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
308  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
309  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
310  d1_calc.evaluate(xp_ref(0), x, dxi);
311  d1_calc.evaluate(xp_ref(1), y, dxi);
312  d2_calc.evaluate(xpp_ref(0), x, dxi);
313  d2_calc.evaluate(xpp_ref(1), y, dxi);
314  xp=af.fp(xi);
315  xpp=af.fpp(xi);
316  if (typeid(data_type)==typeid(float))
317  {
318  TEST_ASSERT((xp-xp_ref).norm()<6e-4);
319  TEST_ASSERT((xpp-xpp_ref).norm()<3e-3);
320  }
321  else
322  {
323  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
324  TEST_ASSERT((xpp-xpp_ref).norm()<8e-6);
325  }
326 
327  // upper surface
328  xi=static_cast<data_type>(0.15);
329  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
330  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
331  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
332  d1_calc.evaluate(xp_ref(0), x, dxi);
333  d1_calc.evaluate(xp_ref(1), y, dxi);
334  d2_calc.evaluate(xpp_ref(0), x, dxi);
335  d2_calc.evaluate(xpp_ref(1), y, dxi);
336  xp=af.fp(xi);
337  xpp=af.fpp(xi);
338  if (typeid(data_type)==typeid(float))
339  {
340  TEST_ASSERT((xp-xp_ref).norm()<2e-3);
341  TEST_ASSERT((xpp-xpp_ref).norm()<2e-2);
342  }
343  else
344  {
345  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
346  TEST_ASSERT((xpp-xpp_ref).norm()<5e-6);
347  }
348 
349  xi=static_cast<data_type>(0.7);
350  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
351  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
352  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
353  d1_calc.evaluate(xp_ref(0), x, dxi);
354  d1_calc.evaluate(xp_ref(1), y, dxi);
355  d2_calc.evaluate(xpp_ref(0), x, dxi);
356  d2_calc.evaluate(xpp_ref(1), y, dxi);
357  xp=af.fp(xi);
358  xpp=af.fpp(xi);
359  if (typeid(data_type)==typeid(float))
360  {
361  TEST_ASSERT((xp-xp_ref).norm()<5e-5);
362  TEST_ASSERT((xpp-xpp_ref).norm()<5e-4);
363  }
364  else
365  {
366  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
367  TEST_ASSERT((xpp-xpp_ref).norm()<8e-5);
368  }
369 // std::cout << "xi=" << xi;
370 // std::cout << "\tdx=" << (x[2]-x[0])/(2*dxi) << "\tdy=" << (y[2]-y[0])/(2*dxi);
371 // std::cout << std::endl << " ";
372 // std::cout << "ddx=" << (x[0]-2*x[1]+x[2])/(dxi*dxi) << "\tddy=" << (y[0]-2*y[1]+y[2])/(dxi*dxi);
373 // std::cout << std::endl << " ";
374 // std::cout << "xp=" << xp;
375 // std::cout << std::endl << " ";
376 // std::cout << "xp_ref=" << xp_ref;
377 // std::cout << std::endl << " ";
378 // std::cout << "dxp=" << xp-xp_ref << "dxp.norm()=" << (xp-xp_ref).norm();
379 // std::cout << std::endl << " ";
380 // std::cout << "xpp= " << xpp;
381 // std::cout << std::endl << " ";
382 // std::cout << "xpp_ref=" << xpp_ref;
383 // std::cout << std::endl << " ";
384 // std::cout << "dxpp=" << xpp-xpp_ref << "\tdxpp.norm()=" << (xpp-xpp_ref).norm();
385 // std::cout << std::endl;
386  }
387  }
388 
389  void camber_test()
390  {
391  airfoil_type af;
392  std::vector<data_type> x(18), y(18);
393  airfoil_point_type xout;
394  const data_type tol(static_cast<data_type>(1e-5));
395 
396  // create 2300 airfoil to test against externally calculated data
397  af.set_camber(2, 3);
398  af.set_thickness(00);
399  x[ 0]=static_cast<data_type>(0); y[ 0]=static_cast<data_type>(0);
400  x[ 1]=static_cast<data_type>(1.25e-2); y[ 1]=static_cast<data_type>(0.1632e-2);
401  x[ 2]=static_cast<data_type>(2.50e-2); y[ 2]=static_cast<data_type>(0.3194e-2);
402  x[ 3]=static_cast<data_type>(5.00e-2); y[ 3]=static_cast<data_type>(0.6111e-2);
403  x[ 4]=static_cast<data_type>(7.50e-2); y[ 4]=static_cast<data_type>(0.8750e-2);
404  x[ 5]=static_cast<data_type>(10.00e-2); y[ 5]=static_cast<data_type>(1.1111e-2);
405  x[ 6]=static_cast<data_type>(15.00e-2); y[ 6]=static_cast<data_type>(1.5000e-2);
406  x[ 7]=static_cast<data_type>(20.00e-2); y[ 7]=static_cast<data_type>(1.7778e-2);
407  x[ 8]=static_cast<data_type>(25.00e-2); y[ 8]=static_cast<data_type>(1.9444e-2);
408  x[ 9]=static_cast<data_type>(30.00e-2); y[ 9]=static_cast<data_type>(2.0000e-2);
409  x[10]=static_cast<data_type>(40.00e-2); y[10]=static_cast<data_type>(1.9592e-2);
410  x[11]=static_cast<data_type>(50.00e-2); y[11]=static_cast<data_type>(1.8367e-2);
411  x[12]=static_cast<data_type>(60.00e-2); y[12]=static_cast<data_type>(1.6327e-2);
412  x[13]=static_cast<data_type>(70.00e-2); y[13]=static_cast<data_type>(1.3469e-2);
413  x[14]=static_cast<data_type>(80.00e-2); y[14]=static_cast<data_type>(0.9796e-2);
414  x[15]=static_cast<data_type>(90.00e-2); y[15]=static_cast<data_type>(0.5306e-2);
415  x[16]=static_cast<data_type>(95.00e-2); y[16]=static_cast<data_type>(0.2755e-2);
416  x[17]=static_cast<data_type>(1); y[17]=static_cast<data_type>(0);
417 
418  // test the airfoil coordinates
419  for (size_t i=0; i<x.size(); ++i)
420  {
421  xout = af.f(x[i]);
422  TEST_ASSERT_DELTA(x[i], xout(0), tol);
423  TEST_ASSERT_DELTA(y[i], xout(1), tol);
424 // std::cout << "i=" << i << "\tycal=" << xout(1) << "\tyref=" << y[i] << "\tdy=" << xout(1)-y[i] << std::endl;
425  }
426  }
427 
429  {
430  // test the known derivatives
431  {
432  airfoil_type af;
433  airfoil_point_type x, xp, x_ref, xp_ref;
434  data_type xi;
435 
436  // configure camber line
437  af.set_camber(3, 2);
438  af.set_thickness(00);
439 
440  // lower surface max camber
441  xi=static_cast<data_type>(-0.2);
442  x_ref << -xi, static_cast<data_type>(0.03);
443  xp_ref << static_cast<data_type>(-1), static_cast<data_type>(0);
444  x=af.f(xi);
445  xp=af.fp(xi);
446  TEST_ASSERT((x-x_ref).norm()<1e-8);
447  TEST_ASSERT((xp-xp_ref).norm()<1e-8);
448 
449  // upper surface max thickness
450  xi=static_cast<data_type>(0.2);
451  x_ref << xi, static_cast<data_type>(0.03);
452  xp_ref << static_cast<data_type>(1), static_cast<data_type>(0);
453  x=af.f(xi);
454  xp=af.fp(xi);
455  TEST_ASSERT((x-x_ref).norm()<1e-8);
456  TEST_ASSERT((xp-xp_ref).norm()<1e-8);
457 // std::cout << "xi=" << xi;
458 // std::cout << "\tx=" << x << "\tx_ref=" << x_ref;
459 // std::cout << "\tdx=" << x-x_ref << "\tdx.norm()=" << (x-x_ref).norm();
460 // std::cout << "\txp=" << xp << "\txp_ref=" << xp_ref;
461 // std::cout << "\tdxp=" << xp-xp_ref << "\tdxp.norm()=" << (xp-xp_ref).norm();
462 // std::cout << std::endl;
463  }
464 
465  // use finite differences to calculate the derivatives
466  {
467  airfoil_type af;
468  airfoil_point_type xr, xp, xpp, xp_ref, xpp_ref;
469  data_type xi, x[3], y[3];
470  const data_type dxi(static_cast<data_type>(1e2)*std::sqrt(std::numeric_limits<data_type>::epsilon()));
473 
474  // configure camber line
475  af.set_camber(3, 2);
476  af.set_thickness(00);
477 
478  // lower surface
479  xi=static_cast<data_type>(-0.8);
480  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
481  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
482  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
483  d1_calc.evaluate(xp_ref(0), x, dxi);
484  d1_calc.evaluate(xp_ref(1), y, dxi);
485  d2_calc.evaluate(xpp_ref(0), x, dxi);
486  d2_calc.evaluate(xpp_ref(1), y, dxi);
487  xp=af.fp(xi);
488  xpp=af.fpp(xi);
489  if (typeid(data_type)==typeid(float))
490  {
491  TEST_ASSERT((xp-xp_ref).norm()<1e-6);
492  TEST_ASSERT((xpp-xpp_ref).norm()<5e-5);
493  }
494  else
495  {
496  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
497  TEST_ASSERT((xpp-xpp_ref).norm()<1e-10);
498  }
499 
500  xi=static_cast<data_type>(-0.25);
501  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
502  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
503  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
504  d1_calc.evaluate(xp_ref(0), x, dxi);
505  d1_calc.evaluate(xp_ref(1), y, dxi);
506  d2_calc.evaluate(xpp_ref(0), x, dxi);
507  d2_calc.evaluate(xpp_ref(1), y, dxi);
508  xp=af.fp(xi);
509  xpp=af.fpp(xi);
510  if (typeid(data_type)==typeid(float))
511  {
512  TEST_ASSERT((xp-xp_ref).norm()<1e-6);
513  TEST_ASSERT((xpp-xpp_ref).norm()<1e-5);
514  }
515  else
516  {
517  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
518  TEST_ASSERT((xpp-xpp_ref).norm()<5e-6);
519  }
520 
521  // upper surface
522  xi=static_cast<data_type>(0.15);
523  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
524  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
525  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
526  d1_calc.evaluate(xp_ref(0), x, dxi);
527  d1_calc.evaluate(xp_ref(1), y, dxi);
528  d2_calc.evaluate(xpp_ref(0), x, dxi);
529  d2_calc.evaluate(xpp_ref(1), y, dxi);
530  xp=af.fp(xi);
531  xpp=af.fpp(xi);
532  if (typeid(data_type)==typeid(float))
533  {
534  TEST_ASSERT((xp-xp_ref).norm()<1e-6);
535  TEST_ASSERT((xpp-xpp_ref).norm()<5e-6);
536  }
537  else
538  {
539  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
540  TEST_ASSERT((xpp-xpp_ref).norm()<2e-6);
541  }
542 
543  xi=static_cast<data_type>(0.7);
544  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
545  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
546  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
547  d1_calc.evaluate(xp_ref(0), x, dxi);
548  d1_calc.evaluate(xp_ref(1), y, dxi);
549  d2_calc.evaluate(xpp_ref(0), x, dxi);
550  d2_calc.evaluate(xpp_ref(1), y, dxi);
551  xp=af.fp(xi);
552  xpp=af.fpp(xi);
553  if (typeid(data_type)==typeid(float))
554  {
555  TEST_ASSERT((xp-xp_ref).norm()<2e-3);
556  TEST_ASSERT((xpp-xpp_ref).norm()<2e-2);
557  }
558  else
559  {
560  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
561  TEST_ASSERT((xpp-xpp_ref).norm()<7e-5);
562  }
563 // std::cout << "xi=" << xi;
564 // std::cout << std::endl << " ";
565 // std::cout << "xp=" << xp;
566 // std::cout << std::endl << " ";
567 // std::cout << "xp_ref=" << xp_ref;
568 // std::cout << std::endl << " ";
569 // std::cout << "dxp=" << xp-xp_ref << "dxp.norm()=" << (xp-xp_ref).norm();
570 // std::cout << std::endl << " ";
571 // std::cout << "xpp= " << xpp;
572 // std::cout << std::endl << " ";
573 // std::cout << "xpp_ref=" << xpp_ref;
574 // std::cout << std::endl << " ";
575 // std::cout << "dxpp=" << xpp-xpp_ref << "\tdxpp.norm()=" << (xpp-xpp_ref).norm();
576 // std::cout << std::endl;
577  }
578  }
579 
581  {
582  airfoil_type af;
583  std::vector<data_type> xi(36), x(36), y(36);
584  airfoil_point_type xout;
585  const data_type tol(static_cast<data_type>(3e-3));
586 
587  // create 2412 airfoil to test against Abbott & von Doenhoff data NACA Report 824 p. 358
588  af.set_camber(2, 4);
589  af.set_thickness(12);
590  xi[ 0] = static_cast<data_type>( 0.0e-2); x[ 0] = static_cast<data_type>( 0.0e-2); y[ 0] = static_cast<data_type>( 0.0e-2);
591  xi[ 1] = static_cast<data_type>( -1.078935e-2); x[ 1] = static_cast<data_type>( 1.25e-2); y[ 1] = static_cast<data_type>(-1.65e-2);
592  xi[ 2] = static_cast<data_type>( -2.265269e-2); x[ 2] = static_cast<data_type>( 2.50e-2); y[ 2] = static_cast<data_type>(-2.27e-2);
593  xi[ 3] = static_cast<data_type>( -4.695760e-2); x[ 3] = static_cast<data_type>( 5.00e-2); y[ 3] = static_cast<data_type>(-3.01e-2);
594  xi[ 4] = static_cast<data_type>( -7.162585e-2); x[ 4] = static_cast<data_type>( 7.50e-2); y[ 4] = static_cast<data_type>(-3.46e-2);
595  xi[ 5] = static_cast<data_type>( -9.650263e-2); x[ 5] = static_cast<data_type>( 10.00e-2); y[ 5] = static_cast<data_type>(-3.75e-2);
596  xi[ 6] = static_cast<data_type>(-14.664316e-2); x[ 6] = static_cast<data_type>( 15.00e-2); y[ 6] = static_cast<data_type>(-4.10e-2);
597  xi[ 7] = static_cast<data_type>(-19.710203e-2); x[ 7] = static_cast<data_type>( 20.00e-2); y[ 7] = static_cast<data_type>(-4.23e-2);
598  xi[ 8] = static_cast<data_type>(-24.774236e-2); x[ 8] = static_cast<data_type>( 25.00e-2); y[ 8] = static_cast<data_type>(-4.22e-2);
599  xi[ 9] = static_cast<data_type>(-29.847722e-2); x[ 9] = static_cast<data_type>( 30.00e-2); y[ 9] = static_cast<data_type>(-4.12e-2);
600  xi[10] = static_cast<data_type>(-40.000000e-2); x[10] = static_cast<data_type>( 40.00e-2); y[10] = static_cast<data_type>(-3.80e-2);
601  xi[11] = static_cast<data_type>(-50.059125e-2); x[11] = static_cast<data_type>( 50.00e-2); y[11] = static_cast<data_type>(-3.34e-2);
602  xi[12] = static_cast<data_type>(-60.101712e-2); x[12] = static_cast<data_type>( 60.00e-2); y[12] = static_cast<data_type>(-2.76e-2);
603  xi[13] = static_cast<data_type>(-70.122161e-2); x[13] = static_cast<data_type>( 70.00e-2); y[13] = static_cast<data_type>(-2.14e-2);
604  xi[14] = static_cast<data_type>(-80.116232e-2); x[14] = static_cast<data_type>( 80.00e-2); y[14] = static_cast<data_type>(-1.50e-2);
605  xi[15] = static_cast<data_type>(-90.079880e-2); x[15] = static_cast<data_type>( 90.00e-2); y[15] = static_cast<data_type>(-0.82e-2);
606  xi[16] = static_cast<data_type>(-95.048848e-2); x[16] = static_cast<data_type>( 95.00e-2); y[16] = static_cast<data_type>(-0.48e-2);
607  xi[17] = static_cast<data_type>( -1.0e+0); x[17] = static_cast<data_type>( 99.99e-2); y[17] = static_cast<data_type>(-0.13e-2);
608  xi[18] = static_cast<data_type>( 0.0e-2); x[18] = static_cast<data_type>( 0.0e-2); y[18] = static_cast<data_type>( 0.0e-2);
609  xi[19] = static_cast<data_type>( 1.444524e-2); x[19] = static_cast<data_type>( 1.25e-2); y[19] = static_cast<data_type>( 2.15e-2);
610  xi[20] = static_cast<data_type>( 2.753309e-2); x[20] = static_cast<data_type>( 2.50e-2); y[20] = static_cast<data_type>( 2.99e-2);
611  xi[21] = static_cast<data_type>( 5.315146e-2); x[21] = static_cast<data_type>( 5.00e-2); y[21] = static_cast<data_type>( 4.13e-2);
612  xi[22] = static_cast<data_type>( 7.842503e-2); x[22] = static_cast<data_type>( 7.50e-2); y[22] = static_cast<data_type>( 4.96e-2);
613  xi[23] = static_cast<data_type>( 10.350449e-2); x[23] = static_cast<data_type>( 10.00e-2); y[23] = static_cast<data_type>( 5.63e-2);
614  xi[24] = static_cast<data_type>( 15.331062e-2); x[24] = static_cast<data_type>( 15.00e-2); y[24] = static_cast<data_type>( 6.61e-2);
615  xi[25] = static_cast<data_type>( 20.283261e-2); x[25] = static_cast<data_type>( 20.00e-2); y[25] = static_cast<data_type>( 7.26e-2);
616  xi[26] = static_cast<data_type>( 25.219585e-2); x[26] = static_cast<data_type>( 25.00e-2); y[26] = static_cast<data_type>( 7.67e-2);
617  xi[27] = static_cast<data_type>( 30.147780e-2); x[27] = static_cast<data_type>( 30.00e-2); y[27] = static_cast<data_type>( 7.88e-2);
618  xi[28] = static_cast<data_type>( 40.000000e-2); x[28] = static_cast<data_type>( 40.00e-2); y[28] = static_cast<data_type>( 7.80e-2);
619  xi[29] = static_cast<data_type>( 49.941485e-2); x[29] = static_cast<data_type>( 50.00e-2); y[29] = static_cast<data_type>( 7.24e-2);
620  xi[30] = static_cast<data_type>( 59.898946e-2); x[30] = static_cast<data_type>( 60.00e-2); y[30] = static_cast<data_type>( 6.36e-2);
621  xi[31] = static_cast<data_type>( 69.878040e-2); x[31] = static_cast<data_type>( 70.00e-2); y[31] = static_cast<data_type>( 5.18e-2);
622  xi[32] = static_cast<data_type>( 79.883299e-2); x[32] = static_cast<data_type>( 80.00e-2); y[32] = static_cast<data_type>( 3.75e-2);
623  xi[33] = static_cast<data_type>( 89.919268e-2); x[33] = static_cast<data_type>( 90.00e-2); y[33] = static_cast<data_type>( 2.08e-2);
624  xi[34] = static_cast<data_type>( 94.950448e-2); x[34] = static_cast<data_type>( 95.00e-2); y[34] = static_cast<data_type>( 1.14e-2);
625  xi[35] = static_cast<data_type>( 1.0e+0); x[35] = static_cast<data_type>(100.01e-2); y[35] = static_cast<data_type>( 0.13e-2);
626 
627  // test the airfoil coordinates
628  for (size_t i=0; i<x.size(); ++i)
629  {
630  xout = af.f(xi[i]);
631  TEST_ASSERT_DELTA(x[i], xout(0), tol);
632  TEST_ASSERT_DELTA(y[i], xout(1), tol);
633 // std::cout << "i=" << i << std::endl;
634 // std::cout << " " << "\txcal=" << xout(0)
635 // << "\txref=" << x[i]
636 // << "\tdx=" << xout(0)-x[i] << std::endl;
637 // std::cout << " " << "\tycal=" << xout(1)
638 // << "\tyref=" << y[i]
639 // << "\tdy=" << xout(1)-y[i] << std::endl;
640  }
641  }
642 
644  {
645  // use finite differences to calculate the derivatives
646  airfoil_type af;
647  airfoil_point_type xr, xp, xpp, xp_ref, xpp_ref;
648  data_type xi, x[3], y[3];
649  const data_type dxi(static_cast<data_type>(1e2)*std::sqrt(std::numeric_limits<data_type>::epsilon()));
652 
653  // configure camber line
654  af.set_camber(3, 2);
655  af.set_thickness(24);
656 
657  // lower surface
658  xi=static_cast<data_type>(-0.8);
659  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
660  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
661  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
662  d1_calc.evaluate(xp_ref(0), x, dxi);
663  d1_calc.evaluate(xp_ref(1), y, dxi);
664  d2_calc.evaluate(xpp_ref(0), x, dxi);
665  d2_calc.evaluate(xpp_ref(1), y, dxi);
666  xp=af.fp(xi);
667  xpp=af.fpp(xi);
668  if (typeid(data_type)==typeid(float))
669  {
670  TEST_ASSERT((xp-xp_ref).norm()<2e-5);
671  TEST_ASSERT((xpp-xpp_ref).norm()<5e-4);
672  }else
673  {
674  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
675  TEST_ASSERT((xpp-xpp_ref).norm()<2e-4);
676  }
677 
678  xi=static_cast<data_type>(-0.25);
679  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
680  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
681  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
682  d1_calc.evaluate(xp_ref(0), x, dxi);
683  d1_calc.evaluate(xp_ref(1), y, dxi);
684  d2_calc.evaluate(xpp_ref(0), x, dxi);
685  d2_calc.evaluate(xpp_ref(1), y, dxi);
686  xp=af.fp(xi);
687  xpp=af.fpp(xi);
688  if (typeid(data_type)==typeid(float))
689  {
690  TEST_ASSERT((xp-xp_ref).norm()<5e-3);
691  TEST_ASSERT((xpp-xpp_ref).norm()<5e-3);
692  }
693  else
694  {
695  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
696  TEST_ASSERT((xpp-xpp_ref).norm()<5e-5);
697  }
698 
699  // upper surface
700  xi=static_cast<data_type>(0.15);
701  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
702  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
703  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
704  d1_calc.evaluate(xp_ref(0), x, dxi);
705  d1_calc.evaluate(xp_ref(1), y, dxi);
706  d2_calc.evaluate(xpp_ref(0), x, dxi);
707  d2_calc.evaluate(xpp_ref(1), y, dxi);
708  xp=af.fp(xi);
709  xpp=af.fpp(xi);
710  if (typeid(data_type)==typeid(float))
711  {
712  TEST_ASSERT((xp-xp_ref).norm()<5e-3);
713  TEST_ASSERT((xpp-xpp_ref).norm()<3e-2);
714  }
715  else
716  {
717  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
718  TEST_ASSERT((xpp-xpp_ref).norm()<6e-3);
719  }
720 
721  xi=static_cast<data_type>(0.7);
722  xr=af.f(xi-dxi); x[0]=xr(0); y[0]=xr(1);
723  xr=af.f(xi); x[1]=xr(0); y[1]=xr(1);
724  xr=af.f(xi+dxi); x[2]=xr(0); y[2]=xr(1);
725  d1_calc.evaluate(xp_ref(0), x, dxi);
726  d1_calc.evaluate(xp_ref(1), y, dxi);
727  d2_calc.evaluate(xpp_ref(0), x, dxi);
728  d2_calc.evaluate(xpp_ref(1), y, dxi);
729  xp=af.fp(xi);
730  xpp=af.fpp(xi);
731  if (typeid(data_type)==typeid(float))
732  {
733  TEST_ASSERT((xp-xp_ref).norm()<1e-4);
734  TEST_ASSERT((xpp-xpp_ref).norm()<5e-4);
735  }
736  else
737  {
738  TEST_ASSERT((xp-xp_ref).norm()<1e-10);
739  TEST_ASSERT((xpp-xpp_ref).norm()<1e-4);
740  }
741 // std::cout << "xi=" << xi;
742 // std::cout << "\tdx=" << (x[2]-x[0])/(2*dxi) << "\tdy=" << (y[2]-y[0])/(2*dxi);
743 // std::cout << std::endl << " ";
744 // std::cout << "ddx=" << (x[0]-2*x[1]+x[2])/(dxi*dxi) << "\tddy=" << (y[0]-2*y[1]+y[2])/(dxi*dxi);
745 // std::cout << std::endl << " ";
746 // std::cout << "xp=" << xp;
747 // std::cout << std::endl << " ";
748 // std::cout << "xp_ref=" << xp_ref;
749 // std::cout << std::endl << " ";
750 // std::cout << "dxp=" << xp-xp_ref << "\tdxp.norm()=" << (xp-xp_ref).norm();
751 // std::cout << std::endl << " ";
752 // std::cout << "xpp= " << xpp;
753 // std::cout << std::endl << " ";
754 // std::cout << "xpp_ref=" << xpp_ref;
755 // std::cout << std::endl << " ";
756 // std::cout << "dxpp=" << xpp-xpp_ref << "\tdxpp.norm()=" << (xpp-xpp_ref).norm();
757 // std::cout << std::endl;
758  }
759 
760 };
761 
762 #endif
763 
~four_digit_test_suite()
Definition: four_digit_test_suite.hpp:81
point_type f(const data_type &xi) const
Definition: four_digit.hpp:130
void airfoil_derivatives_test()
Definition: four_digit_test_suite.hpp:643
void camber_derivatives_test()
Definition: four_digit_test_suite.hpp:428
Definition: four_digit_test_suite.hpp:30
bool set_camber(const data_type &cam, const data_type &cam_loc)
Definition: four_digit.hpp:78
four_digit_test_suite()
Definition: four_digit_test_suite.hpp:77
eli::geom::curve::pseudo::four_digit< data_type > airfoil_type
Definition: four_digit_test_suite.hpp:34
data_type get_thickness() const
Definition: four_digit.hpp:76
int evaluate(data__ &d, __itphi itphi, const data__ &dx) const
Definition: d2o2.hpp:147
Eigen::Matrix< data_type, 5, 1 > coefficient_type
Definition: four_digit.hpp:37
std::string get_name() const
Definition: four_digit.hpp:121
coefficient_type get_thickness_coefficients() const
Definition: four_digit.hpp:54
point_type fpp(const data_type &xi) const
Definition: four_digit.hpp:148
data_type get_maximum_camber() const
Definition: four_digit.hpp:90
void AddTests(const float &)
Definition: four_digit_test_suite.hpp:39
Eigen::Matrix< data_type, 1, 2 > point_type
Definition: four_digit.hpp:36
Definition: d1o2.hpp:27
airfoil_type::coefficient_type airfoil_thickness_coefficient_type
Definition: four_digit_test_suite.hpp:35
Definition: d2o2.hpp:27
int evaluate(data__ &d, itphi__ itphi, const data__ &dx) const
Definition: d1o2.hpp:131
void thickness_test()
Definition: four_digit_test_suite.hpp:132
void airfoil_coefficients_test()
Definition: four_digit_test_suite.hpp:86
void create_airfoil_test()
Definition: four_digit_test_suite.hpp:102
void camber_test()
Definition: four_digit_test_suite.hpp:389
void AddTests(const double &)
Definition: four_digit_test_suite.hpp:51
data__ data_type
Definition: four_digit_test_suite.hpp:33
void set_sharp_trailing_edge(bool fl)
Definition: four_digit.hpp:59
void AddTests(const long double &)
Definition: four_digit_test_suite.hpp:63
void thickness_derivatives_test()
Definition: four_digit_test_suite.hpp:211
data_type get_maximum_camber_location() const
Definition: four_digit.hpp:91
point_type fp(const data_type &xi) const
Definition: four_digit.hpp:139
bool set_thickness(const data_type &t)
Definition: four_digit.hpp:66
void airfoil_test()
Definition: four_digit_test_suite.hpp:580
airfoil_type::point_type airfoil_point_type
Definition: four_digit_test_suite.hpp:36