[petsc-users] Mass matrix with PetscFE

[Julian Andrej]

Any news about this? I've seen you merged the dmforest branch into next.

On 2016-02-26 01:22, Matthew Knepley wrote:
> I am sorry about the delay. I have your example working but it exposed
> a bug in Plex so I need to push the fix first. I should have
> everything for you early next week.
> 
>   Thanks
> 
>      Matt
> 
> On Feb 25, 2016 2:04 AM, "Julian Andrej" <juan at tf.uni-kiel.de> wrote:
> 
>> After a bit of rethinking the problem, the discrepancy between the
>> size of matrix A and the mass matrix M arises because of the
>> Dirichlet boundary conditions. So why aren't the BCs not imposed on
>> the mass matrix? Do I need to handle Dirichlet BCs differently in
>> this context (like zero rows and put one the diagonal?)
>> 
>> On 24.02.2016 20 [1]:54, juan wrote:
>> I attached another example which creates the correct mass matrix
>> but also overwrites the DM for the SNES solve. Somehow i cannot
>> manage
>> to really copy the DM to dm_mass and use that. If i try to do that
>> with
>> DMClone(dm, &dm_mass) i get a smaller mass matrix (which is not of
>> size A).
>> 
>> Maybe this helps in the discussion.
>> 
>> Relevant code starts at line 455.
>> 
>> On 2016-02-24 15:03, Julian Andrej wrote:
>> Thanks Matt,
>> 
>> I attached the modified example.
>> 
>> the corresponding code (and only changes to ex12) is starting at
>> line
>> 832.
>> 
>> It also seems that the mass matrix is of size 169x169 and the
>> stiffness matrix is of dimension 225x225. I'd assume that if i
>> multiply test and trial function i'd get a matrix of same size (if
>> the
>> space/quadrature is the same for the stiffness matrix)
>> 
>> On 24.02.2016 14 [2]:56, Matthew Knepley wrote:
>> On Wed, Feb 24, 2016 at 7:47 AM, Julian Andrej <juan at tf.uni-kiel.de
>> <mailto:juan at tf.uni-kiel.de>> wrote:
>> 
>> 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.
>> 
>> Its difficult to know when I cannot see the whole code. You can
>> always
>> insert
>> 
>> MatViewFromOptions(M, NULL, "-mat_view");
>> 
>> Using
>> MatView(M, PETSC_VIEWER_STDOUT_WORLD);
>> 
>> gives me a reasonable output to stdout.
>> 
>> Good.
>> 
>> 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);
>> 
>> I cannot explain this, but it has to be something like you are
>> viewing
>> the matrix before it is
>> actually assembled. Feel free to send the code. It sounds like it is
>> mostly working.
>> 
>> Matt
>> 
>> Any hints?
>> 
>> On 24.02.2016 13 [3] <tel:24.02.2016%2013>: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>
>> <mailto: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
>> 
>> --
>> 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
> 
> 
> Links:
> ------
> [1] tel:24.02.2016%2020
> [2] tel:24.02.2016%2014
> [3] tel:24.02.2016%2013