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

27 messages
12
Open this post in threaded view
|

## 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
Open this post in threaded view
|

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

 I should add, A only has elements on the main diagonal, and the next upper and lower diagonals. All other entries are zero.
Open this post in threaded view
|

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

 In reply to this post by Joza ----- Original Message ----- > From: Joza <[hidden email]> > To: [hidden email] > 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. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 ----- Original Message ----- > From: Sergei Steshenko <[hidden email]> > To: Joza <[hidden email]>; "[hidden email]" <[hidden email]> > Cc: > Sent: Saturday, November 3, 2012 12:44 AM > Subject: Re: Successive Over-Relaxation ... what is wrong? Improvements? > > > > > > > ----- Original Message ----- >>  From: Joza <[hidden email]> >>  To: [hidden email] >>  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. > _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 In reply to this post by Joza Sorry, I forgot I shouldn't have such numbers there ... I'll replace the likes of j=999:1000 with j = n-1:n. The code works fine for a 10x10 matrix that I tried. But Octave won't churn out an answer with 1000x1000. I'm quite new to sparse matrices so I'm wondering if there are ways to speed up the code or optimize it, perhaps by taking advantage of the sparse matrix's properties?
Open this post in threaded view
|

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

 What about just adding a print statement which shows the value of iter and norm_delta_x and tol + 4.0*(10^-16)*norm_x to see if you ever come near to your convergence criteria? _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 Well, here's the read-out for the last iteration: iters: 20374    norm: 2.0769230769230695 yet iters: 100      norm: 2.0769230769230771
Open this post in threaded view
|

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

 In reply to this post by Joza ----- Original Message ----- > From: Joza <[hidden email]> > To: [hidden email] > Cc: > Sent: Saturday, November 3, 2012 1:34 AM > Subject: Re: Successive Over-Relaxation ... what is wrong? Improvements? > > Sorry, I forgot I shouldn't have such numbers there ... I'll replace the > likes of j=999:1000 with j = n-1:n. > > The code works fine for a 10x10 matrix that I tried. But Octave won't churn > out an answer with 1000x1000. I'm quite new to sparse matrices so I'm > wondering if there are ways to speed up the code or optimize it, perhaps by > taking advantage of the sparse matrix's properties? > > > > -- > View this message in context: > http://octave.1599824.n4.nabble.com/Successive-Over-Relaxation-what-is-wrong-Improvements-tp4646041p4646046.html> Sent from the Octave - General mailing list archive at Nabble.com. Insert some diagnostic print statements allowing you to monitor progress. Depending on how you count, you may say that 100 x 100 in the middle of 10 x 10 and 1000 x 1000 piece. So, try 100 x 100 first. Regards,   Sergei. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 In reply to this post by Joza On Fri, Nov 2, 2012 at 5:13 PM, Joza wrote: Well, here's the read-out for the last iteration: iters: 20374    norm: 2.0769230769230695 yet iters: 100      norm: 2.0769230769230771 -- View this message in context: http://octave.1599824.n4.nabble.com/Successive-Over-Relaxation-what-is-wrong-Improvements-tp4646041p4646048.html Sent from the Octave - General mailing list archive at Nabble.com. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave It looks like it is oscillating about the solution, which probably means you need to fiddle withthe relaxation parameter-- Ed Meyer _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 In reply to this post by Joza --- On Fri, 11/2/12, Ed Meyer <[hidden email]> wrote: From: Ed Meyer <[hidden email]> Subject: Re: Successive Over-Relaxation ... what is wrong? Improvements? To: "Joza" <[hidden email]> Cc: [hidden email] Date: Friday, November 2, 2012, 8:37 PM On Fri, Nov 2, 2012 at 5:13 PM, Joza <[hidden email]> wrote: Well, here's the read-out for the last iteration: iters: 20374    norm: 2.0769230769230695 yet iters: 100      norm: 2.0769230769230771 -- View this message in context: http://octave.1599824.n4.nabble.com/Successive-Over-Relaxation-what-is-wrong-Improvements-tp4646041p4646048.htmlSent from the Octave - General mailing list archive at Nabble.com. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octaveIt looks like it is oscillating about the solution, which probably means you need to fiddle with the relaxation parameter -- Ed Meyer -----Inline Attachment Follows----- _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octaveIIRC, the OP says Octave built-in algorithm converges, so the problem is most likely an error in the OP's implementation. Regards,   Sergei. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 Yea.... I've tried a few values for the relaxation parameter and the matrix size ... it seems that either the norm converges to a constant, or it keeps oscillating around a certain value. As far as I am aware, I've followed the implementation as per textbook - I can't see anything obviously wrong with it. Any idea what could be causing this?
Open this post in threaded view
|

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

 ----- Original Message ----- > From: Joza <[hidden email]> > To: [hidden email] > Cc: > Sent: Sunday, November 4, 2012 3:12 AM > Subject: Re: Successive Over-Relaxation ... what is wrong? Improvements? > > Yea.... > > I've tried a few values for the relaxation parameter and the matrix size ... > > it seems that either the norm converges to a constant, or it keeps > oscillating around a certain value. > > As far as I am aware, I've followed the implementation as per textbook - I > can't see anything obviously wrong with it. Any idea what could be causing > this? > > > Do you have any test case that converges ? --Sergei. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 Yea, for A = 10x10, for example, here is the result, with relaxation parameter = 1.4: ...  x: 4.999998     iters: 710      norm: 0.0000010676089660  x: 4.999998     iters: 711      norm: 0.0000003683371855  x: 4.999998     iters: 712      norm: 0.0000006345346101  x: 4.999998     iters: 713      norm: 0.0000007962736616  x: 4.999998     iters: 714      norm: 0.0000008603756717  x: 4.999998     iters: 715      norm: 0.0000008603756717  x: 4.999998     iters: 716      norm: 0.0000008603756717  x: 4.999998     iters: 717      norm: 0.0000008603756717  x: 4.999998     iters: 718      norm: 0.0000008603756717  x: 4.999998     iters: 719      norm: 0.0000008603756717  x: 4.999998     iters: 720      norm: 0.0000008603756717 ans =     5.0000     9.0000    12.0000    14.0000    15.0000    15.0000    14.0000    12.0000     9.0000     5.0000 So it seems to work here?
Open this post in threaded view
|

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

 Here is my textbook's implementation of the Algorithm which I followed: The sparse matrix structure a contains val, col, row, d ie. 4 arrays containing non-zero values, column indices of those values, where each row starts, and the diagonal values. (I haven't done this part exactly as I figured it wouldn't make a difference since Octave stores the matrix I created anyway as a sparse matrix...right?...) for k := 1 to maxits do       normx := 0; dnormx := 0       for i := 1 to n do             sum : = 0          for j :=rowstart[i ] to rowstart[i + 1] − 1 do             sum := sum + val [ j ] × x [ col [ j ] ]          endfor j          delx := ω × (b[i ] − sum)/d[i ]          x [i ] := x [i ] + delx          normx := max(normx, abs( x [i ]) )          dnormx := max(dnormx, abs(delx ) )       endfor i       if dnormx ≤ tol + 4 × m × normx then return (x, k) endfor k
Open this post in threaded view
|

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

 In reply to this post by Joza On 4 Nov 2012, at 14:05, Joza wrote: > Yea, for A = 10x10, for example, here is the result, with relaxation > parameter = 1.4: > > ... > x: 4.999998     iters: 710      norm: 0.0000010676089660 > x: 4.999998     iters: 711      norm: 0.0000003683371855 > x: 4.999998     iters: 712      norm: 0.0000006345346101 > x: 4.999998     iters: 713      norm: 0.0000007962736616 > x: 4.999998     iters: 714      norm: 0.0000008603756717 > x: 4.999998     iters: 715      norm: 0.0000008603756717 > x: 4.999998     iters: 716      norm: 0.0000008603756717 > x: 4.999998     iters: 717      norm: 0.0000008603756717 > x: 4.999998     iters: 718      norm: 0.0000008603756717 > x: 4.999998     iters: 719      norm: 0.0000008603756717 > x: 4.999998     iters: 720      norm: 0.0000008603756717 > > ans = > >    5.0000 >    9.0000 >   12.0000 >   14.0000 >   15.0000 >   15.0000 >   14.0000 >   12.0000 >    9.0000 >    5.0000 > > So it seems to work here? > not really, it actually seems to be stagnating. I just tried a very simple implementation of the algorithm wiyh your parameters and it converges within a tolerance of 1e-15 in about 170 iterations. c.   _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 Would you mind posting the algorithm code? I'd like to compare it to mine and see exactly what is going on.
Open this post in threaded view
|

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

 ----- Original Message ----- > From: Joza <[hidden email]> > To: [hidden email] > Cc: > Sent: Sunday, November 4, 2012 8:29 PM > Subject: Re: Successive Over-Relaxation ... what is wrong? Improvements? > > Would you mind posting the algorithm code? > > I'd like to compare it to mine and see exactly what is going on. > > > > -- > View this message in context: > http://octave.1599824.n4.nabble.com/Successive-Over-Relaxation-what-is-wrong-Improvements-tp4646041p4646075.html> Sent from the Octave - General mailing list archive at Nabble.com. > _______________________________________________ > Help-octave mailing list > [hidden email] > https://mailman.cae.wisc.edu/listinfo/help-octave> The code you've posted as an example in written in Pascal-like language (if not in Pascal itself) and thus is trivially translatable into any other procedural language. So, show some diligence and perseverance. I.e. translate the code, take a simple case (4 x 4 matrix), do the calculations "manually" (say, using a calculator), compare the results to what's happening in your program. This all is called debugging, and is not Octave-specific. Regards,   Sergei. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave
Open this post in threaded view
|

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

 I understand that much very well, but the code I posted in the OP was written Octave. I have already explained that the code works fine for a 10 by 10 sparse matrix, but not for a 1000 by 1000 sparse matrix. Hence, I have been asking for help with my Octave code, from individuals here, whom, for the most part, are far more experience with the language than I am, and hopefully could shed some light on implementation problems with my code written in Octave.
Open this post in threaded view
|

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

 In reply to this post by Joza On Sat, Nov 3, 2012 at 6:12 PM, Joza wrote: Yea.... I've tried a few values for the relaxation parameter and the matrix size ... it seems that either the norm converges to a constant, or it keeps oscillating around a certain value. As far as I am aware, I've followed the implementation as per textbook - I can't see anything obviously wrong with it. Any idea what could be causing this? -- View this message in context: http://octave.1599824.n4.nabble.com/Successive-Over-Relaxation-what-is-wrong-Improvements-tp4646041p4646059.html Sent from the Octave - General mailing list archive at Nabble.com. _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave It works for me if I use the formula for the overrelaxation factor  w = 2/(1 + pi/n)that is usually recommended-- Ed Meyer _______________________________________________ Help-octave mailing list [hidden email] https://mailman.cae.wisc.edu/listinfo/help-octave