Commit 69bcab53 authored by Christoph Schwarzmeier's avatar Christoph Schwarzmeier

Merge branch 'rename_math-constants' into 'master'

[API] Rename mathematical constants to make them differ from C macros

See merge request !216
parents 57f2e7db 1a95afc8
Pipeline #16667 passed with stages
in 382 minutes and 32 seconds
......@@ -1229,7 +1229,7 @@ int main( int argc, char **argv )
const Vector3<uint_t> domainSize( XBlocks * blockSize, YBlocks * blockSize, ZBlocks * blockSize );
const auto domainVolume = real_t(domainSize[0] * domainSize[1] * domainSize[2]);
const real_t sphereVolume = math::M_PI / real_t(6) * diameter * diameter * diameter;
const real_t sphereVolume = math::pi / real_t(6) * diameter * diameter * diameter;
const uint_t numberOfSediments = uint_c(std::ceil(solidVolumeFraction * domainVolume / sphereVolume));
real_t expectedSedimentVolumeFraction = (useBox||useHopper) ? real_t(0.45) : real_t(0.52);
......
......@@ -463,8 +463,8 @@ int main( int argc, char **argv )
real_t domainHeight = real_c(ZCells) - topWallOffset;
real_t fluidVolume = real_c( XCells * YCells ) * domainHeight;
real_t solidVolume = solidVolumeFraction * fluidVolume;
uint_t numberOfParticles = uint_c(std::ceil(solidVolume / ( math::M_PI / real_t(6) * diameter * diameter * diameter )));
diameter = std::cbrt( solidVolume / ( real_c(numberOfParticles) * math::M_PI / real_t(6) ) );
uint_t numberOfParticles = uint_c(std::ceil(solidVolume / ( math::pi / real_t(6) * diameter * diameter * diameter )));
diameter = std::cbrt( solidVolume / ( real_c(numberOfParticles) * math::pi / real_t(6) ) );
auto densityRatio = real_t(2.5);
......@@ -821,7 +821,7 @@ int main( int argc, char **argv )
}
real_t sphereVolume = diameter * diameter * diameter * math::M_PI / real_t(6);
real_t sphereVolume = diameter * diameter * diameter * math::pi / real_t(6);
Vector3<real_t> gravitationalForce( real_t(0), real_t(0), -(densityRatio - real_t(1)) * gravitationalAcceleration * sphereVolume );
timeloop.addFuncAfterTimeStep(pe_coupling::ForceOnBodiesAdder( blocks, bodyStorageID, gravitationalForce ), "Gravitational force" );
......
......@@ -33,13 +33,13 @@ namespace dem {
real_t calcCoefficientOfRestitution(const real_t k, const real_t gamma, const real_t meff)
{
auto a = real_t(0.5) * gamma / meff;
return std::exp(-a * math::M_PI / std::sqrt(k / meff - a*a));
return std::exp(-a * math::pi / std::sqrt(k / meff - a*a));
}
real_t calcCollisionTime(const real_t k, const real_t gamma, const real_t meff)
{
auto a = real_t(0.5) * gamma / meff;
return math::M_PI / std::sqrt( k/meff - a*a);
return math::pi / std::sqrt( k/meff - a*a);
}
}
......
......@@ -675,7 +675,7 @@ int main( int argc, char **argv )
auto refinementTimestep = lbm::refinement::makeTimeStep< LatticeModel_T, BoundaryHandling_T >( blocks, sweep, pdfFieldID, boundaryHandlingID );
// add force evaluation and logging
real_t normalizationFactor = ( zeroShearTest ) ? real_t(1) : ( math::M_PI / real_t(8) * densityFluid * shearRate * shearRate * wallDistance * wallDistance * diameter * diameter );
real_t normalizationFactor = ( zeroShearTest ) ? real_t(1) : ( math::pi / real_t(8) * densityFluid * shearRate * shearRate * wallDistance * wallDistance * diameter * diameter );
std::string loggingFileName( baseFolderLogging + "/LoggingForcesNearPlane");
loggingFileName += "_lvl" + std::to_string(numberOfLevels);
loggingFileName += "_D" + std::to_string(uint_c(diameter));
......
......@@ -1080,7 +1080,7 @@ int main( int argc, char **argv )
WALBERLA_LOG_INFO_ON_ROOT("Initial simulation has ended.")
//evaluate the gravitational force necessary to keep the sphere at a approximately fixed position
gravity = forceEval->getForce() / ( (densityRatio - real_t(1) ) * diameter * diameter * diameter * math::M_PI / real_t(6) );
gravity = forceEval->getForce() / ( (densityRatio - real_t(1) ) * diameter * diameter * diameter * math::pi / real_t(6) );
GalileoSim = std::sqrt( ( densityRatio - real_t(1) ) * gravity * diameter * diameter * diameter ) / viscosity;
ReynoldsSim = uIn * diameter / viscosity;
u_ref = std::sqrt( std::fabs(densityRatio - real_t(1)) * gravity * diameter );
......@@ -1234,7 +1234,7 @@ int main( int argc, char **argv )
}
// add gravity
Vector3<real_t> extForcesOnSphere( real_t(0), real_t(0), - gravity * ( densityRatio - real_t(1) ) * diameter * diameter * diameter * math::M_PI / real_t(6));
Vector3<real_t> extForcesOnSphere( real_t(0), real_t(0), - gravity * ( densityRatio - real_t(1) ) * diameter * diameter * diameter * math::pi / real_t(6));
timeloop.addFuncAfterTimeStep( pe_coupling::ForceOnBodiesAdder( blocks, bodyStorageID, extForcesOnSphere ), "Add external forces (gravity and buoyancy)" );
// evaluate the sphere properties
......
......@@ -762,7 +762,7 @@ void run( const shared_ptr< Config > & config, const LatticeModel_T & latticeMod
{
setup.maxVelocity_L = ( setup.acceleration_L * setup.radius_L * setup.radius_L ) / ( real_t(4) * setup.viscosity_L );
setup.meanVelocity_L = ( setup.acceleration_L * setup.radius_L * setup.radius_L ) / ( real_t(8) * setup.viscosity_L );
setup.flowRate_L = setup.meanVelocity_L * math::M_PI * setup.radius_L * setup.radius_L;
setup.flowRate_L = setup.meanVelocity_L * math::pi * setup.radius_L * setup.radius_L;
}
else
{
......
......@@ -744,7 +744,7 @@ public:
void operator()( const real_t t )
{
tConstTerm_ = ( sinPeriod_ > real_t(0) ) ? ( std::abs( std::sin( math::M_PI * t / sinPeriod_ ) ) ) : real_t(1);
tConstTerm_ = ( sinPeriod_ > real_t(0) ) ? ( std::abs( std::sin( math::pi * t / sinPeriod_ ) ) ) : real_t(1);
tConstTerm_ *= uTerm_ * HTerm_;
tConstTerm_ *= ( raisingTime_ > real_t(0) ) ? std::min( t / raisingTime_, real_t(1) ) : real_t(1);
}
......
......@@ -371,7 +371,7 @@ uint_t createSpheresRandomly( StructuredBlockForest & forest, pe::BodyStorage &
pe::createSphere( globalBodyStorage, forest.getBlockStorage(), bodyStorageID, 0, Vector3<real_t>( xParticle, yParticle, zParticle ), creationDiameter * real_t(0.5), material );
currentSphereVolume += math::M_PI / real_t(6) * creationDiameter * creationDiameter * creationDiameter;
currentSphereVolume += math::pi / real_t(6) * creationDiameter * creationDiameter * creationDiameter;
++numberOfSpheres;
}
......@@ -851,16 +851,16 @@ int main( int argc, char **argv ) {
const real_t restitutionCoeff = real_t(0.88);
const real_t frictionCoeff = real_t(0.25);
real_t sphereVolume = diameterAvg * diameterAvg * diameterAvg * math::M_PI / real_t(6); // based on avg. diameter
real_t sphereVolume = diameterAvg * diameterAvg * diameterAvg * math::pi / real_t(6); // based on avg. diameter
const real_t particleMass = densityRatio * sphereVolume;
const real_t Mij = particleMass * particleMass / (real_t(2) * particleMass);
const real_t lnDryResCoeff = std::log(restitutionCoeff);
const real_t collisionTime = real_t(0.5);
const real_t stiffnessCoeff = math::M_PI * math::M_PI * Mij / (collisionTime * collisionTime *
(real_t(1) - lnDryResCoeff * lnDryResCoeff / (math::M_PI * math::M_PI + lnDryResCoeff * lnDryResCoeff)));
const real_t stiffnessCoeff = math::pi * math::pi * Mij / (collisionTime * collisionTime *
(real_t(1) - lnDryResCoeff * lnDryResCoeff / (math::pi * math::pi + lnDryResCoeff * lnDryResCoeff)));
const real_t dampingCoeff = -real_t(2) * std::sqrt(Mij * stiffnessCoeff) *
(std::log(restitutionCoeff) / std::sqrt(math::M_PI * math::M_PI + (std::log(restitutionCoeff) * std::log(restitutionCoeff))));
const real_t contactDuration = real_t(2) * math::M_PI * Mij / (std::sqrt(real_t(4) * Mij * stiffnessCoeff - dampingCoeff * dampingCoeff)); //formula from Uhlman
(std::log(restitutionCoeff) / std::sqrt(math::pi * math::pi + (std::log(restitutionCoeff) * std::log(restitutionCoeff))));
const real_t contactDuration = real_t(2) * math::pi * Mij / (std::sqrt(real_t(4) * Mij * stiffnessCoeff - dampingCoeff * dampingCoeff)); //formula from Uhlman
WALBERLA_LOG_INFO_ON_ROOT("Created particle material with:\n"
<< " - coefficient of restitution = " << restitutionCoeff << "\n"
......
......@@ -70,9 +70,9 @@ void initBC( const shared_ptr< StructuredBlockStorage > & blocks, const BlockDat
// obtain the physical coordinate of the center of the current cell
const Vector3< real_t > p = blocks->getBlockLocalCellCenter( *block, *cell );
// set the values of the Dirichlet boundary condition given by the function u(x,1) = sin(2*M_PI*x)*sinh(2*M_PI) in the source and destination field
src->get( *cell ) = std::sin( real_c(2) * math::M_PI * p[0] ) * std::sinh( real_c(2) * math::M_PI );
dst->get( *cell ) = std::sin( real_c(2) * math::M_PI * p[0] ) * std::sinh( real_c(2) * math::M_PI );
// set the values of the Dirichlet boundary condition given by the function u(x,1) = sin(2*PI*x)*sinh(2*PI) in the source and destination field
src->get( *cell ) = std::sin( real_c(2) * math::pi * p[0] ) * std::sinh( real_c(2) * math::pi );
dst->get( *cell ) = std::sin( real_c(2) * math::pi * p[0] ) * std::sinh( real_c(2) * math::pi );
}
}
}
......@@ -98,9 +98,9 @@ void initRHS( const shared_ptr< StructuredBlockStorage > & blocks, const BlockDa
// obtain the physical coordinate of the center of the current cell
const Vector3< real_t > p = blocks->getBlockLocalCellCenter( *block, *cell );
// set the right-hand side, given by the function f(x,y) = 4*M_PI*PI*sin(2*M_PI*x)*sinh(2*M_PI*y)
f->get( *cell ) = real_c(4) * math::M_PI * math::M_PI * std::sin( real_c(2) * math::M_PI * p[0] ) *
std::sinh( real_c(2) * math::M_PI * p[1] );
// set the right-hand side, given by the function f(x,y) = 4*PI*PI*sin(2*PI*x)*sinh(2*PI*y)
f->get( *cell ) = real_c(4) * math::pi * math::pi * std::sin( real_c(2) * math::pi * p[0] ) *
std::sinh( real_c(2) * math::pi * p[1] );
}
}
}
......@@ -147,7 +147,7 @@ void JacobiSweep::operator()( IBlock * const block )
dst->get(x,y,z) += ( real_c(1) / (dx_ * dx_) ) * src->get( x-1, y , z );
dst->get(x,y,z) += ( real_c(1) / (dy_ * dy_) ) * src->get( x , y+1, z );
dst->get(x,y,z) += ( real_c(1) / (dy_ * dy_) ) * src->get( x , y-1, z );
dst->get(x,y,z) /= ( real_c(2) / (dx_ * dx_) + real_c(2)/(dy_ * dy_) + real_c(4) * math::M_PI * math::M_PI );
dst->get(x,y,z) /= ( real_c(2) / (dx_ * dx_) + real_c(2)/(dy_ * dy_) + real_c(4) * math::pi * math::pi );
)
// swap source and destination fields
......@@ -282,7 +282,7 @@ int main( int argc, char ** argv )
// ...or the variant using the stencil concept
// set up the stencil weights
std::vector< real_t > weights( Stencil_T::Size );
weights[ Stencil_T::idx[ stencil::C ] ] = real_c(2) / ( dx * dx ) + real_c(2) / ( dy * dy ) + real_c(4) * math::M_PI * math::M_PI;
weights[ Stencil_T::idx[ stencil::C ] ] = real_c(2) / ( dx * dx ) + real_c(2) / ( dy * dy ) + real_c(4) * math::pi * math::pi;
weights[ Stencil_T::idx[ stencil::E ] ] = real_c(-1) / ( dx * dx );
weights[ Stencil_T::idx[ stencil::W ] ] = real_c(-1) / ( dx * dx );
weights[ Stencil_T::idx[ stencil::N ] ] = real_c(-1) / ( dy * dy );
......
......@@ -71,7 +71,7 @@ void JacobiSweep::operator()( IBlock * const block )
dst->get(x,y,z) += ( real_c(1) / (dx_ * dx_) ) * src->get( x-1, y , z );
dst->get(x,y,z) += ( real_c(1) / (dy_ * dy_) ) * src->get( x , y+1, z );
dst->get(x,y,z) += ( real_c(1) / (dy_ * dy_) ) * src->get( x , y-1, z );
dst->get(x,y,z) /= ( real_c(2) / (dx_ * dx_) + real_c(2)/(dy_ * dy_) + real_c(4) * math::M_PI * math::M_PI );
dst->get(x,y,z) /= ( real_c(2) / (dx_ * dx_) + real_c(2)/(dy_ * dy_) + real_c(4) * math::pi * math::pi );
)
// swap source and destination fields
......@@ -90,7 +90,7 @@ typedef stencil::D2Q5 Stencil_T;
// set up the stencil weights
std::vector< real_t > weights( Stencil_T::Size );
weights[ Stencil_T::idx[ stencil::C ] ] = real_c(2) / ( dx * dx ) + real_c(2) / ( dy * dy ) + real_c(4) * math::M_PI * math::M_PI;
weights[ Stencil_T::idx[ stencil::C ] ] = real_c(2) / ( dx * dx ) + real_c(2) / ( dy * dy ) + real_c(4) * math::pi * math::pi;
weights[ Stencil_T::idx[ stencil::E ] ] = real_c(-1) / ( dx * dx );
weights[ Stencil_T::idx[ stencil::W ] ] = real_c(-1) / ( dx * dx );
weights[ Stencil_T::idx[ stencil::N ] ] = real_c(-1) / ( dy * dy );
......@@ -185,9 +185,9 @@ void initRHS( const shared_ptr< StructuredBlockStorage > & blocks, const BlockDa
// obtain the physical coordinate of the center of the current cell
const Vector3< real_t > p = blocks->getBlockLocalCellCenter( *block, *cell );
// set the right-hand side, given by the function f(x,y) = 4*M_PI*PI*sin(2*M_PI*x)*sinh(2*M_PI*y)
f->get( *cell ) = real_c(4) * math::M_PI * math::M_PI * std::sin( real_c(2) * math::M_PI * p[0] ) *
std::sinh( real_c(2) * math::M_PI * p[1] );
// set the right-hand side, given by the function f(x,y) = 4*PI*PI*sin(2*PI*x)*sinh(2*PI*y)
f->get( *cell ) = real_c(4) * math::pi * math::pi * std::sin( real_c(2) * math::pi * p[0] ) *
std::sinh( real_c(2) * math::pi * p[1] );
}
}
}
......@@ -236,10 +236,10 @@ void initBC( const shared_ptr< StructuredBlockStorage > & blocks, const BlockDat
// obtain the physical coordinate of the center of the current cell
const Vector3< real_t > p = blocks->getBlockLocalCellCenter( *block, *cell );
// set the values of the Dirichlet boundary condition given by the function u(x,1) = sin(2*M_PI*x)*sinh(2*M_PI)
// set the values of the Dirichlet boundary condition given by the function u(x,1) = sin(2*PI*x)*sinh(2*PI)
// in the source and destination field
src->get( *cell ) = std::sin( real_c(2) * math::M_PI * p[0] ) * std::sinh( real_c(2) * math::M_PI );
dst->get( *cell ) = std::sin( real_c(2) * math::M_PI * p[0] ) * std::sinh( real_c(2) * math::M_PI );
src->get( *cell ) = std::sin( real_c(2) * math::pi * p[0] ) * std::sinh( real_c(2) * math::pi );
dst->get( *cell ) = std::sin( real_c(2) * math::pi * p[0] ) * std::sinh( real_c(2) * math::pi );
}
}
}
......
......@@ -66,8 +66,8 @@ void initU( const shared_ptr< StructuredBlockStorage > & blocks, const BlockData
// obtain the physical coordinate of the center of the current cell
const Vector3< real_t > p = blocks->getBlockLocalCellCenter( *block, *cell );
// set the initial condition, given by the function u(x,y,0) = sin(M_PI*x)*sin(M_PI*y)
u->get( *cell ) = std::sin( math::M_PI * p[0] ) * std::sin( math::M_PI * p[1] );
// set the initial condition, given by the function u(x,y,0) = sin(PI*x)*sin(PI*y)
u->get( *cell ) = std::sin( math::pi * p[0] ) * std::sin( math::pi * p[1] );
}
}
}
......
......@@ -67,8 +67,8 @@ void initU( const shared_ptr< StructuredBlockStorage > & blocks, const BlockData
// obtain the physical coordinate of the center of the current cell
const Vector3< real_t > p = blocks->getBlockLocalCellCenter( *block, *cell );
// set the initial condition, given by the function u(x,y,0) = sin(M_PI*x)*sin(M_PI*y)
u->get( *cell ) = std::sin( math::M_PI * p[0] ) * std::sin( math::M_PI * p[1] );
// set the initial condition, given by the function u(x,y,0) = sin(PI*x)*sin(PI*y)
u->get( *cell ) = std::sin( math::pi * p[0] ) * std::sin( math::pi * p[1] );
}
}
}
......
......@@ -90,7 +90,7 @@ public:
const real_t meff)
{
auto a = real_t(0.5) * getDampingN(type1, type2) / meff;
return std::exp(-a * math::M_PI / std::sqrt(getStiffness(type1, type2) / meff - a*a));
return std::exp(-a * math::pi / std::sqrt(getStiffness(type1, type2) / meff - a*a));
}
inline
......@@ -99,7 +99,7 @@ public:
const real_t meff)
{
auto a = real_t(0.5) * getDampingN(type1, type2) / meff;
return math::M_PI / std::sqrt( getStiffness(type1, type2)/meff - a*a);
return math::pi / std::sqrt( getStiffness(type1, type2)/meff - a*a);
}
inline
......@@ -110,8 +110,8 @@ public:
const real_t meff)
{
const real_t lnDryResCoeff = std::log(cor);
setStiffness(type1, type2, math::M_PI * math::M_PI * meff / ( collisionTime * collisionTime * ( real_t(1) - lnDryResCoeff * lnDryResCoeff / ( math::M_PI * math::M_PI + lnDryResCoeff* lnDryResCoeff )) ));
setDampingN( type1, type2, - real_t(2) * std::sqrt( meff * getStiffness(type1, type2) ) * ( lnDryResCoeff / std::sqrt( math::M_PI * math::M_PI + ( lnDryResCoeff * lnDryResCoeff ) ) ));
setStiffness(type1, type2, math::pi * math::pi * meff / ( collisionTime * collisionTime * ( real_t(1) - lnDryResCoeff * lnDryResCoeff / ( math::pi * math::pi + lnDryResCoeff* lnDryResCoeff )) ));
setDampingN( type1, type2, - real_t(2) * std::sqrt( meff * getStiffness(type1, type2) ) * ( lnDryResCoeff / std::sqrt( math::pi * math::pi + ( lnDryResCoeff * lnDryResCoeff ) ) ));
}
private:
uint_t numParticleTypes_;
......
......@@ -16,57 +16,43 @@
//! \file Constants.h
//! \author Klaus Iglberger
//! \author Sebastian Eibl <sebastian.eibl@fau.de>
//! \author Christoph Schwarzmeier <christoph.schwarzmeier@fau.de>
//! \brief Header file for mathematical constants
//
//======================================================================================================================
#pragma once
//*************************************************************************************************
// Includes
//*************************************************************************************************
#include <core/DataTypes.h>
#include <cmath>
#include <core/DataTypes.h>
// Disable false warnings in GCC 5
#if ( defined __GNUC__ ) && ( __GNUC__ == 5 ) && ( __GNUC_MINOR__ == 1 )
# pragma GCC diagnostic push
# pragma GCC diagnostic ignored "-Wunused-variable"
#if (defined __GNUC__) && (__GNUC__ == 5) && (__GNUC_MINOR__ == 1)
# pragma GCC diagnostic push
# pragma GCC diagnostic ignored "-Wunused-variable"
#endif
namespace walberla {
namespace math {
# undef M_E
# undef M_LOG2E
# undef M_LOG10E
# undef M_LN2
# undef M_LN10
# undef M_PI
# undef M_PI_2
# undef M_PI_4
# undef M_1_PI
# undef M_2_PI
# undef M_2_SQRTPI
# undef M_SQRT2
# undef M_SQRT1_2
constexpr real_t M_E = real_t( 2.7182818284590452354 ); /* e */
constexpr real_t M_LOG2E = real_t( 1.4426950408889634074 ); /* log_2 e */
constexpr real_t M_LOG10E = real_t( 0.43429448190325182765 ); /* log_10 e */
constexpr real_t M_LN2 = real_t( 0.69314718055994530942 ); /* log_e 2 */
constexpr real_t M_LN10 = real_t( 2.30258509299404568402 ); /* log_e 10 */
constexpr real_t M_PI = real_t( 3.14159265358979323846 ); /* pi */
constexpr real_t M_PI_2 = real_t( 1.57079632679489661923 ); /* pi/2 */
constexpr real_t M_PI_4 = real_t( 0.78539816339744830962 ); /* pi/4 */
constexpr real_t M_1_PI = real_t( 0.31830988618379067154 ); /* 1/pi */
constexpr real_t M_2_PI = real_t( 0.63661977236758134308 ); /* 2/pi */
constexpr real_t M_2_SQRTPI = real_t( 1.12837916709551257390 ); /* 2/sqrt(pi) */
constexpr real_t M_SQRT2 = real_t( 1.41421356237309504880 ); /* sqrt(2) */
constexpr real_t M_SQRT1_2 = real_t( 0.70710678118654752440 ); /* 1/sqrt(2) */
namespace walberla
{
namespace math
{
constexpr real_t e = real_t(2.7182818284590452354); /* e */
constexpr real_t log2_e = real_t(1.4426950408889634074); /* log_2 e */
constexpr real_t log10_e = real_t(0.43429448190325182765); /* log_10 e */
constexpr real_t ln_two = real_t(0.69314718055994530942); /* log_e 2 */
constexpr real_t ln_ten = real_t(2.30258509299404568402); /* log_e 10 */
constexpr real_t pi = real_t(3.14159265358979323846); /* pi */
constexpr real_t half_pi = real_t(1.57079632679489661923); /* pi/2 */
constexpr real_t fourth_pi = real_t(0.78539816339744830962); /* pi/4 */
constexpr real_t one_div_pi = real_t(0.31830988618379067154); /* 1/pi */
constexpr real_t two_div_pi = real_t(0.63661977236758134308); /* 2/pi */
constexpr real_t two_div_root_pi = real_t(1.12837916709551257390); /* 2/sqrt(pi) */
constexpr real_t root_two = real_t(1.41421356237309504880); /* sqrt(2) */
constexpr real_t one_div_root_two = real_t(0.70710678118654752440); /* 1/sqrt(2) */
} // namespace math
}
} // namespace walberla
\ No newline at end of file
......@@ -111,8 +111,8 @@ public:
inline bool containsClosedInterval( const vector_type & point ) const;
inline bool containsClosedInterval( const vector_type & point, const value_type dx ) const;
inline GenericAABB getExtended( const value_type e ) const;
inline GenericAABB getExtended( const vector_type & e ) const;
inline GenericAABB getExtended( const value_type eps ) const;
inline GenericAABB getExtended( const vector_type & eps ) const;
inline GenericAABB getTranslated( const vector_type & translation ) const;
......@@ -162,8 +162,8 @@ public:
inline void setAxisBounds( const uint_t index, const value_type value1, const value_type value2 );
inline void extend( const value_type e );
inline void extend( const vector_type & e );
inline void extend( const value_type eps );
inline void extend( const vector_type & eps );
inline void setCenter( const vector_type & center );
inline void translate( const vector_type & translation );
......
......@@ -674,19 +674,19 @@ bool GenericAABB< T >::containsClosedInterval( const vector_type & point, const
/**
* \brief Creates a new GenericAABB by extending this one
*
* \param e epsilon by which the bounding box is extended in each direction
* \param eps epsilon by which the bounding box is extended in each direction
*
* \returns The extended GenericAABB
*/
template< typename T >
GenericAABB< T > GenericAABB< T >::getExtended( const value_type e ) const
GenericAABB< T > GenericAABB< T >::getExtended( const value_type eps ) const
{
vector_type newMinCorner( minCorner_[0] - e, minCorner_[1] - e, minCorner_[2] - e );
vector_type newMinCorner( minCorner_[0] - eps, minCorner_[1] - eps, minCorner_[2] - eps );
return createFromMinMaxCorner( newMinCorner[0], newMinCorner[1], newMinCorner[2],
std::max( newMinCorner[0], maxCorner_[0] + e ),
std::max( newMinCorner[1], maxCorner_[1] + e ),
std::max( newMinCorner[2], maxCorner_[2] + e ) );
std::max( newMinCorner[0], maxCorner_[0] + eps ),
std::max( newMinCorner[1], maxCorner_[1] + eps ),
std::max( newMinCorner[2], maxCorner_[2] + eps ) );
}
......@@ -694,19 +694,19 @@ GenericAABB< T > GenericAABB< T >::getExtended( const value_type e ) const
/**
* \brief Creates a new GenericAABB by extending this one
*
* \param e epsilon vector by which the bounding box is extended. The box is extended in each direction by
* \param eps epsilon vector by which the bounding box is extended. The box is extended in each direction by
* the corresponding vector component.
*
* \returns The extended GenericAABB
*/
template< typename T >
GenericAABB< T > GenericAABB< T >::getExtended( const vector_type & e ) const
GenericAABB< T > GenericAABB< T >::getExtended( const vector_type & eps ) const
{
vector_type newMinCorner( minCorner_ - e );
vector_type newMinCorner( minCorner_ - eps );
return createFromMinMaxCorner( newMinCorner[0], newMinCorner[1], newMinCorner[2],
std::max( newMinCorner[0], maxCorner_[0] + e[0] ),
std::max( newMinCorner[1], maxCorner_[1] + e[0] ),
std::max( newMinCorner[2], maxCorner_[2] + e[0] ) );
std::max( newMinCorner[0], maxCorner_[0] + eps[0] ),
std::max( newMinCorner[1], maxCorner_[1] + eps[0] ),
std::max( newMinCorner[2], maxCorner_[2] + eps[0] ) );
}
......@@ -1519,18 +1519,18 @@ void GenericAABB< T >::setAxisBounds( const uint_t index, const value_type value
/**
* \brief Extends this GenericAABB
*
* \param e epsilon by which the bounding box is extended in each direction
* \param eps epsilon by which the bounding box is extended in each direction
*/
template< typename T >
void GenericAABB< T >::extend( const value_type e )
void GenericAABB< T >::extend( const value_type eps )
{
minCorner_[0] -= e;
minCorner_[1] -= e;
minCorner_[2] -= e;
minCorner_[0] -= eps;
minCorner_[1] -= eps;
minCorner_[2] -= eps;
maxCorner_[0] += e;
maxCorner_[1] += e;
maxCorner_[2] += e;
maxCorner_[0] += eps;
maxCorner_[1] += eps;
maxCorner_[2] += eps;
for( uint_t i = 0; i < 3; ++i )
if( minCorner_[i] > maxCorner_[i] )
......@@ -1544,14 +1544,14 @@ void GenericAABB< T >::extend( const value_type e )
/**
* \brief Extends this GenericAABB
*
* \param e epsilon vector by which the bounding box is extended. The box is extended in each direction by
* \param eps epsilon vector by which the bounding box is extended. The box is extended in each direction by
* the corresponding vector component.
*/
template< typename T >
void GenericAABB< T >::extend( const vector_type & e )
void GenericAABB< T >::extend( const vector_type & eps )
{
minCorner_ -= e;
maxCorner_ += e;
minCorner_ -= eps;
maxCorner_ += eps;
for( uint_t i = 0; i < 3; ++i )
if( minCorner_[i] > maxCorner_[i] )
......
......@@ -189,12 +189,12 @@ NodePtr EquationParser::parseFunction( const std::string& str, size_t& index ) c
{
case OP_FUNC_EXP:
funcPtr = std::make_shared<Node>( OP_PROD );
funcPtr->left() = std::make_shared<Node>( M_E );
funcPtr->left() = std::make_shared<Node>( math::e );
funcPtr->right() = nodePtr;
return funcPtr;
case OP_FUNC_LN:
funcPtr = std::make_shared<Node>( OP_LOG );
funcPtr->right() = std::make_shared<Node>( M_E );
funcPtr->right() = std::make_shared<Node>( math::e );
funcPtr->left() = nodePtr;
return funcPtr;
case OP_FUNC_SQRT:
......
......@@ -37,7 +37,7 @@ public:
void updateMassAndInertia(const real_t density) override;
real_t getVolume() const override { return (real_t(4) / real_t(3)) * math::M_PI * getRadius() * getRadius() * getRadius(); }
real_t getVolume() const override { return (real_t(4) / real_t(3)) * math::pi * getRadius() * getRadius() * getRadius(); }
static const int SHAPE_TYPE = 1; ///< Unique shape type identifier for spheres.\ingroup mesa_pd_shape
......@@ -48,7 +48,7 @@ private:
inline
void Sphere::updateMassAndInertia(const real_t density)
{
const real_t m = (real_c(4.0)/real_c(3.0) * math::M_PI) * getRadius() * getRadius() * getRadius() * density;
const real_t m = (real_c(4.0)/real_c(3.0) * math::pi) * getRadius() * getRadius() * getRadius() * density;
const Mat3 I = Mat3::makeDiagonalMatrix( real_c(0.4) * m * getRadius() * getRadius() );
mass_ = m;
......
......@@ -100,7 +100,7 @@ public:
const real_t meff)
{
auto a = real_t(0.5) * getDampingN(type1, type2) / meff;
return std::exp(-a * math::M_PI / std::sqrt(getStiffness(type1, type2) / meff - a*a));
return std::exp(-a * math::pi / std::sqrt(getStiffness(type1, type2) / meff - a*a));
}
inline
......@@ -109,7 +109,7 @@ public:
const real_t meff)
{
auto a = real_t(0.5) * getDampingN(type1, type2) / meff;
return math::M_PI / std::sqrt( getStiffness(type1, type2)/meff - a*a);
return math::pi / std::sqrt( getStiffness(type1, type2)/meff - a*a);
}
inline
......@@ -120,8 +120,8 @@ public:
const real_t meff)
{
const real_t lnDryResCoeff = std::log(cor);
setStiffness(type1, type2, math::M_PI * math::M_PI * meff / ( collisionTime * collisionTime * ( real_t(1) - lnDryResCoeff * lnDryResCoeff / ( math::M_PI * math::M_PI + lnDryResCoeff* lnDryResCoeff )) ));
setDampingN( type1, type2, - real_t(2) * std::sqrt( meff * getStiffness(type1, type2) ) * ( lnDryResCoeff / std::sqrt( math::M_PI * math::M_PI + ( lnDryResCoeff * lnDryResCoeff ) ) ));
setStiffness(type1, type2, math::pi * math::pi * meff / ( collisionTime * collisionTime * ( real_t(1) - lnDryResCoeff * lnDryResCoeff / ( math::pi * math::pi + lnDryResCoeff* lnDryResCoeff )) ));
setDampingN( type1, type2, - real_t(2) * std::sqrt( meff * getStiffness(type1, type2) ) * ( lnDryResCoeff / std::sqrt( math::pi * math::pi + ( lnDryResCoeff * lnDryResCoeff ) ) ));
}
private:
uint_t numParticleTypes_;
......
......@@ -536,7 +536,7 @@ inline void EPA::createInitialSimplex( size_t numPoints, GeomPrimitive &geom1, G
}
Vec3 direction1 = (d % axis).getNormalized();
Quat q(d, (real_t(2.0)/real_t(3.0)) * real_t(walberla::math::M_PI));
Quat q(d, (real_t(2.0)/real_t(3.0)) * real_t(walberla::math::pi));
Mat3 rot = q.toRotationMatrix();
Vec3 direction2 = (rot*direction1).getNormalized();
Vec3 direction3 = (rot*direction2).getNormalized();
......
......@@ -1203,8 +1203,8 @@ inline real_t HardContactSemiImplicitTimesteppingSolvers::relaxInelasticGenerali
alpha_left = -angleI - shiftI;
alpha_right = +angleI - shiftI;
if( alpha_left < 0 ) {
alpha_left += 2 * math::M_PI;
alpha_right += 2 * math::M_PI;
alpha_left += 2 * math::pi;
alpha_right += 2 * math::pi;
}
}
else if( contactCache.diag_nto_[i](0, 0) > contactCache.mu_[i] * a3 ) {
......@@ -1213,8 +1213,8 @@ inline real_t HardContactSemiImplicitTimesteppingSolvers::relaxInelasticGenerali
alpha_left = -angleJ - shiftJ;
alpha_right = +angleJ - shiftJ;
if( alpha_left < 0 ) {
alpha_left += 2 * math::M_PI;
alpha_right += 2 * math::M_PI;
alpha_left += 2 * math::pi;
alpha_right += 2 * math::pi;
}
}
else {
......@@ -1223,15 +1223,15 @@ inline real_t HardContactSemiImplicitTimesteppingSolvers::relaxInelasticGenerali
real_t alpha1_left( -angleJ - shiftJ );
real_t alpha1_right( +angleJ - shiftJ );
if( alpha1_left < 0 ) {
alpha1_left += 2 * math::M_PI;
alpha1_right += 2 * math::M_PI;
alpha1_left += 2 * math::pi;
alpha1_right += 2 * math::pi;
}
const real_t angleI( std::acos( fractionI ) );
real_t alpha2_left( -angleI - shiftI );
real_t alpha2_right( +angleI - shiftI );
if( alpha2_left < 0 ) {
alpha2_left += 2 * math::M_PI;
alpha2_right += 2 * math::M_PI;
alpha2_left += 2 * math::pi;
alpha2_right += 2 * math::pi;
}
// Swap intervals if second interval does not start right of the first interval.
......@@ -1241,7 +1241,7 @@ inline real_t HardContactSemiImplicitTimesteppingSolvers::relaxInelasticGenerali
}
if( alpha2_left > alpha1_right ) {
alpha2_right -= 2*math::M_PI;
alpha2_right -= 2*math::pi;
if( alpha2_right > alpha1_right ) {
// [alpha1_left; alpha1_right] \subset [alpha2_left; alpha2_right]
}
......@@ -1807,7 +1807,7 @@ inline void HardContactSemiImplicitTimesteppingSolvers::integratePositions( Body
v = v * (edge * getSpeedLimitFactor() / dt / speed );
}
const real_t maxPhi = real_t(2) * math::M_PI * getSpeedLimitFactor();
const real_t maxPhi = real_t(2) * math::pi * getSpeedLimitFactor();
const real_t phi = w.length() * dt;
if (phi > maxPhi)
{
......
......@@ -197,7 +197,7 @@ void Raytracer::setupView_() {
// viewing plane setup
d_ = (cameraPosition_ - lookAtPoint_).length();
aspectRatio_ = real_t(pixelsHorizontal_) / real_t(pixelsVertical_);
real_t fov_vertical_rad = fov_vertical_ * math::M_PI / real_t(180.0);
real_t fov_vertical_rad = fov_vertical_ * math::pi / real_t(180.0);
viewingPlaneHeight_ = real_c(tan(fov_vertical_rad/real_t(2.))) * real_t(2.) * d_;
viewingPlaneWidth_ = viewingPlaneHeight_ * aspectRatio_;
viewingPlaneOrigin_ = lookAtPoint_ - u_*viewingPlaneWidth_/real_t(2.) - v_*viewingPlaneHeight_/real_t(2.);
......
......@@ -194,8 +194,8 @@ void Capsule::calcBoundingBox()
*/
Mat3 Capsule::calcInertia( const real_t radius, const real_t length, const real_t density)
{
const real_t sphereMass( real_t (4)/real_t (3) * math::M_PI * radius*radius*radius * density );
const real_t cylinderMass( math::M_PI * radius*radius * length * density );
const real_t sphereMass( real_t (4)/real_t (3) * math::pi * radius*radius*radius * density );
const real_t cylinderMass( math::pi * radius*radius * length * density );
// 'Ia' represent the moment of inertia along the x-axis. 'Ia' contains the following two parts:
// - cylinder : I = (1/2)*mass*radius^2
......
......@@ -223,7 +223,7 @@ inline real_t Capsule::getVolume() const
*/
inline real_t Capsule::calcVolume( real_t radius, real_t length )