[petsc-users] Mass matrix with PetscFE

[Julian Andrej]

I'm now using the petsc git master branch.

I tried adding my code to the ex12

   DM dm_mass;
   PetscDS prob_mass;
   PetscFE fe;
   Mat M;
   PetscFECreateDefault(dm, user.dim, 1, PETSC_TRUE, NULL, -1, &fe);

   DMClone(dm, &dm_mass);
   DMGetDS(dm_mass, &prob_mass);
   PetscDSSetDiscretization(prob_mass, 0, (PetscObject) fe);
   PetscDSSetJacobian(prob_mass, 0, 0, mass_kernel, NULL, NULL, NULL);
   DMCreateMatrix(dm_mass, &M);

   MatSetOptionsPrefix(M, "M_";)

and receive the error on running
./exe -interpolate -refinement_limit 0.0125 -petscspace_order 2 
-M_mat_view binary

WARNING! There are options you set that were not used!
WARNING! could be spelling mistake, etc!
Option left: name:-M_mat_view value: binary

I don't know if the matrix is actually there and assembled or if the 
option is ommitted because something is wrong.

Using
MatView(M, PETSC_VIEWER_STDOUT_WORLD);

gives me a reasonable output to stdout.

But saving the matrix and analysing it in matlab, results in an all zero 
matrix.

PetscViewerBinaryOpen(PETSC_COMM_WORLD, "Mout",FILE_MODE_WRITE, &viewer);
MatView(M, viewer);

Any hints?


On 24.02.2016 13:58, Matthew Knepley wrote:
> On Wed, Feb 24, 2016 at 6:47 AM, Julian Andrej <juan at tf.uni-kiel.de
> <mailto:juan at tf.uni-kiel.de>> wrote:
>
>     Hi,
>
>     i'm trying to assemble a mass matrix with the PetscFE/DMPlex
>     interface. I found something in the examples of TAO
>
>     https://bitbucket.org/petsc/petsc/src/da8116b0e8d067e39fd79740a8a864b0fe207998/src/tao/examples/tutorials/ex3.c?at=master&fileviewer=file-view-default
>
>     but using the lines
>
>     DMClone(dm, &dm_mass);
>     DMSetNumFields(dm_mass, 1);
>     DMPlexCopyCoordinates(dm, dm_mass);
>     DMGetDS(dm_mass, &prob_mass);
>     PetscDSSetJacobian(prob_mass, 0, 0, mass_kernel, NULL, NULL, NULL);
>     PetscDSSetDiscretization(prob_mass, 0, (PetscObject) fe);
>     DMPlexSNESComputeJacobianFEM(dm_mass, u, M, M, NULL);
>     DMCreateMatrix(dm_mass, &M);
>
>     leads to errors in DMPlexSNESComputeJacobianFEM (u is a global vector).
>
>     I don't can understand the necessary commands until
>     DMPlexSNESComputeJacobianFEM. What does it do and why is it
>     necessary? (especially why does the naming involve SNES?)
>
>     Is there another/easier/better way to create a mass matrix (the
>     inner product of the function space and the test space)?
>
>
> 1) That example needs updating. First, look at SNES ex12 which is up to
> date.
>
> 2) I assume you are using 3.6. If you use the development version, you
> can remove DMPlexCopyCoordinates().
>
> 3) You need to create the matrix BEFORE calling the assembly
>
> 4) Always always always send the entire error messge
>
>    Matt
>
>     Regards
>     Julian Andrej
>
>
>
>
> --
> What most experimenters take for granted before they begin their
> experiments is infinitely more interesting than any results to which
> their experiments lead.
> -- Norbert Wiener