[petsc-users] Mass matrix with PetscFE
Matthew Knepley
knepley at gmail.com
Thu Feb 25 18:22:04 CST 2016
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: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: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 <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
>>>>
>>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20160225/4a66858d/attachment-0001.html>
More information about the petsc-users
mailing list