<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;color:rgb(0,0,0);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 class="HOEnZb"><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 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>