Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
solution_transfer.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 2018 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 #ifndef dealii_solution_transfer_h
17 # define dealii_solution_transfer_h
18 
19 
20 /*---------------------------- solutiontransfer.h ----------------------*/
21 
22 
23 # include <deal.II/base/config.h>
24 
25 # include <deal.II/base/exceptions.h>
26 # include <deal.II/base/smartpointer.h>
27 
28 # include <deal.II/dofs/dof_handler.h>
29 
30 # include <deal.II/lac/vector.h>
31 
32 # include <vector>
33 
34 DEAL_II_NAMESPACE_OPEN
35 
332 template <int dim,
333  typename VectorType = Vector<double>,
334  typename DoFHandlerType = DoFHandler<dim>>
336 {
337 # ifndef DEAL_II_MSVC
338  static_assert(dim == DoFHandlerType::dimension,
339  "The dimension explicitly provided as a template "
340  "argument, and the dimension of the DoFHandlerType "
341  "template argument must match.");
342 # endif
343 public:
347  SolutionTransfer(const DoFHandlerType &dof);
348 
353 
358  void
359  clear();
360 
366  void
368 
376  void
377  prepare_for_coarsening_and_refinement(const std::vector<VectorType> &all_in);
378 
383  void
384  prepare_for_coarsening_and_refinement(const VectorType &in);
385 
397  void
398  refine_interpolate(const VectorType &in, VectorType &out) const;
399 
419  void
420  interpolate(const std::vector<VectorType> &all_in,
421  std::vector<VectorType> & all_out) const;
422 
432  void
433  interpolate(const VectorType &in, VectorType &out) const;
434 
439  std::size_t
440  memory_consumption() const;
441 
446  "You are attempting an operation for which this object is "
447  "not prepared. This may be because you either did not call "
448  "one of the prepare_*() functions at all, or because you "
449  "called the wrong one for the operation you are currently "
450  "attempting.");
451 
457  "You are attempting to call one of the prepare_*() functions "
458  "of this object to prepare it for an operation for which it "
459  "is already prepared. Specifically, the object was "
460  "previously prepared for pure refinement.");
461 
467  "You are attempting to call one of the prepare_*() functions "
468  "of this object to prepare it for an operation for which it "
469  "is already prepared. Specifically, the object was "
470  "previously prepared for both coarsening and refinement.");
471 
472 private:
476  SmartPointer<const DoFHandlerType,
479 
484 
491  {
504  };
505 
510 
511 
517  std::vector<std::vector<types::global_dof_index>> indices_on_cell;
518 
531  {
532  Pointerstruct()
533  : indices_ptr(nullptr)
534  , dof_values_ptr(nullptr)
535  , active_fe_index(0)
536  {}
537  Pointerstruct(std::vector<types::global_dof_index> *indices_ptr_in,
538  const unsigned int active_fe_index_in = 0)
539  : indices_ptr(indices_ptr_in)
540  , dof_values_ptr(nullptr)
541  , active_fe_index(active_fe_index_in)
542  {}
544  std::vector<Vector<typename VectorType::value_type>> *dof_values_ptr_in,
545  const unsigned int active_fe_index_in = 0)
546  : indices_ptr(nullptr)
547  , dof_values_ptr(dof_values_ptr_in)
548  , active_fe_index(active_fe_index_in)
549  {}
550  std::size_t
551  memory_consumption() const;
552 
553  std::vector<types::global_dof_index> * indices_ptr;
554  std::vector<Vector<typename VectorType::value_type>> *dof_values_ptr;
555  unsigned int active_fe_index;
556  };
557 
564  std::map<std::pair<unsigned int, unsigned int>, Pointerstruct> cell_map;
565 
570  std::vector<std::vector<Vector<typename VectorType::value_type>>>
572 };
573 
574 
575 DEAL_II_NAMESPACE_CLOSE
576 
577 #endif // dealii_solutiontransfer_h
578 /*---------------------------- solutiontransfer.h ---------------------------*/
types::global_dof_index n_dofs_old
SmartPointer< const DoFHandlerType, SolutionTransfer< dim, VectorType, DoFHandlerType > > dof_handler
static ::ExceptionBase & ExcAlreadyPrepForRef()
static ::ExceptionBase & ExcNotPrepared()
std::map< std::pair< unsigned int, unsigned int >, Pointerstruct > cell_map
std::vector< std::vector< types::global_dof_index > > indices_on_cell
void prepare_for_pure_refinement()
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
static ::ExceptionBase & ExcAlreadyPrepForCoarseAndRef()
SolutionTransfer(const DoFHandlerType &dof)
unsigned int global_dof_index
Definition: types.h:89
void prepare_for_coarsening_and_refinement(const std::vector< VectorType > &all_in)
void refine_interpolate(const VectorType &in, VectorType &out) const
std::vector< std::vector< Vector< typename VectorType::value_type > > > dof_values_on_cell
std::size_t memory_consumption() const
PreparationState prepared_for
void interpolate(const std::vector< VectorType > &all_in, std::vector< VectorType > &all_out) const