<div dir="ltr">On Wed, Aug 7, 2013 at 7:26 AM, Bishesh Khanal <span dir="ltr"><<a href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>></span> wrote:<br><div class="gmail_extra"><div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">On Wed, Aug 7, 2013 at 2:15 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</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 dir="ltr"><div><div>On Wed, Aug 7, 2013 at 7:07 AM, Bishesh Khanal <span dir="ltr"><<a href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>></span> wrote:<br>

</div></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 dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">On Tue, Aug 6, 2013 at 11:34 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</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 dir="ltr"><div><div>On Tue, Aug 6, 2013 at 10:59 AM, Bishesh Khanal <span dir="ltr"><<a href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>></span> wrote:<br>



</div></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 dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">
On Tue, Aug 6, 2013 at 4:40 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</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 dir="ltr"><div><div>On Tue, Aug 6, 2013 at 8:06 AM, Bishesh Khanal <span dir="ltr"><<a href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>></span> wrote:<br>





</div></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 dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">

On Mon, Aug 5, 2013 at 4:14 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</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 dir="ltr"><div><div>On Mon, Aug 5, 2013 at 8:48 AM, Bishesh Khanal <span dir="ltr"><<a href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>></span> wrote:<br>







</div></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 dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">

On Mon, Aug 5, 2013 at 3:17 PM, Matthew Knepley <span dir="ltr"><<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</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 dir="ltr"><div><div>On Mon, Aug 5, 2013 at 7:54 AM, Bishesh Khanal <span dir="ltr"><<a href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>></span> wrote:<br>









</div></div><div class="gmail_extra"><div><div><div class="gmail_quote">
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><br><div class="gmail_extra"><br><br><div class="gmail_quote">

On Wed, Jul 17, 2013 at 9:48 PM, Jed Brown <span dir="ltr"><<a href="mailto:jedbrown@mcs.anl.gov" target="_blank">jedbrown@mcs.anl.gov</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>Bishesh Khanal <<a href="mailto:bisheshkh@gmail.com" target="_blank">bisheshkh@gmail.com</a>> writes:<br>







<br>
> Now, I implemented two different approaches, each for both 2D and 3D, in<br>
> MATLAB. It works for the smaller sizes but I have problems solving it for<br>
> the problem size I need (250^3 grid size).<br>
> I use staggered grid with p on cell centers, and components of v on cell<br>
> faces. Similar split up of K to cell center and faces to account for the<br>
> variable viscosity case)<br>
<br>
</div>Okay, you're using a staggered-grid finite difference discretization of<br>
variable-viscosity Stokes.  This is a common problem and I recommend<br>
starting with PCFieldSplit with Schur complement reduction (make that<br>
work first, then switch to block preconditioner).  You can use PCLSC or<br>
(probably better for you), assemble a preconditioning matrix containing<br>
the inverse viscosity in the pressure-pressure block.  This diagonal<br>
matrix is a spectrally equivalent (or nearly so, depending on<br>
discretization) approximation of the Schur complement.  The velocity<br>
block can be solved with algebraic multigrid.  Read the PCFieldSplit<br>
docs (follow papers as appropriate) and let us know if you get stuck.<br></blockquote><div><br></div><div>I was trying to assemble the inverse viscosity diagonal matrix to use as the preconditioner for the Schur complement solve step as you suggested. I've few questions about the ways to implement this in Petsc:<br>











</div><div>A naive approach that I can think of would be to create a vector with its components as reciprocal viscosities of the cell centers corresponding to the pressure variables, and then create a diagonal matrix from this vector. However I'm not sure about:<br>











How can I make this matrix, (say S_p) compatible to the Petsc distribution of the different rows of the main system matrix over different processors ? The main matrix was created using the DMDA structure with 4 dof as explained before. <br>











</div></div>The main matrix correspond to the DMDA with 4 dofs but for the S_p matrix would correspond to only pressure space. Should the distribution of the rows of S_p among different processor not correspond to the distribution of the rhs vector, say h' if it is solving for p with Sp = h' where S = A11 inv(A00) A01 ? <br>










</div></div></blockquote></div><br clear="all"></div></div><div>PETSc distributed vertices, not dofs, so it never breaks blocks. The P distribution is the same as the entire problem divided by 4.</div></div></div></blockquote>









<div><br></div><div>Thanks Matt. So if I create a new DMDA with same grid size but with dof=1 instead of 4, the vertices for this new DMDA will be identically distributed as for the original DMDA ? Or should I inform PETSc by calling a particular function to make these two DMDA have identical distribution of the vertices ?<br>








</div></div></div></div></blockquote><div><br></div></div></div><div>Yes.</div><div><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 dir="ltr">
<div class="gmail_extra">
<div class="gmail_quote"><div>Even then I think there might be a problem due to the presence of "fictitious pressure vertices". The system matrix (A) contains an identity corresponding to these fictitious pressure nodes, thus when using a -pc_fieldsplit_detect_saddle_point, will detect a A11 zero block of size that correspond to only non-fictitious P-nodes. So the preconditioner S_p for the Schur complement outer solve with Sp = h' will also need to correspond to only the non-fictitious P-nodes. This means its size does not directly correspond to the DMDA grid defined for the original problem. Could you please suggest an efficient way of assembling this S_p matrix ?</div>








</div></div></div></blockquote><div><br></div></div><div>Don't use detect_saddle, but split it by fields -pc_fieldsplit_0_fields 0,1,2 -pc_fieldsplit_1_fields 4</div></div></div></div></blockquote><div><br></div><div>







How can I set this split in the code itself without giving it as a command line option when the system matrix is assembled from the DMDA for the whole system with 4 dofs. (i.e. <i>without</i> using the DMComposite or <i>without</i> using the nested block matrices to assemble different blocks separately and then combine them together). <br>







I need the split to get access to the fieldsplit_1_ksp in my code, because not using detect_saddle_point means I cannot use -fieldsplit_1_ksp_constant_null_space due to the presence of identity for the fictitious pressure nodes present in the fieldsplit_1_ block. I need to use PCFieldSplitGetSubKsp() so that I can set proper null-space basis.<br>






</div></div></div></div></blockquote><div><br></div></div></div><div>This is currently a real problem with the DMDA. In the unstructured case, where we always need specialized spaces, you can</div><div>use something like</div>





<div><br>
</div><div><div>    PetscObject  pressure;</div><div>    MatNullSpace nullSpacePres;</div><div><br></div><div>    ierr = DMGetField(dm, 1, &pressure);CHKERRQ(ierr);</div><div>    ierr = MatNullSpaceCreate(PetscObjectComm(pressure), PETSC_TRUE, 0, NULL, &nullSpacePres);CHKERRQ(ierr);</div>






<div>    ierr = PetscObjectCompose(pressure, "nullspace", (PetscObject) nullSpacePres);CHKERRQ(ierr);</div><div>    ierr = MatNullSpaceDestroy(&nullSpacePres);CHKERRQ(ierr);</div></div><div><br></div><div>and then DMGetSubDM() uses this information to attach the null space to the IS that is created using the information in the PetscSection.</div>






<div>If you use a PetscSection to set the data layout over the DMDA, I think this works correctly, but this has not been tested at all and is very</div><div>new code. Eventually, I think we want all DMs to use this mechanism, but we are still working it out.</div>





</div></div></div></blockquote><div><br>Currently I do not use PetscSection. If this makes a cleaner approach, I'd try it too but may a bit later (right now I'd like test my model with a quickfix even if it means a little dirty code!)<br>





</div><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 dir="ltr"><div class="gmail_extra">
<div class="gmail_quote">
<div><br></div><div>Bottom line: For custom null spaces using the default layout in DMDA, you need to take apart the PCFIELDSPLIT after it has been setup,</div><div>which is somewhat subtle. You need to call KSPSetUp() and then reach in and get the PC, and the subKSPs. I don't like this at all, but we</div>






<div>have not reorganized that code (which could be very simple and inflexible since its very structured).</div></div></div></div></blockquote><div><br></div><div>So I tried to get this approach working but I could not succeed and encountered some errors. Here is a code snippet:<br>





<br>//mDa is the DMDA that describes the whole grid with all 4 dofs (3 velocity components and 1 pressure comp.)<br></div><div></div><div>ierr = DMKSPSetComputeRHS(mDa,computeRHSTaras3D,this);CHKERRQ(ierr);<br>    ierr = DMKSPSetComputeOperators(mDa,computeMatrixTaras3D,this);CHKERRQ(ierr);<br>





    ierr = KSPSetDM(mKsp,mDa);CHKERRQ(ierr);<br>    ierr = KSPSetNullSpace(mKsp,mNullSpaceSystem);CHKERRQ(ierr);   //I've the mNullSpaceSystem based on mDa, that contains a null space basis for the complete system.<br>





    ierr = KSPSetFromOptions(mKsp);CHKERRQ(ierr);                                //This I expect would register these options I give:-pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_0_fields 0,1,2 //-pc_fieldsplit_1_fields 3<br>





<br>    ierr = KSPSetUp(mKsp);CHKERRQ(ierr);<br><br>    ierr = KSPGetPC(mKsp,&mPcOuter);     //Now get the PC that was obtained from the options (fieldsplit)<br><br>    ierr = PCFieldSplitSchurPrecondition(mPcOuter,PC_FIELDSPLIT_SCHUR_PRE_USER,mPcForSc);CHKERRQ(ierr);   //I have created the matrix mPcForSc using a DMDA with identical //size to mDa but with dof=1 corresponding to the pressure nodes (say mDaPressure).<br>





  <br>  ierr = PCSetUp(mPcOuter);CHKERRQ(ierr);  <br><br>    KSP *kspSchur;<br>    PetscInt kspSchurPos = 1;<br>    ierr = PCFieldSplitGetSubKSP(mPcOuter,&kspSchurPos,&kspSchur);CHKERRQ(ierr);<br>    ierr = KSPSetNullSpace(kspSchur[1],mNullSpacePressure);CHKERRQ(ierr);                 //The null space is the one that correspond to only pressure nodes, created using the mDaPressure.<br>





    ierr = PetscFree(kspSchur);CHKERRQ(ierr);<br><br>    ierr = KSPSolve(mKsp,NULL,NULL);CHKERRQ(ierr);<br></div></div></div></div></blockquote><div><br></div></div></div><div>Sorry, you need to return to the old DMDA behavior, so you want</div>




<div><br></div><div>  -pc_fieldsplit_dm_splits 0</div></div></div></div></blockquote><div><br></div><div>Thanks, with this it seems I can attach the null space properly, but I have a question regarding whether the Schur complement ksp solver is actually using the preconditioner matrix I provide. <br>



When using -ksp_view, the outer level pc object of type fieldsplit does report that: "Preconditioner for the Schur complement formed from user provided matrix", but in the KSP solver for Schur complement S, the pc object (fieldsplit_1_) is of type ilu and doesn't say that it is using the matrix I provide. Am I missing something here ? <br>



Below are the relevant commented code snippet and the output of the -ksp_view<br></div><div>(The options I used: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_dm_splits 0 -pc_fieldsplit_0_fields 0,1,2 -pc_fieldsplit_1_fields 3 -ksp_converged_reason -ksp_view )<br>


</div></div></div></div></blockquote><div><br></div></div></div><div>If ILU does not error, it means it is using your matrix, because the Schur complement matrix cannot be factored, and FS says it is using your matrix.</div>

</div></div></div></blockquote><div><br></div><div>Thanks Matt! By the way, what do these statements mean in -ksp_view results:<br>not using I-node routines<br>or<br>using I-node routines: found 729 nodes, limit used is 5</div>
</div></div></div></blockquote><div><br></div><div>This is an optimization in the sparse matrix storage format for rows with identical nonzero structure.</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 dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div class="gmail_extra">
<div class="gmail_quote"><div><br>

</div><div>   Matt</div><div><div><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 dir="ltr"><div class="gmail_extra"><div class="gmail_quote">

<div>Code snippet:<br>
ierr = KSPSetNullSpace(mKsp,mNullSpaceSystem);CHKERRQ(ierr);   //The nullspace for the whole system<br>    ierr = KSPSetFromOptions(mKsp);CHKERRQ(ierr);                    <br>    ierr = KSPSetUp(mKsp);CHKERRQ(ierr);                   //Set up mKsp with the options provided with fieldsplit and the fields associated with the two splits.<br>



<br>    ierr = KSPGetPC(mKsp,&mPcOuter);CHKERRQ(ierr);                //Get the fieldsplit pc set up from the options<br><br>    ierr = PCFieldSplitSchurPrecondition(mPcOuter,PC_FIELDSPLIT_SCHUR_PRE_USER,mPcForSc);CHKERRQ(ierr);  //Use mPcForSc as the preconditioner for Schur Complement<br>



<br>    KSP *kspSchur;<br>    PetscInt kspSchurPos = 1;<br>    ierr = PCFieldSplitGetSubKSP(mPcOuter,&kspSchurPos,&kspSchur);CHKERRQ(ierr);<br>    ierr = KSPSetNullSpace(kspSchur[1],mNullSpacePressure);CHKERRQ(ierr);                                  //Attach the null-space for the Schur complement ksp solver.<br>



    ierr = PetscFree(kspSchur);CHKERRQ(ierr);<br><br>    ierr = KSPSolve(mKsp,NULL,NULL);CHKERRQ(ierr);<br><br><br></div><div><br></div><div>the output of the -ksp_view<br>KSP Object: 1 MPI processes<br>  type: gmres<br>


    GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement<br>
    GMRES: happy breakdown tolerance 1e-30<br>  maximum iterations=10000, initial guess is zero<br>  tolerances:  relative=1e-05, absolute=1e-50, divergence=10000<br>  left preconditioning<br>  has attached null space<br>



  using PRECONDITIONED norm type for convergence test<br>PC Object: 1 MPI processes<br>  type: fieldsplit<br>    FieldSplit with Schur preconditioner, blocksize = 4, factorization FULL<br>    Preconditioner for the Schur complement formed from user provided matrix<br>



    Split info:<br>    Split number 0 Fields  0, 1, 2<br>    Split number 1 Fields  3<br>    KSP solver for A00 block <br>      KSP Object:      (fieldsplit_0_)       1 MPI processes<br>        type: gmres<br>          GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement<br>



          GMRES: happy breakdown tolerance 1e-30<br>        maximum iterations=10000, initial guess is zero<br>        tolerances:  relative=1e-05, absolute=1e-50, divergence=10000<br>        left preconditioning<br>        using PRECONDITIONED norm type for convergence test<br>



      PC Object:      (fieldsplit_0_)       1 MPI processes<br>        type: ilu<br>          ILU: out-of-place factorization<br>          0 levels of fill<br>          tolerance for zero pivot 2.22045e-14<br>          using diagonal shift on blocks to prevent zero pivot<br>



          matrix ordering: natural<br>          factor fill ratio given 1, needed 1<br>            Factored matrix follows:<br>              Matrix Object:               1 MPI processes<br>                type: seqaij<br>



                rows=2187, cols=2187<br>                package used to perform factorization: petsc<br>                total: nonzeros=140625, allocated nonzeros=140625<br>                total number of mallocs used during MatSetValues calls =0<br>



                  using I-node routines: found 729 nodes, limit used is 5<br>        linear system matrix = precond matrix:<br>        Matrix Object:         1 MPI processes<br>          type: seqaij<br>          rows=2187, cols=2187<br>



          total: nonzeros=140625, allocated nonzeros=140625<br>          total number of mallocs used during MatSetValues calls =0<br>            using I-node routines: found 729 nodes, limit used is 5<br>    KSP solver for S = A11 - A10 inv(A00) A01 <br>



      KSP Object:      (fieldsplit_1_)       1 MPI processes<br>        type: gmres<br>          GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement<br>          GMRES: happy breakdown tolerance 1e-30<br>



        maximum iterations=10000, initial guess is zero<br>        tolerances:  relative=1e-05, absolute=1e-50, divergence=10000<br>        left preconditioning<br>        has attached null space<br>        using PRECONDITIONED norm type for convergence test<br>



      PC Object:      (fieldsplit_1_)       1 MPI processes<br>        type: ilu<br>          ILU: out-of-place factorization<br>          0 levels of fill<br>          tolerance for zero pivot 2.22045e-14<br>          using diagonal shift on blocks to prevent zero pivot<br>



          matrix ordering: natural<br>          factor fill ratio given 1, needed 1<br>            Factored matrix follows:<br>              Matrix Object:               1 MPI processes<br>                type: seqaij<br>



                rows=729, cols=729<br>                package used to perform factorization: petsc<br>                total: nonzeros=15625, allocated nonzeros=15625<br>                total number of mallocs used during MatSetValues calls =0<br>



                  not using I-node routines<br>        linear system matrix followed by preconditioner matrix:<br>        Matrix Object:         1 MPI processes<br>          type: schurcomplement<br>          rows=729, cols=729<br>



            Schur complement A11 - A10 inv(A00) A01<br>            A11<br>              Matrix Object:               1 MPI processes<br>                type: seqaij<br>                rows=729, cols=729<br>                total: nonzeros=15625, allocated nonzeros=15625<br>



                total number of mallocs used during MatSetValues calls =0<br>                  not using I-node routines<br>            A10<br>              Matrix Object:               1 MPI processes<br>                type: seqaij<br>



                rows=729, cols=2187<br>                total: nonzeros=46875, allocated nonzeros=46875<br>                total number of mallocs used during MatSetValues calls =0<br>                  not using I-node routines<br>



            KSP of A00<br>              KSP Object:              (fieldsplit_0_)               1 MPI processes<br>                type: gmres<br>                  GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement<br>



                  GMRES: happy breakdown tolerance 1e-30<br>                maximum iterations=10000, initial guess is zero<br>                tolerances:  relative=1e-05, absolute=1e-50, divergence=10000<br>                left preconditioning<br>



                using PRECONDITIONED norm type for convergence test<br>              PC Object:              (fieldsplit_0_)               1 MPI processes<br>                type: ilu<br>                  ILU: out-of-place factorization<br>



                  0 levels of fill<br>                  tolerance for zero pivot 2.22045e-14<br>                  using diagonal shift on blocks to prevent zero pivot<br>                  matrix ordering: natural<br>                  factor fill ratio given 1, needed 1<br>



                    Factored matrix follows:<br>                      Matrix Object:                       1 MPI processes<br>                        type: seqaij<br>                        rows=2187, cols=2187<br>                        package used to perform factorization: petsc<br>



                        total: nonzeros=140625, allocated nonzeros=140625<br>                        total number of mallocs used during MatSetValues calls =0<br>                          using I-node routines: found 729 nodes, limit used is 5<br>



                linear system matrix = precond matrix:<br>                Matrix Object:                 1 MPI processes<br>                  type: seqaij<br>                  rows=2187, cols=2187<br>                  total: nonzeros=140625, allocated nonzeros=140625<br>



                  total number of mallocs used during MatSetValues calls =0<br>                    using I-node routines: found 729 nodes, limit used is 5<br>            A01<br>              Matrix Object:               1 MPI processes<br>



                type: seqaij<br>                rows=2187, cols=729<br>                total: nonzeros=46875, allocated nonzeros=46875<br>                total number of mallocs used during MatSetValues calls =0<br>                  using I-node routines: found 729 nodes, limit used is 5<br>



        Matrix Object:         1 MPI processes<br>          type: seqaij<br>          rows=729, cols=729<br>          total: nonzeros=15625, allocated nonzeros=15625<br>          total number of mallocs used during MatSetValues calls =0<br>



            not using I-node routines<br>  linear system matrix = precond matrix:<br>  Matrix Object:   1 MPI processes<br>    type: seqaij<br>    rows=2916, cols=2916, bs=4<br>    total: nonzeros=250000, allocated nonzeros=250000<br>



    total number of mallocs used during MatSetValues calls =0<br>      using I-node routines: found 729 nodes, limit used is 5<br><br><br><br></div><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 dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><div><br></div><div>or</div><div><br></div><div>  PCFieldSplitSetDMSplits(pc, PETSC_FALSE)</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</div>



<div><div><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 dir="ltr"><div class="gmail_extra"><div class="gmail_quote">
<div></div><div>The errors I get when running with options: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_0_fields 0,1,2 -pc_fieldsplit_1_fields 3<br>
[0]PETSC ERROR: --------------------- Error Message ------------------------------------<br>[0]PETSC ERROR: No support for this operation for this object type!<br>[0]PETSC ERROR: Support only implemented for 2d!<br>[0]PETSC ERROR: ------------------------------------------------------------------------<br>





[0]PETSC ERROR: Petsc Release Version 3.4.2, Jul, 02, 2013 <br>[0]PETSC ERROR: See docs/changes/index.html for recent updates.<br>[0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.<br>[0]PETSC ERROR: See docs/index.html for manual pages.<br>





[0]PETSC ERROR: ------------------------------------------------------------------------<br>[0]PETSC ERROR: src/AdLemMain on a arch-linux2-cxx-debug named edwards by bkhanal Tue Aug  6 17:35:30 2013<br>[0]PETSC ERROR: Libraries linked from /home/bkhanal/Documents/softwares/petsc-3.4.2/arch-linux2-cxx-debug/lib<br>





[0]PETSC ERROR: Configure run at Fri Jul 19 14:25:01 2013<br>[0]PETSC ERROR: Configure options --with-cc=gcc --with-fc=g77 --with-cxx=g++ --download-f-blas-lapack=1 --download-mpich=1 -with-clanguage=cxx --download-hypre=1<br>





[0]PETSC ERROR: ------------------------------------------------------------------------<br>[0]PETSC ERROR: DMCreateSubDM_DA() line 188 in /home/bkhanal/Documents/softwares/petsc-3.4.2/src/dm/impls/da/dacreate.c<br>[0]PETSC ERROR: DMCreateSubDM() line 1267 in /home/bkhanal/Documents/softwares/petsc-3.4.2/src/dm/interface/dm.c<br>





[0]PETSC ERROR: PCFieldSplitSetDefaults() line 337 in /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/impls/fieldsplit/fieldsplit.c<br>[0]PETSC ERROR: PCSetUp_FieldSplit() line 458 in /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/impls/fieldsplit/fieldsplit.c<br>





[0]PETSC ERROR: PCSetUp() line 890 in /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/pc/interface/precon.c<br>[0]PETSC ERROR: KSPSetUp() line 278 in /home/bkhanal/Documents/softwares/petsc-3.4.2/src/ksp/ksp/interface/itfunc.c<br>





</div><div>[0]PETSC ERROR: solveModel() line 181 in "unknowndirectory/"/user/bkhanal/home/works/AdLemModel/src/PetscAdLemTaras3D.cxx<br>WARNING! There are options you set that were not used!<br>WARNING! could be spelling mistake, etc!<br>





Option left: name:-pc_fieldsplit_1_fields value: 3<br><br><br><br></div><div><br></div><div><br></div><div><br></div><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 dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><div><br></div><div>   Matt</div><div><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 dir="ltr"><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 dir="ltr"><div class="gmail_extra"><div class="gmail_quote">

<div><br></div><div>   Matt</div><div>
<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 dir="ltr"><div class="gmail_extra"><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">

<div dir="ltr"><div class="gmail_extra"><div><br></div>

<div>   Matt</div><span><font color="#888888"><span><font color="#888888"><div>
<div><br></div>-- <br>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></font></span></div></div><span><font color="#888888">
</font></span></blockquote></div><span><font color="#888888"><br></font></span></div></div><span><font color="#888888">
</font></span></blockquote></div></div><span><font color="#888888"><div><br><br clear="all"><span><font color="#888888"><div><br></div>-- <br>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
</font></span></div></font></span></div></div><span><font color="#888888">
</font></span></blockquote></div><span><font color="#888888"><br></font></span></div></div><span><font color="#888888">
</font></span></blockquote></div></div><span><font color="#888888"><div><br><br clear="all"><span><font color="#888888"><div><br></div>-- <br>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
</font></span></div></font></span></div></div><span><font color="#888888">
</font></span></blockquote></div><span><font color="#888888"><br></font></span></div></div><span><font color="#888888">
</font></span></blockquote></div></div></div><span><font color="#888888"><div><div><br><br clear="all"><span class="HOEnZb"><font color="#888888"><div><br></div>-- <br>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
</font></span></div></div></font></span></div></div><span class="HOEnZb"><font color="#888888">
</font></span></blockquote></div><span class="HOEnZb"><font color="#888888"><br></font></span></div></div><span class="HOEnZb"><font color="#888888">
</font></span></blockquote></div></div></div><span class="HOEnZb"><font color="#888888"><div><div><br><br clear="all"><div><br></div>-- <br>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></font></span></div></div>
</blockquote></div><br></div></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br>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>