Godunovs metod

Den aktuella versionen av sidan har ännu inte granskats av erfarna bidragsgivare och kan skilja sig väsentligt från versionen som granskades den 9 april 2016; kontroller kräver 4 redigeringar .

Godunovmetoden är en implementering av genomräkningsscheman som kan användas för att beräkna gasdynamiska flöden med diskontinuiteter i parametrar inom beräkningsdomänen. Detta schema föreslogs av S. K. Godunov 1959. Godunov-metoden är en variant av kontrollvolymmetoden . Flödena genom sidoytorna bestäms av lösningen av problemet med sönderfallet av en godtycklig diskontinuitet . Låt oss förklara med ett exempel.

Exempel

Överväg konstruktionen av den numeriska Godunov-metoden av första noggrannhetsordningen med hjälp av exemplet på att lösa ekvationssystemet för endimensionell instabil gasdynamik, skriven i divergerande form :

{\begin{cases}{\begin{array}{lll}{\frac {\partial \rho }{\partial t))+{\frac {\partial \rho u}{\partial x))&=&0 \\{\frac {\partial (\rho u)}{\partial t))+{\frac {\partial (p+\rho u^{2})}{\partial x))&=&0\\{ \frac {\partial E}{\partial t}}+{\frac {\partial u(E+p)}{\partial x}}&=&0\end{array}}\end{cases}}

Här:

$\rho$ - densitet
$sid$ - tryck
$u$ - hastighet
$E$ är den totala energin per volymenhet

Lägg märke till att:

Den första ekvationen uttrycker lagen om massans bevarande
Den andra ekvationen uttrycker momentumbevarandelagen
Den tredje ekvationen uttrycker lagen om energibevarande $E=e+{\frac {\rho u^{2}}{2}}$ $e=\rho \varepsilon ,$ $\varepsilon =\varepsilon (p,\rho ),$ $\varepsilon ={\frac {1}{\gamma -1}}\cdot {\frac {p}{\rho }}$ , var är den adiabatiska exponenten $\gamma$

Differentialform

Det initiala systemet kan skrivas i en mer kompakt form:

{\frac {\partial q}{\partial t))+{\frac {\partial f}{\partial x))=0

var:

$q$ är vektorn för konservativa variabler $q=\left({\begin{array}{c}\rho \\\rho \,u\\E\end{array}}\right)$
$f$ är flödesvektorn $f=\left({\begin{array}{c}\rho \,u\\p+\rho \,u^{2}\\u(E+p)\end{array}}\right)$

Integralform

Istället för ekvationernas differentialform härleder vi en ny integralform av ekvationerna, mer lämpad för att representera en svag lösning . Här är en svag lösning en generaliserad funktion definierad av integrallikheter erhållna från motsvarande differentialekvationer och de initiala förhållandena för problemet. För att göra detta väljer vi en viss kontrollvolym och integrerar ekvationssystemet över denna volym. Vi tillämpar den generaliserade Stokes-satsen på den erhållna integralen av divergens (för två oberoende variabler kommer detta att vara Greens sats och Ostrogradsky-Gauss-formeln i tredimensionellt rymd). I det här fallet introducerar vi riktningen för att kringgå konturen moturs . $\Omega$

Separat, med tanke på kontinuitetsekvationen , får vi:

\iint \limits _{{\Omega }}\left({\frac {\partial \rho }{\partial t}}+{\frac {\partial \rho u}{\partial x}}\right)d \,xd\,t=\oint \limits _{{\partial \Omega }}\rho dx-\rho udt

För hela ekvationssystemet

\iint \limits _{{\Omega }}\left({\frac {\partial q}{\partial t}}+{\frac {\partial f}{\partial x}}\right)d\,xd \,t=0\Leftrightarrow \oint \limits _{{\partial \Omega }}(qdx-fdt)=0

Att skriva systemet i utökad form:

{\begin{cases}{\begin{array}{lll}\oint \limits _{{\partial \Omega }}(\rho \;d\,x-\rho u\;d\,t)&= &0\\\oint \limits _{{\partial \Omega }}(\rho u\;d\,x-(p+\rho u^{2})\;d\,t)&=&0\\\ oint \limits _{{\partial \Omega }}(E\;d\,xu(p+E)\;d\,t)&=&0\\\end{array}}\end{cases}}

Approximation

Övergången från differentialformen för att skriva det ursprungliga ekvationssystemet till integralformen görs. Integralformen skrivs som lika med noll av integralerna över konturen (gränsen för den valda kontrollvolymen) från vektorerna för konservativa variabler och flöden. Vi representerar konturintegralen som summan av integraler över sektioner (intervall) 1-2 , 2-3 , 3-4 , 4-1 av kontrollvolymen i figuren (som ännu inte är tillgänglig) och på varje sektion uppskattar vi integralen med hjälp av metoden med rektanglar som en produkt av integranden i mitten av intervallet med längden på integrationsintervallet:

{\begin{array}{c}q_{{12}}\cdot (x_{{2}}-x_{{1}})-f_{{12}}\cdot (t_{{2}}-t_ {{1}})+\\\qquad q_{{23}}\cdot (x_{{3}}-x_{{2}})-f_{{23}}\cdot (t_{{3}} -t_{{2}})+\\\qquad \qquad q_{{34}}\cdot (x_{{4}}-x_{{3}})-f_{{34}}\cdot (t_{ {4}}-t_{{3}})+\\\qquad \qquad \qquad q_{{41}}\cdot (x_{{1}}-x_{{4}})-f_{{41} }\cdot (t_{{1}}-t_{{4}})=0\end{array}}\Rightarrow q_{{34}}=q_{{12}}-{\frac {\Delta t} {\Delta x}}(f_{{23}}-f_{{41}})

med hänsyn till de likheter som är giltiga för kontrollvolymen byggd på det kartesiska beräkningsnätet:

$x_{3}-x_{2}=0$
$x_{1}-x_{4}=0$
$t_{2}-t_{1}=0$
$t_{4}-t_{3}=0$

Dessutom:

$x_{2}-x_{1}=\Delta x$
$x_{4}-x_{3}=-\Delta x$
$t_{3}-t_{2}=\Delta t$
$t_{1}-t_{4}=-\Delta t$

hitta värdena för vektorn av konservativa variabler på intervallet 3-4 som hör till det nya lagret:

q_{{34}}=q_{{12}}-{\frac {\Delta t}{\Delta x}}(f_{{23}}-f_{{41}})\qquad \Rightarrow \qquad q_ {{j}}^{{n+1}}=q_{{j}}^{{n}}-{\frac {\Delta t}{\Delta x}}(f_{{j+{\frac { 1}{2))))-f_{{j-{\frac {1}{2))))))

I detta fall betecknar värdena med halvheltalsindex flödena av lagrade värden genom gränserna för beräkningscellen under tiden eller flödena genom sidoytorna ( 2-3 och 4-1 ) av kontrollen volym. Om flödeshastigheten är riktad i samma riktning som den yttre normalen till sidoytan, är flödet negativt , det vill säga det rinner ut ur kontrollvolymen och vice versa.

Expanderat:

{\begin{cases}{\begin{array}{lll}\rho _{{j}}^{{n+1}}&=&\rho _{{j}}^{{n}}-{ \frac {\Delta t}{\Delta x))((\rho u)_{{j+{\frac {1}{2))))-(\rho u)_{{j-{\frac { 1}{2}}}})\\(\rho u)_{{j}}^{{n+1}}&=&(\rho u)_{{j}}^{{n}} -{\frac {\Delta t}{\Delta x))((p+\rho u^{2})_{{j+{\frac {1}{2))))-(p+\rho u^{ 2})_{{j-{\frac {1}{2)))))\\E_{{j}}^{{n+1}}&=&E_{{j}}^{{n} }-{\frac {\Delta t}{\Delta x))((u(p+E))_{{j+{\frac {1}{2))))-(u(p+E)) _{{j-{\frac {1}{2}}}})\end{array}}\end{cases}}

Flödena genom sidoytorna och bestäms från lösningen av problemet med förfallet av en godtycklig diskontinuitet . $f_{{j+{\frac {1}{2))))$ $f_{{j-{\frac {1}{2))))$

Angivande av gränsvillkor

Ett kännetecken för inställningen och implementeringen av gränsvillkor i kontrollvolymmetoderna (inklusive Godunov-metoden) är behovet av att specificera eller beräkna flöden genom kontrollvolymytan som sammanfaller med gränsen för beräkningsdomänen. För de första och sista cellerna i beräkningsskiktet är det nödvändigt att bestämma massan, rörelsemängden och energiflödet genom ytorna.

Ofta introduceras "virtuella" beräkningsceller för att ställa in gränsvillkoren. För att göra detta introduceras ytterligare en cell till vänster om den första cellen och till höger om den sista cellen, i var och en av vilka sådana flödesparametrar specificeras så att de erforderliga flödena modelleras på sidoytan när Riemann löser problem .

Typer av gränsprocedurer

Alla antaganden görs i förhållande till den vänstra gränsen

Fast styv vägg

Huvudvillkoret är frånvaron av flöde av gasmassflödet genom gränsen, vilket motsvarar villkoret för nollflödeshastighet på den givna ytan . I den virtuella cellen måste följande flödesparametrar ställas in: $U=0$

{\begin{cases}{\begin{array}{lcl}p_{w}&=&p_{1}\\\rho _{w}&=&\rho _{1}\\u_{w}&= &-u_{1}\\\end{array}}\end{cases}}

"w" - parametrar i en virtuell cell
"1" - parametrar i den första cellen

Flödesparametrarna som erhålls i diskontinuitetsavklingningsproblemet på sidoytan realiserar ett nollmassflöde genom denna yta.

Reservoar med obegränsad kapacitet

Matematiskt motsvarar detta fall att ställa in tryckvärdet på ytan . Inflödeshastigheten kan bestämmas med formeln ${\tilde {P}}$

{\tilde {U}}=u_{1}+{\frac {{\tilde {P}}-p_{1}}{c_{1}}}

Vart i:

om , då ${\tilde {P}}>p_{1}$ $c_{1}={\frac {{\sqrt {\rho _{1}\cdot ((\gamma +1){\tilde {P))+(\gamma +1)p_{1})))} {2}}$
om , då ${\tilde {P}}<p_{1}$ $c_{1}={\frac {a_{1}\rho _{1}(\gamma -1)\left(1-{\frac {P}{p_{1))}\right)}{2\ gamma \left(1-\left({\frac {P}{p_{1))}\right)^{({\frac {\gamma -1}{2\gamma ))))\höger)))$

Inströmmande överljudsflöde

Låt understrecket beteckna parametrarna för överljudsflödet, sedan om , då ${\bar {U}}>{\bar {c}}={\sqrt {\frac {\gamma {\bar {P}}}{\bar {\rho }}}}$

{\begin{cases}{\begin{array}{lcl}p_{w}&=&{\bar {P}}\\\rho _{w}&=&{\bar {\rho }}\\ u_{w}&=&{\bar {U}}\\\end{array}}\end{cases}}

Flytt överljudsflöde

I det här fallet ställs följande flödesparametrar in i den virtuella cellen:

{\begin{cases}{\begin{array}{lcl}p_{w}&=&p_{1}\\\rho _{w}&=&\rho _{1}\\u_{w}&= &u_{1}\\\end{array}}\end{cases}}

Välja mesh-alternativ

Steget för beräkningsrutnätet längs tidskoordinaten i Godunov-metoden kan bestämmas från Courant-Friedrichs-Levy stabilitetskriteriet . När det gäller det aktuella systemet är detta villkor formulerat enligt följande:

Vågorna som uppstår i problemet med förfallet av en godtycklig diskontinuitet vid punkten får inte nå sidoytorna i tid och och förvränga den självliknande lösningen . $j+{\frac {1}{2}}$ $\Delta t$ $j+{\frac {3}{2}}$ $j-{\frac {1}{2}}$

Genomförandet av denna princip leder till följande relationer:

\Delta t_{j}=r{\frac {\Delta x_{j}}{\max \left(|D_{{j+{\frac {1}{2}}}}^{{L}}|, |D_{{j-{\frac {1}{2}}}}^{{L}}|\right)))

var

$D_{{j+{\frac {1}{2}}}}^{{L}}$ är värdet på hastigheten för vågen längst till vänster i diskontinuitetens avklingning;
$D_{{j-{\frac {1}{2}}}}^{{R}}$ är värdet på hastigheten för vågen längst till höger i diskontinuitetens avklingning;

Som ett resultat tar vi:

\Delta t=\min _{j}\Delta t_{j}

Litteratur

Numerisk lösning av flerdimensionella problem med gasdynamik. Album / redaktör Godunov S. K. . — M .: Nauka, 1976. — 400 sid. - 6500 exemplar.
Samarsky A.A., Popov Yu.P. Skillnadsmetoder för att lösa problem med gasdynamik. - M . : Nauka, 1992. - 2470 exemplar.

Länkar

Metoder för att lösa differentialekvationer

Rutnätsmetoder

Finita elementmetoder	Finita elementmetod Galerkin metod Diskontinuerlig Galerkin-metod Flerskalig finita elementmetod Multigrid-metod
Andra metoder	Ändlig skillnadsmetod Metod med ändlig volym Godunovs metod Gränselementmetod

Icke-grid-metoder