From:	martin DOT kahlert AT keksy DOT mchp DOT siemens DOT de (Martin Kahlert)
Newsgroups:	comp.os.msdos.djgpp
Subject:	Re: A useful class I made for producing formulas
Date:	30 Jul 1997 12:21:12 GMT
Organization:	Siemens Inc.
Lines:	47
Message-ID:	<5rnbjo$jbk$1@salomon.mchp.siemens.de>
References:	<01bc9ae9$7bc81960$32f282c1 AT damien>
Reply-To:	Martin DOT Kahlert AT mchp DOT siemens DOT de
NNTP-Posting-Host:	keksy.mchp.siemens.de
Mime-Version:	1.0
To:	djgpp AT delorie DOT com
DJ-Gateway:	from newsgroup comp.os.msdos.djgpp

In article <01bc9ae9$7bc81960$32f282c1 AT damien>,
	"Richard Birch" <richardbirch AT dial DOT pipex DOT com> writes:
> Anyone interested in a class that gives you a formula when you give it a
> set of results. e.g. give it x=1 when y=1, x=2 when y=4, x=3 when y=9, the
> class will give you the formula y=x^2. That was a VERY simple example - it
> will work given *any* real numbers, and produces formulas up to the order
> of 50. It uses Shawn Hargreaves' ALLEGRO for graphics manipulation so it
> needs that to compile. It produces the formula using simultaneous
> equations.
> Other functions are:
> 
> - drawing the graph of the formula on the screen or to a bitmap
> - setting the scale, axis etc to your needs
> - giving a y result given an x by running x through the formula
> 
> If *anyone* is interested, e-mail me and I'll reply attaching the source
> code (one 16KB header). I'm quite proud of this as this is the only code of
> this nature I know of, and would be pleased to receive ideas for
> improvement.

What you do is called polynomial interpolation.
Your way to get the formula is (I guess) instable.
Try to find a formula for
(0 0),(1 1),...,(50,50) (exakt result: y(x)=x)
With simultanious equations, you have a linear set
of equations of the form
a_0                                     =0
a_0+ 1*a_1+...................+   1*a_50=1
a_0+ 2*a_1+4*a_2+.............+2^50*a_50=2
.......................................
.......................................
a_0+50*a_1+50^2+a_2+..........+50^50*a_50=50

50^50 is a number with 85 digits, while your
computer can manage only its first 16. That causes trouble.
(Verify in C: pow(50.0e0,50.0e0)+1.0e0 == pow(50.0e0,50.0e0))

There is a much faster, much more reliable way of obtaining
your formula: Look for Newton interpolation formula in
any beginners book about numerical mathematics. 

The algorithm is about 10 lines of code and handles as much
data you want, and it only needs RAM in the order of N
(N x/y pairs) instead of your N*N.

Hope that helps,
Martin.

- Raw text -