MATLAB: MATLAB program in a script file that determines

while loop

The equation we are determing pi from is pi^3/32=sigma(n=0 to infinity) -1^n/(2n+1)^3 The program must use a function error and must end when the error>= 1.0e-6, i am struggling to find the problem with my code as I either get an infinite loop and n will not move up from 1. Appreciate your help and this must be done with a while loop
error=1.0e-6;
pivalue=0;
n=0;
while pi>=sum(pivalue)
n=n+1;
sigma(n)=(-1).^(n)/(2.*n+1).^3;
pivalue= abs(sum(sigma(n)).*32).^(1/3);
if error<= abs(pi - pivalue)
break
end
end

Best Answer

Ok, you made a decent start, so lets see what went wrong.
First, try not to use variables named error, or any other existing, useful function names. I'd suggest a different name, perhaps errortol, that indicates better what the variable means, here an error tolerance.
Next, you have a test in the while loop definition. But your test is p>=pivalue. Since this is an alternating series, I expect the result to flip-flop over and under pi. So >= is just abad idea. Anyway, you don't need it, as you then test inside the loop, and use a break.
Next, you are trying to stuff the terms into a vector sigma that will grow with each iteration. growing a vector dynamically is usually a bad idea, as that forces MATLAB to reallocate memory for the vector at every iteration. Luckily, we know this cannot take too many iterations. But, better is to preallocate sigma if you want to keep the individual terms in the sum. I'm not sure why you would bother though. Anyway, that way you are forced to re-sum all of the terms at EVERY iteration. So that is quite inefficient.
A big problem was you were computing pivalue wrong. Here is the offending line:
pivalue= abs(sum(sigma(n)).*32).^(1/3);
See that you formed this subexpression:
sum(sigma(n))
n is a SCALAR! so sigma(n) is just the last term in that vector. If you wanted to form the sum, sum(sigma) would have been better.
Finally, you had one other problem, and this is why I wanted to change the name of error to errortol. It makes it clearer to you what the test inside the loop should be.
if errortol >= abs(pi - pivalue)
You used <=, but really, you need to stop only when the tolerance is greater than the error from the true value of pi.
I might use nthroot, instead of taking a cuberoot using .^(1/3). nthroot is a bit safer. The point is, if the number is negative, then the form you used will compute ONE of the nth roots of a negative number, which will just happen to be a complex number.
nthroot(-1,3)
ans =
-1
(-1).^(1/3)
ans =
0.5 + 0.866025403784439i
As I said, nthroot will not surprise you.
oops. One other thing. In that series, the first term starts at n=0. You increment n before computing the term,so you are effectively starting at n=1.
There are many things to think of when you write good code. You will learn to be careful, but that takes some practice.
So here is how I would have changed your code:
errortol=1.0e-6;
pivalue=0;
sigma=0;
n=0;
% the while loop could run forever like this, but we
% know that it won't, because of that test inside.
while true
sigma = sigma + (-1).^(n)/(2.*n+1).^3;
pivalue = nthroot(sigma.*32,3);
if errortol >= abs(pi - pivalue)
break
end
n=n+1;
end
I will point out that my version is very little different from yours. You had the right idea all along. You just needed some minor fixes in the implementation.
When I run this, I see it has stopped at n=40;
n
n =
40
pivalue
pivalue =
3.14159363278731
Yes, there were a bunch of things I had to fix, but they were all subtle changes, worth thinking about why I did it that way.
Having said all of that, how might I have written this in a more classical MATLAB style? Yes, I know that you need to use a loop here. That teaches you much. But you can also learn by comparison.
n = 0:100;
piest = nthroot(cumsum((-1).^n./(2*n+1).^3)*32,3);
semilogy(n,abs(pi-piest))
grid on
Here you can see clearly at what term the approximation crosses below 1e-6.
It is worth looking at what I wrote, and understanding that too.