flexiblesusy is hosted by Hepforge, IPPP Durham
FlexibleSUSY
sfermions.cpp
Go to the documentation of this file.
1// ====================================================================
2// This file is part of FlexibleSUSY.
3//
4// FlexibleSUSY is free software: you can redistribute it and/or modify
5// it under the terms of the GNU General Public License as published
6// by the Free Software Foundation, either version 3 of the License,
7// or (at your option) any later version.
8//
9// FlexibleSUSY is distributed in the hope that it will be useful, but
10// WITHOUT ANY WARRANTY; without even the implied warranty of
11// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12// General Public License for more details.
13//
14// You should have received a copy of the GNU General Public License
15// along with FlexibleSUSY. If not, see
16// <http://www.gnu.org/licenses/>.
17// ====================================================================
18
25#include "sfermions.hpp"
26#include "linalg2.hpp"
27
28#include <cmath>
29
30#include <Eigen/Core>
31
32namespace flexiblesusy {
33
34namespace {
35
36double conj(double x) noexcept { return x; }
37
38double sqr(double x) noexcept { return x*x; }
39
40int sign(double x) noexcept { return x >= 0.0 ? 1 : -1; }
41
42} // anonymous namespace
43
44namespace sfermions {
45
46static constexpr double oneOverRoot2 = 0.7071067811865475; // 1/sqrt(2.)
47
49 0.5, -0.5, 0.5, -0.5
50};
51
53 1./3., 1./3., -1., -1.
54};
55
57 -4./3., 2./3., 0., 2.
58};
59
65double diagonalize_sfermions_2x2(const Mass_data& pars, double& msf1, double& msf2)
66{
67 const double ml2 = pars.ml2;
68 const double mr2 = pars.mr2;
69 const double yf = pars.yf;
70 const double vd = pars.vd;
71 const double vu = pars.vu;
72 const double gY = pars.gY;
73 const double g2 = pars.g2;
74 const double Tyf = pars.Tyf;
75 const double mu = pars.mu;
76 const double T3 = pars.T3;
77 const double Yl = pars.Yl;
78 const double Yr = pars.Yr;
79 const double vev2 = 0.25 * (sqr(vd) - sqr(vu));
80 Eigen::Matrix<double,2,2> mass_matrix;
81
82 // fill sfermion phi in mass matix in basis (phi_L phi_R)
83 if (sign(T3) > 0) {
84 mass_matrix(0,0) = ml2 + 0.5 * sqr(yf) * sqr(vu)
85 + (T3 * sqr(g2) - 0.5 * Yl * sqr(gY)) * vev2;
86 mass_matrix(0,1) = oneOverRoot2 * (vu*conj(Tyf) - vd*conj(yf)*mu);
87 mass_matrix(1,0) = conj(mass_matrix(0,1));
88 mass_matrix(1,1) = mr2 + 0.5 * sqr(yf) * sqr(vu)
89 - 0.5 * Yr * sqr(gY) * vev2;
90 } else {
91 mass_matrix(0,0) = ml2 + 0.5 * sqr(yf) * sqr(vd)
92 + (T3 * sqr(g2) - 0.5 * Yl * sqr(gY)) * vev2;
93 mass_matrix(0,1) = oneOverRoot2 * (vd*conj(Tyf) - vu*conj(yf)*mu);
94 mass_matrix(1,0) = conj(mass_matrix(0,1));
95 mass_matrix(1,1) = mr2 + 0.5 * sqr(yf) * sqr(vd)
96 - 0.5 * Yr * sqr(gY) * vev2;
97 }
98
99 Eigen::Array<double,2,1> msf;
100 Eigen::Matrix<double, 2, 2> Zf;
101 diagonalize_hermitian(mass_matrix, msf, Zf);
102
103 double theta;
104
105 if (sign(Zf(0,0)) == sign(Zf(1,1))) {
106 theta = std::acos(std::abs(Zf(0,0)));
107 } else {
108 theta = std::acos(std::abs(Zf(0,1)));
109 Zf.col(0).swap(Zf.col(1));
110 std::swap(msf(0), msf(1));
111 }
112
113 msf1 = msf(0);
114 msf2 = msf(1);
115 theta = sign(mass_matrix(0,1) / (mass_matrix(0,0) - mass_matrix(1,1)))
116 * std::abs(theta);
117
118 return theta;
119}
120
121} // namespace sfermions
122} // namespace flexiblesusy
double diagonalize_sfermions_2x2(const Mass_data &pars, double &msf1, double &msf2)
returns mixing angle and sets (squared) mass eigenvalues
Definition: sfermions.cpp:65
const double Hypercharge_left[NUMBER_OF_MSSM_SPARTICLES]
Definition: sfermions.cpp:52
static constexpr double oneOverRoot2
Definition: sfermions.cpp:46
const double Isospin[NUMBER_OF_MSSM_SPARTICLES]
Definition: sfermions.cpp:48
const double Hypercharge_right[NUMBER_OF_MSSM_SPARTICLES]
Definition: sfermions.cpp:56
T sqr(T x)
Definition: mathdefs.hpp:62
void diagonalize_hermitian(const Eigen::Matrix< Scalar, N, N > &m, Eigen::Array< Real, N, 1 > &w, Eigen::Matrix< Scalar, N, N > &z)
Definition: linalg2.hpp:390
T sign(T x)
Definition: mathdefs.hpp:65
constexpr Complex< T > conj(const Complex< T > &z) noexcept
Definition: complex.hpp:48
void swap(nlohmann::basic_json< ObjectType, ArrayType, StringType, BooleanType, NumberIntegerType, NumberUnsignedType, NumberFloatType, AllocatorType, JSONSerializer, BinaryType, CustomBaseClass > &j1, nlohmann::basic_json< ObjectType, ArrayType, StringType, BooleanType, NumberIntegerType, NumberUnsignedType, NumberFloatType, AllocatorType, JSONSerializer, BinaryType, CustomBaseClass > &j2) noexcept(//NOLINT(readability-inconsistent-declaration-parameter-name, cert-dcl58-cpp) is_nothrow_move_constructible< nlohmann::basic_json< ObjectType, ArrayType, StringType, BooleanType, NumberIntegerType, NumberUnsignedType, NumberFloatType, AllocatorType, JSONSerializer, BinaryType, CustomBaseClass > >::value &&//NOLINT(misc-redundant-expression, cppcoreguidelines-noexcept-swap, performance-noexcept-swap) is_nothrow_move_assignable< nlohmann::basic_json< ObjectType, ArrayType, StringType, BooleanType, NumberIntegerType, NumberUnsignedType, NumberFloatType, AllocatorType, JSONSerializer, BinaryType, CustomBaseClass > >::value)
exchanges the values of two JSON objects
Definition: json.hpp:24537
double Yr
Hypercharge of right-handed sfermion.
Definition: sfermions.hpp:50
double mu
Superpotential parameter.
Definition: sfermions.hpp:47
double yf
Yukawa coupling.
Definition: sfermions.hpp:43
double Tyf
trilinear coupling
Definition: sfermions.hpp:46
double Yl
Hypercharge of left-handed sfermion.
Definition: sfermions.hpp:49
double mr2
soft mass of right-handed sfermion
Definition: sfermions.hpp:42
double g2
gauge couplings (not GUT normalized)
Definition: sfermions.hpp:45
double ml2
soft mass of left-handed sfermion
Definition: sfermions.hpp:41