<div class="gmail_quote"><div>Dear Matt,</div><div><br></div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">&gt; Currently I am planning to use Petsc Krylov solver for implicit time<br>

&gt; integration of my CFD code. I need some pointer on the implementation. I am<br>
&gt; providing a detailed list of things I am not very clear with in both the<br>
&gt; documentation and the mailing list archives.<br>
&gt;<br>
&gt; --&gt; The equations to be solved is of the form :  [ I / \delta t + J ] \delta<br>
&gt; U^n = - R[ U^n ]<br>
&gt;<br>
&gt; --&gt; The Jacobian J is hard to evaluate exactly and hence we use the finite<br>
&gt; difference method to evaluate J. I want to use a matrix-free approach so<br>
&gt; all I need is the action of this J on some vector v.<br>
&gt;<br>
&gt; --&gt; Jv = \frac{\partial R}{\partial U}v = ( R[ U^n + hv ] - R[ U^] ) / h,<br>
&gt; h is some parameter. Since I is identity matrix and \delta t is fixed this<br>
&gt; diagonal matrix is trivial to evaluate. Thus [ I / \delta t + J ] is<br>
&gt; available as a user define function in my solver.<br>
&gt;<br>
&gt; 1) Vector U^n is from an unstructured node distribution and is partitioned<br>
&gt; outside of Petsc. I have both the global indexing and local indexing<br>
&gt; information along with the ghost node information.<br>
&gt; 2) I have already defined functions to move data across the ghost nodes<br>
&gt; given the local variable address and memory stride/offset information.<br>
&gt; 3) In Jacobian-vector Jv evaluation, one has to exchange the ghost cell<br>
&gt; values of vector v with adjacent processor before the calculation is done.<br>
&gt; 4) I can write the Jacobian-vector Jv evaluation function to perform ghost<br>
&gt; node exchange before evaluation but I am not clear as to how I can<br>
&gt; interface my local arrays and local / global indexing information into<br>
&gt; Petsc.<br>
&gt; 5) I understand that I have to create a Matrix Shell for the matrix free<br>
&gt; method and define my user-defined function to evaluate the matrix-vector<br>
&gt; product given an input vector. But it is not clear as to how Petsc<br>
&gt; understands how I order my vectors locally and how that corresponds to the<br>
&gt; global indexes.<br>
&gt; 6) I cannot use the DA object as it is for structured/contiguous<br>
&gt; partitioning .i.e., the local indices don&#39;t correspond to contiguous global<br>
&gt; ordering/indexing.<br>
&gt; 7) I also don&#39;t want to use the unstructured mesh and decomposition<br>
&gt; routines in Petsc and I already have my own.<br>
&gt;<br>
&gt; Can this be done in Petsc ?<br>
&gt;<br>
<br>
Yes, there are ??? parts:<br>
<br>
  1) You want to manage unstructured data yourself. The easiest way to do<br>
this now is to implement a DM. You just need to provide<br>
      routines for making global and local vectors (you have this), and<br>
routines for mapping between them (you have this). You will have<br>
      to use petsc-dev, but this definitely looks like the easiest thing to<br>
do.<br>
<br></blockquote><div><br></div><div><br></div><div><div>Thanks for the reply. </div><div><br></div><div> </div></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

     a) If you consider using part of our infrastructure, you can build a<br>
PetscSection and a PetscSF, and everything else will be automatically<br>
         calculated. These are very new things, so there is not a lot of<br>
documentation, and thus you might prefer to do the whole DM yourself.<br>
<br></blockquote><div><br></div><div><div><div>Most of the tutorial slides point out that DM is for structured arrays (maybe this mislead me). Can you point me to the functions that i would have to use to create </div><div>
<br></div><div>1) DM array ghost cells </div><div>2) Local/Global indexing</div><div>3) DM array exchange</div></div></div><div> </div><div>It would be nice you can share some code snippet on writing our own DM&#39;s. </div>
<div><br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
  2) You want to do matrix-free action of your Jacobian. This is easy using<br>
MatShell, however the real question is How will you precondition that<br>
       equation?</blockquote><div><br></div><div>Currently I use the Jacobi diagonal preconditioner. I use the largest spectral radius of my Euler equations \lambda = U + a. BTW I actually have implemented my own parallel Krylov solver for my CFD code using the algorithm in Yousuf Saad&#39;s book. But my advisor urged me to use Petsc so that I would be able to use the various pre-conditioners in Petsc and the code would run more faster/optimal. I have not written my Krylov solver optimally.</div>
<div><br></div><div>So is there is no way I can use preconditioning in matrix-free context in Petsc ?</div><div><br></div><div>Regards,</div><div> </div></div>--<div>  <br>Pavanakumar Mohanamuraly<br>
</div>