<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Thu, Jul 9, 2015 at 11:02 PM, Marco Zocca <span dir="ltr"><<a href="mailto:zocca.marco@gmail.com" target="_blank">zocca.marco@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">Hello,<br>
   I'm looking for some non-optimized PETSc code; namely, I struggle a<br>
bit with generalizing the provided examples. On one hand inlined<br>
elementary operations help keep track of the FLOPs but make for<br>
hard-to read code.<br>
<br>
For instance, KSP example 2, (<br>
<a href="http://www.mcs.anl.gov/petsc/petsc-current/src/ksp/ksp/examples/tutorials/ex2.c.html" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/petsc-current/src/ksp/ksp/examples/tutorials/ex2.c.html</a><br>
), isn't there a better abstraction to defining the finite-difference<br>
grid Laplacian than the following?<br>
<br>
for (Ii=Istart; Ii<Iend; Ii++) {<br>
      v = -1.0; i = Ii/n; j = Ii - i*n;<br>
      if (i>0)   {J = Ii - n; MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);}<br>
      if (i<m-1) {J = Ii + n; MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);}<br>
      if (j>0)   {J = Ii - 1; MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);}<br>
      if (j<n-1) {J = Ii + 1; MatSetValues(A,1,&Ii,1,&J,&v,ADD_VALUES);}<br>
      v = 4.0; MatSetValues(A,1,&Ii,1,&Ii,&v,ADD_VALUES);<br>
    }<br>
<br>
I see this approach as very error prone, yet this is one of the<br>
simplest examples.<br></blockquote><div><br></div><div>You can use MatSetValuesStencil() as is done here:</div><div><br></div><div><a href="http://www.mcs.anl.gov/petsc/petsc-current/src/ksp/ksp/examples/tutorials/ex50.c.html">http://www.mcs.anl.gov/petsc/petsc-current/src/ksp/ksp/examples/tutorials/ex50.c.html</a><br></div><div><br></div><div>on line 162.</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">
Another problematic trend I see is hardcoding constants and<br>
assumptions on the functions, e.g. in the following (taken from<br>
<a href="http://www.mcs.anl.gov/petsc/petsc-current/src/snes/examples/tutorials/ex2.c.html" rel="noreferrer" target="_blank">http://www.mcs.anl.gov/petsc/petsc-current/src/snes/examples/tutorials/ex2.c.html</a><br>
, removed comments for space) we are interested in forming a matrix A<br>
from the values contained in a vector x; i.e. some 'map' operation<br>
that sweeps x and prepares rows of A at a time according to some fixed<br>
stencil, with separate treatment of boundary values .<br>
<br>
Isn't there an idiom to abstract out this functionality? It's one of<br>
the most common operations in numerical PDE codes:<br></blockquote><div><br></div><div>Actually, this is only common is toy example. Almost no real PDEs have a fixed stencil</div><div>because you have material coefficients, or nonlinearities, or constraints, etc.</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">
VecGetArrayRead(x,&xx);<br>
<br>
VecGetSize(x,&n);<br>
d    = (PetscReal)(n - 1); d = d*d;<br>
<br>
for (i=1; i<n-1; i++) {<br>
  j[0] = i - 1; j[1] = i; j[2] = i + 1;<br>
  A[0] = A[2] = d; A[1] = -2.0*d + 2.0*xx[i];<br>
  MatSetValues(jac,1,&i,3,j,A,INSERT_VALUES);<br>
}<br>
<br>
i = 0;   A[0] = 1.0;<br>
MatSetValues(jac,1,&i,1,&i,A,INSERT_VALUES);<br>
<br>
i = n-1; A[0] = 1.0;<br>
<br>
MatSetValues(jac,1,&i,1,&i,A,INSERT_VALUES);<br>
<br>
VecRestoreArrayRead(x,&xx);<br>
<br>
<br>
I think it would be extremely helpful if the developers of the<br>
examples spent a little time in adding theory annotations to the code,<br>
and I would be very thankful if anybody could point me to relevant<br>
tutorials and literature.<br>
<br>
Thank you in advance and kind regards,<br>
Marco Zocca<br>
</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>