Pushed. ML never ever worked right. Ugh.<div><br></div><div>   Matt<br><br><div class="gmail_quote">On Wed, Jun 16, 2010 at 7:24 PM, Barry Smith <span dir="ltr"><<a href="mailto:bsmith@mcs.anl.gov">bsmith@mcs.anl.gov</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;"><div style="word-wrap:break-word"><br><div><div class="im"><div>On Jun 16, 2010, at 9:38 AM, Matthew Knepley wrote:</div>
<br><blockquote type="cite">Barry,<div><br></div><div>  You inserted a check in VecAXPY:563 (which Jed corrected)</div><div><br></div><div><span style="font-family:sans-serif;font-size:12px"><table style="padding-top:8px;padding-right:4px;padding-bottom:8px;padding-left:4px">

<tbody><tr style="font-family:monospace"><td style="padding-top:2px;padding-right:5px;padding-bottom:2px;padding-left:5px;font-size:12px;vertical-align:top;color:rgb(153, 153, 153);text-decoration:none;text-align:right">

<a title="c5b7ef1f9e2d: Fix some misuse of SETERRQ, Barry was naughty and didn't compile" style="color:rgb(136, 0, 0)">jed@16195</a></td>
<td style="padding-top:2px;padding-right:5px;padding-bottom:2px;padding-left:5px;font-size:12px;vertical-align:top"><pre style="margin-top:0px;margin-right:0px;margin-bottom:0px;margin-left:0px"><a style="color:rgb(136, 0, 0);text-decoration:none">   563</a></pre>

</td><td style="padding-top:2px;padding-right:5px;padding-bottom:2px;padding-left:5px;font-size:12px;vertical-align:top"><pre style="margin-top:0px;margin-right:0px;margin-bottom:0px;margin-left:0px">  if (x == y) SETERRQ(((PetscObject)x)->comm,PETSC_ERR_ARG_IDN,"x and y cannot be the same vector");
</pre></td></tr><tr style="background-color:rgb(237, 236, 230);font-family:monospace"></tr></tbody></table><br></span></div><div><span style="font-family:sans-serif;font-size:12px">However, this seems to make no sense in light of mg.c:53</span></div>

<div><span style="font-family:sans-serif;font-size:12px"><br></span></div><div><span style="font-family:sans-serif;font-size:12px"><div>    ierr = MatInterpolateAdd(mglevels->interpolate,mgc->x,mglevels->x,mglevels->x);CHKERRQ(ierr);</div>

<div><br></div><div>which will eventually call this with y == w, and has stack</div></span><div><br></div><div>0]PETSC ERROR: VecAXPY() line 563 in src/vec/vec/interface/rvector.c</div><div>[0]PETSC ERROR: MatMultAdd_ML() line 179 in src/ksp/pc/impls/ml/ml.c</div>
</div></blockquote><div><br></div></div>   Is this code not completely wrong?</div><div><br></div><div><div>  ML_Operator_Apply(shell->mlmat,x_length,xarray,y_length,yarray); </div><div>  ierr = VecRestoreArray(x,&xarray);CHKERRQ(ierr); </div>
<div>  ierr = VecRestoreArray(y,&yarray);CHKERRQ(ierr); </div><div>  ierr = VecAXPY(y,1.0,w);CHKERRQ(ierr); </div><div><br></div><div>   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.</div>
<div><br></div><div>   I think it needs to have two cases: </div><div><br></div><div>   if (w != y) {</div><div><div>      ML_Operator_Apply(shell->mlmat,x_length,xarray,y_length,yarray); </div><div>      ierr = VecRestoreArray(x,&xarray);CHKERRQ(ierr); </div>
<div>     ierr = VecRestoreArray(y,&yarray);CHKERRQ(ierr); </div><div>      ierr = VecAXPY(y,1.0,w);CHKERRQ(ierr); </div><div>  } else {</div><div>        copy w into a work vector then do</div><div><div>      ierr = VecRestoreArray(x,&xarray);CHKERRQ(ierr); </div>
<div>     ierr = VecRestoreArray(y,&yarray);CHKERRQ(ierr); </div><div>      ierr = VecAXPY(y,1.0,wworkvector);CHKERRQ(ierr); </div><div>  }</div><div><br></div><div>  Or do I totally misunderstand the code? Since it appears you are using this routine can you fix it?</div>
</div><div><br></div></div><div class="im"><blockquote type="cite"><div><div><br></div><div><br></div><div>Can I remove this check?</div></div></blockquote><div><br></div></div>   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.  </div>
<div><br></div><div>   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.</div>
<div><br></div><font color="#888888"><div>   Barry</div></font><div class="im"><div><br><blockquote type="cite"><div><div>
<br></div><div>  Thanks,</div><div><br></div><div>      Matt</div><div><br></div>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>

-- Norbert Wiener<br>
</div>
</blockquote></div><br></div></div></blockquote></div><br><br clear="all"><br>-- <br>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>
-- Norbert Wiener<br>
</div>