BSMPT 3.0.7
BSMPT - Beyond the Standard Model Phase Transitions : A C++ package for the computation of the EWPT in BSM models
Loading...
Searching...
No Matches
BSMPT Namespace Reference

This classes calculates the Bounce action of the potential with a set temperature. More...

Namespaces

namespace  ModelID
 operator << overload for the model parameter
 

Classes

class  BounceActionInt
 
class  BounceSolution
 BounceSolution class that handles the calculation of the bounce solution as well as the calculation of the charateristic temperature scales. More...
 
class  BSMPTLogger
 The BSMPTLogger class. More...
 
class  Class_Potential_Origin
 The Class_Potential_Origin class Base class for all models. This class contains all numerical calculations on the tensors and the inherited classes only have to set them. More...
 
struct  CoexPhases
 CoexPhases struct to save pair of coexisting phases (false and true phase) More...
 
class  GravitationalWave
 
struct  GravitationalWaveData
 struct to store all calculated GW data More...
 
struct  gw_data
 gravitational wave data struct More...
 
struct  ISMConstants
 The ISMConstants struct containing all necessary SM constants. More...
 
class  Logger
 
struct  Minimum
 struct to store minimum and temperature More...
 
class  MinimumTracer
 
struct  output
 
class  parser
 The parser class provides the argument parser for the CLI and JSON methods. This is case insensitive. More...
 
class  parserException
 
struct  Phase
 Phase object. More...
 
struct  resultErrorPair
 
struct  status_codes
 status codes struct More...
 
struct  transition_data
 transition data struct More...
 
class  TransitionTracer
 
struct  user_input
 user_input struct to store user input and distribute to the classes More...
 
struct  Vacuum
 Complete vacuum structure of the theory for this parameter point. More...
 

Enumerations

enum class  StatusNLOStability { NotSet , Off , Success , NoNLOStability }
 Possible NLO stability status.
 
enum class  StatusEWSR {
  NotSet , Off , Failure , NotBFB ,
  FlatRegion , EWSymNonRes , EWSymRes
}
 Possible electroweak symmetry restoration status.
 
enum class  StatusTracing {
  NotSet , Success , NoCoverage , NoMinsAtBoundaries ,
  NoGlobMinCoverage , Failure
}
 Possible tracing results.
 
enum class  StatusCoexPair { NotSet , Success , NoCoexPairs }
 Possible status for the coex phase.
 
enum class  StatusCrit {
  NotSet , Success , FalseLower , TrueLower ,
  Failure
}
 Possible status for the critical temperature.
 
enum class  StatusTemperature { NotSet , Success , NotMet , NaN }
 Possible status for the approximated nucleation, exact nucleation, percolation and completion temperature.
 
enum class  StatusGW { NotSet , Success , Failure }
 Possible results for the GW and bounce_sol class.
 
enum class  LoggingLevel {
  None , Default , MinimizerDetailed , MinTracerDetailed ,
  TransitionDetailed , GWDetailed , BounceDetailed , ProgDetailed ,
  EWBGDetailed , Debug , Complete
}
 

Functions

struct resultErrorPair Nintegrate_Inner (BounceSolution &obj, const double &Tprime)
 Nintegrate_Inner Numerical integration of inner integral over inverse Hubble rate for the percolation temperature calculation.
 
struct resultErrorPair Nintegrate_Outer (BounceSolution &obj)
 Nintegrate_Outer Numerical integration of outer integral for the percolation temperature calculation.
 
struct resultErrorPair Nderive_BounceRatio (BounceSolution &obj)
 Nderive_BounceRatio Numerical derivative for the inverse time scale calculation.
 
double SIfunc (const double f)
 SIfunc.
 
double Rfunc (const double f)
 Rfunc.
 
double powspec_density (const double f)
 powspec_density
 
double h2OmSens (const double f)
 return the value of LISA mission nominal sensitivity
 
struct resultErrorPair Nintegrate_SNR (GravitationalWave &obj, const double fmin, const double fmax)
 Nintegrate_SNR Numerical integration of SNR integral.
 
double Getkappa_sw (const double &alpha, const double &vwall, const double &Csound)
 Get efficiency factor kappa_sw for sound waves.
 
double GetK_sw (const double &alpha, const double &vwall, const double &Csound)
 Get K for sound waves.
 
double GetHR (const double &invTimeScale, const double &vwall, const double &Csound)
 Get HR for sound waves.
 
double GetK_turb (const double &alpha, const double &kappa)
 Get K for turbulence.
 
bool IsFluidTurnoverApproxOne (const double &HR, const double &K)
 Determine fluid turnover time regime.
 
std::ostream & operator<< (std::ostream &os, const StatusNLOStability &status)
 Override << operator to handle StatusNLOStability.
 
std::ostream & operator<< (std::ostream &os, const StatusEWSR &status)
 Override << operator to handle StatusEWSR.
 
std::ostream & operator<< (std::ostream &os, const StatusTracing &status)
 Override << operator to handle StatusTracing.
 
std::ostream & operator<< (std::ostream &os, const StatusCoexPair &status)
 Override << operator to handle StatusCoexPair.
 
std::ostream & operator<< (std::ostream &os, const StatusCrit &status)
 Override << operator to handle StatusCrit.
 
std::ostream & operator<< (std::ostream &os, const StatusGW &status)
 Override << operator to handle StatusGW.
 
std::ostream & operator<< (std::ostream &os, const StatusTemperature &status)
 Override << operator to handle StatusTemperature.
 
std::vector< std::vector< double > > Create1DimGrid (const std::vector< double > &point, const int k, const double low_value, const double high_value, const int nsteps=100)
 Create1DimGrid creates a 1-dim grid of given size in index-direction.
 
std::vector< std::vector< double > > Create1DimGrid (const std::vector< double > &min_start, const std::vector< double > &min_end, const int npoints=100)
 Create1DimGrid creates a 1-dim grid of given size between two points.
 
bool almost_the_same (const double &a, const double &b, const double &rel_precision=0.01, const double &num_zero=1e-10)
 
bool almost_the_same (const std::vector< double > &a, const std::vector< double > &b, const bool &allow_for_sign_flip=false, const double &rel_precision=0.01, const double &num_zero=1e-10)
 
void ShowInputError ()
 ShowInputError shows all the available models in the terminal.
 
const std::complex< double > II (0, 1)
 imaginary number i
 
const ISMConstants GetSMConstants ()
 GetSMConstants returns a set of SM contants as indicated by the sources described for each parameter.
 
void ShowLoggerHelp ()
 
void SetLogger (const parser &argparser)
 
void SetLogger (const std::vector< std::string > &args)
 
std::vector< double > NablaNumerical (const std::vector< double > &phi, const std::function< double(std::vector< double >)> &f, const double &eps)
 Numerical method to calculate the gradient of a function f using finite differences method.
 
std::vector< std::vector< double > > HessianNumerical (const std::vector< double > &phi, const std::function< double(std::vector< double >)> &V, double eps)
 Numerical method to calculate the potential's (or other functions's) hessian matrix using finite differences method.
 
bool StringStartsWith (const std::string &str, const std::string &prefix)
 StringStartsWith checks if str starts with prefix.
 
bool StringEndsWith (const std::string &str, const std::string &suffix)
 StringEndsWith tests if str ends with suffix.
 
int factorial (const int &a)
 factorial function
 
template<typename T >
std::vector< T > push_back (std::vector< T > &a, const std::vector< T > &b)
 push back vector into vector
 
template<typename T >
std::string vec_to_string (const std::vector< T > &vec)
 vector to_string
 
std::vector< std::string > split (const std::string &str, char delimiter)
 split string separated by delimiter into substrings
 
template<typename T >
std::ostream & operator<< (std::ostream &os, const std::vector< T > &vec)
 Overload to print out vectors with the << operator.
 
template<typename T >
std::vector< T > operator+ (const std::vector< T > &a, const std::vector< T > &b)
 vector addition
 
template<typename T >
std::vector< T > operator- (const std::vector< T > &a, const std::vector< T > &b)
 vector subtraction
 
template<typename T , typename T2 >
std::vector< T > operator* (const T2 &a, const std::vector< T > &b)
 multiplication of vector with scalar
 
template<typename T , typename T2 >
std::vector< T > operator/ (const std::vector< T > &a, const T2 &b)
 division of vector by scalar
 
template<typename T >
operator* (const std::vector< T > &a, const std::vector< T > &b)
 dot product of two vectors
 
template<typename T >
std::vector< T > operator* (const std::vector< std::vector< T > > &a, const std::vector< T > &b)
 multiplication of matrix with vector
 
template<typename T >
std::vector< T > flatten (std::vector< std::vector< T > > const &vec)
 flatten matrix
 
double L2NormVector (const std::vector< double > &vec)
 L2NormVector.
 
std::vector< std::vector< double > > Transpose (const std::vector< std::vector< double > > &A)
 Calculates the tranpose of a matrix.
 
double Li2 (const double &x)
 Dilogarithm of x.
 
double EllipIntSecond (const double &x)
 Incomplete elliptic integral of the second kind of x with a different parameterization and k^2 = -2.
 
std::ostream & operator<< (std::ostream &os, const ModelID::ModelIDs &Model)
 
std::string ModelIDToString (const ModelID::ModelIDs &Model)
 
double inner_integrand (double Temp, void *params)
 
double outer_integrand (double Temp, void *params)
 
double action_ratio (double Temp, void *params)
 
double snr_integrand (double f, void *params)
 

Variables

const double EmptyValue = NAN
 Value to be store in the columns without any values.
 
const std::unordered_map< StatusNLOStability, std::string > StatusNLOStabilityToString
 Map to convert StatusNLOStability to strins.
 
const std::unordered_map< StatusEWSR, std::string > StatusEWSRToString
 Map to convert StatusEWSRToString to strins.
 
const std::unordered_map< StatusTracing, std::string > StatusTracingToString
 Map to convert StatusTracingToString to strins.
 
const std::unordered_map< StatusCoexPair, std::string > StatusCoexPairToString
 Map to convert StatusCoexPairToString to strins.
 
const std::unordered_map< StatusCrit, std::string > StatusCritToString
 Map to convert StatusCritToString to strins.
 
const std::unordered_map< StatusTemperature, std::string > StatusTemperatureToString
 Map to convert StatusTemperature to strins.
 
const std::unordered_map< StatusGW, std::string > StatusGWToString
 Map to convert StatusGWToString to strins.
 
const bool C_UseParwani = false
 C_UseParwani Use the Parwani Method instead of Arnold-Espinosa.
 
const double C_PT = 0
 C_PT Lower threshold to stop the EWPT calculation.
 
const double C_threshold = 1e-4
 C_threshold threshold to check if a mass is numerically zero.
 
const double C_CWcbGB = 5.0 / 6.0
 C_CWcbGB constant used in the CW potential for gauge bosons.
 
const double C_CWcbHiggs = 1.5
 C_CWcbHiggs constant used in the CW potential for Higgs bosons.
 
const std::string sep = "\t"
 seperator used to write into output files
 

Detailed Description

This classes calculates the Bounce action of the potential with a set temperature.

Function Documentation

◆ almost_the_same() [1/2]

bool BSMPT::almost_the_same ( const double &  a,
const double &  b,
const double &  rel_precision = 0.01,
const double &  num_zero = 1e-10 
)

Returns true if two values are the same given some relative precision

◆ almost_the_same() [2/2]

bool BSMPT::almost_the_same ( const std::vector< double > &  a,
const std::vector< double > &  b,
const bool &  allow_for_sign_flip = false,
const double &  rel_precision = 0.01,
const double &  num_zero = 1e-10 
)

Returns true if two vectors are the element-wise the same given some relative precision

◆ Create1DimGrid() [1/2]

std::vector< std::vector< double > > BSMPT::Create1DimGrid ( const std::vector< double > &  min_start,
const std::vector< double > &  min_end,
const int  npoints = 100 
)

Create1DimGrid creates a 1-dim grid of given size between two points.

Parameters
min_start
min_end
npoints
Returns
npoints long vector of steps on connecting line between min_start and min_end

◆ Create1DimGrid() [2/2]

std::vector< std::vector< double > > BSMPT::Create1DimGrid ( const std::vector< double > &  point,
const int  k,
const double  low_value,
const double  high_value,
const int  nsteps = 100 
)

Create1DimGrid creates a 1-dim grid of given size in index-direction.

Parameters
point
kindex direction in which to create grid
low_value
high_value
nsteps
Returns
gridsize-dim vector of vectors

◆ EllipIntSecond()

double BSMPT::EllipIntSecond ( const double &  x)

Incomplete elliptic integral of the second kind of x with a different parameterization and k^2 = -2.

\( \text{EllipIntSecond}(x) = -i E(i \phi, 2) = \int_0^\phi \sqrt{1+2\sinh^2{\theta}}\,d\theta \)

https://en.wikipedia.org/wiki</Elliptic_integralIncomplete_elliptic_integral_of_the_second_kind https://mathworld.wolfram.com/EllipticIntegraloftheSecondKind.html

Parameters
xreal argument
Returns
double

◆ GetHR()

double BSMPT::GetHR ( const double &  invTimeScale,
const double &  vwall,
const double &  Csound 
)

Get HR for sound waves.

Parameters
invTimeScalebeta/H, inverse time scale
vwallbubble wall velocity
Csoundspeed of sound
Returns
double

◆ GetK_sw()

double BSMPT::GetK_sw ( const double &  alpha,
const double &  vwall,
const double &  Csound 
)

Get K for sound waves.

Parameters
alphastrength of the phase transition
vwallbubble wall velocity
Csoundspeed of sound
Returns
double

◆ GetK_turb()

double BSMPT::GetK_turb ( const double &  alpha,
const double &  kappa 
)

Get K for turbulence.

Parameters
alphastrength of the phase transition
kappaturbulence efficiency factor
Returns
double

◆ Getkappa_sw()

double BSMPT::Getkappa_sw ( const double &  alpha,
const double &  vwall,
const double &  Csound 
)

Get efficiency factor kappa_sw for sound waves.

Parameters
alphastrength of the phase transition
vwallbubble wall velocity
Csoundspeed of sound
Returns
double

◆ GetSMConstants()

const ISMConstants BSMPT::GetSMConstants ( )

GetSMConstants returns a set of SM contants as indicated by the sources described for each parameter.

Returns
The SM Constants used by default in BSMPT

◆ h2OmSens()

double BSMPT::h2OmSens ( const double  f)

return the value of LISA mission nominal sensitivity

Parameters
ffrequency
Returns
value of LISA's nominal sensitivity for a given frequency

◆ HessianNumerical()

std::vector< std::vector< double > > BSMPT::HessianNumerical ( const std::vector< double > &  phi,
const std::function< double(std::vector< double >)> &  V,
double  eps 
)

Numerical method to calculate the potential's (or other functions's) hessian matrix using finite differences method.

\(\frac{\partial^2 V}{\partial \phi_i \phi_j} = \frac{1}{4 \epsilon^2}\left(V(\dots, \vec{\phi}_i + \epsilon , \vec{\phi}_j + \epsilon) - V(\dots, \vec{\phi}_i - \epsilon , \vec{\phi}_j + \epsilon) - V(\dots, \vec{\phi}_i + \epsilon , \vec{\phi}_j - \epsilon) + V(\dots, \vec{\phi}_i - \epsilon , \vec{\phi}_j - \epsilon) \right)\)

where \( \epsilon \) is a small step.

Parameters
phiWhere we want to calculate the Hessian matrix
VPotential (or other function)
epsSize of finite differences step
Returns
std::vector<std::vector<double>> The \( dim \times \dim \) hessian matrix of V taken at phi

◆ IsFluidTurnoverApproxOne()

bool BSMPT::IsFluidTurnoverApproxOne ( const double &  HR,
const double &  K 
)

Determine fluid turnover time regime.

Parameters
HRmean bubble separation
Kfraction of the kinetic energy in the fluid
Returns
true if H*tauSH approx. 1, false if smaller than 1

◆ L2NormVector()

double BSMPT::L2NormVector ( const std::vector< double > &  vec)

L2NormVector.

Parameters
vecvector
Returns
L2 norm of vector

◆ Li2()

double BSMPT::Li2 ( const double &  x)

Dilogarithm of x.

https://en.wikipedia.org/wiki/Dilogarithm

Parameters
xreal argument of from \( (-\infty, 1)\)
Returns
double

◆ NablaNumerical()

std::vector< double > BSMPT::NablaNumerical ( const std::vector< double > &  phi,
const std::function< double(std::vector< double >)> &  f,
const double &  eps 
)

Numerical method to calculate the gradient of a function f using finite differences method.

This method is used while BSMPT is not able to calculate the potential derivative analytically. We used the 4th order method

\(\frac{\partial f}{\partial \phi_i} = \frac{1}{12 \epsilon}\left(-f(\dots ,\vec{\phi}_i + 2 \epsilon ) + 8 f(\dots ,\vec{\phi}_i + \epsilon )- 8 f(\dots ,\vec{\phi}_i - \epsilon ) + f(\dots ,\vec{\phi}_i - 2 \epsilon )\right)\)

where \( \epsilon \) is a small step.

Parameters
phiWhere we want to calculate the gradient
ffunction
epsSize of finite differences step
Returns
std::vector<double> The \( dim \times 1 \) gradient of V taken at phi

◆ Nderive_BounceRatio()

struct resultErrorPair BSMPT::Nderive_BounceRatio ( BounceSolution obj)

Nderive_BounceRatio Numerical derivative for the inverse time scale calculation.

Parameters
objClass reference to pass all needed parameters
Returns
Numerical value of derivative and absolute error

◆ Nintegrate_Inner()

struct resultErrorPair BSMPT::Nintegrate_Inner ( BounceSolution obj,
const double &  Tprime 
)

Nintegrate_Inner Numerical integration of inner integral over inverse Hubble rate for the percolation temperature calculation.

Parameters
objClass reference to pass all needed parameters
Tupper integration boundary
Returns
Numerical value of integral and absolute error

◆ Nintegrate_Outer()

struct resultErrorPair BSMPT::Nintegrate_Outer ( BounceSolution obj)

Nintegrate_Outer Numerical integration of outer integral for the percolation temperature calculation.

Parameters
objClass reference to pass all needed parameters
Returns
Numerical value of integral and absolute error

◆ Nintegrate_SNR()

struct resultErrorPair BSMPT::Nintegrate_SNR ( GravitationalWave obj,
const double  fmin,
const double  fmax 
)

Nintegrate_SNR Numerical integration of SNR integral.

Parameters
objClass reference to pass all needed parameters
Returns
Numerical value of integral and absolute error

◆ operator<<() [1/7]

std::ostream & BSMPT::operator<< ( std::ostream &  os,
const StatusCoexPair status 
)

Override << operator to handle StatusCoexPair.

Parameters
osostream buffer
statusstatus to be printed
Returns
std::ostream& buffer

◆ operator<<() [2/7]

std::ostream & BSMPT::operator<< ( std::ostream &  os,
const StatusCrit status 
)

Override << operator to handle StatusCrit.

Parameters
osostream buffer
statusstatus to be printed
Returns
std::ostream& buffer

◆ operator<<() [3/7]

std::ostream & BSMPT::operator<< ( std::ostream &  os,
const StatusEWSR status 
)

Override << operator to handle StatusEWSR.

Parameters
osostream buffer
statusstatus to be printed
Returns
std::ostream& buffer

◆ operator<<() [4/7]

std::ostream & BSMPT::operator<< ( std::ostream &  os,
const StatusGW status 
)

Override << operator to handle StatusGW.

Parameters
osostream buffer
statusstatus to be printed
Returns
std::ostream& buffer

◆ operator<<() [5/7]

std::ostream & BSMPT::operator<< ( std::ostream &  os,
const StatusNLOStability status 
)

Override << operator to handle StatusNLOStability.

Parameters
osostream buffer
statusstatus to be printed
Returns
std::ostream& buffer

◆ operator<<() [6/7]

std::ostream & BSMPT::operator<< ( std::ostream &  os,
const StatusTemperature status 
)

Override << operator to handle StatusTemperature.

Parameters
osostream buffer
statusstatus to be printed
Returns
std::ostream& buffer

◆ operator<<() [7/7]

std::ostream & BSMPT::operator<< ( std::ostream &  os,
const StatusTracing status 
)

Override << operator to handle StatusTracing.

Parameters
osostream buffer
statusstatus to be printed
Returns
std::ostream& buffer

◆ powspec_density()

double BSMPT::powspec_density ( const double  f)

powspec_density

Parameters
ffrequency
Returns
value of power spectral density-function from LISA mission performance requirement

◆ Rfunc()

double BSMPT::Rfunc ( const double  f)

Rfunc.

Parameters
ffrequency
Returns
value of R-function from LISA mission performance requirement for power spectral density

◆ SIfunc()

double BSMPT::SIfunc ( const double  f)

SIfunc.

Parameters
ffrequency
Returns
value of SI-function from LISA mission performance requirement for power spectral density

◆ StringEndsWith()

bool BSMPT::StringEndsWith ( const std::string &  str,
const std::string &  suffix 
)

StringEndsWith tests if str ends with suffix.

Parameters
str
suffix
Returns

◆ Transpose()

std::vector< std::vector< double > > BSMPT::Transpose ( const std::vector< std::vector< double > > &  A)

Calculates the tranpose of a matrix.

Parameters
Amatrix to be transposed
Returns
std::vector<std::vector<double>> transposed matrix

Variable Documentation

◆ StatusCoexPairToString

const std::unordered_map<StatusCoexPair, std::string> BSMPT::StatusCoexPairToString
Initial value:
{
{StatusCoexPair::NotSet, "not_set"},
{StatusCoexPair::Success, "success"},
{StatusCoexPair::NoCoexPairs, "no_coex_pair"}}

Map to convert StatusCoexPairToString to strins.

◆ StatusCritToString

const std::unordered_map<StatusCrit, std::string> BSMPT::StatusCritToString
Initial value:
{
{StatusCrit::NotSet, "not_set"},
{StatusCrit::Success, "success"},
{StatusCrit::FalseLower, "false_lower"},
{StatusCrit::TrueLower, "true_lower"},
{StatusCrit::Failure, "failure"}}

Map to convert StatusCritToString to strins.

◆ StatusEWSRToString

const std::unordered_map<StatusEWSR, std::string> BSMPT::StatusEWSRToString
Initial value:
{
{StatusEWSR::NotSet, "not_set"},
{StatusEWSR::Off, "off"},
{StatusEWSR::Failure, "failure"},
{StatusEWSR::NotBFB, "non_bfb"},
{StatusEWSR::FlatRegion, "flat_region"},
{StatusEWSR::EWSymNonRes, "ew_syum_non_res"},
{StatusEWSR::EWSymRes, "ew_sym_res"}}

Map to convert StatusEWSRToString to strins.

◆ StatusGWToString

const std::unordered_map<StatusGW, std::string> BSMPT::StatusGWToString
Initial value:
{
{StatusGW::NotSet, "not_set"},
{StatusGW::Success, "success"},
{StatusGW::Failure, "failure"}}

Map to convert StatusGWToString to strins.

◆ StatusNLOStabilityToString

const std::unordered_map<StatusNLOStability, std::string> BSMPT::StatusNLOStabilityToString
Initial value:
{
{StatusNLOStability::NotSet, "not_set"},
{StatusNLOStability::Off, "off"},
{StatusNLOStability::Success, "success"},
{StatusNLOStability::NoNLOStability, "no_nlo_stability"}}

Map to convert StatusNLOStability to strins.

◆ StatusTemperatureToString

const std::unordered_map<StatusTemperature, std::string> BSMPT::StatusTemperatureToString
Initial value:
{{StatusTemperature::NotSet, "not_set"},
{StatusTemperature::Success, "success"},
{StatusTemperature::NotMet, "not_met"},
{StatusTemperature::NaN, "nan"}}

Map to convert StatusTemperature to strins.

◆ StatusTracingToString

const std::unordered_map<StatusTracing, std::string> BSMPT::StatusTracingToString
Initial value:
{
{StatusTracing::NotSet, "not_set"},
{StatusTracing::Success, "success"},
{StatusTracing::NoCoverage, "no_coverage"},
{StatusTracing::NoMinsAtBoundaries, "no_mins_at_boundaries"},
{StatusTracing::NoGlobMinCoverage, "no_glob_min_coverage"},
{StatusTracing::Failure, "failure"}}

Map to convert StatusTracingToString to strins.