What are you doing? You used a simple forward difference to compute the derivative, as an approximation to the derivative of sin(x). You did so while decreasing the step size, by a factor of 4 after every iteration.
At some point, the error actually starts to increase. Why? Look at what happened.
Now, I'll just look at some of the steps in that loop.
ind = [1 5 10 15 20 25:30]';
So, for the first 15 iterations, the error is getting smaller. But when h gets smaller than roughly 1e-8 or so, things start to get bad. The error actually gets worse, and steadily so.
What you need to understand is that for small values of h, you are subtracting two values that are VERY close to each other. For example, look t what happens when h is as small as 1e-12.
Now subtract them.
Finally, divide the difference by dx.
(sin(x0 + dx) - sin(x0))/dx
And compare that to cos(x0)
You should see the two results are starting to deviate after around 4 digits or so there.
The problem is, when dx gets very, very, very tiny, that difference between sin(x+dx) and sin(x) starts to get less accurate. Double precision arithmetic starts to run out of room to compute that difference exactly.
In the case of this problem, the optimal value for dx should be roughly sqrt(eps).
So go back, and look very carefully at what happened. These two values:
Are the same, almost identically, out for the first 12 digits or so. So now, when you subtract them, the result has only roughly 4 significant digits that mean anything.
Even though we see more than 4 digits in the result, after roughly 4 digits there, those digits are complete garbage.
This is usually called massive subtractive cancellation. Remember that a double precision number is only capable of carrying roughly 16 decimal digits or so.