3.3.3. Splines#
import numpy as np
np.set_printoptions(precision=2, linewidth=65)
import matplotlib.pyplot as plt
plt.rc('figure', dpi=150)
import seaborn as sns
sns.set()
sns.set_style('ticks')
sns.set_context("notebook", font_scale=1.2, rc={"lines.linewidth": 1.2})
Polynominterpolation versucht eine globale Modellfunktion zu finden. Jedoch eignen sich Polynome mit hohen Graden im Allgemeinen nicht für eine Interpolation vieler Punkte. Einen anderen Ansatz verfolgen Splines. Diese sind Polynomzüge, welche die einzelnen Messpunkte verbinden und deren Grad klein – typischerweise zwischen eins und drei – ist.
Definition#
Für \(\sf n+1\) Messpunkte \(\sf (x_i, y_i)\) kann eine Splinefunktion \(\sf s_k\), hier ein Polynomspline, wie folgt definiert werden.
Vorausgesetzt ist, dass die Messpunkte sortiert sind, d.h. \(\sf x_0 < x_1 < \cdots < x_n\)
für jedes \(\sf i = 0\dots n − 1\) ist \(\sf s_k\) ein Polynom vom Grad \(\sf k\) auf dem Intervall \(\sf \left[x_i , x_{i+1}\right]\)
\(\sf s_k\) ist auf \(\sf \left[x_0 , x_n \right]\) \(\sf (k − 1)\)-mal stetig differenzierbar
Beispiele:
\(\sf k = 1\): Polygonzug
\(\sf k = 3\): kubische Polynomsplines (B-Splines)
Kubische Splines#
Die in der Praxis häufig eingesetzten kubischen Polynomsplines \(\sf s_3\) (\(\sf k= 3\)) haben folgende Eigenschaften:
\(\sf s_3|\left[x_i,x_{i+1}\right] = \beta_0 +\beta_1 x + \beta_2 x^2 + \beta_3 x^3\)
\(\sf s_3\) ist zweimal stetig differenzierbar auf \(\sf \left[x_0,x_n\right]\), also insbesondere an den Stützpunkten \(\sf x_i\) der Messpunkte
Die Koeffizienten \(\sf \beta_i\) werden wie folgt bestimmt
aus den \(\sf n + 1\) Messpunkten ergeben sich \(\sf n\) Intervalle, d.h. mit jeweils vier Koeffizienten sind es insgesamt \(\sf 4n\) Koeffizienten
Exakte Darstellung der Messpunkte (\(\sf n + 1\) Gleichungen), d.h.: \(\sf s_3(x_i) = y_i\)
Glattheitsbedingungen an den inneren Messpunkten (\(\sf i=1\dots n-1\)), mit jeweils (\(\sf n − 1\) Gleichungen): $\(\sf s_3'(x_i)_- =s_3'(x_i)_+\)\( \)\(\sf s_3''(x_i)_- =s_3''(x_i)_+\)\( \)\(\sf s_3'''(x_i)_- =s_3'''(x_i)_+\)$
Damit sind es \(\sf 4n − 2\) Gleichungen für \(\sf 4n\) Koeffizienten
Um die beiden fehlenden Gleichungen zu finden bzw. zu bestimmen werden Randbedingungen oder Abschlussbedingungen benötigt. Die gängigsten Bedingungen sind:
natürliche Splines: die Krümung am Rand verschwindet, d.h.: $\(\sf s_3''(x_0) = s_3''(x_n) = 0 \)$
periodische Splines: die Steigung und Krümung ist an beiden Rändern gleich $\(\sf s_3'(x_0) = s_3'(x_n)\)\( \)\(\sf s_3''(x_0) = s_3''(x_n)\)$
Hermite Splines: die Steigungen am Rand werden explizit vorgegeben (hier durch \(\sf u\) und \(\sf v\)) $\(\sf s_3'(x_0) = u \)\( \)\(\sf s_3'(x_n) = v \)$
Anwendung#
Im Folgenden werden zwei Beispiele, \(\sf s_1\) und \(\sf s_3\), für die Erstellung von Splines mit Python vorgestellt.
# Erzeugung von Messpunkten
n = 7
xi = np.linspace(0, np.pi, n)
yi = np.sin(xi)
Für die \(\sf s_1\) Splines, kann die Funktion np.interp
verwendet werden. Sie führt eine lineare Interpolation zwischen gegebenen Wertepaaren durch.
# Wertebereich für die Visualisierung der Interpolation
x = np.linspace(0, np.pi, n*6)
y = np.sin(x)
# Interpolation
y_s1 = np.interp(x, xi, yi)
plt.plot(x,y, alpha=0.3, color='C2', lw=5,
label='Generierende Funktion')
plt.plot(x, y_s1, color='C0', label='Interpolation')
plt.scatter(x, y_s1, s=3, zorder=3, color='C0')
plt.scatter(xi, yi, color='C1', label='Messpunkte')
plt.xlabel('x')
plt.ylabel('y')
plt.legend();
Die \(\sf s_3\) Splines können mit Funktionen aus dem scipy-Modul berechnet werden. Dazu werden zunächst die Koeffizienten bestimmt (scipy.interpolate.splrep
) und diese ermöglichen die gewünschte Auswertung, welche mit der Funktion scipy.interpolate.splev
vorgenommen werden kann.
import scipy.interpolate as si
s3 = si.splrep(xi, yi)
y_s3 = si.splev(x, s3)
plt.plot(x,y, alpha=0.3, color='C2', lw=5,
label='Generierende Funktion')
plt.plot(x, y_s3, color='C0', label='Interpolation')
plt.scatter(x, y_s3, s=3, zorder=3, color='C0')
plt.scatter(xi, yi, color='C1', label='Messpunkte')
plt.xlabel('x')
plt.ylabel('y')
plt.legend();