MercuryDPM  Trunk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ExtendedMath.h
Go to the documentation of this file.
1 //Copyright (c) 2013-2020, The MercuryDPM Developers Team. All rights reserved.
2 //For the list of developers, see <http://www.MercuryDPM.org/Team>.
3 //
4 //Redistribution and use in source and binary forms, with or without
5 //modification, are permitted provided that the following conditions are met:
6 // * Redistributions of source code must retain the above copyright
7 // notice, this list of conditions and the following disclaimer.
8 // * Redistributions in binary form must reproduce the above copyright
9 // notice, this list of conditions and the following disclaimer in the
10 // documentation and/or other materials provided with the distribution.
11 // * Neither the name MercuryDPM nor the
12 // names of its contributors may be used to endorse or promote products
13 // derived from this software without specific prior written permission.
14 //
15 //THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
16 //ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
17 //WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
18 //DISCLAIMED. IN NO EVENT SHALL THE MERCURYDPM DEVELOPERS TEAM BE LIABLE FOR ANY
19 //DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
20 //(INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
21 //LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
22 //ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
23 //(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24 //SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25 
26 #ifndef EXTENDEDMATH_H
27 #define EXTENDEDMATH_H
28 
29 #include <iostream> //std::istream and std::stringstream
30 #include <fstream> //std::fstream
31 #include <cmath>
32 #include <complex>
33 #include <limits>
34 
35 #include "NumericalVector.h"
36 #include "Vector.h"
37 #include "Quaternion.h"
38 
39 /*
40  * \brief
41  */
42 namespace constants
43 {
44 //Values from WolframAlpha
45 const Mdouble pi = 3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117068;
46 const Mdouble sqrt_pi = 1.772453850905516027298167483341145182797549456122387128213807789852911284591032181374950656738544665;
47 const Mdouble sqr_pi = 9.869604401089358618834490999876151135313699407240790626413349376220044822419205243001773403718552232;
48 const Mdouble sqrt_2 = 1.414213562373095048801688724209698078569671875376948073176679737990732478462107038850387534327641573;
49 const Mdouble sqrt_3 = 1.732050807568877293527446341505872366942805253810380628055806979451933016908800037081146186757248576;
50 const Mdouble R = 8.31446261815324;
51 const std::complex<Mdouble> i = {0.0, 1.0};
52 const Mdouble degree = pi / 180.; // degree-to-radian conversion
53 }
54 
58 namespace mathsFunc
59 {
63 Mdouble gamma(Mdouble gamma_in);
64 
69 
70 
74 Mdouble chi_squared(Mdouble x, unsigned int k);
75 
79 Mdouble chi_squared_prob(Mdouble x, unsigned int k);
80 
90 Mdouble goldenSectionSearch(Mdouble (* function)(const Mdouble), Mdouble min, Mdouble cur, Mdouble max,
91  Mdouble endCondition, Mdouble curVal = std::numeric_limits<Mdouble>::quiet_NaN());
92 
96 template<typename T>
97 int sign(T val)
98 {
99  return (T(0) < val) - (val < T(0));
100 }
101 
105 template<typename T>
106 T square(const T val)
107 {
108  return val * val;
109 }
110 
114 template<typename T>
115 T cubic(const T val)
116 {
117  return val * val * val;
118 }
119 
127 bool isEqual(Mdouble v1, Mdouble v2, Mdouble absError);
128 
136 bool isEqual(Vec3D v1, Vec3D v2, Mdouble absError);
137 
145 bool isEqual(Matrix3D m1, Matrix3D m2, Mdouble absError);
146 
147 bool isEqual(MatrixSymmetric3D m1, MatrixSymmetric3D m2, Mdouble absError);
148 
149 bool isEqual(Quaternion v1, Quaternion v2, double absError);
150 
154 template<typename T>
155 constexpr T factorial(const T t)
156 {
157  return (t == 0) ? 1 : t * factorial(t - 1);
158 }
159 
160 //platform independent implementation of sine and cosine, taken from
161 // http://stackoverflow.com/questions/18662261/fastest-implementation-of-sine-cosine-and-square-root-in-c-doesnt-need-to-b
162 // (cosine was implemented wrongly on the website, here is a corrected version)
163 
164 // sin(x) = x - x^3/3! + x^5/5! - x^7/7! + ...
165 Mdouble sin(Mdouble x);
166 
167 // cos(x) = 1 - x^2/2! + x^4/4! - x^6/6! + ...
168 Mdouble cos(Mdouble x);
169 
170 Mdouble exp(Mdouble Exponent);
171 
172 Mdouble log(Mdouble Power);
173 
174 
176 // tan=sin/cos
177 template<typename T>
178 T tan(T x)
179 {
180  return sin(x) / cos(x);
181 }
182 
183 
187 Mdouble chebyshev(Mdouble x, const Mdouble coef[], int N);
188 
190 
191 Mdouble I0(Mdouble x);
192 
193 }
194 
195 /*
196  * \brief Namespace for functions required to calculate spherical harmonics
197  */
198 
200 {
201 
202 //Compute all the associated LegenderePolynomials up to order n, and only positive order m at location x
204 
205 //Compute all spherical harmonics up to order p, at angles theta and phi
207 
208 //Compute all squaredFactorials (see eqn 5.23 in a short course on fast multipole methods) up to order p
210 }
211 
212 #endif
Implementation of a 3D quaternion (by Vitaliy).
Definition: Quaternion.h:62
double Mdouble
Definition: GeneralDefine.h:34
Mdouble goldenSectionSearch(Mdouble(*function)(const Mdouble), Mdouble min, Mdouble cur, Mdouble max, Mdouble endCondition, Mdouble curVal=std::numeric_limits< Mdouble >::quiet_NaN())
This function performs a golden section search to find the location of the minimum of a function...
Mdouble exp(Mdouble Exponent)
Definition: ExtendedMath.cc:84
Mdouble beta(Mdouble z, Mdouble w)
This is the beta function, returns the approximation based on cmath's implementation of ln(gamma) ...
const std::complex< Mdouble > i
Definition: ExtendedMath.h:51
Mdouble I0(Mdouble x)
const Mdouble sqrt_pi
Definition: ExtendedMath.h:46
int sign(T val)
This is a sign function, it returns -1 for negative numbers, 1 for positive numbers and 0 for 0...
Definition: ExtendedMath.h:97
T cubic(const T val)
calculates the cube of a number
Definition: ExtendedMath.h:115
Mdouble log(Mdouble Power)
bool isEqual(Mdouble v1, Mdouble v2, Mdouble absError)
Compares the difference of two Mdouble with an absolute error, useful in UnitTests.
Mdouble chi_squared_prob(Mdouble x, unsigned int k)
This is the function which actually gives the probability back using a chi squared test...
Mdouble cos(Mdouble x)
Definition: ExtendedMath.cc:64
constexpr T factorial(const T t)
factorial function
Definition: ExtendedMath.h:155
Mdouble sin(Mdouble x)
Definition: ExtendedMath.cc:44
Mdouble chi_squared(Mdouble x, unsigned int k)
This is a chi_squared function return the value x and degrees of freedom k.
const Mdouble pi
Definition: ExtendedMath.h:45
Mdouble I0_exp(Mdouble x)
T tan(T x)
Definition: ExtendedMath.h:178
Mdouble chebyshev(Mdouble x, const Mdouble coef[], int N)
Namespace for evaluating the zeroth modified Bessel function of the first kind, I0(x), required in StatisticsPoint.hcc.
const Mdouble sqrt_2
Definition: ExtendedMath.h:48
NumericalVector< std::complex< Mdouble > > sphericalHarmonics(int p, Mdouble theta, Mdouble phi)
Mdouble gamma(Mdouble gamma_in)
This is the gamma function returns the true value for the half integer value.
const Mdouble R
Definition: ExtendedMath.h:50
NumericalVector associatedLegendrePolynomials(int n, Mdouble x)
const Mdouble degree
Definition: ExtendedMath.h:52
Implementation of a 3D matrix.
Definition: Matrix.h:37
Definition: Vector.h:49
T square(const T val)
squares a number
Definition: ExtendedMath.h:106
const Mdouble sqrt_3
Definition: ExtendedMath.h:49
const Mdouble sqr_pi
Definition: ExtendedMath.h:47
Implementation of a 3D symmetric matrix.
NumericalVector computeSquaredFactorialValues(int p)