Laplace-Eigenmoden in einem halbkreisförmigen Bereich mit Finite-Differenzen-Methode

8

Die Berechnung der Eigenmoden einer halbkreisförmigen Membran reduziert sich auf das folgende Eigenwertproblem

2u=k2u,

wobei der interessierende Bereich ein Halbkreis ist, der durch und φ [ 0 , π ] definiert ist .r[0,1]]φ[0,π]]

Es ist angebracht, in Zylinderkoordinaten zu arbeiten, in denen der Laplace geschrieben ist als

2u=2ur2+1rur+1r22uφ2.

Randbedingungen legen den Wert von an der Grenze des Halbkreises fest, wobei u = 0 ist .uu=0

Zuerst machen wir eine Diskretisierung von mit u i j = u ( r i , φ j ) , wobei r i = ( i + 1)uuichj=u(rich,φj)undφj=(j+1rich=(ich+12)hri,j=0...N-1undhr=1/N,hr=π/N. Dies ist einzentriertesNetz.φj=(j+12)hφ i,j=0N1hr=1/N.hr=π/.N.

Wir verwenden dann eine Finite-Differenzen-Näherung für den Laplace-Wert und erhalten

2uuich+1,j- -2uich,j+uich- -1,jhr2+1ichhruich+1,j- -uich- -1,j2hr++1(ichhr)2uich,j+1- -2uich,j+uich,j- -1hφ2=k2uichj

oder

uich+1,j(1+12ich)+uich- -1,j(1- -12ich)+1ich2hφ2(uich,j+1+uich,j- -1)+uich,j(- -2- -2ich2hφ2- -k2hr2)=0.

Da unser Netz zentriert ist, müssen wir in der obigen Gleichung den folgenden Ersatz vornehmen: . Diese Ersetzung hilft uns auch, die Koordinatensingularität füri=0 zu beseitigen.ichich+12ich=0

Randbedingungen bei und r = 0 , 1 können alle mit demselben Trick behandelt werden , wobei wir an der Grenze setzenφ=0,πr=0,1

u i , j + 1 = - u i , j u i - 1 , j = - u i , j u i + 1 , j = - u i , j .

uich,j- -1=- -uich,j
uich,j+1=- -uich,j
uich- -1,j=- -uich,j
uich+1,j=- -uich,j.

Aus wir dann einen Vektor v und erhalten ein klassisches Eigenwertproblem für eine Matrix A , die sorgfältig aus den obigen Gleichungen A v = k 2 h 2 r v gebildet wirduichjvEIN

EINv=k2hr2v.

Die Matrix ist eine unsymmetrische reelle Matrix und Eigenwerte und Eigenvektoren können mit einer Routine dgeevvon LAPACK erhalten werden.

Analytische Lösungen können leicht durch das Verfahren der Trennung von Variablen erhalten werden

u(r,φ)=R.(r)Φ(φ).

Sie sind

Wobei J n ist eine (zylindrische) BesselFunktion erster Art der Ordnung n und ξ ( m ) n

u(r,φ)nm=Sünde(nφ)J.n(ξn(m)R.r),
J.nnξn(m) ist die - te Nullstelle von J n .mJ.n

Eigenwerte und Frequenzen sind

ωnm=- -k2=ξn(m)R..

r=0

Hier ist die Darstellung der analytischen Lösung für die erste Eigenfunktion:

Analytische Lösung für die erste Eigenfunktion.

Das folgende Diagramm zeigt den Vergleich der numerischen Ergebnisse für drei verschiedene Diskretisierungen, soweit meine Rechenressourcen dies zulassen.

Vergleich numerischer Lösungen bei verschiedenen Diskretisierungen $ N $.

N.2L.2(u- -ueinneinlytichceinl)/.N.20,5N.

Normalisierter absoluter Fehler als Funktion der Anzahl der Gitterpunkte $ N ^ 2

Liberias
quelle
EINhr,hφ0EIN
N.=60r=0r=0u(r=0,φ)=0
1
N.N.
N.
@ David Ist dieser einfache Konvergenztest das, was Sie gemeint haben? Ist diese (lineare) Konvergenz schnell genug? Immerhin haben wir für die erste Ableitung eine Näherung zweiter Ordnung verwendet.
Liberias

Antworten:

2

Um zumindest eine Antwort zu erhalten, scheint die obige Analyse korrekt zu sein. Randbedingungen erster Ordnung beeinflussen die Genauigkeit einer Diskretisierung zweiter Ordnung (wenn ich mich richtig an LeVeques Buch über Näherungen endlicher Differenzen erinnere; DavidKetcheson kann mir dabei helfen), so dass die Lösungen des diskreten Problems konvergieren die Lösung des analytischen Problems mit der entsprechenden Geschwindigkeit.

Geoff Oxberry
quelle
Gilbert Strang erwähnt in seinen Vorlesungen auch, dass sich Randbedingungen erster Ordnung in der gesamten Lösung ausbreiten. Eigentlich frage ich mich, ob diese Art von Verhalten in der Nähe des Ursprungs typisch für Zylinderkoordinaten ist und möglicherweise aus der Koordinatensingularität stammt.
Liberias
1

Nachdem sie den größten Teil der Analyse überflogen und sich die hübschen Bilder angesehen haben, erzählen sie eine überzeugende Geschichte: Die analytische Lösung weist keine Zirkularität auf, Ihre numerische Lösung jedoch; Der zentrale Bereich ist gekrümmt.

Ich würde mich also ganz auf Ihr Finite-Differenzen-Schema konzentrieren. Ich kann noch nicht sehen, was los ist, aber vielleicht ist ein Teil der Winkelabhängigkeit der Koordinaten nicht vollständig beseitigt.

Ich bezweifle, dass es die Randbedingung ist; Beide Lösungen scheinen an allen Grenzen Null zu haben. Es könnte jedoch denkbar sein, dass die Art der Bedingung eher eine Ableitung als eine Wertbedingung ist.

Phil H.
quelle
-3

Ich verstehe die Randbedingungen nicht.

  • Warum das Minus?
  • r=0r=1φ=0φ=π

Könnten Sie bitte mehr Details geben?

Giuseppe
quelle
u(x)u1x=- -ein- -u1x=einx=00rφφ=0,πr=0,1