//this is a program to implement the exercise on pp. 34-35 in Landau
//this program demonstrates substractive cancellation in the quadratic equation

#include <stdio.h>
#include <math.h>

double a,b,c;//parameter for the quadratic equations
double xx1,xx2;//solutions #1
double yy1,yy2;//alternative solutions #2

FILE *outfile1,*outfile2;

main()
{
	double temp;
	double eps;
	int n;

	a=1.0;
	b=1.0;
	c=1.0;

	outfile1=fopen("quad.dat","w");
	outfile2=fopen("quad_err.dat","w");

	for(n=1;n<=20;n++)
	{
		c=c/10.0;

		eps=4.0*a*c;
		temp=sqrt(b*b-eps);

		xx1=(-b+temp)/(2*a);
		xx2=(-b-temp)/(2*a);

		yy1=-2*c/(b+temp);
		yy2=-2*c/(b-temp);

		//printing out the four solutions
		fprintf(outfile1,
			"%22.18lg\t%22.18lg\t%22.18lg\t%22.18lg\n",
			xx1,xx2,yy1,yy2);
		//printing out their percentage errors and the value of 4ac
		//NOTE:  since b>1, subtractractive cancellation occurs
		//	only for xx1 and yy2
		//	So, one can safely consider yy1 and xx2 as accurate solutions
		fprintf(outfile2,
			"%22.18lg\t%22.18lg\t%22.18lg\n",
			(yy1-xx1)/yy1,(xx2-yy2)/xx2,eps);

	}
	return 0;
}