/*
Prove varie di http://maxima.sourceforge.net/
*/

/*
Innanzi tutto ( solo per esercizio dimostrativo ) estendo il path in modo che sia 
accessibile anche la cartella C:\zmaxima\
Notare che in Maxima non si deve usare la controbarra ma la barra...
*/

(if atom(path_iniziale) then path_iniziale:file_search_maxima)$

file_search_maxima : cons( 
sconcat("C:/zmaxima/###.{lisp,mac,mc}"), path_iniziale)$

/*
Uso il pacchetto ctensor per studiare la metrica di Kerr e Newman ossia 
la metrica del buco nero carico e rotante ( la più generale).
La metrica data qui... non è quella ( fisicamente equivalente ) detta di 
Kerr Newman Schild caratterizzata dal fatto di avere il determinante
del tensore metrico costante ed uguale a -1.
*/

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

/*
Uso una delle metriche preconfezionate ossia appunto non quella di
Schwarzschild ma di kerr_newman nella forma di Boyer Lindquist.
*/

ct_coordsys(kerr_newman);

/*
Il tensore metrico covariante pertanto è questo:
*/

lg;

/*
Il tensore metrico covariante lg è memorizzato come una matrice
ovviamente simmetrica. Verifichiamolo:
*/

matrixp(lg);

/*
Tener presente che la numerazione delle matrici parte da 1.
*/

/*
Innanzi tutto si constata che il determinante della matrice lg
è una funzione di r e di theta.
*/

menogdet:-trigsimp(determinant(lg));

/*
Le componenti non nulle di lg sono:
*/

lg[1,1];

lg[2,2];

lg[3,3];

lg[4,4];

/*
Se a non è nulla ossia se la rotazione esiste effettivamente non
è nullo il seguente termine non diagonale:
*/

lg[1,4];

/*
Quali simboli usa come variabili del tensore metrico.... ovviamente:
*/

ct_coords;

/*
Calcolo ug ( il tensore metrico controvariante) 
*/

cmetric();

ug:trigsimp(ug)$

ug;

/*
Il determinante di ug è ovviamente l'inverso di quello di lg.
*/

trigsimp(determinant(ug));

ug[1,1];

ug[2,2];

ug[3,3];

ug[4,4];

ug[1,4];

trigsimp(ug.lg);

/*
I simboli di Christoffel di prima e seconda specie. Specificando all li
visualizzo tutti.
*/

christof(false);

lcs:trigsimp(lcs);

lcs[2,1,4];

/*
I primi due indici sono simmetrici 
*/

lcs[1,2,4]-lcs[2,1,4];

/*
Chiedo informazioni sulla array lcs ( simboli di prima specie)
*/

arrayinfo(lcs);

/*
Ed ora guardo i simboli di Christoffel di seconda specie 
*/

mcs:trigsimp(mcs);

mcs[1,1,1];

trigsimp(mcs[1,1,2]);

mcs[1,1,3];

mcs[1,1,4];

trigsimp(mcs[1,2,1]);

/*
Anche qui constato che i primi due indici sono simmetrici...
Per esempio...
*/

mcs[1,2,1]-mcs[2,1,1];

mcs[1,2,2];

mcs[1,2,3];

trigsimp(mcs[1,2,4]);

mcs[2,1,4]-mcs[1,2,4];

trigsimp(mcs[1,3,1]);

mcs[3,1,1]-mcs[1,3,1];

mcs[1,3,2];

mcs[1,3,3];

mcs[1,3,4];

mcs[3,1,4]-mcs[1,3,4];

mcs[1,4,1];

mcs[1,4,2];

mcs[4,1,2]-mcs[1,4,2];

mcs[1,4,3];

mcs[4,1,3]-mcs[1,4,3];

mcs[1,4,4];

trigsimp(mcs[2,1,1]);

mcs[2,1,1]-mcs[1,2,1];

trigsimp(mcs[2,1,2]);

trigsimp(mcs[2,1,3]);

trigsimp(mcs[2,1,4]);

mcs[1,2,4]-mcs[2,1,4];

mcs[2,2,1];

mcs[2,2,2];

mcs[2,2,3];

mcs[2,2,4];

mcs[2,3,1];

mcs[2,3,2];

mcs[3,2,2]-mcs[2,3,2];

mcs[2,3,3];

mcs[3,2,3]-mcs[2,3,3];

mcs[2,3,4];

trigsimp(mcs[2,4,1]);

mcs[4,2,1]-mcs[2,4,1];

trigsimp(mcs[2,4,2]);

trigsimp(mcs[2,4,3]);

trigsimp(mcs[2,4,4]);

mcs[4,2,4]-mcs[2,4,4];

trigsimp(mcs[3,1,1]);

trigsimp(mcs[3,1,2]);

mcs[3,1,3];

mcs[3,1,4];

mcs[3,2,1];

mcs[3,2,2];

mcs[2,3,2]-mcs[3,2,2];

mcs[3,2,3];

mcs[2,3,3]-mcs[3,2,3];

mcs[3,2,4];

mcs[3,3,1];

mcs[3,3,2];

mcs[3,3,3];

mcs[3,3,4];

trigsimp(mcs[3,4,1]);

mcs[4,3,1]-mcs[3,4,1];

mcs[3,4,2];

mcs[3,4,3];

trigsimp(mcs[3,4,4]);

mcs[4,3,4]-mcs[3,4,4];

mcs[4,1,1];

mcs[4,1,2];

mcs[1,4,2]-mcs[4,1,2];

mcs[4,1,3];

mcs[1,4,3]-mcs[4,1,3];

mcs[4,1,4];

trigsimp(mcs[4,2,1]);

mcs[2,4,1]-mcs[4,2,1];

trigsimp(mcs[4,3,4]);

mcs[3,4,4]-mcs[4,3,4];

mcs[4,4,1];

trigsimp(mcs[4,4,2]);

trigsimp(mcs[4,4,3]);

trigsimp(mcs[4,4,4]);

/*
Chiedo informazioni sulla array mcs ( simboli di Christoffel di seconda specie )
Usa la convenzione di scrivere il termine controvariante come ultimo.
*/

arrayinfo(mcs);

/*
Ora calcolo il tensore di Riemann a quattro indici, visualizzandone solo
quelli di valore non nullo ( scrivendo true ma la stampa è troppo lunga ).
*/

riemann(false);

/*
Calcola il tensore di Ricci che dovrebbe essere identicamente nullo
ma siccome non fa semplificazioni viene così....
*/

ricci(false);

/*
Trasformo l'array ric in una lista per vederla in modo convenzionale
*/

listaric:trigsimp(listarray(ric));

/*
Ricopio l'array in una matrice. Notare che la numerazione degli 
elementi della array potrebbe partire da zero ma gli indici
con zero non vengono usati per cui il tensore di Ricci 
è questo, semplificato con la semplificazione trigonometrica:
*/

matricci: trigsimp(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]]));

/*
Se semplifico le espressioni...
*/

/*
L'array ric è di tipo hashed ossia sparsa perché non vengono definiti
gli indici della riga 0 e della colonna 0.
*/

arrayinfo(ric);

/*
Questa e' la verifica della metrica di Kerr e Newman. 
*/

matricci[1,1];

/*
azzerando a...
*/

aa:subst(0,a,matricci[1,1]);

/*
Se invece si pone nulla la carica e ma si fa la semplificazione
trigonometrica si deve ottenere 0 perché ci si riduce alla metrica
di Kerr dove non c'è energia nello spazio attorno al buco
nero rotante.
*/

/*
Azzero e....
*/

bb:trigsimp(subst(0,e,matricci[1,1]));

matricci[1,2];

matricci[1,3];

/*
Evidentemente se azzero a o e ottengo zero...
*/

matricci[1,4];

matricci[2,1];

aa:subst(0,a,matricci[2,2]);

matricci[2,3];

matricci[2,4];

aa:subst(0,a,matricci[3,3]);

aa:subst(0,a,matricci[4,4]);

/*
Verifiche: debbono risultare uguali:
*/

trigsimp(diff(menogdet,r)/(2*menogdet));

trigsimp(mcs[2,1,1]+mcs[2,2,2]+mcs[2,3,3]+mcs[2,4,4]);

/*
Giusto! Vedere il Landau per questa importante formula di verifica
che deve valere SEMPRE.
*/

/*
Ed anche questa coppia deve risultare uguale anche se Maxima fatica 
un po' ad accorgersene...
*/

va:trigsimp(diff(menogdet,theta)/(2*menogdet));

vb:trigsimp(mcs[3,1,1]+mcs[3,2,2]+mcs[3,3,3]+mcs[3,4,4]);

/*
Possono sembrare due espressioni diverse ma non lo sono...
*/

trigsimp(va-vb);

/*
PROSIT !
*/