<div dir="ltr"><div dir="ltr">On Thu, Aug 20, 2020 at 7:45 AM Thibault Bridel-Bertomeu <<a href="mailto:thibault.bridelbertomeu@gmail.com">thibault.bridelbertomeu@gmail.com</a>> wrote:<br></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"><div>Dear all, <br></div><div><br></div><div>I have a finite volume code inspired from TS example ex11.c (with a riemann solver, etc...).</div><div>So far, I used only explicit time stepping through the TSSSP, and to set the RHS of my hyperbolic system I used :</div><div><pre width="80">TSSetType(ts, TSSSP);<br>DMTSSetRHSFunctionLocal(dm, <a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/DMPlexTSComputeRHSFunctionFVM.html#DMPlexTSComputeRHSFunctionFVM" target="_blank">DMPlexTSComputeRHSFunctionFVM</a>, &ctx);</pre></div><div>after setting the right Riemann solver in the TS associated to the DM.</div><div>Now, in some cases where the physics is stationary, I would like to reach the steady state faster and use an implicit timestepping operator to do so.</div><div>After looking through the examples, especially TS examples ex48.c and ex53.c, I simply tried to set</div><div></div><div></div><div><pre width="80">TSSetType(ts, TSBEULER);<br>DMTSSetIFunctionLocal(dm, <a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/DMPlexTSComputeIFunctionFEM.html#DMPlexTSComputeIFunctionFEM" target="_blank">DMPlexTSComputeIFunctionFEM</a>, &ctx);<br>DMTSSetIJacobianLocal(dm, <a href="https://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/TS/DMPlexTSComputeIJacobianFEM.html#DMPlexTSComputeIJacobianFEM" target="_blank">DMPlexTSComputeIJacobianFEM</a>, &ctx);<br><br></pre><div>instead of the previous calls. It compiles fine, and it runs. However, nothing happens : it behaves like there is no time evolution at all, the solution does not change from its initial state.<br></div><div>From the source code, it is my understanding that the DMPlexTSComputeIFunctionFEM and DMPlexTSComputeIJacobianFEM methods, in spite of their names, call respectively DMPlexComputeResidual_Internal and DMPlexComputeJacobian_Internal that can handle a FVM discretization.</div><div><br></div><div>What am I missing ? Are there other steps to take before I can simply try to run with a finite volume discretization and an implicit time stepping algorithm ?</div><div><br></div><div>Thank you very much for your help in advance !</div></div></div></blockquote><div><br></div><div>I could never get the FVM stuff to make sense to me for implicit methods. Here is my problem understanding. If you have an FVM method, it decides</div><div>to move "stuff" from one cell to its neighboring cells depending on the solution to the Riemann problem on each face, which computed the flux. This is</div><div>fine unless the timestep is so big that material can flow through into the cells beyond the neighbor. Then I should have considered the effect of the</div><div>Riemann problem for those interfaces. That would be in the Jacobian, but I don't know how to compute that Jacobian. I guess you could do everything</div><div>matrix-free, but without a preconditioner it seems hard.</div><div><br></div><div>Operationally, I always mark FVM things explicit, so they are only handled by RHSFunction/Jacobian, except for the u_t term which I stick in the IJacobian.</div><div>I was never successful in running an implicit FVM and still do not understand it. I could not find any literature that I understood either.</div><div><br></div><div>If there is a definite thing you want changed, I might be able to do it.</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:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div dir="ltr"><div dir="ltr"><div><div dir="ltr"><div dir="ltr"><div><div><div><div><div>Thibault Bridel-Bertomeu<br>—<br></div></div></div></div>Eng, MSc, PhD</div><div>Research Engineer</div><div>CEA/CESTA</div><div>33114 LE BARP</div><div>Tel.: (+33)557046924</div><div>Mob.: (+33)611025322<br></div><div>Mail: <a href="mailto:thibault.bridelbertomeu@gmail.com" target="_blank">thibault.bridelbertomeu@gmail.com</a><br></div></div></div></div></div></div></div></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><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.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>