Quadratic Equation from 3 points

theoarc · May 18, 2010

Okay so this seems like a very simple problem to me but I can't figure it out so I definitely need some help. I'm not a mathematician but a programmer and in this case I need to figure out how to get the curve that passes through 3 points. I know this can be done easily with substitution but that tends to requires human logic in order to decide what equations to use in which order. I need another method. I need to be able to do this programmatically.

I know when I put 3 points into the TI-83 and ask it to do a Quadratic Regression it will give me a function and plot a graph in y = ax² + bx + c

That's what I am looking for a method of solving for y = ax² + bx + c from 3 points hopefully without substitution.

I hope I am making sense.
Thanks for any help.

Tedjn · May 18, 2010

You will need to substitute the points somehow in order to find the coefficients. If you have a library available that can handle common matrix operations, this can be made pretty simple, however.

phizo · May 18, 2010

theoarc said:

Okay so this seems like a very simple problem to me but I can't figure it out so I definitely need some help. I'm not a mathematician but a programmer and in this case I need to figure out how to get the curve that passes through 3 points. I know this can be done easily with substitution but that tends to requires human logic in order to decide what equations to use in which order. I need another method. I need to be able to do this programmatically.

I know when I put 3 points into the TI-83 and ask it to do a Quadratic Regression it will give me a function and plot a graph in y = ax² + bx + c

That's what I am looking for a method of solving for y = ax² + bx + c from 3 points hopefully without substitution.

I hope I am making sense.
Thanks for any help.

Hmmm care to give us the three points? or is it hypothetical exercise?

Anyway I do a bit of programing and have a spattering of maths too.

So I assume some points such as (1,6) (7,8)(10,2) or whatever, and you need to find A B &C? And you want to write a progam to do it?

Xitami · May 18, 2010

[itex]ax_1^2+bx_1+c=y_1 [/itex]
[itex]ax_2^2+bx_2+c=y_2[/itex]
[itex]ax_3^2+bx_3+c=y_2[/itex]

[itex]
\left[\begin{array}{ccc}a\\b\\c\end{array}\right] =
\left[\begin{array}{ccc}y_1\\y_2\\y_3\end{array}\right] \left[\begin{array}{ccc}x_1^2&x_1&1\\x_2^2&x_2&1\\x_3^2&x_3&1\end{array}\right] ^{-1}
[/itex]

http://en.wikipedia.org/wiki/Invertible_matrix#Inversion_of_3.C3.973_matrices

phizo · May 18, 2010

I might have a go using the points I have given and see how far I get, or not, as the case maybe. :rofl:

Unit · May 18, 2010

Do you realize that three points may not necessarily all lie on a parabola? For example, there is no parabola that contains (0, 0), (1, 1) and (2, 2).

However, there is an infinite number of parabolae that contain any two distinct points. For example, the points (-1, 0) and (1, 0). For this particular case, the axis of symmetry is at 0, so the parabolae are of the form [itex]y = ax^2 - a[/itex] where a is a real number.

phizo · May 18, 2010

OK I just started out with an idea and coubbled togeather a bit of a program.
It's not completed, just a bit of a framework to get and idea of what I am doing.

I am trying stick in some values for a b and c and work out the values which give the lowest
error, then when I have that I will have to try and narrow the range and step based on the results.

Anyway here is how far I got, no comment, obviously ;)

Code:

#include <stdio.h>

main(argc,argv)
	int  argc;
	char *argv[]; 
{

	printf("started\n");

	// y = ax2 + bx + c
	// (1,6) (7,8)(10,2) 

	float y1,y2,y3,x1,y2,x3, sa,sb,sc, na,nb,nc, step, ea, olodea;

	x1=1;y1=6;x2=7; y2=8;x3=10; y3=2; oldea=100000;

	sa=sb=sc=250;
	step=1;

	for(a=-sa; a<sa; a+=step){
		for(b=-sb; a<sb; b+=step){
			for(c=-sc; a<sc; c+=step){				e1 = a*36 + b*6 +c -1;
				e2 = a*64 + b*8 +c -7;
				e3 = a*4  + b*2 + c-10;
				et = e1 + e2 + e3;
				ea=abs(et)
				if (ea<oldea) { 
					oldea=ea;
					na=a; nb=b; nc=c; //n=new, so i am saving the nearest guess

				}			}
		}
	}

	printf("\n a=%f b-%f c=%c"; na,nb,nc);

	
}

Obviously it won't even compile at the moment!

Might give you an idea though.

phizo · May 18, 2010

So I think after I had got some results I would run the whole thing again using say 30% above and below my range and with and appropiate smaller step.

Then you just repeat the process keeping your fingers crossed! :rofl:

phizo · May 18, 2010

OK so I have version running, I can hear the fan has kicked in and I can smell burning!
250^3 is a big number I think I might reduce that!

Xitami · May 18, 2010

[itex]a_{13}=a_{23}=a_{33}=1[/itex]

phizo · May 18, 2010

OK got a bit further:-

Code:

#include <stdio.h>
#include <math.h>
main(argc,argv)
	int  argc;
	char *argv[]; 
	// y = ax2 + bx + c	// (1,6) (7,8)(10,2) 
	float y1,y2,y3,x1,x2,x3, sa,sb,sc, na,nb,nc, step, ea, oldea, a, b, c, e1,e2,e3,et;
	x1=1;y1=6;x2=7; y2=8;x3=10; y3=2; oldea=10000000;
	int cc=0;
	sa=sb=sc=10;
	step=1;

	for(a=-sa; a<sa; a+=step){
		for(b=-sb; b<sb; b+=step){
			for(c=-sc; c<sc; c+=step){
				e1 = a*36 + b*6 +c -1;
				e2 = a*64 + b*8 +c -7;
				e3 = a*4  + b*2 + c-10;
				et = e1 + e2 + e3;
				ea=fabs(et);
				printf("ea %f et %f oldea %f \n",ea,et,oldea); //fflush(stdout);
				if (ea< oldea) { 
					oldea=ea;
					na=a; nb=b; nc=c; //n=new, so i am saving the nearest guess
				}
				printf("sa %f sb %f sc %f \n",a,b,c); fflush(stdout);
			}
		}
		if (cc>10) break;
	}
	printf("\n a=%f b= %f c=%f", na,nb,nc);
}

Results

Code:

...
sa 9.000000 sb 9.000000 sc 7.000000 
ea 1086.000000 et 1086.000000 oldea 0.000000 
sa 9.000000 sb 9.000000 sc 8.000000 
ea 1089.000000 et 1089.000000 oldea 0.000000 
sa 9.000000 sb 9.000000 sc 9.000000 

 a=-1.000000 b= 8.000000 c=-2.000000

Dickfore · May 18, 2010

The matrix of the system is called Vandermonde matrix.

phizo · May 18, 2010

So for first point, 1,6, which in the program is 6,1 due to error

You have 6 = -1 x 1 + 8x1 -2 = 5, so not far out

Second 7,8 or 8,7 due to error

7 = -1 x 64 + 64 -2 = 5 so 2 out

Third 10,2 but 2,10 in program

10 = -1 x 4 + 16 + 2 = 14 so 4 out.

So it seems to work in a rough manner

Now if I plugged some new values in based on those results and also reduced the step size
I would hopefully home in on the solution. However if I get the estimates wrong I
might miss it.

I will try and modify it later to do that.

phizo · May 18, 2010

Dickfore said:

The matrix of the system is called Vandermonde matrix.

I hate matrices! :rofl:

Xitami · May 18, 2010

Dickfore said:

The matrix of the system is called Vandermonde matrix.

yes!

phizo · May 18, 2010

I need to make some adjustment to my prog, will try that tomorrow and see how I get on.

Martin Rattigan · May 19, 2010

You could use the Lagrange interpolation formula if you can guarantee that the [itex]x[/itex] coördinates of your points will be distinct.

That is, if your points are [itex](x_1,y_1),(x_2,y_2),(x_3,y_3)[/itex] with no two of [itex]x_1,x_2,x_3[/itex] the same

[tex]y=y_1\frac{(x-x_2)(x-x_3)}{(x_1-x_2)(x_1-x_3)}+y_2\frac{(x-x_1)(x-x_3)}{(x_2-x_1)(x_2-x_3)}+y_3\frac{(x-x_1)(x-x_2)}{(x_3-x_1)(x_3-x_2)}[/tex]

In the case [itex](0,0),(1,1),(2,2)[/itex] mentioned by Unit this will just give you a straight line. If you have duplicate [itex]x[/itex] coördinates there won't be an equation of the form you asked for that goes through the three points unless the corresponding [itex]y[/itex] coördinates are also duplicated, in which case you have only one or two points and you can drop the duplicate(s) and include an arbitrary extra point(s) with distinct [itex]x[/itex] coördinate(s).

The formula generalizes in an obvious way to a different number of points.

You're asking for a hyperbola with a vertical axis in general. Is that really what you want? E.g. a unique circle will go thro' three distinct points. Also if you may want to expand it to more points it would be a good idea to think in terms of spline functions. See http://en.wikipedia.org/wiki/Spline_(mathematics ).

theoarc · May 19, 2010

haha Phizo, that looks like a fun experiment. Not exactly sure it's what I am after... not exactly sure I understand it. You should do some speed tests on it to solve something simple like:
(2, 0) , (-2, 0) , (0, 4)

Currently I'm writing the code in AS3; however, in the end it will need to be in C so I'll probably get better speed later. Anyway using matrices as suggested by Xitami I was able to solve this rather easily, that was exactly what I was looking for.

There's a Matrix library for as3 called KMatrix that worked great. I'm sure there are many in C as well.

Anyway it took 5ms to solve for a, b, and c from the above numbers. So I'm happy with that.

I had seen the Lagrange interpolation formula. However, I didn't want to use it because it looks to require a lot more calculations in the long run. I could do a test and see which is faster over a period of time but I'm pretty happy with the Matrix solution at the moment. I haven't tested it in an actual useful situation yet so I'll find out soon.

Xitami · May 19, 2010

(x1,y1), (x2, y2), (x3, y3)
y=a*x*x+b*x+c[tex]a=\frac{ \left(y_{2} - y_{3}\right) x_{1}
+ \left( y_{1} - y_{2}\right) x_{3}
+ \left(y_{3} - y_{1}\right) x_{2}
}{
\left(x_{3}
- x_{2}\right) x_{1}^2
+ \left(-x_{3}^2
+ x_{2}^2\right) x_{1}
+ \left(x_{2} x_{3}^2
- x_{2}^2 x_{3}\right) }[/tex]

[tex]b= \frac{ \left(-y_{2}
+ y_{3}\right) x_{1}^2
+ \left(-y_{1} x_{3}
+ y_{1} x_{2}\right) x_{1}
+ \left( \left(-y_{1}
+ y_{2}\right) x_{3}^2
+ y_{1} x_{2} x_{3}
- y_{3} x_{2}^2\right) }{ \left(x_{3}
- x_{2}\right) x_{1}^2
+ \left(-x_{3}^2
+ x_{2}^2\right) x_{1}
+ \left(x_{2} x_{3}^2
- x_{2}^2 x_{3}\right) } [/tex]

[tex] c=\frac{ \left(y_{2} x_{3}
- y_{3} x_{2}\right) x_{1}^2
+ \left( \left(y_{1}
- y_{2}\right) x_{3}^2
- y_{1} x_{2} x_{3}
+ y_{3} x_{2}^2\right) x_{1}}{ \left(x_{3}
- x_{2}\right) x_{1}^2
+ \left(-x_{3}^2
+ x_{2}^2\right) x_{1}
+ \left(x_{2} x_{3}^2
- x_{2}^2 x_{3}\right) }[/tex]

[tex]\mu S[/tex] i think

theoarc · May 19, 2010

Xitami,

that's perfect!

I've only tested it once but it worked brilliantly.
Also takes less than 1ms to complete. That's at least 5x the speed of using the matrix class.

Thanks a lot :D

phizo · May 19, 2010

theoarc said:

haha Phizo, that looks like a fun experiment. Not exactly sure it's what I am after... not exactly sure I understand it. You should do some speed tests on it to solve something simple like:
(2, 0) , (-2, 0) , (0, 4)

Currently I'm writing the code in AS3; however, in the end it will need to be in C so I'll probably get better speed later. Anyway using matrices as suggested by Xitami I was able to solve this rather easily, that was exactly what I was looking for.

There's a Matrix library for as3 called KMatrix that worked great. I'm sure there are many in C as well.

Anyway it took 5ms to solve for a, b, and c from the above numbers. So I'm happy with that.

I had seen the Lagrange interpolation formula. However, I didn't want to use it because it looks to require a lot more calculations in the long run. I could do a test and see which is faster over a period of time but I'm pretty happy with the Matrix solution at the moment. I haven't tested it in an actual useful situation yet so I'll find out soon.

OK, I have made some modifications, there were some programing and logical error, and a few more remain but I think it came up with the right values for you inputs.

Code:

 best guess  a=-2.00    b= 0.00     c=8.00    ea 21.00    et 21.00    oldea 4.00     ***  sa -1.00    sb 0.00     sc -3.00    
 best guess  a=-2.00    b= 0.00     c=8.00    ea 18.00    et 18.00    oldea 4.00     ***  sa -1.00    sb 0.00     sc -2.00    
 best guess  a=-2.00    b= 0.00     c=8.00    ea 15.00    et 15.00    oldea 4.00     ***  sa -1.00    sb 0.00     sc -1.00    
 best guess  a=-2.00    b= 0.00     c=8.00    ea 12.00    et 12.00    oldea 4.00     ***  sa -1.00    sb 0.00     sc 0.00     
 best guess  a=-2.00    b= 0.00     c=8.00    ea 9.00     et 9.00     oldea 4.00     ***  sa -1.00    sb 0.00     sc 1.00     
 best guess  a=-2.00    b= 0.00     c=8.00    ea 6.00     et 6.00     oldea 4.00     ***  sa -1.00    sb 0.00     sc 2.00     
 best guess  a=-2.00    b= 0.00     c=8.00    ea 3.00     et 3.00     oldea 4.00     ***  sa -1.00    sb 0.00     sc 3.00     
 best guess  a=-1.00    b= 0.00     c=3.00    ea 0.00     et 0.00     oldea 3.00     ***  sa -1.00    sb 0.00     sc 4.00     
 best guess  a=-1.00    b= 0.00     c=4.00    ea 3.00     et 3.00     oldea 0.00     ***  sa -1.00    sb 0.00     sc 5.00     

...
 best guess  a=-1.00    b= 0.00     c=4.00    ea 86.00    et 86.00    oldea 0.00     ***  sa 9.00     sb 9.00     sc 6.00     
 best guess  a=-1.00    b= 0.00     c=4.00    ea 89.00    et 89.00    oldea 0.00     ***  sa 9.00     sb 9.00     sc 7.00     
 best guess  a=-1.00    b= 0.00     c=4.00    ea 92.00    et 92.00    oldea 0.00     ***  sa 9.00     sb 9.00     sc 8.00     
 best guess  a=-1.00    b= 0.00     c=4.00    ea 95.00    et 95.00    oldea 0.00     ***  sa 9.00     sb 9.00     sc 9.00     
 best guess  a=-1.00    b= 0.00     c=4.00    
****************************************** a=-1.000000 b= 0.000000 c=4.000000
****************************************** a=-1.000000 b= 0.000000 c=4.000000
****************************************** a=-1.000000 b= 0.000000 c=4.000000

New prog.

Code:

#include <stdio.h>
#include <math.h>
main(argc,argv)
	int  argc;
	char *argv[]; {
	// y = ax2 + bx + c	// (1,6) (7,8)(10,2) 
	float y1,y2,y3,x1,x2,x3, sa,sb,sc, sae,sbe,sce, na,nb,nc;
	float stepa, stepb, stepc, ea, oldea, a, b, c, e1, e2, e3, et;
	int c1;
	c1=0;
	x1=2;y1=0;
    x2=-2; y2=0;
	x3=0; y3=4; 

	oldea=10000000;
	sa=sb=sc=-10;
	sae=sbe=sce=10;

	stepa=1;
	stepb=1;
	stepc=1;	while(c1++<3){

	oldea=10000000;
		for(a=sa; a<sae; a+=stepa){
			for(b=sb; b<sbe; b+=stepb){
				for(c=sc; c<sce; c+=stepc){
					e1 = a*x1*x1 + b*x1 +c -y1;
					e2 = a*x2*x2 + b*x2 +c -y2;
					e3 = a*x3*x3 + b*x3 +c -y3;
					et = fabs(e1) + fabs(e2) + fabs(e3);
					ea=fabs(et);
					printf("ea %-08.2f et %-08.2f oldea %-08.2f ***  ",ea,et,oldea); //fflush(stdout);
					if (ea< oldea) { 
						oldea=ea;
						na=a; nb=b; nc=c; //n=new, so i am saving the nearest guess
					}
					printf("sa %-08.2f sb %-08.2f sc %-08.2f \n",a,b,c); 
				printf(" best guess  a=%-8.2f b= %-8.2f c=%-8.2f", na,nb,nc);

				}
			}
		}

		printf("\n****************************************** a=%f b= %f c=%f", na,nb,nc);

		sa=na*0.3;
		sb=nb*0.3;
		sc=nc*0.3;
		sae=na*2.5;
		sbe=nb*2.5;
		sce=nc*2.5;
		stepa=(float)(sa-sae)/30;
		stepb=(float)(sb-sbe)/30;
		stepc=(float)(sc-sce)/30;

	}

}

I would need some better test data as part of the program is not used because it
finds the result easilly.

phizo · May 19, 2010

Hence I would appreciate some test data!

phizo · May 19, 2010

I tried it with
x1=2; y1=0;
x2=-2; y2=0;
x3=5; y3=4;

And it came up with a=0.000000 b= 0.000000 c=0.000000

Maybe there is no solution or more likely errors on my part. (things never work first effort).

Well it did also say oldea=4, which means it did not find a solution as the error was 4 not 0.

Xitami · May 19, 2010

phizo said:

I tried it with
x1=2; y1=0;
x2=-2; y2=0;
x3=5; y3=4;
...

y=4/21x^2 +0x -16/21

phizo · May 19, 2010

Yes there is a problem when it's best guess for a b or c is 0, because I then do 2.5 times 0
and 0.3 times zero, which both give me zero as an answer for my test range!

Xitami · May 19, 2010

x1=17; y1=18; x2=19; y2= 9; x3=20; y3= 2; /*-0.83333333 x^2 + 31.500000 x + -294.66667 */
x1= 5; y1=14; x2= 7; y2=18; x3= 8; y3= 4; /*-5.3333333 x^2 + 70.666667 x + -220.00000 */
x1=17; y1= 6; x2=19; y2=12; x3=20; y3= 0; /*-5.0000000 x^2 + 185.00000 x + -1700.0000 */
x1= 7; y1= 8; x2= 9; y2= 9; x3=10; y3=19; /*3.1666667 x^2 + -47.500000 x + 177.33333 */
x1=13; y1= 2; x2=15; y2=12; x3=16; y3=15; /*-0.66666667 x^2 + 24.333333 x + -203.66667 */
x1=14; y1=19; x2=16; y2=13; x3=17; y3=15; /*1.6666667 x^2 + -46.666667 x + 326.66667 */
x1=20; y1=11; x2=22; y2=14; x3=23; y3= 6; /*-3.1666667 x^2 + 138.16667 x + -1496.6667 */
x1= 6; y1=17; x2= 8; y2=12; x3= 9; y3= 8; /*-0.50000000 x^2 + 10.166667 x + -43.000000 */
x1= 8; y1=12; x2=10; y2= 2; x3=11; y3=12; /*5.0000000 x^2 + -91.000000 x + 408.00000 */
x1= 9; y1=16; x2=11; y2=18; x3=12; y3= 3; /*-5.3333333 x^2 + 113.00000 x + -585.00000 */
x1=19; y1= 0; x2=21; y2=18; x3=22; y3=13; /*-4.6666667 x^2 + 195.66667 x + -2033.0000 */
x1= 8; y1=18; x2=10; y2=8; x3=11; y3= 1; /*-0.66666667 x^2 + 13.000000 x + -61.333333 */
x1=12; y1= 8; x2=14; y2= 3; x3=15; y3=10; /*3.1666667 x^2 + -82.166667 x + 530.00000 */
x1=14; y1= 8; x2=16; y2=10; x3=17; y3= 8; /*-1.0000000 x^2 + 33.666667 x + -275.33333 */
x1=14; y1= 1; x2=16; y2=14; x3=17; y3= 4; /*-5.5000000 x^2 + 171.83333 x + -1327.6667 */
x1= 1; y1= 1; x2= 3; y2=19; x3= 4; y3= 6; /*-7.3333333 x^2 + 38.666667 x + -31.333333 */
x1=10; y1= 4; x2=12; y2=10; x3=13; y3= 0; /*-4.3333333 x^2 + 99.666667 x + -563.33333 */
x1= 8; y1=15; x2=10; y2=19; x3=11; y3= 6; /*-5.0000000 x^2 + 97.000000 x + -456.00000 */
x1=16; y1=16; x2=18; y2=6; x3=19; y3= 0; /*-0.33333333 x^2 + 11.666667 x + -101.33333 */
x1=15; y1=14; x2=17; y2= 3; x3=18; y3= 8; /*3.5000000 x^2 + -112.83333 x + 905.00000 */

phizo · May 19, 2010

Xitami said:

y=4/21x^2 +0x -16/21

Well it did seem 0 0 0 was not far off the solution!

Mark44 · May 19, 2010

Xitami said:

y=4/21x^2 +0x -16/21

phizo said:

Well it did seem 0 0 0 was not far off the solution!

Only if you think that 4/21 and -16/21 are both close to zero.

phizo · May 19, 2010

Mark44 said:

Only if you think that 4/21 and -16/21 are both close to zero.

Well when I was using a step of 1 they are not far off, because I would be trying 1, 0 an -1 and seeing which one produced the smallest error.
If think if you do that by hand 0 may well produce the smallest error.

Martin Rattigan · May 19, 2010

theoarc said:

I had seen the Lagrange interpolation formula. However, I didn't want to use it because it looks to require a lot more calculations in the long run.

Assuming you want to calculate many [itex](x,y)[/itex] for each set of three given points you would rearrange the formula as a quadratic, which should then look very much like Xitami's following post. But preferably as [itex](ax+b)x+c[/itex] to reduce the number of calculations per point. (Two multiplications and two additions as opposed to three multiplications and two additions.)

Mark44 · May 19, 2010

phizo said:

Well when I was using a step of 1 they are not far off, because I would be trying 1, 0 an -1 and seeing which one produced the smallest error.
If think if you do that by hand 0 may well produce the smallest error.

I guess your definition of "not far off" and mine are different. Let us know when you get an algorithm with more than 1 decimal place of precision.

phizo · May 19, 2010

Mark44 said:

Only if you think that 4/21 and -16/21 are both close to zero.

My latest effort gets:-

a=0.166590 = 0.070498 c=-0.517241 the actual decimal values are.

a = 0.1904762 b =0 c = -0.7619048

So not right, but not a million miles away either.

Still wrong though :cry:

phizo · May 19, 2010

And my final effort!

A=0.190478 B= 0.000001 B=-0.761904( by calculator a = 0.1904762 b =0 c = -0.7619048 )

Mark44 · May 19, 2010

Much improved!

phizo · May 19, 2010

NEW RUN 1 New sa=-1.000000 b= -1.000000 c=-1.000000 ends sa=1.000000 b= 1.000000 c=1.000000
NEW RUN 2 New sa=0.133333 b= -0.066667 c=-0.866667 ends sa=0.266667 b= 0.066667 c=-0.733333
NEW RUN 3 New sa=0.186667 b= -0.004445 c=-0.768889 ends sa=0.195556 b= 0.004444 c=-0.760001
NEW RUN 4 New sa=0.190222 b= -0.000296 c=-0.762371 ends sa=0.190815 b= 0.000296 c=-0.761778
NEW RUN 5 New sa=0.190459 b= -0.000020 c=-0.761937 ends sa=0.190499 b= 0.000020 c=-0.761897
NEW RUN 6 New sa=0.190475 b= -0.000001 c=-0.761907 ends sa=0.190478 b= 0.000001 c=-0.761904

Quadratic Equation from 3 points

Similar threads

Hot Threads

Recent Insights