1[/ 2Copyright (c) 2012 John Maddock 3Use, modification and distribution are subject to the 4Boost Software License, Version 1.0. (See accompanying file 5LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt) 6] 7 8[section:airy Airy Functions] 9 10[section:ai Airy Ai Function] 11 12[heading Synopsis] 13 14`` 15 #include <boost/math/special_functions/airy.hpp> 16`` 17 18 namespace boost { namespace math { 19 20 template <class T> 21 ``__sf_result`` airy_ai(T x); 22 23 template <class T, class Policy> 24 ``__sf_result`` airy_ai(T x, const Policy&); 25 26 }} // namespaces 27 28[heading Description] 29 30The function __airy_ai calculates the Airy function Ai which is the first solution 31to the differential equation: 32 33[equation airy] 34 35See Weisstein, Eric W. "Airy Functions." From MathWorld--A Wolfram Web Resource. 36[@http://mathworld.wolfram.com/AiryFunctions.html] 37 38and [@https://en.wikipedia.org/wiki/Airy_zeta_function Airy Zeta function]. 39 40[optional_policy] 41 42The following graph illustrates how this function changes as /x/ changes: for negative /x/ 43the function is cyclic, while for positive /x/ the value tends to zero: 44 45[graph airy_ai] 46 47[heading Accuracy] 48 49This function is implemented entirely in terms of the Bessel functions 50__cyl_bessel_j and __cyl_bessel_k - refer to those functions for detailed accuracy information. 51 52In general though, the relative error is low (less than 100 [epsilon]) for /x > 0/ while 53only the absolute error is low for /x < 0/ as the following error plot illustrates: 54 55[graph ai__double] 56 57[heading Testing] 58 59Since this function is implemented in terms of other special functions, there are only a few 60basic sanity checks, using test values from [@http://functions.wolfram.com/ Wolfram Airy Functions]. 61 62[heading Implementation] 63 64This function is implemented in terms of the Bessel functions using the relations: 65 66[equation airy_ai] 67 68[endsect] [/section:ai Airy Ai Function] 69 70[section:bi Airy Bi Function] 71 72[heading Synopsis] 73 74`` 75 #include <boost/math/special_functions/airy.hpp> 76`` 77 78 namespace boost { namespace math { 79 80 template <class T> 81 ``__sf_result`` airy_bi(T x); 82 83 template <class T, class Policy> 84 ``__sf_result`` airy_bi(T x, const Policy&); 85 86 }} // namespaces 87 88[heading Description] 89 90The function __airy_bi calculates the Airy function Bi which is the second solution to the differential equation: 91 92[equation airy] 93 94[optional_policy] 95 96The following graph illustrates how this function changes as /x/ changes: for negative /x/ 97the function is cyclic, while for positive /x/ the value tends to infinity: 98 99[graph airy_bi] 100 101[heading Accuracy] 102 103This function is implemented entirely in terms of the Bessel functions 104__cyl_bessel_i and __cyl_bessel_j - refer to those functions for detailed accuracy information. 105 106In general though, the relative error is low (less than 100 [epsilon]) for /x > 0/ while 107only the absolute error is low for /x < 0/ as the following error plot illustrate: 108 109[graph bi__double] 110 111[heading Testing] 112 113Since this function is implemented in terms of other special functions, there are only a few 114basic sanity checks, using test values from [@http://functions.wolfram.com functions.wolfram.com]. 115 116[heading Implementation] 117 118This function is implemented in terms of the Bessel functions using the relations: 119 120[equation airy_bi] 121 122[endsect] [/section:bi Airy Bi Function] 123 124[section:aip Airy Ai' Function] 125 126[heading Synopsis] 127 128`` 129 #include <boost/math/special_functions/airy.hpp> 130`` 131 132 namespace boost { namespace math { 133 134 template <class T> 135 ``__sf_result`` airy_ai_prime(T x); 136 137 template <class T, class Policy> 138 ``__sf_result`` airy_ai_prime(T x, const Policy&); 139 140 }} // namespaces 141 142[heading Description] 143 144The function __airy_ai_prime calculates the Airy function Ai' which is the derivative of the first solution to the differential equation: 145 146[equation airy] 147 148[optional_policy] 149 150The following graph illustrates how this function changes as /x/ changes: for negative /x/ 151the function is cyclic, while for positive /x/ the value tends to zero: 152 153[graph airy_aip] 154 155[heading Accuracy] 156 157This function is implemented entirely in terms of the Bessel functions 158__cyl_bessel_j and __cyl_bessel_k - refer to those functions for detailed accuracy information. 159 160In general though, the relative error is low (less than 100 [epsilon]) for /x > 0/ while 161only the absolute error is low for /x < 0/ as the following error plot illustrates: 162 163[graph ai_prime__double] 164 165[heading Testing] 166 167Since this function is implemented in terms of other special functions, there are only a few 168basic sanity checks, using test values from [@http://functions.wolfram.com functions.wolfram.com]. 169 170[heading Implementation] 171 172This function is implemented in terms of the Bessel functions using the relations: 173 174[equation airy_aip] 175 176[endsect] [/section:aip Airy Ai' Function] 177 178[section:bip Airy Bi' Function] 179 180[heading Synopsis] 181 182`` 183 #include <boost/math/special_functions/airy.hpp> 184`` 185 186 namespace boost { namespace math { 187 188 template <class T> 189 ``__sf_result`` airy_bi_prime(T x); 190 191 template <class T, class Policy> 192 ``__sf_result`` airy_bi_prime(T x, const Policy&); 193 194 }} // namespaces 195 196[heading Description] 197 198The function __airy_bi_prime calculates the Airy function Bi' which is the derivative of the second solution to the differential equation: 199 200[equation airy] 201 202[optional_policy] 203 204The following graph illustrates how this function changes as /x/ changes: for negative /x/ 205the function is cyclic, while for positive /x/ the value tends to infinity: 206 207[graph airy_bi] 208 209[heading Accuracy] 210 211This function is implemented entirely in terms of the Bessel functions 212__cyl_bessel_i and __cyl_bessel_j - refer to those functions for detailed accuracy information. 213 214In general though, the relative error is low (less than 100 [epsilon]) for /x > 0/ while 215only the absolute error is low for /x < 0/ as the following error plot illustrates: 216 217[graph bi_prime__double] 218 219[heading Testing] 220 221Since this function is implemented in terms of other special functions, there are only a few 222basic sanity checks, using test values from [@http://functions.wolfram.com functions.wolfram.com]. 223 224[heading Implementation] 225 226This function is implemented in terms of the Bessel functions using the relations: 227 228[equation airy_bip] 229 230[endsect] [/section:bip Airy Bi' Function] 231 232[section:airy_root Finding Zeros of Airy Functions] 233 234[h4 Synopsis] 235 236`#include <boost/math/special_functions/airy.hpp>` 237 238Functions for obtaining both a single zero or root of the Airy functions, 239and placing multiple zeros into a container like `std::vector` 240by providing an output iterator. 241 242The signature of the single value functions are: 243 244 template <class T> 245 T airy_ai_zero( 246 int m); // 1-based index of zero. 247 248 template <class T> 249 T airy_bi_zero( 250 int m); // 1-based index of zero. 251 252and for multiple zeros: 253 254 template <class T, class OutputIterator> 255 OutputIterator airy_ai_zero( 256 int start_index, // 1-based index of first zero. 257 unsigned number_of_zeros, // How many zeros to generate. 258 OutputIterator out_it); // Destination for zeros. 259 260 template <class T, class OutputIterator> 261 OutputIterator airy_bi_zero( 262 int start_index, // 1-based index of zero. 263 unsigned number_of_zeros, // How many zeros to generate 264 OutputIterator out_it); // Destination for zeros. 265 266There are also versions which allow control of the __policy_section for error handling and precision. 267 268 template <class T> 269 T airy_ai_zero( 270 int m, // 1-based index of zero. 271 const Policy&); // Policy to use. 272 273 template <class T> 274 T airy_bi_zero( 275 int m, // 1-based index of zero. 276 const Policy&); // Policy to use. 277 278 279 template <class T, class OutputIterator> 280 OutputIterator airy_ai_zero( 281 int start_index, // 1-based index of first zero. 282 unsigned number_of_zeros, // How many zeros to generate. 283 OutputIterator out_it, // Destination for zeros. 284 const Policy& pol); // Policy to use. 285 286 template <class T, class OutputIterator> 287 OutputIterator airy_bi_zero( 288 int start_index, // 1-based index of zero. 289 unsigned number_of_zeros, // How many zeros to generate. 290 OutputIterator out_it, // Destination for zeros. 291 const Policy& pol); // Policy to use. 292 293[h4 Description] 294 295The Airy Ai and Bi functions have an infinite 296number of zeros on the negative real axis. The real zeros on the negative real 297axis can be found by solving for the roots of 298 299[:['Ai(x[sub m]) = 0]] 300 301[:['Bi(y[sub m]) = 0]] 302 303Here, ['x[sub m]] represents the ['m[super th]] 304root of the Airy Ai function, 305and ['y[sub m]] represents the ['m[super th]] 306root of the Airy Bi function. 307 308The zeros or roots (values of `x` where the function crosses the horizontal `y = 0` axis) 309of the Airy Ai and Bi functions are computed by two functions, 310`airy_ai_zero` and `airy_bi_zero`. 311 312In each case the index or rank of the zero 313returned is 1-based, which is to say: 314 315 airy_ai_zero(1); 316 317returns the first zero of Ai. 318 319Passing an `start_index <= 0` results in a __domain_error being raised. 320 321The first few zeros returned by these functions have approximate values as follows: 322 323[table 324[[m][Ai][Bi]] 325[[1][-2.33811...][-1.17371...]] 326[[2][-4.08795...][-3.27109...]] 327[[3][-5.52056...][-4.83074...]] 328[[4][-6.78671...][-6.16985...]] 329[[5][-7.94413...][-7.37676...]] 330[[6][-9.02265...][-8.49195...]] 331] 332 333[graph airy_zeros] 334 335[h4 Examples of finding Airy Zeros] 336 337[import ../../example/airy_zeros_example.cpp] 338 339[airy_zeros_example_1] 340[airy_zeros_example_2] 341 342Produces the program output: 343[pre 344boost::math::airy_ai_zero<double>(1) = -2.33811 345boost::math::airy_ai_zero<double>(2) = -4.08795 346boost::math::airy_bi_zero<double>(3) = -4.83074 347airy_ai_zeros: 348-2.33811 349-4.08795 350-5.52056 351-6.78671 352-7.94413 353 354boost::math::airy_bi_zero<float_type>(1) = -2.3381074104597670384891972524467354406385401456711 355boost::math::airy_bi_zero<float_type>(2) = -4.0879494441309706166369887014573910602247646991085 356boost::math::airy_bi_zero<float_type>(7) = -9.5381943793462388866329885451560196208390720763825 357airy_ai_zeros: 358-2.3381074104597670384891972524467354406385401456711 359-4.0879494441309706166369887014573910602247646991085 360-5.5205598280955510591298555129312935737972142806175 361] 362 363The full code (and output) for this example is at 364[@../../example/airy_zeros_example.cpp airy_zeros_example.cpp], 365 366[h3 Implementation] 367 368Given the following function (A&S 10.4.105): 369 370[equation airy_zero_1] 371 372Then an initial estimate for the n[super th] zero a[sub n] of Ai is given by (A&S 10.4.94): 373 374[equation airy_zero_2] 375 376and an initial estimate for the n[super th] zero b[sub n] of Bi is given by (A&S 10.4.98): 377 378[equation airy_zero_3] 379 380Thereafter the roots are refined using Newton iteration. 381 382[h3 Testing] 383 384The precision of evaluation of zeros was tested at 50 decimal digits using `cpp_dec_float_50` 385and found identical with spot values computed by __WolframAlpha. 386 387[endsect] [/section:airy_root Finding Zeros of Airy Functions] 388 389[endsect] [/section:airy Airy Functions] 390 391