Transcript 1. Interpolācija
Interpolēšanas metožu izmantošana
Delphi un Matlab
vidē.
Praksē nereti rodas nepieciešamība noskaidrot sakarības dažādos procesos un parādībās un izteikt šīs sakarības matemātiskas izteiksmes veidā.
Ir zināmi trīs funkciju uzdošanas veidi: analītiskais, grafiskais un tabulārais.
Interpolācija visbiežāk jāizmanto tad, ja funkcija ir uzdota tabulas veidā.
Funkciju interpolācija
y x 1 2 … i i+1 … n – 1 n
y
a
*
x
b
Lineārā interpolācija - Taisnes vienādojums
y i y i
1
a
*
x i a
*
x i
1
b
b
Funkcijas vērtības kaimiņ mezglu punktos
Atrodam
a
un
b
vērtības:
b b
y i
1
a
y
y i i
a
*
a
*
x i x i
1 1
y i x i
1
x i y i
x i
*
y i
y i
1
x i
1
y x i i a
*
x i
Nosacījums
x i < X t < x i+1
Piemērs.
xt := strtofloat(edit1.text); For i := 1 to n-1 do begin if (xt < x[i+1]) and (xt > x[i]) then begin a:= (y[i+1] – y[i])/(x[i+1] – x[i]); b:= y[i] – x[i] * a; yt:= a*xt + b; end end;
Kvadrātiskā interpolācija a*x 2 + b*x + c = 0; Parabolas vienādojums y[i-1] = x[i-1]*x[i-1]*a + x[i-1]*b + c y[i] = x[i]*x[i]*a + x[i]*b + c y[i+1] = x[i+1]*x[i+1]*a + x[i+1]*b + c
No vienādojumu sistēmas atrodam
a, b
un
c
vērtības: a = ∆ a /∆; b = ∆ b /∆; c = ∆ c /∆; kur ∆, ∆ a , ∆ b , ∆ c – trešās kārtas determinanti.
x i-1 <= x <= x i+1
Splain-interpolācija
m
kārtības splain dēvējas funkcija, kura ir m pakāpes polinoms uz katras atstarpes starp mezglu punktiem un visos iekšējos mezglos apmierina funkcijas un atvasinājumu kontinuitātes nosacījumiem.
Vislielāko izplatīšanu saņēma kuba splains
y=ax
3
+bx
2
+cx+d.
Atradīsim atkarības viņa a, b, c, d koeficentu aprēķināšanai pa funkcijas un tās pirmo atvasinājumu vērtībām kaimiņ mezglu punktos.
Izkravāšanu vienkāršojumam pieņemsim, ka atstarpes kreisa leņķiska punkta
x i ,x i+1
argumenta vērtība, kurai meklējam splain funkciju, vienāds nullei. To vienmēr var izdarīt pārveidi
x=x i +z
, kur
z
jauns arguments.
y=az 3 +bz 2 +cz+d
z = x - x i x i z 1 = 0 x i+1 z 2 = x i+1 - x i
Tad uzdevums ir noformulējams sekojoši: atrast nezināmus a, b, c, d splaina koeficentus, apmierinoši mezglu punktos nosacījumam
z 1 =0; y(z 1 )=f1; y'(z 1 )=f 1 '; y(z 2 )=f 2 ; y' (z 2 )=f 2 '.
Šeit
f 1 ,f 2 ,f 1 ',f 2 '
- attiecīgi interpolēšanas funkcijas un tās atvasinājuma vērtības mezglu punktos
z 1 , z 2
.
Lai atrastu
a,b,c,d
saliksim vienādojumu sistēmu, izmantojot nosacījumus mezglu punktos f 1 = d; f 1 ′ = c; f 2 = az 2 3 + bz 2 2 + cz 2 +d; f 2 ′ = 3az 2 2 + 2bz +c; tātad d = f 1 ; c = f 1 ′
Paliekot atrastās
d, c
vērtības divos pēdējos vienādojumos un atļaujot tos relatīvi
a, b
saņemsim a = (z 2 (f 2 ′ f 1 ′ ) – 2(f 2 - f 1 ′ * z 2 – f 1 ))/z 2 3 b = (3(f 2 - f 1 ′ * z 2 – f 1 ) – z 2 (f 2 ′ f 1 ′ ))/z 2 2
Paliekot atrastos koeficentus kuba splainu un uzdodot argumenta vērtību, aprēķināsim funkcijas vērtību starp mezglu punktiem.
y = az 3 + bz 2 + cz + d
Bieži pie interpolācijas ir esamas tikai funkcijas vērtības mezglu punktos un nav atvasinājumu vērtību. Tad pirmie atvasinājumi iekšējos mezglu punktos aprēķina pa formulām f i ′ = (f i+1 – f i-1 )/(x i+1 – x i-1 ), f i ′ = (f i+1 – f i-1 )/(2h), i = 2, .. N-1 h- attālums starp mezglu punktiem.
N – mezglu punktu skaits
Atvasinājuma vērtību pirmajā un pēdējā mezglu punktos aprēķina pa formulām: f 1 ′ = f 2 ′ f 2 ″ h; f n ′ = f n-1 ′ – f n-1 ″ h; f 2 ″ = (f 1 + f 3 – 2f 2 )/h 2; f n-1 ″ = (f n-2 + f n – 2f n-1 )/h 2; f 2 ″ , f n-1 ″ - otro atvasinājumu pietuvinātās vērtības otrajā un n - 1 punktos.
Programmas fragmenti Delphi vidē
Funkcijas kubisks splains
function
k_spl(xx:Real):Real; var d:Real; // koeficient d { pirmā atvasinājuma aprēķins mezgla punktos (i>= 2) and (i <= n-1)}
function
fi_1(i:Byte):Real; begin fi_1 := (y[i+1] - y[i-1])/(x[i+1]-x[i-1]); end; { otrā atvasinājuma aprēķins punktā i = 2}
function
f22:single; begin f22 := (y[1] + y[3] - 2*y[2])/sqr(x[2]-x[1]); end; { pirmā atvasinājuma aprēķins punktā i = 1}
function
f12:single; begin f12 := fi_1(2) - f22*(x[2]-x[1]); end;
{ otrā atvasinājuma aprēķins punktā i = n-1}
function
fn2:single; begin fn2 := (y[n-2] + y[n] - 2*y[n-1])/sqr(x[n]-x[n-1]); end; { pirmā atvasinājuma aprēķins punktā i = n}
function
f1n:single; begin f1n := fi_1(n-1) - fn2*(x[n]-x[n-1]); end; {koeficienta c aprēķināšana}
function
c(i:byte):single; begin if (i>= 2) and (i < n-1) then c := fi_1(i); if i = 1 then c := f12; if i = n-1 then c := f1n; end;
{koeficienta a aprēķinšāna} function a(i:Byte):Real; begin if (i>= 2) and (i < n-1) then begin a := ((x[i+1]-x[i])*(fi_1(i+1)-fi_1(i)) - 2*(y[i+1] - fi_1(i)*(x[i+1]-x[i]) - y[i]))/power(x[i+1]-x[i],3); end; if i= 1 then begin a := ((x[i+1]-x[i])*(fi_1(i+1)-fi_1(i)) - 2*(y[i+1] - f12*(x[i+1]-x[i]) - y[i]))/power(x[i+1]-x[i],3); end; if i= n-1 then begin a := ((x[i+1]-x[i])*(f1n-fi_1(i)) - 2*(y[i+1] - fi_1(i)*(x[i+1]-x[i]) - y[i]))/power(x[i+1]-x[i],3); end; end;
{koeficienta b aprēķināšana} function b(i:Byte):Real; begin if (i>= 2) and (i < n-1) then begin b := (3*(y[i+1] - fi_1(i)*(x[i+1]-x[i]) - y[i])-(x[i+1]-x[i])* (fi_1(i+1)-fi_1(i)))/sqr(x[i+1]-x[i]); end; if i = 1 then begin b := (3*(y[i+1] - f12*(x[i+1]-x[i]) - y[i])- (x[i+1]-x[i])*(fi_1(i+1) fi_1(i)))/sqr(x[i+1]-x[i]); end; if i= n-1 then begin b := (3*(y[i+1] - fi_1(i)*(x[i+1]-x[i]) - y[i]) (x[i+1]-x[i])*(f1n-fi_1(i)))/sqr(x[i+1]-x[i]); end; end;
{Funkcijas ķermenis} begin for j := 1 to n-1 do if(xx > x[j]) and (xx < x[j+1]) then it := j; {it – kreisa mezgla numurs x[it] - kreisa mezgla koordināta xx – argumenta vērtība} d := y[it]; k_spl := a(it)*power(xx-x[it],3) + b(it)*sqr(xx-x[it])+ c(it)*(xx-x[it]) + d; end;
{Funkcijas y(xt ) vērtības aprēķināšana}
procedure
TForm1.Button1Click(Sender: TObject); begin xt := strtofloat(edit1.Text); { lasīt vektorus x un y no stringgrid komponenta} for i := 1 to n do begin x[i] := strtofloat(stringgrid1.Cells[i-1,0]); y[i] := strtofloat(stringgrid1.Cells[i-1,1]); end; yt := k_spl(xt) ; {vērtības y(xt) izvade} edit2.Text := floattostrf(yt,fffixed,10,3); end;
{Grafika konstruēšana}
procedure
TForm1.Button2Click(Sender: TObject); begin Series1.Clear; Series2.Clear; for i := 1 to n do begin x[i] := strtofloat(stringgrid1.Cells[i-1,0]); y[i] := strtofloat(stringgrid1.Cells[i-1,1]); Series1.AddXY(x[i],y[i]); end; xg := x[1]; yg := y[1]; Series2.AddXY(xg,yg); repeat xg := xg + 0.05; yg := k_spl(xg) ; end; Series2.AddXY(xg,yf); until xg >= x[n];
{Lasīt vektorus x un y no faila, izmantojot OpenDialog komponentu} procedure TForm1.Open1Click(Sender: TObject); begin if OpenDialog1.Execute then begin AssignFile(ff,OpenDialog1.FileName); Reset(ff); i:=0; while not Eof(ff) do begin i := i+1; readln(ff,x[i],y[i]); end; end; CloseFile(ff); end;
{Fails dati.txt Notepadā} -5.0
-4.2
-3.4
-2.6
-1.8
-1.0
-0.2
0.6
1.4
2.2
9.13
1.26
-0.7
1.03
-2.05
-7.23
-6.15
-2.02
-3.02
-5.43