SUSYSpectrum.cpp
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2012 HEPfit Collaboration
3  * All rights reserved.
4  *
5  * For the licensing terms see doc/COPYING.
6  */
7 
8 #include <iostream>
9 #include <stdexcept>
10 #include <limits>
11 #include "SUSYSpectrum.h"
12 #include "SUSY.h"
13 
14 
16 : mySUSY(SUSY_in), Mchargino(2,2,0.), Mneutralino(4,4,0.),
17  U(2,2,0.), V(2,2,0.), N(4,4,0.),
18  Msup2(6,0.), Msdown2(6,0.), Msneutrino2(6,0.), Mselectron2(6,0.),
19  mch(2,0.), mneu(4,0.), m_su2(6,0.), m_sd2(6,0.), m_sn2(6,0.), m_se2(6,0.),
20  Ru(6,6,0.), Rd(6,6,0.), Rn(6,6,0.), Rl(6,6,0.)
21 {
22 }
23 
24 bool SUSYSpectrum::CalcHiggs(double mh[4], gslpp::complex& saeff_i)
25 {
26  double Mw = mySUSY.Mw_tree();
27  double Mz = mySUSY.getMz();
28  double cos2b = 2.0 * mySUSY.getCosb() * mySUSY.getCosb() - 1.0;
29 
30  /* charged Higgs mass */
31  mh[3] = mySUSY.mHptree;
32 
33  /* pseudo-scalar Higgs mass */
34  mh[2] = sqrt(mh[3] * mh[3] - Mw * Mw);
35 
36  double temp = mh[2] * mh[2] + Mz * Mz;
37  double temp1 = 2.0 * mh[2] * Mz * cos2b;
38  double temp2 = sqrt(temp * temp - temp1 * temp1);
39 
40  /* light Higgs mass */
41  mh[0] = sqrt((temp - temp2)/2.0);
42 
43  /* heavy Higgs mass */
44  mh[1] = sqrt((temp + temp2)/2.0);
45 
46  /* CP-even Higgs mixing angle*/
47  saeff_i = sin(atan((mh[2]*mh[2]+Mz*Mz)*sqrt(1-cos2b*cos2b)/((mh[2]*mh[2]-Mz*Mz)*cos2b))/2.0);
48 
49  return true;
50 }
51 
53 {
54  double Mw = mySUSY.Mw_tree();
55 
56  Mchargino.assign(0, 0, mySUSY.getM2());
57  Mchargino.assign(0, 1, sqrt(2) * Mw * mySUSY.getSinb());
58  //
59  Mchargino.assign(1, 0, sqrt(2) * Mw * mySUSY.getCosb());
60  Mchargino.assign(1, 1, mySUSY.getMuH());
61 
62  /*
63  * M.singularvalue(&U, &V, &S) decompose M into
64  * M = U diag(S) V^+.
65  */
66  gslpp::matrix<gslpp::complex> Utmp(2, 2, 0.), Vtmp(2, 2, 0.);
67  Mchargino.singularvalue(Utmp, Vtmp, mch_i);
68 
72  U_i = Utmp.transpose();
73  V_i = Vtmp.hconjugate();
74 
75  return true;
76 }
77 
79 {
80  double Mw = mySUSY.Mw_tree();
81  double Mz = mySUSY.getMz();
82  double cW2 = Mw*Mw/Mz/Mz;
83  double cW = sqrt(cW2);
84  double sW = sqrt(1.0 - cW2);
85  double sb = mySUSY.getSinb();
86  double cb = mySUSY.getCosb();
87 
88  Mneutralino.assign(0, 0, mySUSY.getM1());
89  Mneutralino.assign(0, 2, - Mz * sW * cb);
90  Mneutralino.assign(0, 3, Mz * sW * sb);
91  //
92  Mneutralino.assign(1, 1, mySUSY.getM2());
93  Mneutralino.assign(1, 2, Mz * cW * cb);
94  Mneutralino.assign(1, 3, - Mz * cW * sb);
95  //
96  Mneutralino.assign(2, 0, Mneutralino(0,2));
97  Mneutralino.assign(2, 1, Mneutralino(1,2));
98  Mneutralino.assign(2, 3, - mySUSY.getMuH());
99  //
100  Mneutralino.assign(3, 0, Mneutralino(0,3));
101  Mneutralino.assign(3, 1, Mneutralino(1,3));
102  Mneutralino.assign(3, 2, Mneutralino(2,3));
103 
104  /*
105  * M.singularvalue(&U, &V, &S) decompose M into
106  * M = U diag(S) V^+.
107  */
108  gslpp::matrix<gslpp::complex> Ntemp1(4, 4, 0.), Ntemp2(4, 4, 0.);
109  Mneutralino.singularvalue(Ntemp1, Ntemp2, mneu_i);
110 
111  //gslpp::matrix<gslpp::complex> Nleft = Ntemp1.transpose();
112  gslpp::matrix<gslpp::complex> Nright = Ntemp2.hconjugate();
113  //std::cout << "Nleft = " << Nleft << std::endl;
114  //std::cout << "Nright = " << Nright << std::endl;
115 
116  /* adopt N=Nright as N^* M N^+ = diag, not as N M N^+.
117  * As a result, a phase rotation is required for N. */
118  gslpp::matrix<gslpp::complex> Mdiag_tmp(4, 4, 0.);
119  Mdiag_tmp = Nright.hconjugate().transpose() * Mneutralino * Nright.hconjugate();
120  //std::cout << "Mdiag_tmp = " << Mdiag_tmp << std::endl;
122  for(int i = 0; i < 4; i++)
123  v1.assign(i, gslpp::complex(1., Mdiag_tmp(i,i).arg()/2.0, true));
124  Nright = gslpp::matrix<gslpp::complex>(v1) * Nright;
125 
126  N_i = Nright;
127 
128  return true;
129 }
130 
132 {
133  double Mw = mySUSY.Mw_tree();
134  double Mz2 = mySUSY.getMz()*mySUSY.getMz();
135  double sW2 = 1.0 - Mw*Mw/Mz2;
136  double cos2b = 2.0 * mySUSY.getCosb() * mySUSY.getCosb() - 1.0;
139 
140  /* DRbar up-type quark masses at scale Q */
141  gslpp::matrix<double> Mu(3, 3, 0.);
142  Mu(0, 0) = mySUSY.Mq_Q(mySUSY.UP);
143  Mu(1, 1) = mySUSY.Mq_Q(mySUSY.CHARM);
144  Mu(2, 2) = mySUSY.Mq_Q(mySUSY.TOP);
145 
146  gslpp::matrix<gslpp::complex> uLL( CKM * mySUSY.msQhat2 * CKM.hconjugate() + Mu * Mu
147  + cos2b * Mz2 * (1.0/2.0 - 2.0/3.0 * sW2) * Id3 );
148  gslpp::matrix<gslpp::complex> uLR( mySUSY.v2()/sqrt(2.0) * mySUSY.getTUhat().hconjugate()
149  - mySUSY.getMuH() * Mu / mySUSY.getTanb() );
150  gslpp::matrix<gslpp::complex> uRR( mySUSY.msUhat2 + Mu * Mu + cos2b * Mz2 * 2.0/3.0 * sW2 * Id3 );
151  for(int i = 0; i < 3; i++)
152  for(int j = 0; j < 3; j++) {
153  Msup2.assign(i, j, uLL(i,j));
154  Msup2.assign(i, j+3, uLR(i,j));
155  Msup2.assign(i+3, j, uLR(j,i).conjugate());
156  Msup2.assign(i+3, j+3, uRR(i,j));
157  }
158 
159  /*
160  * M.singularvalue(&U, &V, &S) decompose M into
161  * M = U diag(S) V^+.
162  */
163  gslpp::matrix<gslpp::complex> RuTmp(6, 6, 0.);
164  Msup2.eigensystem(RuTmp, m_su2_i);
165 
166  Ru_i = RuTmp.hconjugate();
167 
168  return true;
169 
170 }
171 
173 {
174  double Mw = mySUSY.Mw_tree();
175  double Mz2 = mySUSY.getMz()*mySUSY.getMz();
176  double sW2 = 1.0 - Mw*Mw/Mz2;
177  double cos2b = 2.0 * mySUSY.getCosb() * mySUSY.getCosb() - 1.0;
179 
180  /* DRbar down-type quark masses at scale Q */
181  gslpp::matrix<double> Md(3, 3, 0.);
182  Md(0, 0) = mySUSY.Mq_Q(mySUSY.DOWN);
183  Md(1, 1) = mySUSY.Mq_Q(mySUSY.STRANGE);
184  Md(2, 2) = mySUSY.Mq_Q(mySUSY.BOTTOM);
185 
187  + cos2b * Mz2 * (- 1.0/2.0 + 1.0/3.0 * sW2) * Id3 );
188  gslpp::matrix<gslpp::complex> dLR( mySUSY.v1()/sqrt(2.0) * mySUSY.getTDhat().hconjugate()
189  - mySUSY.getMuH() * Md * mySUSY.getTanb() );
190  gslpp::matrix<gslpp::complex> dRR( mySUSY.msDhat2 + Md * Md - cos2b * Mz2 /3.0 * sW2 * Id3 );
191  for(int i = 0; i < 3; i++)
192  for(int j = 0; j < 3; j++) {
193  Msdown2.assign(i, j, dLL(i,j));
194  Msdown2.assign(i, j+3, dLR(i,j));
195  Msdown2.assign(i+3, j, dLR(j,i).conjugate());
196  Msdown2.assign(i+3, j+3, dRR(i,j));
197  }
198 
199  /*
200  * M.singularvalue(&U, &V, &S) decompose M into
201  * M = U diag(S) V^+.
202  */
203  gslpp::matrix<gslpp::complex> RdTmp(6, 6, 0.);
204  Msdown2.eigensystem(RdTmp, m_sd2_i);
205 
206  Rd_i = RdTmp.hconjugate();
207 
208  return true;
209 
210 }
211 
213 {
214  double Mz2 = mySUSY.getMz()*mySUSY.getMz();
215  double cos2b = 2.0 * mySUSY.getCosb() * mySUSY.getCosb() - 1.0;
217 
218  // this section is useful to re-define the sneutrino matrix...
219 
220 // double delta12=0.;
221 // double delta13=0.1;
222 // double delta23=0.1;
223 // gslpp::complex sLmass=mySUSY.msLhat2(0,0);
224 // gslpp::matrix<gslpp::complex> msLhat2modified(3,3,0.);
225 // msLhat2modified.assign(0, 0, sLmass);
226 // msLhat2modified.assign(1, 1, sLmass);
227 // msLhat2modified.assign(2, 2, sLmass);
228 // msLhat2modified.assign(0, 1, 0.);
229 // msLhat2modified.assign(1, 0, 0.);
230 // msLhat2modified.assign(0, 2, delta13*sLmass);
231 // msLhat2modified.assign(2, 0, delta13*sLmass);
232 // msLhat2modified.assign(1, 2, delta23*sLmass);
233 // msLhat2modified.assign(2, 1, delta23*sLmass);
234 // gslpp::matrix<gslpp::complex> nLL( msLhat2modified + cos2b * Mz2 /2.0 * Id3 );
235 
236  // ... until here
237 
238  gslpp::matrix<gslpp::complex> nLL( mySUSY.msLhat2 + cos2b * Mz2 /2.0 * Id3 );
239  gslpp::matrix<gslpp::complex> nLR( mySUSY.v2()/sqrt(2.0) * mySUSY.getTNhat().hconjugate() );
241  for(int i = 0; i < 3; i++) {
242  for(int j = 0; j < 3; j++) {
243  Msneutrino2.assign(i, j, nLL(i,j));
244  Msneutrino2.assign(i, j+3, nLR(i,j));
245  Msneutrino2.assign(i+3, j, nLR(j,i).conjugate());
246  Msneutrino2.assign(i+3, j+3, nRR(i,j));
247  }
248  }
249  /*
250  * M.singularvalue(&U, &V, &S) decompose M into
251  * M = U diag(S) V^+.
252  */
253  gslpp::matrix<gslpp::complex> RnTmp(6, 6, 0.);
254  Msneutrino2.eigensystem(RnTmp, m_sn2_i);
255 
256  m_sn2_i(0)=std::numeric_limits<double>::max();
257  m_sn2_i(1)=std::numeric_limits<double>::max();
258  m_sn2_i(2)=std::numeric_limits<double>::max();
259  Rn_i = RnTmp.hconjugate();
260 
261  return true;
262 
263 }
264 
266 {
267  double Mw = mySUSY.Mw_tree();
268  double Mz2 = mySUSY.getMz()*mySUSY.getMz();
269  double sW2 = 1.0 - Mw*Mw/Mz2;
270  double cos2b = 2.0 * mySUSY.getCosb() * mySUSY.getCosb() - 1.0;
272 
273  gslpp::matrix<double> Me(3, 3, 0.);
274  Me(0, 0) = mySUSY.Ml_Q(mySUSY.ELECTRON);
275  Me(1, 1) = mySUSY.Ml_Q(mySUSY.MU);
276  Me(2, 2) = mySUSY.Ml_Q(mySUSY.TAU);
277 
278  // this section is useful to re-define the slepton matrix
279  //if you want to use deltas, use this...
280  double delta12=0.1;
281  //double delta13=0.;
282  //double delta23=0.;
283  gslpp::complex sLmass=mySUSY.msLhat2(0,0);
284  gslpp::matrix<gslpp::complex> msLhat2modified(3,3,0.);
285  msLhat2modified.assign(0, 0, sLmass);
286  msLhat2modified.assign(1, 1, sLmass);
287  msLhat2modified.assign(2, 2, sLmass);
288  msLhat2modified.assign(0, 1, delta12*sLmass);
289  msLhat2modified.assign(1, 0, delta12*sLmass);
290  msLhat2modified.assign(0, 2, 0.);
291  msLhat2modified.assign(2, 0, 0.);
292  msLhat2modified.assign(1, 2, 0.);
293  msLhat2modified.assign(2, 1, 0.);
294  gslpp::matrix<gslpp::complex> msEhat2modified(3,3,0.);
295  msEhat2modified.assign(0, 0, sLmass);
296  msEhat2modified.assign(1, 1, sLmass);
297  msEhat2modified.assign(2, 2, sLmass);
298  msEhat2modified.assign(0, 1, 0.);
299  msEhat2modified.assign(1, 0, 0.);
300  msEhat2modified.assign(0, 2, 0.);
301  msEhat2modified.assign(2, 0, 0.);
302  msEhat2modified.assign(1, 2, 0.);
303  msEhat2modified.assign(2, 1, 0.);
304  gslpp::matrix<gslpp::complex> eLL( msLhat2modified + Me * Me
305  + cos2b * Mz2 * (- 1.0/2.0 + sW2) * Id3 );
306  gslpp::matrix<gslpp::complex> eLR( mySUSY.v1()/sqrt(2.0) * mySUSY.getTEhat().hconjugate()
307  - mySUSY.getMuH() * Me * mySUSY.getTanb() );
308  gslpp::matrix<gslpp::complex> eRR( msEhat2modified + Me * Me - cos2b * Mz2 * sW2 * Id3 );
309 
310  // ... until here
311  //else use the following...
312 // gslpp::matrix<gslpp::complex> eLL( mySUSY.msLhat2 + Me * Me
313 // + cos2b * Mz2 * (- 1.0/2.0 + sW2) * Id3 );
314 // gslpp::matrix<gslpp::complex> eLR( mySUSY.v1()/sqrt(2.0) * mySUSY.getTEhat().hconjugate()
315 // - mySUSY.getMuH() * Me * mySUSY.getTanb() );
316 // gslpp::matrix<gslpp::complex> eRR( mySUSY.msEhat2 + Me * Me - cos2b * Mz2 * sW2 * Id3 );
317  // ... until here
318 
319  for(int i = 0; i < 3; i++)
320  {
321  for(int j = 0; j < 3; j++)
322  {
323  Mselectron2.assign(i, j, eLL(i,j));
324  Mselectron2.assign(i, j+3, eLR(i,j));
325  Mselectron2.assign(i+3, j, eLR(j,i).conjugate());
326  Mselectron2.assign(i+3, j+3, eRR(i,j));
327  }
328  }
329 
330  /*
331  * M.singularvalue(&U, &V, &S) decompose M into
332  * M = U diag(S) V^+.
333  */
334  gslpp::matrix<gslpp::complex> RlTmp(6, 6, 0.);
335  Mselectron2.eigensystem(RlTmp, m_se2_i);
336 
337  Rl_i = RlTmp.hconjugate();
338 
339  return true;
340 
341 }
342 
344 {
345  int newIndex[6];
346  for (int i = 0; i < 6; i++)
347  newIndex[i] = i;
348 
349  /* sort sfermion masses in increasing order */
350  for (int i = 0; i < 5; i++)
351  for (int k = i + 1; k < 6; k++)
352  if (m_sf2(i) > m_sf2(k)) {
353  std::swap(m_sf2(i), m_sf2(k));
354  std::swap(newIndex[i], newIndex[k]);
355  }
356 
357  /* sort the corresponding rotation matrix, where the first(second) index
358  * denotes mass(gauge) eigenstates. */
359  gslpp::matrix<gslpp::complex> myRf(6, 6, 0.);
360  for (int i = 0; i < 6; i++)
361  for (int k = 0; k < 6; k++)
362  myRf.assign(k, i, Rf(newIndex[k], i));
363  Rf = myRf;
364 }
gslpp::matrix< gslpp::complex > msUhat2
Definition: SUSY.h:499
A class for the chargino masses.
Definition: Mchargino.h:22
gslpp::matrix< gslpp::complex > msDhat2
Definition: SUSY.h:499
double Ml_Q(const lepton l) const
Definition: SUSY.h:467
bool CalcChargino(gslpp::matrix< gslpp::complex > &U_i, gslpp::matrix< gslpp::complex > &V_i, gslpp::vector< double > &mch_i)
Computes the chargino spectrum at tree level.
gslpp::matrix< gslpp::complex > getTUhat() const
Gets the trilinear-coupling matrix for up-type squarks.
Definition: SUSY.h:305
double v1() const
Definition: SUSY.cpp:268
virtual double Mw_tree() const
The tree-level mass of the boson, .
bool CalcSelectron(gslpp::matrix< gslpp::complex > &Rl_i, gslpp::vector< double > &m_se2_i)
Computes the charged-slepton spectrum at tree level.
A class for constructing and defining operations on real matrices.
bool CalcSdown(gslpp::matrix< gslpp::complex > &Rd_i, gslpp::vector< double > &m_sd2_i)
Computes the down-type squark spectrum at tree level.
gslpp::matrix< gslpp::complex > msNhat2
Definition: SUSY.h:499
bool CalcHiggs(double mh[4], gslpp::complex &saeff_i)
Computes the Higgs spectrum at tree level.
Definition: QCD.h:731
const SUSY & mySUSY
Definition: SUSYSpectrum.h:300
gslpp::matrix< gslpp::complex > msLhat2
Definition: SUSY.h:499
gslpp::matrix< gslpp::complex > getVCKM() const
A get method to retrieve the CKM matrix.
Definition: QCD.h:735
gslpp::matrix< gslpp::complex > getTEhat() const
Gets the trilinear-coupling matrix for charged sleptons.
Definition: SUSY.h:406
SUSYSpectrum(const SUSY &SUSY_in)
A SUSYSpectrum constructor.
double getTanb() const
Gets .
Definition: SUSY.h:128
gslpp::matrix< gslpp::complex > getTDhat() const
Gets the trilinear-coupling matrix for down-type squarks.
Definition: SUSY.h:314
Definition: QCD.h:732
bool CalcSup(gslpp::matrix< gslpp::complex > &Ru_i, gslpp::vector< double > &m_su2_i)
Computes the up-type squark spectrum at tree level.
void SortSfermionMasses(gslpp::vector< double > &m_sf2, gslpp::matrix< gslpp::complex > &Rf) const
bool CalcNeutralino(gslpp::matrix< gslpp::complex > &N_i, gslpp::vector< double > &mneu_i)
Computes the neutralino spectrum at tree level.
double getSinb() const
Gets .
Definition: SUSY.h:160
gslpp::matrix< gslpp::complex > Msneutrino2
Definition: SUSYSpectrum.h:325
double Mq_Q(const quark q) const
Definition: SUSY.h:451
double getCosb() const
Gets .
Definition: SUSY.h:169
A class for the CKM matrix elements.
Definition: CKM.h:23
gslpp::matrix< gslpp::complex > getTNhat() const
Gets the trilinear-coupling matrix for sneutrinos.
Definition: SUSY.h:397
An observable class for the -boson mass.
Definition: Mw.h:22
double mHptree
Definition: SUSY.h:492
gslpp::complex getM2() const
Gets the wino mass.
Definition: SUSY.h:92
A class for constructing and defining operations on real vectors.
gslpp::matrix< gslpp::complex > Msdown2
Definition: SUSYSpectrum.h:325
gslpp::matrix< gslpp::complex > Msup2
Stores the tree-level Up-squark, Down-squark, Sneutrino, and Slepton mass matrix. ...
Definition: SUSYSpectrum.h:325
gslpp::matrix< gslpp::complex > Mselectron2
Definition: SUSYSpectrum.h:325
bool CalcSneutrino(gslpp::matrix< gslpp::complex > &Rn_i, gslpp::vector< double > &m_sn2_i)
Computes the sneutrino spectrum at tree level.
gslpp::complex getMuH() const
Gets the parameter in the superpotential.
Definition: SUSY.h:119
gslpp::complex getM1() const
Gets the bino mass.
Definition: SUSY.h:83
A class for defining operations on and functions of complex numbers.
Definition: gslpp_complex.h:35
A class for the neutralino masses.
Definition: Mneutralino.h:24
complex sin(const complex &z)
gslpp::matrix< gslpp::complex > msQhat2
Definition: SUSY.h:499
A base class for SUSY models.
Definition: SUSY.h:25
double getMz() const
A get method to access the mass of the boson .
gslpp::matrix< gslpp::complex > Mneutralino
Definition: SUSYSpectrum.h:315
double v2() const
Definition: SUSY.cpp:273
complex sqrt(const complex &z)