Homing in on an answer¶

In the lab for this week, you are using a scheme invented by Newton to find the cube root of a number. It starts off with some guess, then “iterates” (loops) working on a better and better answer until, hopefully, it “homes in” on a good answer.

What is a good answer?¶

In working with computers, you can never expect to get exactly the right answer except in special cases. Cube roots are not one of those cases. So what is good enough. It must be when the answer is as close as you need it to be. How to we figure that out?

In the case of the cube root, let’s say you have a guess and want to know if it is close.

Take the guess and cube it to get a number. Is this number equal to the number you started with? Yes? Great, you accidentally got the right answer. More likely, it is off by some amount. Let’s call that amount the error.

Here is some code to do this:

#include <iostream>
using namespace std;

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

float number = 27.0;
float guess = number / 2.0;
float calculated_number;
float error;

calculated_number = guess*guess*guess;
error = calculated_number - number;

cout << "The cube root of " << number;
cout << " is about " << guess;
cout << "(the error is " << error << ")" << endl;
}

Running this will tell us how good this is:

The cube root of 27 is about 13.5(the error is 2433.38)

Wow, that was not very good. Well, we called it a guess, didn’t we! Apparently it was not good enough. Here is where Newton comes in. His scheme tells us how to get a better guess:

guess = 1.0/3.0 * ( number / (guess*guess) + 2.0 * guess);
calculated_number = guess*guess*guess;
error = calculated_number - number;

cout << "The cube root of " << number;
cout << " is about " << guess;
cout << "(the error is " << error << ")" << endl;

Why is this written the way it is? Leave off the decimal points in the 1/3 and see what happens. Everything in this “expression” is a floating point number, so the math will be right.

Running this give this:

The cube root of 27 is about 13.5(the error is 2433.38)
The cube root of 27 is about 9.04938(the error is 714.066)

Well, this is better, but still off. We need a loop!

Working harder¶

At this point, you should see that we have to do this work more times to get a better and better answer. Hopefully, that error term will keep getting smaller. If so we need to decide when to stop.

How might we do that?

One way would be to make this loop happen a fixed number of times. Can we do that. Sure. All we need is a variable to keep track of the number of times we go through the loop, then we need to add one to that variable each time we go through the loop. We stop when we hit the number we want.

Something like this:

int count = 0;
...
while (count<10) {
//do something
count = count + 1;
}

This will do “something” 10 times. Let’s try this idea:

#include <iostream>
using namespace std;

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

float number = 27.0;
float guess = number / 2.0;
float calculated_number;
float error;
int count = 0;

calculated_number = guess*guess*guess;
error = calculated_number - number;

while (count < 10) {
cout << "The cube root of " << number;
cout << " is about " << guess;
cout << "(the error is " << error << ")" << endl;

guess = 1.0/3.0 * ( number / (guess*guess) + 2.0 * guess);
calculated_number = guess*guess*guess;
error = calculated_number - number;

count = count + 1;
}
}

Look closely at how things changed! We want to see how things are going as the loop spins, so we output something each time we go through the loop body.

Here is the output:

The cube root of 27 is about 13.5(the error is 2433.38)
The cube root of 27 is about 9.04938(the error is 714.066)
The cube root of 27 is about 6.14282(the error is 204.795)
The cube root of 27 is about 4.33373(the error is 54.3925)
The cube root of 27 is about 3.36835(the error is 11.2167)
The cube root of 27 is about 3.03881(the error is 1.06159)
The cube root of 27 is about 3.00049(the error is 0.0133343)
The cube root of 27 is about 3(the error is 0)
The cube root of 27 is about 3(the error is 0)
The cube root of 27 is about 3(the error is 0)

Wow, we seem to have reached the answer in less than 10 trys! Not bad Mr. Newton!

How might we stop this loop as soon as we hit the accuracy we want? There is no reason to expect all numbers will do as well (try 753).

The problem with this is that we want the error to be small. How small is the question. Let’s see if we can get withing 4 decimal places or error < 0.0001:

Let’s see if we can change the loop to stop when the error is small. Something like this:

error = 10;     // we want to get into the loop!
while (error > 0.0001) {
// do something and calculate a new value for error
}

This loop will only stop if the error is less than 0.0001. There is a problem with this. What happens if the error is negative and gets more negative. You will stop with a lousy guess. It is better to use the absolute value in the error calculation, which we can do if we include the math.h library:

error = fabs(guess*guess*guess - number);

That fabs module returns the value of the parameter as a positive number even if it was negative.

Now we do not care if the error is positive or negative.

When I coded this up, I got this for 753:

The cube root of 753 is about 376.5(the error is 5.3369e+07)
The cube root of 753 is about 251.002(the error is 1.58128e+07)
The cube root of 753 is about 167.339(the error is 4.68509e+06)
The cube root of 753 is about 111.568(the error is 1.38798e+06)
The cube root of 753 is about 74.3988(the error is 411058)
The cube root of 753 is about 49.6446(the error is 121600)
The cube root of 753 is about 33.1982(the error is 35835.4)
The cube root of 753 is about 22.3599(the error is 10426.1)
The cube root of 753 is about 15.4086(the error is 2905.4)
The cube root of 753 is about 11.3296(the error is 701.261)
The cube root of 753 is about 9.5085(the error is 106.679)
The cube root of 753 is about 9.11519(the error is 4.35181)
The cube root of 753 is about 9.09773(the error is 0.00836182)

Hmmm, why is that last error not less than 0.0001? Think about that to see why!

Durn, I just gave you most of the solution to Lab9!