<div dir="ltr">On Sat, Dec 29, 2012 at 1:40 AM, amlan barua <span dir="ltr"><<a href="mailto:abarua@iit.edu" target="_blank">abarua@iit.edu</a>></span> wrote:<br><div class="gmail_extra"><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div dir="ltr"><div><div>Hi Barry,<br></div>I wrote the following piece according to your suggestions. Currently it does nothing but creates a vector with 1 at 1th position, 2 at 2th and so on. But I made it serial, i.e. (n+1)th place is computed using the value of nth place. My question, did I do it correctly, i.e. is it safe or results may change depending on problem size? This is much faster than VecSetValues, I believe the communication is minimum here because I take the advantage of ghost points.<br>
</div>Amlan<br><br>PetscInitialize(&argc,&argv,(char *)0,help);<br> ierr = MPI_Comm_size(PETSC_COMM_WORLD, &size); CHKERRQ(ierr); <br> ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank); CHKERRQ(ierr); <br>
ierr = DMDACreate1d(PETSC_COMM_WORLD,DMDA_BOUNDARY_NONE,size*5,1,1,PETSC_NULL,&da); CHKERRQ(ierr);<br> ierr = DMCreateGlobalVector(da,&vec); CHKERRQ(ierr);<br> ierr = VecSet(vec,1.00);<br> ierr = DMCreateLocalVector(da,&local);<br>
ierr = DMDAGetLocalInfo(da,&info);<br> ierr = DMDAVecGetArray(da,vec,&arr);<br> ierr = DMDAVecGetArray(da,local,&array);<br> temp = 1;<br> for (j=0;j<size;j++) {<br></div></blockquote><div>
<br></div><div style>This is needlessly sequential (or should be).</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"> ierr = DMGlobalToLocalBegin(da,vec,INSERT_VALUES,local); CHKERRQ(ierr);<br>
ierr = DMGlobalToLocalEnd(da,vec,INSERT_VALUES,local); CHKERRQ(ierr);<br></div></blockquote><div style><br>You should never use a communication routine while you have access to the array (*VecGetArray()).<br> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"> if (rank==j) {<br> for (i=info.xs;i<info.xs+info.xm;i++) {<br> if ((!i)==0) {<br>
array[i] = array[i] + array[i-1];<br></div></blockquote><div><br></div><div style>What sort of recurrence do you actually want to implement. When possible, it's much better to reorganize so that you can do local work followed by an MPI_Scan followed by more local work. MPI_Scan is fast (logarithmic).</div>
<div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr">
arr[i] = array[i];<br> }<br> }<br> }<br> }<br> ierr = DMDAVecRestoreArray(da,local,&array);<br> ierr = DMDAVecRestoreArray(da,vec,&arr);<br> ierr = VecView(vec,PETSC_VIEWER_STDOUT_WORLD);<br>
PetscFinalize();<br> return 0;<br><br>
</div><div class="HOEnZb"><div class="h5"><div class="gmail_extra"><br><br><div class="gmail_quote">On Thu, Dec 27, 2012 at 10:40 AM, Barry Smith <span dir="ltr"><<a href="mailto:bsmith@mcs.anl.gov" target="_blank">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><br>
On Dec 27, 2012, at 10:34 AM, amlan barua <<a href="mailto:abarua@iit.edu" target="_blank">abarua@iit.edu</a>> wrote:<br>
<br>
> I think I can use VecSetValues, is that right?<br>
<br>
</div> Yes you could do that. But since you are using a DMDA you could also use DMGetLocalVector(), DMGlobalToLocalBegin/End() followed by DMDAVecGetArray() to access the ghost values.<br>
<span><font color="#888888"><br>
Barry<br>
</font></span><div><div><br>
> Amlan<br>
><br>
><br>
> On Thu, Dec 27, 2012 at 9:04 AM, amlan barua <<a href="mailto:abarua@iit.edu" target="_blank">abarua@iit.edu</a>> wrote:<br>
> Hi Barry,<br>
> Is this scattering a very costly operation? I have to compute x[i] = f(x[i-1]) where f is known. Since this operation is strictly sequential, I thought of gathering the entire vector on processor 0, do the sequential operation there and scatter the result back. However this is unnecessary because I only need the bordering x[i] values. What can be a better way?<br>
> Amlan<br>
><br>
><br>
> On Thu, Dec 27, 2012 at 8:18 AM, Barry Smith <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>> wrote:<br>
><br>
> ierr = DMDACreateNaturalVector(da,&natural);CHKERRQ(ierr);<br>
> ierr = DMDAGlobalToNaturalBegin(da,xin,INSERT_VALUES,natural);CHKERRQ(ierr);<br>
> ierr = DMDAGlobalToNaturalEnd(da,xin,INSERT_VALUES,natural);CHKERRQ(ierr);<br>
><br>
> Now do VecScatterCreateToZero() from natural and the vector will be in the natural ordering on process zero with the dof interlaced.<br>
><br>
><br>
> Barry<br>
><br>
> On Dec 27, 2012, at 12:22 AM, amlan barua <<a href="mailto:abarua@iit.edu" target="_blank">abarua@iit.edu</a>> wrote:<br>
><br>
> > Hi,<br>
> > Is there an analogue of VecScatterCreateToZero for DA vectors? The DMDA object has more than one degrees of freedom.<br>
> > If there isn't any, should I use an IS object to do the scattering?<br>
> > Amlan<br>
><br>
><br>
><br>
<br>
</div></div></blockquote></div><br></div>
</div></div></blockquote></div><br></div></div>