[petsc-users] MatSOR and GaussSeidel

[Cyrill Vonplanta]

I am finished with what I set up to do and I just want to leave this note to other potential PETSc-newbies that browse the mailing list:

The important thing to point out here, is that PETSc (with the configuration below) in general does NOT do a Gauss-Seidel step. Instead it might do a block (!) Gauss-Seidel step using the inodes of the matrix. This leads to iterates that look different from what you would expect if the SOR-step is done coordinate-wise and you get a different convergence history.

(It’s mentioned on the documentation page, but one quickly overreads it)

Cyrill

> On 23 Aug 2016, at 16:54, Jed Brown <jed at jedbrown.org> wrote:
> 
> Cyrill Vonplanta <cyrill.von.planta at usi.ch> writes:
> 
>> Dear PETSc-Users,
>> 
>> I am debugging a smoother of ours and i was wondering what settings of MatSOR exactly form one ordinary Gauss Seidel smoothing step. Currently I use:
>> 
>> ierr = MatSOR(A,b,1.0,(MatSORType)(SOR_ZERO_INITIAL_GUESS | SOR_FORWARD_SWEEP),0,1,1,x); CHKERRV(ierr);
> 
> Yes, this is a standard forward sweep of GS.  Note that your code below
> computes half zeros (because the vector starts as 0), but that it
> handles the diagonal incorrectly if you were to use a nonzero initial
> guess.
> 
>> I expect this to be the same as this naïve Gauss-Seidel step:
>> 
>> 
>> for (int i=0;i<m;i++){
>> 
>>      sum_i = 0;
>> 
>>      sum_i += ps_b_values[i];
>> 
>>      for (int j=0;j<m;j++){
>> 
>>        sum_i -= ps_A_values[i+j*m]*ps_x_values[j];
>> 
>>      }
>> 
>>      ps_x_values[i] += sum_i/ps_A_values[i*m +i];
>> 
>> }
>> 
>> The ps_* refer to the data parts of PETSc types (everything is serial and dense in my toy example. Initial x is zero.m is dimension of A). However the convergence history looks different. Am I missing something here?
>> 
>> Best Cyrill