Invar package

Manipulation of the invariants of Riemann tensors

José M. Martín-García     jose@xact.es

David Yllanes     yllanes@xact.es
Universidad Complutense de Madrid

Renato Portugal     
Laboratório Nacional de Computação Científica, Petrópolis, Brazil

(c) 2006-2020, under GPL

http://www.xact.es/

http://groups.google.com/group/xAct

Set up of the system

Date[]

xAct`Invar`$Version={"2.0.5",{2013,7,1}}

xAct`Invar`$xTensorVersionExpected={"1.1.0",{2013,9,1}}

Off[General::spell]
Off[General::spell1]

TODO

- Nondual invariants can be handled in any dimension, but so far dual invariants only work in dimension 4. Generalize that. Answer to Teake in the group (16 June 2011): “And I agree that the codes of Invar2 would require slight generalizations to be able to generate and deal with the databases for dual invariants in dimensions other than 4, for those cases in which it makes sense. I just had a look at Invar.m and I see that the
private function riemannof is the only place where dimension 4 is hardwired. The function InvSimplify has also the error message “Duals
are only handled in dimension 4.”, which would have to be removed. Finally, some bookkeeping functions like MaxDualIndex would have to
take an additional dimension argument. I think that’s all.”

GPL

(* Invar, a free package for simplification of Riemann invariants *)

(* Copyright (C) 2006-2020 Jose M. Martin-Garcia, D. Yllanes and Renato Portugal *)

(* This program is free software; you can redistribute it and/or modify
 it under the terms of the GNU General Public License as published
 by the Free Software Foundation; either version 2 of the License,
  or (at your option) any later version.

This program is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
 General Public License for more details.

You should have received a copy of the GNU General Public License
 along with this program; if not, write to the Free Software
 Foundation, Inc., 59 Temple Place-Suite 330, Boston, MA 02111-1307,
  USA. 
*)

Info package

(* :Title: Invar *)

(* :Authors: Jose M. Martin-Garcia, David Yllanes and Renato Portugal *)

(* :Summary: Free package for simplification of Riemann invariants *)

(* :Brief Discussion:
   - Based on xTensor.
   - There is a companion version for Maple (using Canon) by Renato
     Portugal.
   - Simplifies Riemann invariants (scalars formed from the Riemann
     tensor of a metric-compatible connection and its Levi-Civita
     connection) by using permutation group theory.
   - Multiterm symmetries (cyclic, Bianchi, commutation and
     dimension-dependent identities) have been solved in advance
     and stored in solution-files.
   - Works with Riemann, Ricci, RicciScalar, Weyl and TFRicci.
*)
  
(* :Context: xAct`Invar` *)

(* :Package Version: 2.0.5 *)

(* :Copyright: Jose M. Martin-Garcia and David Yllanes (2006-2020) *)

(* :History: 
	- 22 Nov 2006 - Invar1
	-  8 Feb 2008 - Invar2
*)

(* :Keywords: *)

(* :Source: Invar.nb *)

(* :Warning: *)

(* :Mathematica Version: 6.0 and later *)

(* :Limitations:
   - Handles algebraic invariants with up to seven Riemann tensors
     or dual invariants with up to five Riemann tensors
   - Handles differential invariants with up to twelve metric derivatives
     or dual invariants with up to eight metric derivatives *)

Begin package

Protect against multiple loading of the package:

With[{xAct`Invar`Private`InvarSymbols=DeleteCases[Join[Names["xAct`Invar`*"],Names["xAct`Invar`Private`*"]],"$Version"|"xAct`Invar`$Version"|"$xTensorVersionExpected"|"xAct`Invar`$xTensorVersionExpected"]},
Unprotect/@xAct`Invar`Private`InvarSymbols;
Clear/@xAct`Invar`Private`InvarSymbols;
]

Explicit (not hidden) import of xTensor. QUESTION: Do we really need ExpressionManipulation?

If[Unevaluated[xAct`xCore`Private`$LastPackage]===xAct`xCore`Private`$LastPackage,xAct`xCore`Private`$LastPackage="xAct`Invar`"];

xAct`xCore`Private`$LastPackage

BeginPackage["xAct`Invar`",{"xAct`xTensor`","xAct`xPerm`","xAct`xCore`","xAct`ExpressionManipulation`"}]

$ContextPath

Check version of xTensor. We simply compare dates:

General::xTensorVersion="Loaded xTensor version `1` but expected version `2` at least.";
If[Not@OrderedQ@Map[Last,{$xTensorVersionExpected,xAct`xTensor`$Version}],Message[General::xTensorVersion,xAct`xTensor`$Version,$xTensorVersionExpected];
Abort[]]

Message:

Print[xAct`xCore`Private`bars];
Print["Package xAct`Invar`  version ",$Version[[1]],", ",$Version[[2]]];
Print["CopyRight (C) 2006-2020, J. M. Martin-Garcia, D. Yllanes and R. Portugal, under the General Public License."]

We specify the context xAct`Invar` to avoid overriding the Disclaimer's of xTensor and xPerm. However we need to turn off the message General:shdw temporarily:

Off[General::shdw]
xAct`Invar`Disclaimer[]:=Print["These are points 11 and 12 of the General Public License:\n\nBECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE, THERE IS NO WARRANTY FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES PROVIDE THE PROGRAM `AS IS´ WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY SERVICING, REPAIR OR CORRECTION.\n\nIN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY AND/OR REDISTRIBUTE THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED TO LOSS OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY YOU OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER PROGRAMS), EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGES."]
On[General::shdw]

If Invar is the last package begin read, print the GPL message:

If[xAct`xCore`Private`$LastPackage==="xAct`Invar`",
Unset[xAct`xCore`Private`$LastPackage];
Print[xAct`xCore`Private`bars];
Print["These packages come with ABSOLUTELY NO WARRANTY; for details type Disclaimer[]. This is free software, and you are welcome to redistribute it under certain conditions. See the General Public License for details."];
Print[xAct`xCore`Private`bars]]

We already have these symbols:

?xAct`Invar`*

Public functions

$InvarDirectory::usage="$InvarDirectory is a global variable giving the directory in which the package file and the database can be found.";

InvarCases::usage="InvarCases[s, d] returns the list of cases of order (number of metric derivatives) s and degree (number of Riemanns) d. InvarCases[s] returns the list of cases of order s. InvarCases[] returns the list of all cases handled by Invar (currently the 47 cases up to order 12 plus the algebraic case of degree 7).";
InvarDualCases::usage="InvarDualCases[s, d] returns the list of dual cases of order (number of metric derivatives) s and degree (number of Riemanns) d in dimension 4. InvarDualCases[s] returns the list of dual cases of order s. InvarDualCases[] returns the list of all dual cases handled by Invar (currently the 14 cases up to order 8 plus the algebraic case of degree 5).";

RInv::usage="RInv[g][c, r] denotes the Riemann invariant of case c and index r associated to the metric g.";
DualRInv::usage="DualRInv[g][c, r] denotes the dual Riemann invariant of case c and index r associated to the metric g.";
RPerm::usage="RPerm[g][{c, ed}, perm] denotes a Riemann invariant of case c with ed epsilon tensors (both associated to metric g), after reordering the indices from a canonical order using permutation perm.";

WInv::usage="WInv[g][c, r] denotes the Weyl invariant of case c and index r associated to the metric g.";
DualWInv::usage="DualWInv[g][c, r] denotes the dual Weyl invariant of case c and index r associated to the metric g.";
WPerm::usage="WPermRPerm[g][{Rd, ed}, perm] denotes a Weyl invariant  formed by Rd Weyl tensors and ed epsilon tensors, after reordering the indices from a canonical order using permutation perm. Contracted Weyls encode the Traceless Ricci tensor.";

RToW::usage="RToW[expr] transforms Riemann invariants into Weyl invariants in expression expr.";
WToR::usage="WToR[expr] transforms Weyl invariants into Riemann invariants in expression expr.";
RToWRules::usage="RToWRules[step, degree]";
RToWDualRules::usage="RToWDualRules[step, degree]";
WToRRules::usage="RToWRules[step, degree]";
WToRDualRules::usage="RToWDualRules[step, degree]";

RInvs::usage="RInvs[step, case] gives all independent Riemann invariants at the given step (1, 2, 3, 4, 5 or 6) and the given case. Using RInvs[g] instead of RInvs we get the particularization of the invariants to metric g.";
WInvs::usage="WInvs[step, case] gives all independent Weyl invariants at the given step (1, 2, 3, 4, 5, 6) and the given case. Using WInvs[g] instead of WInvs we get the particularization of the invariants to metric g.";
DualRInvs::usage="DualRInvs[step, case] gives all independent dual Riemann invariants at the given step (1, 2, 3, 4, 5 or 6) and the given case. Using DualRInvs[g] instead of RInvs we get the particularization of the dual invariants to metric g.";
DualWInvs::usage="DualWInvs[step, case] gives all independent dual Weyl invariants at the given step (1, 2, 3, 4, 5 or 6) and the given case. Using DualWInvs[g] instead of WInvs we get the particularization of the dual invariants to metric g.";
RInvRules::usage="RInvRules[step, case] gives the list of rules to be used to simplify Riemann non-dual invariants of the given case at the given step.";
WInvRules::usage="WInvRules[step, case] gives the list of rules to be used to simplify Weyl non-dual invariants of the given case at the given step.";
DualRInvRules::usage="DualRInvRules[step, case] gives the list of rules to be used to simplify Riemann dual invariants of the given case at the given step.";
DualWInvRules::usage="DualWInvRules[step, case] gives the list of rules to be used to simplify Weyl dual invariants of the given case at the given step.";

RemoveRInvRules::usage="RemoveRInvRules[step, case] removes the stored rules for the invariants of the given step and case.";
RemoveDualRInvRules::usage="RemoveDualRInvRules[step, case] removes the stored rules for the dual invariants of the given step and case.";

MaxIndex::usage="MaxIndex[step] gives the largest index of the non-dual invariants at the given step.";
MaxDualIndex::usage="MaxDualIndex[step] gives the largest index of the dual invariants at the given step.";

InvSimplify::usage="InvSimplify[expr, level] simplifies the invariants of the form inv[metric][case,count] in the expression expr, where inv is one of RInv, DualRInv, WInv, DualWInv. There are six possible levels of simplification: 1 means no simplification; 2 uses the cyclic identity; 3 uses the Bianchi identity; 4 commutes covariant derivatives; 5 uses dimensionally dependent (Lovelock type) identities and 6 reduces some non-dual invariants into products of dual invariants. The default level used in InvSimplify[expr] is $InvSimplifyLevel.";
InvToPerm::usage="InvToPerm[expr] converts all invariants in expr into their corresponding permutations with head RPerm or WPerm.";
PermToRiemann::usage="PermToRiemann[expr, cr] converts all permutations into their explicit tensor expressions. With cr=True contracted Riemann tensors are replaced by Ricci tensors. The default value used in PermToRiemann[expr] is $CurvatureRelations.";
InvToRiemann::usage="InvToRiemann[expr, cr] converts all invariants into their explicit tensor expressions. With cr=True contracted Riemann tensors are replaced by Ricci tensors. The default value used in InvToRiemann[expr] is $CurvatureRelations.";
RiemannToPerm::usage="RiemannToPerm[expr, g] converts all Riemann scalars of metric g into their canonical permutations. If g is a list of metrics then the command is folded over the list. RiemannToPerm[expr] is automatically converted into RiemannToPerm[expr, $Metrics].";
PermToInv::usage="PermToInv[expr] converts all permutations in expr into their correspondent invariants. It is assumed that the permutations have been already canonicalized.";
RiemannToInv::usage="RiemannToInv[expr, g] converts all Riemann scalars of metric g into their invariant form. If g is a list of metrics then the command is folded over the list. RiemannToInv[expr] is automatically converted into RiemannToInv[expr, $Metrics].";
RiemannSimplify::usage="RiemannSimplify[expr, level, cr, g] simplifies the Riemann scalars of metric g using relations up to the given level (see usage message for InvSimplify). With cr=True contracted Riemann tensors are replaced by Ricci tensors. If g is a list of metrics then the command is folded over the list. RiemannSimplify[expr] uses the defaults $InvSimplifyLevel, $CurvatureRelations and $Metrics, respectively. See usage messages for those global variables.";

dim::usage="dim is a constant symbol used in the Invar rules to denote the dimension of a manifold.";
sigma::usage="sigma is a constant symbol used in the Invar rules to denote the sign of the determinant of a metric.";

$InvSimplifyLevel::usage="$InvSimplifyLevel is a global variable specifying the default level of simplification in InvSimplify. Initially it is set to 6.";
$CurvatureRelations::usage="$CurvatureRelations is a global variable specifying the default behaviour of contracted Riemann conversion in the functions PermToRiemann, InvToRiemann and RiemannSimplify.";

RandomRiemannMonomial::usage="RandomRiemannMonomial[c, g] generates a random monomial of case c with Riemann tensors of the metric g. If not specified, g is taken to be the first metric in the list $Metrics.";

$ExpandedCommuteOrder12Q::usage="$ExpandedRulesDiagonalQ is a Boolean global variable stating whether we want to use the expanded form for the Commute step of some hard cases of order 12. Its default value is False.";

Begin Private

Begin["`Private`"]

$InvarNames=Names["xAct`Invar`*"];

sigma

?sigma

Block[{$InvarNames={}},DefConstantSymbol[sigma,PrintAs->"σ"]]

sigma/:Power[sigma,_Integer?EvenQ]=1;
sigma/:Power[sigma,_Integer?OddQ]=sigma;
sigma/:Abs[sigma]=1;

?dim

Block[{$InvarNames={}},DefConstantSymbol[dim,PrintAs->"d"]]

Main setup

Options:

$xpermQ

Options[CanonicalPerm]

Options[ToCanonical]

ReportSetOption[ContractMetric,AllowUpperDerivatives->True];

Default values for the global values of the package:

$InvSimplifyLevel=6;
$CurvatureRelations=False;

Nonstandard configuration:

These two are just aesthetical changes:

ReportSet[$PrePrint,ScreenDollarIndices];

ReportSet[$CovDFormat,"Postfix"];

JMM 2011: This is not well designed as $CurvatureRelations is not directly related to CurvatureRelations as other options are. $CurvatureRelations is information for the Riemann simplification functions, while Invar currently requires the option CurvatureRelations to be False. The proper solution is twofold:
	1) Rename $CurvatureRelations to something else.
	2) Make Invar work properly independently of the value of the option CurvatureRelations. Same thing happens with $CommuteCovDsOnScalars.
Currently we just sent a few warning messages when loading the package.

ReportSetOption[DefCovD,CurvatureRelations->False];

ReportSet[$CommuteCovDsOnScalars,False];

Formatting of invariants

tostring[list_List]:=StringJoin[ToString/@list];

tostring[{0,0,0}]

Invariants with metric:

Format[RInv[_][case_,index_]]:=Subscript["I",StringJoin[tostring[case],",",ToString[index]]];
Format[DualRInv[_][case_,index_]]:=Subscript["D",StringJoin[tostring[case],",",ToString[index]]];
Format[WInv[_][case_,index_]]:=Subscript["W",StringJoin[tostring[case],",",ToString[index]]];
Format[DualWInv[_][case_,index_]]:=Subscript["X",StringJoin[tostring[case],",",ToString[index]]];

Invariants without metric :

Format[RInv[case_,index_]]:=Subscript["I",StringJoin[tostring[case],",",ToString[index]]];
Format[DualRInv[case_,index_]]:=Subscript["D",StringJoin[tostring[case],",",ToString[index]]];
Format[WInv[case_,index_]]:=Subscript["W",StringJoin[tostring[case],",",ToString[index]]];
Format[DualWInv[case_,index_]]:=Subscript["X",StringJoin[tostring[case],",",ToString[index]]];

RInv[{0,1,3},100]^2

Permutations RPerm and WPerm are not formatted.

Reading Invar files as rules

Switch to read the Weyl section of the package:

$ReadWeyl=False;

The installation directory of Invar:

$UserBaseDirectory

Default directory for the package :

$InvarDirectory=StringJoin[$xActDirectory,"/Invar"]

$DataBaseDir="Riemann/";

RInvs

Backwards compatibility notation :

(* Invariants *)
RInv[deg_Integer,index_]:=RInv[Table[0,{deg}],index];
DualRInv[deg_Integer,index_]:=DualRInv[Table[0,{deg}],index];
(* Invariants with metric argument *)
RInv[metric_][deg_Integer,index_]:=RInv[metric][Table[0,{deg}],index];
DualRInv[metric_][deg_Integer,index_]:=DualRInv[metric][Table[0,{deg}],index];
(* Independent Riemann invariants *)
RInvs[step_Integer,deg_Integer,dim___]:=RInvs[step,Table[0,{deg}],dim];DualRInvs[step_Integer,deg_Integer,dim___]:=DualRInvs[step,Table[0,{deg}],dim];
(* Lists of rules *)
RInvRules[step_,deg_Integer,dim___]:=RInvRules[step,Table[0,{deg}],dim];
DualRInvRules[step_,deg_Integer,dim___]:=DualRInvRules[step,Table[0,{deg}],dim];

(* From Riemann to Weyl invariants *)
WInvs[step_Integer,deg_,dim___]:=RInvs[step,deg,dim]/.RInv->WInv;DualWInvs[step_Integer,deg_,dim___]:=DualRInvs[step,deg,dim]/.DualRInv->DualWInv;

The lists of independent invariants could be protected. We define this alternative form of Set:

SetAttributes[ProtSet,HoldFirst];
ProtSet[head_[args___],value_]:=Module[{prot=Unprotect[head]},Set[head[args],value];Protect[Evaluate[prot]];value];

File names

intercase[case_List]:=StringJoin@Insert[ToString/@case,"_",List/@Range[2,Length[case]]];

intercase[{0,1,3}]

filename[step:(5|6),case_,dim_]:=filename[ToString[step],case]<>"_"<>ToString[dim];
filename[step_,case_]:=StringJoin["RInv-",intercase[case],"-",ToString[step]];

filename[3,{0,1,3}]

dualfilename[5,case_,dim_]:=dualfilename["5",case]<>"_"<>ToString[dim];
dualfilename[step_,case_]:=StringJoin["DInv-",intercase[case],"-",ToString[step]];

dualfilename[5,{0,1,3},4]

Check the import of the files

Code contributed by Teake:

RInv::nodatabase="Can't find the Invar database. You can download it from http://www.xact.es/Invar/.";
CheckImport[args__]:=Quiet[Check[Import[args],Message[RInv::nodatabase];Abort[]],Import::nffil];

Reading functions. Maple notation

The permutation files are written in Maple notation. We need a collection of functions to translate the data to a format that Mathematica can understand. Note that the files are read without a metric argument.

Intermediate simple functions. This selects the first of a collection of positions:

firstof[{}]:=0;
firstof[{{a_,b_}}]:=a;
firstof[{{a_,b_},__}]:=a;

positionof[string_String,substring_String]:=firstof@StringPosition[string,substring];
movetofirst[0]:=1;
movetofirst[x_]:=x+4;

Translation from brackets to braces for lists, and from equality to rules:

replacebrackets[string_String]:=StringReplace[string,{"["->"{","]"->"}"}]
replaceequal[string_String]:=StringReplace[string,{"="->"->"}]

This is the main conversion function. It only reads what is found to the right of := and removing the last character (always a colon).

readline[string_]:=StringTake[string,{movetofirst@positionof[string," := "],StringLength[string]-1}]

Final reading functions, for Maple notation:

ReadInvarPerms[filename_,"Maple"]:=ToExpression/@replacebrackets/@readline/@CheckImport[StringJoin[$InvarDirectory,"/",filename],"Lines"];

ReadInvarRules[filename_,"Maple"]:=ToExpression/@replaceequal/@readline/@CheckImport[StringJoin[$InvarDirectory,"/",filename],"Lines"];

Reading functions. Mathematica notation

ReadInvarPerms[filename_,"Mathematica"]:=ToExpression/@CheckImport[StringJoin[$InvarDirectory,"/",filename],"Lines"];

ReadInvarRules[filename_,"Mathematica"]:=ToExpression/@CheckImport[StringJoin[$InvarDirectory,"/",filename],"Lines"];

Removing rules

RemoveRInvRules[step_Integer,degree_Integer,dim___]:=RemoveRInvRules[step,Table[0,{degree}]];
RemoveRInvRules[step_Integer,case_List,dim___]:=(
Unset[DispatchRules[RInv,step,case,dim]];
Unprotect[RInvRules];
Unset[RInvRules[step,case,dim]];
Protect[RInvRules];
);
RemoveDualRInvRules[step_Integer,degree_Integer,dim___]:=RemoveDualRInvRules[step,Table[0,{degree}]];
RemoveDualRInvRules[step_Integer,case_List,dim___]:=(
Unset[DispatchRules[DualRInv,step,case,dim]];
Unprotect[DualRInvRules];
Unset[DualRInvRules[step,case,dim]];
Protect[DualRInvRules];
);

Metric argument

We could have invariants associated to different metrics. Therefore the invariants need to carry a third argument specifying the metric:

addRdep[Rule[inv_[case_,count_],rhs_]]:=Rule[inv[metric_][case,count],rhs/.{RInv->RInv[metric],DualRInv->DualRInv[metric],WInv->WInv[metric],DualWInv->DualWInv[metric]}]

The main idea is that we shall use the dispatched rules with the metric option for the actual simplification algorithm, but the normal (Nometric) rules to generate the files for the package.

(* With a metric argument *)
RInvs[metric_][args__]:=RInvs[args]/.RInv->RInv[metric];WInvs[metric_][args__]:=WInvs[args]/.WInv->WInv[metric];
DualRInvs[metric_][args__]:=DualRInvs[args]/.DualRInv->DualRInv[metric];
DualWInvs[metric_][args__]:=DualWInvs[args]/.DualWInv->DualWInv[metric];

Note that there is no confusion between both (metric and nometric) definitions because in one case RInvs has one argument, but in the other case it has two arguments.

$DefaultDim:=DimOfVBundle@VBundleOfMetric@First@$Metrics;

Cases

The function InvarCases contains the list of cases which are currently handled by Invar:

InvarCases[]:=Flatten[Table[InvarCases[order],{order,2,14,2}],1];
InvarCases[order_Integer?EvenQ]:=InvarCases[order,Reverse@Range[order/2]];
InvarCases[order_Integer?EvenQ,list_List]:=Flatten[InvarCases[order,#]&/@list,1];

InvarCases[14]:=InvarCases[14,7];

InvarCases[2,1]={{0}};

InvarCases[4,1]={{2}};
InvarCases[4,2]={{0,0}};

InvarCases[6,1]={{4}};
InvarCases[6,2]={{0,2},{1,1}};
InvarCases[6,3]={{0,0,0}};

InvarCases[8,1]={{6}};
InvarCases[8,2]={{0,4},{1,3},{2,2}};
InvarCases[8,3]={{0,0,2},{0,1,1}};
InvarCases[8,4]={{0,0,0,0}};

InvarCases[10,1]={{8}};
InvarCases[10,2]={{0,6},{1,5},{2,4},{3,3}};
InvarCases[10,3]={{0,0,4},{0,1,3},{0,2,2},{1,1,2}};
InvarCases[10,4]={{0,0,0,2},{0,0,1,1}};
InvarCases[10,5]={{0,0,0,0,0}};

InvarCases[12,1]={{10}};
InvarCases[12,2]={{0,8},{1,7},{2,6},{3,5},{4,4}};
InvarCases[12,3]={{0,0,6},{0,1,5},{0,2,4},{1,1,4},{0,3,3},{1,2,3},{2,2,2}};
InvarCases[12,4]={{0,0,0,4},{0,0,1,3},{0,0,2,2},{0,1,1,2},{1,1,1,1}};
InvarCases[12,5]={{0,0,0,0,2},{0,0,0,1,1}};
InvarCases[12,6]={{0,0,0,0,0,0}};

InvarCases[14,7]={{0,0,0,0,0,0,0}};

InvarCases[x_]:=Throw[Print["Invalid order ",x]];
InvarCases[args__]:=Throw[Print["Invar cannot work with the cases ",{args}]];

InvarCases[]

Length[%]

InvarCases[10]

InvarCases[10,{2,3}]

indicesofcase[case_]:=4Length[case]+Plus@@case

Factorial/@(indicesofcase/@InvarCases[][[Range[47]]])

Plus@@%

N[%]

Almost Avogadro's number!

InvarDualCases[]:=Flatten[Table[InvarDualCases[order],{order,2,10,2}],1];
InvarDualCases[order_Integer?EvenQ]:=InvarDualCases[order,Reverse@Range[order/2]];
InvarDualCases[order_Integer?EvenQ,list_List]:=Flatten[InvarDualCases[order,#]&/@list,1];

InvarDualCases[10]:=InvarDualCases[10,5];

InvarDualCases[2,1]={{0}};

InvarDualCases[4,1]={{2}};
InvarDualCases[4,2]={{0,0}};

InvarDualCases[6,1]={{4}};
InvarDualCases[6,2]={{0,2},{1,1}};
InvarDualCases[6,3]={{0,0,0}};

InvarDualCases[8,1]={{6}};
InvarDualCases[8,2]={{0,4},{1,3},{2,2}};
InvarDualCases[8,3]={{0,0,2},{0,1,1}};
InvarDualCases[8,4]={{0,0,0,0}};

InvarDualCases[10,5]={{0,0,0,0,0}};

InvarDualCases[x_]:=Throw[Print["Invalid order ",x]];
InvarDualCases[args__]:=Throw[Print["Invar cannot work with the dual cases ",{args}]];

InvarDualCases[]

Length[%]

Maximum degree and index

These are the numbers (directly copied from the files) of canonical invariants for each degree (not containing products):

?MaxIndex

MaxIndex[n_Integer?Positive]:=MaxIndex[Table[0,{n}]];

MaxIndex[{}]:=0;

MaxIndex[{0}]:=1;

MaxIndex[{0,0}]:=3;
MaxIndex[{2}]:=2;

MaxIndex[{0,0,0}]:=9;
MaxIndex[{0,2}]:=12;
MaxIndex[{1,1}]:=12;
MaxIndex[{4}]:=12;

MaxIndex[{0,0,0,0}]:=38;
MaxIndex[{0,0,2}]:=99;
MaxIndex[{0,1,1}]:=125;
MaxIndex[{0,4}]:=126;
MaxIndex[{1,3}]:=138;
MaxIndex[{2,2}]:=86;
MaxIndex[{6}]:=105;

MaxIndex[{0,0,0,0,0}]:=204;
MaxIndex[{0,0,0,2}]:=1020;
MaxIndex[{0,0,1,1}]:=1749;
MaxIndex[{0,0,4}]:=1473;
MaxIndex[{0,1,3}]:=3099;
MaxIndex[{0,2,2}]:=1622;
MaxIndex[{1,1,2}]:=1617;
MaxIndex[{0,6}]:=1665;
MaxIndex[{1,5}]:=1770;
MaxIndex[{2,4}]:=1746;
MaxIndex[{3,3}]:=962;
MaxIndex[{8}]:=1155;

MaxIndex[{0,0,0,0,0,0}]:=1613;
MaxIndex[{0,0,0,0,2}]:=12722;
MaxIndex[{0,0,0,1,1}]:=27022;
MaxIndex[{0,0,0,4}]:=19617;
MaxIndex[{0,0,1,3}]:=60984;
MaxIndex[{0,0,2,2}]:=30974;
MaxIndex[{0,1,1,2}]:=62465;
MaxIndex[{1,1,1,1}]:=5606;
MaxIndex[{0,0,6}]:=25590;
MaxIndex[{0,1,5}]:=53160;
MaxIndex[{0,2,4}]:=52764;
MaxIndex[{1,1,4}]:=27396;
MaxIndex[{0,3,3}]:=27024;
MaxIndex[{1,2,3}]:=54654;
MaxIndex[{2,2,2}]:=9104;
MaxIndex[{0,8}]:=25515;
MaxIndex[{1,7}]:=26670;
MaxIndex[{2,6}]:=26460;
MaxIndex[{3,5}]:=26670;
MaxIndex[{4,4}]:=13607;
MaxIndex[{10}]:=15120;

MaxIndex[{0,0,0,0,0,0,0}]:=16532;

MaxIndex[{0,0,0,0,0,0,0,0}]:=217395;

MaxIndex[{0,0,0,0,0,0,0,0,0}]:=3406747;

MaxIndex[case_]:=Throw@Print["Case ",case," not included."];

allnonduals=Apply[Plus,MaxIndex/@InvarCases[]]

Duals.

MaxDualIndex[n_Integer]:=MaxDualIndex[Table[0,{n}]];

MaxDualIndex[{}]=0;

MaxDualIndex[{0}]=1;

MaxDualIndex[{0,0}]=4;
MaxDualIndex[{2}]=3;

MaxDualIndex[{0,0,0}]=27;
MaxDualIndex[{0,2}]=58;
MaxDualIndex[{1,1}]=36;
MaxDualIndex[{4}]=32;

MaxDualIndex[{0,0,0,0}]=232;
MaxDualIndex[{0,0,2}]=967;
MaxDualIndex[{0,1,1}]=1047;
MaxDualIndex[{0,4}]=876;
MaxDualIndex[{1,3}]=920;
MaxDualIndex[{2,2}]=478;
MaxDualIndex[{6}]=435;

MaxDualIndex[{0,0,0,0,0}]=2582;

MaxDualIndex[{0,0,0,0,0,0}]=35090;

MaxDualIndex[{0,0,0,0,0,0,0}]=558323;

MaxDualIndex[case_]:=Throw[Print["Dual case", case," not included."]];

allduals=Apply[Plus,MaxDualIndex/@InvarDualCases[]]

allnonduals+allduals

More on RInvs

Construct Invs lists for level 1. There is a call to MaxIndex and hence a check of maximum degree.

RInvs[1,case_List]:=ProtSet[RInvs[1,case],Map[RInv[case,#]&,Range[MaxIndex[case]]]];
DualRInvs[1,case_List]:=ProtSet[DualRInvs[1,case],Map[DualRInv[case,#]&,Range[MaxDualIndex[case]]]];
WInvs[1,case_List]:=ProtSet[WInvs[1,case],Map[WInv[case,#]&,Range[MaxIndex[case]]]];
DualWInvs[1,case_List]:=ProtSet[DualWInvs[1,case],Map[DualWInv[case,#]&,Range[MaxDualIndex[case]]]];

DispatchRules

Dispatch rules per degree. Again, we do not need separate dispatched rules for the Weyl invariants:

DispatchRules[RInv,step_,case_,dim___]:=Block[{sigma},DispatchRules[RInv,step,case,dim]=Dispatch[RInvRules[step,case,dim]]];
DispatchRules[DualRInv,step_,case_,dim___]:=DispatchRules[DualRInv,step,case,dim]=Dispatch[DualRInvRules[step,case,dim]];
DispatchRules[WInv,step_,case_,dim___]:=Block[{sigma},DispatchRules[WInv,step,case,dim]=Dispatch[WInvRules[step,case,dim]]];
DispatchRules[DualWInv,step_,case_,dim___]:=DispatchRules[DualWInv,step,case,dim]=Dispatch[DualWInvRules[step,case,dim]];

Read the Canon permutations for the invariants: step 1

First we read the correspondence between invariants and permutations. Note that the symbols Invars and DualInvars are never loaded into Mathematica.

RInvRules[1,case_List]:=Module[{result},
MaxIndex[case];
Print["Reading InvRules for step 1 and case ",case];result=Inner[Rule,tmphead@@RInvs[1,case],tmphead@@ReadInvarPerms[$DataBaseDir<>"1/"<>filename[1,case],"Maple"],List];
ProtSet[RInvRules[1,case],result]
];

DualRInvRules[1,case_List]:=Module[{result},
MaxDualIndex[case];
Print["Reading DualInvRules for step 1 and case ",case];result=Inner[Rule,tmphead@@DualRInvs[1,case],tmphead@@ReadInvarPerms[$DataBaseDir<>"1/"<>dualfilename[1,case],"Maple"],List];
ProtSet[DualRInvRules[1,case],result]
];

Dispatch rules per degree:

reversecycles[inv_->{cycles__}]:=Cycles[cycles]->inv;

DispatchRInvToPermRules[case_List]:=DispatchRInvToPermRules[case]=Dispatch[RInvRules[1,case]];
DispatchPermToRInvRules[case_List]:=DispatchPermToRInvRules[case]=Dispatch[reversecycles/@RInvRules[1,case]];
DispatchDualRInvToPermRules[case_List]:=DispatchDualRInvToPermRules[case]=Dispatch[DualRInvRules[1,case]];
DispatchPermToDualRInvRules[case_List]:=DispatchPermToDualRInvRules[case]=Dispatch[reversecycles/@DualRInvRules[1,case]];

The files for WInvariants are the same as those for RInvariants. We do not create dispatched rules for the WInvariants.

WInvRules[1,case_]:=RInvRules[1,case]/.RInv->WInv;
DualWInvRules[1,case_]:=DualRInvRules[1,case]/.DualRInv->DualWInv;

Examples with invariants:

RInvRules[1,3]

Length[RInvRules[1,7]]//Timing

Length[RInvRules[1,7]]//Timing

RInvRules[1,7][[-1]]

Examples with dual invariants:

DualRInvRules[1,3]//Length

Errors:

Catch@RInvRules[1,10]

Catch@DualWInvRules[1,8]

Read cyclic rules: step 2

Now we read the rules given by the cyclic symmetry.

Invariants:

RInvs[2,case_]:=ProtSet[RInvs[2,case],Complement[RInvs[1,case],First/@RInvRules[2,case]]];
DualRInvs[2,case_]:=ProtSet[DualRInvs[2,case],Complement[DualRInvs[1,case],First/@DualRInvRules[2,case]]];

Rules:

RInvRules[2,case_List]:=Module[{result},
MaxIndex[case];
Print["Reading InvRules for step 2 and case ",case];
result=ReadInvarRules[$DataBaseDir<>"2/"<>filename[2,case],"Mathematica"];
ProtSet[RInvRules[2,case],result]
];

DualRInvRules[2,case_List]:=Module[{result},
MaxDualIndex[case];
Print["Reading DualInvRules for step 2 and case ",case];result=ReadInvarRules[$DataBaseDir<>"2/"<>dualfilename[2,case],"Mathematica"];
ProtSet[DualRInvRules[2,case],result]
];

Examples:

Timing[Length[result=RInvRules[2,{0,1,3}]]]

result[[1]]

result[[-1]]

RInvs[2,{0,1,3}]//Length

Timing[Length[result=DualRInvRules[2,4]]]

result[[1]]

Catch@DualRInvRules[2,8]

DualRInvRules[2,{1,3}][[-1]]

Read Bianchi rules: step 3

Invariants:

RInvs[3,case_]:=ProtSet[RInvs[3,case],Complement[RInvs[2,case],First/@RInvRules[3,case]]];
DualRInvs[3,case_]:=ProtSet[DualRInvs[3,case],Complement[DualRInvs[2,case],First/@DualRInvRules[3,case]]];

Rules:

RInvRules[3,case_List]:=Module[{result},
MaxIndex[case];
Print["Reading InvRules for step 3 and case ",case];
result=ReadInvarRules[$DataBaseDir<>"3/"<>filename[3,case],"Mathematica"];
ProtSet[RInvRules[3,case],result]
];

DualRInvRules[3,case_List]:=Module[{result},
MaxDualIndex[case];
Print["Reading DualInvRules for step 3 and case ",case];result=ReadInvarRules[$DataBaseDir<>"3/"<>dualfilename[3,case],"Mathematica"];
ProtSet[DualRInvRules[3,case],result]
];

Examples:

Timing[Length[RInvRules[3,{0,1,3}]]]

Timing[Length[RInvRules[3,{0,1,3}]]]

RInvRules[3,{0,1,3}]//Length

RInvs[3,{0,1,3}]//Length

Duals:

Timing[Length[DualRInvRules[3,{2,2}]]]

Timing[Length[DualRInvRules[3,{2,2}]]]

DualRInvRules[3,{2,2}]//Length

RInvs[3,{2,2}]//Length

Read Commute rules: step 4

Invariants:

RInvs[4,case_]:=ProtSet[RInvs[4,case],Complement[RInvs[3,case],First/@RInvRules[4,case]]];
DualRInvs[4,case_]:=ProtSet[DualRInvs[4,case],Complement[DualRInvs[3,case],First/@DualRInvRules[4,case]]];

Rules:

$ExpandedCommuteOrder12Q=False;

extendedQ[case_]:=$ExpandedCommuteOrder12Q/;MemberQ[InvarCases[12,{1,2,3}],case];
extendedQ[case_]:=True;

extendedQ[{8}]

extendedQ[{1,1,1,1}]

extendedQ[{0,8}]

RInvRules[4,case_List]:=Module[{result,eQ=If[extendedQ[case],"","NE"]},
MaxIndex[case];
Print["Reading "<>eQ<>"InvRules for step 4 and case ",case];
result=ReadInvarRules[$DataBaseDir<>"4/"<>filename[4,case]<>eQ,"Mathematica"];
ProtSet[RInvRules[4,case],result]
];

DualRInvRules[4,case_List]:=Module[{result},
MaxDualIndex[case];
Print["Reading DualInvRules for step 4 and case ",case];result=ReadInvarRules[$DataBaseDir<>"4/"<>dualfilename[4,case],"Mathematica"];
ProtSet[DualRInvRules[4,case],result]
];

Examples:

Timing[Length[RInvRules[4,{0,1,3}]]]

Timing[Length[RInvRules[4,{0,1,3}]]]

RInvRules[4,{0,1,3}]//Length

RInvs[4,{0,1,3}]//Length

Read Dim-dep rules: step 5

Invariants:

RInvs[5,case_]:=RInvs[5,case,$DefaultDim];
RInvs[5,case_,dim_]:=ProtSet[RInvs[5,case,dim],Complement[RInvs[4,case],First/@RInvRules[5,case,dim]]];
DualRInvs[5,case_]:=DualRInvs[5,case,$DefaultDim];
DualRInvs[5,case_,dim_]:=ProtSet[DualRInvs[5,case,dim],Complement[DualRInvs[4,case],First/@DualRInvRules[5,case,dim]]];

Rules:

RInvRules[5,case_List]:=RInvRules[5,case,$DefaultDim];
RInvRules[5,case_List,dim_]:=Module[{result},
MaxIndex[case];
Print["Reading InvRules for step 5, case ",case," and dimension ",dim];
result=ReadInvarRules[$DataBaseDir<>"5_"<>ToString[dim]<>"/"<>filename[5,case,dim],"Mathematica"];
ProtSet[RInvRules[5,case,dim],result]
];

DualRInvRules[5,case_List]:=DualRInvRules[5,case,$DefaultDim];
DualRInvRules[5,case_List,dim_]:=Module[{result},
MaxDualIndex[case];
Print["Reading DualInvRules for step 5 and case ",case];result=ReadInvarRules[$DataBaseDir<>"5_"<>ToString[dim]<>"/"<>dualfilename[5,case,dim],"Mathematica"];
ProtSet[DualRInvRules[5,case,dim],result]
];

Examples:

Timing[Length[RInvRules[5,{0,1,3},4]]]

Timing[Length[RInvRules[5,{0,1,3},4]]]

RInvs[5,{0,1,3},4]//Length

Read metric-signature rules: step 6

Finally we read the metric-signature dependent relations. There are only 25 independent objects:

Invariants:

RInvs[6,case_]:=RInvs[6,case,$DefaultDim];
RInvs[6,case_,dim_]:=ProtSet[RInvs[6,case,dim],Complement[RInvs[5,case,dim],First/@RInvRules[6,case,dim]]];
DualRInvs[6,case_]:=DualRInvs[6,case,$DefaultDim];
DualRInvs[6,case_,dim_]:=ProtSet[DualRInvs[6,case,dim],DualRInvs[5,case,dim]];

Higher degrees:

RInvRules[6,case_List]:=RInvRules[6,case,$DefaultDim];
RInvRules[6,case_List,dim_]:=Block[{sigma},
Module[{result},
MaxIndex[case];
Print["Reading InvRules for step 6, case ",case," and dimension ",dim];
result=ReadInvarRules[$DataBaseDir<>"6_"<>ToString[dim]<>"/"<>filename[6,case,dim],"Mathematica"];
ProtSet[RInvRules[6,case,dim],result]
]];

Numbers

Diagonals 1 to 5

memory=MemoryInUse[]

Read first five diagonals:

$DefaultDim=4;

Table[Length[RInvs[step,InvarCases[][[case]]]],{step,6},{case,34}]

Transpose[Join[{InvarCases[][[Range[34]]]},%]]//MatrixForm

MemoryInUse[]/1024./1000

Reading the first five diagonals takes 50 Mbytes. With eight more cases from the sixth diagonal it takes 350 Mbytes.

(MemoryInUse[]-memory)/1024./1000

Analysis of numbers for the first five diagonals:
	There are originally 18550 invariants at step 1.
	We have 9133 rules at step 2 in terms of 6673 of the 9717 free invariants at step 2.
	We have 5839 rules at step 3 in terms of 3801 of the 3878 free invariants at step 3.
	We have 3255 rules at step 4 in terms of 618 of the 623 free invariants at step 4.
	We have 287 rules at step 5 in terms of 308 of the 336 free invariants at step 5.
	We have 6 rules at step 6 in terms of 57 of the 330 free invariants at step 6 and 6 free dual invariants.

Total[Length[RInvs[1,#]]&/@InvarCases[][[Range[26]]]]

Total[Length[RInvs[2,#]]&/@InvarCases[][[Range[26]]]]

Total[Length[RInvs[3,#]]&/@InvarCases[][[Range[26]]]]

Total[Length[RInvs[4,#]]&/@InvarCases[][[Range[26]]]]

Total[Length[RInvs[5,#]]&/@InvarCases[][[Range[26]]]]

Total[Length[RInvs[6,#]]&/@InvarCases[][[Range[26]]]]

Check that all previous rules have been used:

check[step_,cases_]:=Do[
case=cases[[count]];
rules=RInvRules[step,case];
right=Union@Cases[Last/@rules,RInv[_,_],{0,Infinity}];
inters=Intersection[right,knowns];
If[inters=!={},
Print["There are ",Length[inters]," rules to use in case ",case],
Print["Checked case ",case]];
knowns=Join[knowns,First/@rules],
{count,1,Length[cases]}]

knowns={};

These are the cases to check:

cases=Append[InvarCases[][[Range[34]]],InvarCases[][[48]]]

check[2,cases]

check[3,cases]

check[4,cases]

check[5,cases]

check[6,cases]

Length[knowns]

MemoryInUse[]/1000./1024

Diagonal 6 non-expanded

Now we read the sixth diagonal (Commute non-expanded):

Table[Length[RInvs[step,InvarCases[][[case]]]],{step,5},{case,35,47}]

Transpose[Join[{InvarCases[][[Range[35,47]]]},%]]//MatrixForm

MemoryInUse[]/1024./1000

Analysis of numbers for the first six diagonals:
	There are originally 623587 invariants at step 1.
	We have 365747 rules at step 2 in terms of 212119 of the 257840 free invariants at step 2.
	We have 169204 rules at step 3 in terms of 84409 of the 84756 free invariants at step 3.
	We have ? rules at step 4 in terms of ? of the 6368 free invariants at step 4.
	We have ? rules at step 5 in terms of 2129 of the 2202 free invariants at step 5.

Total[Length[RInvs[1,#]]&/@InvarCases[][[Range[47]]]]

Total[Length[RInvs[2,#]]&/@InvarCases[][[Range[47]]]]

Length@Union@Cases[Map[Last,Flatten[RInvRules[2,#]&/@InvarCases[][[Range[47]]]]],RInv[__],{0,Infinity}]

Total[Length[RInvs[3,#]]&/@InvarCases[][[Range[47]]]]

Length@Union@Cases[Map[Last,Flatten[RInvRules[3,#]&/@InvarCases[][[Range[47]]]]],RInv[__],{0,Infinity}]

Total[Length[RInvs[4,#]]&/@InvarCases[][[Range[47]]]]

Total[Length[RInvs[5,#]]&/@InvarCases[][[Range[47]]]]

Length@Union@Cases[Map[Last,Flatten[RInvRules[5,#]&/@InvarCases[][[Range[47]]]]],RInv[__],{0,Infinity}]

Checks (now we include the algebraic degree 7 case):

knowns={};

check[2,InvarCases[]]

Length[knowns]

check[3,InvarCases[]]

Length[knowns]

The sixth diagonal is not expanded:

check[4,InvarCases[]]

Length[knowns]

check[5,InvarCases[]]

Length[knowns]

MemoryInUse[]/1024./1000

Numbers for the abstract:

Total[Length[RInvs[1,#]]&/@InvarCases[]]

Total[Length[DualRInvs[1,#]]&/@InvarCases[][[Range[15]]]]

r2=Total[Length[RInvRules[2,#]]&/@InvarCases[]]

dr2=Total[Length[DualRInvRules[2,#]]&/@InvarCases[][[Range[15]]]]

r3=Total[Length[RInvRules[3,#]]&/@InvarCases[]]

dr3=Total[Length[DualRInvRules[3,#]]&/@InvarCases[][[Range[15]]]]

r4=Total[Length[RInvRules[4,#]]&/@InvarCases[]]

dr4=Total[Length[DualRInvRules[4,#]]&/@InvarCases[][[Range[15]]]]

r5=Total[Length[RInvRules[5,#]]&/@InvarCases[]]

dr5=Total[Length[DualRInvRules[5,#]]&/@InvarCases[][[Range[15]]]]

{allnonduals,allduals,allnonduals+allduals,r2+r3+r4+r5+51,dr2+dr3+dr4+dr5,r2+r3+r4+r5+51+dr2+dr3+dr4+dr5}

Number of independent invariants:

allnonduals-(r2+r3+r4+r5)

%-51

allduals-(dr2+dr3+dr4+dr5)

Diagonal 6 expanded

Now we read the sixth diagonal (Commute expanded):

$ExpandedCommuteOrder12Q=True;

Table[Length[RInvs[step,InvarCases[][[case]]]],{step,5},{case,35,47}]

MemoryInUse[]

Transpose[Join[{InvarCases[][[Range[35,47]]]},%]]//MatrixForm

MemoryInUse[]/1024./1000

Analysis of numbers for the first six diagonals, only steps 2 and 3:
	There are originally 623587 invariants at step 1.
	We have 365747 rules at step 2 in terms of 212119 of the 257840 free invariants at step 2.
	We have 169204 rules at step 3 in terms of 84409 of the 84756 free invariants at step 3.
	We have ? rules at step 4 in terms of ? of the 6368 free invariants at step 4.
	We have ? rules at step 5 in terms of 2129 of the 2202 free invariants at step 5.

Total[Length[RInvs[1,#]]&/@InvarCases[][[Range[47]]]]

Total[Length[RInvs[2,#]]&/@InvarCases[][[Range[47]]]]

Length@Union@Cases[Map[Last,Flatten[RInvRules[2,#]&/@InvarCases[][[Range[47]]]]],RInv[__],{0,Infinity}]

Total[Length[RInvs[3,#]]&/@InvarCases[][[Range[47]]]]

Length@Union@Cases[Map[Last,Flatten[RInvRules[3,#]&/@InvarCases[][[Range[47]]]]],RInv[__],{0,Infinity}]

Total[Length[RInvs[4,#]]&/@InvarCases[][[Range[47]]]]

Total[Length[RInvs[5,#]]&/@InvarCases[][[Range[47]]]]

Length@Union@Cases[Map[Last,Flatten[RInvRules[5,#]&/@InvarCases[][[Range[47]]]]],RInv[__],{0,Infinity}]

Checks (now we include the algebraic degree 7 case):

knowns={};

check[2,InvarCases[]]

Length[knowns]

check[3,InvarCases[]]

Length[knowns]

The sixth diagonal is now expanded:

check[4,InvarCases[]]

Length[knowns]

check[5,InvarCases[]]

Length[knowns]

MemoryInUse[]/1024./1000

Dual invariants

Now read dual invariants:

dualtable=Table[Length[DualRInvs[step,InvarCases[][[case]]]],{step,5},{case,15}]

Transpose[Join[{InvarCases[][[Range[15]]]},dualtable]]//MatrixForm

Total/@dualtable

First[%]-Last[%]

There are 7664 relations for 7698 canonical dual invariants.

MemoryInUse[]/1024/1000.

Length[knowns]

checkdual[step_,cases_]:=Do[
case=cases[[count]];
rules=DualRInvRules[step,case];
right=Union@Cases[Last/@rules,RInv[_,_]|DualRInv[_,_],{0,Infinity}];
inters=Intersection[right,knowns];
If[inters=!={},
Print["There are ",inters," rules to use in case ",case],
Print["Checked case ",case]];
knowns=Join[knowns,First/@rules],
{count,1,Length[cases]}]

checkdual[2,InvarCases[][[Range[15]]]]

checkdual[3,InvarCases[][[Range[15]]]]

checkdual[4,InvarCases[][[Range[15]]]]

checkdual[5,InvarCases[][[Range[15]]]]

Algebraic invariants

table=Table[Length[RInvs[step,degree]],{step,1,5},{degree,1,7}]

TableForm[Transpose[table]]

Total/@table

First[%]-Last[%]

dualtable=Table[Length[DualRInvs[step,degree]],{step,1,5},{degree,1,5}]

TableForm[Transpose[dualtable]]

Total/@dualtable

First[%]-Last[%]

Number of canonical forms handled: 21246

18400+2846

Number of independent objects: 25. Therefore we have 21221 relations.

Weyl

Riemann <-> Weyl conversion rules. INCOMPLETE. DEACTIVATED

DispatchRToWRules[step_Integer,case_]:=DispatchRToWRules[step,case]=Dispatch[RToWRules[step,case]];
DispatchRToWDualRules[step_Integer,case_]:=DispatchRToWDualRules[step,case]=Dispatch[RToWDualRules[step,case]];
DispatchWToRRules[step_Integer,case_]:=DispatchWToRRules[step,case]=Dispatch[WToRRules[step,case]];
DispatchWToRDualRules[step_Integer,case_]:=DispatchWToRDualRules[step,case]=Dispatch[WToRDualRules[step,case]];

fileread[RToWRules,step_Integer,case_]:=StringJoin["RToW/RToW",ToString[step],"_R",ToString[case],".mpl"];
fileread[RToWDualRules,step_Integer,case_]:=StringJoin["RToW/RToW",ToString[step],"_D",ToString[case],".mpl"];
fileread[WToRRules,step_Integer,case_]:=StringJoin["WToR/WToR",ToString[step],"_R",ToString[case],".mpl"];
fileread[WToRDualRules,step_Integer,case_]:=StringJoin["WToR/WToR",ToString[step],"_D",ToString[case],".mpl"];
read[file_,step_Integer,case_]:=Module[{},
Print["Reading ",file," for step ",step," and case ",case];
ReadInvarRules[fileread[file,step,case],"Mathematica"]]

RToWRules[step_Integer,case_]:=ProtSet[RToWRules[step,case],read[RToWRules,step,case]];
RToWDualRules[step_Integer,case_]:=ProtSet[RToWDualRules[step,case],read[RToWDualRules,step,case]];
WToRRules[step_Integer,case_]:=ProtSet[WToRRules[step,case],read[WToRRules,step,case]];
WToRDualRules[step_Integer,case_]:=ProtSet[WToRDualRules[step,case],read[WToRDualRules,step,case]];

Weyl simplification: levels 3 and 4. INCOMPLETE. DEACTIVATED

These two files are still incomplete for degree 6 and 7:

DispatchWInvRules3[case_]:=DispatchWInvRules3[case]=Dispatch[WInvRules[3,case]];
DispatchDualWInvRules3[case_]:=DispatchDualWInvRules3[case]=Dispatch[DualWInvRules[3,case]];
DispatchWInvRules4[case_]:=DispatchWInvRules4[case]=Dispatch[WInvRules[4,case]];
DispatchDualWInvRules4[case_]:=DispatchDualWInvRules4[case]=Dispatch[DualWInvRules[4,case]];

fileread[WInvRules,level_,case_]:=StringJoin["Weyl/WInvar",ToString[level],"_R",ToString[case],".mpl"];
fileread[DualWInvRules,level_,case_]:=StringJoin["Weyl/WInvar",ToString[level],"_D",ToString[case],".mpl"];
WInvRules[3,case_]:=ProtSet[WInvRules[3,case],read[WInvRules,3,case]];
DualWInvRules[3,case_]:=ProtSet[DualWInvRules[3,case],read[DualWInvRules,3,case]];
WInvRules[4,case_]:=ProtSet[WInvRules[4,case],read[WInvRules,4,case]];
DualWInvRules[4,case]:={};

The ten main functions

InvSimplify

This function takes an expression containing RInv and DualRInv objects (or their Weyl counterparts) and simplify them to the level given by the user. Note that there is no use of permutations.

?InvSimplify

User driver: conversion to InvSimplify1, which knows the dimension of the vbundle of the metric and already acts on a single invariant. Note that the most external function of InvSimplify is Expand, and not Simplify.

InvSimplify[expr_]:=InvSimplify[expr,$InvSimplifyLevel];
InvSimplify[expr_,1]:=expr;
InvSimplify[expr_,step_]:=Expand[expr/.{inv:(RInv|DualRInv|WInv|DualWInv)[metric_][_,_]:>Block[{dim=DimOfVBundle@VBundleOfMetric@metric,sigma=SignDetOfMetric[metric]},InvSimplify1[inv,step,dim]]}];
SetNumberOfArguments[InvSimplify,{1,2}];

InvSimplify1 receives one invariant. Initial checks. If correct, move on to InvSimplify2.

InvSimplify1[inv_,1,_]:=inv;
InvSimplify1[(DualRInv|DualWInv)[_][_,_],_,dim_Integer]:=Print["Duals are only handled in dimension 4."]/;dim=!=4;
InvSimplify1[_,5|6,dim_Integer]:=Print["Levels 5 and 6 of simplication are only possible in dimension 4."]/;dim=!=4;
InvSimplify1[inv:_[metric_][_,_],step:(2|3|4|5|6),dim_]:=Expand[InvSimplify56[InvSimplify234[inv,step],step]/.{RInv->RInv[metric],DualRInv->DualRInv[metric],WInv->WInv[metric],DualWInv->DualWInv[metric]}];
InvSimplify1[inv_,step_,dim_]:=Print["Invalid level of simplification: ",step];

InvSimplify234[RInv[metric_][case_,count_],step_]:=Module[{newexpr=RInv[case,count]},
If[FreeQ[RInvs[step,case],newexpr],
newexpr=newexpr/.DispatchRules[RInv,2,case];
newexpr=If[step≥3,newexpr/.DispatchRules[RInv,3,case],newexpr];
newexpr=If[step≥4,newexpr/.DispatchRules[RInv,4,case],newexpr]];
newexpr
];
InvSimplify234[DualRInv[metric_][case_,count_],step_]:=
Module[{newexpr=DualRInv[case,count]},
If[FreeQ[DualRInvs[step,case],newexpr],
newexpr=newexpr/.DispatchRules[DualRInv,2,case];
newexpr=If[step≥3,newexpr/.DispatchRules[DualRInv,3,case],newexpr];
newexpr=If[step≥4,newexpr/.DispatchRules[DualRInv,4,case],newexpr]];
newexpr
];
InvSimplify234[WInv[metric_][case_,count_],step_]:=
If[MemberQ[RInvs[step,case],RInv[case,count]],
WInv[case,count],
Module[{newexpr=RInv[case,count]/.DispatchRules[RInv,2,case]/.RInv->WInv},
newexpr=If[step≥3,newexpr/.DispatchRules[WInv,3,case],newexpr];
newexpr=If[step≥4,newexpr/.DispatchRules[WInv,4,case],newexpr];
newexpr]
];
InvSimplify234[DualWInv[metric_][case_,count_],step_]:=If[MemberQ[DualRInvs[step,case],DualRInv[case,count]],
DualWInv[case,count],
Module[{newexpr=DualRInv[case,count]/.DispatchRules[DualRInv,2,case]/.DualRInv->DualWInv},
newexpr=If[step≥3,newexpr/.DispatchRules[DualWInv,3,case],newexpr];
newexpr=If[step≥4,newexpr/.DispatchRules[DualWInv,4,case],newexpr];
newexpr]
];

RulesFor[inv_[case_,_],step_]:=DispatchRules[inv,step,case,dim]

InvSimplify56[expr_,step_/;step<5]:=expr;
InvSimplify56[expr_,5]:=expr/.{inv:(RInv|DualRInv|WInv|DualWInv)[case_,_]:>(inv/.RulesFor[inv,5])};
InvSimplify56[expr_,6]:=InvSimplify56[expr,5]/.{inv:(RInv|WInv)[case_,_]:>(inv/.RulesFor[inv,6])};

InvToPerm

?InvToPerm

With the Canon permutations we need to use the Cycles notation.

toCycles[cycles:{___List}]:=Apply[Cycles,cycles];

Focus InvToPerm on single invariants:

InvToPerm[expr_]:=expr/.inv:(RInv|DualRInv|WInv|DualWInv)[_][_,_]:>InvToPerm[inv];
InvToPerm[RInv[metric_][case_,count_]]:=RPerm[metric][{case,0},toCycles[RInv[case,count]/.DispatchRInvToPermRules[case]]];
InvToPerm[DualRInv[metric_][case_,count_]]:=RPerm[metric][{case,1},toCycles[DualRInv[case,count]/.DispatchDualRInvToPermRules[case]]];
InvToPerm[WInv[metric_][case_,count_]]:=WPerm[metric][{case,0},toCycles[RInv[case,count]/.DispatchRInvToPermRules[case]]];
InvToPerm[DualWInv[metric_][case_,count_]]:=WPerm[metric][{case,1},toCycles[DualRInv[case,count]/.DispatchDualRInvToPermRules[case]]];
SetNumberOfArguments[InvToPerm,1];

There is nothing to decide in this function. All invariants will be converted into permutations, unless the degree and index supplied is incorrect.

PermToRiemann

We need to introduce indices. This is done through these two auxiliary functions. Note that they are defined with remembered values, and so they are computed just once. We specify 2*9 because that is the most complicated case we shall need:

upindices[metric_]:=upindices[metric]=GetIndicesOfVBundle[VBundleOfMetric[metric],18];
canon[metric_]:=canon[metric]=Flatten@Transpose[{upindices[metric],Minus/@upindices[metric]}];

?PermToRiemann

Any permutation can be converted into a Riemann expression using PermToRiemann. The first point is translating the permutation into Images notation:

translate[RPerm[metric_][{case_,dege_},perm_]]:=First@TranslatePerm[perm,{Images,4Length[case]+Plus@@case+4dege}];

RPerm[metric_][info_,-perm_]:=-RPerm[metric][info,perm];
WPerm[metric_][info_,-perm_]:=-WPerm[metric][info,perm];

We construct a function which transforms a given permutation into a Riemann (or Weyl) expression:

riemannof[metric_,0,indices_,{other___,n_}]:=With[{CD=CovDOfMetric[metric]},Apply[Composition,CD/@indices[[Range[4+n,5,-1]]]][Apply[Riemann[CD],Take[indices,4]]]riemannof[metric,0,Drop[indices,4+n],{other}]];
riemannof[metric_,1,indices_,case_]:=riemannof[metric,0,Drop[indices,-4],case]Apply[epsilon[metric],Take[indices,-4]];
riemannof[_,_,{},{}]:=1;

PermToRiemann[rperm_]:=PermToRiemann[rperm,$CurvatureRelations];
(* On Riemann permutations *)
PermToRiemann[rperm:RPerm[metric_][{case_,dege_},_],cr_]:=ToRicci[metric,cr][riemannof[metric,dege,canon[metric][[translate[rperm]]],Reverse[case]]];
(* On Weyl permutations *)
PermToRiemann[wperm:WPerm[metric_][{_,_},_],_]:=ToWeyl[PermToRiemann[wperm/.WPerm->RPerm,True],CovDOfMetric[metric]];
(* On general expressions *)
PermToRiemann[expr_,cr_]:=expr/.perm:(RPerm|WPerm)[_][__]:>PermToRiemann[perm,cr];
SetNumberOfArguments[PermToRiemann,{1,2}];

It is returned with canon indices, not with dollar-indices. That is why we need to wrap it with Scalar. There is a second argument in PermToRiemman to specify whether the contracted Riemann must be replaced by Ricci:

ToRicci[metric_,False][expr_]:=Scalar[expr];
ToRicci[metric_,True][expr_]:=Scalar@ReplaceDummies[ContractMetric[expr/.CurvatureRelations[CovDOfMetric[metric],Riemann],metric]/.CurvatureRelations[CovDOfMetric[metric],Ricci],IndexList@@upindices[metric]];
ToWeyl[expr_,cd_]:=expr/.{Riemann[cd]:>Weyl[cd],Ricci[cd]:>TFRicci[cd]};

Only cr (True or False, with default True) needs to be decided in this function.

InvToRiemann

The invariant Inv[metric][n,m] is converted into a Riemann expression using the function InvToRiemann1:

InvToRiemann1[cr_][inv:(RInv|DualRInv|WInv|DualWInv)[_][_,_]]:=PermToRiemann[InvToPerm[inv],cr];

This intermediate function is in charge of detecting these very frequent cases:

InvToRiemann1[False][(RInv|WInv)[metric_][{0},1]]:=Scalar[Riemann[CovDOfMetric[metric]]@@canon[metric][[{1,3,2,4}]]];
InvToRiemann1[True][(RInv|WInv)[metric_][{0},1]]:=RicciScalar[CovDOfMetric[metric]][];

Main translation function. Note that there is no ToCanonical. We assume that the permutations of the invariants have been already canonicalized:

?InvToRiemann

InvToRiemann[expr_]:=InvToRiemann[expr,$CurvatureRelations];
InvToRiemann[expr_,cr_]:=expr/.inv:(RInv|DualRInv|WInv|DualWInv)[_][_,_]:>InvToRiemann1[cr][inv];
SetNumberOfArguments[InvToRiemann,{1,2}];

RiemannToPerm

Given any **monomial** Riemann expression, we need to know the degree of that expression in terms of Riemann's and Ricci's, and the degree in terms of epsilons. Note that the three objects are given by the metric:

addone[head_[n_Integer]]:=head[n+1];
addone[casePlus[{head_,n_Integer}]]:=casePlus[{head,n+1}];

SetAttributes[casePlus,{Flat,Orderless}]
casePlus[{head_,n__},c___,{head_,m__}]:=casePlus[{head,m,n},c];
casePlus[0,c___]:=casePlus[c];
casePlus[ERROR]:=ERROR;
casePlus[ERROR,__]:=ERROR;

casePlus[{Ricci,8},{Riemann,6},0,0,{Riemann,2},{Riemann,4},casePlus[casePlus[],{Ricci,3},0,{Riemann,1},{Weyl,5}]]

degrees[expr_,cd_,metric_]:=With[{
riemann=Riemann[cd],
ricci=Ricci[cd],
ricciscalar=RicciScalar[cd],
weyl=Weyl[cd],
tfricci=TFRicci[cd],
eps=epsilon[metric]},
Module[{deg},
deg[cd[_][expr1_]x_.]:=casePlus[addone[deg[expr1]],deg[x]];
deg[_riemann x_.]:=casePlus[{Riemann,0},deg[x]];
deg[_ricci x_.]:=casePlus[{Ricci,0},deg[x]];
deg[_ricciscalar x_.]:=casePlus[{RicciScalar,0},deg[x]];
deg[_weyl x_.]:=casePlus[{Weyl,0},deg[x]];
deg[_tfricci x_.]:=casePlus[{TFRicci,0},deg[x]];
deg[_eps x_.]:=casePlus[{epsilon,0},deg[x]];
deg[1|-1]:=0;
deg[x_]:=Infinity;
deg[expr]
]
]

The function is essentially based on the use sowing of the new list of indices by ToCanonical (this was its GivePerm option before xTensor 1.1.0). We need the function RiemannToPerm1 to avoid using PutScalar several times, which would be time-consuming. In this case we need to separate the action on the different metrics which have been defined and therefore this function does need a metric argument. It is a different question whether we make this feature public or not. October 2011: We have added a wrapping ContractMetric to fix cases like metric[a,b]Ricci[-a,-b] or metric[a,b]Riemann[-a,-c,-b,c].

(* First expand *)
RiemannToPerm[expr_,args___]:=RiemannToPerm0[ContractMetric@Expand[expr],args];
SetNumberOfArguments[RiemannToPerm,{1,2}];
(* Thread over sums *)
RiemannToPerm0[expr_Plus,args___]:=RiemannToPerm0[#,args]&/@expr;
RiemannToPerm0[expr_,args___]:=RiemannToPerm1[BreakScalars@PutScalar[expr],args];
(* Deal with the second argument *)
RiemannToPerm1[expr_]:=RiemannToPerm1[expr,$Metrics];
RiemannToPerm1[expr_,metrics_List]:=Fold[RiemannToPerm1,expr,metrics];
(* Possibilities for the first argument *)
RiemannToPerm1[Scalar[expr_],metric_Symbol]:=arrange[expr,degrees[expr,CovDOfMetric[metric],metric],metric];
RiemannToPerm1[expr_Times,metric_Symbol]:=RiemannToPerm1[#,metric]&/@expr;
RiemannToPerm1[tensor_Symbol[],metric_Symbol]:=RPerm[metric][{{0},0},xAct`xPerm`Cycles[{2,3}]]/;tensor===RicciScalar[CovDOfMetric[metric]];
RiemannToPerm1[Power[expr_,n_Integer],metric_Symbol]:=Power[RiemannToPerm1[expr,metric],n];
RiemannToPerm1[expr_,x_]:=expr;

Note that the RicciScalar is always converted into RInv[metric][{0},1] and never into WInv[metric][{0},1]. This additional step must be done by the user. The function arrange always receives a monomial.

(* 0. Trivial cases *)
arrange[0,_,_]:=0;
arrange[expr_,Infinity,_]:=expr;
(* 1. There is Ricci or RicciScalar: convert to Riemann *)
HoldPattern[arrange[expr_,casePlus[{Ricci|RicciScalar,m__},c___],metric_]]:=arrange[RicciToRiemann[expr,CovDOfMetric[metric]],casePlus[{Riemann,m},c],metric];
(* 2a. There are only Riemanns: compute non-dual permutation *)
HoldPattern[arrange[expr_,casePlus[{Riemann,Rcase__}],metric_]]:=With[{canonical={#[[1]],#[[2,1,1]]}&@Reap[ToCanonical[expr],"NewIndices"]},arrange1[canonical,{Sort[{Rcase}],0},metric]];
(* 2b. There are only Riemanns and an epsilon tensor: compute dual permutation. Oct 2011: epsilon must be sorted last *)
HoldPattern[arrange[expr_,casePlus[{Riemann,Rcase__},{epsilon,0}],metric_]]:=With[{epsilonname=epsilon[metric]},Module[{result},
TagSet[epsilonname,xSortPrecedence[epsilonname],Infinity];
With[{canonical={#[[1]],#[[2,1,1]]}&@Reap[ToCanonical[expr],"NewIndices"]},
result=arrange1[canonical,{Sort[{Rcase}],1},metric]
];
TagUnset[epsilonname,xSortPrecedence[epsilonname]];
result
]
];
(* 3. There are Weyl/TFRicci only: convert to Riemann. NOT UPDATED *)
arrange[expr_,{0,0,degWeyl_,degTFRic_,degeps_},metric_]:=With[{cd=CovDOfMetric[metric]},arrange[ToRiemann[FromWeyl[expr,cd],cd,metric],{degWeyl+degTFRic,0,0,0,degeps},metric]/.RPerm->WPerm];
(* 4. Final arrangements *)
arrange1[0,_,_]:=0;
arrange1[_Symbol[a_,b_,-a_,-b_],{{0},0},metric_]:=RPerm[metric][{{0},0},xAct`xPerm`Cycles[{2,3}]];
arrange1[{expr_,{sign_,perm_IndexList}},degs:{_,_},metric_]:=RPerm[metric][degs,sign TranslatePerm[PermutationFromTo[List@@perm,List@@IndexSort[perm]],xAct`xPerm`Cycles]];

The function FromWeyl changes temporarily the heads WeylCD and TRicciCD by RiemannCD and RicciCD respectively. This is a trick to avoid duplicated code.

FromWeyl[expr_,cd_]:=expr/.{Weyl[cd]:>Riemann[cd],TFRicci[cd]:>Ricci[cd]};

This function converts any Ricci tensor into its corresponding Riemann tensor (added $RicciSign in version 2.0.4, as pointed out by Teake):

RicciToRiemann[expr_,cd_]:=With[{ricciscalar=RicciScalar[cd],ricci=Ricci[cd],riemann=Riemann[cd],vbundle=First@VBundlesOfCovD[cd]},expr/.{
ricciscalar[]:>Module[{c=DummyIn@vbundle,d=DummyIn@vbundle},$RicciSign riemann[c,d,-c,-d]],ricci[a_,b_]:>Module[{c=DummyIn@vbundle},$RicciSign riemann[a,-c,b,c]]
}
];

PermToInv

Note also that these dispatched rules have no metric argument; it is added at the end. It is important to realize that right now this command will only work if the permutation is already in canonical form.

PermToInv[RPerm[metric_][{case_,0},perm_]]:=perm/.DispatchPermToRInvRules[case]/.RInv->RInv[metric];
PermToInv[RPerm[metric_][{case_,1},perm_]]:=perm/.DispatchPermToDualRInvRules[case]/.DualRInv->DualRInv[metric];
PermToInv[WPerm[metric_][{case_,0},perm_]]:=perm/.DispatchPermToRInvRules[case]/.RInv->WInv[metric];
PermToInv[WPerm[metric_][{case_,1},perm_]]:=perm/.DispatchPermToDualRInvRules[case]/.DualRInv->DualWInv[metric];
PermToInv[expr_]:=expr/.perm:(RPerm|WPerm)[_][__]:>PermToInv[perm];
SetNumberOfArguments[PermToInv,1];

RiemannToInv

Main transformation function. It is just a combination of RiemannToPerm and PermToInv:

RiemannToInv[expr_]:=RiemannToInv[expr,$Metrics];
RiemannToInv[expr_,metrics_List]:=Fold[RiemannToInv,expr,metrics];
RiemannToInv[expr_,metric_Symbol?MetricQ]:=PermToInv[RiemannToPerm[expr,metric]];
RiemannToInv[expr_,x_]:=(Message[RiemannToInv::unknown,"metric or list of metrics",x];expr);
SetNumberOfArguments[RiemannToInv,{1,2}];

RiemannSimplify

This function requires two additional arguments: the level of simplification (to be passed to InvSimplify) and the metric argument (to be passed to RiemannToInv):

?RiemannSimplify

RiemannSimplify[expr_]:=RiemannSimplify[expr,$InvSimplifyLevel,$CurvatureRelations,$Metrics];
RiemannSimplify[expr_,step_Integer]:=RiemannSimplify[expr,step,$CurvatureRelations,$Metrics];
RiemannSimplify[expr_,cr:(True|False)]:=RiemannSimplify[expr,$InvSimplifyLevel,cr,$Metrics];
RiemannSimplify[expr_,metric:(_Symbol?MetricQ|{__Symbol?MetricQ})]:=RiemannSimplify[expr,$InvSimplifyLevel,$CurvatureRelations,metric];
RiemannSimplify[expr_,step_Integer,cr:(True|False),metric_]:=InvToRiemann[InvSimplify[RiemannToInv[expr,metric],step],cr];
SetNumberOfArguments[RiemannSimplify,{1,3}];

RToW and WToR

Transformation between Riemann and Weyl invariants. Internal function:

RToWInv[RInv[metric_][deg_,index_],level_]:=setdim[RInv[deg,index]/.DispatchRToWRules[level,deg]/.WInv->WInv[metric],metric];
RToWInv[DualRInv[metric_][deg_,index_],level_]:=setdim[DualRInv[deg,index]/.DispatchRToWDualRules[level,deg]/.{DualWInv->DualWInv[metric],WInv->WInv[metric]},metric];
WToRInv[WInv[metric_][deg_,index_],level_]:=setdim[WInv[deg,index]/.DispatchWToRRules[level,deg]/.RInv->RInv[metric],metric];
WToRInv[DualWInv[metric_][deg_,index_],level_]:=setdim[DualWInv[deg,index]/.DispatchWToRDualRules[level,deg]/.{DualRInv->DualRInv[metric],RInv->RInv[metric]},metric];

Set correct dimension (replacing the constant symbol dim in the stored rules):

setdim[expr_,metric_]:=With[{dimension=DimOfVBundle@VBundleOfMetric[metric]},expr/.dim->dimension]

User driver:

?RToW

?WToR

RToW[expr_,level_:2]:=expr/.inv:(RInv|DualRInv)[_][_,_]:>RToWInv[inv,level];
WToR[expr_,level_:2]:=expr/.inv:(WInv|DualWInv)[_][_,_]:>WToRInv[inv,level];
SetNumberOfArguments[RToW,{1,2}];
SetNumberOfArguments[WToR,{1,2}];

RandomRiemannMonomial

RandomRiemannMonomial[case_]:=RandomRiemannMonomial[case,First[$Metrics]];
RandomRiemannMonomial[degree_Integer,metric_]:=RandomRiemannMonomial[Table[0,{degree}],metric];
RandomRiemannMonomial[case_List,metric_]:=PermToRiemann[RPerm[metric][{case,0},RandomPerm[4Length@case+Plus@@case,Images]]];
SetNumberOfArguments[RandomRiemannMonomial,{1,2}];

End package

End private part:

End[]

We have defined these new names, and they are now all protected:

Protect["xAct`Invar`*"];

The three user-global-variables and the constant-symbols dim, sigma must be unprotected:

Unprotect[xAct`Invar`dim,xAct`Invar`sigma,xAct`Invar`$InvSimplifyLevel,xAct`Invar`$CurvatureRelations,xAct`Invar`$ExpandedCommuteOrder12Q,xAct`Invar`InvarDirectory];

EndPackage[]

Examples

Setup

Off[General::spell,General::spell1]

DefManifold[M,4,IndexRange[a,q],z]

DefMetric[-1,metric[-a,-b],CD,{";","▽"},CurvatureRelations->False]

PrintAs[RiemannCD]^="R";
PrintAs[RicciCD]^="R";
PrintAs[RicciScalarCD]^="R";
PrintAs[epsilonmetric]^="ϵ";
PrintAs[metric]^="g";

Special cases

RiemannToInv[RicciScalarCD[]]

RiemannToInv[-3RicciScalarCD[]]

RiemannToInv[-3RiemannCD[-a,-b,a,b]]

RiemannToInv[-3CD[-a]@RicciScalarCD[]CD[a]@RicciScalarCD[]]

RiemannToInv[-3CD[-a]@RiemannCD[b,-c,c,-b]CD[a]@RicciScalarCD[]]

RiemannToInv[-3CD[-a]@RiemannCD[b,-c,c,-b]CD[a]@RicciCD[d,-d]]

RiemannToInv[-3CD[-a]@RiemannCD[b,-c,c,-b]CD[a]@RiemannCD[-d,-e,d,e]]

Check of RiemannToInv

checkperms[case_]:=Module[{riemanns,invs},
riemanns=xAct`Invar`Private`ReadInvarPerms[xAct`Invar`Private`$DataBaseDir<>"0/"<>xAct`Invar`Private`filename[0,case],"Mathematica"];
Print[AbsoluteTiming[Length[invs=RiemannToInv/@riemanns]]];
(invs/.RInv[metric]->RInv)===RInvs[1,case]]

First diagonal:

checkperms[{0}]

Second diagonal:

checkperms[{0,0}]

checkperms[{2}]

Third diagonal:

checkperms[{0,0,0}]

checkperms[{0,2}]

checkperms[{1,1}]

checkperms[{4}]

Fourth diagonal:

checkperms[{0,0,0,0}]

checkperms[{0,0,2}]

checkperms[{0,1,1}]

checkperms[{0,4}]

checkperms[{1,3}]

checkperms[{2,2}]

checkperms[{6}]

Fifth diagonal:

checkperms[{0,0,0,0,0}]

checkperms[{0,0,0,2}]

checkperms[{0,0,1,1}]

checkperms[{0,0,4}]

checkperms[{0,1,3}]

checkperms[{0,2,2}]

checkperms[{1,1,2}]

checkperms[{0,6}]

checkperms[{1,5}]

checkperms[{2,4}]

checkperms[{3,3}]

checkperms[{8}]

Sixth diagonal:

checkperms[{0,0,0,0,0,0}]

checkperms[{0,0,0,0,2}]

checkperms[{0,0,0,1,1}]

checkperms[{0,0,0,4}]

checkperms[{0,0,1,3}]

checkperms[{0,0,2,2}]

checkperms[{0,1,1,2}]

checkperms[{1,1,1,1}]

checkperms[{0,0,6}]

checkperms[{0,1,5}]

checkperms[{0,2,4}]

checkperms[{1,1,4}]

checkperms[{0,3,3}]

checkperms[{1,2,3}]

checkperms[{2,2,2}]

checkperms[{0,8}]

checkperms[{1,7}]

checkperms[{2,6}]

checkperms[{3,5}]

checkperms[{4,4}]

checkperms[{10}]

Seventh diagonal:

checkperms[{0,0,0,0,0,0,0}]

Check of InvToRiemann

$CurvatureRelations

checkperms[case_]:=Module[{riemanns,frominvs},
riemanns=Scalar/@xAct`Invar`Private`ReadInvarPerms[xAct`Invar`Private`$DataBaseDir<>"0/"<>xAct`Invar`Private`filename[0,case],"Mathematica"];
Print[AbsoluteTiming[Length[frominvs=InvToRiemann/@(RInvs[1,case]/.RInv->RInv[metric])]]];
riemanns===frominvs]

First diagonal:

checkperms[{0}]

Second diagonal:

checkperms[{0,0}]

checkperms[{2}]

Third diagonal:

checkperms[{0,0,0}]

checkperms[{0,2}]

checkperms[{1,1}]

checkperms[{4}]

Fourth diagonal:

checkperms[{0,0,0,0}]

checkperms[{0,0,2}]

checkperms[{0,1,1}]

checkperms[{0,4}]

checkperms[{1,3}]

checkperms[{2,2}]

checkperms[{6}]

Fifth diagonal:

checkperms[{0,0,0,0,0}]

checkperms[{0,0,0,2}]

checkperms[{0,0,1,1}]

checkperms[{0,0,4}]

checkperms[{0,1,3}]

checkperms[{0,2,2}]

checkperms[{1,1,2}]

checkperms[{0,6}]

checkperms[{1,5}]

checkperms[{2,4}]

checkperms[{3,3}]

checkperms[{8}]

Sixth diagonal:

checkperms[{0,0,0,0,0,0}]

checkperms[{0,0,0,0,2}]

checkperms[{0,0,0,1,1}]

checkperms[{0,0,0,4}]

checkperms[{0,0,1,3}]

checkperms[{0,0,2,2}]

checkperms[{0,1,1,2}]

checkperms[{1,1,1,1}]

checkperms[{0,0,6}]

checkperms[{0,1,5}]

checkperms[{0,2,4}]

checkperms[{1,1,4}]

checkperms[{0,3,3}]

checkperms[{1,2,3}]

checkperms[{2,2,2}]

checkperms[{0,8}]

checkperms[{1,7}]

checkperms[{2,6}]

checkperms[{3,5}]

checkperms[{4,4}]

checkperms[{10}]

Seventh diagonal:

checkperms[{0,0,0,0,0,0,0}]

The invariant in the paper

paperinv=RiemannCD[-a,-b,-c,-d]CD[-e]@RiemannCD[e,c,f,g]CD[a]@CD[-f]@CD[-h]@RiemannCD[b,d,h,-g]

canon=ToCanonical[paperinv]

RiemannToInv[%]

RiemannCD[a,b,c,d]CD[-e]@RiemannCD[-a,e,f,g](2CD[-d]@CD[-g]@CD[-h]@RiemannCD[-b,-c,-f,h]+CD[-b]@CD[-g]@CD[-h]@RiemannCD[-c,-d,-f,h])

Simplification[%]

RiemannToInv[%]

RiemannCD[a,b,c,d]CD[-e]@RiemannCD[-a,e,f,g](CD[-d]@CD[-g]@CD[-h]@RiemannCD[-b,-c,-f,h]-CD[-d]@CD[-g]@CD[-h]@RiemannCD[-b,-f,-c,h]+CD[-d]@CD[-g]@CD[-h]@RiemannCD[-b,h,-c,-f])

Simplification[%]

RiemannToInv[%]

RiemannCD[a,b,c,d](CD[-e]@RiemannCD[-a,e,f,g]CD[-d]@CD[-g]@CD[-h]@RiemannCD[-b,-c,-f,h]+CD[-f]@RiemannCD[-a,e,f,g]CD[-d]@CD[-e]@CD[-h]@RiemannCD[-b,-c,-g,h]-CD[-f]@RiemannCD[-a,e,f,g]CD[-d]@CD[-g]@CD[-h]@RiemannCD[-b,-c,-e,h])

Simplification[%]

RiemannToInv[%]

RiemannCD[a,b,c,d](CD[-e]@RiemannCD[-a,e,f,g]CD[-d]@CD[-g]@CD[-h]@RiemannCD[-b,-c,-f,h]+CD[g]@RiemannCD[-a,e,-e,f]CD[-d]@CD[-g]@CD[-h]@RiemannCD[-b,-c,-f,h]-CD[g]@RiemannCD[-a,e,-e,f]CD[-d]@CD[-f]@CD[-h]@RiemannCD[-b,-c,-g,h])

Simplification[%]

RiemannToInv[%]

InvSimplify[%,2]

canon

canon-CommuteCovDs[canon,CD,{-h,-g}]//RiemannToInv

InvSimplify[%,3]//Expand

InvSimplify[%,4]//Expand

canon-CommuteCovDs[canon,CD,{-g,-d}]//RiemannToInv

InvSimplify[%,3]//Expand

InvSimplify[%,4]//Expand

Nice examples

rules=RInvRules[4,{1,2,3}];

InvToRiemann/@(First/@Cases[rules,HoldPattern[_->0]]/.RInv->RInv[metric])

Block[{$CommuteCovDsOnScalars=True},ToCanonical[%]]

rules=RInvRules[3,{2,2,2}];

rulesb=Cases[rules,HoldPattern[RInv[{2,2,2},_]->_ RInv[{2,2,2},_]]];

Position[Denominator/@Last/@rulesb,8]

Position[rules,rulesb[[446]]]

Scalar[(Subscript[R,       ;a]^(abcd;e  )) (Subscript[R, cfgi;d;h]^(        )) (Subscript[R, be      ]^(  fg;h;i))]-(1/8) Scalar[(Subscript[R,       ;e]^(abcd;e  )) (Subscript[R, cdfg;i;h]^(        )) (Subscript[R, ab      ]^(  fg;h;i))]

%//RiemannToInv//InvSimplify//Simplify

Renato's examples

First example:

invar1=DualRInv[metric][{0,2},1]

invar2=InvSimplify[invar1]

InvToRiemann[invar1]

NoScalar[%]

CD[a]@CD[-a][%]

expr1=RiemannToInv[%]

CD[a]@CD[-a]@NoScalar@InvToRiemann[invar2]

expr2=RiemannToInv[%]

InvSimplify[expr1-expr2]//Simplify

Second example:

$CurvatureRelations=False;

R[a_,b_,c_,d_]:=RiemannCD[a,b,c,d];
R[a_,b_]:=xAct`Invar`Private`RicciToRiemann[RicciCD[a,b]];

CD[a_,b__][exp_]:=CD[a][CD[b][exp]]

IndicesOf[][expr]

expr=HoldComplete[-1/2*CD[a][R[g9,f5,f4,d2]]*CD[-a,-f9,c1][R[-g9,-f4,f9,b4]]*epsilonmetric[-f5,-d2,-b4,-c1]-1/2*R[g9,f5,f4,d2]*CD[a,-a,-f9,c1][R[-g9,-f4,f9,b4]]*epsilonmetric[-f5,-d2,-b4,-c1]-1/2*CD[a,h2][R[d3,f5,a0,z]]*CD[-a,-f9][R[-d3,f9,-h2,c1]]*epsilonmetric[-f5,-a0,-z,-c1]-1/2*CD[h2][R[d3,f5,a0,z]]*CD[a,-a,-f9][R[-d3,f9,-h2,c1]]*epsilonmetric[-f5,-a0,-z,-c1]-2*CD[a,b8][R[d3,f5,b7,d2]]*CD[-a,-b8][R[-d3,g6,-b7,c1]]*epsilonmetric[-f5,-d2,-g6,-c1]-2*CD[b8][R[d3,f5,b7,d2]]*CD[a,-a,-b8][R[-d3,g6,-b7,c1]]*epsilonmetric[-f5,-d2,-g6,-c1]-CD[a][R[g9,f5,h1,d2]]*CD[-a,-e6,e6][R[-g9,b0,-h1,e5]]*epsilonmetric[-f5,-d2,-b0,-e5]-R[g9,f5,h1,d2]*CD[a,-a,-e6,e6][R[-g9,b0,-h1,e5]]*epsilonmetric[-f5,-d2,-b0,-e5]-2*CD[a,c7][R[d3,c5,h9,f2]]*CD[-a,-h9][R[-d3,g6,-c5,c1]]*epsilonmetric[-f2,-c7,-g6,-c1]-2*CD[c7][R[d3,c5,h9,f2]]*CD[a,-a,-h9][R[-d3,g6,-c5,c1]]*epsilonmetric[-f2,-c7,-g6,-c1]+1/2*CD[a,b8][R[f3,f2]]*CD[-a,-b8][R[-f3,g7,b4,c1]]*epsilonmetric[-f2,-g7,-b4,-c1]+1/2*CD[b8][R[f3,f2]]*CD[a,-a,-b8][R[-f3,g7,b4,c1]]*epsilonmetric[-f2,-g7,-b4,-c1]-1/2*CD[a,b8][R[d3,f5,a0,z]]*CD[-a,-b8][R[-d3,c1]]*epsilonmetric[-f5,-a0,-z,-c1]-1/2*CD[b8][R[d3,f5,a0,z]]*CD[a,-a,-b8][R[-d3,c1]]*epsilonmetric[-f5,-a0,-z,-c1]+CD[a][R[g9,f5,i1,d2]]*CD[-a,-f9,-i1][R[-g9,b0,f9,e5]]*epsilonmetric[-f5,-d2,-b0,-e5]+R[g9,f5,i1,d2]*CD[a,-a,-f9,-i1][R[-g9,b0,f9,e5]]*epsilonmetric[-f5,-d2,-b0,-e5]-1/2*CD[a,-g9][R[g9,f3,h9,f2]]*CD[-a,-h9][R[-f3,g7,b4,c1]]*epsilonmetric[-f2,-g7,-b4,-c1]-1/2*CD[-g9][R[g9,f3,h9,f2]]*CD[a,-a,-h9][R[-f3,g7,b4,c1]]*epsilonmetric[-f2,-g7,-b4,-c1]];

Union@Flatten@Cases[expr,epsilonmetric[inds__]->{inds},{0,Infinity},Heads->True]

Union@Flatten@Cases[expr,R[inds__]->{inds},{0,Infinity},Heads->True]

Union@Flatten@Cases[expr,CD[inds__]->{inds},{0,Infinity},Heads->True]

Cases[Union[%,%%,%%%],_Symbol]

Drop[%,1]

Drop[%,-1]

AddIndices[TangentM,%]

IndicesOfVBundle[TangentM]

R[a_,b_,c_,d_]:=RiemannCD[a,b,c,d];
R[a_,b_]:=xAct`Invar`Private`RicciToRiemann[RicciCD[a,b],CD];

expr2=ReleaseHold[expr]

expr2//ToCanonical

(IndexList@@First@IndicesOfVBundle[TangentM])[[Range[15]]]

RiemannToInv[expr2]

InputForm[%]

Length[expr2]

InvSimplify[%%]

Relation in the previous paper

epsilonmetric[a,b,c,d]RiemannCD[-a,-b,e,f]RiemannCD[-c,-e,-f,g]RiemannCD[-d,h,i,j]RicciCD[-g,-i]RicciCD[-h,-j]+1/8epsilonmetric[a,b,c,d]RiemannCD[-a,-b,e,f]RiemannCD[-c,-d,-e,-f]RiemannCD[g,h,i,j]RicciCD[-g,-i]RicciCD[-h,-j]

RiemannToInv[%]

InvSimplify[%]

More examples

(memory=MemoryInUse[])/1000./1024

Examples of RiemannToInv and InvToRiemann:

test[1]=CD[-h]@CD[-a]@RiemannCD[a,b,c,d]RiemannCD[-b,-d,e,f]CD[-g]@CD[-c]@RiemannCD[-e,-f,g,h]+RInv[metric][{0,2,2},1157];
test[2]=RiemannCD[-a,e,f,g]RiemannCD[a,b,c,d]RiemannCD[h,i,j,k]CD[-d][RiemannCD[-h,-i,-j,-k]]CD[-e][RiemannCD[-b,-c,-f,-g]]-RInv[metric][{0,0,0,1,1},27015];
test[3]=RiemannCD[a,b,c,d]CD[-h]@CD[-g]@CD[-c]@CD[h]@CD[-f]@CD[-e]@CD[-d]@CD[g][RiemannCD[-a,e,-b,f]]+RInv[metric][{0,8},20018];
test[4]=epsilonmetric[-b,-d,-g,-h]CD[-c]@CD[-a]@CD[-f][RiemannCD[g,f,-e,h]]CD[e][RiemannCD[a,b,c,d]]-DualRInv[metric][{1,3},914];

RiemannToInv/@Array[test,4]

InvToRiemann/@Array[test,4]//NoScalar//ToCanonical

Shortcuts:

inv=RInv[metric];
dinv=DualRInv[metric];

Examples of relations at step 2:

test[5]=inv[{1,1,4},15073]+inv[{1,1,4},15063]/4-inv[{1,1,4},15072];
test[6]=dinv[{1,3},914]-dinv[{1,3},908]/4;
test[7]=inv[{0,0,0,1,1},27015]+2inv[{0,0,0,1,1},26365];
test[8]=inv[{8},1116]+inv[{8},1092]-inv[{8},1115];

$InvSimplifyLevel=2;
InvSimplify/@Array[test,4,5]

Examples of relations at step 3:

test[9]=inv[{1,1,4},15072]-inv[{1,1,4},1750]+inv[{1,1,4},1751]-2inv[{1,1,4},1753]-inv[{1,1,4},1754];
test[10]=dinv[{1,3},908]+4dinv[{1,3},328]-4dinv[{1,3},332];
test[11]=inv[{0,0,0,1,1},26365]+4inv[{0,0,0,1,1},26327]-4inv[{0,0,0,1,1},26335];
test[12]=inv[{8},1092]+inv[{8},1090]-inv[{8},1091];

$InvSimplifyLevel=3;
InvSimplify/@Array[test,4,9]

Examples of relations at step 4:

test[13]=inv[{0,0,4},70]-inv[{0,0,0,2},50]-inv[{0,0,1,1},100];
test[14]=-inv[{0,0,4},74]+2inv[{0,0,4},18]-2inv[{0,0,4},19]+4inv[{0,0,1,1},71]-4inv[{0,0,1,1},251]-4inv[{0,0,1,1},264]+4inv[{0,0,1,1},267]-2inv[{0,0,0,2},5]+4inv[{0,0,0,2},6]-2inv[{0,0,0,2},8]-2inv[{0,0,0,2},12]-4inv[{0,0,0,2},31]+8inv[{0,0,0,2},32]-4inv[{0,0,0,2},34]+8inv[{0,0,0,2},35]+8inv[{0,0,0,2},37]-4inv[{0,0,0,2},41]-4inv[{0,0,0,2},42]+4inv[{0,0,0,0,0},2]+10inv[{0,0,0,0,0},6]+4inv[{0,0,0,0,0},8]-4inv[{0,0,0,0,0},9]+8inv[{0,0,0,0,0},11]+6inv[{0,0,0,0,0},12]+4inv[{0,0,0,0,0},25]+4inv[{0,0,0,0,0},27]+8inv[{0,0,0,0,0},31]+8inv[{0,0,0,0,0},33]+8inv[{0,0,0,0,0},35]-8inv[{0,0,0,0,0},36]-8inv[{0,0,0,0,0},38]-8inv[{0,0,0,0,0},43]-8inv[{0,0,0,0,0},45];
test[15]=-dinv[{1,3},11]+dinv[{1,3},10]+dinv[{0,1,1},2]-dinv[{0,1,1},4]-2dinv[{0,1,1},6]+dinv[{0,1,1},88]-2dinv[{0,1,1},89];
test[16]=dinv[{1,3},5]-dinv[{0,1,1},1]-dinv[{0,1,1},18]/2;

$InvSimplifyLevel=4;
InvSimplify/@Array[test,4,13]

Examples of relations at step 5:

test[17]=-inv[{1,2,3},5861]+inv[{1,2,3},7]/8+inv[{1,2,3},49]/4+inv[{1,2,3},380]/2-inv[{1,2,3},392]/2-inv[{1,2,3},474]/2+inv[{1,2,3},480]/2-inv[{1,2,3},625]/2+inv[{1,2,3},631]/2+inv[{1,2,3},1356]-2inv[{1,2,3},1362]+inv[{1,2,3},1380]-inv[{1,2,3},1416]+inv[{1,2,3},1422]-inv[{1,2,3},2229]/2+inv[{1,2,3},2626]/2-2inv[{1,2,3},4985]-inv[{1,2,3},5804]+inv[{1,2,3},5816]+inv[{1,2,3},5837];
test[18]=dinv[{0,0,2},95]+dinv[{0,0},2]inv[{2},1]/4;
test[19]=-inv[{0,3,3},3901]+inv[{0},1]inv[{3,3},7]/8-inv[{0},1]inv[{3,3},150]/2+inv[{0},1]inv[{3,3},155]/2-inv[{0,3,3},14]/4-inv[{0,3,3},101]-inv[{0,3,3},111]-2inv[{0,3,3},585]+inv[{0,3,3},590]+2inv[{0,3,3},602]+inv[{0,3,3},718]+inv[{0,3,3},1144]+inv[{0,3,3},3567]-2inv[{0,3,3},3572]+2inv[{0,3,3},3875];
test[20]=dinv[{0,1,1},4]+dinv[{0,1,1},1]/2+dinv[{0,1,1},2]-dinv[{0,1,1},3];

$InvSimplifyLevel=5;
InvSimplify/@Array[test,4,17]

Examples of relations at step 6:

test[21]=-inv[{1,1,1,1},53]-dinv[{1,1},1]^2/4-inv[{1,1},1]^2/64+3inv[{1,1},1]inv[{1,1},4]/8-5inv[{1,1},4]^2/4-3inv[{1,1},1]inv[{1,1},5]/8+5inv[{1,1},4]inv[{1,1},5]/2-5inv[{1,1},5]^2/4-inv[{1,1,1,1},3]/4+inv[{1,1,1,1},5]/4+inv[{1,1,1,1},12]-inv[{1,1,1,1},14]-inv[{1,1,1,1},15]-inv[{1,1,1,1},16]+inv[{1,1,1,1},20]-2inv[{1,1,1,1},22]-inv[{1,1,1,1},23]-3inv[{1,1,1,1},38]+2inv[{1,1,1,1},39]+2inv[{1,1,1,1},40]+6inv[{1,1,1,1},41]-4inv[{1,1,1,1},42]-inv[{1,1,1,1},43]-inv[{1,1,1,1},44]+inv[{1,1,1,1},45]/2+3inv[{1,1,1,1},46]/2;
test[22]=inv[{0,0,1,1},1104]-dinv[{0,0},2]dinv[{1,1},25]/4-inv[{0,0,1,1},1098]/2;
test[23]=-inv[{0,0,1,3},25986]+dinv[{0,0},2]dinv[{1,3},294]/8+dinv[{0,0},2]dinv[{0,1,1},1]/8+dinv[{0,0},2]dinv[{0,1,1},2]/4+3dinv[{0,0},2]dinv[{0,1,1},18]/32+dinv[{0,0},2]dinv[{0,1,1},63]/8-dinv[{0,0},2]dinv[{0,1,1},86]/2+5dinv[{0,0},2]dinv[{0,1,1},87]/8+dinv[{0,0},2]dinv[{0,1,1},97]/16+dinv[{0,0},2]dinv[{0,1,1},117]/32+3dinv[{0,0},2]dinv[{0,1,1},271]/16-dinv[{0,0},2]dinv[{1,1},1]inv[{0},1]/16+inv[{0,0,1,3},25974]/2;
test[24]=-inv[{0,0,2,2},5470]-dinv[{0,0},2]dinv[{2,2},19]/8-dinv[{0,0},2]dinv[{0,0,2},2]/4+dinv[{0,0},2]dinv[{0,0,2},5]/16-dinv[{0,0},2]dinv[{0,0,2},14]/16+dinv[{0,0},2]dinv[{0,0,2},32]/4-dinv[{0,0},2]dinv[{0,2},15]inv[{0},1]/16+dinv[{0,0},2]dinv[{0,0,0},2]inv[{0},1]/16-dinv[{0,0},2]^2inv[{0},1]^2/128+dinv[{0,0},2]^2inv[{0,0},1]/32-inv[{0,0,2,2},5461]+inv[{0,0,2,2},5466];

$InvSimplifyLevel=6;
InvSimplify/@Array[test,4,21]

(MemoryInUse[]-memory)/1000./1024

Dual products

SignDetOfMetric[metric]^=sigma

Read dual invariants

ReadDuals[case_]:=(
Dlist[case]=DualRInvs[5,case]/.DualRInv->DualRInv[metric];
Slist[case]=InvToRiemann/@Dlist[case]);

ReadDuals[{0,0}]

ReadDuals[{0,2}]

ReadDuals[{0,0,2}]

ReadDuals[{0,1,1}]

ReadDuals[{1,1}]

ReadDuals[{0,0,0}]

ReadDuals[{0,0,0,0}]

ReadDuals[{0,0,0,0,0}]

ReadDuals[{0,4}]

ReadDuals[{1,3}]

ReadDuals[{2,2}]

Products

MakeDualProducts[case1_,case2_]:=Module[{
Dtmp=Outer[Times,Dlist[case1],Dlist[case2]],Stmp=Map[RiemannToInv,Outer[Times,Slist[case1],Slist[case2]]//NoScalar//ContractMetric,{2}]},
Union@Flatten[Map[InvSimplify[#,5]&,Stmp-Dtmp,{2}]//Expand]
];

duals[1]=MakeDualProducts[{0,0},{0,0}];

duals[2]=MakeDualProducts[{0,0},{0,0,0}];

duals[3]=MakeDualProducts[{0,0},{0,0,0,0}];

duals[4]=MakeDualProducts[{0,0},{0,0,0,0,0}];

duals[5]=MakeDualProducts[{0,0,0},{0,0,0}];

duals[6]=MakeDualProducts[{0,0,0},{0,0,0,0}];

duals[7]=MakeDualProducts[{0,0},{0,0,2}];

duals[8]=MakeDualProducts[{0,0},{0,1,1}];

duals[9]=MakeDualProducts[{0,0},{0,2}];

duals[10]=MakeDualProducts[{0,0},{1,1}];

duals[11]=MakeDualProducts[{0,0},{0,4}];

duals[12]=MakeDualProducts[{0,0},{1,3}];

duals[13]=MakeDualProducts[{0,0},{2,2}];

duals[14]=MakeDualProducts[{0,0,0},{0,2}];

duals[15]=MakeDualProducts[{0,0,0},{1,1}];

duals[16]=MakeDualProducts[{1,1},{1,1}];

duals[17]=MakeDualProducts[{0,2},{1,1}];

duals[18]=MakeDualProducts[{0,2},{0,2}];

Length[allduals=Flatten@Array[duals,18]]

ByteCount[allduals]

Two simple functions:

invsof[expr_]:=Union[Cases[expr,RInv[_][_,_]|RInv[_,_]|DualRInv[_][_,_]|DualRInv[_,_],{0,Infinity}]]

casesof[expr_]:=Union[First/@invsof[expr]]

We have 1044 invariants on the 51 equations:

Length[invsofallduals=invsof[allduals]]

and these 29 different cases:

listofcases=casesof[allduals]

Length[%]

Let us check whether there are still invariants which can be expanded:

checkinvs[case_]:=Complement[Cases[invsofallduals,RInv[metric][case,_]],RInvs[5,case]/.RInv->RInv[metric]]

checkinvs/@listofcases

Invariants to solve for:

order[case_]:=2Length[case]+Plus@@case;
degree[case_]:=Length[case];

largestRInv[expr_]:=Last@Last@Sort[{order[First[#]],-degree[First[#]],#}&/@Union@Cases[expr,RInv[_][__],{0,Infinity}]]

largests=largestRInv/@allduals

Sort@Split@Sort[largests]

solutions={};
Do[
eq=Expand[allduals[[count]]/.solutions];
sol=Expand/@Solve[eq==0,largestRInv[eq]][[1,1]];
solutions=Join[Map[Expand,solutions/.sol,{2}],{sol}],
{count,1,Length[allduals]}];

solutions[[-1]]

solutions>>"solutions"

Intersection[First/@solutions,invsof[Last/@solutions]]

cases=Union[First/@First/@solutions]

newrules=ReleaseHold/@(Hold[Cases[solutions,HoldPattern[RInv[_][#,_]->_]]]&)/@cases;

Directory[]

Do[
Export[xAct`Invar`Private`filename[6,cases[[count]],4],Sort[newrules[[count]]/.{RInv[metric]->RInv,DualRInv[metric]->DualRInv}],"Lines"],{count,1,Length[cases]}]

othercases=Complement[InvarCases[],cases]

Do[Export[xAct`Invar`Private`filename[6,othercases[[count]]],{},"Lines"],{count,1,Length[othercases]}]

Transpose[{cases,Length/@newrules}]//MatrixForm