<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Thu, Jan 30, 2014 at 8:36 AM, Luc Berger-Vergiat <span dir="ltr"><<a href="mailto:lb2653@columbia.edu" target="_blank">lb2653@columbia.edu</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
  
    
  
  <div bgcolor="#FFFFFF" text="#000000">
    <div>Thank you so much this helps a lot with
      my research. I hope that you can reuse that little example code
      for some documentation on DMShell. The features available thru
      that object are really good, it is just a little hard (at least
      for me) to work around all the functions call to set it up.<br>
      <br>
      For instance I understand that the DM can be used in the KSP to
      work on the Kmat operator but what is the purpose of the PCSetDM?
      Using the DM on the preconditioner operator passed to the KSP?<br></div></div></blockquote><div><br></div><div>KSPSetDM() calls PCSetDM() internally. All the solver objects use the DM for sizes of Vec/Mat objects, and also</div>
<div>for dividing up the problem.</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 bgcolor="#FFFFFF" text="#000000"><div>
      Thanks again for your help, I will try to transfer this to
      FORTRAN!<br>
      <br>
      Cheers!<br>
      <pre cols="72">Best,
Luc</pre>
      On 01/30/2014 08:52 AM, Matthew Knepley wrote:<br>
    </div>
    <blockquote type="cite">
      <div dir="ltr">
        <div class="gmail_extra">
          <div class="gmail_quote">On Wed, Jan 29, 2014 at 7:19 PM, Luc
            Berger-Vergiat <span dir="ltr"><<a href="mailto:lb2653@columbia.edu" target="_blank">lb2653@columbia.edu</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">
              <div bgcolor="#FFFFFF" text="#000000">
                <div>Hi Matt,<br>
                  I thought I would spend the time writing a small c
                  example of what I am doing in my FORTRAN code so that
                  you can more easily see what it is I am doing wrong.<br>
                  <br>
                  I attached the code and its makefile. I run the code
                  the following arguments:<br>
                  <br>
                  -ksp_type gmres -pc_type fieldsplit
                  -fieldsplit_0_pc_type ilu -fieldsplit_0_ksp_type
                  preonly -fieldsplit_1_pc_type jacobi
                  -fieldsplit_1_ksp_type preonly -ksp_view -ksp_monitor<br>
                </div>
              </div>
            </blockquote>
            <div><br>
            </div>
            <div>1) You need to also call PetscSectionSetDof() in order
              to get full sizes. I don't just add up the fields.</div>
            <div><br>
            </div>
            <div>2) You need to set the DM into the KSP using KSPSetDM()</div>
            <div><br>
            </div>
            <div>3) You need to make the DM inactive if you want to
              provide your own Mat/Vec using KSPSetDMActive()</div>
            <div><br>
            </div>
            <div>
              4) Then you can get rid of the special code for FS</div>
            <div><br>
            </div>
            <div>5) Your options were wrong because you named the fields</div>
            <div><br>
            </div>
            <div>I include a revised code sample, and you should run
              with</div>
            <div>
              <br>
            </div>
            <div>  -ksp_type gmres -pc_type fieldsplit
              -fieldsplit_Field_0_pc_type ilu
              -fieldsplit_Field_0_ksp_type preonly
              -fieldsplit_Field_1_pc_type jacobi
              -fieldsplit_Field_1_ksp_type preonly -ksp_view
              -ksp_monitor</div>
            <div><br>
            </div>
            <div> Thanks,</div>
            <div><br>
            </div>
            <div>     Matt</div>
            <div><br>
            </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">
              <div bgcolor="#FFFFFF" text="#000000">
                <div> and get the following output:<br>
                  <br>
                  start=1, End=9.<br>
                  PetscSection with 2 fields<br>
                    field 0 with 1 components<br>
                  Process 0:<br>
                    (   1) dim  1 offset   0<br>
                    (   2) dim  1 offset   0<br>
                    (   3) dim  1 offset   0<br>
                    (   4) dim  1 offset   0<br>
                    (   5) dim  0 offset   0<br>
                    (   6) dim  0 offset   0<br>
                    (   7) dim  0 offset   0<br>
                    (   8) dim  0 offset   0<br>
                    field 1 with 1 components<br>
                  Process 0:<br>
                    (   1) dim  0 offset   1<br>
                    (   2) dim  0 offset   1<br>
                    (   3) dim  0 offset   1<br>
                    (   4) dim  0 offset   1<br>
                    (   5) dim  1 offset   0<br>
                    (   6) dim  1 offset   0<br>
                    (   7) dim  1 offset   0<br>
                    (   8) dim  1 offset   0<br>
                  PC flag for GetDMSplits: true<br>
                    0 KSP Residual norm 0.000000000000e+00 <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>
                    using PRECONDITIONED norm type for convergence test<br>
                  PC Object: 1 MPI processes<br>
                    type: fieldsplit<br>
                      FieldSplit with MULTIPLICATIVE composition: total
                  splits = 2<br>
                      Solver info for each split is in the following KSP
                  objects:<br>
                      Split number 0 Defined by IS<br>
                      KSP Object:    (fieldsplit_Field_0_)     1 MPI
                  processes<br>
                        type: preonly<br>
                        maximum iterations=10000, initial guess is zero<br>
                        tolerances:  relative=1e-05, absolute=1e-50,
                  divergence=10000<br>
                        left preconditioning<br>
                        using NONE norm type for convergence test<br>
                      PC Object:    (fieldsplit_Field_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=0, cols=0<br>
                                package used to perform factorization:
                  petsc<br>
                                total: nonzeros=0, allocated nonzeros=0<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=0, cols=0<br>
                          total: nonzeros=0, allocated nonzeros=0<br>
                          total number of mallocs used during
                  MatSetValues calls =0<br>
                            not using I-node routines<br>
                      Split number 1 Defined by IS<br>
                      KSP Object:    (fieldsplit_Field_1_)     1 MPI
                  processes<br>
                        type: preonly<br>
                        maximum iterations=10000, initial guess is zero<br>
                        tolerances:  relative=1e-05, absolute=1e-50,
                  divergence=10000<br>
                        left preconditioning<br>
                        using NONE norm type for convergence test<br>
                      PC Object:    (fieldsplit_Field_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=0, cols=0<br>
                                package used to perform factorization:
                  petsc<br>
                                total: nonzeros=0, allocated nonzeros=0<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=0, cols=0<br>
                          total: nonzeros=0, allocated nonzeros=0<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=8, cols=8<br>
                      total: nonzeros=64, allocated nonzeros=160<br>
                      total number of mallocs used during MatSetValues
                  calls =8<br>
                        using I-node routines: found 2 nodes, limit used
                  is 5<br>
                  The solution to this problem is:<br>
                  Vector Object: 1 MPI processes<br>
                    type: seq<br>
                  0<br>
                  0<br>
                  0<br>
                  0<br>
                  0<br>
                  0<br>
                  0<br>
                  0<br>
                  WARNING! There are options you set that were not used!<br>
                  WARNING! could be spelling mistake, etc!<br>
                  Option left: name:-fieldsplit_0_ksp_type value:
                  preonly<br>
                  Option left: name:-fieldsplit_0_pc_type value: ilu<br>
                  Option left: name:-fieldsplit_1_ksp_type value:
                  preonly<br>
                  Option left: name:-fieldsplit_1_pc_type value: jacobi<br>
                  <br>
                  when I run it with the more simple arguements:<br>
                  <br>
                  -ksp_type gmres -pc_type jacobi -ksp_view -ksp_monitor<br>
                  <br>
                  I get the following output:<br>
                  <br>
                  start=1, End=9.<br>
                  PetscSection with 2 fields<br>
                    field 0 with 1 components<br>
                  Process 0:<br>
                    (   1) dim  1 offset   0<br>
                    (   2) dim  1 offset   0<br>
                    (   3) dim  1 offset   0<br>
                    (   4) dim  1 offset   0<br>
                    (   5) dim  0 offset   0<br>
                    (   6) dim  0 offset   0<br>
                    (   7) dim  0 offset   0<br>
                    (   8) dim  0 offset   0<br>
                    field 1 with 1 components<br>
                  Process 0:<br>
                    (   1) dim  0 offset   1<br>
                    (   2) dim  0 offset   1<br>
                    (   3) dim  0 offset   1<br>
                    (   4) dim  0 offset   1<br>
                    (   5) dim  1 offset   0<br>
                    (   6) dim  1 offset   0<br>
                    (   7) dim  1 offset   0<br>
                    (   8) dim  1 offset   0<br>
                    0 KSP Residual norm 4.759858191165e+00 <br>
                    1 KSP Residual norm 2.344421567248e+00 <br>
                    2 KSP Residual norm 3.208390394507e-01 <br>
                    3 KSP Residual norm 7.171256210359e-02 <br>
                    4 KSP Residual norm 1.301032901980e-02 <br>
                    5 KSP Residual norm 1.104121978197e-15 <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>
                    using PRECONDITIONED norm type for convergence test<br>
                  PC Object: 1 MPI processes<br>
                    type: jacobi<br>
                    linear system matrix = precond matrix:<br>
                    Matrix Object:   1 MPI processes<br>
                      type: seqaij<br>
                      rows=8, cols=8<br>
                      total: nonzeros=64, allocated nonzeros=160<br>
                      total number of mallocs used during MatSetValues
                  calls =8<br>
                        using I-node routines: found 2 nodes, limit used
                  is 5<br>
                  The solution to this problem is:<br>
                  Vector Object: 1 MPI processes<br>
                    type: seq<br>
                  -0.977778<br>
                  -0.0444444<br>
                  -0.844444<br>
                  -0.327778<br>
                  3<br>
                  1<br>
                  1.5<br>
                  3<br>
                  <br>
                  It is obvious that something goes wrong in the
                  fieldsplit, I'm just clueless to what it is. I think
                  that everything is setup correctly but I'm missing
                  something to get the KSP to recognize my splits...<br>
                  <br>
                  Thanks for your help.<br>
                  <pre cols="72">Best,
Luc</pre>
                  On 01/28/2014 11:20 AM, Matthew Knepley wrote:<br>
                </div>
                <blockquote type="cite">
                  <div dir="ltr">
                    <div class="gmail_extra">
                      <div class="gmail_quote">On Tue, Jan 28, 2014 at
                        10:11 AM, Luc Berger-Vergiat <span dir="ltr"><<a href="mailto:lb2653@columbia.edu" target="_blank">lb2653@columbia.edu</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">
                          <div bgcolor="#FFFFFF" text="#000000">
                            <div>What I don't really understand is why
                              the size of all the sub fields is zero?<br>
                              As you can see all the matrix object in my
                              fieldsplit preconditioner have
                              total:nonzeros=0, allocated nonzeros=0.<br>
                              <br>
                              This seems strange since I issued the
                              following command:<br>
                                   <em></em>call
                              DMSetDefaultSection(FSDM, FSSection, ierr)<br>
                            </div>
                          </div>
                        </blockquote>
                        <div><br>
                        </div>
                        <div>That section looks like you never called
                          PetscSectionSetUp().</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">
                          <div bgcolor="#FFFFFF" text="#000000">
                            <div> with a the following FSSection:<br>
                              PetscSection with 4 fields<br>
                                field 0 with 1 components<br>
                              Process 0:<br>
                                (   0) dim  0 offset   0<br>
                                (   1) dim  1 offset   0<br>
                                (   2) dim  0 offset   0<br>
                                (   3) dim  0 offset   0<br>
                                (   4) dim  1 offset   0<br>
                                (   5) dim  0 offset   0<br>
                                (   6) dim  0 offset   0<br>
                                (   7) dim  0 offset   0<br>
                                (   8) dim  0 offset   0<br>
                                (   9) dim  0 offset   0<br>
                                (  10) dim  0 offset   0<br>
                                (  11) dim  0 offset   0<br>
                                (  12) dim  0 offset   0<br>
                                (  13) dim  0 offset   0<br>
                                (  14) dim  0 offset   0<br>
                                (  15) dim  0 offset   0<br>
                                field 1 with 1 components<br>
                              Process 0:<br>
                                (   0) dim  1 offset   0<br>
                                (   1) dim  0 offset   1<br>
                                (   2) dim  1 offset   0<br>
                                (   3) dim  1 offset   0<br>
                                (   4) dim  0 offset   1<br>
                                (   5) dim  1 offset   0<br>
                                (   6) dim  0 offset   0<br>
                                (   7) dim  0 offset   0<br>
                                (   8) dim  0 offset   0<br>
                                (   9) dim  0 offset   0<br>
                                (  10) dim  0 offset   0<br>
                                (  11) dim  0 offset   0<br>
                                (  12) dim  0 offset   0<br>
                                (  13) dim  0 offset   0<br>
                                (  14) dim  0 offset   0<br>
                                (  15) dim  0 offset   0<br>
                                field 2 with 1 components<br>
                              Process 0:<br>
                                (   0) dim  0 offset   1<br>
                                (   1) dim  0 offset   1<br>
                                (   2) dim  0 offset   1<br>
                                (   3) dim  0 offset   1<br>
                                (   4) dim  0 offset   1<br>
                                (   5) dim  0 offset   1<br>
                                (   6) dim  1 offset   0<br>
                                (   7) dim  1 offset   0<br>
                                (   8) dim  1 offset   0<br>
                                (   9) dim  1 offset   0<br>
                                (  10) dim  1 offset   0<br>
                                (  11) dim  1 offset   0<br>
                                (  12) dim  0 offset   0<br>
                                (  13) dim  0 offset   0<br>
                                (  14) dim  0 offset   0<br>
                                (  15) dim  0 offset   0<br>
                                field 3 with 1 components<br>
                              Process 0:<br>
                                (   0) dim  0 offset   1<br>
                                (   1) dim  0 offset   1<br>
                                (   2) dim  0 offset   1<br>
                                (   3) dim  0 offset   1<br>
                                (   4) dim  0 offset   1<br>
                                (   5) dim  0 offset   1<br>
                                (   6) dim  0 offset   1<br>
                                (   7) dim  0 offset   1<br>
                                (   8) dim  0 offset   1<br>
                                (   9) dim  0 offset   1<br>
                                (  10) dim  0 offset   1<br>
                                (  11) dim  0 offset   1<br>
                                (  12) dim  1 offset   0<br>
                                (  13) dim  1 offset   0<br>
                                (  14) dim  1 offset   0<br>
                                (  15) dim  1 offset   0<br>
                              <br>
                              I thought that by using
                              DMSetDefaultSection() I would be done
                              setting the fields and that the fieldsplit
                              would detect that section and use it to
                              construct the splits.<br>
                              <br>
                              Should I use another command to tell the
                              PC to use the DM section?<br>
                              <pre cols="72">Best,
Luc</pre>
                              On 01/28/2014 10:25 AM, Matthew Knepley
                              wrote:<br>
                            </div>
                            <blockquote type="cite">
                              <div dir="ltr">
                                <div class="gmail_extra">
                                  <div class="gmail_quote">On Mon, Jan
                                    27, 2014 at 1:35 PM, Luc
                                    Berger-Vergiat <span dir="ltr"><<a href="mailto:lb2653@columbia.edu" target="_blank">lb2653@columbia.edu</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">
                                      <div bgcolor="#FFFFFF" text="#000000">
                                        <div>Thanks Matt,<br>
                                          this indeed setting the number
                                          of fields earlier solve my
                                          issue!<br>
                                          <br>
                                          I have now made some progress
                                          getting my DM setup but I am
                                          having troubles setting my
                                          splitfield options.<br>
                                          I ran my problem with the
                                          -ksp_view option to see what
                                          is going on and I guess that
                                          somehow the section that I
                                          define in my DM is not used by
                                          the preconditioner to split
                                          the fields.<br>
                                          Here is the output of PETSc<br>
                                          <br>
                                          <font face="Helvetica, Arial,
                                            sans-serif">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-08,
                                            absolute=1e-16,
                                            divergence=1e+16<br>
                                              left preconditioning<br>
                                              using PRECONDITIONED norm
                                            type for convergence test<br>
                                            PC Object: 1 MPI processes<br>
                                              type: fieldsplit<br>
                                                FieldSplit with
                                            MULTIPLICATIVE composition:
                                            total splits = 4<br>
                                                Solver info for each
                                            split is in the following
                                            KSP objects:<br>
                                                Split number 0 Defined
                                            by IS<br>
                                          </font></div>
                                      </div>
                                    </blockquote>
                                    <div><br>
                                    </div>
                                    <div>There are 4 splits here and
                                      they are defined by an IS. Why do
                                      you think that is not what you
                                      asked for?</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">
                                      <div bgcolor="#FFFFFF" text="#000000">
                                        <div><font face="Helvetica,
                                            Arial, sans-serif">     KSP
                                            Object:   
                                            (fieldsplit_Field_0_)     1
                                            MPI processes<br>
                                                  type: preonly<br>
                                                  maximum
                                            iterations=10000, initial
                                            guess is zero<br>
                                                  tolerances: 
                                            relative=1e-05,
                                            absolute=1e-50,
                                            divergence=10000<br>
                                                  left preconditioning<br>
                                                  using NONE norm type
                                            for convergence test<br>
                                                PC Object:   
                                            (fieldsplit_Field_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=0, cols=0<br>
                                                          package used
                                            to perform factorization:
                                            petsc<br>
                                                          total:
                                            nonzeros=0, allocated
                                            nonzeros=0<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=0, cols=0<br>
                                                    total: nonzeros=0,
                                            allocated nonzeros=0<br>
                                                    total number of
                                            mallocs used during
                                            MatSetValues calls =0<br>
                                                      not using I-node
                                            routines<br>
                                                Split number 1 Defined
                                            by IS<br>
                                                KSP Object:   
                                            (fieldsplit_Field_1_)     1
                                            MPI processes<br>
                                                  type: preonly<br>
                                                  maximum
                                            iterations=10000, initial
                                            guess is zero<br>
                                                  tolerances: 
                                            relative=1e-05,
                                            absolute=1e-50,
                                            divergence=10000<br>
                                                  left preconditioning<br>
                                                  using NONE norm type
                                            for convergence test<br>
                                                PC Object:   
                                            (fieldsplit_Field_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=0, cols=0<br>
                                                          package used
                                            to perform factorization:
                                            petsc<br>
                                                          total:
                                            nonzeros=0, allocated
                                            nonzeros=0<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=0, cols=0<br>
                                                    total: nonzeros=0,
                                            allocated nonzeros=0<br>
                                                    total number of
                                            mallocs used during
                                            MatSetValues calls =0<br>
                                                      not using I-node
                                            routines<br>
                                                Split number 2 Defined
                                            by IS<br>
                                                KSP Object:   
                                            (fieldsplit_Field_2_)     1
                                            MPI processes<br>
                                                  type: preonly<br>
                                                  maximum
                                            iterations=10000, initial
                                            guess is zero<br>
                                                  tolerances: 
                                            relative=1e-05,
                                            absolute=1e-50,
                                            divergence=10000<br>
                                                  left preconditioning<br>
                                                  using NONE norm type
                                            for convergence test<br>
                                                PC Object:   
                                            (fieldsplit_Field_2_)     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=0, cols=0<br>
                                                          package used
                                            to perform factorization:
                                            petsc<br>
                                                          total:
                                            nonzeros=0, allocated
                                            nonzeros=0<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=0, cols=0<br>
                                                    total: nonzeros=0,
                                            allocated nonzeros=0<br>
                                                    total number of
                                            mallocs used during
                                            MatSetValues calls =0<br>
                                                      not using I-node
                                            routines<br>
                                                Split number 3 Defined
                                            by IS<br>
                                                KSP Object:   
                                            (fieldsplit_Field_3_)     1
                                            MPI processes<br>
                                                  type: preonly<br>
                                                  maximum
                                            iterations=10000, initial
                                            guess is zero<br>
                                                  tolerances: 
                                            relative=1e-05,
                                            absolute=1e-50,
                                            divergence=10000<br>
                                                  left preconditioning<br>
                                                  using NONE norm type
                                            for convergence test<br>
                                                PC Object:   
                                            (fieldsplit_Field_3_)     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=0, cols=0<br>
                                                          package used
                                            to perform factorization:
                                            petsc<br>
                                                          total:
                                            nonzeros=0, allocated
                                            nonzeros=0<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=0, cols=0<br>
                                                    total: nonzeros=0,
                                            allocated nonzeros=0<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=16, cols=16<br>
                                                total: nonzeros=256,
                                            allocated nonzeros=256<br>
                                                total number of mallocs
                                            used during MatSetValues
                                            calls =0<br>
                                                  using I-node routines:
                                            found 4 nodes, limit used is
                                            5</font><br>
                                                <br>
                                          I am also attaching part of my
                                          code which I use to generate
                                          the DM, the KSP and the PC
                                          objects.<br>
                                          <pre cols="72">Best,
Luc</pre>
                                          On 01/25/2014 10:31 AM,
                                          Matthew Knepley wrote:<br>
                                        </div>
                                        <blockquote type="cite">
                                          <div dir="ltr">
                                            <div class="gmail_extra">
                                              <div class="gmail_quote">On
                                                Fri, Jan 24, 2014 at
                                                5:10 PM, Luc
                                                Berger-Vergiat <span dir="ltr"><<a href="mailto:lb2653@columbia.edu" target="_blank">lb2653@columbia.edu</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">Hi
                                                  all,<br>
                                                  I want to use PETSc as
                                                  a solver for an FEM
                                                  problem: modelization
                                                  of shear bands using a
                                                  4 fields mixed
                                                  formulation.<br>
                                                  <br>
                                                  So far I am just
                                                  trying to set up a two
                                                  fields in a Petsc
                                                  section, assign dof
                                                  too each field and
                                                  then set up the
                                                  section to pass it to
                                                  a DM. I am taking this
                                                  approach because in
                                                  general I want for
                                                  fieldsplit to work
                                                  purely on the
                                                  algebraic level
                                                  without knowledge of
                                                  boundary conditions or
                                                  geometry.<br>
                                                  <br>
                                                  As of now I have
                                                  issues when I try to
                                                  assign a point and its
                                                  associated degree of
                                                  freedom to a field
                                                  using:
                                                  PetscSectionSetFieldDof.<br>
                                                  Is this the correct
                                                  way to associate a
                                                  dof/point to a field?<br>
                                                </blockquote>
                                                <div><br>
                                                </div>
                                                <div>You have to set the
                                                  number of fields
                                                  before the chart. I am
                                                  updating the docs.</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">
                                                  I attached an example
                                                  code and its makefile.<span><font color="#888888"><span><font color="#888888"><br>
                                                          <br>
                                                          -- <br>
                                                          Best,<br>
                                                          Luc<br>
                                                          <br>
                                                        </font></span></font></span></blockquote>
                                                <span><font color="#888888"> </font></span></div>
                                              <span><font color="#888888"> <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></font></span></div>
                                            <span><font color="#888888">
                                              </font></span></div>
                                          <span><font color="#888888"> </font></span></blockquote>
                                        <span><font color="#888888"> <br>
                                          </font></span></div>
                                      <span><font color="#888888"> </font></span></blockquote>
                                    <span><font color="#888888"> </font></span></div>
                                  <span><font color="#888888"> <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></font></span></div>
                                <span><font color="#888888"> </font></span></div>
                              <span><font color="#888888"> </font></span></blockquote>
                            <span><font color="#888888"> <br>
                              </font></span></div>
                          <span><font color="#888888"> </font></span></blockquote>
                        <span><font color="#888888"> </font></span></div>
                      <span><font color="#888888"> <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></font></span></div><span class="HOEnZb"><font color="#888888">
                  </font></span></div><span class="HOEnZb"><font color="#888888">
                </font></span></blockquote><span class="HOEnZb"><font color="#888888">
                <br>
              </font></span></div><span class="HOEnZb"><font color="#888888">
            </font></span></blockquote><span class="HOEnZb"><font color="#888888">
          </font></span></div><span class="HOEnZb"><font color="#888888">
          <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
        </font></span></div>
      </div>
    </blockquote>
    <br>
  </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>