[Matlab Debutant] Boucle For

[zenip]

Bonjour à tous,

Pour un projet d'école je suis en train de programmer une petite simulation sur Matlab. J'ai fais un premier code qui marche bien et qui permet de résoudre mon problème avec des conditions initiales "fixes".

J'aimerai maintenant calculer les nouvelles valeurs à partir des valeurs précédentes grâce a une boucle for mais j'ai de la peine à écrire un code pour cela.

Mon code est le suivant :

Ta = [8 9 10 10 10 10 11 12 14 16 18 25 22 21 19 28 17 26 15 14 12 11 10 10];                             
v_wind = 1.5;                          
Aire = 1.58;                           
G = [0 0 0 0 0 0 30 100 400 600 700 1000 1100 1100 1050 900 850 700 500 300 0 0 0 0];    
debit=100;
Acell = 0.156*0.156;                   
Ncell = 60;                            
Tcell = ones(1,24)*(25+273);           
Amod = 1.66;                          
PF = (Acell*Ncell)/Amod;               
rend_ech = 0.7;                       
rend_optique = 0.578;
e_g = 0.003;                       
lambda_g = 0.93;                    
tau_g = 0.92;                      
epsilon_g = 0.88;                     
alpha_g = 0.05;                       
Cw=4180;
e_pv = 3e-4;                  
lambda_pv = 150;                       
alpha_pv = 0.75;                       
tau_pv = 0.95;                        
e_ab = 0.005;                        
lambda_ab = 25;                       
epsilon_ab = 0.5;                    
e_bs = 5e-4;
lambda_bs = 0.033;
eff_ab = rend_optique/(alpha_pv*tau_pv);
 
Ta_K = Ta + 273.15;
Tsol = Ta_K + 2;
Tin = 20+273.15;
Tsky = 0.0552*(Ta_K).^ 1.5;
Tg0 = Ta_K + 2;
Tpv0 = Ta_K + 2;
Tab0 = Ta_K + 2;
Tmoy = (Tpv0+Tg0)/2;
 
hc_n = 5.67+3.86*v_wind;              
c_boltz = 5.67e-8;                   
hrad_sky = c_boltz*epsilon_g.*[(Tg0.^2 - Tsky.^2) .* (Tg0.^2 + Tsky.^2)./(Tg0-Ta_K)]; 
hc_g = lambda_g/e_g;                         
hc_air = 2.4;                         
hc_eva = ones(1,24)*30;          
hrad_g_pv = 4*c_boltz*E*Tmoy.^3;       
hw = 232.3;  
hloss = 200;
 
E = (1/e_pv + 1/e_g -1)^-1;
 
%% Connues matrice B :
b1 = G*Aire*alpha_g + Ta_K.*Aire.*hc_n + Tsky.*Aire.*hrad_sky;
b2 = G*Aire*alpha_pv*tau_g - (G*tau_g*PF*Aire*0.8); 
b3 = 0*G ;
b4 = Tin.*(debit.*Cw.*rend_ech + hw*Aire) + Ta_K.*hloss.*Aire;
 
B = [b1 ; b2 ; b3 ; b4];

J'aimerai donc par exemple mettre les températures initiales Tg0 = 15°C, lancé le calcul, avoir un résultat et utiliser ce résultat pour la prochain itération et faire cela pour toutes mes températures et paramètres.

Donc exemple :

Je fixe au temps t=0 :

Tg0 = 15;
Tsky = 12;
Ta_K = 10;

Je calcule mon hrad_sky = c_boltz*epsilon_g.*[(Tg0.^2 - Tsky.^2) .* (Tg0.^2 + Tsky.^2)./(Tg0-Ta_K)]; au temps 0.

Je résous mes équations, je trouve mes Tg, Tsky, Ta puis je remets ces nouvelles valeurs dans ma boucle afin d'obtenir un nouveau hrad_sky..

Pour cela je pensais initialiser ma matrice b à zéro, puis faire une boucle for mais impossible d'écrire quelques chose qui fonctionne..

J'espère avoir été clair, sinon je me re explique volontiers.

Merci d'avance.
-----

[Antoane]

Bonjour,

j’admets ne pas avoir lu en détail ton code, mais si j'ai bien compris, tu cherches à résoudre pas à pas ton équation, le résultat de l'étape k+1 étant fonction des résultats de l'étape k.
Autrement dit, tu cherches à calculer, pas à pas, les termes d'une suite définie par récurrence.
Pour uns suite définie par : T_k+1 = lambda * T_k⁴, T₀=22 ;
Le code pourra donc être du genre :

T(1) = 22 ;
lambda = 8.623 ;
for k=1:beaucoup
        T(k+1) = lambda*T(k).^4 ;
end

figure, 
plot(T) ;

Sachant que T_k+1 peut être un vecteur, à condition que T soit une matrice dans ton code.
Sachant que T_k+1 peut être un matrice, à condition que T soit une matrice 3D dans ton code.

Sachant que lambda peut être fonction de k.

[zenip]

J'ai bien compris ton code, mais cela ne marche pas quand je veux appliquer cela au mien..

Par exemple :

Ta_K = [8 9 10 10 10 10 11 12 14 16 18 25 22 21 19 28 17 26 15 14 12 11 10 10];    
Tsky(1) = 0.0552*(Ta_K).^ 1.5;
Tg0(1) = 288;
for i = 1:24
    
    hrad_sky(i+1) = c_boltz*epsilon_g.*[(Tg0(i).^2 - Tsky(i).^2) .* (Tg0(i).^2 + Tsky(i).^2)./(Tg0(i)-Ta_K)];    % [W/m2.K]
    
end

J'ai l'erreur suivante : In an assignment A(I) = B, the number of elements in B and I must be the same.

Je pense que que cela vient du fait que mon Ta est défini comme un vecteur 1x24 ?

[Antoane]

Matlab renvoie des messages d'erreur assez explicites.
En l'occurrence, il te dit que tu essayes de mettre un truc d'une certaine taille dans un truc d'une autre taille. Forcément ca marche pas : ou bien le contenant est trop petit et ca rentre pas, ou bien il est trop grand et on ne sait pas où mettre le truc dedans, ou bien les deux ( trop grand pour une dimension, trop petit pour une autre).
Trouve quelle ligne pose problème et pourquoi, tu pourras ensuite avancer.

[A voir en vidéo sur Futura]

[zenip]

Salut ! Alors j'ai bien avancé dans mon code, la boucle tourne mais j'ai encore une question.

Pour résumer :

- Je pose mes conditions initiales
- Je calcul mes divers coefficient
- Je résous un système linéaire grâce a des matrices/vecteurs AX=B
- J'obtiens de nouvelles températures que je remplace à la place de mes conditions initiales pour re calculer mes coefficient
- J'itère jusqu’à ce que j'obtiens mes 24 itérations

Je n'arrive pas à bien écrire mes matrices qui évoluent à chaque itération, j'ai essayé le bout de code ci dessous mais cela ne semble pas fonctionner..

b1(i+1) = G(i)*Aire*alpha_g + Ta_K(i).*Aire.*hc_n + Tsky(i).*Aire.*hrad_sky(i);
b2(i+1) = G(i)*Aire*alpha_pv*tau_g - (G(i)*tau_g*PF*Aire*0.8); 
b3(i+1) = 0*G(i) ;
b4(i+1) = Tin.*(debit.*Cw.*rend_ech + hw*Aire) + Ta_K(i).*hloss(i).*Aire;

% B = [b1(i) ; b2(i) ; b3(i) ; b4(i)];

B = [b1 ; b2 ; b3 ; b4];

[Antoane]

Bonjour,
> cela ne semble pas fonctionner..
Peux-tu détailler ?

Avec ça :

b1(i+1) = G(i)*Aire*alpha_g + Ta_K(i).*Aire.*hc_n + Tsky(i).*Aire.*hrad_sky(i);
b2(i+1) = G(i)*Aire*alpha_pv*tau_g - (G(i)*tau_g*PF*Aire*0.8); 
b3(i+1) = 0*G(i) ;
b4(i+1) = Tin.*(debit.*Cw.*rend_ech + hw*Aire) + Ta_K(i).*hloss(i).*Aire;

% B = [b1(i) ; b2(i) ; b3(i) ; b4(i)];

B = [b1 ; b2 ; b3 ; b4];

si b1, b2 b3 et b4 ne sont pas initilisées, elles seront construites peu à peu au cour de l'exécution du programme et ce seront des vecteurs ligne.
Le résultat sera que B sera une matrice de 4 lignes et ayant 1 colonne lors du premier passage dans la boucle (i.e. @i=1), puis 2 colonnes lors du second passage dans la boucle (i.e. @i=2), etc.

[zenip]

J'obtiens des valeurs complètement absurde donc je me demande si c'est l'écriture du code au niveau des matrices ou une erreur dans une de mes équations.. Je vais repasser cela en revu !

Pour la matrice B je l'ai initialisé comme les autres variables.

%% Matrice B
B = zeros(4,24);
b1(1) = G(1)*Aire*alpha_g + Ta_K(1).*Aire.*hc_n + Tsky(1).*Aire.*hrad_sky(1);
b2(1) = G(1)*Aire*alpha_pv*tau_g - (G(1)*tau_g*PF*Aire*0.8); 
b3(1) = 0*G(1) ;
b4(1) = Tin.*(debit.*Cw.*rend_ech + hw*Aire) + Ta_K(1).*hloss(1).*Aire;

Ensuite je fais la même chose pour la matrice A

Puis je calcul mes inconnues X=A\B

Puis je récupère mes valeurs dans la matrice X Tg0(i+1) = X(3,i+1); pour recalculer les températures suivantes.

Donc pour toi mon code est correcte ?

Merci beaucoup pour ton aide !

[zenip]

"cela ne semble pas fonctionner.." mes matrices n'arrive pas à évoluer dans la boucle. En gros elle sont bloqué à la première itération..

[zenip]

J'ai compris mon erreur mais je n'ai aucune idée de comment la résoudre.

Je m'explique :

à i=1 (première itération) j'obtiens le système suivant :

A=[4x8]
B=[4x2]
X=[8x2]

Les températures qui m'intéresse dans la matrice X se trouve sur les lignes 1,3,5,7 je les récupère et je calcul ma boucle à l'itération suivante :

i=2 j'obtiens :

A=[4x12]
B=[4x3]
X=[12x3]

Mes températures se sont décaler de ligne et sont maintenant aux lignes 3,4,7,10.

Et ainsi de suite, mes températures vont se décaler de ligne à chaque itération, comment je peux récupérer ces valeurs aux bonnes lignes correspondantes ?

Car pour l'instant je récupère T1 = X(1,i+1) mais cela marche que pour i=1, pour i=2 je dois récupéré T1 = X(3,2) car si je laisse mon T1 = X(1,i+1) il va récupérer un zéro et fausser le calcul..

[zenip]

Je pense que je fais une erreur d'écriture de mes matrice dans ma boucle

b1(i+1) = G(i)*Aire*alpha_g + Ta_K(i).*Aire.*hc_n + Tsky(i).*Aire.*hrad_sky(i);
b2(i+1) = G(i)*Aire*alpha_pv*tau_g - (G(i)*tau_g*PF*Aire*0.8); 
b3(i+1) = 0*G(i) ;
b4(i+1) = Tin.*(debit.*Cw.*rend_ech + hw*Aire) + Ta_K(i).*hloss(i).*Aire;

% B = [b1(i) ; b2(i) ; b3(i) ; b4(i)];
% B = [b1 ; b2 ; b3 ; b4];

Je n'arrive pas à écrire ma matrice à i=1 B(1) = [.... ..... ..... ...]

Mon erreur vient peut être d'ici..

[invite73192618]

Perso je te suggérerais de décomposer ton problème en plusieurs fonctions simples. Utiliser une structure pour les paramètres peut aider. Un squelette pourrait ressembler à ceci:

function [resultat]=myfonction % mélange français/anglais juste pour faire rager les puristes ^^
% calcule l'évolution d'une température en fonction de conditions initiales
 
	for n=1:maxn % de 1 à maxn paramètres
		[param]=initialisation(n);
		for t=1:beaucoup % de 1 à beaucoup itérations
			[Ts]=suivant(T0, param);
			T0=Ts;
			resultat(n,t)=TS;
		end
	end
end

function [param]=initialisation(n)
	switch n
	case 1
		param.lambda_g = 0.93;                    
		param.tau_g = 0.92;                      
		param.epsilon_g = 0.88;                     
		param.alpha_g = 0.05;                       
		param.Cw = 4180;
		% etc
	case 2
		param.lambda_g = 0.95;                    
		param.tau_g = 0.97;                      
		param.epsilon_g = 0.85;                     
		param.alpha_g = 0.06;                       
		param.Cw = 3180;
		% etc
	% etc 
	otherwise
		disp('set de paramètres non spécifiés dans la fonction initialisation')
	end
end

function [Ts]=suivant(T0, param)
	b1 = param.G*param.Aire*param.alpha_g + ...
end

[zenip]

Merci de ta réponse mais vu mon niveau en matlab je t'avoue que je n'ai pas compris grand chose a ton squelette... Tu n'as pas de solution pour ma question de dessus ?

Peut etre que si tu voyais mon code dans son ensemble ce serait plus claire.

[invite73192618]

L'idée n'est pas de te fournir une solution toute-faite, mais de t'aider à la trouver toi-même. Si mon squelette est encore trop compliqué, ma suggestion est de le simplifier encore plus.

Pour reprendre ton dernier exemple, si tu as un calcul qui "marche pour i=1 mais fait un décalage nuisible pour i=2", alors tu peux faire une fonction comme suit:

function moncalcul(data, nbetapes)

	for i=1:nbetapes
		[resultat]=moncalculpouriegal1(data);
		[data]=boucle(resultat);
	end

end

function [resultat]=moncalculpouriegal1(data)
	% mettre ici le calcul qui marche pour i=1
end

function [data]=boucle(resultat)
	% mettre ici un calcul qui compense le décalage
end

[zenip]

Bon j'ai réussi à faire ce que je voulais sans boucle finalement car ca me compliquait la vie.. Le code est très moche mais au moins ca fonctionne et je comprends ce que je fais !