1 // --------------------------------------------------------------------- 2 // 3 // Copyright (C) 2001 - 2020 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_mg_transfer_h 17 #define dealii_mg_transfer_h 18 19 #include <deal.II/base/config.h> 20 21 #include <deal.II/base/mg_level_object.h> 22 23 #include <deal.II/distributed/tria_base.h> 24 25 #include <deal.II/dofs/dof_handler.h> 26 27 #include <deal.II/lac/affine_constraints.h> 28 #include <deal.II/lac/block_sparsity_pattern.h> 29 #include <deal.II/lac/block_vector.h> 30 #include <deal.II/lac/la_parallel_vector.h> 31 #include <deal.II/lac/petsc_sparse_matrix.h> 32 #include <deal.II/lac/petsc_vector.h> 33 #include <deal.II/lac/sparse_matrix.h> 34 #include <deal.II/lac/trilinos_sparse_matrix.h> 35 #include <deal.II/lac/vector_memory.h> 36 37 #include <deal.II/multigrid/mg_base.h> 38 #include <deal.II/multigrid/mg_constrained_dofs.h> 39 40 #include <memory> 41 42 43 DEAL_II_NAMESPACE_OPEN 44 45 46 namespace internal 47 { 48 template <typename VectorType> 49 struct MatrixSelector 50 { 51 using Sparsity = ::dealii::SparsityPattern; 52 using Matrix = ::dealii::SparseMatrix<typename VectorType::value_type>; 53 54 static const bool requires_distributed_sparsity_pattern = false; 55 56 template <typename SparsityPatternType, int dim, int spacedim> 57 static void reinitMatrixSelector58 reinit(Matrix & matrix, 59 Sparsity & sparsity, 60 int level, 61 const SparsityPatternType &sp, 62 const DoFHandler<dim, spacedim> &) 63 { 64 sparsity.copy_from(sp); 65 (void)level; 66 matrix.reinit(sparsity); 67 } 68 }; 69 70 #ifdef DEAL_II_WITH_TRILINOS 71 template <typename Number> 72 struct MatrixSelector<LinearAlgebra::distributed::Vector<Number>> 73 { 74 using Sparsity = ::dealii::TrilinosWrappers::SparsityPattern; 75 using Matrix = ::dealii::TrilinosWrappers::SparseMatrix; 76 77 static const bool requires_distributed_sparsity_pattern = false; 78 79 template <typename SparsityPatternType, int dim, int spacedim> 80 static void 81 reinit(Matrix &matrix, 82 Sparsity &, 83 int level, 84 const SparsityPatternType & sp, 85 const DoFHandler<dim, spacedim> &dh) 86 { 87 const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria = 88 dynamic_cast< 89 const dealii::parallel::TriangulationBase<dim, spacedim> *>( 90 &(dh.get_triangulation())); 91 MPI_Comm communicator = 92 dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF; 93 94 matrix.reinit(dh.locally_owned_mg_dofs(level + 1), 95 dh.locally_owned_mg_dofs(level), 96 sp, 97 communicator, 98 true); 99 } 100 }; 101 102 template <> 103 struct MatrixSelector<dealii::TrilinosWrappers::MPI::Vector> 104 { 105 using Sparsity = ::dealii::TrilinosWrappers::SparsityPattern; 106 using Matrix = ::dealii::TrilinosWrappers::SparseMatrix; 107 108 static const bool requires_distributed_sparsity_pattern = false; 109 110 template <typename SparsityPatternType, int dim, int spacedim> 111 static void 112 reinit(Matrix &matrix, 113 Sparsity &, 114 int level, 115 const SparsityPatternType & sp, 116 const DoFHandler<dim, spacedim> &dh) 117 { 118 const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria = 119 dynamic_cast< 120 const dealii::parallel::TriangulationBase<dim, spacedim> *>( 121 &(dh.get_triangulation())); 122 MPI_Comm communicator = 123 dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF; 124 matrix.reinit(dh.locally_owned_mg_dofs(level + 1), 125 dh.locally_owned_mg_dofs(level), 126 sp, 127 communicator, 128 true); 129 } 130 }; 131 132 # ifdef DEAL_II_WITH_MPI 133 # ifdef DEAL_II_TRILINOS_WITH_TPETRA 134 template <typename Number> 135 struct MatrixSelector<dealii::LinearAlgebra::TpetraWrappers::Vector<Number>> 136 { 137 using Sparsity = ::dealii::TrilinosWrappers::SparsityPattern; 138 using Matrix = ::dealii::TrilinosWrappers::SparseMatrix; 139 140 static const bool requires_distributed_sparsity_pattern = false; 141 142 template <typename SparsityPatternType, int dim, int spacedim> 143 static void 144 reinit(Matrix &matrix, 145 Sparsity &, 146 int level, 147 const SparsityPatternType & sp, 148 const DoFHandler<dim, spacedim> &dh) 149 { 150 const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria = 151 dynamic_cast< 152 const dealii::parallel::TriangulationBase<dim, spacedim> *>( 153 &(dh.get_triangulation())); 154 MPI_Comm communicator = 155 dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF; 156 matrix.reinit(dh.locally_owned_mg_dofs(level + 1), 157 dh.locally_owned_mg_dofs(level), 158 sp, 159 communicator, 160 true); 161 } 162 }; 163 # endif 164 165 template <> 166 struct MatrixSelector<dealii::LinearAlgebra::EpetraWrappers::Vector> 167 { 168 using Sparsity = ::dealii::TrilinosWrappers::SparsityPattern; 169 using Matrix = ::dealii::TrilinosWrappers::SparseMatrix; 170 171 static const bool requires_distributed_sparsity_pattern = false; 172 173 template <typename SparsityPatternType, int dim, int spacedim> 174 static void 175 reinit(Matrix &matrix, 176 Sparsity &, 177 int level, 178 const SparsityPatternType & sp, 179 const DoFHandler<dim, spacedim> &dh) 180 { 181 const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria = 182 dynamic_cast< 183 const dealii::parallel::TriangulationBase<dim, spacedim> *>( 184 &(dh.get_triangulation())); 185 MPI_Comm communicator = 186 dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF; 187 matrix.reinit(dh.locally_owned_mg_dofs(level + 1), 188 dh.locally_owned_mg_dofs(level), 189 sp, 190 communicator, 191 true); 192 } 193 }; 194 # endif 195 196 #else 197 // ! DEAL_II_WITH_TRILINOS 198 template <typename Number> 199 struct MatrixSelector<LinearAlgebra::distributed::Vector<Number>> 200 { 201 using Sparsity = ::dealii::SparsityPattern; 202 using Matrix = ::dealii::SparseMatrix<Number>; 203 204 static const bool requires_distributed_sparsity_pattern = false; 205 206 template <typename SparsityPatternType, int dim, int spacedim> 207 static void 208 reinit(Matrix &, 209 Sparsity &, 210 int, 211 const SparsityPatternType &, 212 const DoFHandler<dim, spacedim> &) 213 { 214 AssertThrow( 215 false, 216 ExcNotImplemented( 217 "ERROR: MGTransferPrebuilt with LinearAlgebra::distributed::Vector currently " 218 "needs deal.II to be configured with Trilinos.")); 219 } 220 }; 221 222 #endif 223 224 #ifdef DEAL_II_WITH_PETSC 225 template <> 226 struct MatrixSelector<dealii::PETScWrappers::MPI::Vector> 227 { 228 using Sparsity = ::dealii::DynamicSparsityPattern; 229 using Matrix = ::dealii::PETScWrappers::MPI::SparseMatrix; 230 231 static const bool requires_distributed_sparsity_pattern = true; 232 233 template <typename SparsityPatternType, int dim, int spacedim> 234 static void 235 reinit(Matrix &matrix, 236 Sparsity &, 237 int level, 238 const SparsityPatternType & sp, 239 const DoFHandler<dim, spacedim> &dh) 240 { 241 const dealii::parallel::TriangulationBase<dim, spacedim> *dist_tria = 242 dynamic_cast< 243 const dealii::parallel::TriangulationBase<dim, spacedim> *>( 244 &(dh.get_triangulation())); 245 MPI_Comm communicator = 246 dist_tria != nullptr ? dist_tria->get_communicator() : MPI_COMM_SELF; 247 // Reinit PETSc matrix 248 matrix.reinit(dh.locally_owned_mg_dofs(level + 1), 249 dh.locally_owned_mg_dofs(level), 250 sp, 251 communicator); 252 } 253 }; 254 #endif 255 } // namespace internal 256 257 /* 258 * MGTransferBase is defined in mg_base.h 259 */ 260 261 /*!@addtogroup mg */ 262 /*@{*/ 263 264 265 266 /** 267 * Implementation of transfer between the global vectors and the multigrid 268 * levels for use in the derived class MGTransferPrebuilt and other classes. 269 */ 270 template <typename VectorType> 271 class MGLevelGlobalTransfer : public MGTransferBase<VectorType> 272 { 273 public: 274 /** 275 * Reset the object to the state it had right after the default constructor. 276 */ 277 void 278 clear(); 279 280 /** 281 * Transfer from a vector on the global grid to vectors defined on each of 282 * the levels separately for the active degrees of freedom. In particular, 283 * for a globally refined mesh only the finest level in @p dst is filled as a 284 * plain copy of @p src. All the other level objects are left untouched. 285 */ 286 template <int dim, class InVector, int spacedim> 287 void 288 copy_to_mg(const DoFHandler<dim, spacedim> &dof_handler, 289 MGLevelObject<VectorType> & dst, 290 const InVector & src) const; 291 292 /** 293 * Transfer from multi-level vector to normal vector. 294 * 295 * Copies data from active portions of an MGVector into the respective 296 * positions of a <tt>Vector<number></tt>. In order to keep the result 297 * consistent, constrained degrees of freedom are set to zero. 298 */ 299 template <int dim, class OutVector, int spacedim> 300 void 301 copy_from_mg(const DoFHandler<dim, spacedim> &dof_handler, 302 OutVector & dst, 303 const MGLevelObject<VectorType> &src) const; 304 305 /** 306 * Add a multi-level vector to a normal vector. 307 * 308 * Works as the previous function, but probably not for continuous elements. 309 */ 310 template <int dim, class OutVector, int spacedim> 311 void 312 copy_from_mg_add(const DoFHandler<dim, spacedim> &dof_handler, 313 OutVector & dst, 314 const MGLevelObject<VectorType> &src) const; 315 316 /** 317 * If this object operates on BlockVector objects, we need to describe how 318 * the individual vector components are mapped to the blocks of a vector. 319 * For example, for a Stokes system, we have dim+1 vector components for 320 * velocity and pressure, but we may want to use block vectors with only two 321 * blocks for all velocities in one block, and the pressure variables in the 322 * other. 323 * 324 * By default, if this function is not called, block vectors have as many 325 * blocks as the finite element has vector components. However, this can be 326 * changed by calling this function with an array that describes how vector 327 * components are to be grouped into blocks. The meaning of the argument is 328 * the same as the one given to the DoFTools::count_dofs_per_component 329 * function. 330 */ 331 void 332 set_component_to_block_map(const std::vector<unsigned int> &map); 333 334 /** 335 * Memory used by this object. 336 */ 337 std::size_t 338 memory_consumption() const; 339 340 /** 341 * Print the copy index fields for debugging purposes. 342 */ 343 void 344 print_indices(std::ostream &os) const; 345 346 protected: 347 /** 348 * Internal function to @p fill copy_indices*. Called by derived classes. 349 */ 350 template <int dim, int spacedim> 351 void 352 fill_and_communicate_copy_indices( 353 const DoFHandler<dim, spacedim> &dof_handler); 354 355 /** 356 * Sizes of the multi-level vectors. 357 */ 358 std::vector<types::global_dof_index> sizes; 359 360 /** 361 * Mapping for the copy_to_mg() and copy_from_mg() functions. Here only 362 * index pairs locally owned 363 * 364 * The data is organized as follows: one vector per level. Each element of 365 * these vectors contains first the global index, then the level index. 366 */ 367 std::vector< 368 std::vector<std::pair<types::global_dof_index, types::global_dof_index>>> 369 copy_indices; 370 371 /** 372 * Additional degrees of freedom for the copy_to_mg() function. These are 373 * the ones where the global degree of freedom is locally owned and the 374 * level degree of freedom is not. 375 * 376 * Organization of the data is like for @p copy_indices_mine. 377 */ 378 std::vector< 379 std::vector<std::pair<types::global_dof_index, types::global_dof_index>>> 380 copy_indices_global_mine; 381 382 /** 383 * Additional degrees of freedom for the copy_from_mg() function. These are 384 * the ones where the level degree of freedom is locally owned and the 385 * global degree of freedom is not. 386 * 387 * Organization of the data is like for @p copy_indices_mine. 388 */ 389 std::vector< 390 std::vector<std::pair<types::global_dof_index, types::global_dof_index>>> 391 copy_indices_level_mine; 392 393 /** 394 * This variable stores whether the copy operation from the global to the 395 * level vector is actually a plain copy to the finest level. This means that 396 * the grid has no adaptive refinement and the numbering on the finest 397 * multigrid level is the same as in the global case. 398 */ 399 bool perform_plain_copy; 400 401 /** 402 * The vector that stores what has been given to the 403 * set_component_to_block_map() function. 404 */ 405 std::vector<unsigned int> component_to_block_map; 406 407 /** 408 * The mg_constrained_dofs of the level systems. 409 */ 410 SmartPointer<const MGConstrainedDoFs, MGLevelGlobalTransfer<VectorType>> 411 mg_constrained_dofs; 412 413 private: 414 /** 415 * This function is called to make sure that build() has been invoked. 416 */ 417 template <int dim, int spacedim> 418 void 419 assert_built(const DoFHandler<dim, spacedim> &dof_handler) const; 420 }; 421 422 423 424 /** 425 * Implementation of transfer between the global vectors and the multigrid 426 * levels for use in the derived class MGTransferPrebuilt and other classes. 427 * This class is a specialization for the case of 428 * LinearAlgebra::distributed::Vector that requires a few different calling 429 * routines as compared to the %parallel vectors in the PETScWrappers and 430 * TrilinosWrappers namespaces. 431 */ 432 template <typename Number> 433 class MGLevelGlobalTransfer<LinearAlgebra::distributed::Vector<Number>> 434 : public MGTransferBase<LinearAlgebra::distributed::Vector<Number>> 435 { 436 public: 437 /** 438 * Reset the object to the state it had right after the default constructor. 439 */ 440 void 441 clear(); 442 443 /** 444 * Transfer from a vector on the global grid to vectors defined on each of 445 * the levels separately for the active degrees of freedom. In particular, for 446 * a globally refined mesh only the finest level in @p dst is filled as a 447 * plain copy of @p src. All the other level objects are left untouched. 448 */ 449 template <int dim, typename Number2, int spacedim> 450 void 451 copy_to_mg(const DoFHandler<dim, spacedim> &dof_handler, 452 MGLevelObject<LinearAlgebra::distributed::Vector<Number>> &dst, 453 const LinearAlgebra::distributed::Vector<Number2> &src) const; 454 455 /** 456 * Transfer from multi-level vector to normal vector. 457 * 458 * Copies data from active portions of an MGVector into the respective 459 * positions of a <tt>Vector<number></tt>. In order to keep the result 460 * consistent, constrained degrees of freedom are set to zero. 461 */ 462 template <int dim, typename Number2, int spacedim> 463 void 464 copy_from_mg( 465 const DoFHandler<dim, spacedim> & dof_handler, 466 LinearAlgebra::distributed::Vector<Number2> &dst, 467 const MGLevelObject<LinearAlgebra::distributed::Vector<Number>> &src) const; 468 469 /** 470 * Add a multi-level vector to a normal vector. 471 * 472 * Works as the previous function, but probably not for continuous elements. 473 */ 474 template <int dim, typename Number2, int spacedim> 475 void 476 copy_from_mg_add( 477 const DoFHandler<dim, spacedim> & dof_handler, 478 LinearAlgebra::distributed::Vector<Number2> &dst, 479 const MGLevelObject<LinearAlgebra::distributed::Vector<Number>> &src) const; 480 481 /** 482 * If this object operates on BlockVector objects, we need to describe how 483 * the individual vector components are mapped to the blocks of a vector. 484 * For example, for a Stokes system, we have dim+1 vector components for 485 * velocity and pressure, but we may want to use block vectors with only two 486 * blocks for all velocities in one block, and the pressure variables in the 487 * other. 488 * 489 * By default, if this function is not called, block vectors have as many 490 * blocks as the finite element has vector components. However, this can be 491 * changed by calling this function with an array that describes how vector 492 * components are to be grouped into blocks. The meaning of the argument is 493 * the same as the one given to the DoFTools::count_dofs_per_component 494 * function. 495 */ 496 void 497 set_component_to_block_map(const std::vector<unsigned int> &map); 498 499 /** 500 * Memory used by this object. 501 */ 502 std::size_t 503 memory_consumption() const; 504 505 /** 506 * Print the copy index fields for debugging purposes. 507 */ 508 void 509 print_indices(std::ostream &os) const; 510 511 protected: 512 /** 513 * Internal function to perform transfer of residuals or solutions 514 * basesd on the flag @p solution_transfer. 515 */ 516 template <int dim, typename Number2, int spacedim> 517 void 518 copy_to_mg(const DoFHandler<dim, spacedim> &dof_handler, 519 MGLevelObject<LinearAlgebra::distributed::Vector<Number>> &dst, 520 const LinearAlgebra::distributed::Vector<Number2> & src, 521 const bool solution_transfer) const; 522 523 /** 524 * Internal function to @p fill copy_indices*. Called by derived classes. 525 */ 526 template <int dim, int spacedim> 527 void 528 fill_and_communicate_copy_indices( 529 const DoFHandler<dim, spacedim> &dof_handler); 530 531 /** 532 * Sizes of the multi-level vectors. 533 */ 534 std::vector<types::global_dof_index> sizes; 535 536 /** 537 * Mapping for the copy_to_mg() and copy_from_mg() functions. Here only 538 * index pairs locally owned is stored. 539 * 540 * The data is organized as follows: one table per level. This table has two 541 * rows. The first row contains the global index, the second one the level 542 * index. 543 */ 544 std::vector<Table<2, unsigned int>> copy_indices; 545 546 /** 547 * Same as above, but used to transfer solution vectors. 548 */ 549 std::vector<Table<2, unsigned int>> solution_copy_indices; 550 551 /** 552 * Additional degrees of freedom for the copy_to_mg() function. These are 553 * the ones where the global degree of freedom is locally owned and the 554 * level degree of freedom is not. 555 * 556 * Organization of the data is like for @p copy_indices. 557 */ 558 std::vector<Table<2, unsigned int>> copy_indices_global_mine; 559 560 /** 561 * Same as above, but used to transfer solution vectors. 562 */ 563 std::vector<Table<2, unsigned int>> solution_copy_indices_global_mine; 564 565 /** 566 * Additional degrees of freedom for the copy_from_mg() function. These are 567 * the ones where the level degree of freedom is locally owned and the 568 * global degree of freedom is not. 569 * 570 * Organization of the data is like for @p copy_indices. 571 */ 572 std::vector<Table<2, unsigned int>> copy_indices_level_mine; 573 574 /** 575 * Same as above, but used to transfer solution vectors. 576 */ 577 std::vector<Table<2, unsigned int>> solution_copy_indices_level_mine; 578 579 /** 580 * This variable stores whether the copy operation from the global to the 581 * level vector is actually a plain copy to the finest level. This means that 582 * the grid has no adaptive refinement and the numbering on the finest 583 * multigrid level is the same as in the global case. 584 */ 585 bool perform_plain_copy; 586 587 /** 588 * This variable stores whether the copy operation from the global to the 589 * level vector is actually a plain copy to the finest level except for a 590 * renumbering within the finest level of the degrees of freedom. This means 591 * that the grid has no adaptive refinement. 592 */ 593 bool perform_renumbered_plain_copy; 594 595 /** 596 * The vector that stores what has been given to the 597 * set_component_to_block_map() function. 598 */ 599 std::vector<unsigned int> component_to_block_map; 600 601 /** 602 * The mg_constrained_dofs of the level systems. 603 */ 604 SmartPointer< 605 const MGConstrainedDoFs, 606 MGLevelGlobalTransfer<LinearAlgebra::distributed::Vector<Number>>> 607 mg_constrained_dofs; 608 609 /** 610 * In the function copy_to_mg, we need to access ghosted entries of the 611 * global vector for inserting into the level vectors. This vector is 612 * populated with those entries. 613 */ 614 mutable LinearAlgebra::distributed::Vector<Number> ghosted_global_vector; 615 616 /** 617 * Same as above but used when working with solution vectors. 618 */ 619 mutable LinearAlgebra::distributed::Vector<Number> 620 solution_ghosted_global_vector; 621 622 /** 623 * In the function copy_from_mg, we access all level vectors with certain 624 * ghost entries for inserting the result into a global vector. 625 */ 626 mutable MGLevelObject<LinearAlgebra::distributed::Vector<Number>> 627 ghosted_level_vector; 628 629 /** 630 * Same as above but used when working with solution vectors. 631 */ 632 mutable MGLevelObject<LinearAlgebra::distributed::Vector<Number>> 633 solution_ghosted_level_vector; 634 635 private: 636 /** 637 * This function is called to make sure that build() has been invoked. 638 */ 639 template <int dim, int spacedim> 640 void 641 assert_built(const DoFHandler<dim, spacedim> &dof_handler) const; 642 }; 643 644 645 646 /** 647 * Implementation of the MGTransferBase interface for which the transfer 648 * operations are prebuilt upon construction of the object of this class as 649 * matrices. This is the fast way, since it only needs to build the operation 650 * once by looping over all cells and storing the result in a matrix for each 651 * level, but requires additional memory. 652 * 653 * See MGTransferBase to find out which of the transfer classes is best for 654 * your needs. 655 */ 656 template <typename VectorType> 657 class MGTransferPrebuilt : public MGLevelGlobalTransfer<VectorType> 658 { 659 public: 660 /** 661 * Constructor without constraint matrices. Use this constructor only with 662 * discontinuous finite elements or with no local refinement. 663 */ 664 MGTransferPrebuilt() = default; 665 666 /** 667 * Constructor with constraints. Equivalent to the default constructor 668 * followed by initialize_constraints(). 669 */ 670 MGTransferPrebuilt(const MGConstrainedDoFs &mg_constrained_dofs); 671 672 /** 673 * Destructor. 674 */ 675 virtual ~MGTransferPrebuilt() override = default; 676 677 /** 678 * Initialize the constraints to be used in build(). 679 */ 680 void 681 initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs); 682 683 /** 684 * Reset the object to the state it had right after the default constructor. 685 */ 686 void 687 clear(); 688 689 /** 690 * Actually build the information required for the transfer operations. Needs 691 * to be called before prolongate() or restrict_and_add() can be used. 692 */ 693 template <int dim, int spacedim> 694 void 695 build(const DoFHandler<dim, spacedim> &dof_handler); 696 697 /** 698 * Actually build the prolongation matrices for each level. 699 * 700 * @deprecated use build() instead. 701 */ 702 template <int dim, int spacedim> 703 DEAL_II_DEPRECATED void 704 build_matrices(const DoFHandler<dim, spacedim> &dof_handler); 705 706 /** 707 * Prolongate a vector from level <tt>to_level-1</tt> to level 708 * <tt>to_level</tt> using the embedding matrices of the underlying finite 709 * element. The previous content of <tt>dst</tt> is overwritten. 710 * 711 * @arg src is a vector with as many elements as there are degrees of 712 * freedom on the coarser level involved. 713 * 714 * @arg dst has as many elements as there are degrees of freedom on the 715 * finer level. 716 */ 717 virtual void 718 prolongate(const unsigned int to_level, 719 VectorType & dst, 720 const VectorType & src) const override; 721 722 /** 723 * Restrict a vector from level <tt>from_level</tt> to level 724 * <tt>from_level-1</tt> using the transpose operation of the @p prolongate 725 * method. If the region covered by cells on level <tt>from_level</tt> is 726 * smaller than that of level <tt>from_level-1</tt> (local refinement), then 727 * some degrees of freedom in <tt>dst</tt> are active and will not be 728 * altered. For the other degrees of freedom, the result of the restriction 729 * is added. 730 * 731 * @arg src is a vector with as many elements as there are degrees of 732 * freedom on the finer level involved. 733 * 734 * @arg dst has as many elements as there are degrees of freedom on the 735 * coarser level. 736 */ 737 virtual void 738 restrict_and_add(const unsigned int from_level, 739 VectorType & dst, 740 const VectorType & src) const override; 741 742 /** 743 * Finite element does not provide prolongation matrices. 744 */ 745 DeclException0(ExcNoProlongation); 746 747 /** 748 * You have to call build() before using this object. 749 */ 750 DeclException0(ExcMatricesNotBuilt); 751 752 /** 753 * Memory used by this object. 754 */ 755 std::size_t 756 memory_consumption() const; 757 758 /** 759 * Print all the matrices for debugging purposes. 760 */ 761 void 762 print_matrices(std::ostream &os) const; 763 764 private: 765 /** 766 * Sparsity patterns for transfer matrices. 767 */ 768 std::vector< 769 std::shared_ptr<typename internal::MatrixSelector<VectorType>::Sparsity>> 770 prolongation_sparsities; 771 772 /** 773 * The actual prolongation matrix. column indices belong to the dof indices 774 * of the mother cell, i.e. the coarse level. while row indices belong to 775 * the child cell, i.e. the fine level. 776 */ 777 std::vector< 778 std::shared_ptr<typename internal::MatrixSelector<VectorType>::Matrix>> 779 prolongation_matrices; 780 781 /** 782 * Degrees of freedom on the refinement edge excluding those on the 783 * boundary. 784 */ 785 std::vector<std::vector<bool>> interface_dofs; 786 }; 787 788 789 /*@}*/ 790 791 792 DEAL_II_NAMESPACE_CLOSE 793 794 #endif 795