[35] | 1 | /*! |
---|
| 2 | * \file |
---|
| 3 | * \brief Elementary mathematical functions - header file |
---|
| 4 | * \author Tony Ottosson and Adam Piatyszek |
---|
| 5 | * |
---|
| 6 | * ------------------------------------------------------------------------- |
---|
| 7 | * |
---|
| 8 | * IT++ - C++ library of mathematical, signal processing, speech processing, |
---|
| 9 | * and communications classes and functions |
---|
| 10 | * |
---|
| 11 | * Copyright (C) 1995-2007 (see AUTHORS file for a list of contributors) |
---|
| 12 | * |
---|
| 13 | * This program is free software; you can redistribute it and/or modify |
---|
| 14 | * it under the terms of the GNU General Public License as published by |
---|
| 15 | * the Free Software Foundation; either version 2 of the License, or |
---|
| 16 | * (at your option) any later version. |
---|
| 17 | * |
---|
| 18 | * This program is distributed in the hope that it will be useful, |
---|
| 19 | * but WITHOUT ANY WARRANTY; without even the implied warranty of |
---|
| 20 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the |
---|
| 21 | * GNU General Public License for more details. |
---|
| 22 | * |
---|
| 23 | * You should have received a copy of the GNU General Public License |
---|
| 24 | * along with this program; if not, write to the Free Software |
---|
| 25 | * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA |
---|
| 26 | * |
---|
| 27 | * ------------------------------------------------------------------------- |
---|
| 28 | */ |
---|
| 29 | |
---|
| 30 | #ifndef ELEM_MATH_H |
---|
| 31 | #define ELEM_MATH_H |
---|
| 32 | |
---|
| 33 | #ifndef _MSC_VER |
---|
| 34 | # include <itpp/config.h> |
---|
| 35 | #else |
---|
| 36 | # include <itpp/config_msvc.h> |
---|
| 37 | #endif |
---|
| 38 | |
---|
| 39 | #include <itpp/base/help_functions.h> |
---|
| 40 | #include <itpp/base/converters.h> |
---|
| 41 | #include <cstdlib> |
---|
| 42 | |
---|
| 43 | |
---|
| 44 | //!\addtogroup miscfunc |
---|
| 45 | //!@{ |
---|
| 46 | |
---|
| 47 | #ifndef HAVE_TGAMMA |
---|
| 48 | //! True gamma function |
---|
| 49 | double tgamma(double x); |
---|
| 50 | #endif |
---|
| 51 | |
---|
| 52 | #if !defined(HAVE_LGAMMA) || (HAVE_DECL_SIGNGAM != 1) |
---|
| 53 | //! Lograrithm of an absolute gamma function |
---|
| 54 | double lgamma(double x); |
---|
| 55 | //! Global variable needed by \c lgamma function |
---|
| 56 | extern int signgam; |
---|
| 57 | #endif |
---|
| 58 | |
---|
| 59 | #ifndef HAVE_CBRT |
---|
| 60 | //! Cubic root |
---|
| 61 | double cbrt(double x); |
---|
| 62 | #endif |
---|
| 63 | |
---|
| 64 | //!@} |
---|
| 65 | |
---|
| 66 | namespace itpp { |
---|
| 67 | |
---|
| 68 | //!\addtogroup miscfunc |
---|
| 69 | //!@{ |
---|
| 70 | |
---|
| 71 | // -------------------- sqr function -------------------- |
---|
| 72 | |
---|
| 73 | //! Square of x |
---|
| 74 | inline double sqr(double x) { return (x * x); } |
---|
| 75 | //! Square of complex-valued x, ||x||^2 |
---|
| 76 | inline double sqr(const std::complex<double>& x) |
---|
| 77 | { |
---|
| 78 | return (x.real() * x.real() + x.imag() * x.imag()); |
---|
| 79 | } |
---|
| 80 | //! Square of elements |
---|
| 81 | inline vec sqr(const vec &x) { return apply_function<double>(sqr, x); } |
---|
| 82 | //! Square of elements |
---|
| 83 | inline mat sqr(const mat &x) { return apply_function<double>(sqr, x); } |
---|
| 84 | //! Square of elements |
---|
| 85 | vec sqr(const cvec &x); |
---|
| 86 | //! Square of elements |
---|
| 87 | mat sqr(const cmat &x); |
---|
| 88 | |
---|
| 89 | |
---|
| 90 | // -------------------- abs function -------------------- |
---|
| 91 | |
---|
| 92 | //! Absolute value |
---|
| 93 | inline vec abs(const vec &x) { return apply_function<double>(std::fabs, x); } |
---|
| 94 | //! Absolute value |
---|
| 95 | inline mat abs(const mat &x) { return apply_function<double>(std::fabs, x); } |
---|
| 96 | //! Absolute value |
---|
| 97 | inline ivec abs(const ivec &x) { return apply_function<int>(std::abs, x); } |
---|
| 98 | //! Absolute value |
---|
| 99 | inline imat abs(const imat &x) { return apply_function<int>(std::abs, x); } |
---|
| 100 | //! Absolute value |
---|
| 101 | vec abs(const cvec &x); |
---|
| 102 | //! Absolute value |
---|
| 103 | mat abs(const cmat &x); |
---|
| 104 | |
---|
| 105 | |
---|
| 106 | // -------------------- sign/sgn functions -------------------- |
---|
| 107 | |
---|
| 108 | //! Signum function |
---|
| 109 | inline double sign(double x) |
---|
| 110 | { |
---|
| 111 | return (x == 0.0 ? 0.0 : (x < 0.0 ? -1.0 : 1.0)); |
---|
| 112 | } |
---|
| 113 | //! Signum function |
---|
| 114 | inline vec sign(const vec &x) { return apply_function<double>(sign, x); } |
---|
| 115 | //! Signum function |
---|
| 116 | inline mat sign(const mat &x) { return apply_function<double>(sign, x); } |
---|
| 117 | |
---|
| 118 | //! Signum function |
---|
| 119 | inline double sgn(double x) { return sign(x); } |
---|
| 120 | //! Signum function |
---|
| 121 | inline vec sgn(const vec &x) { return apply_function<double>(sign, x); } |
---|
| 122 | //! Signum function |
---|
| 123 | inline mat sgn(const mat &x) { return apply_function<double>(sign, x); } |
---|
| 124 | |
---|
| 125 | //! Signum function |
---|
| 126 | inline int sign_i(int x) |
---|
| 127 | { |
---|
| 128 | return (x == 0 ? 0 : (x < 0 ? -1 : 1)); |
---|
| 129 | } |
---|
| 130 | //! Signum function |
---|
| 131 | inline ivec sign_i(const ivec &x) { return apply_function<int>(sign_i, x); } |
---|
| 132 | //! Signum function |
---|
| 133 | inline imat sign_i(const imat &x) { return apply_function<int>(sign_i, x); } |
---|
| 134 | |
---|
| 135 | //! Signum function |
---|
| 136 | inline int sgn_i(int x) { return sign_i(x); } |
---|
| 137 | //! Signum function |
---|
| 138 | inline ivec sgn_i(const ivec &x) { return apply_function<int>(sign_i, x); } |
---|
| 139 | //! Signum function |
---|
| 140 | inline imat sgn_i(const imat &x) { return apply_function<int>(sign_i, x); } |
---|
| 141 | |
---|
| 142 | //! Signum function |
---|
| 143 | inline int sign_i(double x) |
---|
| 144 | { |
---|
| 145 | return (x == 0.0 ? 0 : (x < 0.0 ? -1 : 1)); |
---|
| 146 | } |
---|
| 147 | |
---|
| 148 | // -------------------- sqrt function -------------------- |
---|
| 149 | |
---|
| 150 | //! Square root of the elements |
---|
| 151 | inline vec sqrt(const vec &x) { return apply_function<double>(std::sqrt, x); } |
---|
| 152 | //! Square root of the elements |
---|
| 153 | inline mat sqrt(const mat &x) { return apply_function<double>(std::sqrt, x); } |
---|
| 154 | |
---|
| 155 | |
---|
| 156 | // -------------------- gamma function -------------------- |
---|
| 157 | |
---|
| 158 | //! Deprecated gamma function - please use tgamma() instead |
---|
| 159 | inline double gamma(double x) { return tgamma(x); } |
---|
| 160 | //! Deprecated gamma function for vectors. Will be changed to tgamma(). |
---|
| 161 | inline vec gamma(const vec &x) { return apply_function<double>(tgamma, x); } |
---|
| 162 | //! Deprecated gamma function for matrices. Will be changed to tgamma(). |
---|
| 163 | inline mat gamma(const mat &x) { return apply_function<double>(tgamma, x); } |
---|
| 164 | |
---|
| 165 | |
---|
| 166 | // -------------------- rem function -------------------- |
---|
| 167 | |
---|
| 168 | //! The reminder of the division x/y |
---|
| 169 | inline double rem(double x, double y) { return fmod(x, y); } |
---|
| 170 | //! Elementwise reminder of the division x/y for vec and double |
---|
| 171 | inline vec rem(const vec &x, double y) |
---|
| 172 | { |
---|
| 173 | return apply_function<double>(rem, x, y); |
---|
| 174 | } |
---|
| 175 | //! Elementwise reminder of the division x/y for double and vec |
---|
| 176 | inline vec rem(double x, const vec &y) |
---|
| 177 | { |
---|
| 178 | return apply_function<double>(rem, x, y); |
---|
| 179 | } |
---|
| 180 | //! Elementwise reminder of the division x/y for mat and double |
---|
| 181 | inline mat rem(const mat &x, double y) |
---|
| 182 | { |
---|
| 183 | return apply_function<double>(rem, x, y); |
---|
| 184 | } |
---|
| 185 | //! Elementwise reminder of the division x/y for double and mat |
---|
| 186 | inline mat rem(double x, const mat &y) |
---|
| 187 | { |
---|
| 188 | return apply_function<double>(rem, x, y); |
---|
| 189 | } |
---|
| 190 | |
---|
| 191 | // -------------------- mod function -------------------- |
---|
| 192 | |
---|
| 193 | //! Calculates the modulus, i.e. the signed reminder after division |
---|
| 194 | inline int mod(int k, int n) |
---|
| 195 | { |
---|
| 196 | return (n == 0) ? k : (k - n * floor_i(static_cast<double>(k) / n )); |
---|
| 197 | } |
---|
| 198 | |
---|
| 199 | |
---|
| 200 | // -------------------- factorial coefficient function -------------------- |
---|
| 201 | |
---|
| 202 | //! Calculates factorial coefficient for index <= 170. |
---|
| 203 | double fact(int index); |
---|
| 204 | |
---|
| 205 | |
---|
| 206 | // -------------------- binomial coefficient function -------------------- |
---|
| 207 | |
---|
| 208 | //! Compute the binomial coefficient "n over k". |
---|
| 209 | double binom(int n, int k); |
---|
| 210 | |
---|
| 211 | //! Compute the binomial coefficient "n over k". |
---|
| 212 | int binom_i(int n, int k); |
---|
| 213 | |
---|
| 214 | //! Compute the base 10 logarithm of the binomial coefficient "n over k". |
---|
| 215 | double log_binom(int n, int k); |
---|
| 216 | |
---|
| 217 | |
---|
| 218 | // -------------------- greatest common divisor function -------------------- |
---|
| 219 | |
---|
| 220 | /*! |
---|
| 221 | * \brief Compute the greatest common divisor (GCD) \a g of the elements |
---|
| 222 | * \a a and \a b. |
---|
| 223 | * |
---|
| 224 | * \a a and \a b must be non-negative integers. \a gdc(0, 0) is 0 by |
---|
| 225 | * convention; all other GCDs are positive integers. |
---|
| 226 | */ |
---|
| 227 | int gcd(int a, int b); |
---|
| 228 | |
---|
| 229 | |
---|
| 230 | // -------------------- complex related functions -------------------- |
---|
| 231 | |
---|
| 232 | //! Real part of complex values |
---|
| 233 | vec real(const cvec &x); |
---|
| 234 | //! Real part of complex values |
---|
| 235 | mat real(const cmat &x); |
---|
| 236 | //! Imaginary part of complex values |
---|
| 237 | vec imag(const cvec &x); |
---|
| 238 | //! Imaginary part of complex values |
---|
| 239 | mat imag(const cmat &x); |
---|
| 240 | |
---|
| 241 | //! Argument (angle) |
---|
| 242 | vec arg(const cvec &x); |
---|
| 243 | //! Argument (angle) |
---|
| 244 | mat arg(const cmat &x); |
---|
| 245 | //! Angle |
---|
| 246 | inline vec angle(const cvec &x) { return arg(x); } |
---|
| 247 | //! Angle |
---|
| 248 | inline mat angle(const cmat &x) { return arg(x); } |
---|
| 249 | |
---|
| 250 | // Added due to a failure in MSVC++ .NET 2005, which crashes on this |
---|
| 251 | // code. |
---|
| 252 | #ifndef _MSC_VER |
---|
| 253 | //! Conjugate of complex value |
---|
| 254 | inline cvec conj(const cvec &x) |
---|
| 255 | { |
---|
| 256 | return apply_function<std::complex<double> >(std::conj, x); |
---|
| 257 | } |
---|
| 258 | //! Conjugate of complex value |
---|
| 259 | inline cmat conj(const cmat &x) |
---|
| 260 | { |
---|
| 261 | return apply_function<std::complex<double> >(std::conj, x); |
---|
| 262 | } |
---|
| 263 | #else |
---|
| 264 | //! Conjugate of complex value |
---|
| 265 | cvec conj(const cvec &x); |
---|
| 266 | |
---|
| 267 | //! Conjugate of complex value |
---|
| 268 | cmat conj(const cmat &x); |
---|
| 269 | #endif |
---|
| 270 | |
---|
| 271 | //!@} |
---|
| 272 | |
---|
| 273 | } // namespace itpp |
---|
| 274 | |
---|
| 275 | #endif // #ifndef ELEM_MATH_H |
---|
| 276 | |
---|
| 277 | |
---|
| 278 | |
---|