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

[Sergei Steshenko]

----- Original Message -----
> From: Sergei Steshenko <address@hidden>
> To: Joza <address@hidden>; "address@hidden" <address@hidden>
> Cc: 
> Sent: Saturday, November 3, 2012 12:44 AM
> Subject: Re: Successive Over-Relaxation ... what is wrong? Improvements?
> 
> 
> 
> 
> 
> 
> ----- Original Message -----
>>  From: Joza <address@hidden>
>>  To: address@hidden
>>  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.


And, by the way, this is completely wrong to debug such problems on data having 
several thousand elements - because one will quickly get lost in this number of 
values.

Unless he/she is a genius making no mistakes.

So, if I understand correctly, you should first reduce your problem to, say, 4 
x 4 sparse matrix and using debugger and/or print statements debug it.

Regards,
  Sergei.


>