[petsc-users] DMSNESSet<Function,Jacobian>

[Matthew Knepley]

On Tue, Apr 2, 2013 at 1:47 PM, Dharmendar Reddy <dharmareddy84 at gmail.com>wrote:

> Hello,
>          I looked at the source code of SNESSetFunction after i sent the
> email. From what i understood, SNESSetFunction calls DMSNESSetFunction
> internally. I was like, why not just use SNESSetFunction with out passing
> DM. But i get it now, not all problems have to use DM. If i am using DM i
> should rather use DMSNESSetFunction. A better way is DMSNESSetFunctionLocal
> right ? DMSNesSetFunction i called inside DMSNESSetFunctionLocal.
>
>
> I did look at DMPlexComputeResidualFEM but that requires userContext to
> have PetscFEM object. Is PetscFEM object accisable from Fortran or is it
> genrated using
>

The object is just a struct, and has only information you must have for
FEM, namely
a quadrature rule, and tabulation of the basis functions are quadrature
points.


> bin/pythonscripts/PetscGenerateFEMQuadrature.py
>
> This script just generates a header with the struct information. You can
make this header
yourself and obviate the need for the script.


> I structured my code such that, user provides the Linear and bilinear form
> subroutines:
>
> subroutine biLinearFrom(u,v,cell,Kelm,eqnCtx)
> type(FunctionSpace) :: u,v
> type(FEMCell) :: cell
> PetscScalar :: Kelm
> class(*),pointer :: eqnCtx
> end subroutine
>
> Function space object provide the routines for evaluation of basis
> functions and gradients, given cell coordinates, and cell degrees of
> freedom.
>

Can you just use these to create the PetscFEM struct?


>
> So i should be able to get the cell geometry information and cell degree
> of freedom data via DMPlexVecGet/SetClosure right ?
>

Yes, although tell me why DMPlexComputeCellGeometry() does not work, and if
not you can look at the code
to see me getting the coordinate information.


> At this point i was able to create a DMPlex object from a mesh created by
> gmsh. I used DMPLexCreateLabel to add the physical region labels i created
> in GMSh . So no i have a mapping from given a region label, list of cells
> in that region.
>
> If i am doing a physical region wise assemble of the problem, i could use
> the DMPlexgetStratumIS for each regionlabel.
>

Yes.


> since i am doing a cell wise assemble, i would be interested in a mapping
> from cell id to  regionId (i can have a local array of labels indexed by
> region Id) so that i can pass the material information to FEMCell.  Right
> now i have a integer array which holds the map from cellId to regionId.
>
> I was wondering if i could add a integer section to DM object ?
>

You use DMLabelGetValue(), and to trade memory for speed
call DMLabelCreateIndex() first if you want.


> I also need to have some other auxiliary scalar data associated with the
> DM. You suggested in one of the earlier email  that i created another DM
> using DMplexClone ? If i clone a object is it creating another copy of the
> mesh in memory or is it just storing pointers ?
>

Just a pointer.


> I was thinking i could create a section using DMPlexCreateSection and pass
> that section to DMPlexVecGetClosure when i have to access auxilary data.
>

This is exactly what I was suggesting, but you have the section be the
"default section" of the cloned DM. That way
it interacts with all PETSc functions nicely and you only have the extra
memory for a DM object which is tiny.

    Matt

Basically the problem is to solve Poisson equation with current continuity
> equation for a semiconductor. I solve the problem by solving each equation
> until the solutions of both are consistent. For Poisson problem, the degree
> freedom is potential and takes charge density as input. For continuation
> equation degree of freedom is charge density and take potential as input.
>
> Do you have alternative approaches to what i am thinking ?
>
> On Mon, Apr 1, 2013 at 8:24 PM, Matthew Knepley <knepley at gmail.com> wrote:
>
>> On Tue, Apr 2, 2013 at 10:16 AM, Dharmendar Reddy <
>> dharmareddy84 at gmail.com> wrote:
>>
>>> Hello,
>>>         I am confused about the usage of SNESSetFunction,
>>> DMSNESSetFunction, DMSNESSetFunctionLocal and the corresponding Jacobin
>>> routines.
>>>
>>> At present my code is set up as follows.
>>>
>>> I initiate SNES with getFunction and getJacobain routines required as
>>> per SNESSet<Function,Jacobian>. I pass a user context  which has the mesh
>>> and field layout information. I assemble the residual function and Jacobin
>>> element wise using the information from user context. My code at this point
>>> is serial, so i run the element loop from 1 to numberofTotalElements.
>>>
>>> Now i want to switch to using DMPlex. How should i change the interface
>>> to my code.
>>>
>>
>> You can always look at SNES ex62
>> http://www.mcs.anl.gov/petsc/petsc-dev/src/snes/examples/tutorials/ex62.c.html
>>
>>
>>> If i set the dm for snes using SNESSetDM, i can see that the DM can be
>>> accessed via SNESGetDM inside getFunction and getJacobian, then i get the
>>> elementStartID and elementEndID to run the loop for assmbler. But now i
>>> will have to assemble to local vector into global vector right ?
>>>
>>
>> No, here is the sequence:
>>
>>   1) SNESSetDM(snes, dm)
>>
>>   2) DMSNESSetFunctionLocal(dm, userResidual, userCtx)
>>
>> where we have
>>
>>   userResidual(DM dm, Vec X, Vec F, void *userCtx)
>>
>> Both X and F are local vectors which I normally interact with using
>> DMPlexVecGet/SetClosure(). If you
>> are using FEM, you can DMPlexComputeResidualFEM() here and
>> use DMPlexSetFEMIntegration() to
>> input point-wise physics functions as is done in ex62.
>>
>>
>>> Looks like the DMSNESSetFunctionLocal will assemble the local
>>> vector into global vector. The function provided to this should just
>>> evaluate the local vector. Did i understand this right ?
>>>
>>> I am confused about passing the DM in DMSNESSet<Function,Jacobian>. When
>>> the DM can be accessed via SNESGetDM, why do we pass it again explicitly ?
>>>
>>
>> You are not passing the SNES, so where would it come from in this call?
>>
>
>
>
>
>
>>
>>    Matt
>>
>>
>>>
>>> --
>>> -----------------------------------------------------
>>> Dharmendar Reddy Palle
>>> Graduate Student
>>> Microelectronics Research center,
>>> University of Texas at Austin,
>>> 10100 Burnet Road, Bldg. 160
>>> MER 2.608F, TX 78758-4445
>>> e-mail: dharmareddy84 at gmail.com
>>> Phone: +1-512-350-9082
>>> United States of America.
>>> Homepage: https://webspace.utexas.edu/~dpr342
>>>
>>
>>
>>
>> --
>> 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
>>
>
>
>
> --
> -----------------------------------------------------
> Dharmendar Reddy Palle
> Graduate Student
> Microelectronics Research center,
> University of Texas at Austin,
> 10100 Burnet Road, Bldg. 160
> MER 2.608F, TX 78758-4445
> e-mail: dharmareddy84 at gmail.com
> Phone: +1-512-350-9082
> United States of America.
> Homepage: https://webspace.utexas.edu/~dpr342
>



-- 
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/20130402/deb7583d/attachment.html>