Discussion :

[laroche1]

Bonjour, merci pour la réponse

En fait je me pose une question.
Est ce qu'il y a des personnes super fortes en algo??

voila j'ai un gros probleme j'ai recuperé ce code or la personne n'a pas indiqué son algo ni meme quelques commentaires donc j'essaye de comprendre ce qu'elle a pu faire .
Or à un endroit je doute , il y a une boucle sur h (for h = 1: 50 ) , et a la fin le vecteur resultat soit R est inclu dans la boucle. Mais n'écrase-t il pas à chaque fois les données??? de plus G' est inclu dedans la boucle h pour X puis on reutilise G' pour definir C hors de la boucle !!! je comprends plus car W c'est ma matrice finale ...donc ça devrait être pareil pour R ?? non???

Merci beaucoup d'avance !

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                                         %
%    Program for ideal orthopyroxene / symmetric clinopyroxene models     %
%                                                                         %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
%clear;
M=DLMREAD('D:/prog/inversion-idsym/data1.txt','\t');			
j=size(M);
m=DLMREAD('D:/prog/inversion-idsym/data1.txt');
m = [m(1,:) m(2,:) m(3,:) m(4,:)]';
X0=[100000 10000 1000 0 0 10 0 0 0 100000   10000 1000 m']';	
e=DLMREAD('D:/prog/inversion-idsym/error1.txt');	
e = [e(1,:) e(2,:) e(3,:) e(4,:)]';
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);	
end    
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('D:/prog/inversion-idsym/lui1.txt',Q,'\t')
%dlmwrite('D:/prog/inversion-idsym/rcomp-idsym1.txt',K,'\t')

[FR119492]

Salut !

Tu peux essayer de retrouver l'algorithme à partir du code, mais c'est souvent difficile, voire impossible si celui-ci n'est pas correctement documenté. Je pense que tu aurais avantage de reprendre ton problème à la base, en appliquant la méthodologie rigoureuse que j'ai toujours essayé d'inculquer à mes étudiants.

Tu connais:

1. la géométrie du dispositif que tu étudies;
2. les caractéristiques pertinentes des matières qui le constituent;
3. les lois physiques (éq. générales) qui régissent les phénomènes étudiés.

A partir de ces trois groupes de données, tu modélise ton dispositif, c'est-à-dire que tu formules mathématiquement (éq. diff. ordinaires, éq. aux dérivées partielles, recherche d'un extrémum) ton problème particulier. Tu obtiens alors un problème mathématique dont la solution correspondra au comportement de ton dispositif.

Quand tout cela est fait, tu choisis les méthodes numériques et les outils logiciels les mieux adaptés. Cela peut te paraître faire preuve d'un excès de pédanterie, mais, durant ma longue carrière, j'ai vu les pires horreurs avec des gens qui choisissaient la méthode, le langage ou le logiciel avant d'avoir analysé le problème à résoudre.

Cela dit, je reste à ta disposition pour t'aider à résoudre ton problème quanttu l'auras correctement formulé.

Bonne chance.
Jean-Marc Blanc

[To-Ma]

Bonjour,

je ne connais pas le langage de ton code ni le contexte dans le détail mais je donne mon humble avis, on ne sait jamais :

Pour moi la boucle for h... est une résolution d'un problème auto-cohérent pour déterminer G (donc G') et X, et si R et G' sont en dehors de la boucle ça ne change rien...

...A moins que pour ce langage la boucle puisse casser en court d'exécution et passer au calcul de C . Dans ce cas, déterminer R et G dans la boucle permet de faire un calcul avec des paramètres cohérents (c.a.d tous les X, F, P, etc... utilisés ont été calculés dans la même itération)....