Merge branch 'development' into feature-SMtower-2L

ChangeLog

@@ -1,4 +1,4 @@
-FlexibleSUSY-1.7.2 [not released yet]
+FlexibleSUSY-1.7.2 [December, 15 2016]
 
  * Feature (commit b052e35): New flag FlexibleSUSY[23] to disable the
    pole mass calculation of the non-SM particles.  This flag is useful
@@ -11,19 +11,20 @@ FlexibleSUSY-1.7.2 [not released yet]
    very large SUSY scales.  In such a case FlexibleSUSY[23] could be
    set to 0 to suppress the calculation of the non-SM pole masses.
 
+ * Feature (commit 998f11e): Slightly improved speed of the RG
+   running.
+
  * Change (commit 189f508): Speed-up the calculation of the 2L Higgs
    mass corrections in the MSSM and NMSSM, if multi-threading is used,
    by locking the mutex only for the O(at*at) corrections.
 
  * Change: The limits sin(2*theta) = 0 and m_stop1 = m_stop2 have been
    implemented for the 2L O(at*as) Higgs pole mass corrections in the
-   MSSM.
-
- * Bugfix (commit c35dcb2): Fix linking problem of LibraryLink on Mac.
+   MSSM to avoid numeric instabilities.
 
  * Bugfix (commit 20f169f): Re-calculate W pole mass in
    FlexibleEFTHiggs.  Before this commit the electroweak gauge
-   couplings in FlexibleEFTHiggs are be wrong in scenarios with very
+   couplings in FlexibleEFTHiggs are wrong in scenarios with very
    small alpha_em(MZ) (< 1/1000) and/or a small Z pole mass (< 10
    GeV).
 
@@ -32,8 +33,8 @@ FlexibleSUSY-1.7.2 [not released yet]
    the running BSM masses (at the SUSY scale) have been used as
    convergence criterion.  However, they tend to converge very fast,
    compared to the running SM masses at the electroweak scale.  For a
-   more reliable convergence criterion, both the running BSM and SM
-   masses have to be used.
+   more reliable convergence criterion, now both the running BSM and
+   SM masses are used.
 
  * Bugfix (commits 5e1b6b3, cc5bfae): Correction of the 2-loop and
    3-loop QCD corrections to the top pole mass in the Standard Model
@@ -43,19 +44,30 @@ FlexibleSUSY-1.7.2 [not released yet]
    top pole mass.  Before these commits, FlexibleSUSY used the
    expressions from theses references, but wrote result in terms of
    Log[Q^2/mt^2], where mt is the MS-bar mass, while not accounting
-   for the difference between Mt and mt in the logarithms.
+   for the difference between Mt and mt in the logarithms.  This
+   bugfix affects the Higgs pole mass at the 3-loop level.
 
  * Bugfix (commit cecff4b): Flag scalar or vector boson gauge singlet
    tachyons.
 
  * Bugfix (commit 4a3fb5b): Input tan(beta) at the SUSY scale, instead
-   of the matching scale.  This matters when the (unphysical) matching
-   is varied through FlexibleSUSY[19].
+   of at the matching scale in the FlexibleEFTHiggs model files.  This
+   difference matters when the (unphysical) matching scale is varied
+   through FlexibleSUSY[19].
+
+ * Bugfix (commit c35dcb2): Fixed linking problem of the LibraryLink
+   on Mac.
 
  * Bugfix (commits a643be5, cc9ebf1): Avoid function call ambiguities
    when multiple LibraryLink libraries are loaded into Mathematica at
    the same time.
 
+ * Bugfix (commit 1f8e135): Correcting FS<model>Set[] function for
+   models with matrix-valued parameters.
+
+ * Bugfix (commit 4097708): The generated LibraryLink files are now
+   added to the model tarball created by `make pack-<model>-src'.
+
 FlexibleSUSY-1.7.1 [October, 15 2016]
 
  * Change (commit b1efa8c): Updated to GM2Calc 1.3.0.

config/config.h.in

@@ -1,6 +1,9 @@
 /* DO NOT EDIT! GENERATED AUTOMATICALLY! */
 /* generated at @DATE@ */
 
+#ifndef CONFIG_H
+#define CONFIG_H
+
 #define PKGNAME "@PKGNAME@"
 
 #define FLEXIBLESUSY_VERSION  "@FLEXIBLESUSY_VERSION@"
@@ -85,3 +88,5 @@
 
 /* Enable statements for addons */
 @DEFINE_ENABLE_ADDONS@
+
+#endif

configure

@@ -6,7 +6,7 @@
 PROGRAM_NAME=FlexibleSUSY
 FLEXIBLESUSY_MAJOR=1
 FLEXIBLESUSY_MINOR=7
-FLEXIBLESUSY_PATCH=1
+FLEXIBLESUSY_PATCH=2
 FLEXIBLESUSY_EXTRA=""
 FLEXIBLESUSY_VERSION="${FLEXIBLESUSY_MAJOR}.${FLEXIBLESUSY_MINOR}.${FLEXIBLESUSY_PATCH}${FLEXIBLESUSY_EXTRA}"
 GIT_COMMIT=$(git rev-parse HEAD 2> /dev/null || echo unknown)

meta/BetaFunction.m

@@ -122,8 +122,11 @@
            expr
           ];
 
+TimeConstrainedSimplify[expr_] :=
+    TimeConstrained[Simplify[expr], FlexibleSUSY`FSSimplifyBetaFunctionsTimeConstraint, expr];
+
 CollectMatMul[expr_] :=
-    TimeConstrained[Collect[expr, SARAH`MatMul[___]],
+    TimeConstrained[Collect[expr, SARAH`MatMul[___], TimeConstrainedSimplify],
                     FlexibleSUSY`FSSimplifyBetaFunctionsTimeConstraint,
                     expr];
 
@@ -236,7 +239,7 @@
             Return[{localDecl, beta1L, beta2L, beta3L}];
           ];
 
-CreateBetaFunction[betaFunctions_List, sarahTraces_List] :=
+CreateBetaFunction[betaFunctions_List] :=
     Module[{def = "",
             localDecl = "", beta1L = "", beta2L = "", beta3L = "",
             allDecl = "", allBeta = "",
@@ -291,7 +294,7 @@
                    (* protect tensor products *)
                    expr = CConversion`ProtectTensorProducts[#, name]& /@ expr;
                    (* simplify expressions *)
-                   expr = TimeConstrained[Simplify[#], FlexibleSUSY`FSSimplifyBetaFunctionsTimeConstraint, #]& /@ expr;
+                   expr = TimeConstrainedSimplify /@ expr;
                    AppendTo[lst, BetaFunction[name, type, expr]];
                   ];
               ];

meta/FlexibleSUSY.m

@@ -582,7 +582,7 @@ FlexibleSUSY model file (FlexibleSUSY.m).
           numberOfParameters = BetaFunction`CountNumberOfParameters[betaFun] + numberOfBaseClassParameters;
           (* create C++ functions and parameter declarations *)
           sarahTraces          = Traces`ConvertSARAHTraces[additionalTraces];
-          beta                 = BetaFunction`CreateBetaFunction[betaFun, sarahTraces];
+          beta                 = BetaFunction`CreateBetaFunction[betaFun];
           setter               = BetaFunction`CreateSetters[betaFun];
           getter               = BetaFunction`CreateGetters[betaFun];
           parameterDef         = BetaFunction`CreateParameterDefinitions[betaFun];

models/SoftsusyFlavourMSSM/flavoursoft.cpp

@@ -212,9 +212,8 @@ const DoubleMatrix FlavourMssmSoftsusy::displayMns() const {
 			      displayThetaB23(), d);
 }
 
-#define HR "----------------------------------------------------------"
-
 ostream & operator <<(ostream & left, const FlavourMssmSoftsusy & m) {
+  const std::string HR = "----------------------------------------------------------";
   left << m.displayMssmSoft();
   left << "Flavour violating parameters:\n";
   left << "CKM data: 12=" << m.displayTheta12() 
@@ -240,8 +239,6 @@ ostream & operator <<(ostream & left, const FlavourMssmSoftsusy & m) {
   return left;
 }
 
-#undef HR
-
 void FlavourMssmSoftsusy::modselSLHA(ostream & out, const char model []) {
   MssmSoftsusy::modselSLHA(out, model);
   out << "     6    " << 1 << "   # flavour violating MSSM\n"; 

models/SoftsusyMSSM/mssmUtils.cpp

@@ -417,6 +417,7 @@ string recogLsp(int temp, int posj) {
 }
 
 ostream & operator <<(ostream &left, const MssmSoftsusy &s) {
+  const std::string HR = "----------------------------------------------------------";
   left << HR << endl;
   left << "Gravitino mass M3/2: " << s.displayGravitino() << endl;
   left << "Msusy: " << s.displayMsusy() << " MW: " << s.displayMw()

models/SoftsusyMSSM/physpars.cpp

@@ -62,8 +62,6 @@ void sPhysical::display(double *a) const {
   a[k++] = thetaA0;
 }
 
-#define HR "---------------------------------------------------------------\n"
-
 std::ostream & operator <<(std::ostream & left, const drBarPars &s) {
   left << s.displaysPhysical();
   left << "BPMZ conventions, N" << s.nBpmz << "U" << s.uBpmz << "V" 
@@ -113,8 +111,6 @@ std::istream & operator >>(std::istream & left, sPhysical &s) {
   return left;
 }
 
-#undef HR
-
 ostream & operator <<(ostream &st, const sProblem & p) {
   if (!p.test()) return st;
   st << "[ ";

models/SoftsusyMSSM/softpars.cpp

@@ -3079,8 +3079,6 @@ void SoftPars<Susy, Brevity>::standardSugra(double m0,  double m12, double a0) {
   universalTrilinears(a0);
 }
 
-#define HR "---------------------------------------------------------------\n"
-
 template<class Susy, class Brevity>
 ostream & operator <<(ostream &left, const SoftPars<Susy, Brevity> &s) {
   left << "SUSY breaking MSSM parameters at Q: " << s.displayMu() << endl;
@@ -3099,8 +3097,6 @@ ostream & operator <<(ostream &left, const SoftPars<Susy, Brevity> &s) {
   return left;
 }
 
-#undef HR
-
 template<class Susy, class Brevity>
 void SoftPars<Susy, Brevity>::inputSoftParsOnly() {
   char c[70];

models/SoftsusyMSSM/softsusy.cpp

@@ -1105,8 +1105,6 @@ void Softsusy<SoftPars>::rewsb(int sgnMu, double mt, const DoubleVector & /* par
   }
 }
 
-#define HR "----------------------------------------------------------"
-
 
 /// Gives a summary of important properties of a SUSY object: mu, m3sq, mH0,
 /// mChi0, lightest stau, lightest mGluino, lightest stop, lightest chargino,

models/SoftsusyMSSM/susy.cpp

@@ -12,8 +12,6 @@
 
 namespace softsusy {
 
-#define HR "---------------------------------------------------------------\n"
-
 const sBrevity & sBrevity::operator=(const sBrevity &s) {
   if (this == &s) return *this;
   dt = s.dt; ut = s.ut; et = s.et; 
@@ -640,8 +638,6 @@ void MssmSusy::getMasses(QedQcd & r, double vev) const {
   r.setMass(mTau, v1 * e1(3, 3));
 }
 
-#undef HR
-
 // Rotates to quark mass basis, returning the mixing matrices defined as 
 // yu_diag = vul yu vur^T  
 // yd_diag = vdl yd vdr^T 

src/rk.cpp

@@ -18,9 +18,11 @@ using namespace Eigen;
 
 namespace runge_kutta {
 
+namespace {
 /// Returns |a| with sign of b in front
 double sign(double a, double b)
 { return b >= 0 ? fabs(a) : -fabs(a); }
+} // anonymous namespace
 
 void integrateOdes(ArrayXd& ystart, double from, double to, double eps,
                    double h1, double hmin, Derivs derivs,
@@ -145,4 +147,4 @@ void rungeKuttaStep(const ArrayXd& y, const ArrayXd& dydx, double x,
 
 } // namespace runge_kutta
 
-}
+} // namespace flexiblesusy

templates/librarylink.m.in

@@ -164,8 +164,67 @@ FS@ModelName@Set[handle_Integer, a___, (fsSettings | fsSMParameters | fsModelPar
 FS@ModelName@Set[handle_Integer, p:OptionsPattern[]] :=
     FS@ModelName@SetLib[
         handle,
-        First /@ Options[FS@ModelName@Set] /.
+        ReleaseHold[Hold[{
+            (* spectrum generator settings *)
+            OptionValue[precisionGoal],
+            OptionValue[maxIterations],
+            OptionValue[calculateStandardModelMasses],
+            OptionValue[poleMassLoopOrder],
+            OptionValue[ewsbLoopOrder],
+            OptionValue[betaFunctionLoopOrder],
+            OptionValue[thresholdCorrectionsLoopOrder],
+            OptionValue[higgs2loopCorrectionAtAs],
+            OptionValue[higgs2loopCorrectionAbAs],
+            OptionValue[higgs2loopCorrectionAtAt],
+            OptionValue[higgs2loopCorrectionAtauAtau],
+            OptionValue[forceOutput],
+            OptionValue[topPoleQCDCorrections],
+            OptionValue[betaZeroThreshold],
+            OptionValue[forcePositiveMasses],
+            OptionValue[poleMassScale],
+            OptionValue[eftPoleMassScale],
+            OptionValue[eftMatchingScale],
+            OptionValue[eftMatchingLoopOrderUp],
+            OptionValue[eftMatchingLoopOrderDown],
+            OptionValue[eftHiggsIndex],
+            OptionValue[calculateBSMMasses],
+            OptionValue[parameterOutputScale],
+
+            (* Standard Model input parameters *)
+            OptionValue[alphaEmMZ],
+            OptionValue[GF],
+            OptionValue[alphaSMZ],
+            OptionValue[MZ],
+            OptionValue[mbmb],
+            OptionValue[Mt],
+            OptionValue[Mtau],
+            OptionValue[Mv3],
+            OptionValue[MW],
+            OptionValue[Me],
+            OptionValue[Mv1],
+            OptionValue[Mm],
+            OptionValue[Mv2],
+            OptionValue[md2GeV],
+            OptionValue[mu2GeV],
+            OptionValue[ms2GeV],
+            OptionValue[mcmc],
+            OptionValue[CKMTheta12],
+            OptionValue[CKMTheta13],
+            OptionValue[CKMTheta23],
+            OptionValue[CKMDelta],
+            OptionValue[PMNSTheta12],
+            OptionValue[PMNSTheta13],
+            OptionValue[PMNSTheta23],
+            OptionValue[PMNSDelta],
+            OptionValue[PMNSAlpha1],
+            OptionValue[PMNSAlpha2],
+            OptionValue[alphaEm0],
+            OptionValue[Mh]
+
+            (* @ModelName@ input parameters *)
+@setInputParameterArguments@
+        }] /. HoldPattern[OptionValue[param_]] :> param /.
         { p } /.
         FS@ModelName@GetSettings[handle] /.
         FS@ModelName@GetSMInputParameters[handle] /.
-        FS@ModelName@GetInputParameters[handle]];
+        FS@ModelName@GetInputParameters[handle]]];

templates/mass_eigenstates.cpp.in

@@ -573,7 +573,7 @@ void CLASSNAME::calculate_DRbar_parameters()
 void CLASSNAME::calculate_pole_masses()
 {
 #ifdef ENABLE_THREADS
-   CLASSNAME* obj_ptr = this;
+   auto obj_ptr = this;
 
 @callAllLoopMassFunctionsInThreads@
 #else

templates/module.mk.in

@@ -236,6 +236,7 @@ pack-$(MODNAME)-src:
 		tar -czf $(@CLASSNAME@_TARBALL) \
 		$(LIB@CLASSNAME@_SRC) $(LIB@CLASSNAME@_HDR) \
 		$(EXE@CLASSNAME@_SRC) \
+		$(LL@CLASSNAME@_SRC) $(LL@CLASSNAME@_MMA) \
 		$(@CLASSNAME@_MK) $(@CLASSNAME@_TWO_SCALE_MK) \
 		$(@CLASSNAME@_SLHA_INPUT) $(@CLASSNAME@_GNUPLOT)
 

test/test_CMSSM_QedQcd_no_convergence.sh

@@ -16,16 +16,18 @@ spinfo_4=$(echo "$slha_out" | \
                   awk -f "$print_block" -v block=SPINFO | \
                   awk '{ if ($1 == 4) { $1 = ""; print $0 } }')
 
-error=1
+echo "SPINFO[4]: $spinfo_4"
+
+error=0
 
 case "$spinfo_4" in
-    *no\ convergence*) error=0 ;;
+    *SM\(5\)*) error=1 ;;
 esac
 
 if [ $error -eq 0 ] ; then
     echo "Test result: OK"
 else
-    echo "Error: point has converged"
+    echo "Error: point failed in determination of SM(5) parameters."
     echo ""
     echo "Test result: FAIL"
 fi

test/test_NMSSM_beta_functions.cpp

@@ -124,7 +124,7 @@ void test_beta_function_equality(const SoftParsNmssm& a, const NMSSM_soft_parame
    TEST_CLOSE(beta_a.displayMh2Squared(), beta_b.get_mHu2(), 2.0e-12);
    TEST_EQUALITY(beta_a.displayMsSquared(), beta_b.get_ms2());
    TEST_EQUALITY(beta_a.displaySoftMassSquared(mQl), beta_b.get_mq2());
-   TEST_EQUALITY(beta_a.displaySoftMassSquared(mUr), beta_b.get_mu2());
+   TEST_CLOSE(beta_a.displaySoftMassSquared(mUr), beta_b.get_mu2(), 2.0e-12);
    TEST_EQUALITY(beta_a.displaySoftMassSquared(mDr), beta_b.get_md2());
    TEST_EQUALITY(beta_a.displaySoftMassSquared(mLl), beta_b.get_ml2());
    TEST_EQUALITY(beta_a.displaySoftMassSquared(mEr), beta_b.get_me2());