Rolle des numerischen Flusses in DG-FEM

Ich lerne die Theorie hinter den DG-FEM-Methoden anhand des Hesthaven / Warburton-Buches und bin etwas verwirrt über die Rolle des 'numerischen Flusses'. Ich entschuldige mich, wenn dies eine grundlegende Frage ist, aber ich habe eine zufriedenstellende Antwort darauf gesucht und nicht gefunden.

Betrachten Sie die lineare Skalarwellengleichung:

$$\frac{\partial u}{\partial t} + \frac{\partial f(u)}{\partial x} = 0$$ wobei der lineare Fluss gegeben ist als$f(u) = au$.

Wie in Hesthavens Buch vorgestellt, erhalten wir für jedes Element N Gleichungen, eine für jede Basisfunktion, die erzwingen, dass der Rest schwach verschwindet:$k$$N$

$$R_h(x,t) = \frac{\partial u_h}{\partial t} + \frac{\partial au_h}{\partial x}$$

$$\int_{D^k} R_h(x,t) \psi_n(x) \, dx = 0$$

Fein. Wir gehen also einmal die Integration nach Teilen durch, um zur „schwachen Form“ (1) zu gelangen, und integrieren zweimal die Integration nach Teilen, um die „starke Form“ (2) zu erhalten. Ich werde Hesthavens Art von Overkill, aber einfach zu verallgemeinernder Oberflächenintegralform in 1D übernehmen:

(1)

$$\int_{D^k} \left( \frac{\partial u_h^k}{\partial t} \psi_n-au_h^k\frac{d \psi_n}{d x} \right)\, dx = - \int_{\partial D^k} \hat{n}\cdot(au_h)^*\psi_n \,dx \qquad 1 \leq n\leq N$$

(2)

$$\int_{D^k} R_h \psi_n\, dx = \int_{\partial D^k} \hat{n}\cdot \left( au_h^k-(au_h)^* \right)\psi_n \,dx \qquad 1 \leq n\leq N$$

Warum wählen wir einen numerischen Fluss? Warum gehen wir nicht den Wert der Verwendung an der Grenze in (1) anstelle ein Flussmittel zu verwenden? Ja, es ist wahr, dass der Wert dieser Größe elementübergreifend multipliziert werden kann, aber jede Gleichung ist nur über 1 Element D k . Warum ist das wichtig?$au_h^k$$D^k$

Weiterhin ergibt die Grenz Term der zweiten partiellen Integration eindeutig eine Menge das zweite Mal in (2), die ich keinen Sinn macht. Wir machen den gleichen Vorgang! Warum werden die beiden Grenzbegriffe nicht einfach aufgehoben, wodurch (2) unbrauchbar wird? Wie haben wir neue Informationen eingeführt?$au_h^k$

Ich vermisse eindeutig etwas, das für die Methode von entscheidender Bedeutung ist, und ich möchte dieses Problem beheben. Ich habe einige reale und funktionale Analysen durchgeführt. Wenn es also eine theoretischere Antwort bezüglich der Formulierung gibt, würde ich gerne wissen!

Der numerische Fluss wird gewählt, um sicherzustellen, dass die Informationen in dem Problem in Richtung der charakteristischen Kurven der Gleichung wandern (Aufwindung). Wie in den Kommentaren erwähnt, ist der numerische Fluss erforderlich, um die für jedes Element definierten Unterprobleme zu koppeln.

Eine Möglichkeit, sich ein Bild über die Rolle des numerischen Flusses zu machen, besteht darin, das folgende einfache Beispiel zu betrachten.

Betrachten Sie die skalare Advektionsgleichung (wobei der Einfachheit halber $a=1$ )

$$\frac{\partial u}{\partial t} + \frac{\partial u}{\partial x} = 0 \qquad\text{on $\Omega$},$$ wobei die Domäne durch$\Omega = [0,1]$ . Da es sich um eine hyperbolische Gleichung handelt und sich die Informationen von links nach rechts ausbreiten, müssen wir eine Randbedingung bei$x = 0$ (aber nicht bei$x = 1$ )erzwingen. Für Konkretion, dass wir die Dirichlet Bedingung erzwingen$u(0,t) = g_D$ für einige gegeben$g_D$ .

Es sei nun angenommen diskretisieren wir diese Gleichung des DG - Verfahren verwendet wird , und wir verwenden , um zwei Elemente, $D_1 = [0,1/2]$ und $D_2 = [1/2,1]$ . Genauso gut könnten wir den folgenden Satz von zwei gekoppelten PDEs (PDE 1) diskretisieren

$$\begin{align*} \text{(PDE 1):}&& \quad v_t + v_x &= 0 \quad\text{on $D_1$},\\ \text{(PDE 2):}&& \quad w_t + w_x &= 0 \quad\text{on $D_2$}, \end{align*}$$ wo wir diese Gleichungen so koppeln, dass sie der ursprünglichen Gleichung entsprechen.

$D_1$$D_2$

$D_1$$v(0,t) = g_D$$D_2$$w(1/2,t) = v(1/2,t)$

The DG method in this case chooses the numerical fluxes precisely to enforce the above boundary conditions. If we multiply by a test function $\psi$ and integrate by parts over each element $D_k$, we obtain boundary terms of the form

$$\begin{align*} \int_{\partial D_1} \hat{n} \cdot v \psi \, dx &= \left[ v\psi \right]_0^{1/2}\\ \int_{\partial D_2} \hat{n} \cdot w \psi \, dx &= \left[ w\psi\right]_{1/2}^1 \\ \end{align*}$$ In order to "weakly" enforce the boundary conditions, we replace $v$ and $w$ with the prescribed values at those points where boundary conditions are specified (i.e. the left endpoints of $D_1$ and $D_2$). This means we replace $v(0,t)$ by $g_D$ and $w(1/2,t)$ by $v(1/2,t)$ in the boundary integrals.

In other words, we define $u_h^* = g_D$ at $x = 0$, and $u_h^* = v(1/2,t)$ at $x=1/2$, and we recover exactly the standard upwind flux that is used in the DG method.

Looking at things this way, we can consider the numerical flux functions as weakly enforcing the boundary conditions on each element that are required to couple the equations in such a way that respects the characteristic structure of the equations.

For equations more complicated than constant-coefficient advection, information may not propagate always in the same direction, and so the numerical flux must be determined by solving (or approximating the solution to) a Riemann problem at the interface. This is discussed for linear problems in Section 2.4 of Hesthaven's book.

Very loosely speaking there are two things most discretization techniques need in order to converge to the actual solution of your PDE as you increase their approximation quality, regardless if you're using DG or not:

1. Consistency (If a function $u$ satisfies the PDE, then it also satisfies your weak formulation)
2. Stability (small changes in data result in small changes in answer)

The first steps of a DG derivation where you integrate by parts on each mesh element preserves (1) because you are starting with the PDE and only applying legal operations from there.

This does not give you (2) though. You can see this yourself by trying to assemble the matrix of a partially formulated DG weak form and looking at its eigenvalues - for time dependent problem we want them all on the left half plane, but without a proper numerical flux they will be everywhere. This leads to a solution that explodes exponentially in time even if the physical problem does not.

Thus you need to add terms to your formulation so that (2) is satisfied, but without hurting (1). This is tricky to do, but not impossible. You can replace function values with cell wall averages without hurting consistency, and you can always add the jump in at cell walls without hurting consistency (because for a solution $u$ with suitable smoothness property, the jump is just 0!)

The trick is to take combinations of jumps and averages and combine them in a way that your scheme is still consistent but also stable. After that a convergence theorem usually reveals itself.

This is the basics, but you can also often bring in additional physics into the numerical flux so that it doesn't simply satisfy these mathematical requirements but also plays nicely with conservation principles.

When you choose the test function equal to the trial function in the DG method you are creating an optimization problem. That is, you have a Galerkin rather than a Petrov-Galerkin method. You are seeking the time derivatives of the trial function amplitudes that will minimise the element residual in the L2 norm, and you make this miminimisation on the assumption of a given flux function at inflow.