On Jun 17, 2010, at 4:07 PM, Matthew Knepley wrote:
> Pushed. ML never ever worked right. Ugh.
Is it now converging better?
Barry
>
> Matt
>
> On Wed, Jun 16, 2010 at 7:24 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
> On Jun 16, 2010, at 9:38 AM, Matthew Knepley wrote:
>
>> Barry,
>>
>> You inserted a check in VecAXPY:563 (which Jed corrected)
>>
>> jed at 16195
>> 563
>> if (x == y) SETERRQ(((PetscObject)x)->comm,PETSC_ERR_ARG_IDN,"x and y
>> cannot be the same vector");
>>
>> However, this seems to make no sense in light of mg.c:53
>>
>> ierr =
>> MatInterpolateAdd(mglevels->interpolate,mgc->x,mglevels->x,mglevels->x);CHKERRQ(ierr);
>>
>> which will eventually call this with y == w, and has stack
>>
>> 0]PETSC ERROR: VecAXPY() line 563 in src/vec/vec/interface/rvector.c
>> [0]PETSC ERROR: MatMultAdd_ML() line 179 in src/ksp/pc/impls/ml/ml.c
>
> Is this code not completely wrong?
>
> ML_Operator_Apply(shell->mlmat,x_length,xarray,y_length,yarray);
> ierr = VecRestoreArray(x,&xarray);CHKERRQ(ierr);
> ierr = VecRestoreArray(y,&yarray);CHKERRQ(ierr);
> ierr = VecAXPY(y,1.0,w);CHKERRQ(ierr);
>
> ML_Operator() overwrites y. Now if w is the same as y that means it
> overwrites w now the VecAXPY(y,1.0,w) does not have the correct w (because y
> values are in w) so it ends up being 2*y which is not waht we want since it
> is suppose to add the original input values of w in.
>
> I think it needs to have two cases:
>
> if (w != y) {
> ML_Operator_Apply(shell->mlmat,x_length,xarray,y_length,yarray);
> ierr = VecRestoreArray(x,&xarray);CHKERRQ(ierr);
> ierr = VecRestoreArray(y,&yarray);CHKERRQ(ierr);
> ierr = VecAXPY(y,1.0,w);CHKERRQ(ierr);
> } else {
> copy w into a work vector then do
> ierr = VecRestoreArray(x,&xarray);CHKERRQ(ierr);
> ierr = VecRestoreArray(y,&yarray);CHKERRQ(ierr);
> ierr = VecAXPY(y,1.0,wworkvector);CHKERRQ(ierr);
> }
>
> Or do I totally misunderstand the code? Since it appears you are using this
> routine can you fix it?
>
>>
>>
>> Can I remove this check?
>
> No, I really like requiring these various arguments to be different with
> the vector routines. It is one consistent model that keeps the code simple
> and easy to understand. And we don't really need to support having it be the
> same vector. If the only place in ALL of PETSc where it is the same vector is
> actually in a bug I think that demonstrates there is no use case for the
> vectors to be the same.
>
> As someone else pointed out the BLAS call cannot have the two arguments be
> the same (because it is Fortran) so we've have to handle that case with
> additional ugly checks if we did support the same vector.
>
> Barry
>
>>
>> Thanks,
>>
>> Matt
>>
>> --
>> What most experimenters take for granted before they begin their experiments
>> is infinitely more interesting than any results to which their experiments
>> lead.
>> -- Norbert Wiener
>
>
>
>
> --
> What most experimenters take for granted before they begin their experiments
> is infinitely more interesting than any results to which their experiments
> lead.
> -- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL:
<http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20100617/17938fd0/attachment.html>
