[petsc-users] DMPlex and Tao

[Justin Chang]

As far as I know this equality constraint should exist outside of the large
square matrix.

What I was thinking of doing is manually creating an *Nele* by *Ndof*
matrix, each row will have the same number of non-zeros. Since the size of
*Ndof* is simply the global number of free dofs formulated by the
PetscSection, I want to know what's the best way to map the free dofs into
the global column indices [0 *Ndof*)

For example, in src/dm/impls/plex/examples/tutorials/ex1.c where a box mesh
with three fields is created, there are a total of 41 dofs to which 8 of
them are constrained. If no constraints were set, then the size of u (i.e.
*Ndof*) would be 41, but with constraints set it is 33. When you view the
constrained global vector u, it lists all the dofs in the order as
specified by the PetscSection but skips all the constrained ones. Would I
have to use something like DMSet/GetDefaultGlobalSection to do the mapping?

Thanks,
Justin


On Fri, Aug 1, 2014 at 7:26 AM, Matthew Knepley <knepley at gmail.com> wrote:

> On Wed, Jul 30, 2014 at 5:01 PM, Justin Chang <jychang48 at gmail.com> wrote:
>
>> That's good to know thanks. I have another question somewhat related.
>>
>> In that quadratic programming problem I have shown earlier, my K matrix
>> is of global size *Ndof* by *Ndof* where *Ndof* is number of free dofs
>> (for instance, if I have a 2x2 quadrilateral mesh and want to evaluate
>> concentration at all the verticies, I will have nine total dofs but if I
>> have Dirichlet constraints all along the boundary, only my center vertex is
>> free so *Ndof* = 1. Thus when I simply call DMCreateMatrix(...) it would
>> give me a 1x1 matrix.)
>>
>> However, I want my M constraint matrix to be of size *Nele* by *Ndof*
>> where *Nele* is the total number of elements. How would I create this
>> matrix? Because ultimately, I am subjecting all of my free u's to *Nele*
>> number of equality constraints, and if possible I want this M matrix to be
>> compatible with the layout from Vec u.
>>
>
> This is not hard in principle. I used to have code that did this, but I
> took it out during the rewrite because no one
> was using it and it was complicated. It just generalizes the code in
> DMPlexPreallocateOperator() to take two
> PetscSections instead of one.
>
> Are you sure that this does not belong inside a large square matrix for
> the entire problem?
>
>   Thanks,
>
>      Matt
>
>
>>
>> On Wed, Jul 30, 2014 at 2:03 PM, Matthew Knepley <knepley at gmail.com>
>> wrote:
>>
>>> On Wed, Jul 30, 2014 at 2:53 PM, Justin Chang <jychang48 at gmail.com>
>>> wrote:
>>>
>>>> Hi,
>>>>
>>>> This might be a silly question, but I am developing a finite element
>>>> code to solve the non-negative diffusion equation. I am using the DMPlex
>>>> structure and plan on subjecting the problem to a bounded constraint
>>>> optimization to ensure non-negative concentrations. Say for instance, I
>>>> want to solve:
>>>>
>>>> min || K*u - F ||^2
>>>> subject to u >= 0
>>>>              M*u = 0
>>>>
>>>> Where K is the Jacobian, F the residual, u the concentration, and M
>>>> some equality constraint matrix. I know how to do this via Tao, but my
>>>> question is can Tao handle Mats and Vecs created via DMPlex? Because from
>>>> the SNES and TS examples I have seen in ex12, ex62, and ex11, it seems
>>>> there are solvers tailored to handle DMPlex created data structures.
>>>>
>>>
>>> Yes, TAO can handle objects created by DMPlex since they are just Vec
>>> and Mat structures. The additions
>>> to ex12, ex62, and ex11 concern the callbacks from the solver to the
>>> constructions routines for residual and
>>> Jacobian. Right now, none of the callbacks are automatic for TAO so you
>>> would have to code them yourself,
>>> probably by just calling the routines from SNES or TS so it should not
>>> be that hard. We are working to get
>>> everything integrated as it is in the other solvers.
>>>
>>>   Thanks,
>>>
>>>      Matt
>>>
>>>
>>>> Thanks,
>>>> Justin
>>>>
>>>
>>>
>>>
>>> --
>>> 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/20140801/631e9eac/attachment.html>