SoftSUSY is hosted by Hepforge, IPPP Durham
SOFTSUSY  4.1
softsusy.h
Go to the documentation of this file.
1 
15 #ifndef SOFTSUSY_H
16 #define SOFTSUSY_H
17 
18 #ifdef HAVE_CONFIG_H
19 #include <config.h>
20 #endif
21 
22 #include <iostream>
23 #include <fstream>
24 #include <sstream>
25 #include <cstring>
26 #include <cstdlib>
27 #include <cmath>
28 #include "def.h"
29 #include "utils.h"
30 #include "numerics.h"
31 #include "physpars.h"
32 #include "lowe.h"
33 #include "softpars.h"
34 #include "twoloophiggs.h"
35 #include "mssmUtils.h"
36 #include "higher_order.h"
37 
38 #define HR "----------------------------------------------------------"
39 
40 #ifdef COMPILE_TWO_LOOP_GAUGE_YUKAWA
41 #include "mssm_twoloop_mt.h"
42 #include "mssm_twoloop_mb.h"
43 #include "mssm_twoloop_mtau.h"
44 #include "mssm_twoloop_as.h"
45 namespace softsusy { class MssmSoftsusy; }
46 #endif
47 
48 #ifdef ENABLE_HIMALAYA
49 #include "HierarchyCalculator.hpp"
50 #endif
51 
52 namespace softsusy {
53 
55  const double mxDefault = 1.9e16;
56 
58  class AltEwsbMssm {
59  private:
62  double mAcond, muCond;
63  public:
64  void setAltEwsb(double ma, double mu) { mAcond = ma; muCond = mu; };
65  void setAltEwsbMssm(const AltEwsbMssm & s) { *this = s; };
66 
67  // int displayConditionStyle() { return conditionStyle; };
68  double displayMaCond() const { return mAcond; };
69  double displayMuCond() const { return muCond; };
70  const AltEwsbMssm & displayAltEwsbMssm() const { return *this; };
71 
72  void setMaCond(double maInput) { mAcond = maInput; };
73  void setMuCond(double muInput) { muCond = muInput; };
74  };
75 
77  class MssmSoftsusy: public MssmSusy, public MssmSoftPars, public AltEwsbMssm,
78  public Approx, public RGE {
88  private:
89  sPhysical physpars;
90  drBarPars forLoops;
91  sProblem problem;
92  double msusy;
93  double minV;
94  double mw;
95  QedQcd dataSet;
96  double fracDiff;
97  bool setTbAtMX;
98  bool altEwsb;
100  double predMzSq;
101  double t1OV1Ms, t2OV2Ms;
102  double t1OV1Ms1loop, t2OV2Ms1loop;
103  double qewsb;
108 
109  protected:
110  void setT1OV1Ms(double t1) { t1OV1Ms = t1; }
111  void setT2OV2Ms(double t2) { t2OV2Ms = t2; }
112  void setT1OV1Ms1loop(double t1) { t1OV1Ms1loop = t1; }
113  void setT2OV2Ms1loop(double t2) { t2OV2Ms1loop = t2; }
114  double mxBC;
115 
116  void check_flags();
117  DoubleMatrix calcHiggs3L(bool is_bottom);
118  public:
121 
122  // void (*boundaryCondition)(Softsusy &, const DoubleVector &);
124  MssmSoftsusy();
126  MssmSoftsusy(const MssmSusyRGE &);
128 
129  MssmSoftsusy(const MssmSoftsusy &);
134  MssmSoftsusy(const MssmSusy &, const MssmSoftPars & s,
135  const sPhysical & sp, double mu, int l,
136  int t, double hv);
137 
139  const MssmSoftsusy & operator=(const MssmSoftsusy & s);
140 
142  virtual const DoubleVector display() const;
143 
145  const MssmSoftsusy & displayMssmSoft() const;
146 
148  const sPhysical & displayPhys() const;
149 
152  const drBarPars & displayDrBarPars() const;
154  const sProblem & displayProblem() const {return problem; };
156  const QedQcd & displayDataSet() const;
158  double displaySoftA(trilinears, int, int) const;
159 
161  double displayFracDiff() const { return fracDiff; };
162  double displayMinpot() const;
163  double displayMsusy() const;
164  double displayMw() const;
165  double displayMwRun() const;
168  double displayMzRun() const;
169  double displayTadpole1Ms() const;
170  double displayTadpole2Ms() const;
171  double displayTadpole1Ms1loop() const;
172  double displayTadpole2Ms1loop() const;
173  double displayMxBC() const { return mxBC; };
174  const MssmSoftsusy & displaySoftsusy() const { return *this; }
177  double displayMz() const { return displayDataSet().displayMu(); }
179  bool displaySetTbAtMX() const { return setTbAtMX; }
180  bool displayAltEwsb() const { return altEwsb; }
181  double displayPredMzSq() const { return predMzSq; }
182  double displayQewsb() const { return qewsb; }
183 
185  void flagMgutOutOfBounds(bool a) { problem.mgutOutOfBounds = a; };
187  void flagIrqfp(bool a) { problem.irqfp = a; };
189  void flagNonperturbative(bool a) { problem.nonperturbative = a; };
191  void flagTachyon(tachyonType a) {
194  if (close(displayMu(), MZ, 1.0e-6)) problem.tachyonWarning = a;
195  else { problem.tachyon = a;
196  if (PRINTOUT > 2) cout << tachyonNames[a] << " tachyon ";
197  }
198  };
199  void flagTachyonWarning(tachyonType a) { problem.tachyonWarning = a; }
201  void flagM3sq(bool a) { problem.m3sq = a; };
203  void flagNoConvergence(bool a) { problem.noConvergence = a; };
205  void flagNoMuConvergence(bool a) { problem.noMuConvergence = a; };
207  void flagNoRhoConvergence(bool a) { problem.noRhoConvergence = a; };
209  void flagMusqwrongsign(bool a) { problem.muSqWrongSign = a; };
212  void flagInaccurateHiggsMass(bool a) { problem.inaccurateHiggsMass = a; };
214  void flagHiggsufb(bool a) { problem.higgsUfb = a; };
216  void flagNotGlobalMin(bool a) { problem.notGlobalMin = a; };
218  void flagAllProblems(bool a, tachyonType b) { problem.irqfp = a;
219  problem.tachyon = b; problem.tachyonWarning = b;
220  problem.m3sq = a; problem.badConvergence = a;
221  problem.noConvergence = a; problem.higgsUfb = a; problem.notGlobalMin = a;
222  problem.nonperturbative = a; problem.noRhoConvergence = a;
223  problem.noMuConvergence = a; problem.muSqWrongSign = a;
224  problem.inaccurateHiggsMass = b; problem.mgutOutOfBounds = a; }
226  void flagProblemThrown(bool a) { problem.problemThrown = a; }
228  void setProblem(const sProblem a) { problem = a; }
229 
231  void setSoftsusy(const MssmSoftsusy & s) { *this = s; };
233  void setData(const QedQcd & r) { dataSet = r; };
235  void setMinpot(double);
238  void setMsusy(double);
240  void setMw(double);
242  void setPhys(const sPhysical & s) { physpars = s; };
244  void setMxBC(double mx) { mxBC = mx; };
246  void setDrBarPars(const drBarPars & s) { forLoops = s; };
248  void setSetTbAtMX(bool a) { setTbAtMX = a; }
250  void useAlternativeEwsb() { altEwsb = true; }
252  void setPredMzSq(double a) { predMzSq = a; }
254  void set(const DoubleVector &);
256  void setQewsb(double q) { qewsb = q; }
257 
260  void setTwoLoopAlphasThresholds(bool sw);
262  void setTwoLoopMtThresholds(bool sw);
268  void setTwoLoopMbMtauThresholds(bool sw);
270  void setAllTwoLoopThresholds(bool sw);
271 
273  DoubleVector beta() const {
274  sBrevity a;
275  DoubleVector y(MssmSusy::beta(a).display());
276  DoubleVector x(MssmSoftPars::beta(displayMssmSusy()));
277  int i; for (i=1; i<=y.displayEnd(); i++) x(i) = y(i);
278  return x;
279  };
280 
282  void setFracDiff(double fD) { fracDiff = fD; };
284  //for use in doCalcTadpole1oneLoop
285  void H1SfSfCouplings(DoubleMatrix & lTS1Lr, DoubleMatrix & lBS1Lr,
286  DoubleMatrix & lTauS1Lr, double gmzOcthW, double mu,
287  double cosb, double v1) const;
289  //for use in doCalcTadpole1oneLoop
290  void H2SfSfCouplings(DoubleMatrix & lTS2Lr, DoubleMatrix & lBS2Lr,
291  DoubleMatrix & lTauS2Lr, double gmzOcthW, double mu,
292  double sinb) const;
294  double doCalcTad1Sfermions(DoubleMatrix lTS1Lr, DoubleMatrix lBS1Lr,
295  DoubleMatrix lTauS1Lr, double costhDRbar) const;
297  double doCalcTad2Sfermions(DoubleMatrix lTS2Lr, DoubleMatrix lBS2Lr,
298  DoubleMatrix lTauS2Lr, double costhDRbar) const;
300  //for use in doCalcTadpole1oneLoop
301  double doCalcTad1fermions(double q, double v1) const;
303  //for use in doCalcTadpole1oneLoop
304  double doCalcTad2fermions(double q) const;
306  // Follwing BPMZ Goldstone bosons are not included in this.
307  double doCalcTad1Higgs(double q, double costhDRbar2, double g,
308  double tanb) const;
309  //one loop H2 tadpole contributions from Higgs bosons in the loops
310  // Follwing BPMZ Goldstone bosons are not included in this.
311  double doCalcTad2Higgs(double q, double costhDRbar2, double g,
312  double tanb) const;
314  double doCalcTad1Neutralinos(double q, double costhDRbar, double g,
315  double cosb) const;
317  double doCalcTad2Neutralinos(double q, double costhDRbar, double g,
318  double sinb) const;
320  double doCalcTad1Charginos(double q, double g, double cosb) const;
321 
322  double doCalcTad2Charginos(double q, double g, double sinb) const;
323 
325  double doCalcTad1GaugeBosons(double q, double costhDRbar2, double g,
326  double tanb) const;
327 
328  double doCalcTad2GaugeBosons(double q, double costhDRbar2,
329  double g, double tanb) const;
331  void doTadpoles(double mt, double sinthDRbar);
333  virtual double doCalcTadpole1oneLoop(double mt, double sinthDRbar);
335  virtual double doCalcTadpole2oneLoop(double mt, double sinthDRbar);
338  virtual void calcTadpole1Ms1loop(double mt, double sinthDRbar);
341  virtual void calcTadpole2Ms1loop(double mt, double sinthDRbar);
345  DoubleMatrix addStopQCD(double p, double mt, DoubleMatrix & strong);
349  DoubleMatrix addStopStop(double p, double mt, DoubleMatrix & stop);
353  DoubleMatrix addStopSbottom(double p, double mt, DoubleMatrix & sbottom);
357  DoubleMatrix addStopHiggs(double p, double mt, DoubleMatrix & higgs);
361  DoubleMatrix addStopEweak(double p, DoubleMatrix & electroweak);
365  DoubleMatrix addStopNeutralino(double p, double mt,
366  DoubleMatrix & neutralino);
370  DoubleMatrix addStopChargino(double p, DoubleMatrix & chargino);
374  virtual void addStopCorrection(double p, DoubleMatrix & mass, double mt);
378  DoubleMatrix addSdownQCD(double p1, double p2, int family,
379  DoubleMatrix & strong);
383  DoubleMatrix addSdownHiggs(double p1, double p2, int family,
384  DoubleMatrix & higgs);
388  DoubleMatrix addSdownEweak(double p1, double p2, int family,
389  DoubleMatrix & electroweak);
393  DoubleMatrix addSdownNeutralino(double p1, double p2, int family,
394  DoubleMatrix & neutralino);
398  DoubleMatrix addSdownChargino(double p1, double p2, int family,
399  DoubleMatrix & chargino);
403  virtual void addSdownCorrection(DoubleMatrix & mass, int family);
407  DoubleMatrix addSbotQCD(double p, double mt, DoubleMatrix & strong);
411  void addSbotSfermion(double p, double mt, DoubleMatrix & stop,
412  DoubleMatrix & sbottom);
416  DoubleMatrix addSbotHiggs(double p, double mt, DoubleMatrix & higgs);
420  DoubleMatrix addSbotEweak(double p, DoubleMatrix & electroweak);
424  DoubleMatrix addSbotNeutralino(double p, double mt,
425  DoubleMatrix & neutralino);
429  DoubleMatrix addSbotChargino(double p, double mt, DoubleMatrix & chargino);
433  virtual void addSbotCorrection(double p, DoubleMatrix & mass, double mb);
437  DoubleMatrix addSlepHiggs(double p1, double p2, int family,
438  DoubleMatrix & higgs);
442  DoubleMatrix addSlepEweak(double p1, double p2, int family,
443  DoubleMatrix & electroweak);
447  void addSlepGaugino(double p1, double p2, int family,
448  DoubleMatrix & chargino, DoubleMatrix & neutralino);
452  virtual void addSlepCorrection(DoubleMatrix & mass, int family);
455  void addSlepCorrection(double p, DoubleMatrix & mass, int family);
459  void addStauSfermion(double p, double mtau, DoubleMatrix & stop,
460  DoubleMatrix & sbottom);
464  void addStauGaugino(double p, double mtau, DoubleMatrix & chargino,
465  DoubleMatrix & neutralino);
469  DoubleMatrix addStauEweak(double p, double mtau, DoubleMatrix & electroweak);
473  DoubleMatrix addStauHiggs(double p, double mtau, DoubleMatrix & higgs);
477  virtual void addStauCorrection(double p, DoubleMatrix & mass, double mtau);
481  DoubleMatrix addSupQCD(double p1, double p2, int family,
482  DoubleMatrix & strong);
486  DoubleMatrix addSupHiggs(double p1, double p2, int family,
487  DoubleMatrix & higgs);
491  DoubleMatrix addSupEweak(double p1, double p2, int family,
492  DoubleMatrix & electroweak);
496  DoubleMatrix addSupNeutralino(double p1, double p2, int family,
497  DoubleMatrix & neutralino);
501  DoubleMatrix addSupChargino(double p1, double p2, int family,
502  DoubleMatrix & chargino);
506  virtual void addSupCorrection(DoubleMatrix & mass, int family);
509  void addSnuTauSfermion(double p, double & stop, double & sbottom);
512  double addSnuTauHiggs(double p, double & higgs);
515  double addSnuTauEweak(double p, double & electroweak);
519  void addSnuTauGaugino(double p, double & chargino, double & neutralino);
523  virtual void addSnuTauCorrection(double & mass);
527  double addSnuHiggs(double p, int family, double & higgs);
531  double addSnuEweak(double p, int family, double & electroweak);
535  void addSnuGaugino(double p, int family, double & chargino,
536  double & neutralino);
540  virtual void addSnuCorrection(double & mass, int family);
545  void addSquarkCorrection(DoubleMatrix & mass);
550  virtual void doUpSquarks(double mt, double pizztMS, double sinthDRbarMS, int
551  accuracy);
556  virtual void doDownSquarks(double mb, double pizztMS, double sinthDRbarMS,
557  int accuracy, double mt);
562  virtual void doChargedSleptons(double mT, double pizztMS,
563  double sinthDRbarMS, int accuracy);
566  virtual void doSnu(double pizztMS, int accuracy = 0);
572  virtual void treeUpSquark(DoubleMatrix & mass, double mtrun, double pizztMS,
573  double sinthDRbarMS, int family);
579  virtual void treeDownSquark(DoubleMatrix & mass, double mbrun, double pizztMS,
580  double sinthDRbarMS, int family);
586  virtual void treeChargedSlepton(DoubleMatrix & mass, double mTrun,
587  double pizztMS, double sinthDRbarMS,
588  int family);
591  void treeSnu(double & mSq, double pizztMS, int family);
592 
594  virtual void calcDrBarPars();
597  void setNeutCurrCouplings(double sinthDRbar, double & sw2, double & guL,
598  double & gdL, double & geL, double & guR,
599  double & gdR, double & geR );
601  void calcDRTrilinears(drBarPars & eg, double vev, double beta);
602 
603  void calcDrBarHiggs(double beta, double mz2, double mw2,
604  double sinthDRbar, drBarPars & eg);
607  void treeCharginos(DoubleMatrix & mCh, double beta, double mw);
610  void treeNeutralinos(DoubleMatrix & mN, double beta, double mz,
611  double mw, double sinth);
612  // calculates the chargino and neutralino DRbar parameters.
613  //It will fill in the chargino and neutralino
615  void calcDrBarGauginos(double beta, double mw, double mz, double sinth,
616  drBarPars & eg);
620  virtual void sparticleThresholdCorrections(double tb);
621  // const DoubleVector & pars);
625  double qedSusythresh(double alphaEm, double Q) const;
629  double qcdSusythresh(double alphasMSbar, double Q);
633  double calcMs() const;
636  virtual void physical(int accuracy);
638  //to the pole mt for use in calcRunningMt
639  virtual double calcRunMtQCD() const;
641  //to the pole mt for use in calcRunningMt
642  virtual double calcRunMtStopGluino() const;
644  //to the pole mt for use in calcRunningMt
645  virtual double calcRunMtHiggs() const;
646 
648  //to the pole mt for use in calcRunningMt
649  virtual double calcRunMtNeutralinos() const;
651  //to the pole mt for use in calcRunningMt
652  virtual double calcRunMtCharginos() const;
653 
656  virtual double calcRunningMt();
658  virtual double calcRunMtauDrBarConv() const;
659  // Obtains (1 / mTAU) times 1-loop squark-chargino corrections
660  //to mtau for use in calcRunningMtau
661  virtual double calcRunMtauCharginos(double mTauSMMZ) const;
662  // Obtains (1 / mTAU) times 1-loop squark-neutralino corrections
663  //to mtau for use in calcRunningMtau
664  virtual double calcRunMtauNeutralinos(double mTauSMMZ) const;
666  // to mtau for use in calcRunningMtau
667  virtual double calcRunMtauHiggs() const;
670  virtual double calcRunningMtau() ;
671 
673  virtual double calcRunMbDrBarConv() const;
675  //to mb for use in calcRunningMb
676  virtual double calcRunMbSquarkGluino() const;
678  //to mb for use in calcRunningMb
679  virtual double calcRunMbChargino() const;
681  // to mb for use in calcRunningMb
682  virtual double calcRunMbHiggs() const;
684  // to mb for use in calcRunningMb
685  virtual double calcRunMbNeutralinos() const;
688  // rruiz: remove const
689  virtual double calcRunningMb();
690  /*
693  double calcHt(double vev);
696  double calcHb(double vev) const;
699  double calcHtau(double vev) const;
700  */
702  double calcSinthdrbar() const;
704  double getVev();
707  double getVev(double pizzt);
711  virtual void charginos(int accuracy, double piwwt);
714  virtual void addChaLoopSfermion(double p, DoubleMatrix & sigmaL, DoubleMatrix & sigmaR, DoubleMatrix & sigmaS) const;
717  virtual void addChaLoopGauge(double p, DoubleMatrix & sigmaL,
718  DoubleMatrix & sigmaR, DoubleMatrix & sigmaS,
719  DoubleMatrix b1pCha, DoubleMatrix b0pCha,
720  DoubleMatrix b1pNeut,
721  DoubleMatrix b0pNeut) const;
724  virtual void addChaLoopHiggs(double p, DoubleMatrix & sigmaL,
725  DoubleMatrix & sigmaR, DoubleMatrix & sigmaS,
726  DoubleMatrix b1pCha, DoubleMatrix b0pCha,
727  DoubleMatrix b1pNeut,
728  DoubleMatrix b0pNeut) const;
730  virtual void addCharginoLoop(double p, DoubleMatrix &);
735  virtual void neutralinos(int accuracy, double piwwt, double pizzt);
739  void addNeutLoopSfermion(double p, DoubleMatrix & sigmaL,
740  DoubleMatrix & sigmaR, DoubleMatrix & sigmaS);
742  void getNeutPassarinoVeltman(double p, double q, DoubleMatrix & b0fn,
743  DoubleMatrix & b1fn);
747  void addNeutLoopGauge(double p, DoubleMatrix & sigmaL,
748  DoubleMatrix & sigmaR, DoubleMatrix & sigmaS);
752  void addNeutLoopHiggs(double p, DoubleMatrix & sigmaL,
753  DoubleMatrix & sigmaR, DoubleMatrix & sigmaS);
755  virtual void addNeutralinoLoop(double p, DoubleMatrix &);
757  virtual void gluino(int accuracy);
760  void calcHiggsAtScale(int accuracy, double & mt, double & sinthDRbar,
761  double & piwwtMS, double & pizztMS, double q = 0.);
768  bool higgs(int accuracy, double piwwt, double pizzt);
774  // virtual void newhiggs(int accuracy, double piwwt, double pizzt);
778  virtual int rewsbMu(int sgnMu, double & mu) const;
781  virtual int rewsbM3sq(double, double &) const;
784  int rewsbM3sqTree(double, double &) const;
786  void alternativeEwsb(double mt);
794  virtual void rewsb(int sgnMu, double mt, const DoubleVector & pars,
795  double muOld = -6.66e66, double eps = 0.);
798  virtual void rewsbTreeLevel(int sgnMu);
801  double treeLevelMuSq();
811  void iterateMu(double & munew, int sgnMu, double mt,
812  int maxTries, double pizztMS, double sinthDRbar, double tol,
813  int & err);
816  double predTanb(double muSusy = -6.66e66) const;
819  double predMzsq(double & tanb, double muOld = -6.66e66, double eps = 0.);
821  double deltaEW() const;
822 
830  DoubleVector fineTune(void (*boundaryCondition)(MssmSoftsusy &,
831  const DoubleVector &),
832  const DoubleVector & bcPars, double MX,
833  bool doTop = false);
837  double it1par(int numPar, const DoubleVector & bcPars);
841  double ufb3sl(double);
844  double realMinMs() const;
845 
847  //for p=external momentum, Q=renormalisation scale
848  virtual double piZZTHiggs(double p, double Q, double thetaWDRbar) const;
850  //for p=external momentum, Q=renormalisation scale
851  virtual double piZZTsfermions(double p, double Q) const;
854  virtual double piZZTfermions(double p, double Q, bool usePoleMt) const;
857  virtual double piZZTNeutralinos(double p, double Q, double thetaWDRbar) const;
860  virtual double piZZTCharginos(double p, double q, double thetaWDRbar) const;
863  virtual double piZZT(double p, double Q, bool usePoleMt = false) const;
866  virtual double piWWTHiggs(double p, double q, double thetaWDRbar) const;
869  virtual double piWWTfermions(double p, double Q, bool usePoleMt) const;
872  virtual double piWWTsfermions(double p, double Q) const;
875  virtual double piWWTgauginos(double p, double Q, double thetaWDRbar) const;
878  virtual double piWWT(double p, double Q, bool usePoleMt = false) const;
881  virtual void getNeutralinoCharginoHpmCoup(ComplexMatrix & apph1,
882  ComplexMatrix & apph2,
883  ComplexMatrix & bpph1,
884  ComplexMatrix & bpph2) const;
887  double piHpHmFermions(double p, double q) const;
891  double piHpHmSfermions(double p, double q, double mu) const;
894  double piHpHmGauge(double p, double q) const;
897  double piHpHmHiggs(double p, double q) const;
900  double piHpHmGauginos(double p, double q) const;
903  virtual double piHpHm(double p, double Q) const;
906  double piZGT(double p, double Q) const;
909  virtual double piAA(double p, double Q) const;
911  //self-energy: for p=external momentum, q=renormalisation scale
912  Complex pis1s1Sfermions(double p, double q, DoubleMatrix ls1tt,
913  DoubleMatrix ls1bb, DoubleMatrix ls1tautau) const;
916  Complex pis1s2Sfermions(double p, double q, DoubleMatrix ls1tt,
917  DoubleMatrix ls1bb, DoubleMatrix ls1tautau,
918  DoubleMatrix ls2tt, DoubleMatrix ls2bb,
919  DoubleMatrix ls2tautau) const;
922  Complex pis2s2Sfermions(double p, double q, DoubleMatrix ls2tt,
923  DoubleMatrix ls2bb, DoubleMatrix ls2tautau) const;
926  Complex pis1s1Fermions(double p, double q) const;
929  Complex pis2s2Fermions(double p, double q) const;
932  Complex pis1s1Higgs(double p, double q) const;
935  Complex pis1s2Higgs(double p, double q) const;
938  Complex pis2s2Higgs(double p, double q) const;
941  Complex pis1s1Neutralinos(double p, double q) const;
944  Complex pis1s2Neutralinos(double p, double q) const;
947  Complex pis2s2Neutralinos(double p, double q) const;
950  Complex pis1s1Charginos(double p, double q) const;
953  Complex pis1s2Charginos(double p, double q) const;
956  Complex pis2s2Charginos(double p, double q) const;
959  Complex pis1s1(double p, double q) const;
962  Complex pis1s2(double p, double q) const;
965  Complex pis2s2(double p, double q) const;
967  double sinSqThetaEff();
969  virtual double h1s2Mix();
978  virtual void rhohat(double & outrho, double & outsin, double alphaMZDRbar,
979  double pizztMZ, double piwwt0, double piwwtMW,
980  double tol, int maxTries);
984  virtual double deltaVb(double outrho, double outsin, double alphaDRbar,
985  double pizztMZ) const;
989  virtual double dRho(double outrho, double outsin, double alphaDRbar,
990  double pizztMZ, double piwwtMW);
994  virtual double dR(double outrho, double outsin, double alphaDRbar,
995  double pizztMZ, double piwwt0);
997  double maxMass() const;
1002  int lsp(double & mass, int & posi, int & posj) const;
1007  int nlsp(double & mass, int & posi, int & posj) const;
1008 
1010  string printShort() const;
1012  string printLong();
1013 
1015  virtual void printObj() {
1016  cout << this->displayMssmSusy() << this->displayMssmSoftPars();
1017  };
1018 
1020  double thet(double a, double b, double c = 0.0);
1021 
1040  void fixedPointIteration(void (*boundaryCondition)
1041  (MssmSoftsusy &, const DoubleVector &),
1042  double mxGuess,
1043  const DoubleVector & pars, int sgnMu, double tanb,
1044  const QedQcd & oneset, bool gaugeUnification,
1045  bool ewsbBCscale = false);
1049  double lowOrg(void (*boundaryCondition)
1050  (MssmSoftsusy &, const DoubleVector &),
1051  double mxGuess,
1052  const DoubleVector & pars, int sgnMu, double tanb,
1053  const QedQcd & oneset, bool gaugeUnification,
1054  bool ewsbBCscale = false) {
1055  fixedPointIteration(boundaryCondition, mxGuess, pars, sgnMu, tanb, oneset,
1056  gaugeUnification, ewsbBCscale);
1057  return displayMxBC();
1058  };
1059 
1069 
1070  void itLowsoft(int maxTries, int sgnMu, double tol,
1071  double tanb, void (*boundaryCondition)(MssmSoftsusy &,
1072  const DoubleVector &),
1073  const DoubleVector & pars, bool gaugeUnification,
1074  bool ewsbBCscale);
1075 
1079  virtual void rpvSet(const DoubleVector & parameters);
1080 
1084  virtual void methodBoundaryCondition(const DoubleVector & pars);
1085 
1102  void isajetNumbers764
1103  (double & mtopPole, double & mGPole, double & smu, double & mA,
1104  double & tanb, double & mq1l, double & mdr, double & mur, double & meL,
1105  double & meR, double & mql3, double & mdr3, double & mur3, double & mtauL,
1106  double & mtauR, double & at, double & ab, double & atau, double & mq2l,
1107  double & msr, double & mcr, double & mmuL, double & mmuR, double & m1,
1108  double & m2) const;
1110  void isajetInterface764(const char fname[80]) const;
1113  void ssrunInterface764Inside(const char fname [80], fstream & ) const;
1117  void ssrunInterface764(const char fname [80], const char softfname [80]) const;
1121  void isawigInterface764(const char fnamein [80], const char fnameout [80],
1122  const char fnamesoft[80]) const;
1135  virtual void lesHouchesAccordOutput(ostream & out, const char model[],
1136  const DoubleVector & pars,
1137  int sgnMu, double tanb, double qMax,
1138  int numPoints,
1139  bool ewsbBCscale);
1140  /* void slha1(ostream & out, const char model[], const DoubleVector & pars,
1141  int sgnMu, double tanb, double qMax, int numPoints,
1142  bool ewsbBCscale);*/
1146  virtual void setEwsbConditions(const DoubleVector & inputs);
1149  void headerSLHA(ostream & out);
1151  virtual void spinfoSLHA(ostream & out);
1153  void modselSLHA(ostream & out, const char model[]);
1155  void sminputsSLHA(ostream & out);
1157  void minparSLHA(ostream & out, const char model [],
1158  const DoubleVector & pars, double tanb, int sgnMu,
1159  bool ewsbBCscale);
1161  virtual void extparSLHA(ostream & out, const DoubleVector & pars,
1162  bool ewsbBCscale);
1164  void massSLHA(ostream & out);
1166  virtual void higgsMSLHA(ostream & out);
1168  virtual void neutralinoCharginoMSLHA(ostream & out);
1170  virtual void sfermionsSLHA(ostream & out);
1172  void sfermionmixSLHA(ostream & out);
1174  virtual void neutralinoMixingSLHA(ostream & out);
1176  void inomixingSLHA(ostream & out);
1178  virtual void softsusySLHA(ostream & out);
1180  void alphaSLHA(ostream & out);
1182  virtual void hmixSLHA(ostream & out);
1184  void gaugeSLHA(ostream & out);
1186  virtual void yukawasSLHA(ostream & out);
1188  virtual void msoftSLHA(ostream & out);
1190  virtual void drbarSLHA(ostream & out, int numPoints, double qMax, int n);
1191 
1195  void doUfb3(double mgut);
1196 
1200  void assignHiggs(DoubleVector & higgsm, DoubleVector & higgsc,
1201  DoubleVector & dnu, DoubleVector & dnd,
1202  DoubleVector & cn, double beta) const;
1203 
1206  void assignHiggs(DoubleVector & higgsm, DoubleVector & higgsc) const;
1207 
1212  void assignHiggsSfermions(DoubleVector & higgsm, DoubleVector & higgsc,
1213  DoubleVector & dnu, DoubleVector & dnd,
1214  DoubleVector & cn, double beta) const;
1215 
1218  double smPredictionMW() const;
1219 
1222  double twoLoopGm2(double amu1Loop) const;
1223 
1226  double twoLpMt() const;
1227  double twoLpMb() const;
1228 
1233  virtual void doQuarkMixing(DoubleMatrix & mDon, DoubleMatrix & mUpq);
1234 
1239  virtual MssmSusyRGE guessAtSusyMt(double tanb, const QedQcd & oneset);
1240  };
1241 
1242  std::istream& operator>>(std::istream& left, MssmSoftsusy& s);
1243 
1245  void printShortInitialise();
1247  double rho2(double r);
1248 
1253  double minimufb3(double);
1254 
1256  double gEff(double x);
1258  double fEff(double x);
1259 
1261  double lep2Likelihood(double mh);
1264  DoubleVector mhIntegrand(double mh, const DoubleVector & y);
1267  double lnLHiggs(double mh);
1268 
1269 }
1270 
1271 #endif
1272 
1273 
1274 
global variable declaration
Definition: def.cpp:13
DoubleVector boundaryCondition(double m0, double m12, double a0, const RpvSusyPars &r)
Definition: rpvsusypars.cpp:373
void flagMgutOutOfBounds(bool a)
Flags weird mgut-type problems.
Definition: softsusy.h:185
DoubleVector is of variable length, and contains double precision.
Definition: linalg.h:35
void setDrBarPars(const drBarPars &s)
Sets tree-level DRbar parameters.
Definition: softsusy.h:246
Masses of the physical particles.
Definition: physpars.h:32
DRbar values of masses and mixings in MSSM.
Definition: physpars.h:114
bool inaccurateHiggsMass
Higgs mass is potentially inaccurate and cant be trusted.
Definition: physpars.h:87
int displayEnd() const
returns end of dimension
Definition: linalg.h:94
void setQewsb(double q)
Sets user-set scale qewsb.
Definition: softsusy.h:256
DoubleVector beta() const
Returns double vector containing numerical beta functions of parameters.
Definition: softsusy.h:273
switches (options) and parameters such as default fermion masses, required accuracy etc ...
double displayMz() const
Returns value of pole MZ being used.
Definition: softsusy.h:177
bool muSqWrongSign
mu^2 came out with wrong sign; no REWSB
Definition: physpars.h:80
double rho2(double r)
Two-loop Standard Model corrections to rho parameter.
Definition: softsusy.cpp:8979
void setProblem(const sProblem a)
Sets all problem flags.
Definition: softsusy.h:228
bool noRhoConvergence
Couldn&#39;t calculate electroweak rho parameter.
Definition: physpars.h:73
bool nonperturbative
Running went non-perturbative.
Definition: physpars.h:84
mass
used to give order of quark masses stored
Definition: lowe.h:46
double minimufb3(double lnH2)
Definition: softsusy.cpp:6728
QedQcd object contains Standard Model quark and lepton masses. It integrates them using 3 loop qcd x ...
bool displaySetTbAtMX() const
Is tan beta set at the user defined SUSY breaking scale?
Definition: softsusy.h:179
tachyonType tachyonWarning
Definition: physpars.h:79
Contains all supersymmetric RPC-MSSM parameters.
Definition: susy.h:75
void flagInaccurateHiggsMass(bool a)
Definition: softsusy.h:212
int PRINTOUT
no verbose debug output
Definition: def.cpp:15
bool irqfp
Infra-red quasi fixed point breached.
Definition: physpars.h:72
DoubleVector beta(const MssmSusy &) const
Returns double vector containing numerical beta functions of parameters.
Definition: softpars.cpp:1031
Two loop corrections to the top mass. This file has been generated at Fri 31 Mar 2017 14:09:29 with t...
Quark and lepton masses and gauge couplings in QEDxQCD effective theory.
Definition: lowe.h:55
bool badConvergence
Nowhere near a decent solution.
Definition: physpars.h:71
void flagIrqfp(bool a)
Flags Infra-red quasi fixed point breach problem.
Definition: softsusy.h:187
void flagTachyon(tachyonType a)
Flags a negative-mass squared scalar (really a CCB problem)
Definition: softsusy.h:191
numerical routines - differential equation solver, differentiator and function minimiser for instance...
A few handy bits and pieces - little mathematical functions and the like.
void setData(const QedQcd &r)
Sets low energy Standard Model fermion mass and gauge coupling data.
Definition: softsusy.h:233
int included_thresholds
Flag allowing to choose which two-loop thresholds have to be included.
Definition: softsusy.h:120
void setPhys(const sPhysical &s)
Sets all physical parameters.
Definition: softsusy.h:242
bool m3sq
m3sq came out with wrong sign; no REWSB
Definition: physpars.h:81
bool noConvergence
Iteration did not converge: not always serious.
Definition: physpars.h:74
void setPredMzSq(double a)
Set MZ^2 predicted after iteration.
Definition: softsusy.h:252
void setMxBC(double mx)
Sets the scale at which high-scale boundary conditions are applied.
Definition: softsusy.h:244
matrix of complex double values dimensions (rows x cols)
Definition: linalg.h:605
tachyonType tachyon
Definition: physpars.h:75
istream & operator>>(istream &left, flavourPhysical &s)
Formatted input of physical parameters.
Definition: flavoursoft.cpp:1220
double fEff(double x)
function used for calculating sin theta_eff
Definition: softsusy.cpp:8999
void flagAllProblems(bool a, tachyonType b)
Sets all problems equal to either true or false (contained in a)
Definition: softsusy.h:218
Contains data needed in beta function calculation to make it faster.
Definition: susy.h:35
double gEff(double x)
function used for calculating sin theta_eff
Definition: softsusy.cpp:9006
const sProblem & displayProblem() const
Returns any problem flags associated with the object.
Definition: softsusy.h:154
void flagNotGlobalMin(bool a)
LCT: Flags problem if not in global minimum of Higgs potential.
Definition: softsusy.h:216
void setFracDiff(double fD)
sets fracDiff, needed for access by NmssmSoftsusy methods
Definition: softsusy.h:282
Various boolean values to flag any problems with the parameter point.
Definition: physpars.h:69
double mxBC
Scale at which SUSY breaking boundary conditions set.
Definition: softsusy.h:114
double sw2
Some parameters used in the computation.
Definition: def.cpp:51
void flagM3sq(bool a)
Flags problem with Higgs potential minimum.
Definition: softsusy.h:201
Describes a set of parameters and RGEs in general.
Definition: rge.h:49
void flagHiggsufb(bool a)
Flags an inconsistent Higgs minimum.
Definition: softsusy.h:214
Matrix from 1..rows, 1..cols of double values.
Definition: linalg.h:211
A different REWSB condition - given by mu and MA instead of mh1,2.
Definition: softsusy.h:58
A few handy routines for the MSSM: SUSY breaking conditions etc.
virtual void printObj()
Prints whole object to standard output.
Definition: softsusy.h:1015
Contains two-loop routines from Slavich et al.
DoubleVector mhIntegrand(double mh, const DoubleVector &)
Definition: softsusy.cpp:10412
double displayMu() const
Return renomalisation scale.
Definition: rge.h:67
Flags for potential problems in sProblem structure, and structure for containing physical MSSM parame...
double displayFracDiff() const
Displays iteration accuracy attained.
Definition: softsusy.h:161
void flagNoMuConvergence(bool a)
Flags fact that mu sub iteration didn&#39;t converge.
Definition: softsusy.h:205
void setSoftsusy(const MssmSoftsusy &s)
Sets whole object equal to another.
Definition: softsusy.h:231
void flagMusqwrongsign(bool a)
Flags point inconsistent with electroweak symmetry breaking.
Definition: softsusy.h:209
Soft SUSY breaking parameters.
double m1
decay global variable declarations
Definition: decays.cpp:14
Definition: rge.h:23
Contains all supersymmetric MSSM parameters, incorporating R_p MSSM.
Definition: softsusy.h:77
bool close(double m1, double m2, double tol)
Definition: utils.cpp:68
bool notGlobalMin
Not in global minimum of Higgs potential.
Definition: physpars.h:83
const double mxDefault
< default SUSY breaking boundary condition scale
Definition: softsusy.h:55
void setSetTbAtMX(bool a)
Sets the setTbAtMX flag.
Definition: softsusy.h:248
double lowOrg(void(*boundaryCondition)(MssmSoftsusy &, const DoubleVector &), double mxGuess, const DoubleVector &pars, int sgnMu, double tanb, const QedQcd &oneset, bool gaugeUnification, bool ewsbBCscale=false)
Definition: softsusy.h:1049
trilinears
SUSY breaking trilinear parameters.
Definition: softpars.h:29
Soft SUSY breaking parameters and beta functions.
Definition: softpars.h:35
void printShortInitialise()
Prints out header line for print-short output.
Definition: softsusy.cpp:1263
MssmSusy beta(sBrevity &) const
Definition: susy.cpp:300
const double accuracy
Approximate accuracy with which 3 body decays are calculated.
Definition: threeBodyDecays.h:38
void flagNonperturbative(bool a)
Flags non-perturbative RG evolution.
Definition: softsusy.h:189
bool higgsUfb
Higgs potential inconsistent with a good minimum.
Definition: physpars.h:82
void useAlternativeEwsb()
Use alernative EWSB conditions: set mu and MA(pole)
Definition: softsusy.h:250
Contains all supersymmetric RPC-MSSM parameters.
Definition: susy.h:222
void flagNoConvergence(bool a)
Flags fact that calculation hasn&#39;t acheived required accuracy.
Definition: softsusy.h:203
double MZ
global pole mass of MZ in GeV - MZCENT is defined in def.h
Definition: def.cpp:31
void flagProblemThrown(bool a)
Flags a numerical exception eg number too big/small.
Definition: softsusy.h:226
void flagNoRhoConvergence(bool a)
Flags fact that rho parameter sub iteration didn&#39;t converge.
Definition: softsusy.h:207
Two loop corrections to the top mass. This file has been generated at Fri 14 Jul 2017 22:12:37 with t...
bool noMuConvergence
Definition: physpars.h:85
double lep2Likelihood(double mh)
Fit to LEP2 Standard Model results.
Definition: softsusy.cpp:10401
double lnLHiggs(double mh)
returns the smeared log likelihood coming from LEP2 Higgs mass bounds
Definition: softsusy.cpp:10421
drop-in replacement for the original home-grown Complex class
Definition: mycomplex.h:17