Kubische Splines LGS Matrix

Einen optimierten Algorithmus zur Splineberechnung kann man z.B. bei Arndt Brünner nachlesen und meine Umsetzung für ggb Javascript verwenden. In diesem Beispiel geht es um die Darstellung der grundsätzlichen Überlegungen einen kubischen Spline über abschnittsweise definierte kubische Polynome zu beschreiben. Im CAS erhalte ich darüberhinaus exakte Ergebnisse! 4 Punkte A0...A3 werden durch 3 kubische Parabeln verbunden Pro Parabel sind 4 Koeffizienten aij zu bestimmen - insgesamt sind 3*4=12 Gleichungen zu finden und das Gleichungssystem zu lösen: 2 Gleichungen (5.2), (5.3) (Zeile.Gleichg) 2 Gleichungen (6.1), (6.2) 2 Gleichungen (7.1), (7.2) Für "glatte" Übergänge an den Stützstellen müssen Steigung und Krümmung der Parabeln gleich sein 2 Gleichungen (5.4), (5.5) 2 Gleichungen (6.3), (6.4) Randbedingungen (für natürliche Splines) 2 Gleichungen (5.1), (7.3) (8) LGS gesamtes Gleichungssystem (9) Lösung LGS (10)(11)(12) Koeffizienten in Parabelgleichung eingesetzt (13) Koeffizientenvektor a (14) Matrix für LGS in Matrixform A a = b

Das LGS nach Funktion geordnet

In dieser Form lässt sich der Algorithmus wohl einfacher verallgemeinern. Die Gleichungen erstelle ich durch User-Def Funktionen: z.B. 1. Polynom zwischen A_0...A_1:
p(k,x_i) Splinepolynom k an der Stützstelle x_i
dp(k,x_i) 1. Ableitung Splinepolynom k an der Stützstelle x_i
ddp(k,x_i) 2. Ableitung Splinepolynom k an der Stützstelle x_i
Anpassung an unterschiedliche viele Stützpunkte A_0...A_9
Lege eine Liste der Stützpunkte an L:={A_0,A_1,A_2,A_3,A_4,A_5};
n Splinepolynomen:=Length(L)-1
Wandle Punkte A_i in Liste S S:=Sequence(Flatten(Vector(Element(L, k))),k,1,n+1);
Abgreifen der Koordinaten der Punkte XY(i,1) = x(A_i) und XY(i,2) = y(A_i) i=0..nXY(ii,jj):=Element(S, ii+1,jj);
Beispiel LGS für 5 Stützstellen A_0...A_4
11 {ddp(1,x(A_0))=0, Randbedingung A_0 2. Ableitung=0A1:={ddp(1,XY(0,1))=0}
12 13 p(1,x(A_0))=y(A_0), p(1,x(A_1))=y(A_1), Spline-Polynome gehen durch die Stützstellen A_0...A_nA2:=Sequence(p(kk,XY(kk-1,1))=XY(kk-1,2),kk,1,n ) A3:=Sequence(p(kk,XY(kk,1))=XY(kk,2),kk,1,n )
21 22p(2,x(A_1))=y(A_1), p(2,x(A_2))=y(A_2)
31 32p(3,x(A_2))=y(A_2), p(3,x(A_3))=y(A_3),
41 42p(4,x(A_3))=y(A_3), p(4,x(A_4))=y(A_4),
14 23 33dp(1,x(A_1))-dp(2,x(A_1))=0, dp(2,x(A_2))-dp(3,x(A_2))=0, dp(3,x(A_3))-dp(4,x(A_3))=0,An den Stützstellen A_1...A_n-1 gleiche Steigung 1. Ableitung gleichA4:=Sequence(dp(kk,XY(kk,1)) ,kk,1,n-1)-Sequence(dp(kk+1,XY(kk,1)),kk,1,n-1) 
15 24 34ddp(1,x(A_1))-ddp(2,x(A_1))=0, ddp(2,x(A_2))-ddp(3,x(A_2))=0, ddp(3,x(A_3))-ddp(4,x(A_3))=0,An den Stützstellen A_1...A_n-1 gleich Krümmung 2. Ableitung gleichA5:=Sequence(ddp(kk,XY(kk,1)) ,kk,1,n-1)-Sequence(ddp(kk+1,XY(kk,1)),kk,1,n-1)
43ddp(4,x(A_4))=0}Randbedingung A_n 2. Ableitung = 0A6:={ddp(n,XY(n,1))=0}
Gleichungssystem lösen Splinepolynome erstellen SPL:=Solve(Join(A1,A2,A3,A4,A5,A6),Take(Pa,1,4n)) f_i(x):=Substitute(p(i, x),SPL)
In Algebra View (Splineabschnitte zusammensetzen) If(x < x(A_1), f1(x), If(x ≤ x(A_2), f2(x), If(x ≤ x(A_3), f3(x), If(x ≤ x(A_4), f4(x), f5(x))))) http://www.tm-mathe.de/Themen/html/funnatsplines.html