/*
Prove di calcolo tensoriale 
*/

/*
Provo a realizzare un documento wxm per il calcolo di tutte
le grandezze derivate da un assegnato tensore metrico.
In questo caso uso una dipendenza funzionale del tutto
generale ossia uso una generica funzione psi(r) dipendente
dalla sola coordinata r. In questo modo posso trattare
sia il buco nero neutro ( Schwarzschild) che quello carico 
( Reissner e Nordstrom ) ma anche un qualsiasi
altro tipo di buco nero che produca il tensore metrico
che ho imposto.
*/

(if atom(lg) then load(ctensor));

/*
Pulisco le variabili se per caso fossero state inizializzate
in modo non corretto.
*/

init_ctensor();

/*
Specifico quali variabili debbono essere considerate variabili
indipendenti. Ovviamente uso i classici nomi in geometria sferica
ovvero t, r, theta e phi.
*/

ct_coords: [t,r,theta,phi];

/*
Dichiaro che il simbolo psi deve essere considerato una funzione della
variabile r.
*/

depends(psi,r);

/*
La libreria ctensor si aspetta di trovare lg definita come una
matrice quadrata di ordine quattro interpretata come il tensore
metrico covariante.
*/

lg:matrix(
[1+psi,0,0,0],
[0,-1/(1+psi),0,0],
[0,0,-r^2,0],
[0,0,0,-(r*sin(theta))^2]);

lg;

/*
Con la funzione cmetric() creo il tensore metrico controvariante
indicato come ug e definito come la matrice inversa di lg.
*/

cmetric();

ug;

/*
A questo punto calcolo i simboli di Christoffel.
Specificando all li visualizzo tutti ossia quelli di prima 
specie contenuti nell'array lcs, sia quelli di seconda
specie contenuti nell'array ucs.
Notare che la numerazione delle array in maxima parte da 0 ma
qui gli indici 0 non sono utilizzati.
Da notare anche che l'ordine e' invertito rispetto a quello
tradizionale per cui sono simmetrici tra loro il primo ed
is secondo indice e non, come si trova su molti libri, il secondo
con il terzo.
*/

christof(all);

/*
Calcolo il tensore di Riemann che opportunamente contratto mi
fornisce il tensore di Ricci.
*/

riemann(true);

/*
Calcolo il tensore di Ricci.
*/

ricci(true);

/*
Lo visualizzo sotto forma di lista semplificandolo.
*/

ratsimp(listarray(ric));

/*
Ora carico l'array ric in una matrice che chiamo mat_ricci
*/

mat_ricci: ratsimp( matrix(
[ ric[1,1],ric[1,2],ric[1,3],ric[1,4]],
[ ric[2,1],ric[2,2],ric[2,3],ric[2,4]],
[ ric[3,1],ric[3,2],ric[3,3],ric[3,4]],
[ ric[4,1],ric[4,2],ric[4,3],ric[4,4]]));

/*
Invece dei simboli di derivata uso u0 per la funzione non derivata,
u1 per la sua derivata prima rispetto ad r e u2 per la sua derivata 
seconda. Faccio questo per usare piu' facilmente la
funzione subst. Chiamo il risultato mu_ricci.
*/

(u0:"u0",u1:"u1",u2:"u2");

mu_ricci: subst(u0,psi,subst(u2,diff(psi,r,2),subst(u1,diff(psi,r),mat_ricci)));

/*
Ora mostro che se assumo che u0 valga -2*m/r ed u1 ed u2 di conseguenza, il
tensore di Ricci si annulla ossia ho la soluzione del buco nero neutro
ovvero la soluzione di Schwarzschild. Se invece q non è nullo e pongo
u0 = -2*m/r+q^2/r^2 ottengo la soluzione di Reissner e Nordstrom.
*/

v0:-2*m/r+q^2/r^2;

v1:diff(v0,r);

v2:diff(v0,r,2);

/*
Sostituisco v0, v1 e v2 ad u0, u1 ed u2.
*/

ratsimp(subst(v0,u0,subst(v1,u1,subst(v2,u2,mu_ricci))));

/*
...etc...
*/