On Wed, Mar 14, 2012 at 1:33 PM, Xiangdong Liang <span dir="ltr">&lt;<a href="mailto:xdliang@gmail.com">xdliang@gmail.com</a>&gt;</span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
On Tue, Mar 13, 2012 at 5:57 PM, Matthew Knepley &lt;<a href="mailto:knepley@gmail.com">knepley@gmail.com</a>&gt; wrote:<br>
&gt; On Tue, Mar 13, 2012 at 4:39 PM, Xiangdong Liang &lt;<a href="mailto:xdliang@gmail.com">xdliang@gmail.com</a>&gt; wrote:<br>
&gt;&gt;<br>
&gt;&gt; Hello everyone,<br>
&gt;&gt;<br>
&gt;&gt; Can someone provide me advice on using OpenMP in PETSc? I am solving a<br>
&gt;&gt; problem like this:<br>
&gt;&gt;<br>
&gt;&gt; int main()<br>
&gt;&gt; {<br>
&gt;&gt;  vec vsum;<br>
&gt;&gt;<br>
&gt;&gt;  for(i=0; i&lt; N; i++)<br>
&gt;&gt;    {<br>
&gt;&gt;      vec vi;<br>
&gt;&gt;      fcomputev(i, vi);<br>
&gt;&gt;      VecAXPY(vsum, 1.0, vi); // vsum +=vi;<br>
&gt;&gt;    }<br>
&gt;&gt;<br>
&gt;&gt; }<br>
&gt;&gt;<br>
&gt;&gt; Can I use OpenMP &quot;omp parallel for&quot; to do the loop in parallel? For<br>
&gt;&gt; example, suppose I have 8 processes.   It would be nice if each petsc<br>
&gt;&gt; subrountine fcomputev uses 2 processes while 4 different i&#39;s are<br>
&gt;&gt; computed in parallel (since different i&#39;s are independent).<br>
&gt;<br>
&gt;<br>
&gt; PETSc is not threadsafe. This is trivial to do in MPI where the comm for eah<br>
&gt; Vec is a group of 2 procs.<br>
<br>
Do you mean the use of MPI_Reduce? If I use MPI_Reduce, should I<br>
convert the Vec objects into regular array by vecgetarray first, then<br>
apply the MPI_Op (MPI_SUM) on these arrays? Is there any way to<br>
circumvent this vec-arrary converting by define MPI_Op and<br>
MPI_datatype on Vec ?<br></blockquote><div><br></div><div>I think you miss my point. This is the difference between writing SPMD programs</div><div>and threaded programs. PETSc is SPMD, and I don&#39;t ever expect that to change.</div>
<div>CUDA, for instance, is also SPMD.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
Another quick question, where can I find the implementation of vec<br>
ops?  For example, In petscvec.h, VecAXPY is implemented like<br>
(*y-&gt;ops-&gt;axpy)(y,alpha,x).  Can you point me to the implementation of<br>
methods ops-&gt;axpy?<br></blockquote><div><br></div><div>src/vec/vec/impls</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">

Thanks.<br>
Xiangdong<br>
<br>
<br>
<br>
&gt;<br>
&gt;   Matt<br>
&gt;<br>
&gt;&gt;<br>
&gt;&gt; Any helps or hints on this would be appreciated.<br>
&gt;&gt;<br>
&gt;&gt; Best,<br>
&gt;&gt; Xiangdong<br>
&gt;<br>
&gt;<br>
&gt;<br>
<span class="HOEnZb"><font color="#888888">&gt;<br>
&gt; --<br>
&gt; What most experimenters take for granted before they begin their experiments<br>
&gt; is infinitely more interesting than any results to which their experiments<br>
&gt; lead.<br>
&gt; -- Norbert Wiener<br>
</font></span></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<br>