Monotónní kubická interpolace - Monotone cubic interpolation

V matematický pole numerická analýza, monotónní kubická interpolace je varianta kubická interpolace který zachovává monotónnost z soubor dat interpolován.

Monotónnost je zachována lineární interpolace ale není zaručeno kubická interpolace.

Monotónní kubická Hermitova interpolace

Příklad zobrazující monotónní kubickou interpolaci (červeně) a monotónní kubickou interpolaci (modře) souboru monotónních dat.

Monotónní interpolaci lze provést pomocí kubický Hermitův spline s tečnami ${ displaystyle m_ {i}}$ upraveno, aby byla zajištěna monotónnost výsledného Hermitova spline.

Algoritmus je k dispozici také pro monotónní kvintický Hermitova interpolace.

Interpolantní výběr

Existuje několik způsobů výběru interpolačních tečen pro každý datový bod. V této části je popsáno použití metody Fritsch-Carlson. Všimněte si, že je vyžadován pouze jeden průchod algoritmu.

Nechť jsou datové body ${ displaystyle (x_ {k}, y_ {k})}$ indexováno v seřazeném pořadí pro ${ displaystyle k = 1, , tečky , n}$.

1. Vypočítejte sklon svahu sekanční linie mezi po sobě následujícími body:

   ${ displaystyle delta _ {k} = { frac {y_ {k + 1} -y_ {k}} {x_ {k + 1} -x_ {k}}}}$

   pro ${ displaystyle k = 1, , tečky , n-1}$.
   
   
2. Tato přiřazení jsou prozatímní a ve zbývajících krocích mohou být nahrazena. Inicializujte tangenty v každém datovém bodě interiéru jako průměr secants,

   ${ displaystyle m_ {k} = { frac { delta _ {k-1} + delta _ {k}} {2}}}$

   pro ${ displaystyle k = 2, , tečky , n-1}$.
   
   Pro koncové body použijte jednostranné rozdíly:

   ${ displaystyle m_ {1} = delta _ {1} quad { text {a}} quad m_ {n} = delta _ {n-1} ,}$.

   Li ${ displaystyle delta _ {k-1}}$ a ${ displaystyle delta _ {k}}$ mít opačné znaky, nastavit ${ displaystyle m_ {k} = 0}$.
   
   
3. Pro ${ displaystyle k = 1, , tečky , n-1}$kde vůbec ${ displaystyle delta _ {k} = 0}$ (kdykoli dva po sobě jdoucí ${ displaystyle y_ {k} = y_ {k + 1}}$ jsou rovny),
   soubor ${ displaystyle m_ {k} = m_ {k + 1} = 0,}$ protože spline spojující tyto body musí být plochá, aby byla zachována monotónnost.
   Pro tyto kroky 4 a 5 ignorujte ${ displaystyle k ,}$.
   
   
4. Nechat

   ${ displaystyle alpha _ {k} = m_ {k} / delta _ {k} quad { text {a}} quad beta _ {k} = m_ {k + 1} / delta _ { k}}$.

   Pokud ano ${ displaystyle alpha _ {k}}$ nebo ${ displaystyle beta _ {k}}$ je záporné, pak vstupní datové body nejsou striktně monotónní a ${ displaystyle (x_ {k}, , y_ {k})}$ je místní extremum. V takových případech, po částech monotónní křivky lze stále generovat výběrem ${ displaystyle m_ {k} = 0 ,}$ -li ${ displaystyle alpha _ {k} <0}$ nebo ${ displaystyle m_ {k + 1} = 0 ,}$ -li ${ displaystyle beta _ {k} <0}$, ačkoli přísná monotónnost není globálně možná.
   
   
5. Aby se zabránilo přestřelení a zajistit monotónnost, musí být splněna alespoň jedna z následujících tří podmínek:

a) funkce

${ displaystyle phi _ {k} = alpha _ {k} - { frac {(2 alpha _ {k} + beta _ {k} -3) ^ {2}} {3 ( alpha _ {k} + beta _ {k} -2)}}> 0 ,}$, nebo

b) ${ displaystyle alpha _ {k} +2 beta _ {k} -3 leq 0 ,}$, nebo

(C) ${ displaystyle 2 alpha _ {k} + beta _ {k} -3 leq 0 ,}$.

K zajištění přísné monotónnosti stačí pouze podmínka (a): ${ displaystyle phi _ {k}}$ musí být pozitivní.

Jedním z jednoduchých způsobů, jak uspokojit toto omezení, je omezit vektor ${ displaystyle ( alpha _ {k}, , beta _ {k})}$ do kruhu o poloměru 3. To znamená, že pokud ${ displaystyle alpha _ {k} ^ {2} + beta _ {k} ^ {2}> 9 ,}$, pak nastavte

${ displaystyle tau _ {k} = { frac {3} { sqrt { alpha _ {k} ^ {2} + beta _ {k} ^ {2}}}} ,}$,

a změňte měřítko tečen pomocí

${ displaystyle m_ {k} = tau _ {k} , alpha _ {k} , delta _ {k} quad { text {a}} quad m_ {k + 1} = tau _ {k} , beta _ {k} , delta _ {k} ,}$.

Jinak stačí omezit ${ displaystyle alpha _ {k} leq 3}$ a ${ displaystyle beta _ {k} leq 3 ,}$. Toho dosáhnout, pokud ${ displaystyle alpha _ {k}> 3 { text {nebo}} beta _ {k}> 3 ,}$, pak nastavte ${ displaystyle m_ {k} = 3 , delta _ {k} ,}$.

Kubická interpolace

Po předzpracování výše je vyhodnocení interpolovaného splajnu ekvivalentní kubický Hermitův spline pomocí údajů ${ displaystyle x_ {k}}$, ${ displaystyle y_ {k}}$, a ${ displaystyle m_ {k}}$ pro ${ displaystyle k = 1, , tečky , n}$.

Vyhodnotit na ${ displaystyle x}$, najděte index ${ displaystyle k}$ v pořadí kde ${ displaystyle x}$, leží mezi ${ displaystyle x_ {k}}$, a ${ displaystyle x_ {k + 1}}$, to znamená: ${ displaystyle x_ {k} leq x leq x_ {k + 1}}$. Vypočítat

${ displaystyle Delta = x_ {k + 1} -x_ {k} quad { text {a}} quad t = { frac {x-x_ {k}} { Delta}}}$

pak je interpolovaná hodnota

${ displaystyle f _ { text {interpolovaný}} (x) = y_ {k} cdot h_ {00} (t) + Delta cdot m_ {k} cdot h_ {10} (t) + y_ {k +1} cdot h_ {01} (t) + Delta cdot m_ {k + 1} cdot h_ {11} (t)}$

kde ${ displaystyle h_ {ii}}$ jsou základní funkce pro kubický Hermitův spline.

Příklad implementace

Následující JavaScript implementace vezme datovou sadu a vytvoří monotónní kubickou spline interpolantovou funkci:

/ * Monotónní kubická spline interpolace   Příklad použití:var f = createInterpolant ([0, 1, 2, 3, 4], [0, 1, 4, 9, 16]);var zpráva = '';pro (var x = 0; x <= 4; x + = 0,5) {var xSquared = f (x);zpráva + = x + 'na druhou je o' + xSquared + ' n';	}výstraha (zpráva);*/var createInterpolant = funkce(xs, ys) {	var i, délka = xs.délka;		// Řešení problémů s délkou	-li (délka != ys.délka) { házet "Potřebujete stejný počet xs a ys."; }	-li (délka === 0) { vrátit se funkce(X) { vrátit se 0; }; }	-li (délka === 1) {		// Impl: Předpočítání výsledku zabrání problémům, pokud je ys později mutován, a umožní smetí ys		// Impl: Unary plus správně převádí hodnoty na čísla		var výsledek = +ys[0];		vrátit se funkce(X) { vrátit se výsledek; };	}		// Uspořádejte xs a ys tak, aby bylo xs tříděno	var indexy = [];	pro (i = 0; i < délka; i++) { indexy.tlačit(i); }	indexy.třídit(funkce(A, b) { vrátit se xs[A] < xs[b] ? -1 : 1; });	var oldXs = xs, oldYs = ys;	// Impl: Vytváření nových polí také zabrání problémům, pokud budou vstupní pole později mutována	xs = []; ys = [];	// Impl: Unary plus správně převádí hodnoty na čísla	pro (i = 0; i < délka; i++) { xs.tlačit(+oldXs[indexy[i]]); ys.tlačit(+oldYs[indexy[i]]); }		// Získejte po sobě jdoucí rozdíly a svahy	var dys = [], dxs = [], slečna = [];	pro (i = 0; i < délka - 1; i++) {		var dx = xs[i + 1] - xs[i], dy = ys[i + 1] - ys[i];		dxs.tlačit(dx); dys.tlačit(dy); slečna.tlačit(dy/dx);	}		// Získejte koeficienty stupně 1	var c1s = [slečna[0]];	pro (i = 0; i < dxs.délka - 1; i++) {		var m = slečna[i], mDalší = slečna[i + 1];		-li (m*mDalší <= 0) {			c1s.tlačit(0);		} jiný {			var dx_ = dxs[i], dxDalší = dxs[i + 1], běžný = dx_ + dxDalší;			c1s.tlačit(3*běžný/((běžný + dxDalší)/m + (běžný + dx_)/mDalší));		}	}	c1s.tlačit(slečna[slečna.délka - 1]);		// Získejte koeficienty stupně 2 a stupně 3	var c2s = [], c3s = [];	pro (i = 0; i < c1s.délka - 1; i++) {		var c1 = c1s[i], m_ = slečna[i], invDx = 1/dxs[i], běžný_ = c1 + c1s[i + 1] - m_ - m_;		c2s.tlačit((m_ - c1 - běžný_)*invDx); c3s.tlačit(běžný_*invDx*invDx);	}		// Návrat interpolantní funkce	vrátit se funkce(X) {		// Bod zcela vpravo v datové sadě by měl poskytnout přesný výsledek		var i = xs.délka - 1;		-li (X == xs[i]) { vrátit se ys[i]; }				// Hledání intervalu x je in, vrací odpovídající y, pokud x je jedním z původních xs		var nízký = 0, střední, vysoký = c3s.délka - 1;		zatímco (nízký <= vysoký) {			střední = Matematika.podlaha(0.5*(nízký + vysoký));			var x tady = xs[střední];			-li (x tady < X) { nízký = střední + 1; }			jiný -li (x tady > X) { vysoký = střední - 1; }			jiný { vrátit se ys[střední]; }		}		i = Matematika.max(0, vysoký);				// Interpolovat		var rozdíl = X - xs[i], diffSq = rozdíl*rozdíl;		vrátit se ys[i] + c1s[i]*rozdíl + c2s[i]*diffSq + c3s[i]*rozdíl*diffSq;	};};

Reference

• Fritsch, F. N .; Carlson, R. E. (1980). "Monotónní po částech kubická interpolace". Časopis SIAM o numerické analýze. SIAM. 17 (2): 238–246. doi:10.1137/0717021.
• Dougherty, R.L .; Edelman, A .; Hyman, J. M. (duben 1989). „Pozitivní, monotonická nebo konvexní kubická a kvintická Hermitova interpolace. Matematika výpočtu. 52 (186): 471–494. doi:10.2307/2008477.

externí odkazy

• GPLv 2 licencované C ++ implementace: MonotCubicInterpolator.cpp MonotCubicInterpolator.hpp