Tracé d'un abaque-formule de Colebrook [Résolu]

[clement2406]

Bonjour à tous,

Tout d'abord, veuillez m'excuser si je ne poste pas dans la zone approprié, mais j'arrive tout juste sur le forum.

Voila mon problème, je souhaite tracer un abaque à partir de la formule de Colebrooke en mécanique du fluide. Mais le programme de tourne pas, et je dois avoir terminé ce programme avant la fin des vacances, seulement je bloque complétement... J'implore donc votre aide.

Voici mon programme actuel :

Code :

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
%----------------------------------------------------
% abaque lambda en fonction du nombre de Reynolds R
% et du rapport Epsilon/D
%----------------------------------------------------

clear all;
echo on;
%---------------------------------------------------
% abaque lambda en fonction du nombre de Reynolds R
% et du rapport Epsilon/D
%---------------------------------------------------
echo off;

Rs=logspace(1,8,10);
K=10^(-5);
for j=1:length(Rs)
    R=Rs(j);
    for i=1:10:10^4
        x=iter3(15,i*K,R,0.001);
        plot(R,x);
        text('String',['\leftarrow' num2str(i*K) ],...
            'Position',[10.0,x(find(R==10))],...
            'HorizontalAlignment','left'); 	% affichage flêche gauche et E/K
        hold on;
    end
end
text(1000,5,'Epsilon/D varie de 10^{-5} à 10^{-1}');
xlabel('Re Reynolds');
ylabel('Lambda');
title('Abaque Lambda en fonction de Re et Epsilon/D');
echo on;
% cliquer un point pour avoir ses coordonnées
echo off;
[xc,yc] = ginput(1);  % attend 1pt de cliquage
text('String',['Re=' num2str(xc,3) ', Lambda=' num2str(yc,4) ],...
    'Position',[60.0,3200],...
    'HorizontalAlignment','left'); %

grid;

avec iter3 qui est le programme suivant :

Code :

1
2
3
4
5
6
7
8
9
10
function [L] = iter3(x,K,R,ep)
if (R>=10 & R<=10^8) & (K>=10^(-5) & K<=10^(-1))
    L=0;
    while x-L>=ep
    L=(-2*log10((K/(3.71)))+(2.51/(R*x^(-1/2))))^(-2);
    x=L; 
    end
else disp('The Colebrook equation is valid for Reynolds numbers 10<=R<=10^(8) & K must be 10^(-5)<=K<=10^(-1)');
end
end

Un très grand merci aux personnes qui m'aideront et me permettront de corriger mon problème.

Cordialement,

Clement

[Dut]

Remplace

Code :

1
2
3
        text('String',['\leftarrow' num2str(i*K) ],...
            'Position',[10.0,x(find(R==10))],...
            'HorizontalAlignment','left'); 	% affichage flêche gauche et E/K

par

Code :

        text(10.0,x(find(R==10)),['\leftarrow' num2str(i*K) ],'HorizontalAlignment','left')

Ajoute aussi drawnow juste après hold on.

Sinon pour le reste, vu le nombre d'itérations à effectuer, je doute que ton code sois très rapide

[clement2406]

Merci beaucoup pour tas réponse Dut, je vais tester cela tout de suite.

Oui j'en conviens, mon programme est lent, mais pour l'instant je veux avoir un programme fonctionnel avant de l'optimiser au mieux. Je ne suis pas encore très bon en programmation, je préfère procéder par étapes...

Si d'autres personnes décèlent d'autres erreurs je vous serais reconnaissant de m'en faire part. Merci encore Dut.

Bonne journée.

[Dut]

Dans la fonction iter3, si la première condition n'est pas remplie, la variable de sortie L n'existera pas, donc la variable x dans le code principale sera vide et cela fera planté l'exécution.

As-tu réellement besoin de tracer l'abaque en temps réel ?

Tu devrais d'abord faire le calcul de toutes les valeurs en les stockant dans une matrice puis faire l'affichage une seule fois, non ?

[clement2406]

Non effectivement, je n'ai pas besoin de tracer l'abaque en temps réel.

Le stockage des valeurs dans une matrices était mon idée initial, mais j'ai rencontré des erreurs plutôt nombreuses lorsque j'ai essayé de faire ça, j'ai donc abandonné pour cette méthode plus "intuitive" pour ma part.

J'avais effectivement décelé le problème avec L et x dont tu parles, mais je n'arrive pas du tout à le résoudre, d'où mon poste sur ce forum.

A vrai dire, on nous a "lâcher" dans ce projet sous Matlab en ayant à peine les bases sous ce logiciel. Travaillant sous Maple avant, ce logiciel travaillant avec une méthode matricielle me trouble un peu.

Quel est le programme pour stocker chacune de mes valeurs dans une matrice ? Mais surtout que le stockage soit tel que le traçage de l'abaque ce fasse aisément par le logiciel, il me semble que mon problème était sur ce point.

Je te remercie de ton aide Dut.

[clement2406]

J'ai lu attentivement ton tuto Dut concernant les matrices sous Matlab, ce qui ma permis de mieux cerner la mise en place de matrice, mais je ne parviens toujours pas à saisir comment je pourrais faire.

Car il faudrait tout d'abord mettre chacune de mes valeurs de lambda (de iter3 donc) dans une matrice (ce que je ne parviens déjà pas à faire..), seulement lambda dépend de 2 variables (R le nombre de Reynolds et K).

Dans mon programme initial, pour faciliter la lecture de l'abaque, j'ai pensé me contenter de 10 valeurs de K environ (un pas pour K tel que il varie 10 fois dans intervalle donné, comme j'ai voulu faire dans mon programme initial), et faire varier R d'un pas de 10 ou moins.

Pour cette méthode, je vais donc avoir 10 matrices (une pour chaque valeurs de K) et je vais devoir tracer sur un même graphique chaque courbes correspondant à chaque valeurs de K. C'est bien cela ?

Sauf que je n'arrive pas la première étapes... Je ne vois pas comment faire varier R n-fois et K 10 fois pour obtenir les 10 matrices; et encore moins les utiliser par la suite.

Quelqu'un pourrait-il m'aider ?

Je vois très bien la procédure à suivre, mais ne pas pouvoir l'appliquer car je ne manipule pas assez bien le logiciel est quelque chose de très frustrant...

Cordialement,

Clément

[Dut]

Sort de ton problème réel et prends un exemple plus simple :

Code :

1
2
3
4
5
6
7
for i = 1:3
    k = 0;
    for j = 0:0.1:0.5
        k = k+1;
        x(i,k) = i*j;        
    end
end

[clement2406]

C'est effectivement ce que j'ai essayé de faire, avec une boucle while.

Voila où j'en suis :

Code :

1
2
3
4
5
6
7
8
9
10
11
function [M]=matrice1(K)

R=10;
j=10;
i=1;
while j*R<10^8  
    x=iter3(150,K,j*R,0.001);
    M[i]=x
    i=i+1;
    j=j*10;
end

Mais ça ne marche pas. Je dois surement rajouter une boucle mais je ne vois pas bien où, vu que ma condition d’arrêt ne porte que sur R.

[Dut]

Je ne comprends pas pourquoi la fonction iter3 devrait renvoyer autre chose qu'un scalaire.

Pourquoi ne pas simplement faire quelque chose comme ceci dans ton script principal ?

Code :

1
2
3
4
5
6
7
for i = 1:3
    k = 0;
    for j = 0:0.1:0.5
        k = k+1;
        x(i,k) = iter3(...);        
    end
end

[Winjerome]

Un retour sur le tutoriel :
III-A. Indexage classique (ligne,colonne)

La méthode d'indexage classique consiste à spécifier la position d'un élément en fonction de l'indice de la ligne et de l'indice de la colonne où il se trouve dans la matrice en prenant comme premier élément, celui situé en haut à gauche. L'indexage s'effectue entre parenthèses avec en premier l'indice de la ligne et en second, l'indice de la colonne, soit : M(idx ligne, idx colonne).

Je te suggère aussi cette lecture.

[FLB]

Je te propose ce bout de code, qui n'optimise certes pas le calcul mais qui évite les 1000 appels à la fonction plot.

Code :

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
%Paramètres
Rs=logspace(1,8,10);
K=10^(-5);
resultats = zeros(2,length(Rs)*length(1:10:10^4));%Préallocation mémoire

%Calculs
for j=1:length(Rs)
    R=Rs(j);
    for i=1:10:10^4
        resultats(1,length(1:10:10^4)*(j-1)+i)=R;
        resultats(2,length(1:10:10^4)*(j-1)+i)=iter3(15,i*K,R,0.001);
    end
end

%Affichage
plot(resultats(1,:),resultats(2,:),'+');

Toutefois une de tes lignes de code me laisse songeur :

Code :

x(find(R==10))

renverra la première valeur de x (mais x est un tableau d'une seule valeur!) lors du premier passage de la boucle j, puis renverra [] au second passage, ce qui fera planter la ligne

Code :

 text(10.0,x(find(R==10)),['\leftarrow' num2str(i*K) ],'HorizontalAlignment','left');

Si tu veux marquer chaque point avec une flèche, je te propose:

Code :

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
%Paramètres
Rs=logspace(1,8,10);
K=10^(-5);
resultats = zeros(3,length(Rs)*length(1:10:10^4));%Préallocation mémoire

%Calculs
for j=1:length(Rs)
    R=Rs(j);
    for i=1:10:10^4
        resultats(1,length(1:10:10^4)*(j-1)+i)=R;
        resultats(2,length(1:10:10^4)*(j-1)+i)=iter3(15,i*K,R,0.001);
        resultats(3,length(1:10:10^4)*(j-1)+i)=i*K;
    end
end

%Affichage
plot(resultats(1,:),resultats(2,:),'+');
hold on
 %Transformation de la legende numerique en caractères
legendes = cellfun(@num2str,resultats(3,:),'UniformOutput',false);
 %Ajout de la flèche
legendes = cellfun(@strcat, repmat({'\leftarrow'},1,length(resultats(3,:))), legendes,'UniformOutput',false);
 %tracé
text(resultats(1,:),resultats(2,:),legendes,'HorizontalAlignment','left')

[clement2406]

Dut : Effectivement Dut, tu as raison... Merci pour ton script suggéré. Seulement maintenant, j'ai une question : quel format pour ma matrice ? Devrait-je mettre mes résultats dans une matrice colonne, carré ou quelconque pour tracé la courbe affiliée ? Je comptais mettre mes résultats dans une matrice colonne unique, est-ce une erreur ?

EDIT : Ce programme fonctionne dorénavant. C'est maintenant à partir de la matrice obtenue que je n'arrive pas à tracer ma courbe correct.

Code :

1
2
3
4
5
6
7
8
9
function [M]=matricetest(K)
R=1;
j=10;
i=1;
while j*R<10^8  
    M(i)=iter3(150,K,j*R,0.001);
    i=i+1;
    j=j*2;
end

Je veux faire :

Code :

plot(M,log(R))

Mais j'ai du mal à définir R car il doit être de même dimension que M.

winjerome : Merci, j'avais effectivement omis cela, c'est une bête erreur de ma part. Comme vous le voyez, je débute encore dans Matlab, veuillez donc m'excusez pour les erreurs grossières que je fais...


FLB : Merci beaucoup. Je vais essayer de me pencher sur le code que m'as fournis. Même si j’avoue ne pas tout saisir étant encore inexpérimenté...


EDIT : J’obtiens bien des tracés correct en rajoutant l’échelle logarithmique, merci beaucoup FLB. Seulement, comment faire pour relier chacun des points construit et affiché les flèches ? Car je ne comprend pas, les commandes m'ont l'air correct, mais l'affichage des flèches et des légendes ne fonctionne pas..

Code :

1
2
3
4
5
6
7
8
??? Error using ==> cellfun
Input #2 expected to be a cell array, was double
instead.

Error in ==> abaque_correction2 at 20
legendes =
cellfun(@num2str,resultats(3,:),'UniformOutput',false); 
>>

J'ai d'ailleurs oublié de le préciser, mais mon objectif final est d'obtenir ceci sous Matlab :

http://www.ac-nancy-metz.fr/enseign/.../colebrook.gif

[FLB]

il y a 2 coquilles dans mon code, il faut corriger par :

Code :

1
2
3
4
5
6
7
8
9
10
simID=1;
for j=1:length(Rs)
    R=Rs(j);
    for i=1:10:10^4
        resultats(1,simID)=R;
        resultats(2,simID)=iter3(15,i*K,R,0.001);
        simID = simID +1;
    end
end

Sous peine d'enregistrer de mauvais résultats
Et j'ai oublié un détail en copiant collant mon code :

Code :

legendes = cellfun(@num2str,num2cell(resultats(3,:)),'UniformOutput',false);

Mais je doute que tu veuilles afficher toutes les flèches de tous les points?

[clement2406]

D'accord, je ne les avais pas vue...

Non, une seule flèche par courbe pour indiquer le K différent sur chaque courbes.

C'est ce que faisait la partie suivante dans mon programme initial :

Code :

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
%----------------------------------------------------
% abaque lambda en fonction du nombre de Reynolds R
% et du rapport Epsilon/D
%----------------------------------------------------

clear all;
echo on;
%---------------------------------------------------
% abaque lambda en fonction du nombre de Reynolds R
% et du rapport Epsilon/D
%---------------------------------------------------
echo off;

Rs=logspace(1,8,10);
K=10^(-5);
for j=1:length(Rs)
    R=Rs(j);
    for i=1:10:10^4
        x=iter3(15,i*K,R,0.001);
        plot(R,x);
        text('String',['\leftarrow' num2str(i*K) ],...
            'Position',[10.0,x(find(R==10))],...
            'HorizontalAlignment','left'); 	% affichage flêche gauche et E/K
        hold on;
    end
end
text(1000,5,'Epsilon/D varie de 10^{-5} à 10^{-1}');
xlabel('Re Reynolds');
ylabel('Lambda');
title('Abaque Lambda en fonction de Re et Epsilon/D');
echo on;
% cliquer un point pour avoir ses coordonnées
echo off;
[xc,yc] = ginput(1);  % attend 1pt de cliquage
text('String',['Re=' num2str(xc,3) ', Lambda=' num2str(yc,4) ],...
    'Position',[60.0,3200],...
    'HorizontalAlignment','left'); %

grid;

Mais j'ai une erreur lorsque j’essaie de rajouter une commande similaire ici..


Et quelle est la commande pour relier les points entre eux ?

Merci.

[FLB]

Tel que tu l'as programmé pour l'instant tu traces 10000 flèches...
du coup tu vas vouloir plutot marquer un truc comme :

Code :

text(resultats(1,indices_des_resultats_où_R_vaut_10^8_cest_a_dire_999_1999_etc),10^8,legendes,'HorizontalAlignment','left');

Tout en t'assurant que tu remplies la variable "legendes" à partir des mêmes indices.

Pour le dessin des courbes, regarde la fonction

Code :

doc plot

En particulier les linestyles, qui permettent de choisir de tracer des points '.', des lignes '-' (par défaut), des lignes en pointillées etc.
Il te faudra ainsi faire un plot pour chaque courbe, avec un "hold on" entre chaque.

[clement2406]

Mais comment faire un plot pour chaque courbe vu qu'elles sont calculées au fur et a mesure ?

J'ai essayé de mettre un :

Code :

plot(resultats(1,:),resultats(2,:),'--r','MarkerSize',10)

Mais il relit les points de courbes différentes, résultat j'ai quelque chose de complétement faux...

[Winjerome]

Bonjour,

Voir Comment relier des points dans une boucle ?

[clement2406]

Merci Winjerome

Je n'arrive à appliquer cela à mon programme actuel... Tant pis je vais rester avec des points je pense.

[FLB]

Tous les calculs sont fait, il faut que tu récupères chaque série de point qui composent une courbe et que tu la traces, par exemple :

Code :

1
2
3
4
5
figure;
for i=0:999
    plot(resultats(1,i+[1:1000:end]),resultats(2,i+[1:1000:end]));
    hold on
end

Ce code dessine 1000 courbes sur une meme figure, en utilisant pour chaque courbe le résultat d'un calcul sur 1000 (une courbe est faite de 10 points).

[clement2406]

J'aurais une autre petite question : comment fait on pour adapter l'échelle ?

Voici mon programme pour l'instant :

Code :

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
%Paramètres
Rs=logspace(1,8,10);
K=10^(-5);
resultats = zeros(3,length(Rs)*length(1:10:10^3));   %Préallocation mémoire

%Calculs
simID=1;
for j=1:length(Rs)
    R=Rs(j);
    for i=1:10:10^2
        resultats(1,simID)=log10(R);
        resultats(2,simID)=iter3(0.06,i*K,R,0.001);
        simID = simID +1;
    end
end

%Affichage
plotyy(resultats(1,:),K,resultats(2,:),'*');
hold on
matricelam();
R2=10:60:2000;
T2=matricelam();
plot(log10(R2),T2);

hold on
text(40,2600,'Epsilon/D varie de 10^-^5 à 10^-^1'); 
xlabel('Re Reynolds');
ylabel('Lambda');
title('Abaque Lambda en fonction de Re et Epsilon/D');
echo on;
% cliquer un point pour avoir ses coordonnées
echo off;
[xc,yc] = ginput(1);  % attend 1pt de cliquage
text('String',['Re=' num2str(xc,3) ', Lambda=' num2str(yc,4) ],...
      'Position',[4,0.003],...
      'HorizontalAlignment','left'); % 

grid;

J'ai rajouté la partie laminaire de la courbe, donc la droite à gauche du graphique sur la photo de l'abaque (http://www.ac-nancy-metz.fr/enseign/.../colebrook.gif)

Sauf que comme vous le voyez, les courbes calculées avant sont dans un intervalle proche de 0, tandis que la droite ,elle, est dans un intervalle beaucoup plus grand. Je voudrais juste avoir le petit bout de courbe de la partie proche de 0 d'affichée. Je ne trouve pas la commande adéquate.