deal.II version GIT relicensing-2865-ga6be64acbe 2025-03-19 21:10:00+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
solution_transfer.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1999 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_solution_transfer_h
16#define dealii_solution_transfer_h
17
18
19/*---------------------------- solutiontransfer.h ----------------------*/
20
21
22#include <deal.II/base/config.h>
23
26
28
29#include <deal.II/lac/vector.h>
30
31#include <boost/range/iterator_range_core.hpp>
32
33#include <vector>
34
36
37
303template <int dim, typename VectorType = Vector<double>, int spacedim = dim>
305{
306public:
320 const bool average_values = false);
321
325 ~SolutionTransfer() = default;
326
333 void
335 const std::vector<const VectorType *> &all_in);
336
340 void
341 prepare_for_coarsening_and_refinement(const std::vector<VectorType> &all_in);
342
347 void
348 prepare_for_coarsening_and_refinement(const VectorType &in);
349
357 void
358 interpolate(std::vector<VectorType *> &all_out);
359
363 void
364 interpolate(std::vector<VectorType> &all_out);
365
375 void
376 interpolate(VectorType &out);
377
384 void
385 prepare_for_serialization(const VectorType &in);
386
390 void
391 prepare_for_serialization(const std::vector<const VectorType *> &all_in);
392
399 void
400 deserialize(VectorType &in);
401
402
406 void
407 deserialize(std::vector<VectorType *> &all_in);
408
413 DEAL_II_DEPRECATED_EARLY void
414 clear();
415
416private:
423
427 const bool average_values;
428
433 std::vector<const VectorType *> input_vectors;
434
439 unsigned int handle;
440
446 std::vector<char>
449 const CellStatus status);
450
456 void
459 const CellStatus status,
460 const boost::iterator_range<std::vector<char>::const_iterator> &data_range,
461 std::vector<VectorType *> &all_out,
462 VectorType &valence);
463
464
469 void
471};
472
473
474
475namespace Legacy
476{
782 template <int dim, typename VectorType = Vector<double>, int spacedim = dim>
784 {
785 public:
790
795
800 void
801 clear();
802
808 void
810
818 void
820 const std::vector<VectorType> &all_in);
821
826 void
827 prepare_for_coarsening_and_refinement(const VectorType &in);
828
840 void
841 refine_interpolate(const VectorType &in, VectorType &out) const;
842
862 void
863 interpolate(const std::vector<VectorType> &all_in,
864 std::vector<VectorType> &all_out) const;
865
875 void
876 interpolate(const VectorType &in, VectorType &out) const;
877
882 std::size_t
883 memory_consumption() const;
884
890 "You are attempting an operation for which this object is "
891 "not prepared. This may be because you either did not call "
892 "one of the prepare_*() functions at all, or because you "
893 "called the wrong one for the operation you are currently "
894 "attempting.");
895
901 "You are attempting to call one of the prepare_*() functions "
902 "of this object to prepare it for an operation for which it "
903 "is already prepared. Specifically, the object was "
904 "previously prepared for pure refinement.");
905
911 "You are attempting to call one of the prepare_*() functions "
912 "of this object to prepare it for an operation for which it "
913 "is already prepared. Specifically, the object was "
914 "previously prepared for both coarsening and refinement.");
915
916 private:
923
928
949
954
955
961 std::vector<std::vector<types::global_dof_index>> indices_on_cell;
962
975 {
977 : indices_ptr(nullptr)
978 , dof_values_ptr(nullptr)
979 , active_fe_index(0)
980 {}
981 Pointerstruct(std::vector<types::global_dof_index> *indices_ptr_in,
982 const unsigned int active_fe_index_in = 0)
983 : indices_ptr(indices_ptr_in)
984 , dof_values_ptr(nullptr)
985 , active_fe_index(active_fe_index_in)
986 {}
988 std::vector<Vector<typename VectorType::value_type>> *dof_values_ptr_in,
989 const unsigned int active_fe_index_in = 0)
990 : indices_ptr(nullptr)
991 , dof_values_ptr(dof_values_ptr_in)
992 , active_fe_index(active_fe_index_in)
993 {}
994 std::size_t
995 memory_consumption() const;
996
997 std::vector<types::global_dof_index> *indices_ptr;
998 std::vector<Vector<typename VectorType::value_type>> *dof_values_ptr;
999 unsigned int active_fe_index;
1000 };
1001
1008 std::map<std::pair<unsigned int, unsigned int>, Pointerstruct> cell_map;
1009
1015 std::vector<std::vector<Vector<typename VectorType::value_type>>>
1017 };
1018} // namespace Legacy
1019
1021
1022#endif
CellStatus
Definition cell_status.h:31
void refine_interpolate(const VectorType &in, VectorType &out) const
void interpolate(const std::vector< VectorType > &all_in, std::vector< VectorType > &all_out) const
std::vector< std::vector< Vector< typename VectorType::value_type > > > dof_values_on_cell
std::size_t memory_consumption() const
std::vector< std::vector< types::global_dof_index > > indices_on_cell
ObserverPointer< const DoFHandler< dim, spacedim >, SolutionTransfer< dim, VectorType, spacedim > > dof_handler
void prepare_for_coarsening_and_refinement(const std::vector< VectorType > &all_in)
types::global_dof_index n_dofs_old
std::map< std::pair< unsigned int, unsigned int >, Pointerstruct > cell_map
std::vector< const VectorType * > input_vectors
std::vector< char > pack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellStatus status)
void unpack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellStatus status, const boost::iterator_range< std::vector< char >::const_iterator > &data_range, std::vector< VectorType * > &all_out, VectorType &valence)
ObserverPointer< const DoFHandler< dim, spacedim >, SolutionTransfer< dim, VectorType, spacedim > > dof_handler
void deserialize(VectorType &in)
void prepare_for_coarsening_and_refinement(const std::vector< const VectorType * > &all_in)
void prepare_for_serialization(const VectorType &in)
void interpolate(std::vector< VectorType * > &all_out)
~SolutionTransfer()=default
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
static ::ExceptionBase & ExcAlreadyPrepForCoarseAndRef()
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcNotPrepared()
static ::ExceptionBase & ExcAlreadyPrepForRef()
std::vector< types::global_dof_index > * indices_ptr
std::vector< Vector< typename VectorType::value_type > > * dof_values_ptr
Pointerstruct(std::vector< types::global_dof_index > *indices_ptr_in, const unsigned int active_fe_index_in=0)
Pointerstruct(std::vector< Vector< typename VectorType::value_type > > *dof_values_ptr_in, const unsigned int active_fe_index_in=0)