Comment tracer une courbe sur Matlab

Résolu/Fermé
hal - 3 mai 2016 à 13:28
 hal - 12 oct. 2016 à 11:26
Bonjour,
s'il vous plait j'ai besoin de tracer une courbe de "D en fonction de h" comme ecrit ci dessous:

B=1.74;
C=10;
A=5.6*C^(-24);
H=0.001;
q=H*A;
for h=0:-0.5:-10
D=(1/B)*log((A+q)/(A*exp(B*h)+q))
end
plot(h,D,'r')

mais toujours j'obtient une fenetre avec les deux axes la courbe n'apparait pas


1 réponse

bjour Messages postés 6071 Date d'inscription lundi 11 janvier 2010 Statut Contributeur Dernière intervention 2 octobre 2024 4 074
9 mai 2016 à 21:35
Bonsoir,

Le problème vient de la définition dans la boucle.
D est une variable 1x1, à chaque itération tu lui assignes une nouvelle valeur.

Il faut donc ajouter un élément qui spécifie la position de ton calcul 1/B*... dans la variable D, qui sera 1x21 au final. Il est souvent préférable de créer un compteur à part pour les boucles for, je l'appelle i ici:
B=1.74; 
C=10;
A=5.6*C^(-24);
H=0.001;
q=H*A;

h=0:-0.5:-10 #liste des abscisses
for i=1:21 # 21 valeurs à calculer
D(i)=(1/B)*log((A+q)/(A*exp(B*h(i))+q)) # attention: mettre D(i) et h(i)
end;
plot(h,D,'r')

Et le tour est joué !
1
Je te suis très reconnaissant pour le temps que tu a pris pour m'aider à résoudre ce problème. Je n'aurais pas avancé aussi rapidement si tu n'avais pas été là pour m'aider, "merci bjour" .
0
bonjour
s'il vous plait est ce que vous pouvez m'aidez à résoudre ce problème car à chaque fois j'obtien une erreur à la ligne 184 et merci

clear all
%% PARAMETRES DE TEMPS
T=20;
n=600;
dt = T/n;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% PARAMETRES D'ESPACE
nx = 100;
nz = 80;
L = 10;%m
Zs = 5;%m
dx =L/nx;
dz =Zs/nz;
x = dx/2 : dx : L;
z = dz/2 : dz : Zs;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PARAMETRES HYDRODYNAMIQUES ET HYDRODISPERSIFS
Ks = 0.275;
b = 14.75;
or = 0.05;
os = 0.46;
qrel = -0.01;
qin = Ks*qrel;
Co = 3;
cin = 1;
D = 0.02;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% CARACTERISTIQUE DE LA NAPPE SUPERFICIELLE DRAINEE
zomax = L*sqrt(abs(qrel));
hgama=(1/b)*log(abs(qrel));
gama=-(1/b)*log((Ks*exp(b*hgama)+abs(qin))/(Ks+abs(qin)))
A=((Ks+abs(qin))/(b*qin))*(1-exp(-b*gama))
og=os*exp(b*hgama); %
I=-sqrt(abs(qrel)./((L^2)-x.^2))
R=I.*A;
zo=sqrt((abs(qrel).*(L^2-x.^2)))
zgama =zo+gama;
I=-sqrt(abs(qrel)./((L^2)-x.^2))
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% MATRICE DES VITESSES EN REGIME UNIFORME
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for j=1:nz
for i=1:nx
xi=x(i);
zj=z(j);
zo(i)=sqrt((abs(qrel).*(L^2-xi.^2)));
zgama(i) =zo(i)+gama;
I(i)=-sqrt(abs(qrel)./((L^2)-(xi.^2)));
I(nx)=I(nx-1);
R(i)=I(i).*A;
zoj (i)=ceil(zo(i)/dz);
if zj < zo(i)
qx(j,i)= sqrt((abs(qin)*Ks)./(((L^2)-(xi.^2)))).*xi-R(i)*qin;
qz(j,i)=0;
Vx(j,i)=qx(j,i)./os;
Vz(j,i)=0;
teta(j,i)=os;

elseif z(j) >= zo (i) & z(j) < zgama(i)
qx(j,i)=-R(i)*qin;
qz(j,i)=(1-R(i))*qin;
h(j,i)=-0.01.*(1/b)*log((1/Ks)*((Ks+qin)*exp(-b*(z(j)- zoj(i)))-qin));
teta(j,i)=os*exp(h(j,i));
Vx(j,i)=qx(j,i)./teta(j,i);
Vz(j,i)=qz(j,i)./teta(j,i);
elseif z(j) >= zgama(i)
qx(j,i)=0;
qz(j,i)=qin;
h(j,i)=hgama;
teta(j,i)=og;
Vx(j,i)=0;
Vz(j,i)=qz(j,i)./og;
end
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% BOUCLE DE TEMPS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for j=1:nz
for i=1:nx
x(i) = i*dx;
z(j)=j*dz;
zo(i)=sqrt((abs(qrel).*(L^2-x(i).^2)));
if z(j) < zo(i)
C(j,i,1)=0;
else
C(j,i,1)=cin;
end
end
end
for k=2:(n/3)
for j=1:nz
for i=1:nx
if i == 1
Fo(j,i,k-1)=0;
else
Fo(j,i,k-1)=((Vx(j,i)+Vx(j,i-1))./2).*C(j,i-1,k-1) - D/dx.*C(j,i,k-1) + D/dx.*C(j,i-1,k-1);
end
if i == nx
Fe(j,i,k-1)=0;
else
Fe(j,i,k-1)=((Vx(j,i)+Vx(j,i+1))./2).*C(j,i,k-1) - D/dx.*C(j,i+1,k-1) + D/dx.* C(j,i,k-1);
end
if j == 1
Fs(j,i,k-1)=0;
else
Fs(j,i,k-1)=((Vz(j,i)+Vz(j-1,i))./2).*C(j,i,k-1) - D/dz.*C(j,i,k-1) + D/dz.*C(j,i,k-1);
end
if j == nz
Fn(j,i,k-1)=-qin*Co*og;
else
Fn(j,i,k-1)=((Vz(j,i)+Vz(j+1,i))./2).*C(j+1,i,k-1) - D/dz.*C(j+1,i,k-1) + D/dz.*C(j,i,k-1);
end
C(j,i,k)=C(j,i,k-1) - (dt/dx).*(Fe(j,i,k-1)-Fo(j,i,k-1)) - (dt/dz).*(Fn(j,i,k-1)-Fs(j,i,k-1));
end
end
end
for k=(n/3)+1:n
for j=1:nz
for i=1:nx
if i==1
Fo(j,i,k-1)=0;
else
Fo(j,i,k-1)=((Vx(j,i)+Vx(j,i-1))./2).*C(j,i-1,k-1) - D/dx.*C(j,i,k-1) + D/dx.*C(j,i-1,k-1);
end
if i==nx
Fe(j,i,k-1)=0;
else
Fe(j,i,k-1)=((Vx(j,i)+Vx(j,i+1))./2).*C(j,i,k-1) - D/dx.*C(j,i+1,k-1) + D/dx.* C(j,i,k-1);
end
if j==1
Fs(j,i,k-1)=0;
else
Fs(j,i,k-1)=((Vz(j,i)+Vz(j-1,i))./2).*C(j,i,k-1) - D/dz.*C(j,i,k-1) + D/dz.*C(j,i,k-1);
end
if j==nz
Fn(j,i,k-1)=0;
else
Fn(j,i,k-1)=((Vz(j,i)+Vz(j+1,i))./2).*C(j+1,i,k-1) - D/dz.*C(j+1,i,k-1) + D/dz.*C(j,i,k-1);
end
C(j,i,k)=C(j,i,k-1) - (dt/dx).*(Fe(j,i,k-1)-Fo(j,i,k-1)) - (dt/dz).*(Fn(j,i,k-1)-Fs(j,i,k-1));
end
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Maillage
dxs2 = dx/2;
dys2 = dz/2;
for i=1:nx
xt(i)=(i-0.5)*dx;
end
for j=1:nz
yt(j)=(j-0.5)*dz;
end
for ix = 1:nx
xx1(2*ix-1) = xt(ix) - dxs2;
xx1(2*ix) = xt(ix) + dxs2;
end
for iy = 1:nz
yy1(2*iy-1) = yt(iy) - dys2;
yy1(2*iy) = yt(iy) + dys2;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% VISUALISATION DES RESULTATS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%HYDRODYNAMIQUE
figure(1)
surf(qx)
xlabel('x')
ylabel('z')
view(0,90)
colorbar
%
figure(2)
plot(x,R)
xlabel('Distance au drain')
ylabel('Ratio')
legend
%figure(3)
subplot(2,2,1)
plot(z,Vx(:,1),'r')
subplot(2,2,2)
plot(z,Vx(:,50),'b')
subplot(2,2,3)
plot(z,Vx(:,120),'r')
subplot(2,2,4)
plot(z,Vx(:,190),'b')
%figure(4)
hold on
plot(z,Vx(:,1),'g')
plot(z,Vx(:,20),'b')
plot(z,Vx(:,50),'k')
plot(z,Vx(:,90),'r')
legend ('x=0.1 m',' x=2m', 'x=5 m', 'x=9 m')
%figure(5)
subplot(2,1,1)
plot(z,Vz(:,1),'r')
subplot(2,1,2)
plot(z,Vz(:,90),'b')
%%%
figure(6)
surf(x,z,Vx)
xlabel('x (m)')
ylabel('z (m)')
zlabel('Flux horizontal (m/h)')
title('flux horizontal')
view (0,90)
colorbar
%
figure(7)
surf(Vz)
xlabel('x')
ylabel('z')
zlabel('Flux vertical')
title('Flux vertical')
view (0,90)
colorbar
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%TRANSPORT
figure(10)
surf (x,z,C(:,:,100))
xlabel('x (m)')
ylabel('z (m)')
zlabel('concentration g/l')
% title('concentration')
view(0,90)
colorbar
%
figure(20)
surf (x,z,C(:,:,200))
xlabel('x (m)')
ylabel('z (m)')
zlabel('concentration g/l')
title('concentration')
view(0,90)
colorbar
%%%
figure(30)
surf (x,z,C(:,:,300))
xlabel('x (m)')
ylabel('z (m)')
zlabel('concentration g/l')
%title('concentration')
view(0,90)
colorbar
%%%
figure(40)
surf (x,z,C(:,:,400))
xlabel('x (m)')
ylabel('z (m)')
zlabel('concentration g/l')
%title('concentration')
view(0,90)
colorbar
%%%
figure(50)
surf (x,z,-C(:,:,500))
xlabel('x (m)')
ylabel('z (m)')
zlabel('concentration g/l')
%title('concentration')
view(0,90)
colorbar
figure(60)
surf (x,z,-C(:,:,n))
xlabel('x (m)')
ylabel('z (m)')
zlabel('concentration g/l')
%title('concentration')
view(0,90)
colorbar
%
C1=C(:,10,600);
C2=C(:,20,600);
C3=C(:,30,600);
C4=C(:,40,600);
C5=C(:,50,600);
C6=C(:,60,600);
C7=C(:,70,600);
C8=C(:,80,600);
C9=C(:,90,600);
C10=C(:,99,600);
%
CA1=C(:,10,200);
CA2=C(:,20,200);
CA3=C(:,30,200);
CA4=C(:,40,200);
CA5=C(:,50,200);
CA6=C(:,60,200);
CA7=C(:,70,200);
CA8=C(:,80,200);
CA9=C(:,90,200);
CA10=C(:,99,200);
%
figure(22)
hold on
plot(z,CA1,'r')
plot(z,CA4,'B')
plot(z,CA6,'k')
plot(z,CA8,'m')
plot(z,CA9,'g')
hold off
xlabel ('z(m)')
ylabel ('concentration g/l')
legend ('x=1 m',' x=4 m', 'x=6 m', 'x=8 m', 'x=9 m')
%
figure(11)
hold on
plot(z,C1,'r')
plot(z,C4,'B')
plot(z,C6,'k')
plot(z,C8,'m')
plot(z,C9,'g')
hold off
xlabel ('z(m)')
ylabel ('concentration g/l')
legend ('x=1 m',' x=4 m', 'x=6 m', 'x=8 m', 'x=9 m')
figure(12)
for j=1:n
t(j)=j.*dt;
C11(j)=C(nz,nx,j);
C12(j)=C(30,30,j);
end
plot(t,C11)
ylabel('Concentration g/l')
xlabel('temps (h)')
%for i=1:5:n
% figure(i)
% surf (C(:,:,i));
%end
for j=1:nz
for i=1:nx
x(i) = i*dx;
z(j)=j*dz;
zo(i)=sqrt((abs(qrel).*(L^2-x(i).^2)));
if z(j) < zo(i)
CC=C(:,:,n);
end
end
end
figure(13)
surf(CC)
view (0,90)
%
C01=C(:,10,10);
C02=C(:,10,100);
C03=C(:,10,200);
C04=C(:,10,300);
C05=C(:,10,400);
C06=C(:,10,500);
C07=C(:,10,600);
figure(16)
hold on
plot(z,C01,'r')
plot(z,C02,'B')
plot(z,C03,'k')
plot(z,C04,'m')
plot(z,C05,'k--')
plot(z,C06,'r--')
plot(z,C07,'g--')
hold off
xlabel ('z(m)')
ylabel ('concentration g/l')
legend ('t= 0.33h','t= 3.33h','t= 6.66 h','t = 10 h','t = 13.33h','t = 16.66 h','t = 20 h')
%
C001=C(:,90,10);
C002=C(:,90,100);
C003=C(:,80,200);
C004=C(:,90,300);
C005=C(:,90,400);
C006=C(:,90,500);
C007=C(:,90,600);
figure(17)
hold on
plot(z,C001,'r')
plot(z,C002,'B')
plot(z,C003,'k')
plot(z,C004,'m')
plot(z,C005,'k--')
plot(z,C006,'r--')
plot(z,C007,'g--')
hold off
xlabel ('z(m)')
ylabel ('concentration g/l')
legend ('t= 0.33h','t= 3.33h','t= 6.66 h','t = 10 h','t = 13.33h','t = 16.66 h','t = 20 h')
0