|
BSMPT 3.1.4
BSMPT - Beyond the Standard Model Phase Transitions : A C++ package for the computation of the EWPT in BSM models
|
This classes calculates the Bounce action of the potential with a set temperature. More...
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... | |
| struct | BPLParameters |
| 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... | |
| struct | DBPLParameters |
| 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 | TransitionTemperature { NotSet , ApproxNucleation , Nucleation , Percolation , Completion } |
| Possible transitions temperatures. | |
| 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 | GetK_sw (const double &alpha, const double &kappa_sw) |
| Get the kinetic energy fraction \( K \). | |
| double | GetHstar0 (const double &temp, const double &gstar) |
| Calculate the Hubble rate at transition time refshifted to today. | |
| double | GetKtilde (const double &alpha) |
| Get \( \tilde{K} = \frac{\alpha}{1+\alpha}\). | |
| double | GetHtauSH (const double HR, const double K_sw) |
| Calculate \( H_*\tau_\text{sh} = H_* R_* / \sqrt{\bar v_f^2} \). | |
| double | GetHtauSW (const double HR, const double K_sw) |
| Calculate \( H_*\tau_\text{sw} = \min(1,H_*\tau_\text{sh}) \). | |
| double | GetYpsilon (const double HR, const double K_sw) |
| Calculate \(
\Upsilon=1-\frac{1}{\sqrt{1+2H_{\ast}\tau_{\mathrm{sw}}}} \) from https://arxiv.org/abs/1903.09642. | |
| 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) |
| 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::ostream & | operator<< (std::ostream &os, const ModelID::ModelIDs &Model) |
| operator << overload for the model parameter | |
| 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. | |
| void | ShowInputError () |
| ShowInputError shows all the available models in the terminal. | |
| template<typename key , typename value > | |
| std::unordered_map< value, key > | InvertMap (const std::unordered_map< key, value > &originalMap, const std::string &errorOnDuplicateValue) |
| Inverts a map. | |
| 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 > | |
| 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. | |
| 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 std::vector< double > | TGstarLowT |
| List of temperatures below T = 214 GeV for the gstar spline construction. | |
| const std::vector< double > | TGstarHighT |
| List of temperatures above T = 214 GeV for the gstar spline construction. | |
| const std::vector< double > | GstarLowT |
| List of gstar below T = 214 GeV for the gstar spline construction. | |
| const std::vector< double > | GstarHighT |
| List of gstar above T = 214 GeV for the gstar spline construction. | |
| 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 strings. | |
| const std::unordered_map< StatusCrit, std::string > | StatusCritToString |
| Map to convert StatusCritToString to strings. | |
| const std::unordered_map< StatusTemperature, std::string > | StatusTemperatureToString |
| Map to convert StatusTemperature to strings. | |
| const std::unordered_map< StatusGW, std::string > | StatusGWToString |
| Map to convert StatusGWToString to strings. | |
| 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 bool | C_UseParwani = false |
| C_UseParwani Use the Parwani Method instead of Arnold-Espinosa. | |
| const std::string | sep = "\t" |
| seperator used to write into output files | |
This classes calculates the Bounce action of the potential with a set temperature.
| 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
| 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
| 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.
| min_start | |
| min_end | |
| npoints |
| 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.
| point | |
| k | index direction in which to create grid |
| low_value | |
| high_value | |
| nsteps |
| 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_integral#Incomplete_elliptic_integral_of_the_second_kind https://mathworld.wolfram.com/EllipticIntegraloftheSecondKind.html
| x | real argument |
| double BSMPT::GetHstar0 | ( | const double & | temp, |
| const double & | gstar | ||
| ) |
Calculate the Hubble rate at transition time refshifted to today.
| temp | transition temperature |
| gstar | effective d.o.f. at transition time |
| double BSMPT::GetHtauSH | ( | const double | HR, |
| const double | K_sw | ||
| ) |
Calculate \( H_*\tau_\text{sh} = H_* R_* / \sqrt{\bar v_f^2} \).
| HR | |
| K_sw |
| double BSMPT::GetHtauSW | ( | const double | HR, |
| const double | K_sw | ||
| ) |
Calculate \( H_*\tau_\text{sw} = \min(1,H_*\tau_\text{sh}) \).
| HR | |
| K_sw |
| double BSMPT::GetK_sw | ( | const double & | alpha, |
| const double & | kappa_sw | ||
| ) |
Get the kinetic energy fraction \( K \).
| alpha | strength of the phase transition |
| kappa_sw | efficiency factor |
| double BSMPT::GetKtilde | ( | const double & | alpha | ) |
Get \( \tilde{K} = \frac{\alpha}{1+\alpha}\).
| alpha | strength of the phase transition |
| const ISMConstants BSMPT::GetSMConstants | ( | ) |
GetSMConstants returns a set of SM contants as indicated by the sources described for each parameter.
| double BSMPT::GetYpsilon | ( | const double | HR, |
| const double | K_sw | ||
| ) |
Calculate \( \Upsilon=1-\frac{1}{\sqrt{1+2H_{\ast}\tau_{\mathrm{sw}}}} \) from https://arxiv.org/abs/1903.09642.
| HR | |
| K_sw |
| double BSMPT::h2OmSens | ( | const double | f | ) |
return the value of LISA mission nominal sensitivity
| f | frequency |
| 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.
| phi | Where we want to calculate the Hessian matrix |
| V | Potential (or other function) |
| eps | Size of finite differences step |
| std::unordered_map< value, key > BSMPT::InvertMap | ( | const std::unordered_map< key, value > & | originalMap, |
| const std::string & | errorOnDuplicateValue | ||
| ) |
Inverts a map.
| originalMap | Map to be inverted |
| errorOnDuplicateValue | Error message on duplicate value |
| std::runtime_error | if the new map would have duplicate keys |
| double BSMPT::L2NormVector | ( | const std::vector< double > & | vec | ) |
L2NormVector.
| vec | vector |
| double BSMPT::Li2 | ( | const double & | x | ) |
Dilogarithm of x.
https://en.wikipedia.org/wiki/Dilogarithm
| x | real argument of from \( (-\infty, 1)\) |
| 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.
| phi | Where we want to calculate the gradient |
| f | function |
| eps | Size of finite differences step |
| struct resultErrorPair BSMPT::Nderive_BounceRatio | ( | BounceSolution & | obj | ) |
Nderive_BounceRatio Numerical derivative for the inverse time scale calculation.
| obj | Class reference to pass all needed parameters |
| 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.
| obj | Class reference to pass all needed parameters |
| T | upper integration boundary |
| struct resultErrorPair BSMPT::Nintegrate_Outer | ( | BounceSolution & | obj | ) |
Nintegrate_Outer Numerical integration of outer integral for the percolation temperature calculation.
| obj | Class reference to pass all needed parameters |
| struct resultErrorPair BSMPT::Nintegrate_SNR | ( | GravitationalWave & | obj, |
| const double | fmin, | ||
| const double | fmax | ||
| ) |
Nintegrate_SNR Numerical integration of SNR integral.
| obj | Class reference to pass all needed parameters |
| std::ostream & BSMPT::operator<< | ( | std::ostream & | os, |
| const StatusCoexPair & | status | ||
| ) |
Override << operator to handle StatusCoexPair.
| os | ostream buffer |
| status | status to be printed |
| std::ostream & BSMPT::operator<< | ( | std::ostream & | os, |
| const StatusCrit & | status | ||
| ) |
Override << operator to handle StatusCrit.
| os | ostream buffer |
| status | status to be printed |
| std::ostream & BSMPT::operator<< | ( | std::ostream & | os, |
| const StatusEWSR & | status | ||
| ) |
Override << operator to handle StatusEWSR.
| os | ostream buffer |
| status | status to be printed |
| std::ostream & BSMPT::operator<< | ( | std::ostream & | os, |
| const StatusGW & | status | ||
| ) |
Override << operator to handle StatusGW.
| os | ostream buffer |
| status | status to be printed |
| std::ostream & BSMPT::operator<< | ( | std::ostream & | os, |
| const StatusNLOStability & | status | ||
| ) |
Override << operator to handle StatusNLOStability.
| os | ostream buffer |
| status | status to be printed |
| std::ostream & BSMPT::operator<< | ( | std::ostream & | os, |
| const StatusTemperature & | status | ||
| ) |
Override << operator to handle StatusTemperature.
| os | ostream buffer |
| status | status to be printed |
| std::ostream & BSMPT::operator<< | ( | std::ostream & | os, |
| const StatusTracing & | status | ||
| ) |
Override << operator to handle StatusTracing.
| os | ostream buffer |
| status | status to be printed |
| double BSMPT::powspec_density | ( | const double | f | ) |
powspec_density
| f | frequency |
| double BSMPT::Rfunc | ( | const double | f | ) |
Rfunc.
| f | frequency |
| double BSMPT::SIfunc | ( | const double | f | ) |
SIfunc.
| f | frequency |
| bool BSMPT::StringEndsWith | ( | const std::string & | str, |
| const std::string & | suffix | ||
| ) |
StringEndsWith tests if str ends with suffix.
| str | |
| suffix |
| std::vector< std::vector< double > > BSMPT::Transpose | ( | const std::vector< std::vector< double > > & | A | ) |
Calculates the tranpose of a matrix.
| A | matrix to be transposed |
| const std::vector< double > BSMPT::GstarHighT |
List of gstar above T = 214 GeV for the gstar spline construction.
| const std::vector< double > BSMPT::GstarLowT |
List of gstar below T = 214 GeV for the gstar spline construction.
| const std::unordered_map<StatusCoexPair, std::string> BSMPT::StatusCoexPairToString |
Map to convert StatusCoexPairToString to strings.
| const std::unordered_map<StatusCrit, std::string> BSMPT::StatusCritToString |
Map to convert StatusCritToString to strings.
| const std::unordered_map<StatusEWSR, std::string> BSMPT::StatusEWSRToString |
Map to convert StatusEWSRToString to strins.
| const std::unordered_map<StatusGW, std::string> BSMPT::StatusGWToString |
Map to convert StatusGWToString to strings.
| const std::unordered_map<StatusNLOStability, std::string> BSMPT::StatusNLOStabilityToString |
Map to convert StatusNLOStability to strins.
| const std::unordered_map<StatusTemperature, std::string> BSMPT::StatusTemperatureToString |
Map to convert StatusTemperature to strings.
| const std::unordered_map<StatusTracing, std::string> BSMPT::StatusTracingToString |
Map to convert StatusTracingToString to strins.
| const std::vector< double > BSMPT::TGstarLowT |
List of temperatures below T = 214 GeV for the gstar spline construction.