| Les deux révisions précédentes Révision précédente | |
| thermoacoustique [2022/06/19 11:18] – physix | thermoacoustique [2022/06/19 11:19] (Version actuelle) – physix |
|---|
| [[https://web.archive.org/web/20090623041338/http://mshades.free.fr/thermoacoustique/thermoacoustique.html|https://web.archive.org/web/20090623041338/http://mshades.free.fr/thermoacoustique/thermoacoustique.html]] | [[https://web.archive.org/web/20090623041338/http://mshades.free.fr/thermoacoustique/thermoacoustique.html|https://web.archive.org/web/20090623041338/http://mshades.free.fr/thermoacoustique/thermoacoustique.html]] |
| |
| |La thermoacoustique s'occupe des relations entre ondes sonores et transfert thermique. Les transferts thermiques peuvent amplifier les ondes sonores et inversement les ondes sonores peuvent occasionner des transferts thermiques contre nature, d'un point froid vers un point chaud (fonctionnement réfrigérateur ou pompe à chaleur). \\ Ceux qui travaillent sur la combustion cherchent à minimiser cet effet, car il déstabilise la flamme, occasionnant des imbrûlés et des vibrations destructrices auxquelles d'éventuelles parties mobiles sont particulièrement sensibles. \\ Mais l'effet thermoacoustique peut être interessant car il produit du travail, de l'exergie, à partir d'un transfert de chaleur (fonctionnement moteur ou machine thermique). \\ Sur cette page, nous nous attacherons à modéliser le cycle thermoacoustique et à optimiser son rendement exergétique, afin d'atteindre le rendement de Carnot. | {{:thermoacoustique.pdf|:thermoacoustique.pdf |
| <code> | |
| | |
| \\ **Les ondes:** \\ Dans une onde sonore les particules de gaz subissent des variations de position, de vitesse et de pression, généralement sinusoidales, en tout cas cycliques. \\ Une onde peut être progressive (en milieu ouvert) ou stationnaire (en milieu clos). Dans un milieu monodimensionnel, une onde stationnaire est la superposition de 2 ondes progressives de sens opposés. Une onde de pression s'exprime **Pa(x,t)=F(x-a.t)+G(x+a.t)**, où **x** est l'abcisse, **t** le temps, **a** la vitesse locale du son, **Pa** la pression acoustique locale instantanée, **F** et **G** des fonctions réelles continues et dérivables presque partout, généralement périodiques de pulsation **w** (rad/s) et de longueur d'onde **l=2ap/w**. Pour une onde stationnaire F=G, et pour une onde progressive F=0 ou G=0. On voit qu'il existe de nombreuses autres possibilités mais la périodicité en un point n'est plus visible: c'est une sorte de clapotis ou de bruit. \\ \\ | | |
| | |
| </code> |[[https://web.archive.org/web/20090623041338/http://enslyon.free.fr/rapports/physique/Jean-Yves_Tinevez_1.ps.gz|{{https://web.archive.org/web/20090623041338im_/http://mshades.free.fr/thermoacoustique/thermoacTinevez2.gif?nolink&320x232|montage thermoacoustique, tiré de la thèse de J.Y. Tinevez}}]][[https://web.archive.org/web/20090623041338/http://mshades.free.fr/echangeurs/generalitesechangeur.html#pericombustion|{{https://web.archive.org/web/20090623041338im_/http://mshades.free.fr/echangeurs/echcomb.gif?nolink&220x162}}]] Voici à gauche le montage thermoacoustique standard et à droite ce qu'on voudrait réaliser. Le **stack** (la pile de plaques) avec son fort gradient thermique amplifie les ondes sonores. Dans le montage de droite ce rôle est assumé par les alvéoles de l'échangeur. En fonctionnement réfrigérateur le montage de gauche fournit un écart de 16°C entre source froide et source chaude pour 1000Pa de pression acoustique. Le montage de droite vise un fonctionnement moteur. Les différences principales: à gauche une onde stationnaire et un stack de longueur quart d'onde, à droite une onde progressive avec flux de circulation et une longueur d'échangeur égale à plusieurs dizaines de longueurs d'onde (fréquence acoustique élevée, rendement exergétique proche de 100%). Pour analyser correctement le fonctionnement de ces systèmes il faut descendre au niveau de la particule fluide. c'est ce que font les 2 schémas ci-dessous.| | |
| |[[https://web.archive.org/web/20090623041338/http://enslyon.free.fr/rapports/physique/Jean-Yves_Tinevez_1.ps.gz|{{https://web.archive.org/web/20090623041338im_/http://mshades.free.fr/thermoacoustique/cyclethermoacTinvz2.gif?nolink&495x265}}]]cycles de la particule en diagramme (P,v)| | | |
| |[[https://web.archive.org/web/20090623041338/http://www.acs.psu.edu/THERMO/frank/thesis.ps|{{https://web.archive.org/web/20090623041338im_/http://mshades.free.fr/thermoacoustique/cyclesfrank.gif?nolink&687x274}}]]cycles de la particule en diagramme (T,x)| | | |
| |{{https://web.archive.org/web/20090623041338im_/http://mshades.free.fr/thermoacoustique/TAcycleTx.gif?nolink&213x238}}Considérons le cycle thermodynamique de la particule fluide de la manière un peu caricaturale représentée ci-contre en diagramme (T,x), c'est à dire température-position. \\ On a 4 évolutions: \\ 0-1. compression adiabatique (sans échange de chaleur) sans frottement, donc isentropique. \\ 1-2. réchauffe à pression constante. \\ 2-3. détente adiabatique isentropique. \\ 3-0. refroidissement à pression constante. \\ Si les compression et détente adiabatiques sont bien synchronisées avec le déplacement de la particule, celle-ci peut suivre exactement le gradient thermique de la paroi. Ensuite, lors de la réchauffe ou du refroidissement isobare, on continue à longer le gradient avec un léger DT pour assurer le transfert thermique, DT qui peut être aussi petit qu'on veut, typiquement qq dixièmes de degré, en jouant sur l'écartement entre plaques ou le diamètre d'alvéole: on l'a prouvé dans l'optimisation de l'échangeur. \\ En diagramme (T,x) le cycle peut donc être quasiment plat et donc les pertes exergétiques par saut de température presque nulles. Cela ne signifie pas que la puissance motrice du cycle est nulle. Pour le voir il faut passer en diagramme (T,S) ci-dessous. Notons toutefois avant de quitter ce diagramme que comme les phases d'échange thermique ne sont pas exactement (et même pas du tout, voir plus bas) superposées dans l'espace, le flux de chaleur moyen sur un cycle n'est pas constant le long de la paroi, ce qui doit pourtant être le cas dans une configuration d'échangeur, il faut le compenser par de multiples cycles décalés en abcisse x correspondant à autant de particules voisines, autrement dit la longueur de déplacement x<sub>2</sub>-x<sub>0</sub> doit être faible devant la longueur totale **L** de l'échangeur. | | | |
| |{{https://web.archive.org/web/20090623041338im_/http://mshades.free.fr/thermoacoustique/TAcycleTS.gif?nolink&182x201}}Voici le même cycle en diagramme (T,S), c'est à dire température-entropie. Toute quantité de chaleur infinitésimale dQ reçue par la particule entraîne une variation d'entropie dS=dQ/T, d'où l'on tire dQ=T*dS, puis Q=integrale<sub>cycle</sub> (T*dS). La surface à l'intérieur du cycle est donc égale à la chaleur **Q** reçue par la particule pendant le cycle. Soit **W** le travail reçu, on a W+Q=0 sur un tour (par conservation de l'énergie). Le travail fourni -W est donc égal à la surface du cycle, qu'on cherche à maximiser. Comme T0 et T2 sont déterminés par les positions x0 et x2 sur la paroi qu'on ne veut pas trop écarter pour éviter les frottements, reste à trouver T1 qui maximise Q à T0 et T2 données. Comme dS=Cp*dT/T-R*dP/P, on a dS=Cp*dT/T sur l'évolution 1-2, donc **S2-S1=Cp*ln(T2/T1)=Cp*ln(T3/T0)**. De plus **Q<sub>1-2</sub> =Cp*(T2-T1)** et **Q<sub>3-0</sub> =Cp*(T0-T3)**. Comme Q=Q<sub>1-2</sub> +Q<sub>3-0</sub>, Q=Cp*(T2-T1+T0-T3). Or T2/T1=T3/T0=**K**. Alors Q=Cp*(T1-T0)*(K-1), Q=Cp*(T1-T0)*(T2/T1-1). \\ Dérivons Q par rapport à T1: dQ/d(T1)=Cp*[(T2/T1-1)+(T0-T1)*T2/T1²]=Cp*[T0*T2/T1²-1]. Le maximum de Q est atteint quand **T1²=T0*T2**, c'est à dire quand T1 est moyenne géométrique de T0 et T2. Alors Q=Cp*T0<sup>1/2</sup> *(T2<sup>1/2</sup> -T0<sup>1/2</sup> )*(T2<sup>1/2</sup> /T0<sup>1/2</sup> *-1), **travail fourni=Q=Cp*(T2<sup>1/2</sup> -T0<sup>1/2</sup> )<sup>2</sup> ** (en J/kg). Quand T0 et T2 sont voisines, posons T0=T et T2=T+dT. On obtient la formule approximative **Q=Cp*dT²/(4*T)**. On voit donc qu'on a intérêt à éviter les dT trop petits. Au passage calculons T3: K=T3/T0=T2/T1=T1/T0 donc **T3=T1**,** x3=x1**; les zones d'échange thermique [x1,x2] et [x3,x0] ne sont pas du tout superposées. La particule qui dépose la chaleur n'est pas la même que celle qui la prend. Corrigez (mentalement) les 2 diagrammes ci-dessus. \\ Ce travail Q est fourni aux particules voisines sous forme de compression, et il faut absolument l'évacuer du système le plus rapidement possible, sans quoi des phénomènes dissipatifs (ondes de choc,frottements intenses,écarts fluide-paroi…) vont se former, réduisant à néant notre optimisation. On peut récupérer ce travail en organisant dans le stack ou l'alvéole d'échangeur thermoacoustique un gradient de pression **dP/dx**, joint à une circulation moyenne **<u>** dans le même sens, et des clapets anti-retour (à temps de réponse rapide) aux bornes du stack. De cette façon le travail produit est emmagasiné au plus près de la source, dans la particule elle-même, sous forme de compression. \\ Alors f*Q*r*f²=f²*<u>*dP/dx (en W/m), où **f** est le diamètre de l'alvéole, **f** la fréquence du cycle (en Hz), **r** la masse volumique (en kg/m3). On simplifie cette formule en **f*Q*r=<u>*dP/dx** (en W/m3). <u> et dP/dx étant déterminés par la chaleur de combustion et l'optimisation de l'échangeur, reste à déterminer la fréquence f. \\ En première approche, le travail fourni à l'abcisse x pendant la phase 0-1, dW/dx=-d(Pu)/dx (en W/m3), est égal à la quantité de chaleur dQ/dx échangée dans la phase isobare 1-2 (valeur appelée compacité à la base du dimensionnement optimisé de l'alvéole, disons 1MW/m3), et la vitesse est de l'ordre de 50cm/s, disons 1m/s. Quand la pression est disons 1bar(=1e5 Pa), Pu varie entre 0 (au point 0) et 1e5 (au point 1). On a donc d(Pu)=1e5, puis dx=1e5/1e6=0.1m=10cm. La phase de compression dure 10cm, la phase de réchauffe 10cm aussi, et pour un gradient de 1500K/m, le dT=T2-T0 est de 300K, donc Q=256 J/kg pour de l'air,Q*r=308 J/m3=308Pa. Si dP/dx est de l'ordre de 2000 Pa/m (du même ordre de grandeur que les frottements sans effet thermoacoustique) et <u>=1m/s, f*Q*r=2000W/m<sup>3</sup> , d'où f=2000/308=6.5Hz et la longueur d'onde 300/6.5=46m bien trop grande. \\ Transformons la formule f*Q*r=<u>*dP/dx. <u>*dP/dx=f*Cp*dT²/(4*T)*r=f*(g/4(g-1))*P*(dT/T)² puis <u>*dP/Pdx=f*(g/4(g-1))*(dT/T)². \\ Il y a une relation entre f,u et dT: dT=dxT*u/f à un facteur près. Donc dP/Pdx=u*dxT²/T²f. Il faut donc maximiser u/f (en m), c'est à dire le déplacement, quand dxT est donné.| | | |
| |{{https://web.archive.org/web/20090623041338im_/http://mshades.free.fr/thermoacoustique/TAcyclePV.gif?nolink&190x169}}Ceux qui ne connaissent pas très bien la thermodynamique préfèrent utiliser le diagramme (P,V) ci-contre, c'est à dire pression-volume, dans lequel le travail W=integrale<sub>cycle</sub> (-P*dV) apparaît directement. Mais comme on travaille à T0 et T2 imposées, son utilisation est moins pratique dans le cas présent, et d'ailleurs maintenant inutile. On peut quand même y lire que le travail est fourni pendant les phases 1-2 et 2-3. | | | |
| | | | | |
| | |
| |équation |unité |coord. euleriennes (x,t) |définitions | | |
| |conservation de la masse |kg/m³ |(dr/dt)<sub>x</sub> +(d(ru)/dx)<sub>t</sub> =0 |u=(dx/dt)<sub>a</sub> =vitesse (faible: <1m/s) | | |
| |conservation de l'impulsion |N/m³ |d(ru)/dt+dP/dx=0 (oublie les frottements) |T<sub>i</sub> =T+u²/(2*cp) | | |
| |travail reçu |W/m³ |dW/dx=-d(Pu)/dx |Pu est le work flow (en W/m²) | | |
| |transfert thermique |W/m³ |dQ/dx=8plDT/f² |DT=Tparoi-T (K) | | |
| |conservation de l'énergie |W/m³ |dH/dt=d(rcp*T)/dt=(dP/dt)*g/(g-1)=(dW/dx)+(dQ/dx) | | | |
| |équation d'état du gaz parfait | |P=rrT | | | |
| |gradient de paroi |K/m |d(Tp)/dx=d<sub>x</sub> T constant, d(Tp)/dt=0 |h=viscosité(Pa.s) | | |
| |frottement (parabolique) | |(dP/dx)<sub>f</sub> =-32uh/f² |grave;s pratiques car les particules successives ne sont pas en début de cycle au même moment. Il vaut mieux remplacer la coordonnée lagrangienne **a** par une coordonnée de phase **j **qui varie de 0 à 1 quand la particule parcourt son cycle. Pour une particule donnée j est adimensionnelle et varie linéairement dans le temps: j=(t-t0)/T, où t est le temps, t0 le temps au début de cycle (point 0) pour la particule considérée et T la période (en secondes). Il faut faire attention à ce que le temps de période T varie éventuellement avec la particule. | | |
| |{{https://web.archive.org/web/20090623041338im_/http://mshades.free.fr/thermoacoustique/gradP.jpg?nolink&330x287|mshades.free.fr_thermoacoustique_gradp.jpg}}**Profil longitudinal instantané des pressions** \\ Sur cette figure on a repéré en rouge les 4 points clés 0,1,2,3 du cycle dans la partie réchauffante de l'échangeur (dans la partie refroidissante le sens de propagation de l'onde sonore est inversé: en fonctionnement moteur, la partie haute pression du cycle doit toujours être aux hautes températures). Au point 0, la vitesse u est nulle, donc dP/dx>0 pour avoir une accélération. Ceci nous donne une relation entre la pression acoustique p=A*sin(w(t-x/a)) où a est la vitesse du son, et le gradient de pression temporel moyen dxP. **-dp/dx=A*(w/a)*cos(w(t-x/a))**>dxP au point 0, donc **A*(w/a)>dxP**. Comme w=2p*f, plus la fréquence est élevée et courte la longueur d'onde, plus on peut se contenter d'une faible amplitude. La vitesse u aussi est sinusoidale, à une constante près: **u=B*sin(w(t-x/a)+D)+<u>**, où D est le déphasage. On en tire **du/dt=B*w*cos(w(t-x/a)+D)** et **du/dx=-B*(w/a)*cos(w(t-x/a)+D)**. En considérant que les efforts visqueux à la paroi compensent en moyenne le gradient de pression dxP, on a l'équation du/dt-u*(du/dx)=-dp/rdx (voir plus bas). Avec des ondes sinusoidales, ça donne: \\ B*w*cos(w(t-x/a)+D) +(B*sin(w(t-x/a)+D)+<u>)*B*(w/a)*cos(w(t-x/a)+D) =A*(w/a)*cos(w(t-x/a))/r. \\ B*w*cos(w(t-x/a)+D) +B²*(w/2a)*sin(2w(t-x/a)+D) +<u>*B*(w/a)*cos(w(t-x/a)+D) =A*(w/a)*cos(w(t-x/a))/r. \\ B*w*(1+<u>/a)*cos(w(t-x/a)+D) +B²*(w/2a)*sin(2w(t-x/a)+D) =A*(w/a)*cos(w(t-x/a))/r. \\ B*w*(1+<u>/a)*[cos(w(t-x/a))cos(D)-sin(w(t-x/a))sin(D)] +B²*(w/2a)*sin(2w(t-x/a)+D) =A*(w/a)*cos(w(t-x/a))/r. \\ w*[B*(1+<u>/a)*cos(D)-A/ar]*cos(w(t-x/a)) -w*B*(1+<u>/a)*sin(D)*sin(w(t-x/a)) +B²*(w/2a)*sin(2w(t-x/a)+D)=0. \\ <u> est petit devant a: \\ w*[B*cos(D)-A/ar]*cos(w(t-x/a)) -w*B*sin(D)*sin(w(t-x/a)) +B²*(w/2a)*sin(2w(t-x/a)+D)=0. \\ B est petit devant a: \\ w*[B*cos(D)-A/ar]*cos(w(t-x/a)) -w*B*sin(D)*sin(w(t-x/a))=0. \\ Ceci est valable en tout (x,t) donc [B*cos(D)-A/ar]=0 et B*sin(D)=0, puis sin(D)=0, D=0 ou p, cos(D)=-1 ou 1. Comme A et B sont positifs, cos(D)=0, **D=0**, puis B=A/ar, B/a=A/a²r. a est la vitesse isotherme du son, a²=sT=P/r, donc **B/a=A/P**, relation qui lie les amplitudes de vitesse et de pression. \\ Calculons maintenant le déplacement Dx d'une particule pendant la phase 0-2. Comme u=(dx/dt) particulaire, dx/dt=B*sin(w(t-x/a))+<u>. \\ x est des 2 cotés de l'équation, c'est ennuyeux mais on peut considérer que pour une particule donnée x/a varie beaucoup moins que t, donc sur un demi-cycle Dx=integrale<sub>[0,p/w]</sub> ((B*sin(wt)+<u>)*dt), Dx=<u>*p/w +B*2p/w, Dx=p/w *(<u>+2B), Dx=(1/2f)*(<u>+2B), Dx=(<u>/2+B)/f Le <u> sert avec le dxP, à récupérer le travail du cycle. Oublions le pour l'instant. On a alors Dx=B/f, puis dT=dxT*B/f. Le travail du cycle est donc Q=rCp*dT²/(4*T) (en W/m3, voir plus haut) puis Q=rCp*(dxT*B/f)²/(4*T)=<u>*dxP. En utilisant B=Aa/P et a²=sT, on obtient Cp*(dxT*A/f)²/4(g-1)P=<u>*dxP. En posant A=K*a*dxP/f (avec K>1/2p), on trouve Cp*(dxT*K/f²)²/4(g-1)r=<u>/dxP. Si <u> est fixé (par l'optimisation de l'échangeur), dxP est donc proportionnel à f 4 . Bon point pour les hautes fréquences. Mais en fait la fréquence f n'est pas choisie librement. Pour un échangeur donné, elle est conditionnée par l'apport de chaleur en bout d'échangeur, apport qu'à l'idéal on souhaiterait variable à loisir.| _FCKG_BLANK_TD_ écoulement laminaire en gradient de pression adverse Une telle configuration du profil des vitesses est instable. Le gradient de pression positif engendre un contre-courant au voisinage de la paroi. On a -(h/r)*(du/dy) paroi @@_SP_@@=t=(f/4)*dP/dx. La pression est homogène dans toute la section car les mécanismes d'équilibre se propagent à la vitesse du son. Le courant principal (central) s'amenuise dans son sens de propagation (S 1 @@_SP_@@se réduit). De la même façon le contre-courant périphérique s'affaiblit en remontant (la section annulaire S 2 @@_SP_@@se réduit). Sans apport de travail l'écoulement central forme donc une sorte de cône. Il serait tentant de modéliser le profil de vitesses par une double parabole (ce qui donne un modèle à 2 paramètres: S1 et (du/dy) paroi ), mais un polynôme unique sur toute la section (forcément une quadrique, car la vitesse s'annule 4 fois, ce qui donne un modèle à 5 paramètres) est plus facile à manipuler et offre plus de possibilités. Soient R=f/2 et r=R-y; u® est une fonction paire qui s'annule en r=R et r=R1, qui possède un maximum local en r=0. u®=A*(1-r/R)*(1+r/R)*(1-r/R1)*(1+r/R1), u®=A*(1-r²/R²)*(1-r²/R1²). Incidemment le modèle n'a plus que 2 paramètres. Le profil parabolique est un cas particulier correspondant à R1=infini. Si R1=R, (du/dy) paroi @@_SP_@@=0. C'est la famille de profils qu'on obtient par dégénerescence du profil parabolique lorsque la pression motrice dP/dx a disparu. Si R1=i*R, u®=A*(1-r 4 /R 4 ) et (du/dy) paroi @@_SP_@@est maximal pour A=Umax donné. C'est la configuration d'entrée dans l'alvéole, la plus proche de la distribution uniforme des vitesses. On peut reformuler ainsi le profil: u®=Umax*(1-p*r²)*(1-b*r²) , où p=1/R² @@_SP_@@est le coefficient de paroi (invariable), b=1/R1² @@_SP_@@le coefficient de rebroussement et Umax la vitesse centrale. u/Umax=1-(p+b)*r²+pb*r 4 @@_SP_@@, (du/dr)/Umax=-2(p+b)*r+4pb*r 3 @@_SP_@@, (du/dr)=Umax*2r*(2pb*r²-(p+b)) @@_SP_@@qui s'annule pour r=0 et r²=(p+b)/(2pb)=R²*R1²*(1/R²+1/R1²)/2=(R²+R1²)/2 donc r(rebroussement)~=(R+R1)/2 ce qui correspond à peu près à la figure. D'autre part le gradient de paroi est du/dy=-(du/dr) R @@_SP_@@=-Umax*2R*(2pb*R²-(p+b))=Umax*2R*(p-b)=-4(f/n)*dP/dx, où n est la viscosité cinématique. d²u/dr²=Umax*(-2(p+b)+12pb*r²) , est la dérivée seconde au point r, à t et x constants. d²u/dr² s'annule quand r²=(p+b)/6pb=(R²+R1²)/6, r(inflexion)=0.58*r(rebroussement). Calculons la vitesse débitante U . U=integrale(2pr*dr*u®)/pR²=integrale(u®*d(r²))/R², U=Umax*(R²-(p+b)R 4 /2+pbR 6 /3)/R², U=Umax*(1-(p+b)R²/2+pbR 4 /3), U=Umax*(1/2-b/2p+b/3p), @@_SP_@@U=Umax*(1-b/3p)/2 Quant au carré moyen de la vitesse <u²> @@_SP_@@=integrale(u²®*d(r²))/R², <u²>=Umax*integrale((1-(p+b)r²+pbr 4 @@_SP_@@)²*d(r²)) /R², <u²>=Umax*integrale((1-2(p+b)r²+((p+b)²+2pb)r 4 -2(p+b)pbr 6 +p²b²r 8 )*d(r²)) /R², <u²>=Umax*(R²-(p+b)R<sup>4</sup> +((p+b)²+2pb)R 6 /3-(p+b)pbR 8 /2+p²b²R 10 /5)/R², <u²>=Umax*(1-(p+b)R 2 +((p+b)²+2pb)R 4 /3-(p+b)pbR 6 /2+p²b²R 8 /5), <u²>=Umax*(-b/p+((1+b/p)²+2b/p)/3-(1+b/p)b/2p+b²/5p²), <u²>=Umax*(1/3+(b/p)*(-1+2/3+2/3-1/ |