Re: implement a SOR like routine efficiently

[Stefan van der Walt]

Hi Ganesh

In Octave it's pretty simple to code this in C++.  Something like

#include <octave/oct.h>

DEFUN_DLD(SOR, args, nargin, "Generate SOR matrix\nB = SOR(A)") {

  octave_value_list retval;

  if (args.length() != 1) {
    print_usage("SOR");
    return retval;
  }

  Matrix A = args(0).matrix_value();
  if (error_state) {
    error("SOR: A not a matrix");
    return retval;
  }

  if (!A.is_square()) {
    error("SOR: A must be square");
    return retval;
  }

  int N = A.rows();

  for (int i = 2; i < N-1; i++) {
    for (int j = 2; j < N-1; j++) {
      A(i, j) = A(i-1,j) + A(i+1,j) + A(i, j+1) + A(i, j-1);
    }
  }

  retval.append(A);
  return retval;

}

Cheers
Stefan

On Mon, Nov 15, 2004 at 04:34:21PM -0600, Ganesh Bikshandi wrote:
> Hi,
> 
> I need to implement the following loop efficiently:
> 
> for i = 2:n-1
>   for j = 2:n-1 
>       a(i,j) = a(i-1,j) + a(i+1,j) + a(i,j+1) + a(i,j-1)
>    end
> end
> 
> The loop is not vectorizable, as it has a carried dependence. If I run
> this loop as it is, it is very slow ( in my program n is  256 and
> there are several iterations of the above two loops).  I am wondering
> if there is a library function in MATLAB or a better way to do this
> efficiently. Sorry, I know this is an octave user group, but I
> couldn't get an answer in MATLAB user group.
> 
> Ganesh
> 
> 
> 
> -------------------------------------------------------------
> Octave is freely available under the terms of the GNU GPL.
> 
> Octave's home on the web:  http://www.octave.org
> How to fund new projects:  http://www.octave.org/funding.html
> Subscription information:  http://www.octave.org/archive.html
> -------------------------------------------------------------
> 



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------