first attempt to generate fermion self-energies

which return L, R and S components
FlexibleSUSY · Mar 28, 2017 · 4a775b4 · 4a775b4
1 parent 90a6c5e
commit 4a775b4
Show file tree

Hide file tree

Showing 4 changed files with 31 additions and 9 deletions.
diff --git a/meta/CConversion.m b/meta/CConversion.m
@@ -6,6 +6,7 @@
 ArrayType::usage="";
 VectorType::usage="";
 ScalarType::usage="";
+ChiralitySum::usage="";
 
 integerScalarCType::usage="represents an integer C type";
 realScalarCType::usage="represents a real scalar C type";
@@ -40,6 +41,7 @@
 UVec::usage="unit vector";
 UMat::usage="matrix projector";
 Eval::usage="calls .eval()";
+PS::usage="scalar projector";
 
 HaveSameDimension::usage = "Checks if given types have same
 dimension";
@@ -151,6 +153,9 @@
 GetElementType[CConversion`MatrixType[type_, __]] := type;
 GetElementType[CConversion`TensorType[type_, __]] := type;
 
+GetScalarElementType[CConversion`ChiralitySum[type_]] :=
+    CConversion`ChiralitySum[GetScalarElementType[type]];
+
 GetScalarElementType[type_] :=
     CConversion`ScalarType[GetElementType[type]];
 
@@ -211,6 +216,9 @@
 CreateCType[CConversion`TensorType[t_, dims__]] :=
     EigenTensor[CreateCType[ScalarType[t]], Sequence @@ (ToString /@ {dims})];
 
+CreateCType[CConversion`ChiralitySum[type_]] :=
+    "LRS_tensor<" <> CreateCType[type] <> ">";
+
 CastTo[expr_String, toType_ /; toType === None] := expr;
 
 CastTo[expr_String, toType_] :=

diff --git a/meta/SelfEnergies.m b/meta/SelfEnergies.m
@@ -207,6 +207,7 @@ therefore not be accessed in the form Glu(gO2).
                AppendTo[result, SelfEnergies`FSSelfEnergy[field[1]       , expr[[1]]]];
                AppendTo[result, SelfEnergies`FSSelfEnergy[field[SARAH`PR], expr[[2]]]];
                AppendTo[result, SelfEnergies`FSSelfEnergy[field[SARAH`PL], expr[[3]]]];
+               AppendTo[result, SelfEnergies`FSSelfEnergy[field          , expr]];
               ];
            result
           ];
@@ -423,14 +424,20 @@ therefore not be accessed in the form Glu(gO2).
 ExpressionToStringSequentially[expr_, heads_, result_String] :=
     result <> " = " <> ExpressionToString[expr, heads] <> ";\n";
 
-PrepareExpr[expr_, vertexRules_] :=
+PrepareExpr[{S_,R_,L_}, CConversion`ChiralitySum[t_], vertexRules_] :=
+    PrepareExpr[S, t, vertexRules] PS +
+    PrepareExpr[R, t, vertexRules] PR +
+    PrepareExpr[L, t, vertexRules] PL;
+
+PrepareExpr[expr_, _, vertexRules_] :=
     expr /.
     vertexRules /.
     a_[List[i__]] :> a[i] /.
     ReplaceGhosts[FlexibleSUSY`FSEigenstates] /.
     C -> 1;
 
 GetNPointFunctionType[n_SelfEnergies`Tadpole] := TreeMasses`GetTadpoleType[GetField[n]];
+GetNPointFunctionType[n:(_[_, __List])]       := CConversion`ChiralitySum[TreeMasses`GetSelfEnergyType[GetField[n]]];
 GetNPointFunctionType[n_]                     := TreeMasses`GetSelfEnergyType[GetField[n]];
 
 DeclareMomentum[_SelfEnergies`Tadpole] := "";
@@ -452,7 +459,7 @@ therefore not be accessed in the form Glu(gO2).
 CreateFunctionPrototype[n_, type_, loops_] :=
     CreateFunctionName[n, loops] <> "(" <> DeclareMomentum[n] <> DeclareFieldIndices[GetField[n], type] <> ") const";
 
-HermitianizeLater[s_] :=
+HermitianizeLater[s_] := False &&
     GetDimension[GetField[s]] > 1 &&
     (IsScalar[GetField[s]] || IsVector[GetField[s]]) &&
     SelfEnergyIsSymmetric[GetField[s]];
@@ -465,17 +472,24 @@ therefore not be accessed in the form Glu(gO2).
            ]& /@ {ex})]
           ];
 
+SumOverTwoExternalFieldIndices[ex_List, self_] :=
+    SumOverTwoExternalFieldIndices[#,self]& /@ ex;
+
 (* sum over external gauge indices *)
-SumOverExternalIndices[self:(s_[f_, ex__])] :=
-    Module[{dim = GetDimension[f]},
-           s[f, Sequence @@ (CConversion`RefactorSums[
+SumOverTwoExternalFieldIndices[ex_, self_] :=
+    Module[{dim = GetDimension[GetField[self]]},
+           CConversion`RefactorSums[
                SARAH`sum[SARAH`gO1, 1, dim,
                          SARAH`sum[SARAH`gO2, If[HermitianizeLater[self], SARAH`gO1 + 1, 1],
-                                   dim, UMat[dim,dim,SARAH`gO1,SARAH`gO2] #]
+                                   dim, UMat[dim,dim,SARAH`gO1,SARAH`gO2] ex]
                         ]
-           ]& /@ {ex})]
+           ]
           ];
 
+(* sum over external gauge indices *)
+SumOverExternalIndices[self:(s_[f_, ex__])] :=
+    s[f, Sequence @@ (SumOverTwoExternalFieldIndices[#,self]& /@ {ex})];
+
 CreateNPointFunction[nPointFunction_, vertexRules_List] :=
     Module[{dim = GetDimension[GetField[nPointFunction]],
             type = GetNPointFunctionType[nPointFunction],
@@ -504,7 +518,7 @@ therefore not be accessed in the form Glu(gO2).
            decl = "\n" <> ctype <> " CLASSNAME::" <> functionName <> "\n{\n" <>
               IndentText[ctype <> " result;\n\n" <>
                          ExpressionToStringSequentially[
-                             PrepareExpr[expr, vertexRules],
+                             PrepareExpr[expr, type, vertexRules],
                              TreeMasses`GetParticles[], "result"]  <>
                          MakeHermitian[nPointFunction, type] <>
                          "\nreturn result * oneOver16PiSqr;"

diff --git a/templates/mass_eigenstates.cpp.in b/templates/mass_eigenstates.cpp.in
@@ -31,7 +31,6 @@
  */
 
 #include "@ModelName@_mass_eigenstates.hpp"
-#include "eigen_utils.hpp"
 #include "wrappers.hpp"
 #include "linalg2.hpp"
 #include "numerics2.hpp"

diff --git a/templates/mass_eigenstates.hpp.in b/templates/mass_eigenstates.hpp.in
@@ -36,6 +36,7 @@
 #include "@ModelName@_physical.hpp"
 #include "@ModelName@_info.hpp"
 #include "two_loop_corrections.hpp"
+#include "eigen_utils.hpp"
 #include "error.hpp"
 #include "problems.hpp"
 #include "config.h"