Skip to content

Commit

Permalink
Calculate the muon pole mass on demand.
Browse files Browse the repository at this point in the history
  • Loading branch information
Jobst Ziebell authored and Jobst Ziebell committed Aug 6, 2015
1 parent 68517cb commit 2c094d4
Show file tree
Hide file tree
Showing 3 changed files with 33 additions and 21 deletions.
42 changes: 27 additions & 15 deletions meta/GMuonMinus2.m
View file
Expand Up @@ -67,26 +67,38 @@
muonFamily = GetMuonFamily[];

prototypes = ("static const unsigned int muonIndex( void );\n" <>
"static const double muonPhysicalMass( const EvaluationContext &context );\n" <>
"static const double muonCharge( const EvaluationContext &context );");
"static const double muonPhysicalMass( EvaluationContext &context );\n" <>
"static const double muonCharge( EvaluationContext &context );");

definitions = ("static const unsigned int muonIndex( void )\n" <>
"{ unsigned int muonIndex" <>
If[muonIndex =!= Null, " = " <> ToString[muonIndex-1], ""] <>
"; return muonIndex; }\n" <>
"static const double muonPhysicalMass( const EvaluationContext &context )\n" <>
"{\n" <>
IndentText @
("if( context.model.do_calculate_sm_pole_masses() )\n" <>
"static const double muonPhysicalMass( EvaluationContext &context )\n" <>
"{\n" <>
IndentText @
("return PHYSICAL(M" <> ParticleToString[muonFamily] <> ")" <>
If[muonIndex =!= Null, "[muonIndex()]", ""] <> ";\n"
) <>
"return context.mass<" <> ParticleToString[muonFamily] <> ">(" <>
If[muonIndex =!= Null, " muonIndex() ", ""] <> ");\n"
) <>
("static double m_muon_pole = 0.0;\n\n" <>

"if( m_muon_pole == 0.0 )\n" <>
"{\n" <>
IndentText @
("m_muon_pole = context.model.get_physical().M" <>
ParticleToString[muonFamily] <> "(" <>
If[muonIndex =!= Null, ToString[muonIndex-1], ""] <> ");\n\n" <>

"if( m_muon_pole == 0.0 )\n" <>
"{\n" <>
IndentText @
("context.model.calculate_M" <> ParticleToString[muonFamily] <> "_pole();\n" <>
"m_muon_pole = context.model.get_physical().M" <>
ParticleToString[muonFamily] <> "(" <>
If[muonIndex =!= Null, ToString[muonIndex-1], ""] <> ");\n") <>
"}\n") <>
"}\n\n" <>

"return m_muon_pole;\n") <>
"}\n" <>
"static const double muonCharge( const EvaluationContext &context )\n" <>
"static const double muonCharge( EvaluationContext &context )\n" <>
"{ return context.model." <>
NameOfCouplingFunction[{GetPhoton[], GetMuonFamily[], SARAH`bar[GetMuonFamily[]]}] <> "PL" <>
If[muonIndex =!= Null,
Expand Down Expand Up @@ -125,7 +137,7 @@
"struct DiagramEvaluator<OneLoopDiagram<" <>
ToString @ type[[1]] <>
">, PhotonEmitter, ExchangeParticle>\n" <>
"{ static double value( const EvaluationContext &context ); };");
"{ static double value( EvaluationContext &context ); };");

calculationCode = Null;
CreateCalculation[] := Module[{code},
Expand Down Expand Up @@ -534,7 +546,7 @@
"{ " <> StringJoin @ Riffle[ToString /@ indexBounds[[2]], ", "] <> " } };"
);];
definition = ("template<> template<> " <> functionClassName <> "::vertex_type\n" <>
functionClassName <> "::vertex( const indices_type &indices, const EvaluationContext &context )\n" <>
functionClassName <> "::vertex( const indices_type &indices, EvaluationContext &context )\n" <>
"{\n" <>
IndentText @ VertexFunctionBody[parsedVertex] <> "\n" <>
"}");
Expand Down
10 changes: 5 additions & 5 deletions templates/g_muon_minus_2.cpp.in
View file
Expand Up @@ -145,7 +145,7 @@ IndexBounds<0>::indices_type IndexBounds<0>::dummyIndex = {};

// Evaluation context struct, needed by DiagramEvaluator<>::value() further below
struct EvaluationContext {
const @ModelName@<Two_scale> &model;
@ModelName@<Two_scale> &model;

template<class P> double mass(void) const;
template<class P> double mass(unsigned int) const;
Expand Down Expand Up @@ -267,7 +267,7 @@ public:
typedef typename Base::vertex_type vertex_type;

template<class Q = Base>
static vertex_type vertex(const indices_type &indices, const EvaluationContext &context)
static vertex_type vertex(const indices_type &indices, EvaluationContext &context)
{
return Q::vertex(indices, context);
}
Expand Down Expand Up @@ -382,7 +382,7 @@ static inline double OneLoopFunctionF2N(double x)
return -3.0 / (y * y * y) * (1.0 - x * x + 2.0 * x * std::log(x));
}

double calculate_amuon(const @ModelName@<Two_scale> &model)
double calculate_amuon(@ModelName@<Two_scale> &model)
{
#ifdef ENABLE_THREADS
@GMuonMinus2_ThreadedCalculation@
Expand All @@ -394,7 +394,7 @@ double calculate_amuon(const @ModelName@<Two_scale> &model)
@GMuonMinus2_Definitions@

template<class PhotonEmitter, class ExchangeParticle>
double DiagramEvaluator<OneLoopDiagram<3>, PhotonEmitter, ExchangeParticle>::value(const EvaluationContext &context)
double DiagramEvaluator<OneLoopDiagram<3>, PhotonEmitter, ExchangeParticle>::value(EvaluationContext &context)
{
double res = 0.0;

Expand Down Expand Up @@ -481,7 +481,7 @@ double DiagramEvaluator<OneLoopDiagram<3>, PhotonEmitter, ExchangeParticle>::val
}

template<class PhotonEmitter, class ExchangeParticle>
double DiagramEvaluator<OneLoopDiagram<4>, PhotonEmitter, ExchangeParticle>::value(const EvaluationContext &context)
double DiagramEvaluator<OneLoopDiagram<4>, PhotonEmitter, ExchangeParticle>::value(EvaluationContext &context)
{
double res = 0.0;

Expand Down
2 changes: 1 addition & 1 deletion templates/g_muon_minus_2.hpp.in
View file
Expand Up @@ -37,7 +37,7 @@ class Two_scale;
template<> class @ModelName@<Two_scale>;

namespace @ModelName@_GMuonMinus2 {
extern double calculate_amuon(const @ModelName@<Two_scale> &model);
extern double calculate_amuon(@ModelName@<Two_scale> &model);
};
} // namespace flexiblesusy

Expand Down

0 comments on commit 2c094d4

Please sign in to comment.