<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Sat, Feb 28, 2015 at 1:56 PM, Gideon Simpson <span dir="ltr"><<a href="mailto:gideon.simpson@gmail.com" target="_blank">gideon.simpson@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word">Wait, then why is my code working when I use the DOF routines?  The number of degrees of freedom per DA vertex is constant across the DA, but that value is not set until run time.  In other words, what’s wrong with the code:<div><br></div><div><div>PetscOptionsGetInt(NULL,"-n_samples",&n_samples,NULL);</div><div>PetscOptionsGetInt(NULL,"-n_points",&n_points,NULL);</div></div><span class=""><div><br></div><div>DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE, n_samples, n_points, 0, NULL, &da);</div></span></div></blockquote><div><br></div><div>Nothing is wrong with that. I meant that you must have the smae  number on every vertex.</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"><div style="word-wrap:break-word"><div><div><div><span class="HOEnZb"><font color="#888888"><div>
<span style="border-collapse:separate;border-spacing:0px">-gideon</span>

</div></font></span><div><div class="h5">
<br><div><blockquote type="cite"><div>On Feb 28, 2015, at 2:54 PM, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>> wrote:</div><br><div><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Sat, Feb 28, 2015 at 1:50 PM, Gideon Simpson <span dir="ltr"><<a href="mailto:gideon.simpson@gmail.com" target="_blank">gideon.simpson@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div style="word-wrap:break-word">Supposing that I do not have a priori information as to the number of degrees of freedom per DA vertex; I want the number of mesh points along each sample path to be variable.  Hence, I can’t really use a statically defined structure as you suggest.  In that case, are the DOF routines the only option?</div></blockquote><div><br></div><div>DMDA just plain does not support that. It is close, but everything is not hooked up right. If you have</div><div>variable numbers of unknowns per site, you can</div><div><br></div><div>  1) Fix DMDA with our help (demands programming)</div><div><br></div><div>  2) Use DMPlex (demands reading and maybe looking at code)</div><div><br></div><div>Thanks,</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"><div style="word-wrap:break-word"><div><div>
<span style="border-collapse:separate;font-family:Helvetica;font-style:normal;font-variant:normal;font-weight:normal;letter-spacing:normal;line-height:normal;text-align:-webkit-auto;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px">-gideon</span>

</div>
<br><div><blockquote type="cite"><div>On Feb 28, 2015, at 2:42 PM, Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>> wrote:</div><br><div><div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Sat, Feb 28, 2015 at 1:35 PM, Gideon Simpson <span dir="ltr"><<a href="mailto:gideon.simpson@gmail.com" target="_blank">gideon.simpson@gmail.com</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">I’m having some trouble understanding what the difference between these two routines are, though I am finding that there certainly is a difference.  I have the following monte carlo problem.  I am generating n_sample paths each of length n_points, and storing them in a 1D DA:<br>
<br>
DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE, n_samples, n_points, 0, NULL, &da);<br>
DMCreateGlobalVector(da,&paths_vec);<br>
<br>
When I then go to access them,<br>
<br>
PetscScalar **u_array;<br>
<br>
I find that:<br>
<br>
DMDAVecGetArrayDOF(da, paths_vec, &u_array);<br>
<br>
works as exepected, in that u_array[i] is a pointer to the first index of the i-th sample path, but if I call:<br>
<br>
DMDAVecGetArray(da, paths_vec, &u_array);<br>
<br>
u_array[i] is something else, and my attempts to manipulate it result in segmentation faults, even though the code compiles and builds.</blockquote><div><br></div><div>Suppose that you have 4 PetscScalar values at each vertex of the 1D DMDA. If you use </div><div><br></div><div>  PetscScalar **u;</div><div><br></div><div>  DMDAVecGetArrayDOF(da, uVec, &u);</div><div><br></div><div>  u[i][2] /* refers to 3rd scalar on vertex i */</div><div><br></div><div>On the other hand you could use</div><div><br></div><div>typedef struct {</div><div>  PetscScalar a, b, c, d;</div><div>} Vals;</div><div><br></div><div> Vals *u;</div><div><br></div><div><div> DMDAVecGetArray(da, uVec, &u);</div></div><div><br></div><div> u[i].c /* refers to the same value as above */</div><div><br></div><div>Basically the DOF version gives you an extra level of indirection for the components.</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</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"><span><font color="#888888"><br>
-gideon<br>
<br>
</font></span></blockquote></div><br><br clear="all"><span><font color="#888888"><div><br></div>-- <br><div>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>
</font></span></div></div>
</div></blockquote></div><br></div></div></blockquote></div><br><br clear="all"><div><br></div>-- <br><div>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>
</div></blockquote></div><br></div></div></div></div></div></div></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>