#include
#include
/* A program to illustrate computer rounding. One annoying thing about C
relative to FORTRAN is the fact that C does not have h^n (power) as a
built in operation. This program is "quick and dirty". I don't think
I will need to look at it later so I haven't written it carefully. */
float summand(int i, int j );
int main() {
double h, r_direct, r_asymp, r_rewrite, d2, d1;
double fact_2, fact_4, fact_6, fact_8, fact_10, fact_12, fact_14;
fact_2 = 2;
fact_4 = fact_2*3*4; // fact_k is k!.
fact_6 = fact_4*5*6; // The computer automatically converts
fact_8 = fact_6*7*8; // the integers 2, 4, ... into doubles
fact_10 = fact_8*9*10; // so this gives the most accurate double
fact_12 = fact_10*11*12; // precision values.
fact_14 = fact_12*13*14;
while(1) { // This is an infinite loop, type ^C to get out.
cout << "Type h: ";
cin >> h;
r_direct = ( 1. - cos(h) ) / (h*h);
r_rewrite = 2.*sin(h/2)*sin(h/2) / (h*h);
r_asymp = 1./fact_2 - h*h/fact_4 + h*h*h*h/fact_6 // This is the
- h*h*h*h*h*h/fact_8 + h*h*h*h*h*h*h*h/fact_10 // Taylor series
- h*h*h*h*h*h*h*h*h*h/fact_12 // approximation.
+ h*h*h*h*h*h*h*h*h*h*h*h/fact_14;
d1 = r_direct - r_asymp;
d2 = r_rewrite - r_asymp;
cout << "For h = " << h << " r_direct = " << r_direct <<
", and r_direct - r_asymp = " << d1 << endl;
cout << "For h = " << h << " r_rewrite = " << r_rewrite <<
", and r_rewrite - r_asymp = " << d2 << endl;}
}