# The scalar backend used for these component tests can be switched
# using either one of
#
# kernel(scalar_backend="sympy")

# Testing a substitution rule in which an abstract index in the 
# expression is matched against a coordinate in the pattern.
# FIXME: this is not supposed to match!
# {m,n,p}::Indices(values={t,r}).
# {r,t}::Coordinate.
# ex:= A_{m n} + B_{m n};
# rl:= A_{t r} + B_{t r} = 3;
# substitute(ex, rl)

def test01():
    __cdbkernel__=create_scope()
    print(__cdbkernel__.scalar_backend)
    {r,t}::Coordinate.
    {m,n}::Indices(values={t,r}, position=fixed).
    ex:= A_{m};
    rl:= A_{t} = 3;
    evaluate(ex, rl, simplify=True)
    tst:= \components_{m}( \comma{ \comma{t} = 3  } ) - @(ex);
    collect_components(_)
    assert(tst==0)
    print("Test 01 passed")

test01()

def test02():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate.
    {m,n,p,s}::Indices(values={t,r}, position=free).
    C::Depends(r,t).
    ex:= A_{n m} B_{m n p} ( C_{p s} + D_{s p} );
    rl:= [ A_{r t} = 3, B_{t r t} = 2, B_{t r r} = 5, C_{t r} = 1, D_{r t} = r**2*t, D_{t r}=t**2 ];
    evaluate(ex, rl)
    tst:= \components_{s}{\comma{r}= 6 r**2 t + 6, \comma{t}=15 t**2} - @(ex):
    collect_components(_)
    assert(tst==0)
    print("Test 02 passed")

test02()


# Riemann curvature example.

def test03():
    __cdbkernel__=create_scope()
    {r,t,\phi,\theta}::Coordinate.
    {\mu,\nu,\rho,\sigma,\lambda,\kappa,\gamma,\chi}::Indices(values={t,r,\phi,\theta}, position=fixed).
    \partial{#}::PartialDerivative.
    g_{\mu\nu}::Metric.
    g^{\mu\nu}::InverseMetric.         
    ss:= { g_{t t} = -(1-2 M/r),  g_{r r} = 1/(1-2 M/r),
           g_{\theta\theta} = r**2,  g_{\phi\phi}=r**2 \sin(\theta)**2 };
    complete(ss, $g^{\mu\nu}$)
    ch:= \Gamma^{\mu}_{\nu\rho} = 1/2 g^{\mu\sigma} ( \partial_{\rho}{g_{\nu\sigma}}
                                                +\partial_{\nu}{g_{\rho\sigma}}
                                                -\partial_{\sigma}{g_{\nu\rho}} );
    evaluate(ch, ss, rhsonly=True)
    chtst=rhs(ch)
        # Check against e.g. arxiv:0904.4184.
    tst:=\components^{\mu}_{\nu\rho}{ {t,r,t}          =  M/r**2 (1-2 M/r)**(-1),
                                {t,t,r}            =  M/r**2 (1-2 M/r)**(-1),
                                {r,r,r}            = -M/r**2 (1-2 M/r)**(-1),
                                {r,t,t}            =  M/r**2 (1-2 M/r),
                                {r,\theta,\theta}  = -r (1-2M/r),
                                {r,\phi,\phi}      = - (r-2M) \sin(\theta)**2,
                                {\theta,\theta,r}  = 1/r,
                                {\theta,r,\theta}  = 1/r,
                                {\phi,\phi,r}      = 1/r,
                                {\phi,r,\phi}      = 1/r,
                                {\theta,\phi,\phi} = -\sin(\theta)\cos(\theta),
                                {\phi,\theta,\phi} = \cos(\theta)/\sin(\theta),
                                {\phi,\phi,\theta} = \cos(\theta)/\sin(\theta)
                            } - @(chtst);
    collect_components(_)
    assert(tst==0)
    print("Test 03a passed")
    
    rm:= R^{\rho}_{\sigma\mu\nu} = \partial_{\mu}{\Gamma^{\rho}_{\nu\sigma}}
                              -\partial_{\nu}{\Gamma^{\rho}_{\mu\sigma}}
                              +\Gamma^{\rho}_{\mu\lambda} \Gamma^{\lambda}_{\nu\sigma}
                              -\Gamma^{\rho}_{\nu\lambda} \Gamma^{\lambda}_{\mu\sigma};

    substitute(rm, ch)
    evaluate(rm, ss, rhsonly=True)
    
    # The Ricci tensor should be zero.
    rc:= R_{\sigma\nu} = R^{\rho}_{\sigma\rho\nu};
    substitute(rc, rm)
    evaluate(rc, ss, rhsonly=True)
    assert(rhs(rc)==0)
    print("Test 03b passed")
    
    # Finally, the Kretschmann scalar.
    K:= K = R^{\mu}_{\nu\rho\sigma} R^{\lambda}_{\kappa\gamma\chi} g_{\mu\lambda} g^{\nu\kappa} g^{\rho\gamma} g^{\sigma\chi};
    substitute(K, rm)
    evaluate(K, ss)
    tmp=rhs(K)
    tst:= 48 M**2/r**6 - @(tmp);
    simplify(tst)
    assert(tst==0)
    print("Test 03 passed")
    
test03()
    
def test04():
    __cdbkernel__=create_scope()
    {r,t,\theta,\phi}::Coordinate.
    {m,n,p,q}::Indices(position=fixed,values={t,r,\theta,\phi}).
    rl:= { X^{t} -> X^{t}, g_{t t}=-1, g_{r r}=3 };
    ex:= g_{m n} X^{n};
    evaluate(ex, rl)
    tst:= \components_{m}{ \comma{ \comma{t} = -X^{t} } } - @(ex);
    collect_components(_)
    assert(tst==0)
    print("Test 04 passed")

test04()

# Component computations leading to a scalar.

def test05():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate.
    {m,n,p}::Indices(values={r,t}).
    ex:= A_{m n}*A_{m n};
    ss:= { A_{r r} -> a+b, A_{r t} -> b, A_{t t} -> c, A_{r r t}->3, A_{t t t}->c };
    evaluate(ex, ss)
    tst:= (a+b)**2+b**2+c**2  - @(ex);
    assert(tst==0)
    print("Test 05 passed")

test05()

# The following test components of tensor outer products. The first one
# leads to a scalar multiplying the components tensor.

def test06():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate.
    {m,n}::Indices(values={t,r}).
    ex:= g_{m n} X^{n} A;
    rl:= { X^{t} -> X^{t}, g_{t t}=1, g_{r r}=3 };
    evaluate(ex,rl)
    tst:= \components_{m}{ \comma{ \comma{t}= A X^{t} } } - @(ex);
    collect_components(_)
    assert(tst==0)
    print("Test 06 passed")

test06()

 
# The second one is manifestly an outer product.

def test07():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate.
    {m,n,p}::Indices(values={t,r}).
    ex:= X^{n} X^{m} X^{p};
    rl:= { X^{t} -> T, X^{r} -> R };
    evaluate(ex,rl)
    tst:=\components^{n m p}{ {t,t,t}=T**3, {t,t,r}=R T**2, {t,r,t}=R T**2,
                            {t,r,r}=T R**2, {r,t,t}=R T**2, {r,t,r}=T R**2,
                            {r,r,t}=R**2 T, {r,r,r}=R**3 } - @(ex);
    collect_components(tst)
    assert(tst==0)
    print("Test 07 passed")

test07()                          

# More GR stuff.
# {\mu,\nu,\rho,\sigma}::Indices(position=fixed).
# \partial{#}::PartialDerivative.
# ex:= R^{\mu}_{\nu\rho\sigma};
# rl:= { R^{\mu}_{\nu\rho\sigma} -> \partial_{\rho}{\Gamma^{\mu}_{\nu\sigma}},
#         \Gamma^{\mu}_{\nu\rho} = 1/2 g^{\mu\sigma} ( \partial_{\rho}{g_{\nu\sigma}}
#                                                  +\partial_{\nu}{g_{\rho\sigma}}
#                                                  -\partial_{\sigma}{g_{\nu\rho}} )
#                                                   };
# substitute(ex, rl, repeat=True);                                                  


# -------------------------------------------

# This one does not yet work, and might be better served with a split-index
# first.
# {r,t,\phi,\theta}::Coordinate.
# {m,n}::Indices("full",values={t,r,\phi,\theta}).
# {i,j}::Indices("3d", parent="full").
# ex:= A_{m n}+A_{n m};
# rl:= { A_{t t} -> 3, A_{i j} -> \delta_{i j} };
# evaluate(ex, rl);

# Self-contracted tensors.

def test08():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate.
    {m,n,p}::Indices(values={r,t})
    ex:= A_{n}^{n} A_{m}^{m};
    rl:= { A_{t}^{r} -> S, A_{t}^{t} -> 3, A_{r}^{r} -> r };
    evaluate(ex, rl)
    tst:= (r+3)**2 - @(ex);
    simplify(tst)
    assert(tst==0)
    print("Test 08 passed")

test08()

def test09():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate.
    {m,n,p}::Indices(values={r,t})
    ex:= A_{m}^{m};
    rl:= { A_{t}^{r} -> S, A_{t}^{t} -> 3, A_{r}^{r} -> r };
    evaluate(ex, rl)
    tst:= r+3 - @(ex);
    assert(tst==0)
    print("Test 09 passed")

test09()

def test10():
    __cdbkernel__=create_scope()
    {r,t,\phi}::Coordinate.
    {m,n,p}::Indices(values={r,t,\phi}).
    ex:= A_{n}^{n p} B_{p};
    rl:= A_{r}^{r t} -> a, A_{r}^{r r} -> b, A_{r}^{t t} -> Q1, A_{t}^{r t} -> Q2, A_{t}^{r r} -> Q3,
         A_{t}^{t t}-> c, A_{t}^{t r} -> d, B_{r} -> Br, B_{t} -> Bt;
    evaluate(ex, rl)
    tst:= a Bt + c Bt + b Br + d Br - @(ex);
    assert(sympy.simplify(tst)==0)
    print("Test 10 passed")

test10()

# This does not work when 'Q' is removed.
def test11():
    __cdbkernel__=create_scope()
    {r,t,\phi}::Coordinate;
    {m,n,p}::Indices(values={r,t,\phi});
    ex:= T_{m n} = Q A_{m n};
    rl:= A_{r r} -> r, A_{t t} -> t, A_{r t} -> t r;
    evaluate(ex, rl, rhsonly=True)
    ex2:= T_{m m};
    substitute(ex2, ex)
    evaluate(ex2, rl)
    tst:= Q (r+t) - @(ex2);
    assert(tst==0)
    print("Test 11 passed")

test11()

def test12():
    __cdbkernel__=create_scope()
    {\mu,\nu,\rho,\sigma}::Indices(values={0,1,2,3});
    ex:=A_{\mu} A_{\mu};
    rl:= A_{0} -> -\xi_{0}, A_{1} -> \xi_{1}, A_{2} -> -\xi_{2}, A_{3} -> \xi_{3};
    evaluate(ex, rl)
    tst:= \xi_{0}\xi_{0} + \xi_{1}\xi_{1} + \xi_{2}\xi_{2} + \xi_{3}\xi_{3} - @(ex);
    assert(tst==0)
    print("Test 12 passed")

test12()

def test13():
    __cdbkernel__=create_scope()
    {i,j,k}::Indices(values={1,2,3});
    \epsilon^{i j k}::EpsilonTensor;
    ex:=\epsilon^{i j k};
    evaluate(ex)
    tst:= \components^{i j k}{ {1,2,3}=1, {1,3,2}=-1, {2,1,3}=-1, {2,3,1}=1, {3,1,2}=1, {3,2,1}=-1 } - @(ex);    
    collect_components(tst)
    assert(tst==0)
    print("Test 13 passed")

test13()
        
def test14():        
    __cdbkernel__=create_scope()
    {t,x,y,z}::Coordinate;        
    {\mu,\nu,\rho,\sigma}::Indices(values={t,x,y,z});
    f::Depends({t,x,y,z});
    \partial{#}::PartialDerivative;
    ex:= F^{\mu} \partial_{\mu}{ f };
    rl:= F^{t} = A, F^{x} = B, F^{z} = C;
    evaluate(ex, rl)
    canonicalise(ex)
    tst:= A \partial_{t}{f} + B \partial_{x}{f} + C \partial_{z}{f} - @(ex);
    assert(tst==0)
    print("Test 14 passed")

test14()

def test15():        
    __cdbkernel__=create_scope()
    {t,x,y,z}::Coordinate;        
    {\mu,\nu,\rho,\sigma}::Indices(values={t,x,y,z});
    f::Depends({t,x,y,z});
    \partial{#}::PartialDerivative;
    ex:= F^{\mu\nu} \partial_{\mu \nu}{ f };
    rl:= F^{t x} = A, F^{x y} = B, F^{z x} = C;
    evaluate(ex, rl)
    canonicalise(ex)
    tst:= A \partial_{t x}{f} + B \partial_{x y}{f} + C \partial_{x z}{f} - @(ex);
    assert(tst==0)
    print("Test 15 passed")

test15()

def test16():        
    __cdbkernel__=create_scope()
    {x,y,z}::Coordinate;
    f::Depends({x,y,z});
    \partial{#}::PartialDerivative;
    {i,j,k,l}::Indices(values={x,y,z}, position=free);
    ex:= 2 \partial_{i i}{3 f} + m**2 f;
    evaluate(ex)
    tst:= 6\partial_{x x}{f} + 6\partial_{y y}{f} + 6\partial_{z z}{f} + m**2 f - @(ex);
    assert(tst==0)
    print("Test 16 passed")

test16()
        
def test17():
    __cdbkernel__=create_scope()
    {x,y,z}::Coordinate;
    {i,j,k,l}::Indices(values={x,y,z}, position=free);
    Ax::Depends(y);
    Ay::Depends(z);
    Az::Depends(x);        
    \partial{#}::PartialDerivative;
    ex:= \partial_{i j}{A_{i}} + m**2 A_{j};
    evaluate(ex, $A_x=Ax, A_y=Ay, A_z=Az$)
    tst:= \components_{j}{ \comma{x}=m**2 Ax, \comma{y}=m**2 Ay, \comma{z}=m**2 Az } - @(ex);
    collect_components(tst)
    assert(tst==0)
    print("Test 17 passed")

test17()

def test18():
    __cdbkernel__=create_scope()
    {x,y,z}::Coordinate;
    {i,j,k,l}::Indices(values={x,y,z}, position=free);
    A_{i}::Depends({x,y,z});
    \partial{#}::PartialDerivative;
    ex:= \partial_{i}{A_{j}};
    evaluate(ex, $A_{z}=z$)
    tst:= \components_{j i}( \comma{ {z,z}=1 } ) - @(ex);
    collect_components(tst)
    assert(tst==0)
    print("Test 18 passed")

test18()

def test19():
    __cdbkernel__=create_scope()
    {x,y,z}::Coordinate;
    {i,j,k,l}::Indices(values={x,y,z}, position=free);
    A_{i}::Depends({x,y,z});
    \partial{#}::PartialDerivative;
    ex:= \partial_{i}{A_{i}} + m**2 f;
    evaluate(ex, $A_{z}=\sin{z}$)
    tst:= m**2 f + \cos{z} - @(ex);
    sort_product(_)
    assert(tst==0)
    print("Test 19 passed")

test19()

def test20():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate;
    {a,b,c}::Indices(values={r,t});
    \partial{#}::PartialDerivative;
    ex:= A_{a} B_{b} C_{a b};
    rl:= [ A_{r} = 1, A_{t} = q, B_{t} = 2, B_{r} = w, C_{r t} = 3, C_{t r}=5 ];
    evaluate(ex, rl)
    tst:= 6 + 5 q w - @(ex);
    assert(tst==0)
    print("Test 20 passed")

test20()

def test21():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate;
    {a,b,c}::Indices(values={r,t});
    \partial{#}::PartialDerivative;
    {f,g}::Depends(r);    
    ex:= A_{a} \partial_{a}{ f / g };
    rl:= { A_{r} = q };
    evaluate(ex, rl)
    tst:= q (- f \partial_{r}{g} + g \partial_{r}{f} )/g**2 - @(ex);
    assert(tst==0)
    print("Test 21 passed")

test21()

def test22():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate;
    {a,b,c}::Indices(values={r,t});
    \partial{#}::PartialDerivative;
    {f,g}::Depends(r);    
    ex:= A_{a} \partial_{a}{ f * g };
    rl:= { A_{r} = q };
    evaluate(ex, rl)
    tst:= q ( f \partial_{r}{g} + g\partial_{r}{f} ) - @(ex);
    assert(tst==0)
    print("Test 22 passed")

test22()

def test23():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate;        
    {a,b,c}::Indices(values={r,t}, position=fixed);
    \partial{#}::PartialDerivative;
    ex:= \partial^{a}{ A_{a} };
    rl:= A_{r} = r;
    try:
        evaluate(ex, rl);
        assert(1==0)
    except:
        print("Test 23a passed")
    ex:= \partial_{a}{ A^{a} };
    rl:= A^{r} = r;
    evaluate(ex, rl)
    tst:= 1 - @(ex);
    assert(tst==0)
    print("Test 23b passed")
        
test23()
        

# 
# {r,t}::Coordinate.    
# {a,b,c}::Indices(values={r,t});
# \partial{#}::PartialDerivative;
# A_{a}::Depends(r);
# ex:= \partial_{a}{ A_{b} };
# rl:= \partial_{r}{ A_{t} } = 1;
# evaluate(ex, rl);

def test24():
    __cdbkernel__=create_scope()
    {\theta, \varphi}::Coordinate;
    {\alpha, \beta, \gamma, \delta, \rho, \sigma, \mu, \nu, \lambda}::Indices(values={\varphi, \theta}, position=fixed);
    \partial{#}::PartialDerivative;
    g_{\alpha\beta}::Metric.
    g^{\alpha\beta}::InverseMetric.

    sphe:={ g_{\theta\theta} = r**2,
        g_{\varphi\varphi} = r**2 \sin(\theta)**2 }.
    complete(sphe, $g^{\alpha\beta}$)

    ch:= \Gamma^{\alpha}_{\mu\nu} = 1/2 g^{\alpha\beta} ( \partial_{\nu}{g_{\beta\mu}}
                                                        +\partial_{\mu}{g_{\beta\nu}}
                                                        -\partial_{\beta}{g_{\mu\nu}} );

    evaluate(ch, sphe, rhsonly=True)

    rm:= R^{\rho}_{\sigma\mu\nu} = +\partial_{\mu}{\Gamma^{\rho}_{\sigma\nu}}
           -\partial_{\nu}{\Gamma^{\rho}_{\sigma\mu}}
           +\Gamma^{\rho}_{\beta\mu} \Gamma^{\beta}_{\sigma\nu}
           -\Gamma^{\rho}_{\beta\nu} \Gamma^{\beta}_{\sigma\mu};

    substitute(rm, ch)
    evaluate(rm, sphe, rhsonly=True)

    ricci:= R_{\sigma\nu} = R^{\rho}_{\sigma\rho\nu};
    substitute(ricci, rm)
    evaluate(ricci, sphe, rhsonly=True)

    R:= R = R_{\sigma\nu} g^{\sigma\nu};
    substitute(R, ricci)
    evaluate(R, sphe, rhsonly=True)
    val = rhs(R)
    tst:= @(val) - 2/r**2;
    assert(tst==0)
    print("Test 24 passed")
        
test24()

# Evaluating with completely unspecified tensors.
#
def test25():        
    __cdbkernel__=create_scope()
    {a,b,c}::Indices(values={0,1,2,3});
    ex:= X^{a} Y_{a};
    rl:= { Y_{0} = 3, Y_{2} = -5 };
    evaluate(ex, rl)
    tst:= 3 X^{0} - 5 X^{2} - @(ex);
    assert(tst==0)
    print("Test 25 passed")

test25()
        
def test26():
    __cdbkernel__=create_scope()
    {a,b,c}::Indices(values={0,1,2});
    ex:= X^{a b} Y_{a b};
    evaluate(ex)
    tst:= X^{0 0} Y_{0 0} + X^{0 1} Y_{0 1} + X^{0 2} Y_{0 2}
        + X^{1 0} Y_{1 0} + X^{1 1} Y_{1 1} + X^{1 2} Y_{1 2}
        + X^{2 0} Y_{2 0} + X^{2 1} Y_{2 1} + X^{2 2} Y_{2 2} - @(ex):
    sort_product(tst)
    assert(tst==0)
    print("Test 26 passed")

test26()

def test27():
    __cdbkernel__=create_scope()        
    {r,t}::Coordinate;        
    {m,n}::Indices(position=fixed, values={r,t});
    {\sigma,\tau}::Coordinate;        
    X^{m}::Depends(\sigma,\tau);
    \partial{#}::PartialDerivative;
    ex:=\partial_{\tau}{X^{m}};
    evaluate(ex)
    tst:=\components^{m}{\comma{ \comma{t} = \partial_{\tau}{X^{t}} }{ \comma{r} = \partial_{\tau}{X^{r}} } } - @(ex);
    collect_components(_)
    assert(tst==0)
    print("Test 27 passed")

test27()
    
        
import cdb.relativity as rel
# the line below does not work on Travis, so we state the metric here explicitly.
#import cdb.relativity.schwarzschild as ss

def test28():
    g:={g_{t t} = -1+2M/r, g_{r r} = 1/(1-2M/r), g_{\theta\theta} = r**2, g_{\phi\phi} = r**2*sin(\theta)**2};
    complete(g,$g^{\mu\nu}$);
    ch=rel.christoffel_from_metric();
    evaluate(ch, g, rhsonly=True);
    \ddot{#}::Accent;
    \dot{#}::Accent;
    geodesic:= \ddot{X^{\mu}} + \Gamma^{\mu}_{\nu\rho} \dot{X^{\nu}} \dot{X^{\rho}} = 0;
    # Doing evaluate here messes up the components node of the \Gamma factor,
    # there's a weird \components inside a \components. Can probably be tested
    # with a simpler example.
    substitute(geodesic, ch);
    evaluate(geodesic);

#test27()
        
    
def test28():
    __cdbkernel__=create_scope()
    {t,x}::Coordinate:
    {m,n}::Indices(position=fixed, values={t,x});
    a::Depends(t,x);
    \partial{#}::PartialDerivative;
    ex:= (\partial_{t}{a})**2/a**3+a;
    evaluate(ex)
    tst:= (\partial_{t}{a})**2/a**3+a - @(ex);
    assert(tst==0)
    print("Test 28 passed")

test28()


def test29():
    __cdbkernel__=create_scope()
    {t,r}::Coordinate;
    {\mu,\nu,\rho}::Indices(values={r,t});
    ex:= T_{\mu\nu} = T_{\mu\nu};
    rl:= { T_{r t} = Q, T_{t t} = N M };
    evaluate(ex, rl, rhsonly=True)
    # Double evaluation should not introduce additional \components wrappers
    # around scalars.
    evaluate(ex, rl, rhsonly=True)
    ex2=ex[1]
    tst:= \components_{\mu\nu}{\comma{ \comma{r}{t} = Q }{ \comma{t}{t} = N M } } - @(ex2);
    collect_components(tst)        
    assert(tst==0)
    print("Test 29 passed")

test29()

def test30():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate;
    A::Depends(t);        
    \partial{#}::PartialDerivative;
    ex:=\partial_{t}{A};
    evaluate(ex)
    tst:=\partial_{t}{A} - @(ex);
    assert(tst==0)
    print("Test 30a passed")
    ex2:= \partial_{t}{ A**2 };
    evaluate(ex2)
    tst2:= 2A\partial_{t}{A} - @(ex2);
    assert(tst2==0)
    print("Test 30b passed")

test30()

def test31():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate;
    A::Depends(t);
    \partial{#}::PartialDerivative;
    {\mu,\nu,\rho}::Indices(values={r,t});
    ex:= T_{\mu\nu} = T_{\mu\nu};
    rl:= { T_{r t} = Q, T_{t t} = M };
    glrl:= { g_{r r} = N, g_{t t} = B };
    Rrl:= R = \partial_{t}{A};
    evaluate(ex, rl, rhsonly=True)
    ex2:= T_{\mu\nu} - g_{\mu\nu} R;
    substitute(ex2, ex)
    substitute(ex2, Rrl)
    evaluate(ex2, glrl)
    tst:= \components_{\mu\nu}{\comma{ \comma{r}{t} = Q }{ \comma{t}{t} = -B \partial_{t}{A} + M }{ \comma{r}{r} = - N \partial_{t}{A} }} - @(ex2);
    assert(tst!=0)
    collect_components(tst)
    assert(tst==0)
    print("Test 31 passed")

test31()

        

# def test20():
# __cdbkernel__=create_scope()
# {a,b,c,d}::Indices(values={0,1,2}, position=fixed);
# {a,b,c,d}::Integer(0..2);
# \Gamma{#}::GammaMatrix(metric=\eta);
# C_{a b c}::AntiSymmetric;
# ex:= C_{a b c} \Gamma^{a b c};
# rl:= [ C_{0 1 2} = \alpha ];
# evaluate(ex, rl);
# tst:= 6 \Gamma^{0 1 2} - @(ex);
# assert(tst)

# def test32():
#     # https://github.com/kpeeters/cadabra2/issues/86    
#     __cdbkernel__=create_scope()
#     {r,t}::Coordinate.
#     {m,n}::Indices(values={t,r}, position=fixed).
#     ex:= A_{m} = A_{m};
#     rl:= A_{t} = 3;
#     evaluate(ex, rl, rhsonly=True, simplify=True)
#     ex2:= A_{m};
#     evaluate(ex2, ex);
#     tst:= \components_{m}( \comma{ \comma{t} = 3  } ) - @(ex2);
#     collect_components(_)
#     assert(tst==0)
#     print("Test 32 passed")

def test32():
    __cdbkernel__=create_scope()
#kernel(scalar_backend="mathematica")
    {x,y,\xi}::Coordinate;
    {i,j,k}::Indices(values={x,y,\xi});
    \partial{#}::PartialDerivative;
    F_{i j}::Depends(x,\xi);
    ex:=\partial_{i}{F_{i k}} + m**2 A_{k};
    evaluate(ex)
    tst:=\components_{k}{ \comma{ \comma{\xi} = m**2 A_{\xi} + \partial_{x}{F_{x \xi}} + \partial_{\xi}{F_{\xi \xi}} }{ \comma{y}=m**2 A_{y} + \partial_{x}{F_{x y}} + \partial_{\xi}{F_{\xi y}}}{ \comma{x} = m**2 A_{x} + \partial_{x}{F_{x x}} + \partial_{\xi}{F_{\xi x}} } - @(ex);
    collect_components(tst)
    assert(tst==0)
    print("Test 31 passed")

test32()

def test33():
    __cdbkernel__=create_scope()
    {\sigma,\tau}::Coordinate;                          
    {a,b,c}::Indices(worldsheet, values={\sigma, \tau});
    #{\mu,\nu,\rho}::Indices(targetspace, values={0,1,2,3} );
    #{\mu,\nu,\rho}::Indices(targetspace);
    \partial{#}::PartialDerivative;
    X^{a?}::Depends(\sigma, \tau);
    X^{\mu}::Depends(\sigma, \tau);    
    ex:= \eta^{a b} \partial_{a}{ X^{\mu} } \partial_{b}{ X^{\nu} } g_{\mu\nu};
    mink:= { \eta^{\tau\tau}=-1, \eta^{\sigma\sigma}=1 };
    evaluate(ex, mink)
    tst:= - g_{\mu\nu}\partial_{\tau}{X^{\mu}}\partial_{\tau}{X^{\nu}} + g_{\mu\nu}\partial_{\sigma}{X^{\mu}}\partial_{\sigma}{X^{\nu}} - @(ex);
    assert(tst==0)
    print("Test 33 passed")

test33()
        
def test34():
    __cdbkernel__=create_scope()
    {\sigma,\tau}::Coordinate;                          
    {a,b,c}::Indices(worldsheet, values={\sigma, \tau});
    \partial{#}::PartialDerivative;
    X^{a?}::Depends(\sigma, \tau);
    X^{\mu}::Depends(\sigma, \tau);    
    ex:= \eta^{a b} \partial_{a}{ X^{\mu} } \partial_{b}{ X^{\nu} };
    mink:= { \eta^{\tau\tau}=-1, \eta^{\sigma\sigma}=1 };
    evaluate(ex, mink)
    tst:= - \partial_{\tau}{X^{\nu}}\partial_{\tau}{X^{\mu}} + \partial_{\sigma}{X^{\nu}}\partial_{\sigma}{X^{\mu}} - @(ex);
    assert(tst==0)
    print("Test 34 passed")

test34()

def test35():
    __cdbkernel__=create_scope()
    {r,t,\phi,\theta}::Coordinate.
    {\mu,\nu,\rho,\sigma,\lambda,\kappa,\gamma,\chi}::Indices(values={t,r,\phi,\theta}, position=fixed).
    \partial{#}::PartialDerivative.
    g_{\mu\nu}::Metric.
    g^{\mu\nu}::InverseMetric.         
    ss:= { g_{t t} = -(1-2 M/r),  g_{r r} = 1/(1-2 M/r),
           g_{\theta\theta} = r**2,  g_{\phi\phi}=r**2 \sin(\theta)**2 };
    complete(ss, $g^{\mu\nu}$)
    ch:= \Gamma^{\mu}_{\nu\rho} = 1/2 g^{\mu\sigma} ( \partial_{\rho}{g_{\nu\sigma}}
                                                +\partial_{\nu}{g_{\rho\sigma}}
                                                -\partial_{\sigma}{g_{\nu\rho}} );
    evaluate(ch, ss, rhsonly=True)
    comp=substitute($\Gamma^{r}_{\phi\phi}$, ch)
    tst:= (2M-r) \sin(\theta)**2 - @(comp);
    assert(tst==0)
    print("Test 35 passed")

test35()

def test36():
    __cdbkernel__=create_scope()
    {r,t}::Coordinate;
    {i,j}::Indices(values={r,t});
    F::Depends(r);
    \partial{#}::PartialDerivative;
    ex:=\partial_{i}{\exp(F)};
    evaluate(ex)
    tst:=\components_{i}( \comma{ \comma{r} = \exp(F) \partial_{r}{F}  } ) - @(ex);
    collect_components(tst)
    assert(tst==0)
    print("Test 36 passed")

test36()

def test37():
    __cdbkernel__=create_scope()
    {r,t,\phi,\theta}::Coordinate.
    {\mu,\nu,\rho,\sigma,\lambda,\kappa,\gamma,\chi}::Indices(values={t,r,\phi,\theta}, position=fixed).
    \partial{#}::PartialDerivative.
    g_{\mu\nu}::Metric.
    g^{\mu\nu}::InverseMetric.
    g::Determinant(g_{\mu\nu}).
    gt::Trace(g_{\mu\nu}).
    ss:= { g_{t t} = -(1-2 M/r),  g_{r r} = 1/(1-2 M/r),
           g_{\theta\theta} = r**2,  g_{\phi\phi}=r**2 \sin(\theta)**2 };
    complete(ss, $g$)
    complete(ss, $gt$)
    ex:= g;
    substitute(ex, ss)
    simplify(ex)
    tst:= -r**4 \sin(\theta)**2 - @(ex);
    assert(tst==0)
    print("Test 37 passed")
    ex:= gt;
    substitute(ex, ss)

test37()


def test38():
    __cdbkernel__=create_scope()
    {\mu,\nu}::Indices(values={1,2,3,4}, position=fixed);
    \partial{#}::PartialDerivative;
    {A^{1},A^{2},A^{3},A^{4}}::Depends(\partial{#});
    ex:=\partial_{\mu}{A^{\mu}};
    evaluate(ex)
    tst:= \partial_{1}{A^{1}} + \partial_{2}{A^{2}} + \partial_{3}{A^{3}} + \partial_{4}{A^{4}} - @(ex);
    assert(tst==0)
    print("Test 38 passed")

#     So this fails, because we do not check the dependence correctly.
#     But the general logic also fails if you have a \partial_{\mu}{A} \partial_{\mu}{A} node with
#     a generic dependence A::Depends(\partial{#}), because neither of the two factors will then
#         know over which values of mu to sum, and handle_derivative only sums over index values
#         for which we have a dependence.
# 
#         Whenever we encounter a depencence on \partial{#}, we need to replace that with
#             a dependence on all index values which the unique_indices->values contains.
#             Also do that when the dep is an index itself.
    

        
test38()

def test39():
    __cdbkernel__=create_scope()
    {i,j,k,l}::Indices(vector,fourD,values={0..3},position=fixed);
    {i,j,k,l}::Integer(0..3);
    \partial{#}::Derivative;
    {\kappa,\lambda,\mu,\nu}::Indices(vector, threeD,parent=fourD,values={1..3},position=fixed);
    {\kappa,\lambda,\mu,\nu}::Integer(1..3);
    \epsilon^{\lambda \mu \nu }::EpsilonTensor;
    ex:=\epsilon^{\lambda \mu \nu };
    evaluate(ex, $$);        
    assert( ex[3][1][1] == -1 )
    print("Test 39 passed")

test39()
        
print(list(map(lambda x: str(x), server.totals())))