IndexCoefficient looks formidable, so I ask of those who have
experience: is it possible to rewrite IndexSolve to use
IndexCoefficient instead of Coefficient, and then check that the
coefficient is indeed a Scalar?
Second question: IndexCoefficient does not seem to work for tensors
which have dummy indices, that is, internal contractions. Is this
true? Should it? It would be nice to be able to take IndexCoefficient
[scalar CD[-A][v[A]],CD[-B][v[B]]] and return scalar. If this could
work, then IndexCoefficient could be substituted for Coefficient in
IndexSolve.
Cheers
Leo
(1) I thought that this would be simple, using SymmetryOf, but this
includes the symmetries on contracted dummy indices. I notice that
there is a function SymmetryEquivalentsOf which has the machinery
needed to only find the symmetry group of a specific set of indices,
e.g. the free indices, but I'm not sure how to take the machinery and
turn it into a symmetry to use when creating the "shadow" tensor. This
method of course also has the disadvantage of having to create a
temporary tensor and then destroy it.
> Still working up the learning curve, so please excuse my ignorance.
This is not an easy question precisely...
> I'm trying to automatically generate rules from equations of motion.
> Naturally I'd like to use IndexSolve to solve for some term which will
> be replaced, and this returns a rule. However, IndexSolve can not
> solve for objects which have internally contracted indices, e.g. CD[-A]
> [v[A]].
Yes. IndexSolve is currently rather trivial, but I find difficult to
generalize it. Solving tensor equations is very hard in general
because you need to worry about many things at the same time:
symmetries, metrics, projections, invertibility, dummies, It is better
to solve case by case by hand.
> One solution to this shortcoming is to create a temporary "shadow"
> tensor whose index structure is identical to the free indices of the
> object you're trying to solve for. Ideally, if one was generalizing
> this, one would like to preserve the symmetries on those free indices,
> if there are any symmetries(1)
> (1) I thought that this would be simple, using SymmetryOf, but this
> includes the symmetries on contracted dummy indices.
Yes. SymmetryOf gives the symmetry of the tensor, independently of the
indices it has. The result is given in terms of slots numbers. Perhaps
it should be called SlotSymmetryOf. This function is designed to be
used from ToCanonical, which canonicalizes whole expressions, all at
once, and hence needs to know all symmetries. Then it reads the
indices and analyzes which ones are repeated or contracted.
> I notice that
> there is a function SymmetryEquivalentsOf which has the machinery
> needed to only find the symmetry group of a specific set of indices,
> e.g. the free indices, but I'm not sure how to take the machinery and
> turn it into a symmetry to use when creating the "shadow" tensor. This
> method of course also has the disadvantage of having to create a
> temporary tensor and then destroy it.
SymmetryEquivalentsOf would not do what you want. The concept of
symmetry of a contraction of tensors with symmetry is not easy, and is
not yet programmed in xAct. It requires some group-theory capabilities
which are not present in xPerm yet.
>Then one generates a replacement rule
> to replace the object with the "shadow", and then proceed how
> IndexSolve does.
I don't think this will work in general, because of what we just
discussed and because I'm not sure the concept itself makes sense. For
example, it is clear what should be returned from
IndexCoefficient[ 3 u[a] u[-a] + T[], u[c] u[-c] ]
But what should we return in the following expression ?
IndexCoefficient[ u[a] u[b], u[c] u[-c] ]
The fact that we can contract u[c] u[-c] indicates that there is a
metric, and perhaps we could decompose the first argument as a pure
trace and a traceless parts, and declare that the result is
proportional to the inverse metric. This is like asking what is the
coefficient of x^2+y^2 in x*y ? We can obviously define the
coefficient of a monomial in a polynomial, but the coefficient of a
polynomial in another polynomial ? We would need some sort of
arbitrary hierarchical structure among monomials, like in Groebner
Basis theory... I haven't thought much about it.
> While in the process of modifying IndexSolve to do the above, I
> learned that there is a function IndexCoefficient. It seems that the
> shortcoming of IndexSolve comes from using Mathematica's built-in
> Coefficient to extract the scalar prefactor, since it does not know
> that two expressions with different dummy indices (or in different
> positions) are equivalent.
>
> IndexCoefficient looks formidable, so I ask of those who have
> experience: is it possible to rewrite IndexSolve to use
> IndexCoefficient instead of Coefficient, and then check that the
> coefficient is indeed a Scalar?
IndexCoefficient was redesigned and improved to work with variational
derivatives, in which there are similar problems to those commented
above, but simpler. As you mention below, it does not work well when
there are dummy indices either, exactly for the same reasons. But this
is not a problem for variational derivation because there we always
differentiate with respect to a tensor with free indices.
> Second question: IndexCoefficient does not seem to work for tensors
> which have dummy indices, that is, internal contractions. Is this
> true? Should it? It would be nice to be able to take IndexCoefficient
> [scalar CD[-A][v[A]],CD[-B][v[B]]] and return scalar. If this could
> work, then IndexCoefficient could be substituted for Coefficient in
> IndexSolve.
It could be possible to make IndexCoefficient work for a restricted
family of cases, but I do not see how to make it work in general, or
how to determine which is the class of cases for which we want it to
work.
Cheers,
Jose.
IndexSolve::free = "Cannot solve for tensor with unknown free indices: `1`.";
IndexSolve[False, object_, options___] := {}
IndexSolve[True, object_, options___] := {{}}
IndexSolve[lhs_ == rhs_, object_, options___] :=
Module[{seq = Simplification[Expand[lhs - rhs]], frees, inds, dummies,
shadow, shadrule,(*syms,*)
scalar, righths,
verb = Verbose /. CheckOptions[options] /. Options[MakeRule]},
frees = FindFreeIndices[Evaluate[seq]];
inds = FindFreeIndices[object];
dummies = FindDummyIndices[object];
If[verb,
Print["frees: ", frees, ", inds: ", inds, ", dummies: ", dummies]];
If[Sort@frees =!= Sort@inds, Throw[Message[IndexSolve::free, inds]]];
(*
Need to do something with symmetries here!
*)
DefTensor[shadow @@ frees, ManifoldsOf[object](*,syms*)]; (* I am only dealing with Tensors here! This is wrong. *)
shadrule =
MakeRule[{object, shadow @@ frees}, MetricOn -> List @@ dummies,
PatternIndices -> List @@ dummies]; (* this should only match those cases with the same free indices as specified in object *)
If[verb, Print["shadow: ", shadow@@frees, ", shadrule: ", shadrule]];
scalar =
ToCanonical@
Scalar@Coefficient[Simplification[seq /. shadrule], shadow @@ frees];
If[verb, Print["scalar: ", scalar]];
If[! ScalarQ[scalar], Throw[Message[IndexSolve::error, ""]]];
righths = Simplification[NoScalar[seq - scalar object]];
If[verb, Print["righths: ", righths]];
Undef[shadow]; (* would prefer if DefTensor and Undef were quiet ... don't know how to do that *)
MakeRule[{Evaluate[object], Evaluate[-righths/scalar]},
PatternIndices -> List @@ inds, options] // Flatten
]
SetNumberOfArguments[IndexSolve, {2, Infinity}];
Protect[IndexSolve];
Thank you for your code. I like it, and works well. Just note that,
when constructing shadrule, you have to wrap the first argument of
MakeRule with Evaluate, because MakeRule has argument HoldFirst.
> I didn't realize that the code to find
> the subgroup of symmetries on a subset of indices was not present (I guess I
> didn't read the documentation carefully enough).
It is possible (and easy) to find the symmetries of a subset of *free*
indices, but not when there are contractions, because then we need to
worry about metrics and whether the metric is symmetric or
antisymmetric, and other things. I know the required algorithm, but it
is not needed for standard canonicalization, and hence never
programmed it.
> I see your point about the ambiguity of finding the coefficient of
> expressions like u[c] u[-c] in something which can be decomposed into trace
> and tracefree parts. It seems that solving for a tensor in an equation in
> the way that IndexSolve does is different though -- the coefficient of the
> tensor can only be a scalar.
True. I've just reread the code of IndexSolve, to compare with yours,
and there I say explicitly that we will only try to address situations
of the form IndexSolve[ scalar tensor1 == tensor2, tensor1 ], and
where tensor1 does not have dummies. Very limited. You want to
generalize this to have dummies in tensor1.
> This means that you can not solve for u[c] u[-c]
> in an expression like T[-a,-b] + c u[-a] u[-b] == 0, since any trace
> decomposition of u[-a] u[-b] has a metric tensor in front of the trace part
> (in reality one would be able to solve for u[c] u[-c] by tracing the whole
> equation, but that's not what we're asking). I'm probably not thinking
> broadly enough, but is there an example you have for solving an equation
> (where one looks for a scalar coefficient, rather than just finding the
> tensor coefficient) where you would find this ambiguity?
Suppose for example the equation
u[b] u[-b] u[a] == T[a]
If you use my IndexSolve or your IndexSolve, we get u[a] -> T[a] /
Scalar[ u[b] u[-b] ], which is rather silly. What we should get is u
[a] -> T[a] / Power[ Scalar[ T[b] T[-b] ], 1/3 ], but that is not
direct. Imagine having a more complicated function of u[b] u[-b] in
front of u[a]. So, in this sense, we also need to restrict that the
scalar and tensor2 in scalar tensor1 == tensor2 do not contain
tensor1.
> Below is my example of how IndexSolve works on a slightly larger class of
> cases (some cases with internal contractions). It does not take symmetries
> into account and therefore should not work in all cases.
Thanks again!
Cheers,
Jose.