Pages : 1
#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.beTo 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.032Insert 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.032Prior 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-05Insert 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-05How 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:50Save 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é
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 "*".
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
#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
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.
Hors ligne
Pages : 1