Reference documentation for deal.II version 8.4.1
dof_tools_sparsity.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2015 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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #include <deal.II/base/thread_management.h>
17 #include <deal.II/base/quadrature_lib.h>
18 #include <deal.II/base/table.h>
19 #include <deal.II/base/template_constraints.h>
20 #include <deal.II/base/utilities.h>
21 #include <deal.II/lac/sparsity_pattern.h>
22 #include <deal.II/lac/dynamic_sparsity_pattern.h>
23 #include <deal.II/lac/trilinos_sparsity_pattern.h>
24 #include <deal.II/lac/block_sparsity_pattern.h>
25 #include <deal.II/lac/vector.h>
26 #include <deal.II/lac/constraint_matrix.h>
27 #include <deal.II/grid/tria.h>
28 #include <deal.II/grid/tria_iterator.h>
29 #include <deal.II/grid/intergrid_map.h>
30 #include <deal.II/grid/grid_tools.h>
31 #include <deal.II/dofs/dof_handler.h>
32 #include <deal.II/dofs/dof_accessor.h>
33 #include <deal.II/fe/fe.h>
34 #include <deal.II/fe/fe_values.h>
35 #include <deal.II/fe/fe_tools.h>
36 #include <deal.II/hp/fe_collection.h>
37 #include <deal.II/hp/q_collection.h>
38 #include <deal.II/hp/fe_values.h>
39 #include <deal.II/dofs/dof_tools.h>
40 #include <deal.II/numerics/vector_tools.h>
41 
42 
43 #include <algorithm>
44 #include <numeric>
45 
46 DEAL_II_NAMESPACE_OPEN
47 
48 
49 
50 namespace DoFTools
51 {
52 
53  template <typename DoFHandlerType, typename SparsityPatternType>
54  void
55  make_sparsity_pattern (const DoFHandlerType &dof,
56  SparsityPatternType &sparsity,
57  const ConstraintMatrix &constraints,
58  const bool keep_constrained_dofs,
59  const types::subdomain_id subdomain_id)
60  {
61  const types::global_dof_index n_dofs = dof.n_dofs();
62  (void)n_dofs;
63 
64  Assert (sparsity.n_rows() == n_dofs,
65  ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
66  Assert (sparsity.n_cols() == n_dofs,
67  ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
68 
69  // If we have a distributed::Triangulation only allow locally_owned
70  // subdomain. Not setting a subdomain is also okay, because we skip
71  // ghost cells in the loop below.
72  Assert (
73  (dof.get_triangulation().locally_owned_subdomain() == numbers::invalid_subdomain_id)
74  ||
75  (subdomain_id == numbers::invalid_subdomain_id)
76  ||
77  (subdomain_id == dof.get_triangulation().locally_owned_subdomain()),
78  ExcMessage ("For parallel::distributed::Triangulation objects and "
79  "associated DoF handler objects, asking for any subdomain other "
80  "than the locally owned one does not make sense."));
81 
82  std::vector<types::global_dof_index> dofs_on_this_cell;
83  dofs_on_this_cell.reserve (max_dofs_per_cell(dof));
84  typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(),
85  endc = dof.end();
86 
87  // In case we work with a distributed sparsity pattern of Trilinos
88  // type, we only have to do the work if the current cell is owned by
89  // the calling processor. Otherwise, just continue.
90  for (; cell!=endc; ++cell)
91  if (((subdomain_id == numbers::invalid_subdomain_id)
92  ||
93  (subdomain_id == cell->subdomain_id()))
94  &&
95  cell->is_locally_owned())
96  {
97  const unsigned int dofs_per_cell = cell->get_fe().dofs_per_cell;
98  dofs_on_this_cell.resize (dofs_per_cell);
99  cell->get_dof_indices (dofs_on_this_cell);
100 
101  // make sparsity pattern for this cell. if no constraints pattern
102  // was given, then the following call acts as if simply no
103  // constraints existed
104  constraints.add_entries_local_to_global (dofs_on_this_cell,
105  sparsity,
106  keep_constrained_dofs);
107  }
108  }
109 
110 
111 
112  template <typename DoFHandlerType, typename SparsityPatternType>
113  void
114  make_sparsity_pattern (const DoFHandlerType &dof,
115  const Table<2,Coupling> &couplings,
116  SparsityPatternType &sparsity,
117  const ConstraintMatrix &constraints,
118  const bool keep_constrained_dofs,
119  const types::subdomain_id subdomain_id)
120  {
121  const types::global_dof_index n_dofs = dof.n_dofs();
122  (void)n_dofs;
123 
124  Assert (sparsity.n_rows() == n_dofs,
125  ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
126  Assert (sparsity.n_cols() == n_dofs,
127  ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
128  Assert (couplings.n_rows() == dof.get_fe().n_components(),
129  ExcDimensionMismatch(couplings.n_rows(), dof.get_fe().n_components()));
130  Assert (couplings.n_cols() == dof.get_fe().n_components(),
131  ExcDimensionMismatch(couplings.n_cols(), dof.get_fe().n_components()));
132 
133  // If we have a distributed::Triangulation only allow locally_owned
134  // subdomain. Not setting a subdomain is also okay, because we skip
135  // ghost cells in the loop below.
136  Assert (
137  (dof.get_triangulation().locally_owned_subdomain() == numbers::invalid_subdomain_id)
138  ||
139  (subdomain_id == numbers::invalid_subdomain_id)
140  ||
141  (subdomain_id == dof.get_triangulation().locally_owned_subdomain()),
142  ExcMessage ("For parallel::distributed::Triangulation objects and "
143  "associated DoF handler objects, asking for any subdomain other "
144  "than the locally owned one does not make sense."));
145 
147 
148  // first, for each finite element, build a mask for each dof, not like
149  // the one given which represents components. make sure we do the right
150  // thing also with respect to non-primitive shape functions, which
151  // takes some additional thought
152  std::vector<Table<2,bool> > dof_mask(fe_collection.size());
153 
154  // check whether the table of couplings contains only true arguments,
155  // i.e., we do not exclude any index. that is the easy case, since we
156  // don't have to set up the tables
157  bool need_dof_mask = false;
158  for (unsigned int i=0; i<couplings.n_rows(); ++i)
159  for (unsigned int j=0; j<couplings.n_cols(); ++j)
160  if (couplings(i,j) == none)
161  need_dof_mask = true;
162 
163  if (need_dof_mask == true)
164  for (unsigned int f=0; f<fe_collection.size(); ++f)
165  {
166  const unsigned int dofs_per_cell = fe_collection[f].dofs_per_cell;
167 
168  dof_mask[f].reinit (dofs_per_cell, dofs_per_cell);
169 
170  for (unsigned int i=0; i<dofs_per_cell; ++i)
171  for (unsigned int j=0; j<dofs_per_cell; ++j)
172  if (fe_collection[f].is_primitive(i) &&
173  fe_collection[f].is_primitive(j))
174  dof_mask[f](i,j)
175  = (couplings(fe_collection[f].system_to_component_index(i).first,
176  fe_collection[f].system_to_component_index(j).first) != none);
177  else
178  {
179  const unsigned int first_nonzero_comp_i
180  = fe_collection[f].get_nonzero_components(i).first_selected_component();
181  const unsigned int first_nonzero_comp_j
182  = fe_collection[f].get_nonzero_components(j).first_selected_component();
183  Assert (first_nonzero_comp_i < fe_collection[f].n_components(),
184  ExcInternalError());
185  Assert (first_nonzero_comp_j < fe_collection[f].n_components(),
186  ExcInternalError());
187 
188  dof_mask[f](i,j)
189  = (couplings(first_nonzero_comp_i,first_nonzero_comp_j) != none);
190  }
191  }
192 
193 
194  std::vector<types::global_dof_index> dofs_on_this_cell(fe_collection.max_dofs_per_cell());
195  typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(),
196  endc = dof.end();
197 
198  // In case we work with a distributed sparsity pattern of Trilinos
199  // type, we only have to do the work if the current cell is owned by
200  // the calling processor. Otherwise, just continue.
201  for (; cell!=endc; ++cell)
202  if (((subdomain_id == numbers::invalid_subdomain_id)
203  ||
204  (subdomain_id == cell->subdomain_id()))
205  &&
206  cell->is_locally_owned())
207  {
208  const unsigned int fe_index = cell->active_fe_index();
209  const unsigned int dofs_per_cell =fe_collection[fe_index].dofs_per_cell;
210 
211  dofs_on_this_cell.resize (dofs_per_cell);
212  cell->get_dof_indices (dofs_on_this_cell);
213 
214 
215  // make sparsity pattern for this cell. if no constraints pattern
216  // was given, then the following call acts as if simply no
217  // constraints existed
218  constraints.add_entries_local_to_global (dofs_on_this_cell,
219  sparsity,
220  keep_constrained_dofs,
221  dof_mask[fe_index]);
222  }
223  }
224 
225 
226 
227  template <typename DoFHandlerType, typename SparsityPatternType>
228  void
229  make_sparsity_pattern (const DoFHandlerType &dof_row,
230  const DoFHandlerType &dof_col,
231  SparsityPatternType &sparsity)
232  {
233  const types::global_dof_index n_dofs_row = dof_row.n_dofs();
234  const types::global_dof_index n_dofs_col = dof_col.n_dofs();
235  (void)n_dofs_row;
236  (void)n_dofs_col;
237 
238  Assert (sparsity.n_rows() == n_dofs_row,
239  ExcDimensionMismatch (sparsity.n_rows(), n_dofs_row));
240  Assert (sparsity.n_cols() == n_dofs_col,
241  ExcDimensionMismatch (sparsity.n_cols(), n_dofs_col));
242 
243 //TODO: Looks like wasteful memory management here
244 
245  const std::list<std::pair<typename DoFHandlerType::cell_iterator,
246  typename DoFHandlerType::cell_iterator> >
247  cell_list
248  = GridTools::get_finest_common_cells (dof_row, dof_col);
249 
250 
251  typename std::list<std::pair<typename DoFHandlerType::cell_iterator,
252  typename DoFHandlerType::cell_iterator> >::const_iterator
253  cell_iter = cell_list.begin();
254 
255  for (; cell_iter!=cell_list.end(); ++cell_iter)
256  {
257  const typename DoFHandlerType::cell_iterator cell_row = cell_iter->first;
258  const typename DoFHandlerType::cell_iterator cell_col = cell_iter->second;
259 
260  if (!cell_row->has_children() && !cell_col->has_children())
261  {
262  const unsigned int dofs_per_cell_row =
263  cell_row->get_fe().dofs_per_cell;
264  const unsigned int dofs_per_cell_col =
265  cell_col->get_fe().dofs_per_cell;
266  std::vector<types::global_dof_index>
267  local_dof_indices_row(dofs_per_cell_row);
268  std::vector<types::global_dof_index>
269  local_dof_indices_col(dofs_per_cell_col);
270  cell_row->get_dof_indices (local_dof_indices_row);
271  cell_col->get_dof_indices (local_dof_indices_col);
272  for (unsigned int i=0; i<dofs_per_cell_row; ++i)
273  sparsity.add_entries (local_dof_indices_row[i],
274  local_dof_indices_col.begin(),
275  local_dof_indices_col.end());
276  }
277  else if (cell_row->has_children())
278  {
279  const std::vector<typename DoFHandlerType::active_cell_iterator >
280  child_cells = GridTools::get_active_child_cells<DoFHandlerType> (cell_row);
281  for (unsigned int i=0; i<child_cells.size(); i++)
282  {
283  const typename DoFHandlerType::cell_iterator
284  cell_row_child = child_cells[i];
285  const unsigned int dofs_per_cell_row =
286  cell_row_child->get_fe().dofs_per_cell;
287  const unsigned int dofs_per_cell_col =
288  cell_col->get_fe().dofs_per_cell;
289  std::vector<types::global_dof_index>
290  local_dof_indices_row(dofs_per_cell_row);
291  std::vector<types::global_dof_index>
292  local_dof_indices_col(dofs_per_cell_col);
293  cell_row_child->get_dof_indices (local_dof_indices_row);
294  cell_col->get_dof_indices (local_dof_indices_col);
295  for (unsigned int r=0; r<dofs_per_cell_row; ++r)
296  sparsity.add_entries (local_dof_indices_row[r],
297  local_dof_indices_col.begin(),
298  local_dof_indices_col.end());
299  }
300  }
301  else
302  {
303  std::vector<typename DoFHandlerType::active_cell_iterator>
304  child_cells = GridTools::get_active_child_cells<DoFHandlerType> (cell_col);
305  for (unsigned int i=0; i<child_cells.size(); i++)
306  {
307  const typename DoFHandlerType::active_cell_iterator
308  cell_col_child = child_cells[i];
309  const unsigned int dofs_per_cell_row =
310  cell_row->get_fe().dofs_per_cell;
311  const unsigned int dofs_per_cell_col =
312  cell_col_child->get_fe().dofs_per_cell;
313  std::vector<types::global_dof_index>
314  local_dof_indices_row(dofs_per_cell_row);
315  std::vector<types::global_dof_index>
316  local_dof_indices_col(dofs_per_cell_col);
317  cell_row->get_dof_indices (local_dof_indices_row);
318  cell_col_child->get_dof_indices (local_dof_indices_col);
319  for (unsigned int r=0; r<dofs_per_cell_row; ++r)
320  sparsity.add_entries (local_dof_indices_row[r],
321  local_dof_indices_col.begin(),
322  local_dof_indices_col.end());
323  }
324  }
325  }
326  }
327 
328 
329 
330  template <typename DoFHandlerType, typename SparsityPatternType>
331  void
333  (const DoFHandlerType &dof,
334  const std::vector<types::global_dof_index> &dof_to_boundary_mapping,
335  SparsityPatternType &sparsity)
336  {
337  if (DoFHandlerType::dimension == 1)
338  {
339  // there are only 2 boundary indicators in 1d, so it is no
340  // performance problem to call the other function
341  typename DoFHandlerType::FunctionMap boundary_ids;
342  boundary_ids[0] = 0;
343  boundary_ids[1] = 0;
344  make_boundary_sparsity_pattern<DoFHandlerType, SparsityPatternType>
345  (dof,
346  boundary_ids,
347  dof_to_boundary_mapping,
348  sparsity);
349  return;
350  }
351 
352  const types::global_dof_index n_dofs = dof.n_dofs();
353  (void)n_dofs;
354 
355  AssertDimension (dof_to_boundary_mapping.size(), n_dofs);
356  AssertDimension (sparsity.n_rows(), dof.n_boundary_dofs());
357  AssertDimension (sparsity.n_cols(), dof.n_boundary_dofs());
358 #ifdef DEBUG
359  if (sparsity.n_rows() != 0)
360  {
361  types::global_dof_index max_element = 0;
362  for (std::vector<types::global_dof_index>::const_iterator i=dof_to_boundary_mapping.begin();
363  i!=dof_to_boundary_mapping.end(); ++i)
364  if ((*i != DoFHandlerType::invalid_dof_index) &&
365  (*i > max_element))
366  max_element = *i;
367  AssertDimension (max_element, sparsity.n_rows()-1);
368  };
369 #endif
370 
371  std::vector<types::global_dof_index> dofs_on_this_face;
372  dofs_on_this_face.reserve (max_dofs_per_face(dof));
373 
374  // loop over all faces to check whether they are at a boundary. note
375  // that we need not take special care of single lines (using
376  // @p{cell->has_boundary_lines}), since we do not support boundaries of
377  // dimension dim-2, and so every boundary line is also part of a
378  // boundary face.
379  typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(),
380  endc = dof.end();
381  for (; cell!=endc; ++cell)
382  for (unsigned int f=0; f<GeometryInfo<DoFHandlerType::dimension>::faces_per_cell; ++f)
383  if (cell->at_boundary(f))
384  {
385  const unsigned int dofs_per_face = cell->get_fe().dofs_per_face;
386  dofs_on_this_face.resize (dofs_per_face);
387  cell->face(f)->get_dof_indices (dofs_on_this_face,
388  cell->active_fe_index());
389 
390  // make sparsity pattern for this cell
391  for (unsigned int i=0; i<dofs_per_face; ++i)
392  for (unsigned int j=0; j<dofs_per_face; ++j)
393  sparsity.add (dof_to_boundary_mapping[dofs_on_this_face[i]],
394  dof_to_boundary_mapping[dofs_on_this_face[j]]);
395  }
396  }
397 
398 
399 
400  template <typename DoFHandlerType, typename SparsityPatternType>
402  (const DoFHandlerType &dof,
403  const typename FunctionMap<DoFHandlerType::space_dimension>::type &boundary_ids,
404  const std::vector<types::global_dof_index> &dof_to_boundary_mapping,
405  SparsityPatternType &sparsity)
406  {
407  if (DoFHandlerType::dimension == 1)
408  {
409  // first check left, then right boundary point
410  for (unsigned int direction=0; direction<2; ++direction)
411  {
412  // if this boundary is not requested, then go on with next one
413  if (boundary_ids.find(direction) ==
414  boundary_ids.end())
415  continue;
416 
417  // find active cell at that boundary: first go to left/right,
418  // then to children
419  typename DoFHandlerType::level_cell_iterator cell = dof.begin(0);
420  while (!cell->at_boundary(direction))
421  cell = cell->neighbor(direction);
422  while (!cell->active())
423  cell = cell->child(direction);
424 
425  const unsigned int dofs_per_vertex = cell->get_fe().dofs_per_vertex;
426  std::vector<types::global_dof_index> boundary_dof_boundary_indices (dofs_per_vertex);
427 
428  // next get boundary mapped dof indices of boundary dofs
429  for (unsigned int i=0; i<dofs_per_vertex; ++i)
430  boundary_dof_boundary_indices[i]
431  = dof_to_boundary_mapping[cell->vertex_dof_index(direction,i)];
432 
433  for (unsigned int i=0; i<dofs_per_vertex; ++i)
434  sparsity.add_entries (boundary_dof_boundary_indices[i],
435  boundary_dof_boundary_indices.begin(),
436  boundary_dof_boundary_indices.end());
437  };
438  return;
439  }
440 
441  const types::global_dof_index n_dofs = dof.n_dofs();
442  (void)n_dofs;
443 
444  AssertDimension (dof_to_boundary_mapping.size(), n_dofs);
445  Assert (boundary_ids.find(numbers::internal_face_boundary_id) == boundary_ids.end(),
446  typename DoFHandlerType::ExcInvalidBoundaryIndicator());
447  Assert (sparsity.n_rows() == dof.n_boundary_dofs (boundary_ids),
448  ExcDimensionMismatch (sparsity.n_rows(), dof.n_boundary_dofs (boundary_ids)));
449  Assert (sparsity.n_cols() == dof.n_boundary_dofs (boundary_ids),
450  ExcDimensionMismatch (sparsity.n_cols(), dof.n_boundary_dofs (boundary_ids)));
451 #ifdef DEBUG
452  if (sparsity.n_rows() != 0)
453  {
454  types::global_dof_index max_element = 0;
455  for (std::vector<types::global_dof_index>::const_iterator i=dof_to_boundary_mapping.begin();
456  i!=dof_to_boundary_mapping.end(); ++i)
457  if ((*i != DoFHandlerType::invalid_dof_index) &&
458  (*i > max_element))
459  max_element = *i;
460  AssertDimension (max_element, sparsity.n_rows()-1);
461  };
462 #endif
463 
464  std::vector<types::global_dof_index> dofs_on_this_face;
465  dofs_on_this_face.reserve (max_dofs_per_face(dof));
466  typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(),
467  endc = dof.end();
468  for (; cell!=endc; ++cell)
469  for (unsigned int f=0; f<GeometryInfo<DoFHandlerType::dimension>::faces_per_cell; ++f)
470  if (boundary_ids.find(cell->face(f)->boundary_id()) !=
471  boundary_ids.end())
472  {
473  const unsigned int dofs_per_face = cell->get_fe().dofs_per_face;
474  dofs_on_this_face.resize (dofs_per_face);
475  cell->face(f)->get_dof_indices (dofs_on_this_face,
476  cell->active_fe_index());
477 
478  // make sparsity pattern for this cell
479  for (unsigned int i=0; i<dofs_per_face; ++i)
480  for (unsigned int j=0; j<dofs_per_face; ++j)
481  sparsity.add (dof_to_boundary_mapping[dofs_on_this_face[i]],
482  dof_to_boundary_mapping[dofs_on_this_face[j]]);
483  }
484  }
485 
486 
487 
488  template <typename DoFHandlerType, typename SparsityPatternType>
489  void
490  make_flux_sparsity_pattern (const DoFHandlerType &dof,
491  SparsityPatternType &sparsity,
492  const ConstraintMatrix &constraints,
493  const bool keep_constrained_dofs,
494  const types::subdomain_id subdomain_id)
495 
496  // TODO: QA: reduce the indentation level of this method..., Maier 2012
497 
498  {
499  const types::global_dof_index n_dofs = dof.n_dofs();
500  (void)n_dofs;
501 
502  AssertDimension (sparsity.n_rows(), n_dofs);
503  AssertDimension (sparsity.n_cols(), n_dofs);
504 
505  // If we have a distributed::Triangulation only allow locally_owned
506  // subdomain. Not setting a subdomain is also okay, because we skip
507  // ghost cells in the loop below.
508  Assert (
509  (dof.get_triangulation().locally_owned_subdomain() == numbers::invalid_subdomain_id)
510  ||
511  (subdomain_id == numbers::invalid_subdomain_id)
512  ||
513  (subdomain_id == dof.get_triangulation().locally_owned_subdomain()),
514  ExcMessage ("For parallel::distributed::Triangulation objects and "
515  "associated DoF handler objects, asking for any subdomain other "
516  "than the locally owned one does not make sense."));
517 
518  std::vector<types::global_dof_index> dofs_on_this_cell;
519  std::vector<types::global_dof_index> dofs_on_other_cell;
520  dofs_on_this_cell.reserve (max_dofs_per_cell(dof));
521  dofs_on_other_cell.reserve (max_dofs_per_cell(dof));
522  typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(),
523  endc = dof.end();
524 
525  // TODO: in an old implementation, we used user flags before to tag
526  // faces that were already touched. this way, we could reduce the work
527  // a little bit. now, we instead add only data from one side. this
528  // should be OK, but we need to actually verify it.
529 
530  // In case we work with a distributed sparsity pattern of Trilinos
531  // type, we only have to do the work if the current cell is owned by
532  // the calling processor. Otherwise, just continue.
533  for (; cell!=endc; ++cell)
534  if (((subdomain_id == numbers::invalid_subdomain_id)
535  ||
536  (subdomain_id == cell->subdomain_id()))
537  &&
538  cell->is_locally_owned())
539  {
540  const unsigned int n_dofs_on_this_cell = cell->get_fe().dofs_per_cell;
541  dofs_on_this_cell.resize (n_dofs_on_this_cell);
542  cell->get_dof_indices (dofs_on_this_cell);
543 
544  // make sparsity pattern for this cell. if no constraints pattern
545  // was given, then the following call acts as if simply no
546  // constraints existed
547  constraints.add_entries_local_to_global (dofs_on_this_cell,
548  sparsity,
549  keep_constrained_dofs);
550 
551  for (unsigned int face = 0;
552  face < GeometryInfo<DoFHandlerType::dimension>::faces_per_cell;
553  ++face)
554  {
555  typename DoFHandlerType::face_iterator cell_face = cell->face(face);
556  if (! cell->at_boundary(face) )
557  {
558  typename DoFHandlerType::level_cell_iterator neighbor = cell->neighbor(face);
559 
560  // in 1d, we do not need to worry whether the neighbor
561  // might have children and then loop over those children.
562  // rather, we may as well go straight to to cell behind
563  // this particular cell's most terminal child
564  if (DoFHandlerType::dimension==1)
565  while (neighbor->has_children())
566  neighbor = neighbor->child(face==0 ? 1 : 0);
567 
568  if (neighbor->has_children())
569  {
570  for (unsigned int sub_nr = 0;
571  sub_nr != cell_face->number_of_children();
572  ++sub_nr)
573  {
574  const typename DoFHandlerType::level_cell_iterator
575  sub_neighbor
576  = cell->neighbor_child_on_subface (face, sub_nr);
577 
578  const unsigned int n_dofs_on_neighbor
579  = sub_neighbor->get_fe().dofs_per_cell;
580  dofs_on_other_cell.resize (n_dofs_on_neighbor);
581  sub_neighbor->get_dof_indices (dofs_on_other_cell);
582 
583  constraints.add_entries_local_to_global
584  (dofs_on_this_cell, dofs_on_other_cell,
585  sparsity, keep_constrained_dofs);
586  constraints.add_entries_local_to_global
587  (dofs_on_other_cell, dofs_on_this_cell,
588  sparsity, keep_constrained_dofs);
589  // only need to add this when the neighbor is not
590  // owned by the current processor, otherwise we add
591  // the entries for the neighbor there
592  if (sub_neighbor->subdomain_id() != cell->subdomain_id())
593  constraints.add_entries_local_to_global
594  (dofs_on_other_cell, sparsity, keep_constrained_dofs);
595  }
596  }
597  else
598  {
599  // Refinement edges are taken care of by coarser
600  // cells
601  if (cell->neighbor_is_coarser(face) &&
602  neighbor->subdomain_id() == cell->subdomain_id())
603  continue;
604 
605  const unsigned int n_dofs_on_neighbor
606  = neighbor->get_fe().dofs_per_cell;
607  dofs_on_other_cell.resize (n_dofs_on_neighbor);
608 
609  neighbor->get_dof_indices (dofs_on_other_cell);
610 
611  constraints.add_entries_local_to_global
612  (dofs_on_this_cell, dofs_on_other_cell,
613  sparsity, keep_constrained_dofs);
614 
615  // only need to add these in case the neighbor cell
616  // is not locally owned - otherwise, we touch each
617  // face twice and hence put the indices the other way
618  // around
619  if (!cell->neighbor(face)->active()
620  ||
621  (neighbor->subdomain_id() != cell->subdomain_id()))
622  {
623  constraints.add_entries_local_to_global
624  (dofs_on_other_cell, dofs_on_this_cell,
625  sparsity, keep_constrained_dofs);
626  if (neighbor->subdomain_id() != cell->subdomain_id())
627  constraints.add_entries_local_to_global
628  (dofs_on_other_cell, sparsity, keep_constrained_dofs);
629  }
630  }
631  }
632  }
633  }
634  }
635 
636 
637 
638  template <typename DoFHandlerType, typename SparsityPatternType>
639  void
640  make_flux_sparsity_pattern (const DoFHandlerType &dof,
641  SparsityPatternType &sparsity)
642  {
643  ConstraintMatrix constraints;
644  make_flux_sparsity_pattern (dof, sparsity, constraints);
645  }
646 
647  template <int dim, int spacedim>
650  const Table<2,Coupling> &component_couplings)
651  {
652  Assert(component_couplings.n_rows() == fe.n_components(),
653  ExcDimensionMismatch(component_couplings.n_rows(),
654  fe.n_components()));
655  Assert(component_couplings.n_cols() == fe.n_components(),
656  ExcDimensionMismatch(component_couplings.n_cols(),
657  fe.n_components()));
658 
659  const unsigned int n_dofs = fe.dofs_per_cell;
660 
661  Table<2,Coupling> dof_couplings (n_dofs, n_dofs);
662 
663  for (unsigned int i=0; i<n_dofs; ++i)
664  {
665  const unsigned int ii
666  = (fe.is_primitive(i) ?
667  fe.system_to_component_index(i).first
668  :
670  );
671  Assert (ii < fe.n_components(), ExcInternalError());
672 
673  for (unsigned int j=0; j<n_dofs; ++j)
674  {
675  const unsigned int jj
676  = (fe.is_primitive(j) ?
677  fe.system_to_component_index(j).first
678  :
680  );
681  Assert (jj < fe.n_components(), ExcInternalError());
682 
683  dof_couplings(i,j) = component_couplings(ii,jj);
684  }
685  }
686  return dof_couplings;
687  }
688 
689 
690 
691  template <int dim, int spacedim>
692  std::vector<Table<2,Coupling> >
695  const Table<2,Coupling> &component_couplings)
696  {
697  std::vector<Table<2,Coupling> > return_value (fe.size());
698  for (unsigned int i=0; i<fe.size(); ++i)
699  return_value[i]
700  = dof_couplings_from_component_couplings(fe[i], component_couplings);
701 
702  return return_value;
703  }
704 
705 
706 
707  namespace internal
708  {
709  namespace
710  {
711 
712  // implementation of the same function in namespace DoFTools for
713  // non-hp DoFHandlers
714  template <typename DoFHandlerType, typename SparsityPatternType>
715  void
716  make_flux_sparsity_pattern (const DoFHandlerType &dof,
717  SparsityPatternType &sparsity,
718  const Table<2,Coupling> &int_mask,
719  const Table<2,Coupling> &flux_mask)
720  {
722 
723  std::vector<types::global_dof_index> dofs_on_this_cell(fe.dofs_per_cell);
724  std::vector<types::global_dof_index> dofs_on_other_cell(fe.dofs_per_cell);
725 
726  const Table<2,Coupling>
727  int_dof_mask = dof_couplings_from_component_couplings(fe, int_mask),
728  flux_dof_mask = dof_couplings_from_component_couplings(fe, flux_mask);
729 
730  Table<2,bool> support_on_face(fe.dofs_per_cell,
732  for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
733  for (unsigned int f=0; f<GeometryInfo<DoFHandlerType::dimension>::faces_per_cell; ++f)
734  support_on_face(i,f) = fe.has_support_on_face(i,f);
735 
736  typename DoFHandlerType::active_cell_iterator cell = dof.begin_active(),
737  endc = dof.end();
738  for (; cell!=endc; ++cell)
739  if (cell->is_locally_owned())
740  {
741  cell->get_dof_indices (dofs_on_this_cell);
742  // make sparsity pattern for this cell
743  for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
744  for (unsigned int j=0; j<fe.dofs_per_cell; ++j)
745  if (int_dof_mask(i,j) != none)
746  sparsity.add (dofs_on_this_cell[i],
747  dofs_on_this_cell[j]);
748 
749  // Loop over all interior neighbors
750  for (unsigned int face = 0;
751  face < GeometryInfo<DoFHandlerType::dimension>::faces_per_cell;
752  ++face)
753  {
754  const typename DoFHandlerType::face_iterator
755  cell_face = cell->face(face);
756  if (cell_face->user_flag_set ())
757  continue;
758 
759  if (cell->at_boundary (face) )
760  {
761  for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
762  {
763  const bool i_non_zero_i = support_on_face (i, face);
764  for (unsigned int j=0; j<fe.dofs_per_cell; ++j)
765  {
766  const bool j_non_zero_i = support_on_face (j, face);
767 
768  if ((flux_dof_mask(i,j) == always)
769  ||
770  (flux_dof_mask(i,j) == nonzero
771  &&
772  i_non_zero_i
773  &&
774  j_non_zero_i))
775  sparsity.add (dofs_on_this_cell[i],
776  dofs_on_this_cell[j]);
777  }
778  }
779  }
780  else
781  {
782  typename DoFHandlerType::level_cell_iterator
783  neighbor = cell->neighbor(face);
784  // Refinement edges are taken care of by coarser
785  // cells
786  if (cell->neighbor_is_coarser(face))
787  continue;
788 
789  const unsigned int
790  neighbor_face = cell->neighbor_of_neighbor(face);
791 
792  if (cell_face->has_children())
793  {
794  for (unsigned int sub_nr = 0;
795  sub_nr != cell_face->n_children();
796  ++sub_nr)
797  {
798  const typename DoFHandlerType::level_cell_iterator
799  sub_neighbor
800  = cell->neighbor_child_on_subface (face, sub_nr);
801 
802  sub_neighbor->get_dof_indices (dofs_on_other_cell);
803  for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
804  {
805  const bool i_non_zero_i = support_on_face (i, face);
806  const bool i_non_zero_e = support_on_face (i, neighbor_face);
807  for (unsigned int j=0; j<fe.dofs_per_cell; ++j)
808  {
809  const bool j_non_zero_i = support_on_face (j, face);
810  const bool j_non_zero_e = support_on_face (j, neighbor_face);
811 
812  if (flux_dof_mask(i,j) == always)
813  {
814  sparsity.add (dofs_on_this_cell[i],
815  dofs_on_other_cell[j]);
816  sparsity.add (dofs_on_other_cell[i],
817  dofs_on_this_cell[j]);
818  sparsity.add (dofs_on_this_cell[i],
819  dofs_on_this_cell[j]);
820  sparsity.add (dofs_on_other_cell[i],
821  dofs_on_other_cell[j]);
822  }
823  else if (flux_dof_mask(i,j) == nonzero)
824  {
825  if (i_non_zero_i && j_non_zero_e)
826  sparsity.add (dofs_on_this_cell[i],
827  dofs_on_other_cell[j]);
828  if (i_non_zero_e && j_non_zero_i)
829  sparsity.add (dofs_on_other_cell[i],
830  dofs_on_this_cell[j]);
831  if (i_non_zero_i && j_non_zero_i)
832  sparsity.add (dofs_on_this_cell[i],
833  dofs_on_this_cell[j]);
834  if (i_non_zero_e && j_non_zero_e)
835  sparsity.add (dofs_on_other_cell[i],
836  dofs_on_other_cell[j]);
837  }
838 
839  if (flux_dof_mask(j,i) == always)
840  {
841  sparsity.add (dofs_on_this_cell[j],
842  dofs_on_other_cell[i]);
843  sparsity.add (dofs_on_other_cell[j],
844  dofs_on_this_cell[i]);
845  sparsity.add (dofs_on_this_cell[j],
846  dofs_on_this_cell[i]);
847  sparsity.add (dofs_on_other_cell[j],
848  dofs_on_other_cell[i]);
849  }
850  else if (flux_dof_mask(j,i) == nonzero)
851  {
852  if (j_non_zero_i && i_non_zero_e)
853  sparsity.add (dofs_on_this_cell[j],
854  dofs_on_other_cell[i]);
855  if (j_non_zero_e && i_non_zero_i)
856  sparsity.add (dofs_on_other_cell[j],
857  dofs_on_this_cell[i]);
858  if (j_non_zero_i && i_non_zero_i)
859  sparsity.add (dofs_on_this_cell[j],
860  dofs_on_this_cell[i]);
861  if (j_non_zero_e && i_non_zero_e)
862  sparsity.add (dofs_on_other_cell[j],
863  dofs_on_other_cell[i]);
864  }
865  }
866  }
867  sub_neighbor->face(neighbor_face)->set_user_flag ();
868  }
869  }
870  else
871  {
872  neighbor->get_dof_indices (dofs_on_other_cell);
873  for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
874  {
875  const bool i_non_zero_i = support_on_face (i, face);
876  const bool i_non_zero_e = support_on_face (i, neighbor_face);
877  for (unsigned int j=0; j<fe.dofs_per_cell; ++j)
878  {
879  const bool j_non_zero_i = support_on_face (j, face);
880  const bool j_non_zero_e = support_on_face (j, neighbor_face);
881  if (flux_dof_mask(i,j) == always)
882  {
883  sparsity.add (dofs_on_this_cell[i],
884  dofs_on_other_cell[j]);
885  sparsity.add (dofs_on_other_cell[i],
886  dofs_on_this_cell[j]);
887  sparsity.add (dofs_on_this_cell[i],
888  dofs_on_this_cell[j]);
889  sparsity.add (dofs_on_other_cell[i],
890  dofs_on_other_cell[j]);
891  }
892  if (flux_dof_mask(i,j) == nonzero)
893  {
894  if (i_non_zero_i && j_non_zero_e)
895  sparsity.add (dofs_on_this_cell[i],
896  dofs_on_other_cell[j]);
897  if (i_non_zero_e && j_non_zero_i)
898  sparsity.add (dofs_on_other_cell[i],
899  dofs_on_this_cell[j]);
900  if (i_non_zero_i && j_non_zero_i)
901  sparsity.add (dofs_on_this_cell[i],
902  dofs_on_this_cell[j]);
903  if (i_non_zero_e && j_non_zero_e)
904  sparsity.add (dofs_on_other_cell[i],
905  dofs_on_other_cell[j]);
906  }
907 
908  if (flux_dof_mask(j,i) == always)
909  {
910  sparsity.add (dofs_on_this_cell[j],
911  dofs_on_other_cell[i]);
912  sparsity.add (dofs_on_other_cell[j],
913  dofs_on_this_cell[i]);
914  sparsity.add (dofs_on_this_cell[j],
915  dofs_on_this_cell[i]);
916  sparsity.add (dofs_on_other_cell[j],
917  dofs_on_other_cell[i]);
918  }
919  if (flux_dof_mask(j,i) == nonzero)
920  {
921  if (j_non_zero_i && i_non_zero_e)
922  sparsity.add (dofs_on_this_cell[j],
923  dofs_on_other_cell[i]);
924  if (j_non_zero_e && i_non_zero_i)
925  sparsity.add (dofs_on_other_cell[j],
926  dofs_on_this_cell[i]);
927  if (j_non_zero_i && i_non_zero_i)
928  sparsity.add (dofs_on_this_cell[j],
929  dofs_on_this_cell[i]);
930  if (j_non_zero_e && i_non_zero_e)
931  sparsity.add (dofs_on_other_cell[j],
932  dofs_on_other_cell[i]);
933  }
934  }
935  }
936  neighbor->face(neighbor_face)->set_user_flag ();
937  }
938  }
939  }
940  }
941  }
942 
943 
944  // implementation of the same function in namespace DoFTools for
945  // non-hp DoFHandlers
946  template <int dim, int spacedim, typename SparsityPatternType>
947  void
948  make_flux_sparsity_pattern (const ::hp::DoFHandler<dim,spacedim> &dof,
949  SparsityPatternType &sparsity,
950  const Table<2,Coupling> &int_mask,
951  const Table<2,Coupling> &flux_mask)
952  {
953  // while the implementation above is quite optimized and caches a
954  // lot of data (see e.g. the int/flux_dof_mask tables), this is no
955  // longer practical for the hp version since we would have to have
956  // it for all combinations of elements in the hp::FECollection.
957  // consequently, the implementation here is simpler and probably
958  // less efficient but at least readable...
959 
960  const ::hp::FECollection<dim,spacedim> &fe = dof.get_fe();
961 
962  std::vector<types::global_dof_index> dofs_on_this_cell(DoFTools::max_dofs_per_cell(dof));
963  std::vector<types::global_dof_index> dofs_on_other_cell(DoFTools::max_dofs_per_cell(dof));
964 
965  const std::vector<Table<2,Coupling> >
966  int_dof_mask
968 
969  typename ::hp::DoFHandler<dim,spacedim>::active_cell_iterator
970  cell = dof.begin_active(),
971  endc = dof.end();
972  for (; cell!=endc; ++cell)
973  {
974  dofs_on_this_cell.resize (cell->get_fe().dofs_per_cell);
975  cell->get_dof_indices (dofs_on_this_cell);
976 
977  // make sparsity pattern for this cell
978  for (unsigned int i=0; i<cell->get_fe().dofs_per_cell; ++i)
979  for (unsigned int j=0; j<cell->get_fe().dofs_per_cell; ++j)
980  if (int_dof_mask[cell->active_fe_index()](i,j) != none)
981  sparsity.add (dofs_on_this_cell[i],
982  dofs_on_this_cell[j]);
983 
984  // Loop over all interior neighbors
985  for (unsigned int face = 0;
986  face < GeometryInfo<dim>::faces_per_cell;
987  ++face)
988  {
989  const typename ::hp::DoFHandler<dim,spacedim>::face_iterator
990  cell_face = cell->face(face);
991  if (cell_face->user_flag_set ())
992  continue;
993 
994  if (cell->at_boundary (face) )
995  {
996  for (unsigned int i=0; i<cell->get_fe().dofs_per_cell; ++i)
997  for (unsigned int j=0; j<cell->get_fe().dofs_per_cell; ++j)
998  if ((flux_mask(cell->get_fe().system_to_component_index(i).first,
999  cell->get_fe().system_to_component_index(j).first)
1000  == always)
1001  ||
1002  (flux_mask(cell->get_fe().system_to_component_index(i).first,
1003  cell->get_fe().system_to_component_index(j).first)
1004  == nonzero))
1005  sparsity.add (dofs_on_this_cell[i],
1006  dofs_on_this_cell[j]);
1007  }
1008  else
1009  {
1010  typename ::hp::DoFHandler<dim,spacedim>::level_cell_iterator
1011  neighbor = cell->neighbor(face);
1012 
1013  // Refinement edges are taken care of by coarser cells
1014  if (cell->neighbor_is_coarser(face))
1015  continue;
1016 
1017  const unsigned int
1018  neighbor_face = cell->neighbor_of_neighbor(face);
1019 
1020  if (cell_face->has_children())
1021  {
1022  for (unsigned int sub_nr = 0;
1023  sub_nr != cell_face->n_children();
1024  ++sub_nr)
1025  {
1026  const typename ::hp::DoFHandler<dim,spacedim>::level_cell_iterator
1027  sub_neighbor
1028  = cell->neighbor_child_on_subface (face, sub_nr);
1029 
1030  dofs_on_other_cell.resize (sub_neighbor->get_fe().dofs_per_cell);
1031  sub_neighbor->get_dof_indices (dofs_on_other_cell);
1032  for (unsigned int i=0; i<cell->get_fe().dofs_per_cell; ++i)
1033  {
1034  for (unsigned int j=0; j<sub_neighbor->get_fe().dofs_per_cell;
1035  ++j)
1036  {
1037  if ((flux_mask(cell->get_fe().system_to_component_index(i).first,
1038  sub_neighbor->get_fe().system_to_component_index(j).first)
1039  == always)
1040  ||
1041  (flux_mask(cell->get_fe().system_to_component_index(i).first,
1042  sub_neighbor->get_fe().system_to_component_index(j).first)
1043  == nonzero))
1044  {
1045  sparsity.add (dofs_on_this_cell[i],
1046  dofs_on_other_cell[j]);
1047  sparsity.add (dofs_on_other_cell[i],
1048  dofs_on_this_cell[j]);
1049  sparsity.add (dofs_on_this_cell[i],
1050  dofs_on_this_cell[j]);
1051  sparsity.add (dofs_on_other_cell[i],
1052  dofs_on_other_cell[j]);
1053  }
1054 
1055  if ((flux_mask(sub_neighbor->get_fe().system_to_component_index(j).first,
1056  cell->get_fe().system_to_component_index(i).first)
1057  == always)
1058  ||
1059  (flux_mask(sub_neighbor->get_fe().system_to_component_index(j).first,
1060  cell->get_fe().system_to_component_index(i).first)
1061  == nonzero))
1062  {
1063  sparsity.add (dofs_on_this_cell[j],
1064  dofs_on_other_cell[i]);
1065  sparsity.add (dofs_on_other_cell[j],
1066  dofs_on_this_cell[i]);
1067  sparsity.add (dofs_on_this_cell[j],
1068  dofs_on_this_cell[i]);
1069  sparsity.add (dofs_on_other_cell[j],
1070  dofs_on_other_cell[i]);
1071  }
1072  }
1073  }
1074  sub_neighbor->face(neighbor_face)->set_user_flag ();
1075  }
1076  }
1077  else
1078  {
1079  dofs_on_other_cell.resize (neighbor->get_fe().dofs_per_cell);
1080  neighbor->get_dof_indices (dofs_on_other_cell);
1081  for (unsigned int i=0; i<cell->get_fe().dofs_per_cell; ++i)
1082  {
1083  for (unsigned int j=0; j<neighbor->get_fe().dofs_per_cell; ++j)
1084  {
1085  if ((flux_mask(cell->get_fe().system_to_component_index(i).first,
1086  neighbor->get_fe().system_to_component_index(j).first)
1087  == always)
1088  ||
1089  (flux_mask(cell->get_fe().system_to_component_index(i).first,
1090  neighbor->get_fe().system_to_component_index(j).first)
1091  == nonzero))
1092  {
1093  sparsity.add (dofs_on_this_cell[i],
1094  dofs_on_other_cell[j]);
1095  sparsity.add (dofs_on_other_cell[i],
1096  dofs_on_this_cell[j]);
1097  sparsity.add (dofs_on_this_cell[i],
1098  dofs_on_this_cell[j]);
1099  sparsity.add (dofs_on_other_cell[i],
1100  dofs_on_other_cell[j]);
1101  }
1102 
1103  if ((flux_mask(neighbor->get_fe().system_to_component_index(j).first,
1104  cell->get_fe().system_to_component_index(i).first)
1105  == always)
1106  ||
1107  (flux_mask(neighbor->get_fe().system_to_component_index(j).first,
1108  cell->get_fe().system_to_component_index(i).first)
1109  == nonzero))
1110  {
1111  sparsity.add (dofs_on_this_cell[j],
1112  dofs_on_other_cell[i]);
1113  sparsity.add (dofs_on_other_cell[j],
1114  dofs_on_this_cell[i]);
1115  sparsity.add (dofs_on_this_cell[j],
1116  dofs_on_this_cell[i]);
1117  sparsity.add (dofs_on_other_cell[j],
1118  dofs_on_other_cell[i]);
1119  }
1120  }
1121  }
1122  neighbor->face(neighbor_face)->set_user_flag ();
1123  }
1124  }
1125  }
1126  }
1127  }
1128  }
1129 
1130  }
1131 
1132 
1133 
1134 
1135  template <typename DoFHandlerType, typename SparsityPatternType>
1136  void
1137  make_flux_sparsity_pattern (const DoFHandlerType &dof,
1138  SparsityPatternType &sparsity,
1139  const Table<2,Coupling> &int_mask,
1140  const Table<2,Coupling> &flux_mask)
1141  {
1142  // do the error checking and frame code here, and then pass on to more
1143  // specialized functions in the internal namespace
1144  const types::global_dof_index n_dofs = dof.n_dofs();
1145  (void)n_dofs;
1146  const unsigned int n_comp = dof.get_fe().n_components();
1147  (void)n_comp;
1148 
1149  Assert (sparsity.n_rows() == n_dofs,
1150  ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
1151  Assert (sparsity.n_cols() == n_dofs,
1152  ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
1153  Assert (int_mask.n_rows() == n_comp,
1154  ExcDimensionMismatch (int_mask.n_rows(), n_comp));
1155  Assert (int_mask.n_cols() == n_comp,
1156  ExcDimensionMismatch (int_mask.n_cols(), n_comp));
1157  Assert (flux_mask.n_rows() == n_comp,
1158  ExcDimensionMismatch (flux_mask.n_rows(), n_comp));
1159  Assert (flux_mask.n_cols() == n_comp,
1160  ExcDimensionMismatch (flux_mask.n_cols(), n_comp));
1161 
1162  // Clear user flags because we will need them. But first we save them
1163  // and make sure that we restore them later such that at the end of
1164  // this function the Triangulation will be in the same state as it was
1165  // at the beginning of this function.
1166  std::vector<bool> user_flags;
1167  dof.get_triangulation().save_user_flags(user_flags);
1169  (dof.get_triangulation()).clear_user_flags ();
1170 
1171  internal::make_flux_sparsity_pattern (dof, sparsity,
1172  int_mask, flux_mask);
1173 
1174  // finally restore the user flags
1176  (dof.get_triangulation()).load_user_flags(user_flags);
1177  }
1178 
1179 
1180 } // end of namespace DoFTools
1181 
1182 
1183 // --------------------------------------------------- explicit instantiations
1184 
1185 #include "dof_tools_sparsity.inst"
1186 
1187 
1188 
1189 DEAL_II_NAMESPACE_CLOSE
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1052
const types::subdomain_id invalid_subdomain_id
Definition: types.h:239
void make_boundary_sparsity_pattern(const DoFHandlerType &dof, const std::vector< types::global_dof_index > &dof_to_boundary_mapping, SparsityPatternType &sparsity_pattern)
::ExceptionBase & ExcMessage(std::string arg1)
unsigned int size() const
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
Definition: fe.cc:1080
bool is_primitive(const unsigned int i) const
Definition: fe.h:2937
Table< 2, Coupling > dof_couplings_from_component_couplings(const FiniteElement< dim, spacedim > &fe, const Table< 2, Coupling > &component_couplings)
unsigned int global_dof_index
Definition: types.h:88
#define Assert(cond, exc)
Definition: exceptions.h:294
void make_flux_sparsity_pattern(const DoFHandlerType &dof_handler, SparsityPatternType &sparsity_pattern)
const ComponentMask & get_nonzero_components(const unsigned int i) const
Definition: fe.h:2915
std::map< types::boundary_id, const Function< dim, Number > * > type
Definition: function_map.h:81
unsigned int n_components() const
unsigned int subdomain_id
Definition: types.h:42
unsigned int first_selected_component(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
const unsigned int dofs_per_cell
Definition: fe_base.h:283
void make_sparsity_pattern(const DoFHandlerType &dof_handler, SparsityPatternType &sparsity_pattern, const ConstraintMatrix &constraints=ConstraintMatrix(), const bool keep_constrained_dofs=true, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id)
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
Definition: fe.h:2742
Definition: table.h:33
const types::boundary_id internal_face_boundary_id
Definition: types.h:210
void add_entries_local_to_global(const std::vector< size_type > &local_dof_indices, SparsityPatternType &sparsity_pattern, const bool keep_constrained_entries=true, const Table< 2, bool > &dof_mask=default_empty_table) const
std::list< std::pair< typename MeshType::cell_iterator, typename MeshType::cell_iterator > > get_finest_common_cells(const MeshType &mesh_1, const MeshType &mesh_2)
Definition: grid_tools.cc:2105