Hi Wolfgang,
Thanks a lot for clarifying! I decided to modify the code to adopt a p::shared::T model with METIS and realized anisotropic refinement doesn’t seem to work for this case either. The error comes from AffineConstraints’ unacceptance of any refinement cases that are not isotropic. The same error can be reproduced for my serial code as well. If I change the RefinementCase, when it’s, say, cut_yz, to cut_xyz, hence going back to isotropic refinement, both codes run without a problem. My code is implementing a space-time IP-HDG method. It would be nice to have anisotropic refinements in order to separate local time stepping and spatial refinement within a space-time slab. Here is the part of the code where AffineConstraints gets involved and it’s pretty much copied verbatim from step-51: constraints.clear(); DoFTools::make_hanging_node_constraints(dof_handler, constraints); std::map<types::boundary_id, const Function<dim> *> boundary_functions; Solution<dim> solution_function(nu); boundary_functions[0] = &solution_function; VectorTools::interpolate_boundary_values(dof_handler, boundary_functions, constraints); constraints.close(); And here is the runtime error message: An error occurred in line <1102> of file <path/to/deal.II/9.4.1-petsc-hypre/source/source/dofs/dof_tools_constraints.cc> in function void dealii::DoFTools::internal::make_hp_hanging_node_constraints(const dealii::DoFHandler&, dealii::AffineConstraints&) [with int dim = 3; int spacedim = 3; number = double] The violated condition was: cell->face(face)->refinement_case() == RefinementCase::isotropic_refinement Additional information: You are trying to use functionality in deal.II that is currently not implemented. In many cases, this indicates that there simply didn’t appear much of a need for it, or that the author of the original code did not have the time to implement a particular case. If you hit this exception, it is therefore worth the time to look into the code to find out whether you may be able to implement the missing functionality. If you do, please consider providing a patch to the deal.II development sources (see the deal.II website on how to contribute). I’m wondering if it would be ill-advised to simply remove the assertion and re-compile the library. If so, I’m thinking about going a bit deeper and see if I can come up with a patch. In that case, I would really appreciate some insights in terms of AffineConstraints’ incompatibilities, if any, with hanging nodes created by anisotropic refinement. Thanks, Greg On Tuesday, April 4, 2023 at 8:18:58 PM UTC Wolfgang Bangerth wrote: > > Greg: > > > We want to construct a 3D triangulation by extruding a 2D triangulation > > (one that potentially contains hanging nodes) and we only want one > > slice/layer of mesh on the extrusion direction. > > > > Looking around in the GridGenerator namespace led me to the > > extrude_triangulation function. It’s doing everything we desire except > > for that (a) the slices/layers on the extrusion direction has to be at > > least two > > I think this is poorly described in the documentation. The number of > slices = the number of cell layers plus one. Two slices => one layer of > cells. I've fixed this here: > https://github.com/dealii/dealii/pull/15028 > > > and (b) the 2D mesh must be a coarse mesh. I’m wondering if > > there are tips on getting around these two restrictions. > > This you can't get around. > > > > Originally I was using GridGenerator::hyper_rectangle teamed with > > anisotropic refinement cut_xy. But this ceases to work after the code > > was re-implemented with parallel::distributed::Triangulation because > > “this class does not support anisotropic refinement, because it relies > > on the p4est library that does not support this” [1]. > > Yes, and this is true regardless of how you generate the mesh. You can't > create an anisotropically refined mesh for p::d::T. This also means that > you cannot extrude a refined 2d mesh into 3d that has only one layer of > cells. You can extrude a coarse mesh, though, and then refine the > resulting mesh -- it will then have more than one layer in z-direction > in some locations, however. > > Best > W. > > -- > ------------------------------------------------------------------------ > Wolfgang Bangerth email: [email protected] > www: http://www.math.colostate.edu/~bangerth/ > -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to [email protected]. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/8470ea57-c826-45fa-8784-3e8bc4b4a581n%40googlegroups.com.
