<html>
  <head>
    <meta http-equiv="Content-Type" content="text/html; charset=utf-8">
  </head>
  <body text="#000000" bgcolor="#FFFFFF">
    <p>Thank you, it clarifies a lot :)<br>
    </p>
    <blockquote type="cite">
      <div>Ignore is the wrong word. How do you form a preconditioner
        from an object with no values in it?</div>
      <div> </div>
      <blockquote class="gmail_quote" style="margin:0 0 0
        .8ex;border-left:1px #ccc solid;padding-left:1ex"> Can I still
        use something like PCSHELL?<br>
      </blockquote>
      <div><br>
      </div>
      <div>What exactly would you do in your shell? If you are only
        using the action of an operator, it</div>
      <div>usually equivalent to some Krylov method.</div>
    </blockquote>
    <br>
    Yes, I was thinking to use FD matrix free Jacobian and for
    preconditioning to use Krylov method to invert the linearized
    Jacobian which also would be FD matrix free. Sounds bad? <br>
    <br>
    Well, from your advices I convinced to try constructing actual
    Jacobian and use it as Jacobian and matrix from which petsc would
    construct preconditioner. <br>
    <br>
    p.s. Do I reply to emails correctly? or should I replay to
    <a class="moz-txt-link-abbreviated" href="mailto:petsc-users@mcs.anl.gov">petsc-users@mcs.anl.gov</a> only?  <br>
    <br>
    Thank you,<br>
    Oleksandr.<br>
    <br>
    <div class="moz-cite-prefix">On 04/29/2018 06:03 PM, Matthew Knepley
      wrote:<br>
    </div>
    <blockquote type="cite"
cite="mid:CAMYG4GmxrmV3VwiWU5DM8fN9t0GR54bGfQGtAjY8p2KFfouqjQ@mail.gmail.com">
      <meta http-equiv="Content-Type" content="text/html; charset=utf-8">
      <div dir="ltr">
        <div class="gmail_extra">
          <div class="gmail_quote">On Sun, Apr 29, 2018 at 7:44 PM,
            Oleksandr Koshkarov <span dir="ltr"><<a
                href="mailto:olk548@mail.usask.ca" target="_blank"
                moz-do-not-send="true">olk548@mail.usask.ca</a>></span>
            wrote:<br>
            <blockquote class="gmail_quote" style="margin:0 0 0
              .8ex;border-left:1px #ccc solid;padding-left:1ex">Thank
              you. A little clarification:<br>
              <br>
              <blockquote class="gmail_quote" style="margin:0 0 0
                .8ex;border-left:1px #ccc solid;padding-left:1ex">
                    The above uses matrix-free to do matrix-vector
                products for linear system but constructs the
                preconditioner by building the Jacobian via differencing
                and then using that matrix to build the preconditioner.<br>
              </blockquote>
              So SNESComputeJacobianDefaultColo<wbr>r will use memory
              (construct Jacobian) even if preconditioning is set to
              PCNONE and J=NULL? (that is what I saw in my example)<br>
            </blockquote>
            <div><br>
            </div>
            <div><span
style="color:rgb(34,34,34);font-family:arial,sans-serif;font-size:small;font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">1)
                SNESComputeJacobianDefaultColo</span><wbr
style="color:rgb(34,34,34);font-family:arial,sans-serif;font-size:small;font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial"><span
style="color:rgb(34,34,34);font-family:arial,sans-serif;font-size:small;font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">r()
                creates the Jacobian, fullstop.</span><br>
            </div>
            <div><span
style="color:rgb(34,34,34);font-family:arial,sans-serif;font-size:small;font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline"><br>
              </span></div>
            <div><span
style="color:rgb(34,34,34);font-family:arial,sans-serif;font-size:small;font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">2)
                When you pass J=NULL, the Jacobian is created
                automatically by DMDA</span></div>
            <div><span
style="color:rgb(34,34,34);font-family:arial,sans-serif;font-size:small;font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline"><br>
              </span></div>
            <div><span
style="color:rgb(34,34,34);font-family:arial,sans-serif;font-size:small;font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;white-space:normal;word-spacing:0px;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">3)
                The PC has no influence on the assembly process.</span></div>
            <div> </div>
            <blockquote class="gmail_quote" style="margin:0 0 0
              .8ex;border-left:1px #ccc solid;padding-left:1ex">
              <blockquote class="gmail_quote" style="margin:0 0 0
                .8ex;border-left:1px #ccc solid;padding-left:1ex">
                    This does not build the Jacobian (only has matrix
                free matrix vector products) so requires much less
                memory but likely for large problem sizes the linear
                system solve will be slow (require many iterations) or
                won't converge at all.   The conditioning of the linear
                system depends on the exact problem you are solving and
                the type of discretization you are using. There is no
                easy rules that always apply but for most
                discretizations of PDEs the number of iterations needed
                by the linear solver increases with the problem sizes.
                This means for most large problems matrix-free (without
                building any sort of jacobean and preconditioner) is
                impractical and one needs to pay the price of using more
                memory to get reasonable convergence.<br>
              </blockquote>
              So the following command  ignores precoditioning?</blockquote>
            <div><br>
            </div>
            <div>Ignore is the wrong word. How do you form a
              preconditioner from an object with no values in it?</div>
            <div> </div>
            <blockquote class="gmail_quote" style="margin:0 0 0
              .8ex;border-left:1px #ccc solid;padding-left:1ex"> Can I
              still use something like PCSHELL?<br>
            </blockquote>
            <div><br>
            </div>
            <div>What exactly would you do in your shell? If you are
              only using the action of an operator, it</div>
            <div>usually equivalent to some Krylov method.</div>
            <div> </div>
            <blockquote class="gmail_quote" style="margin:0 0 0
              .8ex;border-left:1px #ccc solid;padding-left:1ex">
              SNESSetJacobian(snes,J,J,MatMF<wbr>FDComputeJacobian,0);<br>
              <br>
              p.s. for my problems if probably be unrealistic to
              construct Jacobian (state size will start from N >
              1000^4).<br>
            </blockquote>
            <div><br>
            </div>
            <div>The number of dofs does not tell us anything. You would
              need to know the sparsity. People regularly solve</div>
            <div>problems with billions of unknowns.</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">
              Thank you,<br>
              <br>
              Oleksandr.<br>
              <br>
              <br>
              On 04/29/2018 05:26 PM, Smith, Barry F. wrote:<br>
              <blockquote class="gmail_quote" style="margin:0 0 0
                .8ex;border-left:1px #ccc solid;padding-left:1ex">
                <br>
                <blockquote class="gmail_quote" style="margin:0 0 0
                  .8ex;border-left:1px #ccc solid;padding-left:1ex">
                  On Apr 29, 2018, at 5:40 PM, Oleksandr Koshkarov <<a
                    href="mailto:olk548@mail.usask.ca" target="_blank"
                    moz-do-not-send="true">olk548@mail.usask.ca</a>>
                  wrote:<br>
                  <br>
                  Dear All,<br>
                  <br>
                  sorry for spam because of my poor PETSc knowledge (I
                  am just starting with this nice framework).<br>
                  <br>
                  I think, I figured part of it out. However, I want to
                  point that src/ts/examples/tutorials/ex15<wbr>.c is
                  misleading. (or maybe it is a bug?)<br>
                  <br>
                  in this example we have<br>
                  <br>
                  TSGetSNES(ts,&snes);<br>
                  MatCreateSNESMF(snes,&Jmf);<br>
                  SNESSetJacobian(snes,Jmf,J,SNE<wbr>SComputeJacobianDefault,NULL);<br>
                  <br>
                  // or this: SNESSetJacobian(snes,Jmf,J,SNE<wbr>SComputeJacobianDefaultColor,<wbr>0);<br>
                  <br>
                  which implies (I think) that Jacobian would be matrix
                  free. And if one would use PCNONE for preconditioning
                  the matrix would never be allocated. However, it seems
                  in reality it allocates matrix.<br>
                </blockquote>
                    The above uses matrix-free to do matrix-vector
                products for linear system but constructs the
                preconditioner by building the Jacobian via differencing
                and then using that matrix to build the preconditioner.<br>
                <br>
                <blockquote class="gmail_quote" style="margin:0 0 0
                  .8ex;border-left:1px #ccc solid;padding-left:1ex">
                  To avoid it, I used<br>
                  <br>
                  MatCreateSNESMF(snes,&J);<br>
                  <br>
                  SNESSetJacobian(snes,J,J,MatMF<wbr>FDComputeJacobian,0);<br>
                  <br>
                  which seems to work fine. I am not sure I fully
                  understand the difference and i have zero intuition
                  and I also have no idea what happens with
                  preconditioning in this case. If someone have some
                  useful comets, please share :) (I read the relevant
                  section in PETSc manual, but still not fully
                  understanding what I should use when)<br>
                </blockquote>
                     This does not build the Jacobian (only has matrix
                free matrix vector products) so requires much less
                memory but likely for large problem sizes the linear
                system solve will be slow (require many iterations) or
                won't converge at all.   The conditioning of the linear
                system depends on the exact problem you are solving and
                the type of discretization you are using. There is no
                easy rules that always apply but for most
                discretizations of PDEs the number of iterations needed
                by the linear solver increases with the problem sizes.
                This means for most large problems matrix-free (without
                building any sort of jacobean and preconditioner) is
                impractical and one needs to pay the price of using more
                memory to get reasonable convergence.<br>
                <br>
                   Barry<br>
                <br>
                <blockquote class="gmail_quote" style="margin:0 0 0
                  .8ex;border-left:1px #ccc solid;padding-left:1ex">
                  Thank you and again sorry for the spam,<br>
                  Oleksandr.<br>
                  <br>
                  On 04/28/2018 07:20 PM, Smith, Barry F. wrote:<br>
                  <blockquote class="gmail_quote" style="margin:0 0 0
                    .8ex;border-left:1px #ccc solid;padding-left:1ex">
                    ~/Src/petsc/src/ts/examples/tu<wbr>torials<br>
                    $ grep SNESComputeJacobianDefaultColo<wbr>r *.c<br>
                    ex10.c:    ierr = SNESSetJacobian(snes,A,B,SNESC<wbr>omputeJacobianDefaultColor,0);<wbr>CHKERRQ(ierr);<br>
                    ex15.c:      ierr = SNESSetJacobian(snes,Jmf,J,SNE<wbr>SComputeJacobianDefaultColor,<wbr>0);CHKERRQ(ierr);<br>
                    ex17.c:    ierr = SNESSetJacobian(snes,J,J,SNESC<wbr>omputeJacobianDefaultColor,0);<wbr>CHKERRQ(ierr);<br>
                    <br>
                        I don't think you need to explicitly create the
                    MatFDColoring object.<br>
                    <br>
                        Please take a look at ex15.c and see if you can
                    get it working like that example. If that doesn't
                    work let us know and we can take a closer look at
                    it.<br>
                    <br>
                    <br>
                        Barry<br>
                    <br>
                    <br>
                    <blockquote class="gmail_quote" style="margin:0 0 0
                      .8ex;border-left:1px #ccc solid;padding-left:1ex">
                      On Apr 28, 2018, at 8:05 PM, Oleksandr Koshkarov
                      <<a href="mailto:olk548@mail.usask.ca"
                        target="_blank" moz-do-not-send="true">olk548@mail.usask.ca</a>>
                      wrote:<br>
                      <br>
                      Hello All,<br>
                      <br>
                      I hope someone can help :) I think I am doing
                      something wrong, but cannot understand what. I
                      have a huge time dependent system with 3d DMDA
                      data structure and I am evolving it with explicit
                      Runge-Kutta by using TS and basically only using
                      "TSSetRHSFunction". Now I want to repeat it with
                      implicit time stepper (for now Crank-Nicolson) and
                      I am trying to provide finite difference Jacobian
                      and I am failing miserably. I also cannot find
                      appropriate example in PETSc tutorial (if you can
                      point me to working example, it would be great).<br>
                      <br>
                      Here is my best attempt (what wrong with it?):<br>
                      <br>
                         DMDACreate3d(PETSC_COMM_<wbr>WORLD,
                      DM_BOUNDARY_PERIODIC, DM_BOUNDARY_PERIODIC,
                      DM_BOUNDARY_PERIODIC,<br>
                                  DMDA_STENCIL_STAR,<br>
                                  NX, NY, NZ,<br>
                                  PETSC_DECIDE, PETSC_DECIDE,
                      PETSC_DECIDE,<br>
                                  2*3+NC*NS,<br>
                                  1,<br>
                                  NULL,  NULL, NULL,  &da);<br>
                         DMSetUp(da);<br>
                         DMCreateGlobalVector(da,&x);<br>
                         TSCreate(PETSC_COMM_WORLD,&<wbr>ts);<br>
                         TSSetProblemType(ts,TS_NONLIN<wbr>EAR);<br>
                         TSSetRHSFunction(ts,NULL,comp<wbr>ute_RHS,NULL);<br>
                         TSSetMaxTime(ts,T_FINAL);<br>
                         TSSetExactFinalTime(ts,TS_EXA<wbr>CTFINALTIME_STEPOVER);<br>
                         TSSetDM(ts,da);<br>
                         TSSetType(ts,TSCN); //it works with:
                      TSSetType(ts,TSRK);<br>
                         set_IC(da,x);<br>
                         TSSetTimeStep(ts,DT);<br>
                         TSSetSolution(ts,x);<br>
                      <br>
                         TSGetSNES(ts,&snes);<br>
                         SNESGetKSP(snes,&ksp);<br>
                         KSPGetPC(ksp,&pc);<br>
                         PCSetType(pc,PCNONE);<br>
                      <br>
                         DMSetMatType(da,MATAIJ);<br>
                         DMCreateMatrix(da,&J);<br>
                         ISColoring iscoloring;<br>
                         MatFDColoring  matfdcoloring;<br>
                         DMCreateColoring(da,IS_COLORI<wbr>NG_GLOBAL,&iscoloring);<br>
                         MatFDColoringCreate(J,iscolor<wbr>ing,&matfdcoloring);<br>
                      <br>
                         MatFDColoringSetType(matfdcol<wbr>oring,MATMFFD_DS);<br>
                      <br>
                      // I think I do something wrong in the following 3
                      lines<br>
                      <br>
                         PetscErrorCode (*temp_f)(SNES,Vec,Vec,void*);<br>
                         SNESGetFunction(snes,NULL,&te<wbr>mp_f,NULL);<br>
                         MatFDColoringSetFunction(matf<wbr>dcoloring,(PetscErrorCode
                      (*)(void))temp_f,NULL);<br>
                      <br>
                         MatFDColoringSetUp(J,iscolori<wbr>ng,matfdcoloring);<br>
                      <br>
                      SNESSetJacobian(snes,J,J,SNESC<wbr>omputeJacobianDefaultColor,<wbr>matfdcoloring);<br>
                         ISColoringDestroy(&<wbr>iscoloring);<br>
                      <br>
                         TSSolve(ts,x);<br>
                      <br>
                      Thank you,<br>
                      <br>
                      Oleksandr Koshkarov.<br>
                      <br>
                    </blockquote>
                  </blockquote>
                </blockquote>
              </blockquote>
              <br>
            </blockquote>
          </div>
          <br>
          <br clear="all">
          <div><br>
          </div>
          -- <br>
          <div class="gmail_signature" data-smartmail="gmail_signature">
            <div dir="ltr">
              <div>
                <div dir="ltr">
                  <div>What most experimenters take for granted before
                    they begin their experiments is infinitely more
                    interesting than any results to which their
                    experiments lead.<br>
                    -- Norbert Wiener</div>
                  <div><br>
                  </div>
                  <div><a href="http://www.caam.rice.edu/%7Emk51/"
                      target="_blank" moz-do-not-send="true">https://www.cse.buffalo.edu/~knepley/</a><br>
                  </div>
                </div>
              </div>
            </div>
          </div>
        </div>
      </div>
    </blockquote>
    <br>
  </body>
</html>