1 /* 2 * Scythe Statistical Library Copyright (C) 2000-2002 Andrew D. Martin 3 * and Kevin M. Quinn; 2002-present Andrew D. Martin, Kevin M. Quinn, 4 * and Daniel Pemstein. All Rights Reserved. 5 * 6 * This program is free software; you can redistribute it and/or 7 * modify under the terms of the GNU General Public License as 8 * published by Free Software Foundation; either version 2 of the 9 * License, or (at your option) any later version. See the text files 10 * COPYING and LICENSE, distributed with this source code, for further 11 * information. 12 * -------------------------------------------------------------------- 13 * scythestat/smath.h 14 * 15 */ 16 17 /*! 18 * \file smath.h 19 * \brief Definitions for functions that perform common mathematical 20 * operations on every element of a Matrix. 21 * 22 * \note As is the case throughout the library, we provide both 23 * general and default template definitions of the Matrix-returning 24 * functions in this file, explicitly providing documentation for only 25 * the general template versions. As is also often the case, Doxygen 26 * does not always correctly add the default template definition to 27 * the function list below; there is always a default template 28 * definition available for every function. 29 * 30 */ 31 32 #ifndef SCYTHE_MATH_H 33 #define SCYTHE_MATH_H 34 35 #ifdef SCYTHE_COMPILE_DIRECT 36 #include "matrix.h" 37 #include "algorithm.h" 38 #include "error.h" 39 #else 40 #include "scythestat/matrix.h" 41 #include "scythestat/algorithm.h" 42 #include "scythestat/error.h" 43 #endif 44 45 #include <cmath> 46 #include <numeric> 47 #include <set> 48 49 namespace scythe { 50 51 namespace { 52 typedef unsigned int uint; 53 } 54 55 /* Almost every function in this file follows one of the two patterns 56 * described by these macros. The first macro handles single-argument 57 * functions. The second handles two-matrix-argument functions (or 58 * scalar-matrix, matrix-scalar. The second macro also permits 59 * cross-type operations (these are limited only by the capabilities 60 * of the underlying functions). 61 */ 62 #define SCYTHE_MATH_OP(NAME, OP) \ 63 template <matrix_order RO, matrix_style RS, typename T, \ 64 matrix_order PO, matrix_style PS> \ 65 Matrix<T,RO,RS> \ 66 NAME (const Matrix<T,PO,PS>& A) \ 67 { \ 68 Matrix<T,RO,RS> res(A.rows(), A.cols(), false); \ 69 std::transform(A.begin_f(), A.end_f(), res.begin_f(), (T (*) (T))OP); \ 70 return res; \ 71 } \ 72 \ 73 template <typename T, matrix_order O, matrix_style S> \ 74 Matrix<T,O,Concrete> \ 75 NAME (const Matrix<T,O,S>& A) \ 76 { \ 77 return NAME<O,Concrete>(A); \ 78 } 79 80 #define SCYTHE_MATH_OP_2ARG(NAME, OP) \ 81 template <matrix_order RO, matrix_style RS, typename T, \ 82 matrix_order PO1, matrix_style PS1, \ 83 matrix_order PO2, matrix_style PS2, typename S> \ 84 Matrix<T,RO,RS> \ 85 NAME (const Matrix<T,PO1,PS1>& A, const Matrix<S,PO2,PS2>& B) \ 86 { \ 87 SCYTHE_CHECK_10 (A.size() != 1 && B.size() != 1 && \ 88 A.size() != B.size(), scythe_conformation_error, \ 89 "Matrices with dimensions (" << A.rows() \ 90 << ", " << A.cols() \ 91 << ") and (" << B.rows() << ", " << B.cols() \ 92 << ") are not conformable"); \ 93 \ 94 Matrix<T,RO,RS> res; \ 95 \ 96 using std::placeholders::_1; \ 97 \ 98 if (A.size() == 1) { \ 99 res.resize2Match(B); \ 100 std::transform(B.template begin_f<RO>(), B.template end_f<RO>(),\ 101 res.begin_f(), std::bind(static_cast<T(*)(T,S)>(&OP), A(0), _1)); \ 102 } else if (B.size() == 1) { \ 103 res.resize2Match(A); \ 104 std::transform(A.template begin_f<RO>(), A.template end_f<RO>(),\ 105 res.begin_f(), std::bind(static_cast<T(*)(T,S)>(&OP), _1, B(0))); \ 106 } else { \ 107 res.resize2Match(A); \ 108 std::transform(A.template begin_f<RO>(), A.template end_f<RO>(),\ 109 B.template begin_f<RO>(), res.begin_f(), (T (*) (T, S))OP); \ 110 } \ 111 \ 112 return res; \ 113 } \ 114 \ 115 template <typename T, matrix_order PO1, matrix_style PS1, \ 116 matrix_order PO2, matrix_style PS2, \ 117 typename S> \ 118 Matrix<T,PO1,Concrete> \ 119 NAME (const Matrix<T,PO1,PS1>& A, const Matrix<S,PO2,PS2>& B) \ 120 { \ 121 return NAME<PO1,Concrete>(A, B); \ 122 } \ 123 \ 124 template<matrix_order RO, matrix_style RS, typename T, \ 125 matrix_order PO, matrix_style PS, typename S> \ 126 Matrix<T,RO,RS> \ 127 NAME (const Matrix<T,PO,PS>& A, S b) \ 128 { \ 129 return NAME<RO,RS>(A, Matrix<S,RO,Concrete>(b)); \ 130 } \ 131 \ 132 template <typename T, typename S, matrix_order PO, matrix_style PS> \ 133 Matrix<T,PO,Concrete> \ 134 NAME (const Matrix<T,PO,PS>& A, S b) \ 135 { \ 136 return NAME<PO,Concrete>(A, Matrix<S,PO,Concrete>(b)); \ 137 } \ 138 \ 139 template<matrix_order RO, matrix_style RS, typename T, \ 140 matrix_order PO, matrix_style PS, typename S> \ 141 Matrix<T,RO,RS> \ 142 NAME (T a, const Matrix<S,PO,PS>& B) \ 143 { \ 144 return NAME<RO,RS>(Matrix<S, RO,Concrete>(a), B); \ 145 } \ 146 \ 147 template <typename T, typename S, matrix_order PO, matrix_style PS> \ 148 Matrix<T,PO,Concrete> \ 149 NAME (T a, const Matrix<S,PO,PS>& B) \ 150 { \ 151 return NAME<PO,Concrete>(Matrix<S,PO,Concrete>(a), B); \ 152 } 153 154 155 /* calc the inverse cosine of each element of a Matrix */ 156 157 /*! 158 * \brief Calculate the inverse cosine of each element of a Matrix 159 * 160 * This function calculates the inverse cosine of each element in a Matrix 161 * 162 * \param A The matrix whose inverse cosines are of interest. 163 * 164 * \see tan() 165 * \see tanh() 166 * \see sin() 167 * \see sinh() 168 * \see cos() 169 * \see cosh() 170 * \see acosh() 171 * \see asin() 172 * \see asinh() 173 * \see atan() 174 * \see atanh() 175 * \see atan2() 176 */ 177 SCYTHE_MATH_OP(acos,::acos)178 SCYTHE_MATH_OP(acos, ::acos) 179 180 /* calc the inverse hyperbolic cosine of each element of a Matrix */ 181 /*! 182 * \brief Calculate the inverse hyperbolic cosine of each element of a Matrix 183 * 184 * This function calculates the inverse hyperbolic cosine of each element 185 * in a Matrix 186 * 187 * \param A The matrix whose inverse hyperbolic cosines are of interest. 188 * 189 * \see tan() 190 * \see tanh() 191 * \see sin() 192 * \see sinh() 193 * \see cos() 194 * \see cosh() 195 * \see acos() 196 * \see asin() 197 * \see asinh() 198 * \see atan() 199 * \see atanh() 200 * \see atan2() 201 */ 202 203 SCYTHE_MATH_OP(acosh, ::acosh) 204 205 /* calc the inverse sine of each element of a Matrix */ 206 207 /*! 208 * \brief Calculate the inverse sine of each element of a Matrix 209 * 210 * This function calculates the inverse sine of each element 211 * in a Matrix 212 * 213 * \param A The matrix whose inverse sines are of interest. 214 * 215 * \see tan() 216 * \see tanh() 217 * \see sin() 218 * \see sinh() 219 * \see cos() 220 * \see cosh() 221 * \see acos() 222 * \see acosh() 223 * \see asinh() 224 * \see atan() 225 * \see atanh() 226 * \see atan2() 227 */ 228 229 SCYTHE_MATH_OP(asin, ::asin) 230 231 /* calc the inverse hyperbolic sine of each element of a Matrix */ 232 233 /*! 234 * \brief Calculate the inverse hyperbolic sine of each element of a Matrix 235 * 236 * This function calculates the inverse hyperbolic sine of each element 237 * in a Matrix 238 * 239 * \param A The matrix whose inverse hyperbolic sines are of interest. 240 * 241 * \see tan() 242 * \see tanh() 243 * \see sin() 244 * \see sinh() 245 * \see cos() 246 * \see cosh() 247 * \see acos() 248 * \see acosh() 249 * \see asin() 250 * \see atan() 251 * \see atanh() 252 * \see atan2() 253 */ 254 255 SCYTHE_MATH_OP(asinh, ::asinh) 256 257 /* calc the inverse tangent of each element of a Matrix */ 258 259 /*! 260 * \brief Calculate the inverse tangent of each element of a Matrix 261 * 262 * This function calculates the inverse tangent of each element 263 * in a Matrix 264 * 265 * \param A The matrix whose inverse tangents are of interest. 266 * 267 * \see tan() 268 * \see tanh() 269 * \see sin() 270 * \see sinh() 271 * \see cos() 272 * \see cosh() 273 * \see acos() 274 * \see acosh() 275 * \see asin() 276 * \see asin() 277 * \see atanh() 278 * \see atan2() 279 */ 280 281 SCYTHE_MATH_OP(atan, ::atan) 282 283 /* calc the inverse hyperbolic tangent of each element of a Matrix */ 284 /*! 285 * \brief Calculate the inverse hyperbolic tangent of each element of a Matrix 286 * 287 * This function calculates the inverse hyperbolic tangent of each element 288 * in a Matrix 289 * 290 * \param A The matrix whose inverse hyperbolic tangents are of interest. 291 * 292 * \see tan() 293 * \see tanh() 294 * \see sin() 295 * \see sinh() 296 * \see cos() 297 * \see cosh() 298 * \see acos() 299 * \see acosh() 300 * \see asin() 301 * \see asinh() 302 * \see atan() 303 * \see atan2() 304 */ 305 306 SCYTHE_MATH_OP(atanh, ::atanh) 307 308 /* calc the angle whose tangent is y/x */ 309 310 /*! 311 * \brief Calculate the angle whose tangent is y/x 312 * 313 * This function calculates the angle whose tangent is y/x, given two 314 * matrices A and B (where y is the ith element of A, and x is the jth element 315 * of matrix B). 316 * 317 * \param A The matrix of y values 318 * \param B The matrix of x values 319 * 320 * \see tan() 321 * \see tanh() 322 * \see sin() 323 * \see sinh() 324 * \see cos() 325 * \see cosh() 326 * \see acos() 327 * \see acosh() 328 * \see asin() 329 * \see asinh() 330 * \see atan() 331 * \see atanh() 332 */ 333 334 SCYTHE_MATH_OP_2ARG(atan2, ::atan2) 335 336 /* calc the cube root of each element of a Matrix */ 337 /*! 338 * \brief Calculate the cube root of each element of a Matrix 339 * 340 * This function calculates the cube root of each element 341 * in a Matrix 342 * 343 * \param A The matrix whose cube roots are of interest. 344 * 345 * \see sqrt() 346 */ 347 348 SCYTHE_MATH_OP(cbrt, ::cbrt) 349 350 /* calc the ceil of each element of a Matrix */ 351 /*! 352 * \brief Calculate the ceiling of each element of a Matrix 353 * 354 * This function calculates the ceiling of each element 355 * in a Matrix 356 * 357 * \param A The matrix whose ceilings are of interest. 358 * 359 * \see floor() 360 */ 361 362 SCYTHE_MATH_OP(ceil, ::ceil) 363 364 /* create a matrix containing the absval of the first input and the 365 * sign of the second 366 */ 367 /*! 368 * \brief Create a matrix containing the absolute value of the first input 369 * and the sign of the second input 370 * 371 * This function creates a matrix containing the absolute value of the first 372 * input, a matrix called A, and the sign of the second input, matrix B. 373 * 374 * \param A The matrix whose absolute values will comprise the resultant matrix. 375 * \param B The matrix whose signs will comprise the resultant matrix 376 */ 377 378 SCYTHE_MATH_OP_2ARG(copysign, ::copysign) 379 380 /* calc the cosine of each element of a Matrix */ 381 382 /*! 383 * \brief Calculate the cosine of each element of a Matrix 384 * 385 * This function calculates the cosine of each element in a Matrix 386 * 387 * \param A The matrix whose cosines are of interest. 388 * 389 * \see tan() 390 * \see tanh() 391 * \see sin() 392 * \see sinh() 393 * \see cosh() 394 * \see acos() 395 * \see acosh() 396 * \see asin() 397 * \see asinh() 398 * \see atan() 399 * \see atanh() 400 * \see atan2() 401 */ 402 403 SCYTHE_MATH_OP(cos, ::cos) 404 405 /* calc the hyperbolic cosine of each element of a Matrix */ 406 /*! 407 * \brief Calculate the hyperbolic cosine of each element of a Matrix 408 * 409 * This function calculates the hyperbolic cosine of each element in a Matrix 410 * 411 * \param A The matrix whose hyperbolic cosines are of interest. 412 * 413 * \see tan() 414 * \see tanh() 415 * \see sin() 416 * \see sinh() 417 * \see cos() 418 * \see acos() 419 * \see acosh() 420 * \see asin() 421 * \see asinh() 422 * \see atan() 423 * \see atanh() 424 * \see atan2() 425 */ 426 427 SCYTHE_MATH_OP(cosh, ::cosh) 428 429 /* calc the error function of each element of a Matrix */ 430 /*! 431 * \brief Calculate the error function of each element of a Matrix 432 * 433 * This function calculates the error function of each element in a Matrix 434 * 435 * \param A The matrix whose error functions are of interest. 436 * 437 * \see erfc() 438 */ 439 440 SCYTHE_MATH_OP(erf, ::erf) 441 442 /* calc the complementary error function of each element of a Matrix */ 443 /*! 444 * \brief Calculate the complementary error function of each element of a Matrix 445 * 446 * This function calculates the complemenatry error function of each 447 * element in a Matrix 448 * 449 * \param A The matrix whose complementary error functions are of interest. 450 * 451 * \see erf() 452 */ 453 454 SCYTHE_MATH_OP(erfc, ::erfc) 455 456 /* calc the vaue e^x of each element of a Matrix */ 457 /*! 458 * \brief Calculate the value e^x for each element of a Matrix 459 * 460 * This function calculates the value e^x for each element of a matrix, where 461 * x is the ith element of the matrix A 462 * 463 * \param A The matrix whose elements are to be exponentiated. 464 * 465 * \see expm1() 466 */ 467 468 SCYTHE_MATH_OP(exp, ::exp) 469 470 /* calc the exponent - 1 of each element of a Matrix */ 471 /*! 472 * \brief Calculate the value e^(x-1) for each element of a Matrix 473 * 474 * This function calculates the value e^(x-1) for each element of a matrix, where 475 * x is the ith element of the matrix A 476 * 477 * \param A The matrix whose elements are to be exponentiated. 478 * 479 * \see exp() 480 */ 481 482 SCYTHE_MATH_OP(expm1, ::expm1) 483 484 /* calc the absval of each element of a Matrix */ 485 /*! 486 * \brief Calculate the absolute value of each element of a Matrix 487 * 488 * This function calculates the absolute value of each element in a Matrix 489 * 490 * \param A The matrix whose absolute values are to be taken. 491 */ 492 493 SCYTHE_MATH_OP(fabs, (T (*) (T))::fabs) 494 495 /* calc the floor of each element of a Matrix */ 496 /*! 497 * \brief Calculate the floor of each element of a Matrix 498 * 499 * This function calculates the floor of each element 500 * in a Matrix 501 * 502 * \param A The matrix whose floors are of interest. 503 * 504 * \see ceil() 505 */ 506 507 SCYTHE_MATH_OP(floor, ::floor) 508 509 /* calc the remainder of the division of each matrix element */ 510 /*! 511 * \brief Calculate the remainder of the division of each matrix element 512 * 513 * This function calculates the remainder when the elements of Matrix A are 514 * divided by the elements of Matrix B. 515 * 516 * \param A The matrix to serve as dividend 517 * \param B the matrix to serve as divisor 518 */ 519 520 SCYTHE_MATH_OP_2ARG(fmod, ::fmod) 521 522 /* calc the fractional val of input and return exponents in int 523 * matrix reference 524 */ 525 526 /*! 527 */ 528 template <matrix_order RO, matrix_style RS, typename T, 529 matrix_order PO1, matrix_style PS1, 530 matrix_order PO2, matrix_style PS2> 531 Matrix<T,RO,RS> 532 frexp (const Matrix<T,PO1,PS1>& A, Matrix<int,PO2,PS2>& ex) 533 { 534 SCYTHE_CHECK_10(A.size() != ex.size(), scythe_conformation_error, 535 "The input matrix sizes do not match"); 536 Matrix<T,PO1,Concrete> res(A.rows(), A.cols()); 537 538 typename Matrix<T,PO1,PS1>::const_forward_iterator it; 539 typename Matrix<T,PO1,Concrete>::forward_iterator rit 540 = res.begin_f(); 541 typename Matrix<int,PO2,PS2>::const_forward_iterator it2 542 = ex.begin_f(); 543 for (it = A.begin_f(); it != A.end_f(); ++it) { 544 *rit = ::frexp(*it, &(*it2)); 545 ++it2; ++rit; 546 } 547 548 return res; 549 } 550 551 template <typename T, matrix_order PO1, matrix_style PS1, 552 matrix_order PO2, matrix_style PS2> 553 Matrix<T,PO1,Concrete> frexp(Matrix<T,PO1,PS1> & A,Matrix<int,PO2,PS2> & ex)554 frexp (Matrix<T,PO1,PS1>& A, Matrix<int,PO2,PS2>& ex) 555 { 556 return frexp<PO1,Concrete>(A,ex); 557 } 558 559 /* calc the euclidean distance between the two inputs */ 560 /*! 561 * \brief Calculate the euclidean distance between two inputs 562 * 563 * This function calculates the euclidean distance between the elements of Matrix 564 * A and the elements of Matrix B. 565 * 566 * \param A Input matrix 567 * \param B Input matrix 568 */ 569 SCYTHE_MATH_OP_2ARG(hypot,::hypot)570 SCYTHE_MATH_OP_2ARG(hypot, ::hypot) 571 572 /* return (int) logb */ 573 SCYTHE_MATH_OP(ilogb, ::ilogb) 574 575 /* compute the bessel func of the first kind of the order 0 */ 576 /*! 577 * \brief Compute the Bessel function of the first kind of the order 0 578 * 579 * This function computes the Bessel function of the first kind of order 0 580 * for each element in the input matrix, A. 581 * 582 * \param A Matrix for which the Bessel function is of interest 583 * 584 * \see j1() 585 * \see jn() 586 * \see y0() 587 * \see y1() 588 * \see yn() 589 */ 590 591 SCYTHE_MATH_OP(j0, ::j0) 592 593 /* compute the bessel func of the first kind of the order 1 */ 594 /*! 595 * \brief Compute the Bessel function of the first kind of the order 1 596 * 597 * This function computes the Bessel function of the first kind of order 1 598 * for each element in the input matrix, A. 599 * 600 * \param A Matrix for which the Bessel function is of interest 601 * 602 * \see j0() 603 * \see jn() 604 * \see y0() 605 * \see y1() 606 * \see yn() 607 */ 608 609 SCYTHE_MATH_OP(j1, ::j1) 610 611 /* compute the bessel func of the first kind of the order n 612 * TODO: This definition causes the compiler to issue some warnings. 613 * Fix 614 */ 615 /*! 616 * \brief Compute the Bessel function of the first kind of the order n 617 * 618 * This function computes the Bessel function of the first kind of order n 619 * for each element in the input matrix, A. 620 * 621 * \param n Order of the Bessel function 622 * \param A Matrix for which the Bessel function is of interest 623 * 624 * \see j0() 625 * \see j1() 626 * \see y0() 627 * \see y1() 628 * \see yn() 629 */ 630 631 SCYTHE_MATH_OP_2ARG(jn, ::jn) 632 633 /* calc x * 2 ^ex */ 634 /*! 635 * \brief Compute x * 2^ex 636 * 637 * This function computes the value of x * 2^ex, where x is the ith element of 638 * the input matrix A, and ex is the desired value of the exponent. 639 * 640 * \param A Matrix whose elements are to be multiplied 641 * \param ex Matrix of powers to which 2 will be raised. 642 */ 643 SCYTHE_MATH_OP_2ARG(ldexp, ::ldexp) 644 645 /* compute the natural log of the absval of gamma function */ 646 647 /*! 648 * \brief Compute the natural log of the absolute value of the gamma function 649 * 650 * This function computes the absolute value of the Gamma Function, evaluated at 651 * each element of the input matrix A. 652 * 653 * \param A Matrix whose elements will serve as inputs for the Gamma Function 654 * 655 * \see log() 656 */ 657 658 SCYTHE_MATH_OP(lgamma, ::lgamma) 659 660 /* calc the natural log of each element of a Matrix */ 661 /*! 662 * \brief Compute the natural log of each element of a Matrix 663 * 664 * This function computes the natural log of each element in a matrix, A. 665 * 666 * \param A Matrix whose natural logs are of interest 667 * 668 * \see log10() 669 * \see log1p() 670 * \see logb() 671 */ 672 673 SCYTHE_MATH_OP(log, (T (*)(T))::log) 674 675 /* calc the base-10 log of each element of a Matrix */ 676 /*! 677 * \brief Compute the log base 10 of each element of a Matrix 678 * 679 * This function computes the log base 10 of each element in a matrix, A. 680 * 681 * \param A Matrix whose logs are of interest 682 * 683 * \see log() 684 * \see log1p() 685 * \see logb() 686 */ 687 688 SCYTHE_MATH_OP(log10, ::log10) 689 690 /* calc the natural log of 1 + each element of a Matrix */ 691 /*! 692 * \brief Compute the natural log of 1 + each element of a Matrix 693 * 694 * This function computes the natural log of 1 + each element of a Matrix. 695 * 696 * \param A Matrix whose logs are of interest 697 * 698 * \see log() 699 * \see log10() 700 * \see logb() 701 */ 702 703 SCYTHE_MATH_OP(log1p, ::log1p) 704 705 /* calc the logb of each element of a Matrix */ 706 /*! 707 * \brief Compute the logb each element of a Matrix 708 * 709 * This function computes the log base b of each element of a Matrix. 710 * 711 * \param A Matrix whose logs are of interest 712 * 713 * \see log() 714 * \see log10() 715 * \see log1p() 716 */ 717 718 SCYTHE_MATH_OP(logb, ::logb) 719 720 /* x = frac + i, return matrix of frac and place i in 2nd matrix 721 */ 722 template <matrix_order RO, matrix_style RS, typename T, 723 matrix_order PO1, matrix_style PS1, 724 matrix_order PO2, matrix_style PS2> 725 Matrix<T,RO,RS> 726 modf (const Matrix<T,PO1,PS1>& A, Matrix<double,PO2,PS2>& ipart) 727 { 728 SCYTHE_CHECK_10(A.size() != ipart.size(), scythe_conformation_error, 729 "The input matrix sizes do not match"); 730 Matrix<T,PO1,Concrete> res(A.rows(), A.cols()); 731 732 typename Matrix<T,PO1,PS1>::const_forward_iterator it; 733 typename Matrix<T,PO1,Concrete>::forward_iterator rit 734 = res.begin_f(); 735 typename Matrix<double,PO2,PS2>::const_forward_iterator it2 736 = ipart.begin_f(); 737 for (it = A.begin_f(); it != A.end_f(); ++it) { 738 *rit = ::modf(*it, &(*it2)); 739 ++it2; ++rit; 740 } 741 742 return res; 743 } 744 745 template <typename T, matrix_order PO1, matrix_style PS1, 746 matrix_order PO2, matrix_style PS2> 747 Matrix<T,PO1,Concrete> modf(Matrix<T,PO1,PS1> & A,Matrix<double,PO2,PS2> & ipart)748 modf (Matrix<T,PO1,PS1>& A, Matrix<double,PO2,PS2>& ipart) 749 { 750 return modf<PO1,Concrete>(A,ipart); 751 } 752 753 /* calc x^ex of each element of a Matrix */ 754 755 /*! 756 * \brief Compute x^ex for each element of a matrix 757 * 758 * This function computes x^ex, where x is the ith element of the matrix A, 759 * and ex is the desired exponent. 760 * 761 * \param A Matrix to be exponentiated 762 * \param ex Desired exponent 763 */ 764 SCYTHE_MATH_OP_2ARG(pow, std::pow) 765 766 /* calc rem == x - n * y */ 767 SCYTHE_MATH_OP_2ARG(remainder, ::remainder) 768 769 /* return x rounded to nearest int */ 770 771 /*! 772 * \brief Return x rounded to the nearest integer 773 * 774 * This function returns x, where x is the ith element of the Matrix A, 775 * rounded to the nearest integer. 776 * 777 * \param A Matrix whose elements are to be rounded 778 */ 779 780 SCYTHE_MATH_OP(rint, ::rint) 781 782 /* returns x * FLT_RADIX^ex */ 783 SCYTHE_MATH_OP_2ARG(scalbn, ::scalbn) 784 785 /* calc the sine of x */ 786 787 /*! 788 * \brief Calculate the sine of each element of a Matrix 789 * 790 * This function calculates the sine of each element in a Matrix 791 * 792 * \param A The matrix whose sines are of interest. 793 * 794 * \see tan() 795 * \see tanh() 796 * \see sinh() 797 * \see cos() 798 * \see cosh() 799 * \see acos() 800 * \see acosh() 801 * \see asin() 802 * \see asinh() 803 * \see atan() 804 * \see atanh() 805 * \see atan2() 806 */ 807 808 SCYTHE_MATH_OP(sin, ::sin) 809 810 /* calc the hyperbolic sine of x */ 811 /*! 812 * \brief Calculate the hyperbolic sine of each element of a Matrix 813 * 814 * This function calculates the hyperbolic sine of each element in a Matrix 815 * 816 * \param A The matrix whose hyperbolic sines are of interest. 817 * 818 * \see tan() 819 * \see tanh() 820 * \see sin() 821 * \see cos() 822 * \see cosh() 823 * \see acos() 824 * \see acosh() 825 * \see asin() 826 * \see asinh() 827 * \see atan() 828 * \see atanh() 829 * \see atan2() 830 */ 831 832 SCYTHE_MATH_OP(sinh, ::sinh) 833 834 /* calc the sqrt of x */ 835 /*! 836 * \brief Calculate the square root of each element in a matrix 837 * 838 * This function calculates the square root of each element in a Matrix 839 * 840 * \param A The matrix whose roots are of interest. 841 * 842 * \see cbrt() 843 844 */ 845 846 SCYTHE_MATH_OP(sqrt, (T (*)(T))::sqrt) 847 848 /* calc the tangent of x */ 849 850 /*! 851 * \brief Calculate the tangent of each element of a Matrix 852 * 853 * This function calculates the tangent of each element in a Matrix 854 * 855 * \param A The matrix whose tangents are of interest. 856 * 857 * \see sinh() 858 * \see tanh() 859 * \see sin() 860 * \see cos() 861 * \see cosh() 862 * \see acos() 863 * \see acosh() 864 * \see asin() 865 * \see asinh() 866 * \see atan() 867 * \see atanh() 868 * \see atan2() 869 */ 870 871 SCYTHE_MATH_OP(tan, ::tan) 872 873 /* calc the hyperbolic tangent of x */ 874 /*! 875 * \brief Calculate the hyperbolic tangent of each element of a Matrix 876 * 877 * This function calculates the hyperbolic tangent of each element in a Matrix 878 * 879 * \param A The matrix whose hyperbolic tangents are of interest. 880 * 881 * \see sinh() 882 * \see tan() 883 * \see sin() 884 * \see cos() 885 * \see cosh() 886 * \see acos() 887 * \see acosh() 888 * \see asin() 889 * \see asinh() 890 * \see atan() 891 * \see atanh() 892 * \see atan2() 893 */ 894 895 SCYTHE_MATH_OP(tanh, ::tanh) 896 897 /* bessel function of the second kind of order 0*/ 898 /*! 899 * \brief Compute the Bessel function of the second kind of order 0 900 * 901 * This function computes the Bessel function of the second kind of order 0 902 * for each element in the input matrix, A. 903 * 904 * \param A Matrix for which the Bessel function is of interest 905 * 906 * \see j0() 907 * \see j1() 908 * \see jn() 909 * \see y1() 910 * \see yn() 911 */ 912 913 SCYTHE_MATH_OP(y0, ::y0) 914 915 /* bessel function of the second kind of order 1*/ 916 /*! 917 * \brief Compute the Bessel function of the second kind of order 1 918 * 919 * This function computes the Bessel function of the second kind of order 1 920 * for each element in the input matrix, A. 921 * 922 * \param A Matrix for which the Bessel function is of interest 923 * 924 * \see j0() 925 * \see j1() 926 * \see jn() 927 * \see y0() 928 * \see yn() 929 */ 930 931 SCYTHE_MATH_OP(y1, ::y1) 932 933 /* bessel function of the second kind of order n 934 * TODO: This definition causes the compiler to issue some warnings. 935 * Fix 936 */ 937 /*! 938 * \brief Compute the Bessel function of the second kind of order n 939 * 940 * This function computes the Bessel function of the second kind of order n 941 * for each element in the input matrix, A. 942 * 943 * \param n Order of the Bessel function 944 * \param A Matrix for which the Bessel function is of interest 945 * 946 * \see j0() 947 * \see j1() 948 * \see jn() 949 * \see y0() 950 * \see y1() 951 */ 952 953 SCYTHE_MATH_OP_2ARG(yn, ::yn) 954 955 } // end namespace scythe 956 957 #endif /* SCYTHE_MATH_H */ 958