Dear Daniel, For non-primitive elements I believe its is necessary to use the project_boundary_values_<xyz>() functions. So, maybe project_boundary_values_curl_conforming() <https://dealii.org/current/doxygen/deal.II/group__constraints.html#ga1c6685360c01c9c46eeb7575e8ef68ac> or project_boundary_values_curl_conforming_l2() <https://dealii.org/current/doxygen/deal.II/group__constraints.html#ga6fca3672ae63b249402460a6ed4538b4> might be one of the functions that you are looking for.
Best, Jean-Paul > On 14 Jun 2019, at 14:54, Daniel Garcia-Sanchez <[email protected]> > wrote: > > Hi Dhananjay, > > Thanks for your message. I removed the dof_handler by mistake when I did the > copy/paste. > > I still get the error when I run the following code: > > -------------------------------------------------------- > constraints.clear(); > > > constraints.reinit(locally_relevant_dofs); > > > DoFTools::make_hanging_node_constraints(dof_handler, constraints); > VectorTools::interpolate_boundary_values(dof_handler, 1, ZeroFunction<dim, > number>(components), constraints, ComponentMask({false, true, true})); > -------------------------------------------------------- > > On Friday, June 14, 2019 at 2:36:13 PM UTC+2, Dhananjay Phansalkar wrote: > Hello Daniel, > I think you are missing first argument for the function > "VectorTools::interpolate_boundary_values" I guess it requires dofhandler > object . Have a look at step 6 > (https://www.dealii.org/current/doxygen/deal.II/step_6.html#Step6setup_system > <https://www.dealii.org/current/doxygen/deal.II/step_6.html#Step6setup_system>) > > Cheers > > Dhananjay > > On Friday, June 14, 2019 at 1:57:42 PM UTC+2, Daniel Garcia-Sanchez wrote: > Hi, > > I'm starting to use the NedelecSZ element to simulate the maxwell equations > in the frequency/harmonic domain. I started to simulate the resonant modes of > a cube that has metallic boundaries. The boundary condition for a metal is n > x E = 0. That means that the tangential components of the electrical field > are zero. > > If I try to use the following code: > > -------------------------------------------------------- > constraints.clear(); > > > constraints.reinit(locally_relevant_dofs); > > > DoFTools::make_hanging_node_constraints(dof_handler, constraints); > VectorTools::interpolate_boundary_values(1, ZeroFunction<dim, > number>(components), constraints, ComponentMask({false, true, true})); > -------------------------------------------------------- > > I obtain the following error: > > -------------------------------------------------------- > > > An error occurred in line <2967> of file > </opt/src/dealii/include/deal.II/numerics/vector_tools.templates.h> in > function > > void dealii::VectorTools::internal::do_interpolate_boundary_values(const > M_or_MC<dim, spacedim> &, const DoFHandlerType<dim, spacedim> &, const > std::map<types::boundary_id, const Function<spacedim, number> *> > &, std::map<types::global_dof_index, number> &, const dealii::ComponentMask > &) [dim = 3, spacedim = 3, number = std::complex<double>, DoFHandlerType = > DoFHandler, M_or_MC = Mapping] > The violated condition was: > > > cell->get_fe().is_primitive(i) > > > Additional information: > > > This function can only deal with requested boundary values that > correspond to primitive (scalar) base elements. You may want to look up in > the deal.II glossary what the term 'primitive' means. > > > > There are alternative boundary value interpolation functions in namespace > 'VectorTools' that you can use for non-primitive finite elements. > -------------------------------------------------------- > > How can I set the boundary condition n x E = 0 for the NedelecSZ element in > 3D? > > Thanks! > Daniel > > -- > The deal.II project is located at http://www.dealii.org/ > <http://www.dealii.org/> > For mailing list/forum options, see > https://groups.google.com/d/forum/dealii?hl=en > <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] > <mailto:[email protected]>. > To view this discussion on the web visit > https://groups.google.com/d/msgid/dealii/b88aaeb0-a754-4d50-87db-4ade7f8f8b47%40googlegroups.com > > <https://groups.google.com/d/msgid/dealii/b88aaeb0-a754-4d50-87db-4ade7f8f8b47%40googlegroups.com?utm_medium=email&utm_source=footer>. > For more options, visit https://groups.google.com/d/optout > <https://groups.google.com/d/optout>. -- 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/B6EAF4BF-36E4-4C42-AFC7-76B06AEFF0B1%40gmail.com. For more options, visit https://groups.google.com/d/optout.
