<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Thu, May 29, 2014 at 10:03 AM, Giacomo Mulas <span dir="ltr"><<a href="mailto:gmulas@oa-cagliari.inaf.it" target="_blank">gmulas@oa-cagliari.inaf.it</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">Hi Jose, and list.<br>
<br>
I am in the process of writing the function to use with<br>
EPSSetArbitrarySelection().<br>
<br>
Inside it, I will need to take some given component (which one is included<br>
in the info passed via ctx) of the eigenvector and square it. To do this,<br>
since the eigenvector is not necessarily local, I will need to first do a<br>
scatter to a local 1-component vector. So this would be like:<br>
<br>
... some omitted machinery to cast the info from *ctx to more easily<br>
accessible form...<br></blockquote><div><br></div><div>There might be an easier way to do this:</div><div> PetscScalar val = 0.0, gval;</div><div><br></div><div> <span style="color:rgb(0,0,0)">VecGetOwnershipRange(xr, &low, &high);</span></div>
<div><span style="color:rgb(0,0,0)"> if ((myindex >= low) && (myindex < high)) {</span></div> VecGetArray(localx1,&a);<br> val = a[myindex-low];<br> VecRestoreArray(localx1, &a);<br><div><span style="color:rgb(0,0,0)"> }</span></div>
<div><span style="color:rgb(0,0,0)"> MPI_Allreduce(&val, &gval, 1, MPIU_SCALAR, MPI_SUM, PETSC_COMM_WORLD);</span></div><div><span style="color:rgb(0,0,0)"><br></span></div><div><span style="color:rgb(0,0,0)">Now everyone has the value at myindex.</span></div>
<div><span style="color:rgb(0,0,0)"><br></span></div><div><span style="color:rgb(0,0,0)"> Matt</span></div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
ierr = ISCreateStride(PETXC_COMM_<u></u>WORLD,1,myindex,1,&is1_from);<u></u>CHKERRQ(ierr);<br>
ierr = ISCreateStride(PETSC_COMM_<u></u>WORLD,1,0,1,&is1_to);CHKERRQ(<u></u>ierr);<br>
ierr = VecCreateSeq(PETSC_COMM_SELF, 1, &localx1);CHKERRQ(ierr);<br>
ierr = VecScatterCreate(xr,is1_from,<u></u>localx1,is1_to,&scatter1); CHKERRQ(ierr);<br>
ierr = VecScatterBegin(scatter1,xr,<u></u>localx1,INSERT_VALUES,<br>
SCATTER_FORWARD);<br>
ierr = VecScatterEnd(scatter1,xr,<u></u>localx1,INSERT_VALUES,<br>
SCATTER_FORWARD);<br>
ierr = VecGetArray(localx1,&comp);<br>
*rr = comp*comp;<br>
ierr = VecRestoreArray(localx1, &comp);<br>
ierr = VecDestroy(localx1);<br>
ierr = VecScatterDestroy(&scatter1);<br>
ierr = ISDestroy(&is1_from);<br>
ierr = ISDestroy(&is1_to);<br>
*ri = 0;<br>
<br>
... some internal housekeeping omitted<br>
<br>
return 0;<br>
<br>
The questions are:<br>
<br>
1) when the arbitrary function is called, is it called on all nodes<br>
simultaneously, so that collective functions can be expected to work<br>
properly, being called on all involved nodes at the same time? Should all<br>
processes compute the *rr and *ri to be returned, and return the same value?<br>
would it be more efficient to create a unit vector uv containing only one<br>
nonzero component, and then use VecDot(xr, uv, &comp), instead of pulling<br>
the component I need and squaring it as I did above?<br>
<br>
<br>
2) since the stride, the 1-component vector, the scatter are presumably the same through all calls within one EPSSolve, can I take them out of the arbitrary function, and make them available to it through *ctx? For this<br>
to work, the structure of xr, the eigenvector passed to the arbitrary<br>
function, must be known outside of EPSSolve.<br>
<br>
<br>
Thanks, bye<br>
Giacomo<br>
<br>
-- <br>
______________________________<u></u>______________________________<u></u>_____<br>
<br>
Giacomo Mulas <<a href="mailto:gmulas@oa-cagliari.inaf.it" target="_blank">gmulas@oa-cagliari.inaf.it</a>><br>
______________________________<u></u>______________________________<u></u>_____<br>
<br>
INAF - Osservatorio Astronomico di Cagliari<br>
via della scienza 5 - 09047 Selargius (CA)<br>
<br>
tel. <a href="tel:%2B39%20070%2071180244" value="+3907071180244" target="_blank">+39 070 71180244</a><br>
mob. : <a href="tel:%2B39%20329%20%206603810" value="+393296603810" target="_blank">+39 329 6603810</a><br>
______________________________<u></u>______________________________<u></u>_____<br>
<br>
"When the storms are raging around you, stay right where you are"<br>
(Freddy Mercury)<br>
______________________________<u></u>______________________________<u></u>_____<br>
</blockquote></div><br><br clear="all"><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
</div></div>