[petsc-dev] KSP in DMPlex

[Smith, Barry F.]

> On Nov 14, 2019, at 12:30 PM, Faibussowitsch, Jacob via petsc-dev <petsc-dev at mcs.anl.gov> wrote:
> 
> Hello,
> 
> So I am trying to make a simple 5 pt stencil finite difference laplacian operator (https://akmitra.public.iastate.edu/aero361/design_web/Laplace.pdf) work on a dmplex mesh, but all of the examples that I have found in ksp/examples/tutorials all use DMDA. The only way that I have managed to get the KSP to even run is by creating the operator matrix myself by hand. The workflow for the working implementation is as follows:
> 
> DMPlexMakeBoxMesh -> KSPCreate -> KSPsetDM -> KSPDeactiveDM -> Make the matrix by hand ->  KSPSetOperator(matrix) -> KSPSetComputeRHS -> Solve
> 
> Given that I tether KSP to the DM by KSPSetDM, but then immediately deactivate it (because then it complains about mallocs), this doesn’t work when reading in an actual unstructured mesh. Does KSP just not support dmplex?

   Nothing to do with KSP supporting DMPLEX. 

   For a given DMDA the matrix structure is completely determined when using finite differences, thus DMDA can always create an empty sparse matrix ready for KSP. Then one can either fill in the sparse matrix (as do many of the examples) with values for a particular PDE or have the MatFDColoring code compute the matrix entries from the provided RHS function.

   For DMPLEX the flex object itself does not have enough information to create a sparse matrix; it would be ill-defined since the matrix depends on the finite element space, boundary conditions etc. So it would really be DMPLEX + PETSCFE and a bunch of stuff to fully define the matrix. 

   See for example snes/examples/tutorials/ex12.c this is a way over complicated with lots of options, but it is one of the most basic for use DMPLEX for a PDE.  See SetupDiscretization() for how it creates the FE and DS objects.


> 
> Best,
> 
> Jacob
>