Code for numerical differentiation doesn't work?

[AyarasSivara]

I am having some problems with my code on numerical differentiation. We are supposed to use the two point, three point, four point, and five point formulas to find the first derivative of x^3 - sin x.

I made the following code:

#include <iostream>
#include <iomanip>
#include <cmath>

using namespace std;

double h1 = 0.01, h2 = 0.001, ActualValue;

double f (double x) //stores the function
{
    double r;
    r = pow (x,3) - sin (x);
    return r;
}



int main ()
{
    ActualValue = 3.0 - cos (1.0);
    cout << "Table of Errors for the First Derivative of x^3 - sin x" << endl;
    cout << setprecision (10);
    cout << "Actual Value = " << ActualValue << endl;
    cout << " " << endl;
    cout << " \t \t \t h = 0.01 \t h = 0.001" << endl;
    cout << "2-Point Left \t \t" << ActualValue - (f(1.0) - f(1.0 - h1))/h1 << "\t" << ActualValue - (f(1.0) - f(1.0 - h2))/h2 <<endl; //Two Point Left formula
    cout << "2-Point Right \t \t" << fabs(ActualValue - (f(1.0 + h1) - f(1.0))/h1) << "\t" << fabs(ActualValue - (f(1.0 + h2) - f(1.0))/h2) << endl; //Two Point Right formula
    cout << "3-Point       \t \t" << ActualValue - (f(1.0 + h1) - f (1.0 - h1))/2 << "\t" << ActualValue - (f(1.0 + h2) - f (1.0 - h2))/2 << endl;
    cout << "4-Point Positive  \t" << ActualValue - (((-2 * f(1.0 - h1)) - 3 * f(1.0) + 6 * f(1.0 + h1) - f(1.0 + (2.0*h1)))/6)<< "\t" <<  ActualValue - (((-2 * f(1.0 - h2)) - 3 * f(1.0) + 6 * f(1.0 + h2) - f(1.0 + (2.0*h2)))/6)<< endl;
    cout << "4-Point Negative  \t" << - ((-2 * f(1.0 + h1)) - 3 * f(1.0) + 6 * f(1.0 - h1) - f(1.0 - (2.0*h1)))/6 << "\t" << - ((-2 * f(1.0 + h2)) - 3 * f(1.0) + 6 * f(1.0 - h2) - f(1.0 - (2.0*h2)))/6 << endl;
    cout << "5-Point           \t" << f(1.0 - (2 * h1)) - 8 * f(1.0 - h1) + 8 * f(1.0 + h1) - f(1.0 - 2 * h1) << "\t" << f(1.0 - (2 * h1)) - 8 * f(1.0 - h2) + 8 * f(1.0 + h2) - f(1.0 - 2 * h2) << endl;
}

It works for the first two formulas, but doesn't work after that. THere doesn't seem to be anything wrong with the formulas I input, but it won't run properly. There is a given value for errors, but it doesn't match up. For example, for the three point formula, it must be 0.0001. However, my code is giving 2.43.

[Ganado]

f(1.0 + h2) - f(1.0 - 2 * h2)endl;

Your code doesn't even compile, so it's hard to tell what you're actually doing wrong. My guess is that there's some typo in all the math you have, perhaps a wrong sign or misplaced parenthesis.

Please post your actual code, and use code formatting tags.
Add [code] and [/code] around your code.

What does "it works"/"doesn't work" mean? Be specific.

[AyarasSivara]

Oh sorry! This is my first time here so I'm sorry if I am doing things wrong. I'll edit it right now!

[dhayden]

Reformatting part of the code:

    cout << "3-Point       \t \t"
         << ActualValue - (f(1.0 + h1) - f(1.0 - h1)) / 2 << "\t"
         << ActualValue - (f(1.0 + h2) - f(1.0 - h2)) / 2 << endl;

Shouldn't that denominator be (2*h)?

Please post a link to the formulas that you're trying to code. Meanwhile, here is your code formatted for readability.

#include <iostream>
#include <iomanip>
#include <cmath>

using namespace std;

double h1 = 0.01, h2 = 0.001, ActualValue;

double
f(double x)			//stores the function
{
    double r;
    r = pow(x, 3) - sin(x);
    return r;
}

int
main()
{
    ActualValue = 3.0 - cos(1.0);
    cout << "Table of Errors for the First Derivative of x^3 - sin x" << endl;
    cout << setprecision(10);
    cout << "Actual Value = " << ActualValue << endl;
    cout << " " << endl;
    cout << " \t \t \t h = 0.01 \t h = 0.001" << endl;
    cout << "2-Point Left \t \t"
	 << ActualValue - (f(1.0) - f(1.0 - h1)) / h1 << "\t"
	 << ActualValue - (f(1.0) - f(1.0 - h2)) / h2 << endl;	//Two Point Left formula
    cout << "2-Point Right \t \t"
	 << fabs(ActualValue - (f(1.0 + h1) - f(1.0)) / h1) << "\t"
	 << fabs(ActualValue - (f(1.0 + h2) - f(1.0)) / h2) << endl;	//Two Point Right formula
    cout << "3-Point       \t \t"
	 << ActualValue - (f(1.0 + h1) - f(1.0 - h1)) / 2 << "\t"
	 << ActualValue - (f(1.0 + h2) - f(1.0 - h2)) / 2 << endl;
    cout << "4-Point Positive  \t"
	 << ActualValue - (((-2 * f(1.0 - h1)) - 3 * f(1.0) + 6 * f(1.0 + h1) - f(1.0 + (2.0 * h1))) / 6) << "\t"
	 << ActualValue - (((-2 * f(1.0 - h2)) - 3 * f(1.0) + 6 * f(1.0 + h2) - f(1.0 + (2.0 * h2))) / 6) << endl;
    cout << "4-Point Negative  \t"
	 << -((-2 * f(1.0 + h1)) - 3 * f(1.0) + 6 * f(1.0 - h1) - f(1.0 - (2.0 * h1))) / 6 << "\t"
	 << -((-2 * f(1.0 + h2)) - 3 * f(1.0) + 6 * f(1.0 - h2) - f(1.0 - (2.0 * h2))) / 6 << endl;
    cout << "5-Point           \t"
	 << f(1.0 - (2 * h1)) - 8 * f(1.0 - h1) + 8 * f(1.0 + h1) - f(1.0 - 2 * h1) << "\t"
	 << f(1.0 - (2 * h1)) - 8 * f(1.0 - h2) + 8 * f(1.0 + h2) - f(1.0 - 2 * h2) << endl;
}

[lastchance]

Agree with @dhayden: h1 and h2 are missing from the denominators in the 3, 4 and 5-point schemes.

You are also calculating the error (i.e. subtracting from the actual derivatives) in the first four schemes, and then just displaying the estimate of the derivative in the last two.

Why don't you have functions like

double approx3pt( double x, double h )
{
    return ( f(x+h) - f(x-h) ) / ( 2.0 * h );
}

and similarly for the other approximations. Then you won't be calculating essentially the same thing twice, and you won't be tied to x=1.

Your 3-point schemes should be divided by 2h, not 2.

Your 4-point schemes should be divided by 6h, not 6.

Your 5-point scheme has an error of sign in the last term on each line and you have put h1 where it should be h2 in one part of the last expression. It should also be divided by 12h (I think). Something like

double approx5pt( double x, double h )
{
    return ( f(x-2*h) - 8*f(x-h) +8*f(x+h) - f(x+2*h) ) / ( 12.0 * h );
}

Actual derivative: 2.4597

Errors in approximations:
2-point right: -3.431632e-02,  -3.421826e-03
2-point left : +3.409831e-02,  +3.419645e-03
3-point      : -1.090050e-04,  -1.090050e-06
4-point right: -7.030152e-08,  -6.998530e-11
4-point left : +6.994131e-08,  +7.010634e-11
5-point      : -1.801046e-10,  +6.977318e-14