Je me dis qu'une étape supplémentaire pour avoir la figure en 3D est d'avoir aussi les isobares. Voilà ce que j'ai fait. Ça marche... presque
Code : Tout sélectionner
unitsize(4cm,8cm);
import graph;
// t est la température réduite
real t;
// C_isotherme_VdW : isotherme de Van der Waals (conduisant pour t<1
// à des propriétés non physiques)
path C_isotherme_VdW;
path C_isobare_VdW;
// C_isotherme_Andrews : isotherme de Van der Waals corrigées
// (palier de liquéfaction)
path C_isotherme_Andrews;
path C_isobare_Andrews;
// Cg et C_palier sont des courbes auxiliaires nécessaires au calcul
path Cg,C_palier;
path Cg_isobare,C_palier_isobare;
// vmin et vmax sont les bornes du volume réduit utilisées pour la
// détermination du palier de liquéfaction (ne correspond pas forcément
// aux bornes de l'affichage définitif)
// On prend vmin très légèrement supérieur à 1/3 car sinon problème de
// définition de la pression réduite p (voir infra)
real vmin=1/3+0.000001;
real vmax=30;
// precision est la précision du calcul d'aire (Maxwell)
real precision=.0001;
// Variable auxiliaire pour le calcul d'aire
real diff_aires;
// nb_points est le nombre de points servant à obtenir C_isotherme_VdW
int nb_points=1000;
// Variables auxiliaires dichotomie
real a,b,c;
// Guide (un peu un équivalent de path) pour la courbe de saturation
// construite point par point dans la boucle suivante
guide C_saturation_liquide,C_saturation_vapeur;
guide C_saturation_liquide_isobare,C_saturation_vapeur_isobare;
// Style du tracé des isothermes
pen couleur_hypercritique = red;
pen couleur_liquide = .5green;
pen couleur_vapeur = purple;
pen couleur_liquide_vapeur = couleur_liquide + couleur_vapeur;
pen style_isotherme;
// Un tableau pour conserver les pressions de changement d'état calculées
// pour les isothermes. Cela peut être utile pour construire les isobares
// "qui vont bien".
real[] pression_palier_isotherme;
real[] vmin_palier,vmax_palier;
//
int imin=-10,imax=10;
// Boucle sur différentes valeurs de la température réduite
for(int i=imin; i<imax; ++i){
t=1+i/100;
write(t);
// Équation de Van der Waals en coordonnées réduites
real p(real v){return 8*t/(3*v-1)-3/v^2;}
// La courbe représentant l'isotherme de Van der Waals est calculée
C_isotherme_VdW=graph(p,vmin,vmax,n=nb_points,Spline);
if(t>=1){
// Si la température critique est supérieure à 1, le fluide de Van der
// Waals est supercritique. Il ne subit pas de changement d'état.
// L'isotherme d'Andrews est identique à l'isotherme de Van der Waals.
C_isotherme_Andrews=C_isotherme_VdW;
pression_palier_isotherme[i-imin]=t;
vmin_palier[i-imin]=-1;
vmax_palier[i-imin]=-1;
// draw(C_isotherme_Andrews,couleur_hypercritique);
}
else{
// Sinon, il faut déterminer la position du palier de liquéfaction.
// Ici, on utilise la construction de Maxwell : voir par exemple
// www.lerepairedessciences.fr/sciences/agregation_fichiers/thermo/maxw.pdf
//
// Soit v1 le volume réduit corespondant au minimum
// (tangente horizontale) de l'isotherme de Van der Waals.
// Cg est l'isotherme de Van der Waals pour v > v1.
Cg=subpath(C_isotherme_VdW,
dirtime(C_isotherme_VdW,(1,0)),
reltime(C_isotherme_VdW,1));
// Soit v2 le volume réduit correspondant au maximum de l'isotherme de
// Van der Waals. Cg est l'isotherme de Van der Waals pour v1 < v < v2.
Cg=subpath(Cg,0,maxtimes(Cg)[1]);
// a est la pression réduite pour le volume v1 (le minimum de p)
a=max(min(Cg).y,.02);
// b est la pression réduite pour le volume v2 (la maximum de p)
b=max(Cg).y;
// On initialise diff_aires à une valeur supérieure pour rentrer dans
// la boucle while
diff_aires=2*precision;
while(abs(diff_aires)>precision)
{
// On amorce une boucle pour procéder par dichotomie afin
// de trouver la pression réduite du palier de liquéfaction
// donnant des aires égales compte tenu de la précision
// demandée.
c=(a+b)/2;
// v_palier_inter_isotherme_VdW est un tableau donnant généralement
// 3 valeurs de v correspondant à la valeur de la pression d'essai
// du palier stockée dans c
real[] v_palier_inter_isotherme_VdW=times(C_isotherme_VdW,(0,c));
if(v_palier_inter_isotherme_VdW.length!=3) {
// Problème...
write("attention");
break;
}
// Pour calculer l'aire de Maxwell, il faut recentrer les courbes
// sur l'axe des abscisses. On retranche dans à la pression de Van
// der Waals p(v) la pression d'essai du palier c.
real g(real v){return p(v)-c;}
// Méthode de simpson pour le calcul d'intégrale
diff_aires=simpson(g,
point(C_isotherme_VdW,
v_palier_inter_isotherme_VdW[0]).x,
point(C_isotherme_VdW,
v_palier_inter_isotherme_VdW[2]).x);
// Dichotomie
if(diff_aires<0) b=c; else a=c;
}
// Construction du palier (la pression réduite de liquéfaction vaut c)
real h(real v){return c;}
pression_palier_isotherme[i-imin]=c;
// Deux points suffisent pour construire une droite
C_palier=graph(h,vmin,vmax,n=2);
real[] v_palier_inter_isotherme_VdW=times(C_isotherme_VdW,(0,c));
if(v_palier_inter_isotherme_VdW.length==3){
C_saturation_liquide=C_saturation_liquide..
point(C_isotherme_VdW,v_palier_inter_isotherme_VdW[0]);
C_saturation_vapeur=C_saturation_vapeur..
point(C_isotherme_VdW,v_palier_inter_isotherme_VdW[2]);
C_isotherme_Andrews=firstcut(C_isotherme_VdW,C_palier).before
-- lastcut(C_isotherme_VdW,C_palier).after;
// draw(firstcut(C_isotherme_VdW,C_palier).before,couleur_liquide);
// draw(point(C_isotherme_VdW,v_palier_inter_isotherme_VdW[0]) --
// point(C_isotherme_VdW,v_palier_inter_isotherme_VdW[2]),
// couleur_liquide_vapeur);
// draw(lastcut(C_isotherme_VdW,C_palier).after,couleur_vapeur);
vmin_palier[i-imin]=v_palier_inter_isotherme_VdW[0];
vmax_palier[i-imin]=v_palier_inter_isotherme_VdW[2];}
write(v_palier_inter_isotherme_VdW);
write("Vmin "+string(vmin_palier[i-imin]));
}
// draw(C_isotherme_Andrews,style_isotherme);
}
//draw(C_saturation_liquide..reverse(C_saturation_vapeur),blue);
write(pression_palier_isotherme);
real p;
// Boucle sur différentes valeurs de la pression réduite
for(int i=imin; i<imax; ++i){
p=pression_palier_isotherme[i-imin];
write(p);
// Équation de Van der Waals en coordonnées réduites
real t(real v){return (p+(3/v^2))*(3*v-1)/8;}
// La courbe représentant l'isotherme de Van der Waals est calculée
C_isobare_VdW=graph(t,vmin,vmax,n=nb_points,Spline);
if(p>=1){
// Si la pression critique est supérieure à 1, le fluide de Van der
// Waals est supercritique. Il ne subit pas de changement d'état.
// L'isobare d'Andrews est identique à l'isobare de Van der Waals.
C_isobare_Andrews=C_isobare_VdW;
draw(C_isobare_Andrews,couleur_hypercritique);}
else{
// Sinon, il faut déterminer la position du palier de liquéfaction.
real t_palier(real v){return t(vmin_palier[i-imin]);}
// Deux points suffisent pour construire une droite
write(t_palier(vmin));
C_palier_isobare=
graph(t_palier,vmin,vmax,n=2);
C_saturation_liquide_isobare=C_saturation_liquide_isobare..
point(C_isobare_VdW,vmin_palier[i-imin]);
C_saturation_vapeur_isobare=C_saturation_vapeur_isobare..
point(C_isobare_VdW,vmax_palier[i-imin]);
C_isobare_Andrews=firstcut(C_isobare_VdW,C_palier_isobare).before
-- lastcut(C_isobare_VdW,C_palier_isobare).after;
draw(firstcut(C_isobare_VdW,C_palier_isobare).before,couleur_liquide);
draw(point(C_isobare_VdW,vmin_palier[i-imin]) --
point(C_isobare_VdW,vmax_palier[i-imin]),
couleur_liquide_vapeur);
draw(lastcut(C_isobare_VdW,C_palier_isobare).after,couleur_vapeur);
}
}
xlimits(-.5,4.5,Crop);
ylimits(-.2,1.5,Crop);
xaxis(Label("$v_r$",position=EndPoint, align=SE),
xmin=0, xmax=4.6,
Ticks(scale(.7)*Label(),
NoZero,
Step=.5,step=.1,
Size=2mm, size=1mm,
pTick=black,ptick=gray),
Arrow);
yaxis(Label("$t_r$",position=EndPoint, align=NE),
ymin=-.1,ymax=1.55,
Ticks(scale(.7)*Label(),
NoZero,
Step=0.5,step=.1,
Size=2mm, size=1mm,
pTick=black,ptick=gray),
Arrow);
write("patienter encore un peu...");
