<div dir="ltr"><div class="gmail_quote"><div dir="ltr">On Tue, Sep 4, 2018 at 7:34 PM Rochan Upadhyay <<a href="mailto:u.rochan@gmail.com">u.rochan@gmail.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr">Dear Petsc Users and Developers,<div><br></div><div>Upon looking at the function listings and their descriptions in the manual suggests that the only finite element that can be automatically used is the continuous Lagrange. Is it possible to "easily" use vector valued Hcurl and Hdiv elements (like RT, NED, BDM and other more exotic types) ?</div></div></blockquote><div><br></div><div>The short answer is no.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div> It seems that there are two methods for generating finite elements : (1) Compute nodal basis functions from orthogonal polynomials with DOFs defined as pointwise evaluations (PETSCSPACEPOLYNOMIAL) (2) Compute basis functions at quadrature points yourself and input them directly (PETSCSPACEPOINT). Even if you input the trial and test basis functions this way, how would we input the various types of Piola transforms needed to map real to reference space. It seems in the implementation at "EvaluateFieldJets" only the covariant Piola transform (using invJ) is hard coded.</div></div></blockquote><div><br></div><div>True. Toby has a branch that fixes this and will support Hdiv and Hcurl. However, I do not know how many weeks away we</div><div>are on this project since the semester has now started.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div>Also I cannot find a good documentation of Petsc FE capabilities (unlike the superb presentations describing the DM routines).</div></div></blockquote><div><br></div><div>Yes, this is somewhat on purpose. PetscFE/FV exist principally to generate examples for me, and to help us</div><div>understand how discretizations should interact with meshes and solvers. We have not yet decided whether</div><div>PETSc will officially support discretizations, or just use something like Firedrake or Libmesh. However, I think</div><div>pushing on this has made those packages better, so its not a waste of time.</div><div><br></div><div>So, I am willing to explain anything you want on the list, and try and translate those explanations to documentation.</div><div>Also, we accept documentation Pull Requests (even a few words long), and Patrick Sanan did a nice presentation</div><div>at PETSc 2018 on how to easily submit them through the web interface. I have since tried it and its easy.</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir="ltr"><div> For me, it is preferable to use Petsc directly for FEM as other packages like fenics, libmesh etc.are too bulky, require too much effort to learn and not so easy to experiment with.</div><div><br></div><div>Regards,</div><div>Rochan</div></div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>