Merge branch 'master' into hotfix-stable

FeynCalc/Dirac/Anti5.m

@@ -54,10 +54,16 @@
 	FCVerbose		-> False
 };
 
-Anti5[expr_, opts:OptionsPattern[]] :=
+Anti5[a_ == b_, rest___] :=
+	Anti5[a,rest] == Anti5[b,rest];
+
+Anti5[expr_List, rest___]:=
+	Anti5[#, rest]&/@expr;
+
+Anti5[expr_/; !MemberQ[{List,Equal},expr], opts:OptionsPattern[]] :=
 	Anti5[expr, 1, opts];
 
-Anti5[expr_, n_/; !OptionQ[n] && MatchQ[n, Infinity | -Infinity | _Integer], OptionsPattern[]] :=
+Anti5[expr_/; !MemberQ[{List,Equal},expr], n_/; !OptionQ[n] && MatchQ[n, Infinity | -Infinity | _Integer], OptionsPattern[]] :=
 	Block[ {new = 0,ex,terms,rest,res, eps, freePart, dsPart, diracObjects,
 			diracObjectsEval, null1, null2, dsHead, time, repRule},
 
@@ -187,15 +193,15 @@
 anti5MoveRight[r_. ch:holdDOT[___, (DiracGamma[5|6|7] | (_. DiracGamma[6] + _. DiracGamma[7])), d__], Infinity] :=
 	r anti5MoveRight[FixedPoint[(# /. chiralTrickAnticommutingDDimRight -> commonGamma5Properties /. commonGamma5Properties -> chiralTrickAnticommutingDDimRight)&,
 	chiralTrickAnticommutingDDimRight@@ch]/. chiralTrickAnticommutingDDimRight->holdDOT, Infinity]/;
-	MemberQ[{"NDR","NDR-Drop"},FeynCalc`Package`DiracGammaScheme] && FreeQ[r,holdDOT] &&
+	MemberQ[{"NDR","NDR-Discard"},FeynCalc`Package`DiracGammaScheme] && FreeQ[r,holdDOT] &&
 	(MatchQ[FCGetDimensions[{d}, FreeQ->{DiracGamma[5],DiracGamma[6],DiracGamma[7]}, ChangeDimension->True],{_Symbol}]);
 
 (*			-> Limited number of rearrangements	*)
 anti5MoveRight[r_. holdDOT[b___, ch: (DiracGamma[5|6|7] | (_. DiracGamma[6] + _. DiracGamma[7])), c_, d___], n_Integer?NonNegative] :=
 	r anti5MoveRight[FixedPoint[(# /. chiralTrickAnticommutingDDimRight -> commonGamma5Properties /. commonGamma5Properties -> chiralTrickAnticommutingDDimRight)&,
 	chiralTrickAnticommutingDDimRight@@Join[{b}/. DiracGamma->diga,{ch,c}, {d} /. DiracGamma->diga]]/.
 	chiralTrickAnticommutingDDimRight->holdDOT /.diga->DiracGamma, n+1]/;
-	MemberQ[{"NDR","NDR-Drop"},FeynCalc`Package`DiracGammaScheme] && FreeQ[r,holdDOT] && (n+1<=max) &&
+	MemberQ[{"NDR","NDR-Discard"},FeynCalc`Package`DiracGammaScheme] && FreeQ[r,holdDOT] && (n+1<=max) &&
 	MatchQ[FCGetDimensions[{c}, FreeQ->{DiracGamma[5],DiracGamma[6],DiracGamma[7]}, ChangeDimension->True],{_Symbol}];
 
 (*	------------------------------------------------------------------------------	*)
@@ -248,15 +254,15 @@
 anti5MoveLeft[r_. ch:holdDOT[d__, (DiracGamma[5|6|7] | (_. DiracGamma[6] + _. DiracGamma[7])), ___], -Infinity] :=
 	r anti5MoveLeft[FixedPoint[(# /. chiralTrickAnticommutingDDimLeft -> commonGamma5Properties /. commonGamma5Properties -> chiralTrickAnticommutingDDimLeft)&,
 							chiralTrickAnticommutingDDimLeft@@ch]/. chiralTrickAnticommutingDDimLeft->holdDOT, -Infinity]/;
-	MemberQ[{"NDR","NDR-Drop"},FeynCalc`Package`DiracGammaScheme] && FreeQ[r,holdDOT] &&
+	MemberQ[{"NDR","NDR-Discard"},FeynCalc`Package`DiracGammaScheme] && FreeQ[r,holdDOT] &&
 	(MatchQ[FCGetDimensions[{d}, FreeQ->{DiracGamma[5],DiracGamma[6],DiracGamma[7]}, ChangeDimension->True],{_Symbol}]);
 
 (*			-> Limited number of rearrangements	*)
 anti5MoveLeft[r_. holdDOT[b___, c_, ch:(DiracGamma[5|6|7] | (_. DiracGamma[6] + _. DiracGamma[7])), d___], n_Integer?NonPositive] :=
 	r anti5MoveLeft[FixedPoint[(# /. chiralTrickAnticommutingDDimLeft -> commonGamma5Properties /. commonGamma5Properties -> chiralTrickAnticommutingDDimLeft)&,
 							chiralTrickAnticommutingDDimLeft@@Join[{b}/. DiracGamma->diga,{c,ch}, {d} /. DiracGamma->diga]]/.
 							chiralTrickAnticommutingDDimLeft->holdDOT /.diga->DiracGamma, n-1]/;
-		MemberQ[{"NDR","NDR-Drop"},FeynCalc`Package`DiracGammaScheme] && FreeQ[r,holdDOT] && (max<=n-1) &&
+		MemberQ[{"NDR","NDR-Discard"},FeynCalc`Package`DiracGammaScheme] && FreeQ[r,holdDOT] && (max<=n-1) &&
 		MatchQ[FCGetDimensions[{c}, FreeQ->{DiracGamma[5],DiracGamma[6],DiracGamma[7]}, ChangeDimension->True],{_Symbol}];
 
 (*	------------------------------------------------------------------------------	*)

FeynCalc/Dirac/Chisholm.m

@@ -49,7 +49,13 @@
 	NonCommutative 		-> True
 };
 
-Chisholm[expr_, OptionsPattern[]] :=
+Chisholm[a_ == b_, opts:OptionsPattern[]] :=
+	Chisholm[a,opts] == Chisholm[b,opts];
+
+Chisholm[expr_List, opts:OptionsPattern[]]:=
+	Chisholm[#, opts]&/@expr;
+
+Chisholm[expr_/; !MemberQ[{List,Equal},expr], OptionsPattern[]] :=
 	Block[{ex, tmp, tmpli, tmpli1, tmpli2, terms,rest,res, holdDOT, eps, freePart, dsPart, diracObjects,
 			diracObjectsEval, null1, null2, dsHead, time, repRule, mode, chisholmRule1, chisholmRule2, maxIterations,
 			nonComm, chisholmRuleInsideTrace},

FeynCalc/Dirac/DiracChainCombine.m

@@ -54,7 +54,13 @@
 	FCI 				-> False
 };
 
-DiracChainCombine[expr_, OptionsPattern[]] :=
+DiracChainCombine[a_ == b_, opts:OptionsPattern[]] :=
+	DiracChainCombine[a,opts] == DiracChainCombine[b,opts];
+
+DiracChainCombine[expr_List, opts:OptionsPattern[]]:=
+	DiracChainCombine[#, opts]&/@expr;
+
+DiracChainCombine[expr_/; !MemberQ[{List,Equal},expr], OptionsPattern[]] :=
 	Block[{ex, res, rules=rulesMain},
 
 		optDiracGammaCombine = OptionValue[DiracGammaCombine];

FeynCalc/Dirac/DiracChainExpand.m

@@ -50,7 +50,13 @@
 	Momentum 			-> All
 };
 
-DiracChainExpand[expr_, OptionsPattern[]] :=
+DiracChainExpand[a_ == b_, opts:OptionsPattern[]] :=
+	DiracChainExpand[a,opts] == DiracChainExpand[b,opts];
+
+DiracChainExpand[expr_List, opts:OptionsPattern[]]:=
+	DiracChainExpand[#, opts]&/@expr;
+
+DiracChainExpand[expr_/; !MemberQ[{List,Equal},expr], OptionsPattern[]] :=
 	Block[{ex, null1, null2, uniqList, solsList, repRule, ruleDelayed, res, rhs},
 
 	momList 			= OptionValue[Momentum];

FeynCalc/Dirac/DiracChainFactor.m

@@ -32,7 +32,13 @@
 	FCE -> False
 };
 
-DiracChainFactor[expr_, OptionsPattern[]] :=
+DiracChainFactor[a_ == b_, opts:OptionsPattern[]] :=
+	DiracChainFactor[a,opts] == DiracChainFactor[b,opts];
+
+DiracChainFactor[expr_List, opts:OptionsPattern[]]:=
+	DiracChainFactor[#, opts]&/@expr;
+
+DiracChainFactor[expr_/; !MemberQ[{List,Equal},expr], OptionsPattern[]] :=
 	Block[ {ex, moms,res, diracChains},
 
 		If[ OptionValue[FCI],

FeynCalc/Dirac/DiracChainJoin.m

@@ -50,11 +50,17 @@
 	FCE 			-> False,
 	FCI 			-> False,
 	FCVerbose		-> False,
-	Factoring		-> True,
+	Factoring		-> {Factor2, 5000},
 	TraceOfOne		-> 4
 };
 
-DiracChainJoin[expr_, OptionsPattern[]] :=
+DiracChainJoin[a_ == b_, opts:OptionsPattern[]] :=
+	DiracChainJoin[a,opts] == DiracChainJoin[b,opts];
+
+DiracChainJoin[expr_List, opts:OptionsPattern[]]:=
+	DiracChainJoin[#, opts]&/@expr;
+
+DiracChainJoin[expr_/; !MemberQ[{List,Equal},expr], OptionsPattern[]] :=
 	Block[{	ex, tmp,  res, diracObjects, diracObjectsEval, null1, null2,
 			dsHead, time, repRule},
 
@@ -120,9 +126,9 @@
 
 			FCPrint[1, "DiracChainJoin: Inserting Dirac objects back.", FCDoControl->dchjVerbose];
 			time=AbsoluteTime[];
-			repRule = MapThread[Rule[#1,#2]&,{diracObjects,diracObjectsEval}];
+			repRule = Thread[Rule[diracObjects,diracObjectsEval]];
 			FCPrint[3,"DiracChainJoin: repRule: ",repRule , FCDoControl->dchjVerbose];
-			res =  ( tmp/.repRule);
+			res =  (tmp/. Dispatch[repRule]);
 			FCPrint[1, "DiracChainJoin: Done inserting Dirac objects back, timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->dchjVerbose],
 
 			(*Fast mode*)

FeynCalc/Dirac/DiracEquation.m

@@ -35,20 +35,26 @@
 	FCE -> False
 };
 
-DiracEquation[ex_, OptionsPattern[]] :=
-	Block[{expr,res},
+DiracEquation[a_ == b_, opts:OptionsPattern[]] :=
+	DiracEquation[a,opts] == DiracEquation[b,opts];
+
+DiracEquation[expr_List, opts:OptionsPattern[]]:=
+	DiracEquation[#, opts]&/@expr;
+
+DiracEquation[expr_/; !MemberQ[{List,Equal},expr], OptionsPattern[]] :=
+	Block[{ex,res},
 
 		If[	!OptionValue[FCI],
-			expr  = FCI[ex],
-			expr  = ex
+			ex  = FCI[expr],
+			ex  = expr
 		];
 
 		(* If there are no spinors or no Dirac matrices, then we can't apply the Dirac equation	*)
-		If[	FreeQ2[expr, {Spinor,DiracGamma}],
-			Return[expr];
+		If[	FreeQ2[ex, {Spinor,DiracGamma}],
+			Return[ex];
 		];
 
-		res  = DotSimplify[FixedPoint[diraceq,expr,5]/.PairContract->Pair, Expanding->False];
+		res  = DotSimplify[FixedPoint[diraceq,ex,5]/.PairContract->Pair, Expanding->False];
 
 		If[ OptionValue[FCE],
 			res = FCE[res]

FeynCalc/Dirac/DiracGammaCombine.m

@@ -32,7 +32,13 @@
 	FCI -> False
 };
 
-DiracGammaCombine[expr_, OptionsPattern[]] :=
+DiracGammaCombine[a_ == b_, opts:OptionsPattern[]] :=
+	DiracGammaCombine[a,opts] == DiracGammaCombine[b,opts];
+
+DiracGammaCombine[expr_List, opts:OptionsPattern[]]:=
+	DiracGammaCombine[#, opts]&/@expr;
+
+DiracGammaCombine[expr_/; !MemberQ[{List,Equal},expr], OptionsPattern[]] :=
 	Block[{	ex, tmp, res, holdDOT, freePart, diracPart, holdPlus, diracList, diracListEval, repRule,
 			null1, null2},
 

FeynCalc/Dirac/DiracGammaExpand.m

@@ -37,7 +37,13 @@
 	Momentum -> All
 };
 
-DiracGammaExpand[expr_, OptionsPattern[]] :=
+DiracGammaExpand[a_ == b_, opts:OptionsPattern[]] :=
+	DiracGammaExpand[a,opts] == DiracGammaExpand[b,opts];
+
+DiracGammaExpand[expr_List, opts:OptionsPattern[]]:=
+	DiracGammaExpand[#, opts]&/@expr;
+
+DiracGammaExpand[expr_/; !MemberQ[{List,Equal},expr], OptionsPattern[]] :=
 	Block[{ex, null1, null2, uniqList, solsList, repRule, momList, ruleDelayed, res},
 
 	momList = OptionValue[Momentum];

FeynCalc/Dirac/DiracOrder.m

@@ -11,17 +11,17 @@
 	Copyright (C) 2014-2020 Vladyslav Shtabovenko
 *)
 
-(* :Summary:  Lexicographic ordering of Dirac matrices						*)
+(* :Summary:  Canonical ordering of Dirac matrices							*)
 
 (* ------------------------------------------------------------------------ *)
 
 DiracOrder::usage =
-"DiracOrder[exp] orders the Dirac matrices in exp lexicographically. \
+"DiracOrder[exp] orders the Dirac matrices in exp canonically. \
 DiracOrder[exp, orderlist] orders the Dirac matrices in exp according \
 to orderlist.\n
 DiracOrder is also an option of DiracSimplify and some other functions dealing \
 with Dirac algebra. If set to True, the function DiracOrder will be applied to \
-the intermediate result to reorder the Dirac matrices lexicographically.
+the intermediate result to reorder the Dirac matrices canonically.
 ";
 
 DiracOrder::failmsg =
@@ -42,6 +42,7 @@
 Options[DiracOrder] = {
 	DiracGammaCombine	-> False,
 	DiracTrick 			-> True,
+	DotSimplify			-> True,
 	FCDiracIsolate		-> True,
 	FCE					-> False,
 	FCI					-> False,
@@ -50,11 +51,18 @@
 	MaxIterations		-> Infinity
 };
 
+DiracOrder[a_ == b_, rest___] :=
+	DiracOrder[a,rest] == DiracOrder[b,rest];
+
+DiracOrder[expr_List, rest___]:=
+	DiracOrder[#, rest]&/@expr;
+
 DiracOrder[expr_, (opts:OptionsPattern[])/;opts=!={}] :=
 	DiracOrder[expr, {}, opts];
 
-DiracOrder[expr_, orderList_List/; (!OptionQ[orderList] || orderList==={}), OptionsPattern[]]:=
-	Block[{ex,res,dsHead,dsPart,freePart,null1,null2,diracObjects,tmp, maxIterations, diracObjectsEval, repRule,time},
+DiracOrder[expr_/; !MemberQ[{List,Equal},expr], orderList_List/; (!OptionQ[orderList] || orderList==={}), OptionsPattern[]]:=
+	Block[{	ex,res,dsHead,dsPart,freePart,null1,null2,diracObjects,
+			tmp, maxIterations, diracObjectsEval, repRule,time},
 
 		maxIterations = OptionValue[MaxIterations];
 
@@ -78,12 +86,14 @@
 		FCPrint[3, "DiracOrder: Entering with ", ex, FCDoControl->doVerbose];
 
 		If[	OptionValue[FCDiracIsolate],
+
 			(* This is the normal mode which works well both for large and small expressions *)
 			FCPrint[1, "DiracOrder: Normal mode.", FCDoControl->doVerbose];
 			time=AbsoluteTime[];
 			FCPrint[1, "DiracOrder: Extracting Dirac objects.", FCDoControl->doVerbose];
 			ex = FCDiracIsolate[ex,FCI->True,Head->dsHead, DotSimplify->True, DiracChain->True,
-				DiracGammaCombine->OptionValue[DiracGammaCombine], FCJoinDOTs->OptionValue[FCJoinDOTs]];
+				DiracGammaCombine->OptionValue[DiracGammaCombine], FCJoinDOTs->OptionValue[FCJoinDOTs], Split->False,
+				LorentzIndex->True, CartesianIndex->True];
 
 
 			{freePart,dsPart} = FCSplit[ex,{dsHead}];
@@ -93,26 +103,56 @@
 			diracObjects = Cases[dsPart+null1+null2, dsHead[_], Infinity]//DeleteDuplicates//Sort;
 			FCPrint[1, "DiracOrder: Done extracting Dirac objects, timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose];
 
+			diracObjectsEval = diracObjects;
+
+			If[	OptionValue[DiracTrick]===True,
+				time=AbsoluteTime[];
+				FCPrint[1, "DiracOrder: Applying DiracTrick.", FCDoControl->doVerbose];
+				diracObjectsEval = DiracTrick[#,FCI->True]&/@diracObjectsEval;
+				diracObjectsEval =FeynCalc`Package`diracChainContract /@ diracObjectsEval;
+				FCPrint[1, "DiracOrder: Done applying DiracTrick,timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose];
+				FCPrint[3, "DiracOrder: After DiracTrick: ", diracObjectsEval, FCDoControl->doVerbose];
+			];
+
 			time=AbsoluteTime[];
 			If[orderList=!={},
 				FCPrint[1, "DiracOrder: Ordering according to: ", orderList, FCDoControl->doVerbose];
-				diracObjectsEval = Map[diracOrderCustom[#,orderList]&, (diracObjects/. DOT -> holdDOT/.dsHead->Identity)];
+				diracObjectsEval = Map[diracOrderCustom[#,orderList]&, (diracObjectsEval/. DOT -> holdDOT/.dsHead->Identity)];
 				FCPrint[1, "DiracOrder: Ordering done, timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose],
 
-				FCPrint[1, "DiracOrder. Using lexicographic ordering.", FCDoControl->doVerbose];
-				diracObjectsEval = Map[diracOrderLex[#,maxIterations]&, (diracObjects/. DOT -> holdDOT/.dsHead->Identity)];
+				FCPrint[1, "DiracOrder. Using canonical ordering.", FCDoControl->doVerbose];
+				diracObjectsEval = Map[diracOrderLex[#,maxIterations]&, (diracObjectsEval/. DOT -> holdDOT/.dsHead->Identity)];
 				FCPrint[1, "DiracOrder: Ordering done, timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose]
 			];
 
+			FCPrint[3, "DiracOrder: Intermediate result: ", diracObjectsEval, FCDoControl->doVerbose];
+
+
 			time=AbsoluteTime[];
 			FCPrint[1, "DiracOrder: Inserting Dirac objects back.", FCDoControl->doVerbose];
+			diracObjectsEval = diracObjectsEval /. holdDOT[]->1 /.holdDOT->DOT /. PairContract->Pair;
+
+			If[	OptionValue[DotSimplify],
+				time=AbsoluteTime[];
+				FCPrint[1, "DiracOrder: Applying DotSimplify.", FCDoControl->doVerbose];
+				diracObjectsEval = DotSimplify[#,FCI->True, Expanding->True]&/@diracObjectsEval;
+				FCPrint[1, "DiracOrder: Done applying DotSimplify,timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose];
+				FCPrint[3, "DiracOrder: After DotSimplify: ", diracObjectsEval, FCDoControl->doVerbose];
+			];
 
-			diracObjectsEval = diracObjectsEval /. holdDOT[]->1 /.holdDOT->DOT /. PairContract -> Pair;
-			If[	!FreeQ[diracObjectsEval,CartesianPair],
-				diracObjectsEval = diracObjectsEval/. CartesianPair->CartesianPairContract /. CartesianPairContract->CartesianPair
+
+			If[	(OptionValue[DiracTrick]===True || OptionValue[DiracTrick]===Last)  && FreeQ[orderList,ExplicitLorentzIndex[0]],
+				time=AbsoluteTime[];
+				FCPrint[1, "DiracOrder: Applying DiracTrick.", FCDoControl->doVerbose];
+				diracObjectsEval = DiracTrick[#,FCI->True, FCDiracIsolate -> {DotSimplify -> False, Expanding -> False,
+					FCJoinDOTs -> False, DiracGammaCombine -> False, Factoring -> False}]&/@diracObjectsEval;
+				diracObjectsEval =FeynCalc`Package`diracChainContract /@ diracObjectsEval;
+				FCPrint[1, "DiracOrder: Done applying DiracTrick,timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose];
+				FCPrint[3, "DiracOrder: After DiracTrick: ", diracObjectsEval, FCDoControl->doVerbose];
 			];
 
 			repRule = Thread[Rule[diracObjects, diracObjectsEval]];
+
 			FCPrint[3,"DiracOrder: repRule: ",repRule , FCDoControl->doVerbose];
 			tmp = freePart + (dsPart/. Dispatch[repRule]);
 			FCPrint[1, "DiracOrder: Done inserting Dirac objects back, timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose];
@@ -133,24 +173,22 @@
 				FCPrint[1, "DiracOrder: Ordering done, timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose]
 
 			];
-			tmp = tmp/. holdDOT[]->1 /.holdDOT->DOT /. PairContract -> Pair;
-			If[	!FreeQ[diracObjectsEval,CartesianPair],
-				diracObjectsEval = diracObjectsEval/. CartesianPair->CartesianPairContract /. CartesianPairContract->CartesianPair
-			];
-
-		];
+			tmp = tmp/. holdDOT[]->1 /.holdDOT->DOT /. PairContract->Pair;
+			tmp = FeynCalc`Package`diracChainContract[tmp];
 
-
-		res = tmp;
-
-		(* If orderingList contains an ExplicitLorentzIndex[0], do not apply DiracTrick, since it would mess up the ordering!	*)
-		If[	OptionValue[DiracTrick] && FreeQ[orderList,ExplicitLorentzIndex[0]],
+			If[	OptionValue[DiracTrick] && FreeQ[orderList,ExplicitLorentzIndex[0]],
 				time=AbsoluteTime[];
 				FCPrint[1, "DiracOrder: Applying DiracTrick.", FCDoControl->doVerbose];
-				res = DiracTrick[res, FCI->True];
+				tmp = DiracTrick[tmp,FCI->True];
+				tmp = FeynCalc`Package`diracChainContract[tmp];
 				FCPrint[1, "DiracOrder: Done applying DiracTrick,timing: ", N[AbsoluteTime[] - time, 4], FCDoControl->doVerbose]
+			]
+
 		];
 
+
+		res = tmp;
+
 		FCPrint[1, "DiracOrder: Leaving.", FCDoControl->doVerbose];
 		FCPrint[3, "DiracOrder: Leaving with ", res, FCDoControl->doVerbose];
 
@@ -163,7 +201,7 @@
 	];
 
 diracOrderLex[x_, maxIterations_]:=
-	FixedPoint[(# /. {
+	FixedPoint[(FeynCalc`Package`diracChainContract[#] /. {
 		holdDOT[a___,DiracGamma[(h1:LorentzIndex|Momentum|CartesianIndex|CartesianMomentum)[ar1__], dim1_:4],DiracGamma[(h2:LorentzIndex|Momentum|CartesianIndex|CartesianMomentum)[ar2__], dim2_:4],b___]/;
 			!OrderedQ[{h1[First[{ar1}],dim1],h2[First[{ar2}],dim2]}] && h1[First[{ar1}],dim1]=!=h2[First[{ar2}],dim2] :>
 			-holdDOT[a, DiracGamma[h2[ar2],dim2], DiracGamma[h1[ar1],dim1] ,b] +
@@ -182,7 +220,7 @@
 	})&, x, maxIterations]
 
 customOrdering[x_, currentElement_]:=
-	x //. {
+	FeynCalc`Package`diracChainContract[x] //. {
 		holdDOT[a___,DiracGamma[(h1:LorentzIndex|Momentum|CartesianIndex|CartesianMomentum)[ar1__], dim1_:4],
 			DiracGamma[(h2:LorentzIndex|Momentum|CartesianIndex|CartesianMomentum)[ar2__], dim2_:4],b___]/;
 			!FreeQ[h2[First[{ar2}],dim2],currentElement] && h1[First[{ar1}],dim1]=!=h2[First[{ar2}],dim2] :>