The answer I'm posting (after two years) is based on an idea by @dominicprice (former student of Kasper, and developer of many aspects of cadabra
).
The code
def expand_nabla(ex, *, index_name=None, index_set=None):
if index_name is None:
index_name = "expandnablaidx"
Indices($expandnablaidx, expandnablaidx#$, $expandnablaindices$)
index_counter = 1
for nabla in ex[r'\nabla']:
nabla.name=r'\partial'
dindex = nabla.indices().__next__()
for arg in nabla.args():
ret:=0;
for index in arg.free_indices():
t2:= @(arg);
newidx = Ex(index_name + str(index_counter))
index_counter += 1
if index.parent_rel==sub:
t1:= -\Gamma^{@(newidx)}_{@(dindex) @(index)};
t2[index]:= _{@(newidx)};
else:
t1:= \Gamma^{@(index)}_{@(dindex) @(newidx)};
t2[index]:= ^{@(newidx)};
ret += Ex(str(nabla.multiplier)) * t1 * t2
nabla += ret
distribute(ex, repeat=True)
if index_set is not None:
rename_dummies(ex, "expandnablaindices", index_set)
return ex
How does it work?
The problem I found with the algorithm expand_nabla
reported in the cadabra book https://cadabra.science/notebooks/ref_programming.html, was that when considering higher order derivatives, the name of the dummy index wouldn't change, resulting in a "repeated indices" error.
It is desirable that the algorithm could use one of the declared indices!
Therefore, we have to give a name to the set of indices:
\nabla{#}::Derivative.
\partial{#}::PartialDerivative.
{a, b, c, d, e, f, g, h, i, j, k ,l}::Indices(space, position=fixed).
{\mu,\nu,\lambda,\rho}::Indices(spacetime, position=fixed).
NOTE: there are two types of indices: space
and spacetime
.
Now, we declare an expression and expand the nabla operator (covariant derivative).
foo := \nabla_{c}{ h^{a}_{b} };
expand_nabla(_, index_set="space");
Since we have asked to expand with dummy indices of type space
, the result is
$$\partial{c}{h^{a},{b}}+\Gamma^{a},{c d} h^{d},{b}-\Gamma^{d},{c b} h^{a},{d}$$
However, if we expand indices of type spacetime
, we get
$$\partial{c}{h^{a}\,{b}}+\Gamma^{a}\,{c \mu} h^{\mu}\,{b}-\Gamma^{\mu}\,{c b} h^{a}\,{\mu}$$
Higher orders
Higher orders (I was able to calculate up to third order) can be computed with ease, for example:
foo := \nabla_{e}{\nabla_{d}{\nabla_{c}{ h^{a}_{b} }}};
#expand_nabla(_, index_set="vector");
expand_nabla(_, index_name=r'\lambda');
The index_name=r'\lambda'
ask the algorithm to use the name \lambda
to expand the dummy indices.