Contenu | Rechercher | Menus

Annonce

Si vous avez des soucis pour rester connecté, déconnectez-vous puis reconnectez-vous depuis ce lien en cochant la case
Me connecter automatiquement lors de mes prochaines visites.

À propos de l'équipe du forum.

#1 Le 23/10/2007, à 13:44

laroche1

aide programmation matlab

Bonjour à tous,
Est ce que parmis vous, il y aurait un pro ou plusieurs de la programmation matlab??
Merci

Hors ligne

#2 Le 23/10/2007, à 14:04

toufalk

Re : aide programmation matlab

pro, je sais pas, mais tu peux toujours poser ta question...
Et pourquoi utiliser un truc propriétaire pas donné quand il y a scilab, un équivalent "semi-libre" aussi performant??

Hors ligne

#3 Le 23/10/2007, à 14:12

snigit

Re : aide programmation matlab

Parfois on n'a pas le choix, quand les instances supérieures l'imposent...

Hors ligne

#4 Le 23/10/2007, à 14:41

laroche1

Re : aide programmation matlab

tout a fait snigit je n'ai pas le choix...

Voila j'explique donc mon probleme qui risque d'être long et compliqué.
J'ai recupéré le programme de quelq'un pour etoffer les résultats qu'ils n'avait pas eu le temps de finir. Le probleme c'est que cette personne n'est plus joignable et le programme ne fonctionne plus. En fait cette personne a retranscrit son programme sous word mais a fait des erreurs de frappe et le programme matlab a disparu !!!

Je suis debutante, je ne comprends pas tout le programme. A la base c'est un informaticien qui a fait ce programme

Je vous joins tout ce qu'on m'a donné. 

Il y a aussi une notice qui explique en principe le programme (je dis bien en principe)

INVERSION PROGRAM (MATLAB© 5.3.0.1)

All the files must be within the WORK folder under MATLAB
For additional information on the program or on its application, please contact ncoussae@africamuseum.be



To prepare the Input Files

Prior data file.
Create a table where each line corresponds to the values of a specific datum type and each column to the set of values for a given experiment.  On each line, individual data are separated by a space or a tab.
Example:
P    0.2    0.5    1    1    1.5
T    1473    1473    1373    1473    1273
Ca cpx    0.802    0.779    0.872    0.806    0.894
Ca opx    0.080    0.085    0.056    0.075    0.032

Insert this table into a data.txt file (discard the first column):

  0.2    0.5    1    1    1.5
1473    1473    1373    1473    1273
0.802    0.779    0.872    0.806    0.894   
0.080    0.085    0.056    0.075    0.032

Prior error file.
Create a table where each line corresponds to the variances of a specific datum type and each column to the set of variances for a given experiment.  As above, the individual data are separated by a space or a tab on each line.

P    4.444E-07    2.778E-06    1.111E-05    1.111E-05    2.500E-05
T    1.361E+00    1.361E+00    1.361E+00    1.361E+00    1.361E+00
Ca cpx    4.391E-05    2.323E-05    1.533E-05    6.254E-05    2.237E-05
Ca opx    5.082E-05    4.799E-05    9.072E-06    2.622E-05    3.906E-05

Insert this table into a error.txt file (discard the first column):

4.444E-07    2.778E-06    1.111E-05    1.111E-05    2.500E-05
1.361E+00    1.361E+00    1.361E+00    1.361E+00    1.361E+00
4.391E-05    2.323E-05    1.533E-05    6.254E-05    2.237E-05
5.082E-05    4.799E-05    9.072E-06    2.622E-05    3.906E-05




How to use the program

Paste the program onto the programming MATLAB® window
Select the prior parameter values through the line :

X0=[100000 10000 1000 0 0 10 0 0 0 100000 10000 1000 m']';

The values must follow the order specified in the comment :
  %line matrix of the posterior parameter values (  deltaU....                ) to which is added the line matrix (m’) of the posterior data.

Select the prior variances on the parameters through line E, e.g.,

E=[1E9;1E8;1E5;0;0;1000;0;0;0;1E9;1E8;1E5;e];

The values must follow the order specified in the comment :
  %column matrix of the prior errors on the parameters (  deltaU...                ) to which is added the column matrix (e) of the prior errors on the data.

Select the number of iterations through the line:
for h=1:50

Save the program.
Start the program within the MATLAB® application window.





Results

When the program has ended, it creates two files within the MATLAB® work folder

reslt.txt : result file for the parameters;  the first line recalls the prior variances, the second gives the posterior values and the others are the posterior covariance matrix.

rcomp.txt : line matrix of the recalculated experimental data


Program:

% M = input matrix of the prior data (format as shown in the .xls file; orthopyroxene data follow clinopyroxene data, unlike in the text or titles)
% m = column matrix of the prior data
% X0 = line matrix of the prior parameters and data
% e = column matrix of the errors
% E = column matrix of the prior errors on the parameters and data
% C0 = creation of a nth order diagonal matrix, where n is the total number of parameters and data
% placement of the errors onto the diagonal to create the prior covariance matrix
% X = matrix of the prior data as used in iteration
% a = number of data
% start of inversion
% h = number of iterations
% p = form of the thermodynamic law used for the En-transfer reaction
% application of p to the whole set of data
% P1 = p derivative for the thermodynamic parameters of the En-reaction
% P2 = p derivative for the thermodynamic parameters of the Di-reaction
% P3 = p derivative for the various pressures
% P4 = p derivative for the various temperatures
% P5 = p derivative for the various Ca contents in clinopyroxene
% P6 = p derivative for the various Ca contents in orthopyroxne
% P whole set of p derivatives
% f = form of the thermodynamic law used for the Di-transfer reaction
% application of p to the whole set of data
% F1 = f derivative for the thermodynamic parameters of the En-reaction
% F2 = f derivative for the thermodynamic parameters of the Di-reaction
% F3 = f derivative for the various pressures
% F4 = f derivative for the various temperatures
% F5 = f derivative for the various Ca contents in clinopyroxene
% F6 = f derivative for the various Ca contents in orthopyroxne
% F whole set of f derivatives
% g = set of the two thermodynamic laws describing both reactions
% G= set of the derivatives for the two reactions
% mathematical form of the Tarantola et Valette's theorem
% parameter intermediate values during inversion
% end of inversion 
% C = computation of the final covariance matrix
% W = final covariance matrix for the parameters
% Z = prior variance matrix before inversion
% Q = set of results and variances for the parameters
% k = matrix of the recalculated data after inversion

voici le code
C'est un code qui permet de faire une inversion generale à partir de données experimentales pour trouver des paramètres thermodynamiques tel que deltaU, deltaV, deltaS

clear
M=DLMREAD('/data.txt','\t');			
j=size(M);
m=DLMREAD('/data.txt');					
X0=[100000	10000 1000 0 0 10 0 0 0 100000	10000	1000 m']';	
e=DLMREAD('/error.txt');		
E=[1E9;1E8;1E5;0;0;1000;0;0;0;1E9;1E8;1E5;e];				
D=size(X0);								
C0=eye(D(1));								
for i=1:D(1)								
   C0(i,i)= E(i);							
X=X0;											
a=[j(1,2)];
for h=1:50									
   p=[X(1)+X(12+1)*X(2)-X(12+a+1)*X(3)+8.314*X(12+a+1)*log((1-X(12+2*a+1))/(1-X(12+3*a+1)))+((X(12+2*a+1))^2)*(X(4)+X(12+1)*X(5)-X(6)*X(12+a+1))-((X(12+3*a+1))^2)*(X(7)+X(8)*X(12+1)-X(9)*X(12+a+1))];
    		for b=14:a+12
   		p=[p;X(1)+X(b)*X(2)-X(a+b)*X(3)+8.314*X(a+b)*log((1-X(2*a+b))/(1-X(3*a+b)))+((X(2*a+b))^2)*(X(4)+X(b)*X(5)-X(6)*X(a+b))-((X(3*a+b))^2)*(X(7)+X(8)*X(b)-X(9)*X(b+a))];
      	end
   P1=[1 X(13) -X(12+a+1) (X(12+2*a+1))^2 (X(12+2*a+1))^2*X(13) -(X(12+2*a+1))^2*X(12+a+1) -(X(12+3*a+1))^2 -(X(12+3*a+1))^2*X(13) (X(12+3*a+1))^2*X(12+a+1)];
   	for z=14:a+12
         P1=[P1;1 X(z) -X(z+a) (X(z+2*a))^2 (X(z+2*a))^2*X(z) -(X(z+2*a))^2*X(z+a) -(X(z+3*a))^2 -(X(z+3*a))^2*X(z) (X(z+3*a))^2*X(z+a)];
      end
   P2=zeros(a,3);   
   P3=eye(a);
   P4=eye(a);
   P5=eye(a);
   P6=eye(a);
     	for y=1:a
         P3(y,y)=X(2)+(X(12+y+2*a))^2*X(5)-(X(12+y+3*a))^2*X(8);
   		P4(y,y)=-X(3)+8.314*log((1-X(2*a+12+y))/(1-X(3*a+12+y)))-X(6)*(X(2*a+12+y))^2+X(9)*(X(3*a+12+y))^2;
         P5(y,y)=-8.314*X(12+a+y)/(1-X(2*a+12+y))+2*X(12+2*a+y)*(X(4)+X(12+y)*X(5)-X(12+a+y)*X(6));
         P6(y,y)=8.314*X(12+a+y)/(1-X(3*a+12+y))-2*X(12+3*a+y)*(X(7)+X(8)*X(12+y)-X(12+a+y)*X(9));
      end 
   P=[P1 P2 P3 P4 P5 P6];   
   f=[X(10)+X(12+1)*X(11)-X(a+12+1)*X(12)+8.314*X(12+a+1)*log(X(12+3*a+1)/X(12+2*a+1))-(1-X(12+2*a+1))^2*(X(4)+X(12+1)*X(5)-X(6)*X(12+1+a))+(1-X(12+3*a+1))^2*(X(7)+X(8)*X(12+1)-X(9)*X(12+a+1))];
    	for b=14:a+12
   	f=[f;X(10)+X(b)*X(11)-X(a+b)*X(12)+8.314*X(a+b)*log(X(3*a+b)/X(2*a+b))-(1-X(2*a+b))^2*(X(4)+X(b)*X(5)-X(6)*X(a+b))+(1-X(3*a+b))^2*(X(7)+X(8)*X(b)-X(9)*X(b+a))];
      end
   F2=[-(1-X(12+2*a+1))^2 -(1-X(12+2*a+1))^2*X(12+1) (1-X(12+1+2*a))^2*X(12+a+1) (1-X(12+3*a+1))^2 (1-X(12+3*a+1))^2*X(12+1) -(1-X(12+3*a+1))^2*X(12+a+1) 1 X(12+1) -X(12+a+1)];
   	for z=14:a+12
         F2=[F2;-(1-X(z+2*a))^2 -(1-X(z+2*a))^2*X(z) (1-X(z+2*a))^2*X(z+a) (1-X(z+3*a)^2) (1-X(z+3*a))^2*X(z) -(1-X(z+3*a)^2)*X(z+a) 1 X(z) -X(z+a)];
      end
   F1=zeros(a,3);   
   F3=eye(a);
   F4=eye(a);
   F5=eye(a);
   F6=eye(a);
   	for y=1:a
         F3(y,y)=X(11)-(1-X(12+y+2*a))^2*X(5)+(1-X(12+y+3*a))^2*X(8);
   		F4(y,y)=-X(12)+8.314*log(X(3*a+12+y)/X(2*a+12+y))+X(6)*(1-X(2*a+12+y))^2-X(9)*(1-X(3*a+12+y))^2;
         F5(y,y)=-8.314*X(12+a+y)/X(2*a+12+y)+2*(1-X(12+2*a+y))*(X(4)+X(12+y)*X(5)-X(12+a+y)*X(6));
         F6(y,y)=8.314*X(12+a+y)/X(3*a+12+y)+2*(X(12+3*a+y)-1)*(X(7)+X(8)*X(12+y)-X(12+y+a)*X(9));
    	end
  	F=[F1 F2 F3 F4 F5 F6];
   g=[p;f];  
   G=[P;F];    
   X=X0+C0*G'*inv(G*C0*G')*(G*(X-X0)-g);
   R=[X(1) X(2) X(3) X(4) X(5) X(6) X(7) X(8) X(9) X(10) X(11) X(12)]
  
end
C=C0-C0*G'*inv(G*C0*G')*G*C0;
W=[C(1,1) C(1,2) C(1,3) C(1,4) C(1,5) C(1,6) C(1,7) C(1,8) C(1,9) C(1,10) C(1,11) C(1,12)];
for v=2:12
   W=[W; C(v,1) C(v,2) C(v,3) C(v,4) C(v,5) C(v,6) C(v,7) C(v,8) C(v,9) C(v,10) C(v,11) C(v,12)];
end
D=[E(1) E(2) E(3) E(4) E(5) E(6) E(7) E(8) E(9) E(10) E(11) E(12)];
Q=[D;R;W];
for i=1:a
   ka(i)=[X(i+12)];
   kb(i)=[X(i+12+a)];
   kc(i)=[X(i+12+2*a)];
   kd(i)=[X(i+12+3*a)];
end
K=[ka;kb;kc;kd];

dlmwrite('/reslt.txt',Q,'\t')
dlmwrite('/rcomp.txt',K,'\t')

les données de data.txt :
0.2    0.5    1     1     1.5
1473     1473     1373     1473     1273
0.802     0.779     0.872     0.806    0.894
0.080     0.085     0.056     0.075    0.032

les données error.txt
4.444E-07 2.778E-06 1.111E-05 1.111E-05    2.500E-05
1.361E+00 1.361E+00 1.361E+00 1.361E+00    1.361E+00
4.391E-05 2.323E-05 1.533E-05 6.254E-05    2.237E-05
5.082E-05 4.799E-05 9.072E-06 2.622E-05    3.906E-05

Voila j'espere avoir été claire, je m'en excuse si c'est pas le cas.


Le message erreur :

??? Error using ==> horzcat
All rows in the bracketed expression must have the same 
number of columns.

Error in ==> W:\matlab\code.m
On line 12  ==> X0=[100000;10000;1000;0;0;10;0;0;0;100000;10000;1000;m']

Merci

Dernière modification par laroche1 (Le 23/10/2007, à 14:43)

Hors ligne

#5 Le 23/10/2007, à 15:14

toufalk

Re : aide programmation matlab

Je comprends bien qu'il râle... mais comme je ne comprends pas le programme, je ne peux que te livrer mon analyse.

Dans la notice il est ecrit

Select the prior parameter values through the line :

X0=[100000 10000 1000 0 0 10 0 0 0 100000 10000 1000 m']';

The values must follow the order specified in the comment :
  %line matrix of the posterior parameter values (  deltaU....                ) to which is added the line matrix (m’) of the posterior data.

Il faut donc modifier X0 pour qu'il corresponde à ce que tu cherche. La syntaxe est expliquée dans la ligne de commentaire que tu ne nous a pas donné hmm
De même, tu aura à modifier E et h.
Tel quel, ton programme ne peux pas fonctionner car il ne peux pas construire le fameux XO tel qu'il est écrit (problème de taille de matrices)...

Dernière modification par toufalk (Le 23/10/2007, à 15:20)

Hors ligne

#6 Le 23/10/2007, à 15:43

laroche1

Re : aide programmation matlab

OK merci deja de m'eclaircir le soucis

Moi aussi j'ai du mal à comprendre pourquoi il a fait tel ou tel truc...

Excuse en pour la syntaxe mais je n'arrive pas à le mettre en fait c'est la signification des valeurs dans X0 (delta U delta V ...) et dans E ce sont les erreurs de ces données

Je pense qu'il a voulu ajouter les données (data.txt ) dans le vecteur Xo mais je ne sais pas comment il a pu faire

Comment puis je faire pour joindre un document pour te montrer??
Merci


j'avais essayer ceci sur un conseil pour info

m=m.';
X0=[100000;10000;1000;0;0;10;0;0;0;100000;10000;1000;m(:)]

ps : ce forum est mon dernier recours ......snif

Dernière modification par laroche1 (Le 23/10/2007, à 15:47)

Hors ligne

#7 Le 24/10/2007, à 09:08

laroche1

Re : aide programmation matlab

Bonjour après avoir modifié les erreurs sur Xo et E

j'ai une erreur à cette ligne

??? Error using ==> *
Inner matrix dimensions must agree.

Error in ==> W:\matlab\code.m
On line 77  ==>    X=X0+C0*G'*inv(G*C0*G')*(G*(X-X0)-g);

Si vous pouviez m'aider à la rectifier ça serait sympa
Merci

Hors ligne

#8 Le 24/10/2007, à 10:16

vinc-mai

Re : aide programmation matlab

J'ai pas regardé ton programme mais il doit s'agir d'une erreur de taille de matrices:

pour que le produit matriciel soit défini entre 2 matrices de taille (m1,n1) et (m2,n2) il faut que n1=m2.

Si tu désires construire une matrice dont chaque terme est le produit des termes des 2 matrices de départ, il faut:
que les matrices aient la même taille,
utiliser lopérateur ".*" et non "*".


ZiK un lecteur audio et son blog.

Hors ligne

#9 Le 24/10/2007, à 10:30

laroche1

Re : aide programmation matlab

Merci pour le truc mais ça ne marche pas ,c'est sur que c'est un soucis de taille
J'ai verifié les tailles ça me parait correct..
je mets le code rectifié au cas ou

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Program for ideal orthopyroxene / symmetric clinopyroxene models

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear
M=dlmread('W:/matlab/data.txt','\t');			
j=size(M);
m=dlmread('W:/matlab/data.txt');	

%if all(m(:,end)==0)
  % m(:,end)=[];   
  %end
m=m.';
X0=[100000;10000;1000;0;0;10;0;0;0;100000;10000;1000;m(:)]'
%E=dlmread('W:/matlab/error.txt','\t');
e=dlmread('W:/matlab/error.txt');
e=e.';
E=[1E9;1E8;1E5;0;0;1000;0;0;0;1E9;1E8;1E5;e(:)]
D=size(X0);								
C0=eye(D(1));								
    for i=1:D(1)								
      C0(i,i)= E(i);							
        X=X0;											
        a=[j(1,2)];
        
        for h=1:50									
            p=[X(1)+X(12+1)*X(2)-X(12+a+1)*X(3)+8.314*X(12+a+1)*log((1-X(12+2*a+1))/(1-X(12+3*a+1)))+((X(12+2*a+1))^2)*(X(4)+X(12+1)*X(5)-X(6)*X(12+a+1))-((X(12+3*a+1))^2)*(X(7)+X(8)*X(12+1)-X(9)*X(12+a+1))];
    		    
           
            
            
            for b=14:a+12
   		            p=[p;X(1)+X(b)*X(2)-X(a+b)*X(3)+8.314*X(a+b)*log((1-X(2*a+b))/(1-X(3*a+b)))+((X(2*a+b))^2)*(X(4)+X(b)*X(5)-X(6)*X(a+b))-((X(3*a+b))^2)*(X(7)+X(8)*X(b)-X(9)*X(b+a))];
            	end
            P1=[1 X(13) -X(12+a+1) (X(12+2*a+1))^2 (X(12+2*a+1))^2*X(13) -(X(12+2*a+1))^2*X(12+a+1) -(X(12+3*a+1))^2 -(X(12+3*a+1))^2*X(13) (X(12+3*a+1))^2*X(12+a+1)];
            
                for z=14:a+12
                     P1=[P1;1 X(z) -X(z+a) (X(z+2*a))^2 (X(z+2*a))^2*X(z) -(X(z+2*a))^2*X(z+a) -(X(z+3*a))^2 -(X(z+3*a))^2*X(z) (X(z+3*a))^2*X(z+a)];
                end
   P2=zeros(a,3);   
   P3=eye(a);
   P4=eye(a);
   P5=eye(a);
   P6=eye(a);
     	       for y=1:a
                     P3(y,y)=X(2)+(X(12+y+2*a))^2*X(5)-(X(12+y+3*a))^2*X(8);
   		            P4(y,y)=-X(3)+8.314*log((1-X(2*a+12+y))/(1-X(3*a+12+y)))-X(6)*(X(2*a+12+y))^2+X(9)*(X(3*a+12+y))^2;
                   P5(y,y)=-8.314*X(12+a+y)/(1-X(2*a+12+y))+2*X(12+2*a+y)*(X(4)+X(12+y)*X(5)-X(12+a+y)*X(6));
                  P6(y,y)=8.314*X(12+a+y)/(1-X(3*a+12+y))-2*X(12+3*a+y)*(X(7)+X(8)*X(12+y)-X(12+a+y)*X(9));
               end 
    P=[P1 P2 P3 P4 P5 P6];   
   f=[X(10)+X(12+1)*X(11)-X(a+12+1)*X(12)+8.314*X(12+a+1)*log(X(12+3*a+1)/X(12+2*a+1))-(1-X(12+2*a+1))^2*(X(4)+X(12+1)*X(5)-X(6)*X(12+1+a))+(1-X(12+3*a+1))^2*(X(7)+X(8)*X(12+1)-X(9)*X(12+a+1))];
        	for b=14:a+12
         	    f=[f;X(10)+X(b)*X(11)-X(a+b)*X(12)+8.314*X(a+b)*log(X(3*a+b)/X(2*a+b))-(1-X(2*a+b))^2*(X(4)+X(b)*X(5)-X(6)*X(a+b))+(1-X(3*a+b))^2*(X(7)+X(8)*X(b)-X(9)*X(b+a))];
             end
   F2=[-(1-X(12+2*a+1))^2 -(1-X(12+2*a+1))^2*X(12+1) (1-X(12+1+2*a))^2*X(12+a+1) (1-X(12+3*a+1))^2 (1-X(12+3*a+1))^2*X(12+1) -(1-X(12+3*a+1))^2*X(12+a+1) 1 X(12+1) -X(12+a+1)];
   
          for z=14:a+12
            F2=[F2;-(1-X(z+2*a))^2 -(1-X(z+2*a))^2*X(z) (1-X(z+2*a))^2*X(z+a) (1-X(z+3*a)^2) (1-X(z+3*a))^2*X(z) -(1-X(z+3*a)^2)*X(z+a) 1 X(z) -X(z+a)];
          end
   F1=zeros(a,3);   
   F3=eye(a);
   F4=eye(a);
   F5=eye(a);
   F6=eye(a);
         	for y=1:a
              F3(y,y)=X(11)-(1-X(12+y+2*a))^2*X(5)+(1-X(12+y+3*a))^2*X(8);
   		    F4(y,y)=-X(12)+8.314*log(X(3*a+12+y)/X(2*a+12+y))+X(6)*(1-X(2*a+12+y))^2-X(9)*(1-X(3*a+12+y))^2;
             F5(y,y)=-8.314*X(12+a+y)/X(2*a+12+y)+2*(1-X(12+2*a+y))*(X(4)+X(12+y)*X(5)-X(12+a+y)*X(6));
              F6(y,y)=8.314*X(12+a+y)/X(3*a+12+y)+2*(X(12+3*a+y)-1)*(X(7)+X(8)*X(12+y)-X(12+y+a)*X(9));
    	    end
  	F=[F1 F2 F3 F4 F5 F6];
   g=[p;f];  
   G=[P;F];    
 X=X0+C0*G'*inv(G*C0*G')*(G*(X-X0)-g);
   R=[X(1) X(2) X(3) X(4) X(5) X(6) X(7) X(8) X(9) X(10) X(11) X(12)]
  
    end
C=C0-C0*G'*inv(G*C0*G')*G*C0;
W=[C(1,1) C(1,2) C(1,3) C(1,4) C(1,5) C(1,6) C(1,7) C(1,8) C(1,9) C(1,10) C(1,11) C(1,12)];
        for v=2:12
            W=[W; C(v,1) C(v,2) C(v,3) C(v,4) C(v,5) C(v,6) C(v,7) C(v,8) C(v,9) C(v,10) C(v,11) C(v,12)];
        end
D=[E(1) E(2) E(3) E(4) E(5) E(6) E(7) E(8) E(9) E(10) E(11) E(12)];
Q=[D;R;W];
        for i=1:a
       ka(i)=[X(i+12)];
      kb(i)=[X(i+12+a)];
      kc(i)=[X(i+12+2*a)];
      kd(i)=[X(i+12+3*a)];
        end
K=[ka;kb;kc;kd];
end
dlmwrite('W:/matlab/reslt.txt',Q,'\t')
dlmwrite('W:/matlab/rcomp.txt',K,'\t')

Hors ligne

#10 Le 25/10/2007, à 09:24

vinc-mai

Re : aide programmation matlab

Sans les fichiers data.txt et error.txt pour essayer, çà semble difficile de trouver l'erreur!


ZiK un lecteur audio et son blog.

Hors ligne

#11 Le 25/10/2007, à 10:10

laroche1

Re : aide programmation matlab

Bonjour alors je mets tout codifié

J'ai un soucis je ne retrouve pas le meme resultat , Matlab ferait-il des arrondis??


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%Program for symmetric orthopyroxene / asymmetric clinopyroxene models

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear
M=DLMREAD('W:/matlab/data.txt','\t');			
j=size(M);
m=DLMREAD('W:/matlab/data.txt');		
m=m.';
X0=[100000;10000;1000;1;1;0;0;1;100000;10000;1000;m(:)];			
e=DLMREAD('W:/matlab/error.txt');	
e=e.';
	E=[1E9;1E8;1E5;1000;1000;0;0;1000;1E9;1E8;1E5;e(:)];			
D=size(X0);									
C0=eye(D(1));								
for i=1:D(1)								
   C0(i,i)= E(i);							
end
X=X0;											
a=[j(1,2)];
for h=1:50
p=[X(1)+X(11+1)*X(2)-X(3)*X(11+a+1)+8.314*X(11+a+1)*log((1-X(11+2*a+1))/(1-X(11+3*a+1)))-X(11+2*a+1)^2*X(11+a+1)*X(4)-2*X(11+2*a+1)^2*(1-X(11+2*a+1))*X(11+a+1)*(X(5)-X(4))-X(11+3*a+1)^2*(X(6)+X(7)*X(11+1)-X(8)*X(11+a+1))];
    		for b=2:a
   		p=[p;X(1)+X(11+b)*X(2)-X(3)*X(11+a+b)+8.314*X(11+a+b)*log((1-X(11+2*a+b))/(1-X(11+3*a+b)))-X(11+2*a+b)^2*X(11+a+b)*X(4)-2*X(11+2*a+b)^2*(1-X(11+2*a+b))*X(11+a+b)*(X(5)-X(4))-X(11+3*a+b)^2*(X(6)+X(7)*X(11+b)-X(8)*X(11+a+b))];
      	end
   P1=[1 X(11+1) -X(11+a+1) -X(11+2*a+1)^2*X(11+a+1)+2*X(11+2*a+1)^2*(1-X(11+2*a+1))*X(11+a+1) -2*X(11+2*a+1)^2*(1-X(11+2*a+1))*X(11+a+1) -X(11+3*a+1)^2 -X(11+3*a+1)^2*X(11+1) X(11+3*a+1)^2*X(11+a+1)];
   	for b=2:a
         P1=[P1;1 X(11+b) -X(11+a+b) -X(11+2*a+b)^2*X(11+a+b)+2*X(11+2*a+b)^2*(1-X(11+2*a+b))*X(11+a+b) -2*X(11+2*a+b)^2*(1-X(11+2*a+b))*X(11+a+b) -X(11+3*a+b)^2 -X(11+3*a+b)^2*X(11+b) X(11+3*a+b)^2*X(11+a+b)];
      end
   P2=zeros(a,3);   
   P3=eye(a);
   P4=eye(a);
   P5=eye(a);
   P6=eye(a);
     	for y=1:a
         P3(y,y)=X(2)-X(11+3*a+y)^2*X(7);
   		P4(y,y)=-X(3)+8.314*log((1-X(11+2*a+y))/(1-X(11+3*a+y)))-X(11+2*a+y)^2*X(4)-2*X(11+2*a+y)^2*(1-X(11+2*a+y))*(X(5)-X(4))+X(11+3*a+y)^2*X(8);
         P5(y,y)=-8.314*X(11+a+y)/(1-X(11+2*a+y))-2*X(11+2*a+y)*X(11+a+y)*X(4)-4*X(11+2*a+y)*(1-X(11+2*a+y))*X(11+a+y)*(X(5)-X(4))+2*X(11+2*a+y)^2*X(11+a+y)*(X(5)-X(4));
         P6(y,y)=8.314*X(11+a+y)/(1-X(11+3*a+y))-2*X(11+3*a+y)*(X(6)+X(7)*X(11+y)-X(8)*X(11+a+y));
      end 
      P=[P1 P2 P3 P4 P5 P6];  
      f=[X(9)+X(10)*X(11+1)-X(11)*X(11+a+1)+8.314*X(11+a+1)*log(X(11+3*a+1)/X(11+2*a+1))+(1-X(11+2*a+1))^2*X(11+a+1)*X(5)+2*(1-X(11+2*a+1))^2*X(11+2*a+1)*X(11+a+1)*(X(4)-X(5))+(1-X(11+3*a+1))^2*(X(6)+X(7)*X(11+1)-X(8)*X(11+a+1))];
    	for b=2:a
   	f=[f;X(9)+X(10)*X(11+b)-X(11)*X(11+a+b)+8.314*X(11+a+b)*log(X(11+3*a+b)/X(11+2*a+b))+(1-X(11+2*a+b))^2*X(11+a+b)*X(5)+2*(1-X(11+2*a+b))^2*X(11+2*a+b)*X(11+a+b)*(X(4)-X(5))+(1-X(11+3*a+b))^2*(X(6)+X(7)*X(11+b)-X(8)*X(11+a+b))];
      end
   F2=[2*(1-X(11+2*a+1))^2*X(11+2*a+1)*X(11+a+1) (1-X(11+2*a+1))^2*X(11+a+1)-2*(1-X(11+2*a+1))^2*X(11+2*a+1)*X(11+a+1) (1-X(11+3*a+1))^2 (1-X(11+3*a+1))^2*X(11+1) -(1-X(11+3*a+1))^2*X(11+a+1) 1 X(11+1) -X(11+a+1)];
   	for b=2:a
         F2=[F2;2*(1-X(11+2*a+b))^2*X(11+2*a+b)*X(11+a+b) (1-X(11+2*a+b))^2*X(11+a+b)-2*(1-X(11+2*a+b))^2*X(11+2*a+b)*X(11+a+b) (1-X(11+3*a+b))^2 (1-X(11+3*a+b))^2*X(11+b) -(1-X(11+3*a+b))^2*X(11+a+b) 1 X(11+b) -X(11+a+b)];
      end
   F1=zeros(a,3);   
   F3=eye(a);
   F4=eye(a);
   F5=eye(a);
   F6=eye(a);
   	for y=1:a
         F3(y,y)=X(10)+(1-X(11+3*a+y))^2*X(7);
   		F4(y,y)=-X(11)+8.314*log(X(11+3*a+y)/X(11+2*a+y))+(1-X(11+2*a+y))^2*X(5)+2*(1-X(11+2*a+y))^2*X(11+2*a+y)*(X(4)-X(5))-(1-X(11+3*a+y))^2*X(8);
         F5(y,y)=-8.314*X(11+a+y)/X(11+2*a+y)-2*(1-X(11+2*a+y))*X(11+a+y)*X(5)-4*(1-X(11+2*a+y))*X(11+2*a+y)*X(11+a+y)*(X(4)-X(5))+2*(1-X(11+2*a+y))^2*X(11+a+y)*(X(4)-X(5));
         F6(y,y)=8.314*X(11+a+y)/X(11+3*a+y)-2*(1-X(11+3*a+y))*(X(6)+X(7)*X(11+y)-X(8)*X(11+a+y));
    	end
   F=[F1 F2 F3 F4 F5 F6];
   g=[p;f];  
   G=[P;F];    
   X=X0+C0*G'*inv(G*C0*G')*(G*(X-X0)-g);
   R=[X(1) X(2) X(3) X(4) X(5) X(6) X(7) X(8) X(9) X(10) X(11)]
  
end
C=C0-C0*G'*inv(G*C0*G')*G*C0;
W=[C(1,1) C(1,2) C(1,3) C(1,4) C(1,5) C(1,6) C(1,7) C(1,8) C(1,9) C(1,10) C(1,11)];
for v=2:11
   W=[W; C(v,1) C(v,2) C(v,3) C(v,4) C(v,5) C(v,6) C(v,7) C(v,8) C(v,9) C(v,10) C(v,11)];
end
D=[E(1) E(2) E(3) E(4) E(5) E(6) E(7) E(8) E(9) E(10) E(11)];
Q=[D;R;W];
for i=1:a											
   ka(i)=[X(i+11)];
   kb(i)=[X(i+11+a)];
   kc(i)=[X(i+11+2*a)];
   kd(i)=[X(i+11+3*a)];
end
K=[ka;kb;kc;kd];

dlmwrite('W:/matlab/reslt3.txt',Q,'\t')
dlmwrite('W:/matlab/rcomp3.txt',K,'\t')

data:

0.2	0.5	1.0	1.0	1.5	1.5	1.5	1.5	1.5	2.0	2.0	2.0	2.0	2.0	2.0	2.0	2.0	2.0	2.0	2.0	2.0	2.0	2.0	2.5	2.5	3.0	3.0	3.0	3.0	3.0	4.0	4.0	4.0	4.0	5.0	5.0	5.0	5.0	6.0
1473	1473	1373	1473	1273	1373	1473	1573	1673	1173	1273	1373	1473	1573	1673	1273	1373	1473	1573	1673	1173	1273	1673	1176	1273	1264	1173	1473	1773	1773	1373	1473	1573	1773	1373	1573	1673	1773	1573
0.8023	0.7790	0.8725	0.8065	0.8940	0.8655	0.8355	0.7550	0.6240	0.9320	0.9075	0.8810	0.8355	0.7765	0.6445	0.9090	0.8710	0.8320	0.7665	0.6300	0.9380	0.9120	0.6275	0.9420	0.9130	0.9100	0.9340	0.8590	0.5635	0.5550	0.8710	0.8570	0.7870	0.5260	0.8830	0.8000	0.7380	0.5970	0.8080
0.0805	0.0855	0.0560	0.0755	0.0320	0.0450	0.0610	0.0865	0.1175	0.0240	0.0360	0.0440	0.0550	0.0720	0.0970	0.0360	0.0430	0.0545	0.0795	0.1055	0.0210	0.0300	0.1025	0.0180	0.0270	0.0250	0.0180	0.0530	0.0985	0.1090	0.0340	0.0360	0.0520	0.0960	0.0275	0.0470	0.0630	0.0860	0.0370

error

0.00000044444444444445	0.00000277777777777778	0.00001111111111111110	0.00001111111111111110	0.00002500000000000000	0.00002500000000000000	0.00002500000000000000	0.00002500000000000000	0.00002500000000000000	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00004444444444444450	0.00006944444444444400	0.00006944444444444400	0.00010000000000000000	0.00010000000000000000	0.00010000000000000000	0.00010000000000000000	0.00010000000000000000	0.00017777777777777800	0.00017777777777777800	0.00017777777777777800	0.00017777777777777800	0.00027777777777777600	0.00027777777777777600	0.00027777777777777600	0.00027777777777777600	0.00040000000000000100
1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000	1.36111111111111000000
0.00004391058208738580	0.00002322543184787350	0.00001533241399332300	0.00006253904623073680	0.00002237034331628930	0.00004067315701790890	0.00003495366143170180	0.00007112788503411250	0.00018976468138275800	0.00000145158949049209	0.00000787432319866708	0.00004120252008765530	0.00007971938775510170	0.00010939452311536800	0.00011324986303871500	0.00001275510204081640	0.00002939468718246480	0.00005225722165771530	0.00012808733564645700	0.00006132965597329090	0.00000730460189919631	0.00000145158949049209	0.00012420162578022900	0.00000416493127863391	0.00000326607635360721	0.00000937109537692630	0.00001306430541442830	0.00005524105186267360	0.00011113732036580000	0.00111863130021914000	0.00005033327246165090	0.00005856934610578940	0.00001778197125852810	0.00026296566837107400	0.00001181640625000000	0.00009290172739149390	0.00007112788503411250	0.00003149729279466830	0.00010412328196584800
0.00005081571264979620	0.00004799317752939470	0.00000907243431557555	0.00002621933517201340	0.00003906250000000010	0.00000580635796196835	0.00000326607635360719	0.00000444549281463199	0.00000486591695501722	0.00000326607635360720	0.00001306430541442880	0.00001778197125852810	0.00000326607635360721	0.00002322543184787350	0.00003628973726230230	0.00001306430541442880	0.00002939468718246480	0.00002041297721004500	0.00000734867179561621	0.00002621933517201340	0.00002939468718246480	0.00000036289737262302	0.00000976562500000056	0.00000145158949049209	0.00000326607635360720	0.00001778197125852810	0.00003906249999999980	0.00001306430541442880	0.00000486591695501731	0.00000216262975778547	0.00001306430541442880	0.00000104123281965848	0.00003906250000000010	0.00000416493127863391	0.00000527124114952103	0.00000216262975778547	0.00003906250000000010	0.00003906249999999900	0.00000216262975778547

resutat que j'obtiens

27132.315	483.1681	12.1576	-13.1959	-5.7039	0	0	-23.6157	40952.7182	1659.7886	28.337

les résultats que je devrais obtenir

27383,78	487,51	12,33	-13,14	-5,77	-23,62	41033,52	1660,56	28,4

Je n'utilise pas la meme version que celui qui a fait les calculs avant moi , est ce que ça pourrait venir de la ??


Merci

Hors ligne

#12 Le 25/10/2007, à 10:35

snigit

Re : aide programmation matlab

Je pense qu'aucune théorie mathématique n'ait changé entre les deux versions smile

Hors ligne

#13 Le 25/10/2007, à 10:41

laroche1

Re : aide programmation matlab

Ok

par contre quand je lui demande de ressortir les valeurs de data.txt à l'écran j'obtiens ceci

est ce que ça ne jouerait pas dans les calculs car là il arrondit??

M =

  1.0e+003 *

  Columns 1 through 13 

    0.0002    0.0005    0.0010    0.0010    0.0015    0.0015    0.0015    0.0015    0.0015    0.0020    0.0020    0.0020    0.0020
    1.4730    1.4730    1.3730    1.4730    1.2730    1.3730    1.4730    1.5730    1.6730    1.1730    1.2730    1.3730    1.4730
    0.0008    0.0008    0.0009    0.0008    0.0009    0.0009    0.0008    0.0008    0.0006    0.0009    0.0009    0.0009    0.0008
    0.0001    0.0001    0.0001    0.0001    0.0000    0.0000    0.0001    0.0001    0.0001    0.0000    0.0000    0.0000    0.0001

  Columns 14 through 26 

    0.0020    0.0020    0.0020    0.0020    0.0020    0.0020    0.0020    0.0020    0.0020    0.0020    0.0025    0.0025    0.0030
    1.5730    1.6730    1.2730    1.3730    1.4730    1.5730    1.6730    1.1730    1.2730    1.6730    1.1760    1.2730    1.2640
    0.0008    0.0006    0.0009    0.0009    0.0008    0.0008    0.0006    0.0009    0.0009    0.0006    0.0009    0.0009    0.0009
    0.0001    0.0001    0.0000    0.0000    0.0001    0.0001    0.0001    0.0000    0.0000    0.0001    0.0000    0.0000    0.0000

  Columns 27 through 39 

    0.0030    0.0030    0.0030    0.0030    0.0040    0.0040    0.0040    0.0040    0.0050    0.0050    0.0050    0.0050    0.0060
    1.1730    1.4730    1.7730    1.7730    1.3730    1.4730    1.5730    1.7730    1.3730    1.5730    1.6730    1.7730    1.5730
    0.0009    0.0009    0.0006    0.0006    0.0009    0.0009    0.0008    0.0005    0.0009    0.0008    0.0007    0.0006    0.0008
    0.0000    0.0001    0.0001    0.0001    0.0000    0.0000    0.0001    0.0001    0.0000    0.0000    0.0001    0.0001    0.0000

Hors ligne

#14 Le 25/10/2007, à 15:36

vinc-mai

Re : aide programmation matlab

Il n'arrondit que pour l'affichage.

Ton erreur

??? Error using ==> *
Inner matrix dimensions must agree.

Error in ==> W:\matlab\code.m
On line 77  ==>    X=X0+C0*G'*inv(G*C0*G')*(G*(X-X0)-g);

provient du fait que G est une matrice dmension 2*15; tu ne peux donc le multiplier avec C0 de taille 167².
Tu devrais faire un essai avec très peu de point de mesure pour voir d'où vient le problème. Créées un petit fichier data.txt pour tes tests.


ZiK un lecteur audio et son blog.

Hors ligne