Die Community zu .NET und Classic VB.
Menü

Der Hobby-Algorithmus

 von 

Einführung 

In der zweidimensionalen Computergrafik stehen wir häufig vor dem Problem, durch eine gegebene Anzahl von Punkten eine glatte Kurve zu legen. Meist werden hierfür Splines verwendet. Dabei handelt es sich um miteinander verbundene Kurvenstücke. Besonders beliebt sind kubische Bézierkurven, da sie durch vier Kontrollpunkte mit einfacher geometrischer Bedeutung definiert sind und mit Hilfe des de-Casteljau-Algorithmus effizient und numerisch stabil berechnet werden können. In Abbildung 1 sehen wir zwei verbundene Bézierkurven, die beispielsweise Teil eines Béziersplines sein könnten. Die Bézierkurve auf der rechten Seite schneidet die beiden Knotenpunkte Latex: z_i und Latex: z_{i+1}. Die Tangenten an den Schnittpunkten schneiden die zugehörigen Kontrollpunkte Latex: p_i und Latex: q_i. Der Abstand von Latex: p_i und Latex: q_i von Latex: z_i bzw. Latex: z_{i+1} gibt grob gesprochen an, wie sehr sich die Kurve an die Tangenten anschmiegt.


Abbildung 1: Winkel und Kontrollpunkte von Bézierkurven

Manchmal wollen wir Béziersplines durch eine Anzahl Knotenpunkte zeichnen, ohne die anderen Kontrollpunkte explizit gegeben zu haben. Vielmehr soll das resultierende Spline die Knotenpunkte so natürlich wie möglich, d. h. mit möglichst geringer Gesamtkrümmung, interpolieren. Dies mathematisch zu quantifizieren, ist ziemlich kompliziert. Ich möchte hier lediglich als Ergebnis den Hobby-Algorithmus vorstellen. Dieser wurde von dem amerikanischen Mathematiker John D. Hobby im Jahre 1985 entwickelt und von Donald E. Knuth für sein Schriftdefinitionssystem METAFONT, welches in TeX und LaTeX benutzt wird, eingesetzt. Verwendung findet der Algorithmus ebenfalls in den auf METAFONT basierenden Vektorgrafiksprachen METAPOST und Asymptote.

Definitionen  

Sei Latex: n \in \mathbb N mit Latex: n > 1. Gegeben seien Latex: n+1 Knotenpunkte Latex: (x_0, y_0), \ldots, (x_n, y_n), durch die wir ein kubisches Bézierspline mit durch den Hobby-Algorithmus ermittelten Kontrollpunkten legen wollen. Zur notationellen Vereinfachung definieren wir:

Latex: \vec z_i : \! {{=}} \, (x_i, y_i)^T \text{ fuer } 0 \le i \le n

Latex: \delta x_i : \! {{=}} \, x_{i+1} - x_i \text{ fuer } 0 \le i < n

Latex: \delta y_i : \! {{=}} \, y_{i+1} - y_i \text{ fuer } 0 \le i < n

Latex: \vec\delta_i : \! {{=}} \, \vec z_{i+1} - \vec z_i \text{ fuer } 0 \le i < n

Latex: d_i : \! {{=}} \, \| \vec \delta_i \| \text{ fuer } 0 \le i < n

Dabei bezeichnen Latex: (\vec v)^T die Transposition des Vektors Latex: \vec v (sodass die Latex: \vec z_i Spaltenvektoren im Latex: \mathbb R^2 sind) und Latex: \|\!\cdot\!\| die übliche euklidische Norm.
Grundsätzlich lassen sich zwei Typen von Splines unterscheiden: offene und geschlossene. Geschlossene Splines zeichnen sich gegenüber offenen dadurch aus, dass der letzte Punkt Latex: \vec z_n durch ein weiteres Kurvenstück wieder mit dem ersten Punkt Latex: \vec z_0 verbunden ist. Geschlossene Splines heißen auch zyklisch. Bezüglich des Hobby-Algorithmus unterscheiden sich offene und geschlossene Splines in einigen wichtigen Details. Das Verhalten an den inneren Knotenpunkten Latex: \vec z_1, \ldots, \vec z_{n-1} ist jedoch genau gleich, sodass wir zunächst diese Punkte betrachten und erst später zwischen den beiden Splinetypen differenzieren.

Die inneren Knotenpunkte  

Zunächst führen wir für jeden inneren Knotenpunkt die in Abbildung 1 eingezeichneten Winkel ein. Latex: \psi_i ist der Winkel zwischen der Verbindungsstrecke zwischen Latex: \vec z_i und Latex: \vec z_{i+1} sowie der Verbindungsstrecke zwischen Latex: \vec z_{i-1} und Latex: \vec z_i, d. h. zwischen den Vektoren Latex: \vec\delta_i und Latex: \vec\delta_{i-1}. Die Winkel Latex: \psi_i werden dabei so gewählt, dass Latex: -\pi < \psi_i \le \pi gilt. Latex: \theta_i ist der Winkel zwischen Latex: \vec\delta_i und der Verbindungsstrecke von Latex: \vec z_i und Latex: \vec p_i. Schließlich ist Latex: \phi_i der Winkel zwischen der Verbindungsstrecke von Latex: \vec z_{i-1} und Latex: \vec z_i und der Verbindungsstrecke von Latex: \vec z_i und Latex: \vec q_{i-1}. Zu beachten ist, dass auf diese Weise auch Latex: \theta_0 und Latex: \phi_n, nicht aber Latex: \theta_n, Latex: \phi_0, Latex: \psi_0 und Latex: \psi_n definiert werden können.
Damit der Spline bei Latex: \vec z_i stetig differenzierbar ist, müssen Latex: \vec q_{i-1}, Latex: \vec z_i und Latex: \vec p_i auf einer Geraden liegen, und zwar derart, dass Latex: \vec z_i zwischen Latex: \vec q_{i-1} und Latex: \vec p_i liegt. Daraus folgt die Bedingung, dass Latex: \phi_i, Latex: \theta_i und Latex: \psi_i sich zum Vollwinkel bzw. Nullwinkel ergänzen. Diese Bedingung werden wir später benötigen.
Durch die Winkel Latex: \phi_i und Latex: \theta_i ist zwar die Richtung der Kontrollpunkte Latex: \vec p_i und Latex: \vec q_{i-1} bezüglich Latex: \vec z_i gegeben, nicht aber ihr Abstand von Latex: \vec z_i. Wir benötigen also pro innerem Knotenpunkt noch zwei weitere Parameter. Diese sollen mit der lokalen Spannungsenergie des Splines an den Kontrollpunkten in Verbindung stehen und werden deshalb Spannungsparameter genannt. Latex: \tau_i soll hierbei den Spannungsparameter des bei Latex: \vec z_i auslaufenden Kurvensegments bezeichnen, Latex: \bar\tau_i entsprechend den des einlaufenden Segments. Im Hobby-Algorithmus werden die Spannungsparameter so definiert, dass sich minimale Gesamtspannungsenergie (und damit natürlichstes Aussehen) für Latex: \tau_i = \bar\tau_i = 1 ergibt.
Zunächst betrachten wir lediglich eine einzelne Bézierkurve zwischen den Punkten Latex: \vec z_i und Latex: \vec z_{i+1}. Wir können das Problem vereinfachen, indem wir ausnutzen, dass bei einer gegebenen affinen Transformation Latex: \mathcal T diejenige Bézierkurve, die durch die mit Latex: \mathcal T transformierten Kontrollpunkte definiert ist, identisch mit der durch Latex: \mathcal T transformierten Bézierkurve ist. Mathematisch formuliert bedeutet dies, dass für die Bézierkurve Latex: \mathfrak B(\vec a, \vec b, \vec c, \vec d) mit den Kontrollpunkten Latex: \vec a, Latex: \vec b, Latex: \vec c und Latex: \vec d gilt:

Latex: \mathcal T(\mathfrak B(\vec a, \vec b, \vec c, \vec d))= \mathfrak B(\mathcal T(\vec a), \mathcal T(\vec b), \mathcal T(\vec c), \mathcal T(\vec d))

Daraus folgern wir, dass wir eine Bézierkurve mit den Knotenpunkten Latex: \vec{\tilde z}_i = (0, 0)^T und Latex: \vec{\tilde z}_{i+1} = (1, 0)^T betrachten und die ermittelten Kontrollpunkte mit der affinen winkeltreuen Transformation Latex: \mathcal T_i, die Latex: \vec{\tilde z}_i auf Latex: \vec z_i und Latex: \vec{\tilde z}_{i+1} auf Latex: \vec z_{i+1} abbildet, transformieren können. Latex: \mathcal T_i ist eindeutig bestimmt durch

Latex: \mathcal T_i(\vec v) = \vec z_i + \begin{pmatrix} \delta x_i & -\delta y_i \\ \delta y_i & \delta x_i \end{pmatrix} \vec v.

Für die zwei untransformierten Kontrollpunkte Latex: \vec{\tilde p}_i und Latex: \vec{\tilde q}_i gilt dann:

Latex: \vec{\tilde p}_i = \vec{\tilde z}_i + \xi_i \cdot \begin{pmatrix} \cos\theta_i \\ \sin\theta_i \end{pmatrix} = \begin{pmatrix} \xi_i \cos\theta_i \\ \xi_i \sin\theta_i \end{pmatrix},

Latex: \vec{\tilde q}_i = \vec{\tilde z}_{i+1} + \zeta_i \cdot \begin{pmatrix} -\cos\phi_{i+1} \\ \sin\phi_{i+1} \end{pmatrix} = \begin{pmatrix} 1 - \zeta_i \cos\phi_{i+1} \\ \zeta_i \sin\phi_{i+1} \end{pmatrix}.

Dabei sind Latex: \xi_i der Abstand zwischen Latex: \vec z_i und Latex: \vec p_i sowie Latex: \zeta_i der Abstand zwischen Latex: \vec z_{i+1} und Latex: \vec q_i. Hat man Latex: \vec{\tilde p}_i und Latex: \vec{\tilde q}_i erst einmal gefunden, ergeben sich die gesuchten Kontrollpunkte Latex: \vec p_i und Latex: \vec q_i einfach durch Anwendung von Latex: \mathcal T_i:

Latex: \vec p_i = \mathcal T_i(\vec{\tilde p}_i),

Latex: \vec q_i = \mathcal T_i(\vec{\tilde q}_i).

Die eigentliche Schwierigkeit am Hobby-Algorithmus ist die Ermittlung von Latex: \xi_i, Latex: \zeta_i, Latex: \phi_i und Latex: \theta_i. Hobby entwickelt in seiner Publikation zunächst eine Methode, um Latex: \xi_i und Latex: \zeta_i aus den Winkeln zu berechnen. Es gilt:

Latex: \xi_i :\!= \frac{\rho(\theta_i, \phi_{i+1})}{3 \tau_i},

Latex: \zeta_i :\!= \frac{\sigma(\theta_i, \phi_{i+1})}{3 \bar\tau_{i+1}}.

Entscheidend für die Qualität des Ergebnisses ist die Auswahl der Funktionen Latex: \rho und Latex: \sigma.
Hobby wählt:

Latex: \rho(\theta, \phi) :\!= \frac{2 + \alpha(\theta, \phi)}{1 + (1-c) \cos\theta + c \cos\phi},

Latex: \sigma(\theta, \phi) :\!= \frac{2 - \alpha(\theta, \phi)}{1 + (1-c) \cos\phi + c \cos\theta}

mit

Latex: \alpha(\theta, \phi) :\!= a \cdot (\sin\theta - b \sin\phi) \cdot (\sin\phi - b \cos\theta) \cdot (\cos\theta - \cos\phi).

Hierbei sind Latex: a, Latex: b und Latex: c experimentell zu ermittelnde Konstanten; METAFONT und Asymptote wählen:

Latex: a :\!= \sqrt2,

Latex: b :\!= \frac{1}{16},

Latex: c :\!= \frac{3 - \sqrt5}{2}.

Damit ist das Problem zumindest für die inneren Knotenpunkte bis auf die Bestimmung der Winkel Latex: phi_i und Latex: theta_i gelöst. Für letzteres führt Hobby eine Krümmungsfunktion k ein:

Latex: k(\theta, \phi, \tau, \bar\tau) :\!= \frac{2 \tau^2}{\rho(\theta, \phi)^2} \left (\frac{\sigma(\theta, \phi) \sin(\theta + \phi)}{\bar\tau} - 3 \sin\theta \right )

Diese Funktion ist jedoch nichtlinear in den Variablen Latex: \theta und Latex: \phi, was auf ein schwierig zu lösendes nichtlineares Gleichungssystem für Latex: \theta_i und Latex: \phi_i hinausliefe. Um dies zu umgehen, wird Latex: k nach Latex: (\theta, \phi) entwickelt und nur der lineare Term Latex: \hat k(\theta, \phi, \tau, \bar\tau) benutzt. Man erhält unter Benutzung der linearen Näherungen Latex: \sin\theta \approx \theta und Latex: \cos\theta \approx 1:

Latex: \alpha(\theta, \phi) \approx a \cdot (\theta - b \phi) \cdot (\phi - b) \cdot (1 - 1) = 0,

Latex: \rho(\theta, \phi) \approx \frac{2}{1 + 1 - c + c} = 1,

Latex: \sigma(\theta, \phi) \approx \frac{2}{1 + 1 - c + c} = 1,

Latex: \hat k(\theta, \phi, \tau, \bar\tau) = 2 \tau^2 \left (\frac{\theta + \phi}{\bar\tau} - 3 \theta \right ).

Diese Linearisierung führt dazu, dass die Interpolation bei sehr großen Winkeln Latex: \theta, Latex: \phi nicht mehr optimal ist. Derartig große Winkel kommen jedoch in der Praxis relativ selten vor, sodass die linearisierte Krümmungsfunktion für unsere Zwecke eine mehr als ausreichende Näherung darstellt.
Damit die Krümmung an den Knotenpunkten übereinstimmt, muss

Latex: \frac{k(\theta_i, \phi_{i+1}, \tau_i, \bar\tau_{i+1})}{d_i} = \frac{k(\phi_i, \theta_{i-1}, \bar\tau_i, \tau_{i-1})}{d_{i-1}}

gelten. Aus der Stetigkeit der linearisierten Krümmungsfunktion ergibt sich wiederum:

Latex: \frac{\hat k(\theta_i, \phi_{i+1}, \tau_i, \bar\tau_{i+1})}{d_i} = \frac{\hat k(\phi_i, \theta_{i-1}, \bar\tau_i, \tau_{i-1})}{d_{i-1}}

Latex: \frac{\tau_i^2}{d_i} \left (\frac{\theta_i + \phi_{i+1}}{\bar\tau_{i+1}} - 3 \theta_i \right ) = \frac{\bar\tau_i^2}{d_{i-1}} \left (\frac{\phi_i + \theta_{i-1}}{\tau_{i-1}} - 3 \phi_i \right )

Latex: \frac{\tau_i^2}{d_i \bar\tau_{i+1}} \phi_{i+1} + \frac{\tau_i^2}{d_i} \left (\frac{1}{\bar\tau_{i+1}} - 3 \right ) \theta_i = \frac{\bar\tau_i^2}{d_{i-1} \tau_{i-1}} \theta_{i-1} + \frac{\bar\tau_i^2}{d_{i-1}} \left (\frac{1}{\tau_{i-1}} - 3 \right ) \phi_i

Latex: d_{i-1} \tau_{i-1} \tau_i^2 \phi_{i+1} + d_{i-1} \tau_{i-1} \tau_i^2 \left (1 - 3 \bar\tau_{i+1} \right ) \theta_i = d_i \bar\tau_{i+1} \bar\tau_i^2 \theta_{i-1} + d_i \bar\tau_{i+1} \bar\tau_i^2 \left (1 - 3 \tau_{i-1} \right ) \phi_i

Mit den Abkürzungen

Latex: \lambda_i :\!= 3 \tau_i - 1,

Latex: \mu_i :\!= 3 \bar\tau_i - 1,

Latex: C_i :\!= d_{i-1} \tau_i^2 \tau_{i-1},

Latex: A_i :\!= d_i \bar\tau_i^2 \bar\tau_{i+1}

vereinfacht sich die letzte Gleichung zu:

Latex: A_i \theta_{i-1} + C_i \mu_i \theta_i - A_i \lambda_i \phi_i - C_i \phi_{i+1} = 0

Weiterhin gilt, wie oben erwähnt, aufgrund der stetigen Differenzierbarkeit an den Knotenpunkten:

Latex: \phi_i = -\psi_i - \theta_i

Dies alles gilt natürlich nur an den inneren Knotenpunkten. Wir haben aus der Stetigkeit der Krümmungsfunktion Latex: n-1 und aus der stetigen Differenzierbarkeit des Splines ebenfalls Latex: n-1, insgesamt also Latex: 2n-2 Bedingungen erhalten. Für die Latex: n Kurvensegmente benötigen wir je zwei zusätzliche Kontrollpunkte. Die Abstände Latex: \xi_i und Latex: \zeta_i der Kontrollpunkte von den Knotenpunkten lässt sich wie oben gezeigt aus den Winkeln Latex: \theta_i und Latex: \phi_i berechnen, somit benötigen wir pro Kontrollpunkt nur eine Gleichung (für den jeweils anliegenden Winkel). Dennoch fehlen immer noch zwei Bedingungen, für die wir die äußeren Knotenpunkte Latex: \vec z_0 und Latex: \vec z_n betrachten und erstmals zwischen offenen und geschlossenen Splines unterscheiden müssen.

Offene Splines  

Hier müssen wir zwei weitere Gleichungen künstlich einführen, da wir keine periodischen Randbedingungen haben. Analog zu den Spannungsparametern definiert Hobby die Wirbelparameter Latex: chi_0 und Latex: \chi_n derart, dass sich für Latex: \chi_0 = \chi_n = 1 ein natürliches Aussehen ergibt. Die Wirbelparameter beeinflussen das Aussehen des Splines an den Endpunkten. Für Latex: \vec z_0 gibt es darüberhinaus einen weiteren Ausgangsspannungsparameter Latex: \tau_0; analog existiert für Latex: \vec z_n ein Eingangsspannungsparameter Latex: \bar\tau_n. Man wählt die folgenden Randbedingungen:

Latex: \hat k(\theta_0, \phi_1, \tau_0, \bar\tau_1) = \chi_0 \,\hat k(\phi_1, \theta_0, \bar\tau_1, \tau_0) \iff \theta_0 = F_1 \phi_1,

Latex: \hat k(\phi_n, \theta_{n-1}, \bar\tau_n, \tau_{n-1}) = \chi_n \,\hat k(\theta_{n-1}, \phi_n, \tau_{n-1}, \bar\tau_n) \iff \phi_n = G_{n-1} \theta_{n-1}.

mit

Latex: F_1 :\!= \frac{\tau_0^3 + \lambda_0 \chi_0 \bar\tau_1^3}{\mu_1 \tau_0^3 + \chi_0 \bar\tau_1^3},

Latex: G_{n-1} :\!= \frac{\bar\tau_n^3 + \mu_n \chi_n \tau_{n-1}^3}{\lambda_{n-1} \bar\tau_n^3 + \chi_n \tau_{n-1}^3}.

Damit lassen sich Latex: \theta_0 und Latex: \phi_n aus den anderen Gleichungen eliminieren:

Latex: 0 = C_1 \mu_1 \theta_1 + A_1 \left (F_1 - \lambda_1 \right ) \phi_1 - C_1 \phi_2 Latex: = \left (A_1 \left (\lambda_1 - F_1 \right ) + C_1 \mu_1 \right ) \theta_1 - C_1 \phi_2 + A_1 \left (\lambda_1 - F_1 \right ) \psi_1,

Latex: 0 = A_{n-1} \theta_{n-2} + C_{n-1} \left (\mu_{n-1} - G_{n-1} \right ) \theta_{n-1} - A_{n-1} \lambda_{n-1} \phi_{n-1} Latex: = A_{n-1} \theta_{n-2} + \left (C_{n-1} \left (\mu_{n-1} - G_{n-1} \right ) + A_{n-1} \lambda_{n-1} \right ) \theta_{n-1} + A_{n-1} \lambda_{n-1} \psi_{n-1}.

Für Latex: n>2 lässt sich Latex: \phi_2 ebenfalls durch Latex: \theta_2 ersetzen:

Latex: \left (C_1 \mu_1 + A_1 \left (\lambda_1 - F_1 \right ) \right ) \theta_1 + C_1 \theta_2 + A_1 \left (\lambda_1 - F_1 \right ) \psi_1 + C_1 \psi_2 = 0

Weiterhin gilt für Latex: 1 < i < n-1:

Latex: A_i \theta_{i-1} + (A_i \lambda_i + C_i \mu_i) \theta_i + C_i \theta_{i+1} + A_i \lambda_i \psi_i + C_i \psi_{i+1} = 0.

Definiert man nun

Latex: F_i :\!= 0 \qquad (1 < i < n),

Latex: G_i :\!= 0 \qquad (0 < i < n-1),

Latex: \psi_n :\!= 0,

Latex: Q_i :\!= A_i \left (\lambda_i - F_i \right ),

Latex: R_i :\!= C_i \left (\mu_i - G_i \right ),

Latex: B_i :\!= Q_i + R_i,

Latex: D_i :\!= -Q_i \psi_i - C_i \psi_{i+1},

so gilt:

Latex: B_1 \theta_1 + C_1 \theta_2 = D_1,

Latex: A_i \theta_{i-1} + B_i \theta_i + C_i \theta_{i+1} = D_i \qquad (i < 1 < n-1),

Latex: A_n \theta_{n-1} + B_n \theta_n = D_n.

Dies ist ein ein lineares tridiagonales Gleichungssystem der Form

Latex: \begin{pmatrix} B_1 & C_1 & & \cdots & & & \\ A_2 & B_2 & C_2 & \cdots & & & \\ & & & \ddots & & & \\ & & & \cdots & A_{n-2} & B_{n-2} & C_{n-2} \\ & & & \cdots & & A_{n-1} & B_{n-1} \end{pmatrix} Latex: \begin{pmatrix} \theta_1 \\ \theta_2 \\ \vdots \\ \theta_{n-2} \\ \theta_{n-1} \end{pmatrix} Latex: = \begin{pmatrix} D_1 \\ D_2 \\ \vdots \\ D_{n-2} \\ D_{n-1} \end{pmatrix}

mit Subdiagonale Latex: (A_i), Diagonale Latex: (B_i), Superdiagonale Latex: (C_i) und rechtsseitigem Vektor Latex: (D_i), welches mit dem Gaußalgorithmus in linearer Zeit gelöst werden kann.
Sind alle Latex: \tau_i und Latex: \bar\tau_i größer als Latex: \frac23, so gilt wegen Latex: C_i > 0, Latex: A_i > 0, Latex: \lambda_i > 1 und Latex: \mu_i > 1 an den inneren Knotenpunkten:

Latex: \left | B_i \right | = \lambda_i A_i + \mu_i C_i > A_i + C_i = \left | A_i \right | + \left | C_i \right |.

An den Randpunkten folgt analog, falls zusätzlich Latex: \lambda_1 > F_1 und Latex: \mu_{n-1} > G_{n-1} gelten:

Latex: \left | B_1 \right | = A_1 \left (\lambda_1 - F_1 \right ) + C_1 \mu_1 > \left | C_1 \right |,

Latex: \left | B_{n-1} \right | = A_{n-1} \lambda_{n-1} + C_{n-1} \left (\mu_{n-1} - G_{n-1} \right ) > \left | A_{n-1} \right |.

Somit ist die Systemmatrix diagonaldominant. Das führt dazu, dass in diesem (in der Realität meistens angenommenen) Fall keine Pivotisierung nötig ist.

Geschlossene Splines  

Hier haben wir periodische Randbedingungen vorliegen, weshalb es sich anbietet, die Indizes zyklisch zu definieren:

Latex: \vec\delta_n :\!= z_0 - z_n,

Latex: d_n :\!= \| \vec\delta_n \| ,

Latex: d_{-1} :\!= d_n,

Latex: d_{n+1} :\!= d_0,

Latex: \theta_{-1} :\!= \theta_n,

Latex: \theta_{n+1} :\!= \theta_0,

Latex: \phi_{-1} :\!= \phi_n,

Latex: \phi_{n+1} :\!= \phi_0,

Latex: \tau_{-1} :\!= \tau_n,

Latex: \bar\tau_{n+1} :\!= \bar\tau_0.

Weiterhin definieren wir Latex: \psi_0 und Latex: \psi_{n+1} als den Winkel zwischen Latex: \vec\delta_n und Latex: \vec\delta_0 sowie Latex: \psi_n als den Winkel zwischen Latex: \vec\delta_{n-1} und Latex: \vec\delta_n. Damit können wir die Gleichungen Kruemmung und Differenzierbarkeit auch für Latex: i = 0 und Latex: i = n anwenden. Mit den Definitionen

Latex: Q_i :\!= A_i \lambda_i,

Latex: R_i :\!= C_i \mu_i,

Latex: B_i :\!= Q_i + R_i,

Latex: D_i :\!= -Q_i \psi_i - C_i \psi_{i+1}

erhalten wir wieder

Latex: A_i \theta_{i-1} + B_i \theta_i + C_i \theta_{i+1} = D_i.

Für Latex: i = 0 und Latex: i = n gelten mit obigen Definitionen die Spezialfälle:

Latex: B_0 \theta_0 + C_0 \theta_1 + A_0 \theta_n = D_0,

Latex: C_n \theta_0 + A_n \theta_{n-1} + B_n \theta_n = D_n.

Aufgrund der Periodizität des Problems ist die Gleichung für Latex: i = n+1 identisch zu der Gleichung für Latex: i = 0. Damit haben wir also für die Latex: n+1 Winkel Latex: \theta_0, \ldots, \theta_n (ein Winkel pro Knotenpunkt) Latex: n+1 Gleichungen gefunden, die nun aber kein einfaches tridiagonales Gleichungssystem mehr, sondern ein zyklisches tridiagonales Gleichungssystem der Form

Latex: \begin{pmatrix} B_0 & C_0 & & \cdots & & & A_0 \\ A_1 & B_1 & C_1 & \cdots & & & \\ & & & \ddots & & & \\ & & & \cdots & A_{n-1} & B_{n-1} & C_{n-1} \\ C_n & & & \cdots & & A_n & B_n \end{pmatrix} Latex: \begin{pmatrix} \theta_0 \\ \theta_1 \\ \vdots \\ \theta_{n-1} \\ \theta_n \end{pmatrix} Latex: = \begin{pmatrix} D_0 \\ D_1 \\ \vdots \\ D_{n-1} \\ D_n \end{pmatrix}

bilden. Unter Zuhilfenahme der Sherman-Morrison-Formel kann auch dieses System mit linearer Zeitkomplexität gelöst werden.

Quellen  

  1. John D. Hobby: Smooth, Easy to Compute Interpolating Splines, 1985
    ftp://db.stanford.edu/pub/cstr/reports/cs/tr/85/1047/
  2. Martin Kerscher, Joachim Puls, Sigmung Stintzing: Numerik für Physiker, 2007
    http://www.mathematik.uni-muenchen.de/~kerscher/vorlesungen/numeriksose08/doc/script.pdf
  3. William H. Press, Saul A. Teukolsky, William T. Vetterling, Brian P. Flannery:
    Numerical Recipes in C++, Second Edition, Cambridge University Press

Downloads  

Beispielprojekt [282559 Bytes]


Artikel als PDF [189103 Bytes]

Ihre Meinung  

Falls Sie Fragen zu diesem Tutorial haben oder Ihre Erfahrung mit anderen Nutzern austauschen möchten, dann teilen Sie uns diese bitte in einem der unten vorhandenen Themen oder über einen neuen Beitrag mit. Hierzu können sie einfach einen Beitrag in einem zum Thema passenden Forum anlegen, welcher automatisch mit dieser Seite verknüpft wird.