Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
dof_accessor_set.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2019 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #include <deal.II/dofs/dof_accessor.h>
17 #include <deal.II/dofs/dof_handler.h>
18 #include <deal.II/dofs/dof_levels.h>
19 
20 #include <deal.II/fe/fe.h>
21 
22 #include <deal.II/grid/tria_iterator.h>
23 #include <deal.II/grid/tria_iterator.templates.h>
24 
25 #include <deal.II/hp/dof_handler.h>
26 
27 #include <deal.II/lac/block_vector.h>
28 #include <deal.II/lac/la_parallel_block_vector.h>
29 #include <deal.II/lac/la_parallel_vector.h>
30 #include <deal.II/lac/la_vector.h>
31 #include <deal.II/lac/petsc_block_vector.h>
32 #include <deal.II/lac/petsc_vector.h>
33 #include <deal.II/lac/sparse_matrix.h>
34 #include <deal.II/lac/trilinos_parallel_block_vector.h>
35 #include <deal.II/lac/trilinos_tpetra_vector.h>
36 #include <deal.II/lac/trilinos_vector.h>
37 #include <deal.II/lac/vector.h>
38 
39 #include <vector>
40 
41 DEAL_II_NAMESPACE_OPEN
42 
43 
44 template <typename DoFHandlerType, bool lda>
45 template <class OutputVector, typename number>
46 void
48  const Vector<number> &local_values,
49  OutputVector & values,
50  const unsigned int fe_index) const
51 {
52  if (!this->has_children() && !this->is_artificial())
53  {
54  if ((dynamic_cast<DoFHandler<DoFHandlerType::dimension,
55  DoFHandlerType::space_dimension> *>(
56  this->dof_handler) != nullptr) ||
57  // for hp-DoFHandlers, we need to require that on
58  // active cells, you either don't specify an fe_index,
59  // or that you specify the correct one
60  (fe_index == this->active_fe_index()) ||
61  (fe_index == DoFHandlerType::default_fe_index))
62  // simply set the values on this cell
63  this->set_dof_values(local_values, values);
64  else
65  {
66  Assert(local_values.size() ==
67  this->dof_handler->get_fe(fe_index).dofs_per_cell,
68  ExcMessage("Incorrect size of local_values vector."));
69 
70  FullMatrix<double> interpolation(
71  this->get_fe().dofs_per_cell,
72  this->dof_handler->get_fe(fe_index).dofs_per_cell);
73 
74  this->get_fe().get_interpolation_matrix(
75  this->dof_handler->get_fe(fe_index), interpolation);
76 
77  // do the interpolation to the target space. for historical
78  // reasons, matrices are set to size 0x0 internally even
79  // we reinit as 4x0, so we have to treat this case specially
80  Vector<number> tmp(this->get_fe().dofs_per_cell);
81  if ((tmp.size() > 0) && (local_values.size() > 0))
82  interpolation.vmult(tmp, local_values);
83 
84  // now set the dof values in the global vector
85  this->set_dof_values(tmp, values);
86  }
87  }
88  else
89  // otherwise distribute them to the children
90  {
91  Assert((dynamic_cast<DoFHandler<DoFHandlerType::dimension,
92  DoFHandlerType::space_dimension> *>(
93  this->dof_handler) != nullptr) ||
94  (fe_index != DoFHandlerType::default_fe_index),
95  ExcMessage(
96  "You cannot call this function on non-active cells "
97  "of hp::DoFHandler objects unless you provide an explicit "
98  "finite element index because they do not have naturally "
99  "associated finite element spaces associated: degrees "
100  "of freedom are only distributed on active cells for which "
101  "the active_fe_index has been set."));
102 
103  const FiniteElement<dim, spacedim> &fe =
104  this->get_dof_handler().get_fe(fe_index);
105  const unsigned int dofs_per_cell = fe.dofs_per_cell;
106 
107  Assert(this->dof_handler != nullptr,
108  typename BaseClass::ExcInvalidObject());
109  Assert(local_values.size() == dofs_per_cell,
110  typename BaseClass::ExcVectorDoesNotMatch());
111  Assert(values.size() == this->dof_handler->n_dofs(),
112  typename BaseClass::ExcVectorDoesNotMatch());
113 
114  Vector<number> tmp(dofs_per_cell);
115 
116  for (unsigned int child = 0; child < this->n_children(); ++child)
117  {
118  if (tmp.size() > 0)
119  fe.get_prolongation_matrix(child, this->refinement_case())
120  .vmult(tmp, local_values);
121  this->child(child)->set_dof_values_by_interpolation(tmp,
122  values,
123  fe_index);
124  }
125  }
126 }
127 
128 
129 // --------------------------------------------------------------------------
130 // explicit instantiations
131 #include "dof_accessor_set.inst"
132 
133 DEAL_II_NAMESPACE_CLOSE
void set_dof_values_by_interpolation(const Vector< number > &local_values, OutputVector &values, const unsigned int fe_index=DoFHandlerType::default_fe_index) const
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
Definition: fe.cc:335
static ::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1407
const unsigned int dofs_per_cell
Definition: fe_base.h:282
size_type size() const