On Sat, Nov 1, 2014 at 12:20 PM, Justin Chang <[email protected]> wrote:
> Thanks for the response. On a related note, > > In the source code of SNES ex12, where is the PetscSection actually being > created? I can't seem to find anywhere in the code or its routines where > the creation has been explicitly called. Because when the > discretization/problem is define, I imagine a PetscSection for the fields > and constraints has to be invoked somewhere. It seems to me as soon as > DMPlexAddBoundary() is called it creates the PetscSectionConstraints but I > can't seem to really confirm this within the source of that function. The > reason I want to know this is because if I really were to declare my own > constraints, I want to make sure that the PetscSection for the field > variables has already been set up. > Much of the mechanism has been moved into the library. The reason I did this was to make things like nonlinear multigrid automatic, rather than making the user define a different section on each level. Let me explain what is happening. When we call DMGetDefaultSection(), if the section is not present it will call the createdefaultsection() member function: https://bitbucket.org/petsc/petsc/src/ea1a9a653238fa98fb68e87b49145608ab6e5301/src/dm/interface/dm.c?at=master#cl-2977 which for Plex is here https://bitbucket.org/petsc/petsc/src/ea1a9a653238fa98fb68e87b49145608ab6e5301/src/dm/impls/plex/plex.c?at=master#cl-5601 It uses the discretization information from PetscDS and the BC information from DMPlexAddBoundary() to build the inputs for DMPlexCreateSection() https://bitbucket.org/petsc/petsc/src/ea1a9a653238fa98fb68e87b49145608ab6e5301/src/dm/impls/plex/plex.c?at=master#cl-5739 I realize that call is still limited. For example, it constrains all dofs on a point, but sometimes you only want to constrain some. In the most general case, you would just build the Section manually, using the same tools I use in DMPlexCreateSection(): https://bitbucket.org/petsc/petsc/src/ea1a9a653238fa98fb68e87b49145608ab6e5301/src/dm/impls/plex/plex.c?at=master#cl-3204 Thanks, Matt > Thanks, > Justin > > > Sent from my iPhone > > On Nov 1, 2014, at 11:40 AM, Matthew Knepley <[email protected]> wrote: > > On Fri, Oct 31, 2014 at 4:15 PM, Justin Chang <[email protected]> wrote: > >> Matt, >> >> Thanks for the response. One more question: >> >> If I go with approach 2, will manually setting the constraints/indices >> and its values be compatible with the DMPlexSNESComputeResidual/JacobianFEM >> routines? When I look at the source code of those routines it seems the >> constrained values are added to the local solution vectors via >> DMPlexInsertBoundaryValuesFEM. If I choose not to use DMPlexAddBoundary and >> wish to manually declare my boundary values, i assume I should call >> DMPlexProjectField with mode INSERT_BC_VALUES? >> > > Yes, exactly. I use DMPlexProjectFunctionLabelLocal(), > > > https://bitbucket.org/petsc/petsc/src/ea1a9a653238fa98fb68e87b49145608ab6e5301/src/dm/impls/plex/plexfem.c?at=master#cl-576 > > Thanks, > > Matt > > >> Thanks, >> Justin >> >> Sent from my iPhone >> >> On Oct 30, 2014, at 4:24 PM, Matthew Knepley <[email protected]> wrote: >> >> On Thu, Oct 30, 2014 at 4:16 PM, Justin Chang <[email protected]> >> wrote: >> >>> Matt, thanks for the quick response. >>> >>> What about for (dirichlet) boundary conditions? Would it be possible to >>> do something similar for those, like using those PetscSectionSetConstraint >>> functions? >>> >> >> Yes. There are generally two ways of handling Dirichlet conditions: >> >> 1) Replace those rows of the Jacobian with the identity and put the >> boundary value in the rhs. I find >> this cumbersome for nonlinear solves. >> >> 2) Remove these unknowns from the system using >> PetscSectionSetConstraintDof() and ConstratinIndices(). >> >> The DMPlexAddBoundary() functions are automating this processing by >> marking boundaries using DMLabels, >> and then constraining dofs on those boundaries. >> >> Thanks, >> >> Matt >> >> >>> Thanks, >>> Justin >>> >>> On Thu, Oct 30, 2014 at 3:35 PM, Matthew Knepley <[email protected]> >>> wrote: >>> >>>> On Thu, Oct 30, 2014 at 3:29 PM, Justin Chang <[email protected]> >>>> wrote: >>>> >>>>> Hi all, >>>>> >>>>> So I am writing an FEM code where it reads input data (e.g., auxiliary >>>>> coefficients, source terms, etc) from a text file. I have preprocessed the >>>>> data text file so that each vertex point has its corresponding data. For >>>>> instance, if my source term for a diffusion problem has a sin or cos >>>>> function of the coordinates, then this data is already tabulated and >>>>> simply >>>>> needs to be fed into my PETSc program. The data text file containing both >>>>> the auxiliary data and the coordinate/connectivities will also be used to >>>>> provide the input arguments for the DMPlexCreateFromDAG() function. >>>>> >>>>> Normally, I would hard-code the sin/cos functions into the source code >>>>> itself, but i want my program to be "universal" in that it can take as >>>>> input any user-defined function for not only these auxiliaries but for >>>>> other things like boundary conditions. I see that PETSc has a lot of >>>>> examples on how to read data into vectors and matrices, but I guess my >>>>> question is is there a way to project data from a text file into a vector >>>>> that has already been created with a defined DM structure? >>>>> >>>> >>>> If you have the functions available, I think it is far more universal >>>> to use the function itself, since then you can be independent >>>> of mesh and discretization when specifying input and BC. >>>> >>>> However, if you want to read it in, and you guarantee that it matches >>>> the mesh and discretization, I think the easiest thing to do is >>>> demand that it come in the same order as the vertices and use >>>> >>>> VecGetArray(V, &a); >>>> for (v = vStart, i = 0; v < vEnd; ++v) { >>>> PetscSectionGetDof(s, v, &dof); >>>> PetscSectionGetOffset(s, v, &off); >>>> for (d = 0; d < dof; ++d) a[off+d] = text_data[i++]; >>>> } >>>> VecRestoreArray(V, &a); >>>> >>>> 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 >> >> > > > -- > 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
