dolfin team mailing list archive

["Garth N. Wells" <gnw20@xxxxxxxxx>]

On 12/01/11 10:21, Anders Logg wrote:
> I was browsing through PETScVector.cpp and found this very strange
> piece of code:
> 
> void PETScVector::set(const double* block, uint m, const uint* rows)
> {
>   assert(x);
>   int _m =  static_cast<int>(m);
>   const int* _rows = reinterpret_cast<const int*>(rows);
>   if (m == 0)
>   {
>     _rows = &_m;
>     double tmp = 0.0;
>     block = &tmp;
>   }
>   VecSetValues(*x, _m, _rows, block, INSERT_VALUES);
> }
> 
> What's the deal with the m == 0 special case?
> 

It can happen in parallel.

> I would think that it would be natural to have
> 
>   assert(m > 0);
> 
> in place of that code.
> 
> What happens now when m = 0 is that the first row of the vector is set
> to 0.0 (and the rows argument is ignored). Why should m = 0 be allowed
> and why should element 0 be set to 0.0?
> 

When m = 0, PETSc will do nothing, since m is the number of entries to set.

Anyway, the code can be removed because

    VecSetValues

is not collective.  For collective functions (that must be called by
each process), we need to deal with the m = 0 case, otherwise MPI
programs will hang.

Garth

> --
> Anders
> 
> _______________________________________________
> Mailing list: https://launchpad.net/~dolfin
> Post to     : dolfin@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~dolfin
> More help   : https://help.launchpad.net/ListHelp