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.
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
