Thanks for the help! I printed out the number of constraints after each of your suggestions, and after setting the hanging node constraints the number of constraints is 0 (as it should be, I'm working with a 1D mesh that has been globally refined some number of times), but after interpolating the boundary values the number of constraints is only 1. Given that dim = 1 and I'm trying to enforce Dirichlet boundaries, I reckon there should be 2 constraints. Running the std::map version of interpolate_boundary_values, I still only get 1 constraint.
So my question is: is this a bug? I'm trying to step through it in GDB, but it hangs when I try to step into this function <https://github.com/dealii/dealii/blob/5ade4a135e8b8d73985a2cace4d41d010257c8e3/include/deal.II/numerics/vector_tools_boundary.templates.h#L410>, so I'm a little bit at a loss as to how to continue debugging. - Lucas On Tuesday, August 30, 2022 at 3:48:06 PM UTC-5 Wolfgang Bangerth wrote: > On 8/30/22 13:56, Lucas Myers wrote: > > > > To enforce these boundary conditions, I initialize constraints with: > > ``` > > constraints.clear(); > > dealii::DoFTools::make_hanging_node_constraints(dof_handler, > constraints); > > > > dealii::VectorTools:: > > interpolate_boundary_values(dof_handler, > > 0, > > > > dealii::Functions::ZeroFunction<dim>(1), > > constraints); > > constraints.close(); > > ``` > > I use the `distribute_local_to_global` in the assembly function, and > then once > > I've solved the linear system (directly, with UMFPACK), I use > > `constraints.distribute`. Not sure what else I need to do to enforce > those > > boundary conditions. > > > > I'm attaching the code (<300 lines including white space) in case that > is > > helpful. Any help is appreciated! > > Rather than trying to debug this myself, let me ask some questions: > * Right after the call you show above, are the necessary constraints in > your > 'constraint' object? > * Right before you call 'distribute_local_to_blobal()' are the necessary > constraints still in the 'constraints' object? > * Right before you call 'distribute()' after you solve the linear system, > are > the necessary constraints still in the 'constraints' object? > * Right after you call 'distribute()', are the entries in the solution > vector > that correspond to boundary nodes correct? > > If the answer to any of the questions above is 'no', then you've got a > spot to > further debug. (All of this might be easiest to check on a very coarse > mesh, > say a 2x2 mesh.) Right now, all you see is that the *end result* is wrong. > You've got to narrow things down to a smaller part of the code, and the > questions above should help you with it! > > Best > W. > > -- > ------------------------------------------------------------------------ > Wolfgang Bangerth email: bang...@colostate.edu > 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 dealii+unsubscr...@googlegroups.com. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/1ebbbc73-dd45-41e0-894f-3543f3ded1adn%40googlegroups.com.
