MercuryDPM  Trunk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
mathsFunc Namespace Reference

Namespace for some extra maths function that are often needed. More...

Functions

Mdouble gamma (Mdouble gamma_in)
 This is the gamma function returns the true value for the half integer value. More...
 
Mdouble beta (Mdouble z, Mdouble w)
 This is the beta function, returns the approximation based on cmath's implementation of ln(gamma) More...
 
Mdouble chi_squared (Mdouble x, unsigned int k)
 This is a chi_squared function return the value x and degrees of freedom k. More...
 
Mdouble chi_squared_prob (Mdouble x, unsigned int k)
 This is the function which actually gives the probability back using a chi squared test. More...
 
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. More...
 
template<typename T >
int sign (T val)
 This is a sign function, it returns -1 for negative numbers, 1 for positive numbers and 0 for 0. More...
 
template<typename T >
square (const T val)
 squares a number More...
 
template<typename T >
cubic (const T val)
 calculates the cube of a number More...
 
bool isEqual (Mdouble v1, Mdouble v2, Mdouble absError)
 Compares the difference of two Mdouble with an absolute error, useful in UnitTests. More...
 
bool isEqual (Vec3D v1, Vec3D v2, Mdouble absError)
 Compares the difference of two Vec3D with an absolute error, useful in UnitTests. More...
 
bool isEqual (Matrix3D m1, Matrix3D m2, Mdouble absError)
 Compares the difference of two Vec3D with an absolute error, useful in UnitTests. More...
 
bool isEqual (MatrixSymmetric3D m1, MatrixSymmetric3D m2, Mdouble absError)
 
bool isEqual (Quaternion v1, Quaternion v2, double absError)
 
template<typename T >
constexpr T factorial (const T t)
 factorial function More...
 
Mdouble sin (Mdouble x)
 
Mdouble cos (Mdouble x)
 
Mdouble exp (Mdouble Exponent)
 
Mdouble log (Mdouble Power)
 
template<typename T >
tan (T x)
 
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. More...
 
Mdouble I0_exp (Mdouble x)
 
Mdouble I0 (Mdouble x)
 

Detailed Description

Namespace for some extra maths function that are often needed.

Function Documentation

Mdouble mathsFunc::beta ( Mdouble  z,
Mdouble  w 
)

This is the beta function, returns the approximation based on cmath's implementation of ln(gamma)

Computes the beta function for real numbers, based on cmath's approximation of the ln of the gamma function. See https://en.wikipedia.org/wiki/Beta_function for more details

Parameters
zfirst Mdouble argument to compute the beta function of
wsecond Mdouble argument to compute the beta function of
Returns
the value of beta(z,w) as an Mdouble

Definition at line 164 of file ExtendedMath.cc.

References exp().

Referenced by SuperQuadricParticle::computeMass(), Multipole::convertMultipoleToLocal(), InfiniteWall::getFurthestPointSuperQuadric(), SuperQuadricParticle::getVolume(), ParhamiMcMeekingSinterSpecies::set(), SuperQuadricParticle::setBoundingRadius(), SuperQuadricParticle::setInertia(), LocalExpansion::translateLocalExpansion(), and Multipole::TranslateMultipoleExpansionTo().

165 {
166  return std::exp(std::lgamma(z) + std::lgamma(w) - std::lgamma(z + w));
167 }
Mdouble exp(Mdouble Exponent)
Definition: ExtendedMath.cc:84
Mdouble mathsFunc::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.

Definition at line 293 of file ExtendedMath.cc.

References constants::i, and logger.

Referenced by I0_exp().

294 {
295  const Mdouble* p = coef;
296  Mdouble b0 = *p++;
297  Mdouble b1 = 0, b2;
298  int i = N - 1;
299 
300  logger.assert(i > 0, "i is greater than 0");
301  do
302  {
303  b2 = b1;
304  b1 = b0;
305  b0 = x * b1 - b2 + *p++;
306  } while (--i);
307 
308  return (0.5 * (b0 - b2));
309 }
Logger< MERCURY_LOGLEVEL > logger("MercuryKernel")
double Mdouble
Definition: GeneralDefine.h:34
const std::complex< Mdouble > i
Definition: ExtendedMath.h:50
Mdouble mathsFunc::chi_squared ( Mdouble  x,
unsigned int  k 
)

This is a chi_squared function return the value x and degrees of freedom k.

Definition at line 172 of file ExtendedMath.cc.

References exp(), and gamma().

Referenced by chi_squared_prob().

173 {
174 
175  Mdouble prefactor = pow(2, k / 2.0) * gamma(k / 2.0);
176  Mdouble mainfactor = pow(x, k / 2.0 - 1) * exp(x / -2.0);
177 
178  return mainfactor / prefactor;
179 
180 }
double Mdouble
Definition: GeneralDefine.h:34
Mdouble exp(Mdouble Exponent)
Definition: ExtendedMath.cc:84
Mdouble gamma(Mdouble gamma_in)
This is the gamma function returns the true value for the half integer value.
Mdouble mathsFunc::chi_squared_prob ( Mdouble  x,
unsigned int  k 
)

This is the function which actually gives the probability back using a chi squared test.

This calculates the probability based on a chi squared test First we calculated the cumulative chi_squared function.

This is the function which actually gives the probability back It is calculated by calling the normal chi_squared function and using the trapezoidal rule. The final results is 1-the cumulative chi_squared function

Definition at line 188 of file ExtendedMath.cc.

References chi_squared(), and constants::i.

Referenced by RNG::test().

189 {
190 
191 //The current value was picked by tried were it stopped effect the 4 d.p.
192  const int num_steps_per_unit = 100;
193  Mdouble sum = 0;
194  Mdouble x = 0;
195  long int num_steps = static_cast<int>(num_steps_per_unit * x_max);
196 //Use trapezional rule, but ignoring the ends
197  for (int i = 0; i < num_steps; i++)
198  {
199  x = x_max / num_steps * (i + 0.5);
200  sum = sum + chi_squared(x, k);
201  }
202  return 1.0 - sum * x_max / num_steps;
203 
204 }
double Mdouble
Definition: GeneralDefine.h:34
const std::complex< Mdouble > i
Definition: ExtendedMath.h:50
Mdouble chi_squared(Mdouble x, unsigned int k)
This is a chi_squared function return the value x and degrees of freedom k.
Mdouble mathsFunc::cos ( Mdouble  x)

Definition at line 64 of file ExtendedMath.cc.

References constants::i, and constants::pi.

Referenced by ChuteWithHopper::addHopper(), LiquidBridgeWilletInteraction::computeAdhesionForce(), LiquidMigrationWilletInteraction::computeAdhesionForce(), BaseCluster::computeInternalStructure(), CircularPeriodicBoundary::createPeriodicParticle(), LevelSetWall::createVTKSphere(), VChute::getDistanceAndNormal(), SineWall::getDistanceAndNormal(), HorizontalScrew::getDistanceAndNormal(), HorizontalBaseScrew::getDistanceAndNormal(), Coil::getDistanceAndNormal(), InfiniteWallWithHole::getDistanceAndNormal(), Screw::getDistanceAndNormalLabCoordinates(), Vec3D::getFromCylindricalCoordinates(), RNG::getNormalVariate(), helpers::objectivenessTest(), BaseCluster::particleInsertionSuccessful(), HopperInsertionBoundary::placeParticle(), CircularPeriodicBoundary::rotateParticle(), Quaternion::setAngleZ(), Chute::setChuteAngleAndMagnitudeOfGravity(), Quaternion::setEuler(), sphericalHarmonics::sphericalHarmonics(), tan(), HorizontalBaseScrew::writeVTK(), HorizontalScrew::writeVTK(), Screw::writeVTK(), AxisymmetricIntersectionOfWalls::writeVTK(), and ScrewsymmetricIntersectionOfWalls::writeVTK().

65 {
66  Mdouble N = floor(x / (2.0 * constants::pi) + 0.5);
67  x -= N * 2.0 * constants::pi;
68  Mdouble Sum = 1.0;
69  Mdouble Power = 1;
70  Mdouble Sign = 1;
71  const Mdouble x2 = x * x;
72  Mdouble Fact = 1.0;
73  for (unsigned int i = 2; i < 25; i += 2)
74  {
75  Power *= x2;
76  Fact *= i * (i - 1);
77  Sign *= -1.0;
78  Sum += Sign * Power / Fact;
79  }
80  return Sum;
81 }
double Mdouble
Definition: GeneralDefine.h:34
const std::complex< Mdouble > i
Definition: ExtendedMath.h:50
const Mdouble pi
Definition: ExtendedMath.h:45
template<typename T >
T mathsFunc::cubic ( const T  val)
template<typename T >
constexpr T mathsFunc::factorial ( const T  t)

factorial function

Definition at line 153 of file ExtendedMath.h.

Referenced by sphericalHarmonics::computeSquaredFactorialValues(), HGridOptimiser::histNumberParticlesPerCell(), and sphericalHarmonics::sphericalHarmonics().

154 {
155  return (t == 0) ? 1 : t * factorial(t - 1);
156 }
constexpr T factorial(const T t)
factorial function
Definition: ExtendedMath.h:153
Mdouble mathsFunc::gamma ( Mdouble  gamma_in)

This is the gamma function returns the true value for the half integer value.

This is the gamma function, gives 'exact' answers for the half integer values This is done using the recussion relation and the known values for 1 and 0.5 Note, return NaN for non-half integer values.

Definition at line 137 of file ExtendedMath.cc.

References constants::sqrt_pi.

Referenced by chi_squared(), InfiniteWall::getFurthestPointSuperQuadric(), HopperInsertionBoundary::placeParticle(), and SuperQuadricParticle::setBoundingRadius().

138 {
139  const Mdouble ep = 1e-5;
140 
141  if (gamma_in > 1.0 + ep)
142  {
143  return ((gamma_in - 1) * gamma(gamma_in - 1));
144  }
145  else
146  {
147 
148  if ((gamma_in - ep < 1.0) && (gamma_in + ep > 1.0))
149  return 1.0;
150  else if ((gamma_in - ep < 0.5) && (gamma_in + ep > 0.5))
151  return constants::sqrt_pi;
152  else
153  return std::numeric_limits<Mdouble>::quiet_NaN();
154  }
155 } //end func gamma
double Mdouble
Definition: GeneralDefine.h:34
const Mdouble sqrt_pi
Definition: ExtendedMath.h:46
Mdouble gamma(Mdouble gamma_in)
This is the gamma function returns the true value for the half integer value.
Mdouble mathsFunc::goldenSectionSearch ( Mdouble(*)(const Mdouble function,
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.

Parameters
[in]functionA function pointer to the function of which you want to calculate the location of its minimum.
[in]minThe minimum location
[in]curThe current location
[in]maxThe maximum location
[in]endConditionThe algorithm terminates when abs(max - min) < endCondition
[in]curValThe value of the function at the current location (on default this value is calculated internally)

Definition at line 206 of file ExtendedMath.cc.

208 {
209  if (std::abs(max - min) < endCondition)
210  {
211  return 0.5 * (min + max);
212  }
213  std::cout << "Min=" << min << " Max=" << max << " diff=" << max - min << std::endl;
214  Mdouble resphi = 2 - 0.5 * (1 + std::sqrt(5));
215  Mdouble x;
216  if (max - cur > cur - min)
217  {
218  x = cur + resphi * (max - cur);
219  }
220  else
221  {
222  x = cur - resphi * (cur - min);
223  }
224  if (std::isnan(curVal))
225  curVal = function(cur);
226  Mdouble xVal = function(x);
227  if (xVal < curVal)
228  {
229  if (max - cur > cur - min)
230  {
231  return goldenSectionSearch(function, cur, x, max, endCondition, xVal);
232  }
233  else
234  {
235  return goldenSectionSearch(function, min, x, cur, endCondition, xVal);
236  }
237  }
238  else
239  {
240  if (max - cur > cur - min)
241  {
242  return goldenSectionSearch(function, min, cur, x, endCondition, curVal);
243  }
244  else
245  {
246  return goldenSectionSearch(function, x, cur, max, endCondition, curVal);
247  }
248  }
249 }
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 mathsFunc::I0 ( Mdouble  x)

Definition at line 391 of file ExtendedMath.cc.

References exp(), and I0_exp().

392 {
393  if (x < 0)
394  x = -x;
395  return exp(x) * I0_exp(x);
396 }
Mdouble exp(Mdouble Exponent)
Definition: ExtendedMath.cc:84
Mdouble I0_exp(Mdouble x)
Mdouble mathsFunc::I0_exp ( Mdouble  x)

Definition at line 311 of file ExtendedMath.cc.

References A, and chebyshev().

Referenced by I0().

312 {
313  // Coefficients for [0..8]
314  const Mdouble A[] =
315  {
316  -4.415341646479339379501E-18,
317  3.330794518822238097831E-17,
318  -2.431279846547954693591E-16,
319  1.715391285555133030611E-15,
320  -1.168533287799345168081E-14,
321  7.676185498604935616881E-14,
322  -4.856446783111929460901E-13,
323  2.955052663129639834611E-12,
324  -1.726826291441555707231E-11,
325  9.675809035373236912241E-11,
326  -5.189795601635262906661E-10,
327  2.659823724682386650351E-9,
328  -1.300025009986248042121E-8,
329  6.046995022541918949321E-8,
330  -2.670793853940611733911E-7,
331  1.117387539120103718151E-6,
332  -4.416738358458750563591E-6,
333  1.644844807072889708931E-5,
334  -5.754195010082103703981E-5,
335  1.885028850958416557291E-4,
336  -5.763755745385823658851E-4,
337  1.639475616941335798421E-3,
338  -4.324309995050575944301E-3,
339  1.054646039459499831831E-2,
340  -2.373741480589946881561E-2,
341  4.930528423967070848781E-2,
342  -9.490109704804764442101E-2,
343  1.716209015222087753491E-1,
344  -3.046826723431983986831E-1,
345  6.767952744094760849951E-1
346  };
347 
348  // Coefficients for [8..infinity]
349  const Mdouble B[] =
350  {
351  -7.233180487874753954561E-18,
352  -4.830504485944182071261E-18,
353  4.465621420296759999011E-17,
354  3.461222867697461093101E-17,
355  -2.827623980516583484941E-16,
356  -3.425485619677219134621E-16,
357  1.772560133056526383601E-15,
358  3.811680669352622420751E-15,
359  -9.554846698828307648701E-15,
360  -4.150569347287222086631E-14,
361  1.540086217521409826911E-14,
362  3.852778382742142701141E-13,
363  7.180124451383666233671E-13,
364  -1.794178531506806117781E-12,
365  -1.321581184044771311881E-11,
366  -3.149916527963241364541E-11,
367  1.188914710784643834241E-11,
368  4.940602388224969589101E-10,
369  3.396232025708386345151E-9,
370  2.266668990498178064591E-8,
371  2.048918589469063741831E-7,
372  2.891370520834756482971E-6,
373  6.889758346916823984261E-5,
374  3.369116478255694089901E-3,
375  8.044904110141088316081E-1
376  };
377 
378  if (x < 0)
379  x = -x;
380 
381  if (x <= 8.0)
382  {
383  Mdouble y = (x / 2.0) - 2.0;
384  return (chebyshev(y, A, 30));
385  }
386 
387  return (chebyshev(32.0 / x - 2.0, B, 25) / sqrt(x));
388 
389 }
double Mdouble
Definition: GeneralDefine.h:34
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.
bool mathsFunc::isEqual ( Mdouble  v1,
Mdouble  v2,
Mdouble  absError 
)

Compares the difference of two Mdouble with an absolute error, useful in UnitTests.

Parameters
[in]v1The first Mdouble
[in]v2The second Mdouble
[in]absErrorThe allowed maximum absolute error
Returns
True if the two Mdouble are equal

Definition at line 251 of file ExtendedMath.cc.

Referenced by checkTemplate(), isEqual(), ChuteInsertionBoundary::placeParticle(), SuperQuadricParticle::setBoundingRadius(), and MercuryBase::setHGridCellOverSizeRatio().

252 {
253  return std::abs(v1 - v2) <= absError;
254 }
bool mathsFunc::isEqual ( Vec3D  v1,
Vec3D  v2,
Mdouble  absError 
)

Compares the difference of two Vec3D with an absolute error, useful in UnitTests.

Parameters
[in]v1The first Vec3D
[in]v2The second Vec3D
[in]absErrorThe allowed maximum absolute error
Returns
true if the two Vec3D are equal

Definition at line 257 of file ExtendedMath.cc.

References isEqual(), Vec3D::X, Vec3D::Y, and Vec3D::Z.

258 {
259  return isEqual(v1.X, v2.X, absError) && isEqual(v1.Y, v2.Y, absError) && isEqual(v1.Z, v2.Z, absError);
260 }
Mdouble X
the vector components
Definition: Vector.h:65
bool isEqual(Mdouble v1, Mdouble v2, Mdouble absError)
Compares the difference of two Mdouble with an absolute error, useful in UnitTests.
Mdouble Y
Definition: Vector.h:65
Mdouble Z
Definition: Vector.h:65
bool mathsFunc::isEqual ( Matrix3D  m1,
Matrix3D  m2,
Mdouble  absError 
)

Compares the difference of two Vec3D with an absolute error, useful in UnitTests.

Parameters
[in]m1The first Vec3D
[in]m2The second Vec3D
[in]absErrorThe allowed maximum absolute error
Returns
true if the two Vec3D are equal

Definition at line 262 of file ExtendedMath.cc.

References isEqual(), Matrix3D::XX, Matrix3D::XY, Matrix3D::XZ, Matrix3D::YX, Matrix3D::YY, Matrix3D::YZ, Matrix3D::ZX, Matrix3D::ZY, and Matrix3D::ZZ.

263 {
264  return (isEqual(m1.XX, m2.XX, absError)
265  && isEqual(m1.XY, m2.XY, absError)
266  && isEqual(m1.XZ, m2.XZ, absError)
267  && isEqual(m1.YX, m2.YX, absError)
268  && isEqual(m1.YY, m2.YY, absError)
269  && isEqual(m1.YZ, m2.YZ, absError)
270  && isEqual(m1.ZX, m2.ZX, absError)
271  && isEqual(m1.ZY, m2.ZY, absError)
272  && isEqual(m1.ZZ, m2.ZZ, absError));
273 }
Mdouble XY
Definition: Matrix.h:43
Mdouble ZY
Definition: Matrix.h:43
Mdouble XZ
Definition: Matrix.h:43
Mdouble ZX
Definition: Matrix.h:43
Mdouble ZZ
Definition: Matrix.h:43
bool isEqual(Mdouble v1, Mdouble v2, Mdouble absError)
Compares the difference of two Mdouble with an absolute error, useful in UnitTests.
Mdouble YZ
Definition: Matrix.h:43
Mdouble YX
Definition: Matrix.h:43
Mdouble XX
all nine matrix elements
Definition: Matrix.h:43
Mdouble YY
Definition: Matrix.h:43
bool mathsFunc::isEqual ( MatrixSymmetric3D  m1,
MatrixSymmetric3D  m2,
Mdouble  absError 
)

Definition at line 275 of file ExtendedMath.cc.

References isEqual(), MatrixSymmetric3D::XX, MatrixSymmetric3D::XY, MatrixSymmetric3D::XZ, MatrixSymmetric3D::YY, MatrixSymmetric3D::YZ, and MatrixSymmetric3D::ZZ.

276 {
277  return (isEqual(m1.XX, m2.XX, absError)
278  && isEqual(m1.XY, m2.XY, absError)
279  && isEqual(m1.XZ, m2.XZ, absError)
280  && isEqual(m1.YY, m2.YY, absError)
281  && isEqual(m1.YZ, m2.YZ, absError)
282  && isEqual(m1.ZZ, m2.ZZ, absError));
283 }
bool isEqual(Mdouble v1, Mdouble v2, Mdouble absError)
Compares the difference of two Mdouble with an absolute error, useful in UnitTests.
Mdouble XX
The six distinctive matrix elements.
bool mathsFunc::isEqual ( Quaternion  v1,
Quaternion  v2,
double  absError 
)

Definition at line 285 of file ExtendedMath.cc.

References Quaternion::getComponent(), and isEqual().

286 {
287  return isEqual(v1.getComponent(0), v2.getComponent(0), absError) &&
288  isEqual(v1.getComponent(1), v2.getComponent(1), absError) &&
289  isEqual(v1.getComponent(2), v2.getComponent(2), absError) &&
290  isEqual(v1.getComponent(3), v2.getComponent(3), absError);
291 }
Mdouble getComponent(int index) const
Returns the requested component of this Quaternion.
Definition: Quaternion.cc:232
bool isEqual(Mdouble v1, Mdouble v2, Mdouble absError)
Compares the difference of two Mdouble with an absolute error, useful in UnitTests.
Mdouble mathsFunc::log ( Mdouble  Power)
Todo:
check if this function works

Definition at line 104 of file ExtendedMath.cc.

References A, exp(), and R.

Referenced by CircularPeriodicBoundary::checkBoundaryAfterParticleMoved(), helpers::computeDisptFromCollisionTimeAndRestitutionCoefficientAndTangentialRestitutionCoefficientAndEffectiveMass(), CircularPeriodicBoundary::createPeriodicParticle(), RNG::getNormalVariate(), SlidingFrictionSpecies::setCollisionTimeAndNormalAndTangentialRestitutionCoefficient(), LinearPlasticViscoelasticNormalSpecies::setCollisionTimeAndRestitutionCoefficient(), LinearViscoelasticNormalSpecies::setCollisionTimeAndRestitutionCoefficient(), HertzianViscoelasticNormalSpecies::setEffectiveElasticModulusAndRestitutionCoefficient(), LinearViscoelasticNormalSpecies::setRestitutionCoefficient(), LinearPlasticViscoelasticNormalSpecies::setRestitutionCoefficient(), and SinterNormalSpecies::setStiffnessAndRestitutionCoefficient().

105 {
106  Mdouble N, P, L, R, A, E;
107  E = 2.71828182845905;
108  P = Power;
109  N = 0.0;
110 
111  // This speeds up the convergence by calculating the integral
112  while (P >= E)
113  {
114  P /= E;
115  N++;
116  }
117  N += (P / E);
118  P = Power;
119  do
120  {
121  A = N;
122  L = (P / (exp(N - 1.0)));
123  R = ((N - 1.0) * E);
124  N = ((L + R) / E);
125  } while (N < A);
126  //} while (N != A);
127 
128  return N;
129 }
double Mdouble
Definition: GeneralDefine.h:34
Mdouble exp(Mdouble Exponent)
Definition: ExtendedMath.cc:84
template<typename T >
int mathsFunc::sign ( val)
Mdouble mathsFunc::sin ( Mdouble  x)

Definition at line 44 of file ExtendedMath.cc.

References constants::i, and constants::pi.

Referenced by ChuteWithHopper::addHopper(), BaseCluster::computeInternalStructure(), LevelSetWall::createVTKSphere(), VChute::getDistanceAndNormal(), SineWall::getDistanceAndNormal(), HorizontalScrew::getDistanceAndNormal(), HorizontalBaseScrew::getDistanceAndNormal(), Coil::getDistanceAndNormal(), InfiniteWallWithHole::getDistanceAndNormal(), Screw::getDistanceAndNormalLabCoordinates(), Vec3D::getFromCylindricalCoordinates(), ChuteWithHopper::getMaximumVelocityInducedByGravity(), RNG::getNormalVariate(), helpers::objectivenessTest(), BaseCluster::particleInsertionSuccessful(), HopperInsertionBoundary::placeParticle(), CircularPeriodicBoundary::rotateParticle(), Quaternion::setAngleZ(), Chute::setChuteAngleAndMagnitudeOfGravity(), Quaternion::setEuler(), tan(), HorizontalBaseScrew::writeVTK(), HorizontalScrew::writeVTK(), Screw::writeVTK(), AxisymmetricIntersectionOfWalls::writeVTK(), and ScrewsymmetricIntersectionOfWalls::writeVTK().

45 {
46  Mdouble N = floor(x / (2.0 * constants::pi) + 0.5);
47  x -= N * 2.0 * constants::pi;
48  Mdouble Sum = 0;
49  Mdouble Power = x;
50  Mdouble Sign = 1;
51  const Mdouble x2 = x * x;
52  Mdouble Fact = 1.0;
53  for (unsigned int i = 1; i < 25; i += 2)
54  {
55  Sum += Sign * Power / Fact;
56  Power *= x2;
57  Fact *= (i + 1) * (i + 2);
58  Sign *= -1.0;
59  }
60  return Sum;
61 }
double Mdouble
Definition: GeneralDefine.h:34
const std::complex< Mdouble > i
Definition: ExtendedMath.h:50
const Mdouble pi
Definition: ExtendedMath.h:45
template<typename T >
T mathsFunc::square ( const T  val)

squares a number

Definition at line 104 of file ExtendedMath.h.

Referenced by DPMBase::checkParticleForInteractionLocal(), helpers::computeDisptFromCollisionTimeAndRestitutionCoefficientAndTangentialRestitutionCoefficientAndEffectiveMass(), SlidingFrictionInteraction::computeFrictionForce(), MindlinRollingTorsionInteraction::computeFrictionForce(), FrictionInteraction::computeFrictionForce(), BaseParticle::computeMass(), SinterInteraction::computeNormalForce(), HertzianSinterInteraction::computeSinterForce(), LinearPlasticViscoelasticNormalSpecies::computeTimeStep(), SinterNormalSpecies::computeTimeStep(), PSD::convertProbabilityDensityNumberDistributionToProbabilityDensityVolumeDistribution(), PSD::convertProbabilityDensityToProbabilityDensityNumberDistribution(), Matrix3D::deviator(), HertzianViscoelasticNormalSpecies::getCollisionTime(), LinearPlasticViscoelasticNormalSpecies::getCollisionTime(), SinterNormalSpecies::getCollisionTime(), LinearViscoelasticNormalSpecies::getCollisionTime(), StatisticsVector< T >::getCutoff2(), ParabolaChute::getDistanceAndNormal(), SineWall::getDistanceAndNormal(), HorizontalScrew::getDistanceAndNormal(), Coil::getDistanceAndNormal(), Screw::getDistanceAndNormalLabCoordinates(), CGCoordinates::XY::getDistanceSquared(), CGCoordinates::Y::getDistanceSquared(), CGCoordinates::YZ::getDistanceSquared(), CGCoordinates::Z::getDistanceSquared(), CGCoordinates::RZ::getDistanceSquared(), CGCoordinates::X::getDistanceSquared(), CGCoordinates::XZ::getDistanceSquared(), CGCoordinates::R::getDistanceSquared(), ChargedBondedInteraction::getElasticEnergy(), SinterInteraction::getElasticEnergy(), LinearViscoelasticInteraction::getElasticEnergy(), HertzianSinterInteraction::getElasticEnergy(), LinearPlasticViscoelasticInteraction::getElasticEnergy(), HertzianViscoelasticInteraction::getElasticEnergy(), IrreversibleAdhesiveInteraction::getElasticEnergy(), HertzianViscoelasticInteraction::getElasticEnergyAtEquilibrium(), CGCoordinates::Base_X_Y_Z::getGaussIntegralPrefactor(), CGCoordinates::Base_XY_XZ_YZ::getGaussIntegralPrefactor(), CGCoordinates::XYZ::getGaussIntegralPrefactor(), CGCoordinates::Base_XY_XZ_YZ::getGaussPrefactor(), CGCoordinates::XYZ::getGaussPrefactor(), helpers::getMaximumVelocity(), CGFields::GradVelocityField::getSquared(), CGFields::LiquidMigrationFields::getSquared(), CGFields::StandardFields::getSquared(), CGCoordinates::YZ::getTangentialSquared(), CGCoordinates::XY::getTangentialSquared(), CGCoordinates::RZ::getTangentialSquared(), CGCoordinates::XZ::getTangentialSquared(), CGCoordinates::XYZ::getTangentialSquared(), Mercury3D::hGridFindContactsWithTargetCell(), MercuryBase::hGridUpdateMove(), BaseParticle::isInContactWith(), CGCoordinates::Base_XY_XZ_YZ::normalisePolynomialCoefficients(), SuperQuadricParticle::setBoundingRadius(), StatisticsVector< T >::setCGWidth(), SlidingFrictionSpecies::setCollisionTimeAndNormalAndTangentialRestitutionCoefficient(), SlidingFrictionSpecies::setCollisionTimeAndNormalAndTangentialRestitutionCoefficientNoDispt(), LinearPlasticViscoelasticNormalSpecies::setCollisionTimeAndRestitutionCoefficient(), LinearViscoelasticNormalSpecies::setCollisionTimeAndRestitutionCoefficient(), LinearViscoelasticNormalSpecies::setRestitutionCoefficient(), LinearPlasticViscoelasticNormalSpecies::setRestitutionCoefficient(), SinterNormalSpecies::setStiffnessAndRestitutionCoefficient(), MatrixSymmetric3D::square(), and Matrix3D::square().

105 {
106  return val * val;
107 }
template<typename T >
T mathsFunc::tan ( x)
Todo:
should be properly computed

Definition at line 176 of file ExtendedMath.h.

References cos(), and sin().

Referenced by ChuteWithHopper::addHopper(), VChute::getDistanceAndNormal(), and ChuteWithHopper::setHopper().

177 {
178  return sin(x) / cos(x);
179 }
Mdouble cos(Mdouble x)
Definition: ExtendedMath.cc:64
Mdouble sin(Mdouble x)
Definition: ExtendedMath.cc:44