Nugroho's blog.: Newton Polinomial.

Tuesday, April 11, 2017

Newton Polinomial.


Here's code for Newton's divided differences interpolation polynomial (quite mouthful huh, :) ).

The purpose of this method is to create a function (polynomial) that passes through given set of data points.




I read data point from several edit box.

procedure TForm1.Button3Click(Sender: TObject);
var i:integer;
begin
for i:=0 to n do begin
x[i]:=strToFloat(kx[i].Text);
y[i]:=strToFloat(ky[i].Text);
end;
xc:=strToFloat(kxc.Text);
yc:=fn(xc);
kyc.Text:=floatToStr(yc);
gambarNewton;
end;

.
kx and ky is tEdit created when button1 is clicked

procedure TForm1.Button1Click(Sender: TObject);
var i:integer;
begin
button2.Enabled:=true;
button3.Enabled:=true;
button4.Enabled:=true;
button5.Enabled:=true;

n:=strToInt(edit1.Text);
kxc:=tEdit.Create(form1); kyc:=tEdit.Create(form1);
kxc.Parent:=form1; kyc.Parent:=form1;
kxc.Left:=36; kyc.Left:=72;
kxc.Width:=36; kyc.Width:=36;
kxc.Text:='0,5';

for i:=0 to n do begin
kx[i]:=tEdit.Create(form1); ky[i]:=tEdit.Create(form1);
kx[i].Parent:=form1; ky[i].Parent:=form1;
kx[i].Top:=36+36*i; ky[i].Top:=36+36*i;
kx[i].Left:=36; ky[i].Left:=72;
kx[i].Width:=36; ky[i].Width:=36;
kx[i].Text:=intToStr(i); ky[i].Text:=intToStr(i);
end;
end;


xc is x coordinate where the corresponding y (yc) is obtained using Newton method by calling it

yc=fn(xc)

function tform1.fn(xs:real):real;
var i:integer;fs:real;
begin
fs:=0;
for i:=0 to n do begin
fs:=fs+b(i,0)*c(xs,i);
end;
fn:=fs;
end;



the fn function call the two other function. The b function, a recursive contain divided difference like this

function tform1.b(i,j:integer):real;
begin
if i=0 then b:=y[0]
else if (i-j)=1 then
b:=(y[i]-y[j])/(x[i]-x[j])
else
b:=(b(i,j+1)-b(i-1,j))/(x[i]-x[j]);
end;

and c function, a recursive function (or you could rewrite it using simple for command)

function tform1.c(xs:real;i:integer):real;
begin
if i=0 then c:=1
else c:=(xs-x[i-1])*c(xs,i-1);
end;

and finally, draw the data and the function on image1

fprocedure tform1.gambarNewton;
var i,x0,y0:integer;px,py:real;
begin
x0:=image1.Width div 2; y0:=image1.Height div 2;
image1.Canvas.Brush.Color:=clLime;
image1.Canvas.Rectangle(0,0,image1.Width,image1.Height);
image1.Canvas.Brush.Color:=clWhite;
image1.Canvas.Pen.Color:=clBlack;
image1.Canvas.MoveTo(0,y0); image1.Canvas.LineTo(image1.Width,y0);
image1.Canvas.MoveTo(x0,0); image1.Canvas.LineTo(x0,image1.Height);
for i:=-300 to 300 do begin
px:=i/skala; py:=skala*fn(px);
image1.Canvas.Pixels[x0+i,y0-round(py)]:=clGreen;
end;
for i:=0 to n do begin
px:=x0+skala*x[i]; py:=y0-skala*y[i];
image1.Canvas.Ellipse(round(px)-7,round(py)-7,round(px)+7,round(py)+7);
end;
px:=x0+skala*xc; py:=y0-skala*yc;
image1.Canvas.Brush.Color:=clred;
image1.Canvas.Ellipse(round(px)-7,round(py)-7,round(px)+7,round(py)+7);
image1.Canvas.Brush.Color:=clwhite;
end;






No comments:

323f (5) amp (1) android (12) apple (7) arduino (18) art (1) assembler (21) astina (4) ATTiny (23) blackberry (4) camera (3) canon (2) cerita (2) computer (106) crazyness (11) debian (1) delphi (39) diary (286) flash (8) fortran (6) freebsd (6) google apps script (8) guitar (2) HTML5 (10) IFTTT (7) Instagram (7) internet (12) iOS (5) iPad (6) iPhone (5) java (1) javascript (1) keynote (2) LaTeX (6) lazarus (1) linux (29) lion (15) mac (28) macbook air (8) macbook pro (3) macOS (1) Math (3) mathematica (1) maverick (6) mazda (4) microcontroler (35) mountain lion (2) music (37) netbook (1) nugnux (6) os x (36) php (1) Physicist (29) Picture (3) programming (189) Python (109) S2 (13) software (7) Soliloquy (125) Ubuntu (5) unix (4) Video (8) wayang (3) yosemite (3)