John — I just added TRI3SUBDIVISION to quadrature_trap_2D.C as a case, and I think it’s working now. I am still debugging things one step at a time, and one place that I got a complaint from is fe_lagrange.C — I may be confusing myself right now, but I think it’s complaining because TRI3SUBDIVISION is listed under FIRST order along with TRI3, but I think TRI3SUBDIVISION should be FOURTH order. I think I am confused because TRI3SUBDIVISION is listed in fe_lagrange.C.
Thank you, Mike On Aug 28, 2019, at 12:05 PM, John Peterson <jwpeter...@gmail.com<mailto:jwpeter...@gmail.com>> wrote: Hi Mike, I think this (using QTrap with TRI3SUBDIVISION elements) should work just fine. Did you try it yet? Let us know and we can update the case statements for quadrature. There are many places like this in the library where TRI3 and TRI3SUBDIVISION elements can be used more or less interchangeably, so there's probably a bunch we haven't found yet. -- John On Wed, Aug 28, 2019 at 12:33 AM Lee, Jae Ho <jaeho...@live.unc.edu<mailto:jaeho...@live.unc.edu>> wrote: Hello, I am currently implementing something that makes use of TRI3SUBDIVISION elements, and one of the initial things that I want to do is to compute the nodal coordinates for each deformed element. I am trying to avoid the use of MeshFunction to do this. Of course since subdivision shape functions are not interpolatory, I can’t just use the solution of the coordinate system directly as the coordinates of the element vertices. So what I am trying to use QTrap that uses nodal quadrature points and evaluating the shape functions at the quadrature points for a given element to get something like std::unique_ptr<libMesh::FEBase> fe(libMesh::FEBase::build(2, libMesh::SUBDIVISION)); std::unique_ptr<libMesh::QBase> qrule = libMesh::QBase::build(libMesh::QTRAP, 2, libMesh::FIRST); fe->attach_quadrature_rule(qrule.get()); fe->reinit(elem); for (unsigned int qp = 0; qp < qrule->n_points(); ++qp) { for (unsigned int d = 0; d < 3; ++d) { for (unsigned int i = 0; i < phi.size(); ++i) { X_node[qp][d] += X_solution[i][d] * phi[i][qp]; } } } X_solution is where the solution of the coordinate system has been stored. Does this seem reasonable? But quadrature_trap_2D.C complained because TRI3SUBDIVISION is not supported. Is it safe to hardcode TRI3SUBDIVISION as a case in the code? Or is there a different way to go about it (e.g. is there a different nodal quadrature rule that supports TRI3SUBDIVISION)? Thank you, Mike _______________________________________________ Libmesh-users mailing list Libmesh-users@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/libmesh-users
