/////////////////////////////////////////////////////////////////////////////////
//////////////   Evolution du spectre de matrices symétriques et gaussiennes  ///
//////////////   (Gaussian Orthogonal Ensemble - GOE  THEOREME DE WIGNER)     ///
/////////////////////////////////////////////////////////////////////////////////



clf();
stacksize("max");
T=1;
n=100;  //Pas de temps T/n
N=400;  //Nombre de browniens = taille des matrices 
       //Exemples pour voir le demi-cercle :  
       //N=500, n=100 => 1 minute) et (N=1000, n=100 => 2minutes

///////////////////////////////////////////////////////////////////////////////////
//Construction d'une matrice symétrique au temps n T/n basée sur les N(N+1)/2//////
//browniens indépdendants                                                    //////
///////////////////////////////////////////////////////////////////////////////////

// Simulation de N trajectoires brownienne sur n pas de temps
// (pour le calcul de la première ligne de la matrice)


W=[];
Spectre=[];



function y=sym(M)
    y=triu(M)+triu(M,+1)'+(sqrt(2)-1)*diag(diag(M));
endfunction



for  k=1:n
dW=sym(sqrt(T/n)*grand(N,N,"nor",0,1));  
W=W+(1/sqrt(N))*dW;
Spectre=[Spectre;spec(W)'];
end


/// Affichage des résultats

clf(0);xset("window",0);show_window();// Evolution temporelle des valeurs propres


    for i=1:N
        plot2d((T/n)*[1:n]',Spectre(:,i),i);
    end




clf(1);xset("window",1);show_window(); //Histogramme des valeurs propres au dernier instant


histplot(linspace(-2*sqrt(T),2*sqrt(T),20),Spectre(n,:),style=2);
x= linspace(-2*sqrt(T),2*sqrt(T))';
plot2d(x,(1/(2*%pi*T))*sqrt(4*T-x^2),style=3)