help-octave
[Top][All Lists]

## Re: Successive Over-Relaxation ... what is wrong? Improvements?

 From: Sergei Steshenko Subject: Re: Successive Over-Relaxation ... what is wrong? Improvements? Date: Fri, 2 Nov 2012 15:44:53 -0700 (PDT)

```

----- Original Message -----
> Cc:
> Sent: Saturday, November 3, 2012 12:25 AM
> Subject: Successive Over-Relaxation ... what is wrong? Improvements?
>
> Below is an implementation I have written of the Successive Over-Relaxation
> Method (for pedagogical purposes). I must use a starting guess of zeros and
> b in Ax=b is all ones. A is sparse.
>
> I've run Octave's CGS and it solved the problem in about 1 second ... my
> code has been running for well over 5 minutes and nothing ...
>
> I'm sure there is something I'm missing here. Any ideas? Feel free to
> run
> the code yourself ... w is the relaxation parameter ( try 1.35) and tol is
> the convergence tolerance.
>
> This is the second major algorithm I've had trouble with recently. Must be
> something up with mu approach! Anyway, I'd really like to hear all
> critiques!
> ************************************************************************************************************************************************
>
> function [x, iters] = SOR(w, tol)
>
> r = sparse([1 1], [1 2], [2 -1], 1, 1000);
> A = toeplitz(r);
> b = ones(1000,1);
> x = zeros(1000,1);
>
> iters = 0;
> converged = 0;
>
> while !converged
>
>     norm_x = 0.0;
>     norm_delta_x = 0.0;
>     for i=1:1:1000;
>         sum = 0.0;
>
>         if i==1
>             for j=1:2
>                 sum = sum + A(i,j)*x(j);
>             end
>
>                 elseif i==1000
>             for j=999:1000
>                 sum = sum + A(i,j)*x(j);
>             end
>
>                 else
>             for j=i-1:i+1
>                 sum = sum + A(i,j)*x(j);
>             end
>         end
>
>         delta_x = w*(b(i) - sum)/A(i,i);
>         x(i) = x(i) + delta_x
>         iters = iters + 1;
>
>         norm_x = max( norm_x, abs(x(i)) );
>         norm_delta_x = max(norm_delta_x, abs(delta_x) )
>     end
>
>     if norm_delta_x <= tol + 4.0*(10^-16)*norm_x
>         converged = 1;
>     end
> end
>

Without going into details I can definitely tell you shouldn't be writing code
this way.

In any decent code review lines like:

for i=1:1:1000;
for j=999:1000

will be ruthlessly ruled out. This is because you are using constants, though
you most likely mean to use something which is a function of matrices sizes.

I.e. your code is error prone because you need to change it whenever matrices
sizes change.

And, by the way, I do not know ';' in

for i=1:1:1000;

does.

Regards,
Sergei.

```