[petsc-dev] finite element energy functions

[Matthew Knepley]

On Mon, Dec 2, 2013 at 4:29 PM, Jed Brown <jedbrown at mcs.anl.gov> wrote:

> Geoffrey Irving <irving at naml.us> writes:
>
> > For weak form PDEs which arise as the gradient of an integrated energy
> > function, it would be great to have a routine to perform the
> > appropriate integral.  It would be best to structure this as a generic
> > "perform an integral" function, evaluating
> >
> >     E(f) = sum_i w_i f(u(x_i), grad u(x_i))
>
> There is a SNES interface for quadratic functionals, SNESSetObjective(),
> and it has been plumbed down to DMDASNESSetObjectiveLocal(), but there
> is no DMSNESSetObjectiveLocal.  I would add that, and a
> DMPlexComputeResidualFEM().
>
> Actually, I would have made DMPlexSNESSetFunctionFEM() and now
> DMPlexSNESSetObjectiveFEM(), which configures the SNES appropriately,
> rather than calling SNES methods and passing in pointers to library
> functions.  All the other DMs work this way.  Matt, why did you do it
> differently with DMPlex?
>

When you and Peter changed the way SNES used to work, I asked you how I
should
change Plex. You told me to do this:

    ierr = DMSNESSetFunctionLocal(dm,  (PetscErrorCode
(*)(DM,Vec,Vec,void*)) DMPlexComputeResidualFEM, &user);CHKERRQ(ierr);
    ierr = DMSNESSetJacobianLocal(dm,  (PetscErrorCode
(*)(DM,Vec,Mat,Mat,MatStructure*,void*)) DMPlexComputeJacobianFEM,
&user);CHKERRQ(ierr);

which I did. Moreover, this is how it always worked, namely specifying the
local function.
What exactly are you proposing to replace this?


>
> > for a set of FE fields u (plus auxiliary fields) and some quadrature
> > rule.  The quadrature rule might need to be passed in independently of
> > the PetscFE objects, since the integral might tie together several
> > different fields connected to several different PetscFE's, conceivably
> > with distinct quadrature rules.  The integral is only valid as a
> > discrete energy if all the quadrature rules match, but there might be
> > other integral application where this doesn't hold.
>
> How is any of this different from evaluating residuals?  Each element
> computes a scalar now, all of which are summed together, but otherwise I
> think it's the same mechanics.


It looks the same as DMPlexComputeL2Diff().


> > Relatedly, how does one ensure that a PetscFE object is Galerkin
> > (matching primal and dual spaces) so that the energy to be discretely
> > valid?
>
> How does one ensure that the weak form itself is symmetric and how do
> you ensure that it is not indefinite (like a saddle; the Lagrangian is
> not the energy)?  We can do consistency checks, like verifying that the
> provided residual function is the gradient of the objective (energy).
>

We do not have any way to check that the spaces are the same, and I think
its not trivial. Here I think we fall back on the user giving us the
correct input.

   Matt

-- 
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-dev/attachments/20131202/9bb38ad5/attachment.html>