<div dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">On 1 May 2014 18:25, Anton Popov <span dir="ltr"><<a href="mailto:popov@uni-mainz.de" target="_blank">popov@uni-mainz.de</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div text="#000000" bgcolor="#FFFFFF"><div class="">
<div>On 5/1/14 10:39 PM, Anush Krishnan
wrote:<br>
</div>
<blockquote type="cite">
<div dir="ltr">Hi Anton,<br>
<div>
<div class="gmail_extra"><br>
<div class="gmail_quote">On 29 April 2014 05:14, Anton Popov
<span dir="ltr"><<a href="mailto:popov@uni-mainz.de" target="_blank">popov@uni-mainz.de</a>></span>
wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div text="#000000" bgcolor="#FFFFFF">
<div>
<div><br>
</div>
</div>
You can do the whole thing much easier (to my
opinion).<br>
Since you created two DMDA anyway, just do:<br>
<br>
- find first index on every processor using MPI_Scan<br>
- create two global vectors (no ghosts)<br>
- put proper global indicies to global vectors<br>
- create two local vectors (with ghosts) and set ALL
entries to -1 (to have what you need in boundary
ghosts)<br>
- call global-to-local scatter<br>
<br>
Done!<br>
</div>
</blockquote>
<div><br>
</div>
<div>Won't the vectors contain floating point values? Are
you storing your indices as real numbers?<br>
</div>
</div>
</div>
</div>
</div>
</blockquote></div>
YES, exactly. And then I cast them to PetscInt when I compose
stencils.<br>
<br>
Something like this:<br>
idx[0] = (PetscInt) ivx[k][j][i];<br>
idx[1] = (PetscInt) ivx[k][j][i+1];<br>
idx[2] = (PetscInt) ivy[k][j][i];<br>
... and so on, where ivx, ivy, ... are the index arrays in x, y
.. directions<br>
<br>
Then I insert (actually add) stencils using MatSetValues.<br>
<br>
By the way, you can ideally preallocate in parallel with
MatMPIAIJSetPreallocation. To count precisely number entries in the
diagonal & off-diagonal blocks use the same mechanism to easily
access global indices, and then compare them with the local row
range, which is also known:<br>
- within the range -> d_nnz[i]++;<br>
- outside the range -> o_nnz[i]++;<br></div></blockquote><div><br>Thanks a lot for the help. I did exactly that and it worked perfectly.<br><br>But
just to clarify: if I was using 32-bit floats, would I start having
trouble when my matrix size reaches ~10 million due to the floating
point precision?<br> <br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div text="#000000" bgcolor="#FFFFFF">
<br>
Anton<div class=""><br>
<blockquote type="cite">
<div dir="ltr">
<div>
<div class="gmail_extra">
<div class="gmail_quote">
<div> </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div text="#000000" bgcolor="#FFFFFF"> <br>
The advantage is that you can access global indices
(including ghosts) in every block using i-j-k indexing
scheme.<br>
I personally find this way quite easy to implement
with PETSc<span><font color="#888888"><br>
<br>
Anton<br>
</font></span></div>
</blockquote>
</div>
<br>
</div>
<div class="gmail_extra">Thank you,<br>
Anush<br>
</div>
</div>
</div>
</blockquote>
<br>
</div></div>
</blockquote></div><br></div></div>