[scilab-Users] Problem with Scilab

[Calixte.Denizet at ac-rennes.fr]

Hello Guillaume,

Works fine for me under windows vista and scilab 5.2.1. Maybe you should load your file with the menu "Load into Scilab Ctrl+L" in the text editor... ;)

Calixte

----- Message d'origine -----
De: Guillaume GRAMPEIX <grampeix at crans.org>
Date: Mardi, Mars 2, 2010 9:52 pm
Objet: [scilab-Users] Problem with Scilab
À: users at lists.scilab.org


> Hello,
> 
> I've been worked with Scilab since few days and this evening i got a 
> problem that i don't understant. Scilab says !--error 2, Invalid 
> factor. 
> Do you know why i have this message? I send to you my file.
> 
> Cordialy
> -- 
> Guillaume Grampeix
> French student form ENPC school.
> stacksize('max');
> clear all;
> //dt utile =>0.1 à 1,7
> //-----------------------
> //  POTENTIEL UTILISE
> //-----------------------
> 
> function a = V(x)
>   a = h*((x-1).^2).*(x+1).^2;
> endfunction
> 
> function a = nablaV(x)
>   a = 2*h*( (x-1).*(x+1).^2 + (x-1).^2*(x+1) );
> endfunction
> 
> //------------------------
> //  VARIABLES A DECLARER
> //------------------------
> h = 1;
> beta = 1;
> //------------------
> // INITIALISATIONS
> //-----------------
> moy_var_v=zeros();
> var_m=zeros();
> moy_var=zeros();
> var_val=zeros();
> u=0;                           //pour obtenir un vecteur des varleurs 
> de la variance
> dt_v=zeros();
> 
>    //boucle pour different dt
>  It1=0.1;
>  for dt= 0.1 :0.2 :It1
>   u=u+1;
>   dt_v(u)=dt;                  //pour obtenir un vecteur des valeurs 
> de dt
>                    
>     //boucle pour moyenne sur N représentatif
> It2 = 20;
> for t=1:It2
>   
>   N1=14; N=2^N1;  
>   x = zeros();
>   X = x;
>   x(1) = 1;
>   V_old = V(x(1));
>   V_new = V_old; 
>   E = zeros();
>   E(1) = V_new;
>  //----------------------------------------
> //  GENERATION DE L'ENSEMBLE STATISTIQUE
> //----------------------------------------
> rejet = 0;
> for i=2:N
> 
>   V_old = V(x(i-1));
> 
>    //--------------------
>   //   PROPOSITION    
>   //--------------------
>   bruit = rand(1,'normal');
>   x_propose = x(i-1) - nablaV(x(i-1))*dt + sqrt(2*dt/beta)*bruit;
>   X(i) = x_propose;
>   V_new = V(x_propose);
>   pxy = exp(-0.5*bruit^2);
>   Bruit = ( x(i-1) - x_propose + nablaV(x_propose)*dt )/sqrt(2*dt/beta);
>   pyx = exp(-0.5*Bruit^2);
>   //-----------------------
>   //   ACCEPTION / REJET
>   //------------------------
>   p = exp(-beta*(V_new - V_old))*pyx/pxy;
>   if (p >= 1)
>     x(i) = x_propose;
>   else 
>     if (rand(1,'uniform') < p)
>       x(i) = x_propose;
>     else
>       x(i) = x(i-1);
>       rejet = rejet + 1;
>     end
>   end
> 
>   // mise a jour des observables
>   E(i) = E(i-1) + V(x(i));
>  end
> 
> // convergence de l'énergie
> //energie = E' ./ (1:N);
> //figure (1)
> //plot(energie);
> //xtitle('Energie de la particule)
> 
> //----------------------------
> //    TAUX DE REJET
> //----------------------------
> 
> rejet = rejet/N;        // valeur instantanée du rejet
> rejet_v(t)=rejet;       // valeurs du rejet pour It2 itérations
> moy_rej=mean(rejet_v);  // moyenne du rejet sur It2 itérations
> //___________________________________________________
> //         block averaging
> //___________________________________________________
> 
> //---------------
> // Initialisation
> //---------------
> 
> //n = taille du bloc
> //m = nombre de block
> 
> energie_bloc=zeros(N,1);
> var=zeros(N,1);
> moy_var=0;
> m=2;
> p=0;
> vect_m = zeros(N,1);
> vect_n = zeros(N,1);
> 
> //----------------
> // boucle
> //----------------
> 
> while m <= N;
>     n=N/m;
>     p=p+1;
> sigma_carre = 0;    
> 
>     for k=1:m
>         energie_bloc(k)= mean(x((k-1)*n+1:k*n));
>     end
>     
>     sigma_carre =n/(m-1)*sum((energie_bloc(1:m)-mean(x)).^2);
>     
>     vect_m(p)=m;
>     vect_n(p)=n;
>     var(p)=sqrt(sigma_carre);  
>     
>   m=2*m;
> end
> moy_var=sum(var)/N1;              //moyenne intantanée de la variance 
> de l'énergie
> moy_var_v(t)=moy_var;             //moyenne de la variance de 
> l'énergie pour t itérations
> //var_m(1:N1,t)=var(1:N1);          //sauvegarde des valeurs de la 
> variance de l'érnergie
> var_moy=mean(moy_var_v);          //moyenne sur It2 itérarions de la 
> variance de l'energie pour une valeur de dt
> 
>   //-------------
>   // Courbes
>   //-------------
>   
> //Energie
> set("figure_style","new");
> f = get("current_figure");
> set("current_figure",1);
> xtitle('Variance en fonction de la taille des bloc','Log de la taille 
> des blocs','Variance de l''énergie')
> plot(log(vect_n(1:p)),var(1:p),'k-+');  
> 
> //Rejet
> set("figure_style","new");
> f = get("current_figure");
> set("current_figure",2);
> xtitle('Rejet en fonction des iterations','itération It2','Rejet')
> plot(rejet_v(1:t),'-+');
> end
> 
> //------------------------------
> //    Sauvegardes de valeurs
> //------------------------------
> 
> disp(var_moy,'Variance moyenne');
> disp(moy_rej,'Rejet moyen');
> disp(dt,'dt=');
> var_m1(1:It2,u)=moy_var_v(1:It2);      // valeur de la variance de 
> l'énergie suivant It2 et dt
> var_val(u)=var_moy;                    // valeur de la variance pour 
> différent dt
> moy_rej_v(u)=moy_rej;                  // valeur du rejet pour 
> différent dt
> 
> fprintfMat('C:\\mala\\var_m1.txt',var_m1);
> fprintfMat('C:\\mala\\var_moy.txt',var_val);
> fprintfMat('C:\\mala\\rejet.txt',moy_rej_v);
> end
> 
> //-----------------------------
> //      Courbes finales
> //-----------------------------
> 
> //set("figure_style","new");
> //f = get("current_figure");
> //set("current_figure",3);
> //subplot(1,2,1);
> //xtitle('Variance de l''énergie en fonction de dt','dt','Variance de 
> l''énergie')
> //plot(dt_v,var_val(1:u),'b-+');
> //subplot(1,2,2)
> //xtitle('Taux acceptation rejet en fonction de dt','dt','Rejet')//plot(dt_v,moy_rej_v(1:u),'r-+');