<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Tue, Mar 1, 2016 at 9:50 AM, Jed Brown <span dir="ltr"><<a href="mailto:jed@jedbrown.org" target="_blank">jed@jedbrown.org</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><span class="">Matthew Knepley <<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>> writes:<br>
> I am trying to understand what you are saying here. I have a vector on 2<br>
> procs, divided as<br>
><br>
> p0  |----interior---|--ghost--|<br>
><br>
> p1  |----interior---|--ghost--|<br>
<br>
</span>Remember that the global Vec consists of just the interiors.  The local<br>
form is a sequential Vec that includes ghosts.<br>
<span class=""><br>
> Now I calculate some values in the interior on both procs using<br>
> VecGetArray(),<br>
><br>
> p0  |xxxxinteriorxxx|--ghost--|<br>
><br>
> p1  |xxxxinteriorxxx|--ghost--|<br>
<br>
</span>If you do it this way, then the global Vec is complete without needing<br>
to communicate anything.<br>
<span class=""><br>
> What is wrong with calling VecGhostUpdate() here to fill in the<br>
> ghosts?<br>
<br>
</span>The ghost values aren't part of the global Vec, so it's just needless<br>
communication.  But more insidiously, it encourages an interface where<br>
the input Vec is assumed to have updated ghost values (a "hidden"<br>
property).<br>
<span class=""><br>
> Is your contention that I should have also calculated the ghosts?<br>
<br>
</span>No.  But if you have a FE discretization then you're presumably using<br>
ADD_VALUES.  With VecGhost, the natural thing is to put your entries in<br>
the local form and use VecGhostUpdate(x,ADD_VALUES,SCATTER_REVERSE).<br></blockquote><div><br></div><div>Okay, this shows that we always need to specify exactly which vector, local or global,</div><div>we are talking about (and also why I don't use VecGhost).</div><div><br></div><div>I was confusing the VecGhostUpdate() above which puts local values in the global vector</div><div>with the update which would refresh the local portion. The former update should be called</div><div>automatically during assembly, which already happens in the FormFunctionLocals, and</div><div>would make sense to have in VecAssembly() for VecGhost I think. However, my data</div><div>arrows, and those of the questioner, were loading data into the ghost portion, which does</div><div>not make sense for VecAssembly().</div><div><br></div><div>So I think, on balance, we should not do this because the relevant assembly is done by</div><div>a better, more transparent mechanism already.</div><div><br></div><div>  Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Note that you would need a subsequent<br>
VecGhostUpdate(x,INSERT_VALUES,SCATTER_FORWARD) if you wanted the local<br>
form to be consistent with the now-assembled global form.  But the<br>
algebraic solvers don't need that, so it's only relevant if user code<br>
will be called on that Vec prior to modification of any type by the<br>
algebraic solver, and having user code assume the ghost values are<br>
current is a dangerous assumption that will make your code fragile and<br>
disgusting.<br>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature">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</div>
</div></div>