Reference documentation for deal.II version 8.1.0
matrix_tools.h
1 // ---------------------------------------------------------------------
2 // @f$Id: matrix_tools.h 31118 2013-10-04 13:05:10Z bangerth @f$
3 //
4 // Copyright (C) 1998 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__matrix_tools_h
18 #define __deal2__matrix_tools_h
19 
20 
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/thread_management.h>
24 #include <deal.II/lac/constraint_matrix.h>
25 #include <deal.II/dofs/function_map.h>
26 
27 #include <map>
28 
30 
31 
32 // forward declarations
33 template <int dim> class Quadrature;
34 
35 
36 template<typename number> class Vector;
37 template<typename number> class FullMatrix;
38 template<typename number> class SparseMatrix;
39 
40 template <typename number> class BlockSparseMatrix;
41 template <typename Number> class BlockVector;
42 
43 template <int dim, int spacedim> class Mapping;
44 template <int dim, int spacedim> class DoFHandler;
45 template <int dim, int spacedim> class MGDoFHandler;
46 template <int dim, int spacedim> class FEValues;
47 
48 namespace hp
49 {
50  template <int> class QCollection;
51  template <int, int> class MappingCollection;
52  template <int, int> class DoFHandler;
53 }
54 
55 
56 #ifdef DEAL_II_WITH_PETSC
57 namespace PETScWrappers
58 {
59  class SparseMatrix;
60  class Vector;
61  namespace MPI
62  {
63  class SparseMatrix;
64  class BlockSparseMatrix;
65  class Vector;
66  class BlockVector;
67  }
68 }
69 #endif
70 
71 #ifdef DEAL_II_WITH_TRILINOS
73 {
74  class SparseMatrix;
75  class Vector;
76  class BlockSparseMatrix;
77  class BlockVector;
78  namespace MPI
79  {
80  class Vector;
81  class BlockVector;
82  }
83 }
84 #endif
85 
86 
227 namespace MatrixCreator
228 {
242  template <int dim, typename number, int spacedim>
243  void create_mass_matrix (const Mapping<dim, spacedim> &mapping,
244  const DoFHandler<dim,spacedim> &dof,
245  const Quadrature<dim> &q,
246  SparseMatrix<number> &matrix,
247  const Function<spacedim> *const a = 0,
248  const ConstraintMatrix &constraints = ConstraintMatrix());
249 
254  template <int dim, typename number, int spacedim>
256  const Quadrature<dim> &q,
257  SparseMatrix<number> &matrix,
258  const Function<spacedim> *const a = 0,
259  const ConstraintMatrix &constraints = ConstraintMatrix());
260 
274  template <int dim, typename number, int spacedim>
275  void create_mass_matrix (const Mapping<dim, spacedim> &mapping,
276  const DoFHandler<dim,spacedim> &dof,
277  const Quadrature<dim> &q,
278  SparseMatrix<number> &matrix,
279  const Function<spacedim> &rhs,
280  Vector<double> &rhs_vector,
281  const Function<spacedim> *const a = 0,
282  const ConstraintMatrix &constraints = ConstraintMatrix());
283 
288  template <int dim, typename number, int spacedim>
290  const Quadrature<dim> &q,
291  SparseMatrix<number> &matrix,
292  const Function<spacedim> &rhs,
293  Vector<double> &rhs_vector,
294  const Function<spacedim> *const a = 0,
295  const ConstraintMatrix &constraints = ConstraintMatrix());
296 
300  template <int dim, typename number, int spacedim>
302  const hp::DoFHandler<dim,spacedim> &dof,
303  const hp::QCollection<dim> &q,
304  SparseMatrix<number> &matrix,
305  const Function<spacedim> *const a = 0,
306  const ConstraintMatrix &constraints = ConstraintMatrix());
307 
311  template <int dim, typename number, int spacedim>
313  const hp::QCollection<dim> &q,
314  SparseMatrix<number> &matrix,
315  const Function<spacedim> *const a = 0,
316  const ConstraintMatrix &constraints = ConstraintMatrix());
317 
321  template <int dim, typename number, int spacedim>
323  const hp::DoFHandler<dim,spacedim> &dof,
324  const hp::QCollection<dim> &q,
325  SparseMatrix<number> &matrix,
326  const Function<spacedim> &rhs,
327  Vector<double> &rhs_vector,
328  const Function<spacedim> *const a = 0,
329  const ConstraintMatrix &constraints = ConstraintMatrix());
330 
334  template <int dim, typename number, int spacedim>
336  const hp::QCollection<dim> &q,
337  SparseMatrix<number> &matrix,
338  const Function<spacedim> &rhs,
339  Vector<double> &rhs_vector,
340  const Function<spacedim> *const a = 0,
341  const ConstraintMatrix &constraints = ConstraintMatrix());
342 
343 
366  template <int dim, int spacedim>
367 
369  const DoFHandler<dim,spacedim> &dof,
370  const Quadrature<dim-1> &q,
371  SparseMatrix<double> &matrix,
372  const typename FunctionMap<spacedim>::type &boundary_functions,
373  Vector<double> &rhs_vector,
374  std::vector<types::global_dof_index> &dof_to_boundary_mapping,
375  const Function<spacedim> *const weight = 0,
376  std::vector<unsigned int> component_mapping = std::vector<unsigned int>());
377 
378 // * Same function, but for 1d.
379 // */
380 //
381 // void create_boundary_mass_matrix (const Mapping<1,1> &mapping,
382 // const DoFHandler<1,1> &dof,
383 // const Quadrature<0> &q,
384 // SparseMatrix<double> &matrix,
385 // const FunctionMap<1>::type &boundary_functions,
386 // Vector<double> &rhs_vector,
387 // std::vector<types::global_dof_index>&dof_to_boundary_mapping,
388 // const Function<1> * const a = 0);
389 // //codimension 1
390 //
391 // void create_boundary_mass_matrix (const Mapping<1,2> &mapping,
392 // const DoFHandler<1,2> &dof,
393 // const Quadrature<0> &q,
394 // SparseMatrix<double> &matrix,
395 // const FunctionMap<2>::type &boundary_functions,
396 // Vector<double> &rhs_vector,
397 // std::vector<types::global_dof_index>&dof_to_boundary_mapping,
398 // const Function<2> * const a = 0);
399 
400 
401 
406  template <int dim, int spacedim>
407 
409  const Quadrature<dim-1> &q,
410  SparseMatrix<double> &matrix,
411  const typename FunctionMap<spacedim>::type &boundary_functions,
412  Vector<double> &rhs_vector,
413  std::vector<types::global_dof_index> &dof_to_boundary_mapping,
414  const Function<spacedim> *const a = 0,
415  std::vector<unsigned int> component_mapping = std::vector<unsigned int>());
416 
420  template <int dim, int spacedim>
421 
423  const hp::DoFHandler<dim,spacedim> &dof,
424  const hp::QCollection<dim-1> &q,
425  SparseMatrix<double> &matrix,
426  const typename FunctionMap<spacedim>::type &boundary_functions,
427  Vector<double> &rhs_vector,
428  std::vector<types::global_dof_index> &dof_to_boundary_mapping,
429  const Function<spacedim> *const a = 0,
430  std::vector<unsigned int> component_mapping = std::vector<unsigned int>());
431 
435 //
436 // void create_boundary_mass_matrix (const hp::MappingCollection<1,1> &mapping,
437 // const hp::DoFHandler<1,1> &dof,
438 // const hp::QCollection<0> &q,
439 // SparseMatrix<double> &matrix,
440 // const FunctionMap<1>::type &boundary_functions,
441 // Vector<double> &rhs_vector,
442 // std::vector<types::global_dof_index>&dof_to_boundary_mapping,
443 // const Function<1> * const a = 0);
444 
448  template <int dim, int spacedim>
449 
451  const hp::QCollection<dim-1> &q,
452  SparseMatrix<double> &matrix,
453  const typename FunctionMap<spacedim>::type &boundary_functions,
454  Vector<double> &rhs_vector,
455  std::vector<types::global_dof_index> &dof_to_boundary_mapping,
456  const Function<spacedim> *const a = 0,
457  std::vector<unsigned int> component_mapping = std::vector<unsigned int>());
458 
472  template <int dim, int spacedim>
473  void create_laplace_matrix (const Mapping<dim, spacedim> &mapping,
474  const DoFHandler<dim,spacedim> &dof,
475  const Quadrature<dim> &q,
476  SparseMatrix<double> &matrix,
477  const Function<spacedim> *const a = 0,
478  const ConstraintMatrix &constraints = ConstraintMatrix());
479 
484  template <int dim, int spacedim>
486  const Quadrature<dim> &q,
487  SparseMatrix<double> &matrix,
488  const Function<spacedim> *const a = 0,
489  const ConstraintMatrix &constraints = ConstraintMatrix());
490 
504  template <int dim, int spacedim>
505  void create_laplace_matrix (const Mapping<dim, spacedim> &mapping,
506  const DoFHandler<dim,spacedim> &dof,
507  const Quadrature<dim> &q,
508  SparseMatrix<double> &matrix,
509  const Function<spacedim> &rhs,
510  Vector<double> &rhs_vector,
511  const Function<spacedim> *const a = 0,
512  const ConstraintMatrix &constraints = ConstraintMatrix());
513 
518  template <int dim, int spacedim>
520  const Quadrature<dim> &q,
521  SparseMatrix<double> &matrix,
522  const Function<spacedim> &rhs,
523  Vector<double> &rhs_vector,
524  const Function<spacedim> *const a = 0,
525  const ConstraintMatrix &constraints = ConstraintMatrix());
526 
531  template <int dim, int spacedim>
533  const hp::DoFHandler<dim,spacedim> &dof,
534  const hp::QCollection<dim> &q,
535  SparseMatrix<double> &matrix,
536  const Function<spacedim> *const a = 0,
537  const ConstraintMatrix &constraints = ConstraintMatrix());
538 
543  template <int dim, int spacedim>
545  const hp::QCollection<dim> &q,
546  SparseMatrix<double> &matrix,
547  const Function<spacedim> *const a = 0,
548  const ConstraintMatrix &constraints = ConstraintMatrix());
549 
554  template <int dim, int spacedim>
556  const hp::DoFHandler<dim,spacedim> &dof,
557  const hp::QCollection<dim> &q,
558  SparseMatrix<double> &matrix,
559  const Function<spacedim> &rhs,
560  Vector<double> &rhs_vector,
561  const Function<spacedim> *const a = 0,
562  const ConstraintMatrix &constraints = ConstraintMatrix());
563 
568  template <int dim, int spacedim>
570  const hp::QCollection<dim> &q,
571  SparseMatrix<double> &matrix,
572  const Function<spacedim> &rhs,
573  Vector<double> &rhs_vector,
574  const Function<spacedim> *const a = 0,
575  const ConstraintMatrix &constraints = ConstraintMatrix());
576 
580  DeclException0 (ExcComponentMismatch);
581 }
582 
583 
584 
786 namespace MatrixTools
787 {
796  using namespace MatrixCreator;
797 
804  template <typename number>
805  void
806  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
807  SparseMatrix<number> &matrix,
808  Vector<number> &solution,
809  Vector<number> &right_hand_side,
810  const bool eliminate_columns = true);
811 
821  template <typename number>
822  void
823  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
825  BlockVector<number> &solution,
826  BlockVector<number> &right_hand_side,
827  const bool eliminate_columns = true);
828 
829 #ifdef DEAL_II_WITH_PETSC
830 
867  void
868  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
870  PETScWrappers::Vector &solution,
871  PETScWrappers::Vector &right_hand_side,
872  const bool eliminate_columns = true);
873 
878  void
879  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
881  PETScWrappers::MPI::Vector &solution,
882  PETScWrappers::MPI::Vector &right_hand_side,
883  const bool eliminate_columns = true);
884 
910  void
911  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
913  PETScWrappers::Vector &solution,
914  PETScWrappers::MPI::Vector &right_hand_side,
915  const bool eliminate_columns = true);
916 
920  void
921  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
924  PETScWrappers::MPI::BlockVector &right_hand_side,
925  const bool eliminate_columns = true);
926 
927 #endif
928 
929 #ifdef DEAL_II_WITH_TRILINOS
930 
968  void
969  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
971  TrilinosWrappers::Vector &solution,
972  TrilinosWrappers::Vector &right_hand_side,
973  const bool eliminate_columns = true);
974 
980  void
981  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
984  TrilinosWrappers::BlockVector &right_hand_side,
985  const bool eliminate_columns = true);
986 
1023  void
1024  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
1027  TrilinosWrappers::MPI::Vector &right_hand_side,
1028  const bool eliminate_columns = true);
1029 
1035  void
1036  apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
1039  TrilinosWrappers::MPI::BlockVector &right_hand_side,
1040  const bool eliminate_columns = true);
1041 #endif
1042 
1079  void
1080  local_apply_boundary_values (const std::map<types::global_dof_index,double> &boundary_values,
1081  const std::vector<types::global_dof_index> &local_dof_indices,
1082  FullMatrix<double> &local_matrix,
1083  Vector<double> &local_rhs,
1084  const bool eliminate_columns);
1085 
1089  DeclException0 (ExcBlocksDontMatch);
1090 }
1091 
1092 
1093 
1094 DEAL_II_NAMESPACE_CLOSE
1095 
1096 #endif
std::map< types::boundary_id, const Function< dim > * > type
Definition: function_map.h:56
void create_laplace_matrix(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< double > &matrix, const Function< spacedim > *const a=0, const ConstraintMatrix &constraints=ConstraintMatrix())
DeclException0(ExcComponentMismatch)
DeclException0(ExcBlocksDontMatch)
void local_apply_boundary_values(const std::map< types::global_dof_index, double > &boundary_values, const std::vector< types::global_dof_index > &local_dof_indices, FullMatrix< double > &local_matrix, Vector< double > &local_rhs, const bool eliminate_columns)
Definition: hp.h:103
void create_boundary_mass_matrix(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim-1 > &q, SparseMatrix< double > &matrix, const typename FunctionMap< spacedim >::type &boundary_functions, Vector< double > &rhs_vector, std::vector< types::global_dof_index > &dof_to_boundary_mapping, const Function< spacedim > *const weight=0, std::vector< unsigned int > component_mapping=std::vector< unsigned int >())
void create_mass_matrix(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const Quadrature< dim > &q, SparseMatrix< number > &matrix, const Function< spacedim > *const a=0, const ConstraintMatrix &constraints=ConstraintMatrix())
void apply_boundary_values(const std::map< types::global_dof_index, double > &boundary_values, SparseMatrix< number > &matrix, Vector< number > &solution, Vector< number > &right_hand_side, const bool eliminate_columns=true)