Calcolo numerico/Interpolazione polinomiale a tratti (spline)

Definizione di spline

$s$ è definita come

$s(x)=\{\begin{array}{cc}{s}_0(x)& ⟺x\in [x_0,x_1]\\ s_1(x)& ⟺x\in [x_1,x_2]\\ s_{n-1}(x)& ⟺x\in [x_{n-1},x_n]\end{array}$

dove le $s_i$ sono polinomi di grado $k$ per $i=0,\dots ,n-1$, localmente $C^{\mathrm{\infty }}$, e tali che globalmente $s_x\in C_{k-1}$.

Nel caso $k=1$, consideriamo l'intervallo $[a,b]$ e una partizione, in ciascun sottointervallo si ha un polinomio di grado 1, e le $s_i$ sono localmente polinomi di grado 1 e $s$ è di classe $C^0$ sull'intervallo; si ha l'interpolazione lineare a tratti.

Nel caso $k=3$ si ha l'interpolazione spline cubica, si considera una partizione di $[a,b]$, in ciascun sottointervallo definisco un polinomio di grado 3. $s$ è di classe $C^2$ sull'intervallo.

Sistema per la determinazione della spline

$s_i(x)$ sono polinomi di grado 3, con $i=0,\dots ,n-1$. Abbiamo $4n$ coefficienti da determinare.

Dobbiamo imporre le seguenti condizioni:

• $n+1$ condizioni di interpolazione, perché la spline deve valere $y_i$ in $x_i$;
• $n-1$ condizioni di raccordo con continuità, cioè $s_{i-1}(x_i)=s_i(x_i)$ per $i=1,\dots n-1$;
• $2n-2$ condizioni sulle derivate, cioè $s_{i^{\prime }-1}(x_i)=s_{i^{\prime }}(x_i)$ e $s_{i^{″}-1}(x_i)=s_{i^{″}}(x_i)$ perché $s$ dev'essere di classe $C^2$.

Quindi mancano due condizioni per avere un sistema determinato di $4n$ equazioni in $4n$ incognite.

In base alle due condizioni aggiuntive che si impongono si hanno diversi tipi di spline:

Spline naturale

imponiamo che $s^{″}(a)=s^{″}(b)=0$.

Spline completa e vincolata

$s_{0^{\prime }}(a)=f^{\prime }(a)$ e $s_{n^{\prime }-1}(b)=f^{\prime }(b)$.

Spline periodica

supponiamo $T=b-a$, allora $s_{0^{\prime }}(a)=s_{n^{\prime }-1}(b)$ e analogamente $s_{0^{″}}(a)=s_{n^{″}-1}(b)$.

Abbiamo un sistema di $4n$ equazioni in $4n$ incognite, che però ha grandi dimensioni ed è spesso mal condizionato.

Metodo dei momenti

Definizione della funzione: sia $M_i$ la derivata seconda di $s$ in $x_i$, per $i=1,\dots ,n$ (incognite da calcolare), dove $M_i$ sono detti momenti. Siccome $s_i$ ha grado 3, la sua derivata è un polinomio di grado 1, inoltre, imponendo le condizioni:

$\{\begin{array}{c}{s}_{i^{″}}(x_i)=M_i\\ s_{i^{″}+1}(x_{i+1})=M_{i+1}\end{array}$

si ottiene

$s_{i^{″}}(x)=M_{i+1}(x-x_i)/h_i-M_i(x-x_{i+1})/h_i$

dove $h_i=x_{i+1}-x_i$.

Integrando due volte otteniamo:

$s_{i^{\prime }}(x)=\frac{M_{i+1}(x-x_i{)}^2}{2h_i}-M_i\frac{(x-x_{i+1}{)}^2}{2H_i}+{\alpha }_i$

Integrando nuovamente

$s_i(x)=M_{i+1}\frac{(x-x_i{)}^3}{6H_i}-M_i\frac{(x-x_{i+1}{)}^3}{6H_i}+\alpha \ast (x-x_i)+{\beta }_i\text{equazione }\ast$

Imposizione delle condizioni di interpolazione: per determinare i nuovi parametri ${\alpha }_i,{\beta }_i$, imponiamo le condizioni di interpolazione e otteniamo:

$s_i(x_i)=-M_i\frac{(x_i-x_{i-1}{)}^3}{6H_i}+{\beta }_i=f(x_i)$

Inoltre

$s_i(x_{i+1})=M_{i+1}\ast (x_{i+1}-x_i{)}^3/(6H_i)+{\alpha }_i(x_{i+1}-x_i)+{\beta }_i=f(x_{i+1})$

dove l'ultimo passaggio vale per la condizione di raccordo con continuità $s_i(x_{i+1}=s_{i+1}(x_{i+1})$. Quindi, per $i=0,\dots ,n-1$ otteniamo:

$\{\begin{array}{c}{M}_i\frac{H_i^2}{6}+{\beta }_i=f(x_i)\\ M_{i+1}\frac{H_i^2}{6}+{\alpha }_i{H}_i+{\beta }_i=f(x_{i+1})\end{array}$

Risolviamo il sistema:

$\{\begin{array}{c}{M}_i\frac{H_i^2}{6}-f(x_i)=-{\beta }_i\\ M_{i+1}\frac{H_i^2}{6}+{\alpha }_i{H}_i-f(x_{i+1})=-{\beta }_i\end{array}$

ed eguagliando le due equazioni otteniamo:

$M_i\frac{H_i^2}{6}-f(x_i)=M_{i+1}\frac{H_i^2}{6}+{\alpha }_i{H}_i-f(x_{i+1})$

$f(x_{i+1})-f(x_i)+(M_i-M_{i+1})\frac{H_i^2}{6}={\alpha }_i{H}_i$

${\alpha }_i=\frac{f(x_{i+1})-f(x_i)}{H_i}+(M_i-M_{i+1})\frac{H_i}{6}$

invece ricaviamo ${\beta }_i$ dalla prima equazione:

${\beta }_i=-M_i\frac{H_i^2}{6}+f(x_i)$

e come soluzione del sistema otteniamo

$\{\begin{array}{c}{\beta }_i=f(x_i)-M_i{H}_i^2/6\\ {\alpha }_i=\frac{f(x_{i+1})-f(x_i)}{H_i}-H_i/6(M_{i+1}-M_i)\end{array}$

Imposizione delle condizioni di raccordo con continuità delle derivate:

$s_{i^{\prime }-1}(x)=\frac{M_i(x-x_{i-1}{)}^2}{2h_{i-1}}-M_{i-1}\frac{(x-x_i{)}^2}{2H_{i-1}}+{\alpha }_i$

$=\frac{M_i(x-x_{i-1}{)}^2}{2h_{i-1}}-M_{i-1}\frac{(x-x_i{)}^2}{2H_{i-1}}+\frac{f(x_i)-f(x_{i-1})}{H_{i-1}}-H_{i-1}/6(M_i-M_{i-1})$

Quindi, valutando in $x_i$ si ottiene:

$s_{i^{\prime }-1}(x_i)=\frac{M_i{h}_{i-1}}{2}+\frac{f(x_i)-f(x_{i-1})}{H_{i-1}}-H_{i-1}/6(M_i-M_{i-1})$

analogamente

$s_{i^{\prime }}(x_i)=\frac{M_{i+1}{h}_i}{2}+\frac{f(x_{i+1})-f(x_i)}{H_i}-H_i/6(M_{i+1}-M_i)$

ed eguagliando le due espressioni:

$M_{i-1}\ast H_{i-1}/6+M_i(H_{i-1}+H_i)/3+M_{i+1}{H}_i/6=\frac{f(x_{i+1})-f(x_i)}{H_i}-\frac{f(x_i)-f(x_{i-1})}{H_{i-1}}$

e inserendo l'espressione della differenza divisa:

$M_{i-1}\ast H_{i-1}/6+M_i(H_{i-1}+H_i)/3+M_{i+1}{H}_i/6=f[x_{i+1},x_i,x_{i-1}]\ast (H_i-H_{i-1})$

e dividendo ambo i membri per $H_i-H_{i-1}$ otteniamo:

$M_{i-1}\ast \frac{H_{i-1}}{6(h_i-h_{i-1})}+M_i/3+M_{i+1}\frac{H_i}{6(h_i-h_{i-1})}=f[x_{i+1},x_i,x_{i-1}]$

e moltiplicando per 6

$M_{i-1}\ast \frac{H_{i-1}}{H_i-H_{i-1}}+2M_i+M_{i+1}\frac{H_i}{H_i-H_{i-1}}=6f[x_{i+1},x_i,x_{i-1}]$

Chiamiamo

${\gamma }_i=\frac{H_{i-1}}{H_i-H_{i-1}}$

e simmetricamente

${\delta }_i=\frac{H_i}{H_i-H_{i-1}}$

e sono due quantità positive, che hanno la proprietà

${\gamma }_i+{\delta }_i=1,\mathrm{\forall }i$

per qualsiasi tipo di nodi.

Abbiamo $n-1$ equazioni in $n+1$ incognite, date da:

$M_{i-1}\ast {\gamma }_i+2M_i+M_{i+1}{\delta }_i=6f[x_{i+1},x_i,x_{i-1}]$

Imposizione delle condizioni iniziali: la situazione dipende dalle condizioni imposte.

Spline naturale: imponiamo che $M_0=s^{″}(x_0)=0$, e $M_n=s^{″}(x_n)=0$. Allora otteniamo un sistema $A𝒙=𝒃$ dove $A$ è una matrice tridiagonale della forma:

$A=\begin{array}{ccccc}2& {\delta }_1& 0& 0& 0\\ {\gamma }_2& 2& {\delta }_2& 0& 0\\ 0& {\gamma }_3& 2& \dots & 0\\ 0& 0& \dots & 2& {\delta }_{i-1}\\ 0& 0& 0& {\gamma }_{n-1}& 2\end{array}$

$𝒙=\begin{array}{c}{M}_1\\ M_2\\ \dots \\ M_{n-2}{M}_{n-1}\end{array}$

è il vettore delle incognite, mentre il termine noto è

$𝒙=\begin{array}{c}f[x_0,x_1,x_2]\\ f[x_1,x_2,x_3]\\ f[x_2,x_3,x_4]\\ \dots \\ f[x_{n-2},x_{n-1},x_n]\end{array}$

La matrice è diagonalmente dominante sulle righe.

Se i nodi sono equispaziati, cioè se $H_i=H\mathrm{\forall }i$, la matrice è simmetrica. Per risolvere il sistema si può usare l'eliminazione gaussiana senza pivot, oppure metodi diretti che convergono per la dominanza diagonale stretta.

Spline completa: imponiamo le condizioni

$s^{\prime }(x_0)=f^{\prime }(a),s^{\prime }(x_n)=f^{\prime }(b)$

Tenendo conto che

$s^{\prime }(x)=\frac{M_i(x-x_{i-1}{)}^2}{2h_{i-1}}-M_{i-1}\frac{(x-x_i{)}^2}{2H_{i-1}}+\frac{f(x_i)-f(x_{i-1})}{H_{i-1}}-H_{i-1}/6(M_i-M_{i-1})$

si ha

$s^{\prime }(x_0)=-M_0{H}_0/2+f[x_0,x_1]+H_0/6(M_0-M_1)$

$f^{\prime }(x_0)=H_{0}f[x_0,x_0]=s^{\prime }(x_0)=-M_0{H}_0/3-M_1{H}_0/6+H_{0}f[x_0,x_1]$

che si riscrive come

$H_0/3M_0-H_0/6M_1=(f[x_0,x_1]-f[x_0,x_0])\ast H_0=f[x_0,x_0,x_1]\ast H_0$

dividendo per $H_0$ e moltiplicando per 6:

$2M_0+M_1=6f[x_0,x_0,x_1]$

Analogamente imponendo $f^{\prime }(x_n)=s^{\prime }(x_n)$ si ottene:

$M_{n-1}+2M_n=6f[x_{n-1},x_n,x_n]$

In questo caso abbiamo un sistema di $n+1$ equazioni della forma $A𝒙=𝒃$ con:

$A=\begin{array}{ccccc}2& 1& 0& 0& 0\\ {\gamma }_1& 2& {\delta }_1& 0& 0\\ 0& {\gamma }_2& 2& {\delta }_2& 0\\ 0& 0& {\gamma }_{n-1}& 2& {\delta }_{n-1}\\ 0& 0& 0& 1& 2\end{array}$

Il vettore delle incognite è

$\begin{array}{c}{M}_0\\ M_1\\ \dots \\ M_{n-1}\\ M_n\end{array}$

Il termine noto è:

$6\ast \begin{array}{c}f[x_0,x_0,x_1]\\ f[x_0,x_1,x_2]\\ \dots \\ f[x_{n-2},x_{n-1},x_n]\\ f[x_{n-1}{x}_n,x_n]\end{array}$

Abbiamo una matrice diagonalmente dominante in senso stretto, non è simmetrica neanche nel caso dei nodi equispaziati.

Spline periodica: $s^{\prime }(x_0)=s^{\prime }(x_n)$ e analogamente $s^{″}(x_0)=s^{″}(x_n)$.

La seconda si traduce nella condizione $M_0=M_n$. Prendendo l'espressione di $s_{i^{\prime }}$, si impone l'uguaglianza e si ha $f(x_n)=f(x_0)$. Introduciamo un nodo fittizio $x_{n+1}$ esterno all'intervallo, imponiamo che $H_n=H_0$. Imponendo tutte le condizioni si ottiene

${\delta }_n{M}_1+{\gamma }_n{M}_{n-1}+2M_n=6f[x_{n-1},x_n,x_{n+1}]$

${\gamma }_1{M}_n+2M_1+{\delta }_1{M}_2=6f[x_0,x_1,x_2]$

Quindi il sistema $A𝒙=𝒃$ è tale che:

$A=\begin{array}{cccccc}2& {\delta }_1& 0& 0& 0& {\gamma }_1\\ {\gamma }_2& 2& {\delta }_2& 0& 0& 0\\ 0& {\gamma }_3& 2& {\delta }_3& 0& 0\\ 0& 0& 0& {\gamma }_{n-1}& 2& {\delta }_{n-1}\\ {\gamma }_n& 0& 0& 0& {\delta }_n& 2\end{array}$

Il vettore delle incognite è

$𝒙=\begin{array}{c}{M}_1\\ M_2\\ \dots \\ M_{n-1}\\ M_n\end{array}$

mentre il vettore del termine noto è

$\begin{array}{c}f[x_0,x_1,x_2]\\ f[x_1,x_2,x_3]\\ \dots \\ f[x_{n-2},x_{n-1},x_n]\\ f[x_{n-1},x_n,x_{n+1}]\end{array}$

Questa matrice $n\times n$ non è tridiagonale a causa di ${\delta }_n$ in posizione $A_{n1}$ e ${\gamma }_1$ in posizione $A_{1n}$. La matrice in ogni caso è diagonalmente dominante in senso stretto. Se i nodi sono equispaziati si ha la simmetria.

Conclusione

Nei tre casi, la matrice è diagonalmente dominante in senso stretto, e per il primo teorema di Gersch-Gorin è non singolare.

Teorema 5.7

Se i nodi $A=x_0<x_1<\dots <x_{n-1}<x_n=b$ sono distinti, allora esiste ed è unica la spline cubica che interpola $f$ funzione assegnata.

Teorema 5.8

Supponiamo che $f$ sia $C_4([a,b])$, consideriamo un partizionamento di $[a,b]$ in sottointervalli $H_i$, poniamo $H=\underset{i}{max}{H}_i$ e chiamiamo $\beta =\frac{H}{\underset{i}{min}{H}_i}$, e sia $s$ la spline cubica interpolante.

$|f_r-s_r{|}^{\mathrm{\infty }}\le c_r^{4-r}|f^{(4)}|$