Stabilitätskriterium für Wellen in anisotropen Festkörpern

Die Bewegungsgleichungen für einen elastischen Festkörper sind gegeben durch

\begin{aligned} \nabla \cdot σ + f = ρ \ddot{u} \\ σ = C ε \\ ε = \frac{1}{2} (\nabla u + [\nabla u]^{T}) \end{aligned}

$\begin{align} &\nabla \cdot \boldsymbol{\sigma} + \mathbf{f} = \rho \ddot{\mathbf{u}}\\ &\boldsymbol{\sigma} = \mathbb{C}\boldsymbol{\varepsilon}\\ &\boldsymbol{\varepsilon} = \frac{1}{2}\left(\nabla \mathbf{u} + [\nabla\mathbf{u}]^T\right) \end{align}$

oder in Indexnotation

\begin{aligned} σ_{i j, j} + f_{i} = ρ \ddot{u_{i}} \\ σ_{i j} = C_{i j k l} ε_{k l} \\ ε = \frac{1}{2} (u_{i, j} + u_{j, i}) \end{aligned}

$\begin{align} &\sigma_{ij,j} + f_i = \rho \ddot{u_i}\\ &\sigma_{ij} = C_{ijkl}\varepsilon_{kl}\\ &\varepsilon = \frac{1}{2}(u_{i,j} + u_{j,i}) \end{align}$

$\mathbf{u}$ ist der Verschiebungsvektor, $\mathbf{f}$ ist die Körperkraft (Quellterm), $\boldsymbol{\sigma}$ ist der Spannungstensor, $\mathbf{\varepsilon}$ ist der Dehnungstensor und $\mathbb{C}$ ist der Steifheitstensor. Im Fall von isotropen Festkörpern wird der Steifheitstensor in Form von zwei verschiedenen Konstanten geschrieben . Für unbegrenzte Domänen lässt die Gleichung zwei Arten von Wellen zu, die entkoppelt sind, und das Kriterium für die Stabilität wird durch den schlimmsten Fall der beiden verschiedenen Fälle angegeben (dh der mit der höheren Geschwindigkeit).

Für transversale isotrope Materialien gibt es 5 unabhängige Parameter, die den Tensor definieren, und 3 Arten von Wellen (2 davon sind gekoppelt). Im allgemeineren Fall beträgt die Anzahl der Parameter 21 und die Welle wird gekoppelt.

Frage: Wie finden Sie das Stabilitätskriterium in einem Zeitmarschalgorithmus für den allgemeinen Fall?

pde numerical-analysis finite-difference cfl solid-mechanics nicoguaro
quelle

Wie definieren Sie "Stabilität" und nach welchen Kriterien suchen Sie?

Wolfgang Bangerth

Ich suche Stabilität in einem expliziten numerischen Schema. Das Äquivalent der CFL-Bedingung.

Nicoguaro

Antworten:

Wellengleichungen wie diese können als hyperbolisches System von Erhaltungsgesetzen erster Ordnung umgeschrieben werden:

q_{t} + \nabla \cdot F (q) = 0.

$q_t + \nabla \cdot F(q) = 0.$

Der stabile Zeitschritt für eine explizite numerische Diskretisierung hängt von der CFL-Nummer ab, die proportional zur maximalen Wellengeschwindigkeit ist, die im Problem auftritt. Diese Geschwindigkeit kann ermittelt werden, indem die Eigenwerte des Jacobi der Flussfunktion berechnet und der größte (in absoluten Werten) genommen werden. $F$

In Dimensionen, hat - Komponenten und die genaue Analyse erfordert die Eigenwerte von beliebigen linearen Kombinationen dieser Komponenten zu finden. Für die meisten Systeme, einschließlich der Elastizität, reicht es jedoch aus, die Eigenwerte jeder der Komponenten von . $d$ $F$ $d$ $F$

Eine schöne Referenz für diese Theorie ist LeVeques Text über Methoden mit endlichem Volumen . Die Elastizität wird in Kapitel 22 ausführlich behandelt.

Es gelten alle üblichen Vorbehalte bezüglich der CFL-Bedingung: Es ist eine notwendige, aber normalerweise nicht ausreichende Bedingung für die Stabilität. Die ausreichende Bedingung für die Stabilität einer gegebenen Diskretisierung ist jedoch im Allgemeinen durch die CFL-Bedingung multipliziert mit einer Konstanten gegeben. Um Ihr Stabilitätskriterium zu finden, müssen Sie sowohl die maximale Wellengeschwindigkeit (basierend auf den Gleichungen, die Sie lösen) als auch diese Konstante (basierend auf der von Ihnen verwendeten Diskretisierung) kennen.

David Ketcheson
quelle

Danke für deine Antwort. Ich habe das Kapitel in Leveques Buch studiert und die Form der Matrizen , und . Und ich habe die charakteristischen Polynome auf kubische Gleichungen reduziert, deren Lösung im Allgemeinen umständlich ist. Um Ihren Vorschlag zu überprüfen, vereinfache ich die Gleichungen für den transversalen isotropen Fall, und die Eigenwerte jedes einzelnen sind nicht die Maxima der Geschwindigkeiten. Ich habe dann die Eigenwerte für lineare Kombinationen berechnet und sie ändern sich je nach Richtung (wie im Buch beschrieben). Ist es dann notwendig, das Maximum für alle möglichen Richtungen zu finden?

A

$A$

B

$B$

C

$C$

Nicoguaro

Im Fall eines anisotropen Materials werden die Phasengeschwindigkeiten von Wellen, die sich durch dieses Material , durch die Christoffelsche Gleichung bestimmt: Hier ist die Dichte des Materials, ist die Ausbreitungsgeschwindigkeit der Welle, ist das Kronecker-Delta, ist der Einheitsvektor in der Richtung .

[ρ c^{2} δ_{i j} - C_{i j k l} n_{j} n_{l}] [u_{k}] = 0

$\begin{equation} [\rho c^2\delta_{ij} - C_{ijkl}n_jn_l][u_k]=0 \end{equation}$

ρ

$\rho$

c

$c$

δ

$\delta$

n_{j}

$n_j$

j

$j$

Es ist ersichtlich, dass die Lösung der Determinante der Klammern auf der linken Seite entspricht, was eine Eigenwertgleichung in den Eigenwerten ergibt . $\rho c^2$

Im dreidimensionalen anisotropen Fall haben wir also immer noch drei verschiedene Phasengeschwindigkeiten für eine bestimmte Ausbreitungsrichtung, von denen die größte für die CFL-Analyse verwendet werden muss, ähnlich wie die Längsgeschwindigkeit in einer isotropes Problem.

DanielRch
quelle

Danke für deine Antwort. Mir ist bekannt, dass die Christoffel-Gleichung die Informationen der Wellengeschwindigkeit enthält, obwohl ich diese hinzufügen werde, um das Maximum zu finden, das Sie für alle Richtungen lösen müssen (unter Verwendung von sphärischen Koordinaten) mit Eigenwerten der Christoffel-Gleichung. Meine Frage bezog sich eher auf die CFL-Bedingung für diesen Fall (ich war mir nicht sicher, ob nur die maximale Geschwindigkeit für das Material berechnet wurde).

max_{θ, ϕ} max_{i} λ_{i} (θ, ϕ)

$\max_{\theta, \phi} \max_i{\lambda_i(\theta, \phi)}\\$

λ_{i}

$\lambda_i$

Nicoguaro

Ich werde die Antwort von @DavidKetcheson erweitern. Zunächst werden die Gleichungen als hyperbolisches System von Erhaltungsgesetzen erster Ordnung umgeschrieben:

q_{t} + \nabla \cdot F (q) = 0

$q_t + \nabla \cdot F(q) = 0$

oder

q_{t} + A q_{x} + B q_{y} + C q_{z} = 0

$q_t + A q_x + B q_y + C q_z = 0$

Wobei ein Zustandsvektor ist, der mit den Komponenten des Spannungstensors gebildet wird und Komponenten des Geschwindigkeitsvektors . $q$ $(\sigma_{11}, \sigma_{22}, \sigma_{33}, \sigma_{12}, \sigma_{23}, \sigma_{13})$ $(u, v, w)$

q = (\begin{matrix} σ_{11} \\ σ_{22} \\ σ_{33} \\ σ_{12} \\ σ_{23} \\ σ_{13} \\ u \\ v \\ w \end{matrix}),

$q = \begin{pmatrix} \sigma_{11}\\ \sigma_{22}\\ \sigma_{33}\\ \sigma_{12}\\ \sigma_{23}\\ \sigma_{13}\\ u\\ v\\ w \end{pmatrix} \enspace ,$

A = (\begin{matrix} 0 & 0 & 0 & 0 & 0 & 0 & c_{11} & c_{16} & c_{15} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{12} & c_{26} & c_{25} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{13} & c_{36} & c_{35} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{14} & c_{46} & c_{45} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{15} & c_{56} & c_{55} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{16} & c_{66} & c_{56} \\ \frac{1}{ρ} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & \frac{1}{ρ} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{ρ} & 0 & 0 & 0 \end{matrix}),

$A = \begin{pmatrix} 0 & 0 & 0 & 0 & 0 & 0 & {c}_{11} & {c}_{16} & {c}_{15}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{12} & {c}_{26} & {c}_{25}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{13} & {c}_{36} & {c}_{35}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{14} & {c}_{46} & {c}_{45}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{15} & {c}_{56} & {c}_{55}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{16} & {c}_{66} & {c}_{56}\\ \frac{1}{\rho} & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & \frac{1}{\rho} & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & \frac{1}{\rho} & 0 & 0 & 0 \end{pmatrix} \enspace ,$

B = (\begin{matrix} 0 & 0 & 0 & 0 & 0 & 0 & c_{16} & c_{12} & c_{14} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{26} & c_{22} & c_{24} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{36} & c_{23} & c_{34} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{46} & c_{24} & c_{44} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{56} & c_{25} & c_{45} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{66} & c_{26} & c_{46} \\ 0 & 0 & 0 & \frac{1}{ρ} & 0 & 0 & 0 & 0 & 0 \\ 0 & \frac{1}{ρ} & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{ρ} & 0 & 0 & 0 & 0 \end{matrix}),

$B = \begin{pmatrix} 0 & 0 & 0 & 0 & 0 & 0 & {c}_{16} & {c}_{12} & {c}_{14}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{26} & {c}_{22} & {c}_{24}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{36} & {c}_{23} & {c}_{34}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{46} & {c}_{24} & {c}_{44}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{56} & {c}_{25} & {c}_{45}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{66} & {c}_{26} & {c}_{46}\\ 0 & 0 & 0 & \frac{1}{\rho} & 0 & 0 & 0 & 0 & 0\\ 0 & \frac{1}{\rho} & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & \frac{1}{\rho} & 0 & 0 & 0 & 0 \end{pmatrix} \enspace ,$

C = (\begin{matrix} 0 & 0 & 0 & 0 & 0 & 0 & c_{15} & c_{14} & c_{13} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{25} & c_{24} & c_{23} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{35} & c_{34} & c_{33} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{45} & c_{44} & c_{34} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{55} & c_{45} & c_{35} \\ 0 & 0 & 0 & 0 & 0 & 0 & c_{56} & c_{46} & c_{36} \\ 0 & 0 & 0 & 0 & 0 & \frac{1}{ρ} & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & \frac{1}{ρ} & 0 & 0 & 0 & 0 \\ 0 & 0 & \frac{1}{ρ} & 0 & 0 & 0 & 0 & 0 & 0 \end{matrix}) .

$C = \begin{pmatrix} 0 & 0 & 0 & 0 & 0 & 0 & {c}_{15} & {c}_{14} & {c}_{13}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{25} & {c}_{24} & {c}_{23}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{35} & {c}_{34} & {c}_{33}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{45} & {c}_{44} & {c}_{34}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{55} & {c}_{45} & {c}_{35}\\ 0 & 0 & 0 & 0 & 0 & 0 & {c}_{56} & {c}_{46} & {c}_{36}\\ 0 & 0 & 0 & 0 & 0 & \frac{1}{\rho} & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & \frac{1}{\rho} & 0 & 0 & 0 & 0\\ 0 & 0 & \frac{1}{\rho} & 0 & 0 & 0 & 0 & 0 & 0 \end{pmatrix} \enspace .$

Um die Geschwindigkeiten des Problems zu berechnen (wie oben beschrieben), müssen wir die Matrix , wobei ist ein Einheitsvektor und bestimmt die Ausbreitungsrichtung. Um die CFL-Bedingung zu finden, muss sie gelöst werden $\hat{A}(n_1, n_2, n_3) = n_1 A + n_2 B + n_3 C$ $\mathbb{n}=(n_1, n_2, n_3)$

max_{(θ, ϕ)} max_{i} γ_{i} (θ, ϕ)

$\max_{(\theta, \phi)} \max_i \gamma_i(\theta, \phi)$

Dabei sind sphärische Winkel und die Eigenwerte der Matrix . $(\theta, \phi)$ $\gamma_i$ $\hat{A}(\theta, \phi)$

Basierend darauf und der Antwort von @DavidKetcheson ist es einfacher, die Eigenwerte der Christoffel-Gleichung zu berechnen und das Optimierungsproblem zu lösen

max_{(θ, ϕ)} max_{i} λ_{i} (θ, ϕ)

$\max_{(\theta, \phi)} \max_i \lambda_i(\theta, \phi)$

mit Eigenwerten der Christoffel-Gleichung. Und die Geschwindigkeit ist nur . $\lambda_i$ $c = \sqrt{\lambda_i/\rho}$

nicoguaro
quelle