diff --git a/meta/GMuonMinus2.m b/meta/GMuonMinus2.m index 6405b293c..ff0d35d7e 100644 --- a/meta/GMuonMinus2.m +++ b/meta/GMuonMinus2.m @@ -357,20 +357,12 @@ If you add new kinds of vertices (e.g for new diagram types): (* Use of memoization gives ~30% speedup for the MSSM! *) memoizedVertices = {}; MemoizingVertex[particles_List, options : OptionsPattern[SARAH`Vertex]] := - Module[{memo, ordering, orderedParticles}, - (* First we sort the particles *) - ordering = Ordering[particles]; - orderedParticles = particles[[ordering]]; - - memo = Select[memoizedVertices, MatchesMemoizedVertex[orderedParticles], 1]; - If[memo =!= {}, memo = memo[[1]], + Module[{memo = Select[memoizedVertices, MatchesMemoizedVertex[particles], 1]}, + If[memo =!= {}, memo[[1]], (* Create a new entry *) - memo = SARAH`Vertex[orderedParticles, options]; - AppendTo[memoizedVertices, memo];]; - - (* Now return the particles to their original order *) - memo = OrderVertex[memo, Ordering[ordering]]; - Return[memo]]; + memo = SARAH`Vertex[particles, options]; + AppendTo[memoizedVertices, memo]; + memo]]; MatchesMemoizedVertex[particles_List][vertex_] := MatchQ[particles, Vertices`StripFieldIndices /@ vertex[[1]]]; @@ -567,7 +559,7 @@ If you add new kinds of vertices (e.g for new diagram types): ">\n" <> "{\n" <> IndentText @ - ("static const bool is_permutation = true;\n" <> + ("static const int permutation = " <> ToString[PermutationFactor[particles, ordering]] <> ";\n" <> "typedef " <> orderedVertexFunction <> " orig_type;\n" <> "typedef boost::mpl::vector_c StringJoin @ Riffle[ToString /@ (Ordering[ordering] - 1), ", "] <> @@ -579,6 +571,15 @@ If you add new kinds of vertices (e.g for new diagram types): Return[{prototypes, definitions}]; ]; +PermutationFactor[particles_List, ordering_List] := Module[{ + particleTypes = SARAH`getType[#, False, FlexibleSUSY`FSEigenstates]& /@ + particles + }, + Which[SARAH`VType @@ particleTypes =!= FFV, 1, + OrderedQ @ Select[ordering, particleTypes[[#]] === F &], 1, + True, -1] +]; + (* Creates local declarations of field indices, whose values are taken from the elements of `arrayName'. *) @@ -628,7 +629,7 @@ If you add new kinds of vertices (e.g for new diagram types): Parameters`CreateLocalConstRefs[expr] <> "\n" <> "const " <> GetComplexScalarCType[] <> " result = " <> Parameters`ExpressionToString[expr] <> ";\n\n" <> - "return vertex_type(result);", + "return vertex_type(result * permutationFactor);", "LeftAndRightComponentedVertex", exprL = Vertices`SortCp @ SARAH`Cp[Sequence @@ indexedParticles][SARAH`PL] /. vertexRules; @@ -641,7 +642,7 @@ If you add new kinds of vertices (e.g for new diagram types): Parameters`ExpressionToString[exprL] <> ";\n\n" <> "const " <> GetComplexScalarCType[] <> " right = " <> Parameters`ExpressionToString[exprR] <> ";\n\n" <> - "return vertex_type(left, right);"]; + "return vertex_type(left * permutationFactor, right * permutationFactor);"]; sarahParticles = SARAH`getParticleName /@ particles; particleInfo = Flatten[(Cases[SARAH`Particles[FlexibleSUSY`FSEigenstates], {#, ___}] &) /@ @@ -699,7 +700,7 @@ If you add new kinds of vertices (e.g for new diagram types): prototype = ("template<> struct " <> dataClassName <> "\n" <> "{\n" <> IndentText @ - ("static const bool is_permutation = false;\n" <> + ("static const int permutation = 0;\n" <> "typedef IndexBounds<" <> ToString @ NumberOfIndices[parsedVertex] <> "> index_bounds;\n" <> "typedef " <> VertexClassName[parsedVertex] <> " vertex_type;\n" <> "typedef boost::mpl::vector_c @@ -723,7 +724,7 @@ If you add new kinds of vertices (e.g for new diagram types): ); ]; definition = ("template<> template<> " <> functionClassName <> "::vertex_type\n" <> - functionClassName <> "::vertex(const indices_type &indices, EvaluationContext &context)\n" <> + functionClassName <> "::vertex(const indices_type &indices, EvaluationContext &context, int permutationFactor)\n" <> "{\n" <> IndentText @ VertexFunctionBody[parsedVertex] <> "\n" <> "}"); diff --git a/model_files/MRSSM2/FlexibleSUSY.m.in b/model_files/MRSSM2/FlexibleSUSY.m.in index c4988b50f..135796081 100644 --- a/model_files/MRSSM2/FlexibleSUSY.m.in +++ b/model_files/MRSSM2/FlexibleSUSY.m.in @@ -87,6 +87,11 @@ InitialGuessAtLowScale = { {Ye, Automatic} }; +ExtraSLHAOutputBlocks = { + {FlexibleSUSYLowEnergy, + {{21, FlexibleSUSYObservable`aMuon} } } +}; + DefaultPoleMassPrecision = HighPrecision; HighPoleMassPrecision = {hh, Ah, Hpm}; MediumPoleMassPrecision = {}; diff --git a/templates/a_muon.cpp.in b/templates/a_muon.cpp.in index af359396e..a4e9b4c7b 100644 --- a/templates/a_muon.cpp.in +++ b/templates/a_muon.cpp.in @@ -296,15 +296,15 @@ public: template struct VertexFunctionData; /** - * @class VertexFunctionHelper + * @class VertexFunctionHelper * @brief A helper template that is used by VertexFunction * * You should never have to use this class directly. It is just a * helper class. */ -template class VertexFunctionHelper; +template class VertexFunctionHelper; -template class VertexFunctionHelper +template class VertexFunctionHelper<0, P...> : public VertexFunctionData { private: typedef VertexFunctionData Base; @@ -326,9 +326,11 @@ public: auto end = indices.begin() + boost::mpl::at_c::type::value; return std::vector(begin, end); } + + static constexpr int permutationFactor = 1; }; -template class VertexFunctionHelper +template class VertexFunctionHelper : public VertexFunctionData::orig_type { private: typedef typename VertexFunctionData::orig_type orig_type; @@ -345,6 +347,8 @@ public: boost::mpl::at_c::type::value > (indices); } + + static constexpr int permutationFactor = permutation; }; /** @@ -359,9 +363,9 @@ public: * possible particle permutations. */ template class VertexFunction - : public VertexFunctionHelper::is_permutation, P...> { + : public VertexFunctionHelper::permutation, P...> { private: - typedef VertexFunctionHelper::is_permutation, P...> Base; + typedef VertexFunctionHelper::permutation, P...> Base; public: typedef typename Base::indices_type indices_type; typedef typename Base::index_bounds index_bounds; @@ -395,9 +399,9 @@ public: * For canonical P... the definitions have to be supplied separately. */ template - static vertex_type vertex(const indices_type &indices, EvaluationContext &context) + static vertex_type vertex(const indices_type &indices, EvaluationContext &context, int permutationFactor = 1) { - return Q::vertex(indices, context); + return Q::vertex(indices, context, Base::permutationFactor); } }; diff --git a/test/module.mk b/test/module.mk index d32af4246..791d2e2ae 100644 --- a/test/module.mk +++ b/test/module.mk @@ -180,6 +180,11 @@ TEST_SRC += \ $(DIR)/test_CMSSM_database.cpp endif +ifeq ($(WITH_MRSSM2),yes) +TEST_SRC += \ + $(DIR)/test_MRSSM2_gmm2.cpp +endif + endif # ifneq ($(findstring two_scale,$(SOLVERS)),) ifneq ($(findstring semi_analytic,$(SOLVERS)),) @@ -800,6 +805,8 @@ $(DIR)/test_sfermions.x: $(LIBSoftsusyMSSM) $(LIBCMSSM) $(LIBFLEXI) $(LIBTEST) $ $(DIR)/test_CMSSM_database.x: $(DIR)/test_CMSSM_database.o $(LIBCMSSM) $(LIBFLEXI) $(LIBTEST) $(filter-out -%,$(LOOPFUNCLIBS)) $(CXX) $(CXXFLAGS) $(CPPFLAGS) -o $@ $(call abspathx,$^) $(BOOSTTESTLIBS) $(BOOSTTHREADLIBS) $(GSLLIBS) $(FLIBS) $(SQLITELIBS) $(THREADLIBS) +$(DIR)/test_MRSSM2_gmm2.x: $(LIBMRSSM2) $(LIBFLEXI) $(filter-out -%,$(LOOPFUNCLIBS)) + $(DIR)/test_CMSSM_model.x: $(LIBSoftsusyMSSM) $(LIBCMSSM) $(LIBFLEXI) $(LIBTEST) $(filter-out -%,$(LOOPFUNCLIBS)) $(DIR)/test_CMSSM_info.x: $(LIBCMSSM) $(LIBFLEXI) $(filter-out -%,$(LOOPFUNCLIBS)) diff --git a/test/test_MRSSM2.hpp b/test/test_MRSSM2.hpp new file mode 100644 index 000000000..491a54922 --- /dev/null +++ b/test/test_MRSSM2.hpp @@ -0,0 +1,45 @@ +// ==================================================================== +// This file is part of FlexibleSUSY. +// +// FlexibleSUSY is free software: you can redistribute it and/or modify +// it under the terms of the GNU General Public License as published +// by the Free Software Foundation, either version 3 of the License, +// or (at your option) any later version. +// +// FlexibleSUSY is distributed in the hope that it will be useful, but +// WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// General Public License for more details. +// +// You should have received a copy of the GNU General Public License +// along with FlexibleSUSY. If not, see +// . +// ==================================================================== + +#ifndef TEST_MRSSM2_H +#define TEST_MRSSM2_H + +#include + +#include "lowe.h" +#include "MRSSM2_two_scale_spectrum_generator.hpp" + +using namespace flexiblesusy; + +MRSSM2_slha> +setup_MRSSM2(const MRSSM2_input_parameters& input) +{ + softsusy::QedQcd qedqcd; + + Spectrum_generator_settings settings; + settings.set(Spectrum_generator_settings::calculate_sm_masses, 0); + settings.set(Spectrum_generator_settings::calculate_bsm_masses, 0); + + MRSSM2_spectrum_generator spectrum_generator; + spectrum_generator.set_settings(settings); + spectrum_generator.run(qedqcd, input); + + return std::get<0>(spectrum_generator.get_models_slha()); +} + +#endif diff --git a/test/test_MRSSM2_gmm2.cpp b/test/test_MRSSM2_gmm2.cpp new file mode 100644 index 000000000..dd8c7cc4c --- /dev/null +++ b/test/test_MRSSM2_gmm2.cpp @@ -0,0 +1,82 @@ +// ==================================================================== +// This file is part of FlexibleSUSY. +// +// FlexibleSUSY is free software: you can redistribute it and/or modify +// it under the terms of the GNU General Public License as published +// by the Free Software Foundation, either version 3 of the License, +// or (at your option) any later version. +// +// FlexibleSUSY is distributed in the hope that it will be useful, but +// WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// General Public License for more details. +// +// You should have received a copy of the GNU General Public License +// along with FlexibleSUSY. If not, see +// . +// ==================================================================== + +#define BOOST_TEST_DYN_LINK +#define BOOST_TEST_MODULE test_MRSSM2_gmm2 + +#include + +#include "test_MRSSM2.hpp" + +#include "wrappers.hpp" +#include "MRSSM2_a_muon.hpp" + +using namespace flexiblesusy; + +BOOST_AUTO_TEST_CASE( test_amu ) +{ + MRSSM2_input_parameters input; + + // chargino dominance + input.TanBeta = 10; + input.Ms = 1000; + input.LamTDInput = -1.0; + input.LamTUInput = -1.0; + input.LamSDInput = 1.1; + input.LamSUInput = -1.1; + input.MuDInput = 400; + input.MuUInput = 400; + input.BMuInput = Sqr(300); + input.mq2Input << Sqr(1000), 0, 0, + 0, Sqr(1000), 0, + 0, 0, Sqr(1000); + input.ml2Input << Sqr(500), 0, 0, + 0, Sqr(500), 0, + 0, 0, Sqr(500); + input.md2Input << Sqr(1000), 0, 0, + 0, Sqr(1000), 0, + 0, 0, Sqr(1000); + input.mu2Input << Sqr(1000), 0, 0, + 0, Sqr(1000), 0, + 0, 0, Sqr(1000); + input.me2Input << Sqr(500), 0, 0, + 0, Sqr(500), 0, + 0, 0, Sqr(500); + input.mS2Input = Sqr(2000); + input.mT2Input = Sqr(3000); + input.moc2Input = Sqr(1000); + input.mRd2Input = Sqr(700); + input.mRu2Input = Sqr(1000); + input.MDBSInput = 1000; + input.MDWBTInput = 500; + input.MDGocInput = 1500; + + MRSSM2_slha> m = setup_MRSSM2(input); + + auto amu = MRSSM2_a_muon::calculate_a_muon(m); + BOOST_CHECK_CLOSE(amu, -8.30e-11, 1.0); + + // neutralino dominance + input.ml2Input << Sqr(8000), 0, 0, + 0, Sqr(8000), 0, + 0, 0, Sqr(8000); + m = setup_MRSSM2(input); + + amu = MRSSM2_a_muon::calculate_a_muon(m); + BOOST_CHECK_CLOSE(amu, 6.33e-12, 1.0); +}