Reference documentation for deal.II version GIT 3779fa9eb4 2023-09-28 13:00:02+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\}}\)
scalapack.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2023 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_scalapack_h
17 #define dealii_scalapack_h
18 
19 #include <deal.II/base/config.h>
20 
21 #ifdef DEAL_II_WITH_SCALAPACK
22 
23 # include <deal.II/base/exceptions.h>
24 # include <deal.II/base/mpi.h>
25 # include <deal.II/base/mpi_stub.h>
26 # include <deal.II/base/mutex.h>
28 
29 # include <deal.II/lac/full_matrix.h>
32 
33 # include <limits>
34 # include <memory>
35 
37 
75 template <typename NumberType>
76 class ScaLAPACKMatrix : protected TransposeTable<NumberType>
77 {
78 public:
82  using size_type = unsigned int;
83 
93  const size_type n_rows,
94  const size_type n_columns,
95  const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
96  const size_type row_block_size = 32,
97  const size_type column_block_size = 32,
99 
109  const size_type size,
110  const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
111  const size_type block_size = 32,
114 
128  const std::string &filename,
129  const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
130  const size_type row_block_size = 32,
131  const size_type column_block_size = 32);
132 
136  ~ScaLAPACKMatrix() override = default;
137 
146  void
147  reinit(
148  const size_type n_rows,
149  const size_type n_columns,
150  const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
151  const size_type row_block_size = 32,
152  const size_type column_block_size = 32,
154 
162  void
163  reinit(const size_type size,
164  const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
165  const size_type block_size = 32,
168 
172  void
174 
179  get_property() const;
180 
185  get_state() const;
186 
195 
205  void
207  const unsigned int rank);
208 
215  void
217 
227  void
228  copy_to(LAPACKFullMatrix<NumberType> &matrix, const unsigned int rank) const;
229 
235  void
237 
262  void
264  const std::pair<unsigned int, unsigned int> &offset_A,
265  const std::pair<unsigned int, unsigned int> &offset_B,
266  const std::pair<unsigned int, unsigned int> &submatrix_size) const;
267 
276  void
278 
290  void
292  const NumberType a = 0.,
293  const NumberType b = 1.,
294  const bool transpose_B = false);
295 
305  void
306  add(const NumberType b, const ScaLAPACKMatrix<NumberType> &B);
307 
317  void
318  Tadd(const NumberType b, const ScaLAPACKMatrix<NumberType> &B);
319 
340  void
341  mult(const NumberType b,
343  const NumberType c,
345  const bool transpose_A = false,
346  const bool transpose_B = false) const;
347 
365  void
368  const bool adding = false) const;
369 
387  void
390  const bool adding = false) const;
391 
409  void
412  const bool adding = false) const;
413 
432  void
435  const bool adding = false) const;
436 
457  void
458  save(const std::string &filename,
459  const std::pair<unsigned int, unsigned int> &chunk_size =
460  std::make_pair(numbers::invalid_unsigned_int,
462 
474  void
475  load(const std::string &filename);
476 
482  void
484 
490  void
492 
505  void
506  invert();
507 
523  std::vector<NumberType>
525  const std::pair<unsigned int, unsigned int> &index_limits,
526  const bool compute_eigenvectors);
527 
536  std::vector<NumberType>
538  const std::pair<NumberType, NumberType> &value_limits,
539  const bool compute_eigenvectors);
540 
557  std::vector<NumberType>
559  const std::pair<unsigned int, unsigned int> &index_limits,
560  const bool compute_eigenvectors);
561 
572  std::vector<NumberType>
574  const std::pair<NumberType, NumberType> &value_limits,
575  const bool compute_eigenvectors);
576 
605  std::vector<NumberType>
607  ScaLAPACKMatrix<NumberType> *VT = nullptr);
608 
649  void
650  least_squares(ScaLAPACKMatrix<NumberType> &B, const bool transpose = false);
651 
675  unsigned int
676  pseudoinverse(const NumberType ratio);
677 
692  NumberType
693  reciprocal_condition_number(const NumberType a_norm) const;
694 
698  NumberType
699  l1_norm() const;
700 
704  NumberType
705  linfty_norm() const;
706 
710  NumberType
711  frobenius_norm() const;
712 
716  size_type
717  m() const;
718 
722  size_type
723  n() const;
724 
728  unsigned int
729  local_m() const;
730 
734  unsigned int
735  local_n() const;
736 
740  unsigned int
741  global_row(const unsigned int loc_row) const;
742 
746  unsigned int
747  global_column(const unsigned int loc_column) const;
748 
752  NumberType
753  local_el(const unsigned int loc_row, const unsigned int loc_column) const;
754 
758  NumberType &
759  local_el(const unsigned int loc_row, const unsigned int loc_column);
760 
771  template <class InputVector>
772  void
773  scale_columns(const InputVector &factors);
774 
785  template <class InputVector>
786  void
787  scale_rows(const InputVector &factors);
788 
789 private:
794  NumberType
795  norm_symmetric(const char type) const;
796 
801  NumberType
802  norm_general(const char type) const;
803 
814  std::vector<NumberType>
816  const bool compute_eigenvectors,
817  const std::pair<unsigned int, unsigned int> &index_limits =
818  std::make_pair(numbers::invalid_unsigned_int,
820  const std::pair<NumberType, NumberType> &value_limits =
821  std::make_pair(std::numeric_limits<NumberType>::quiet_NaN(),
822  std::numeric_limits<NumberType>::quiet_NaN()));
823 
843  std::vector<NumberType>
845  const bool compute_eigenvectors,
846  const std::pair<unsigned int, unsigned int> &index_limits =
847  std::make_pair(numbers::invalid_unsigned_int,
849  const std::pair<NumberType, NumberType> &value_limits =
850  std::make_pair(std::numeric_limits<NumberType>::quiet_NaN(),
851  std::numeric_limits<NumberType>::quiet_NaN()));
852 
853  /*
854  * Stores the distributed matrix in @p filename
855  * using serial routines
856  */
857  void
858  save_serial(const std::string &filename,
859  const std::pair<unsigned int, unsigned int> &chunk_size) const;
860 
861  /*
862  * Loads the distributed matrix from file @p filename
863  * using serial routines
864  */
865  void
866  load_serial(const std::string &filename);
867 
868  /*
869  * Stores the distributed matrix in @p filename
870  * using parallel routines
871  */
872  void
873  save_parallel(const std::string &filename,
874  const std::pair<unsigned int, unsigned int> &chunk_size) const;
875 
876  /*
877  * Loads the distributed matrix from file @p filename
878  * using parallel routines
879  */
880  void
881  load_parallel(const std::string &filename);
882 
888 
894 
900  std::shared_ptr<const Utilities::MPI::ProcessGrid> grid;
901 
905  int n_rows;
906 
911 
916 
921 
926 
931 
935  int descriptor[9];
936 
940  mutable std::vector<NumberType> work;
941 
945  mutable std::vector<int> iwork;
946 
951  std::vector<int> ipiv;
952 
957  const char uplo;
958 
963  const int first_process_row;
964 
970 
975  const int submatrix_row;
976 
981  const int submatrix_column;
982 
987 };
988 
989 // ----------------------- inline functions ----------------------------
990 
991 # ifndef DOXYGEN
992 
993 template <typename NumberType>
994 inline NumberType
995 ScaLAPACKMatrix<NumberType>::local_el(const unsigned int loc_row,
996  const unsigned int loc_column) const
997 {
998  return (*this)(loc_row, loc_column);
999 }
1000 
1001 
1002 
1003 template <typename NumberType>
1004 inline NumberType &
1005 ScaLAPACKMatrix<NumberType>::local_el(const unsigned int loc_row,
1006  const unsigned int loc_column)
1007 {
1008  return (*this)(loc_row, loc_column);
1009 }
1010 
1011 
1012 template <typename NumberType>
1013 inline unsigned int
1015 {
1016  return n_rows;
1017 }
1018 
1019 
1020 
1021 template <typename NumberType>
1022 inline unsigned int
1024 {
1025  return n_columns;
1026 }
1027 
1028 
1029 
1030 template <typename NumberType>
1031 unsigned int
1033 {
1034  return n_local_rows;
1035 }
1036 
1037 
1038 
1039 template <typename NumberType>
1040 unsigned int
1042 {
1043  return n_local_columns;
1044 }
1045 
1046 
1047 # endif // DOXYGEN
1048 
1050 
1051 #endif // DEAL_II_WITH_SCALAPACK
1052 
1053 #endif
int descriptor[9]
Definition: scalapack.h:935
std::vector< NumberType > compute_SVD(ScaLAPACKMatrix< NumberType > *U=nullptr, ScaLAPACKMatrix< NumberType > *VT=nullptr)
Definition: scalapack.cc:2024
NumberType & local_el(const unsigned int loc_row, const unsigned int loc_column)
std::vector< NumberType > eigenpairs_symmetric_by_value(const std::pair< NumberType, NumberType > &value_limits, const bool compute_eigenvectors)
Definition: scalapack.cc:1453
unsigned int local_n() const
NumberType frobenius_norm() const
Definition: scalapack.cc:2413
unsigned int pseudoinverse(const NumberType ratio)
Definition: scalapack.cc:2236
std::vector< NumberType > eigenpairs_symmetric_by_value_MRRR(const std::pair< NumberType, NumberType > &value_limits, const bool compute_eigenvectors)
Definition: scalapack.cc:1796
const int first_process_column
Definition: scalapack.h:969
void copy_from(const LAPACKFullMatrix< NumberType > &matrix, const unsigned int rank)
Definition: scalapack.cc:363
void save_parallel(const std::string &filename, const std::pair< unsigned int, unsigned int > &chunk_size) const
Definition: scalapack.cc:2817
void least_squares(ScaLAPACKMatrix< NumberType > &B, const bool transpose=false)
Definition: scalapack.cc:2144
NumberType local_el(const unsigned int loc_row, const unsigned int loc_column) const
ScaLAPACKMatrix< NumberType > & operator=(const FullMatrix< NumberType > &)
Definition: scalapack.cc:332
void Tadd(const NumberType b, const ScaLAPACKMatrix< NumberType > &B)
Definition: scalapack.cc:1057
void mTmult(ScaLAPACKMatrix< NumberType > &C, const ScaLAPACKMatrix< NumberType > &B, const bool adding=false) const
Definition: scalapack.cc:1212
std::vector< int > iwork
Definition: scalapack.h:945
void add(const ScaLAPACKMatrix< NumberType > &B, const NumberType a=0., const NumberType b=1., const bool transpose_B=false)
Definition: scalapack.cc:991
size_type n() const
LAPACKSupport::State get_state() const
Definition: scalapack.cc:323
LAPACKSupport::Property get_property() const
Definition: scalapack.cc:314
int column_block_size
Definition: scalapack.h:920
void Tmmult(ScaLAPACKMatrix< NumberType > &C, const ScaLAPACKMatrix< NumberType > &B, const bool adding=false) const
Definition: scalapack.cc:1198
std::vector< NumberType > eigenpairs_symmetric_MRRR(const bool compute_eigenvectors, const std::pair< unsigned int, unsigned int > &index_limits=std::make_pair(numbers::invalid_unsigned_int, numbers::invalid_unsigned_int), const std::pair< NumberType, NumberType > &value_limits=std::make_pair(std::numeric_limits< NumberType >::quiet_NaN(), std::numeric_limits< NumberType >::quiet_NaN()))
Definition: scalapack.cc:1814
void scale_rows(const InputVector &factors)
Definition: scalapack.cc:3517
const int first_process_row
Definition: scalapack.h:963
std::shared_ptr< const Utilities::MPI::ProcessGrid > grid
Definition: scalapack.h:900
ScaLAPACKMatrix(const size_type n_rows, const size_type n_columns, const std::shared_ptr< const Utilities::MPI::ProcessGrid > &process_grid, const size_type row_block_size=32, const size_type column_block_size=32, const LAPACKSupport::Property property=LAPACKSupport::Property::general)
Definition: scalapack.cc:81
void load(const std::string &filename)
Definition: scalapack.cc:3053
const int submatrix_column
Definition: scalapack.h:981
const int submatrix_row
Definition: scalapack.h:975
void mmult(ScaLAPACKMatrix< NumberType > &C, const ScaLAPACKMatrix< NumberType > &B, const bool adding=false) const
Definition: scalapack.cc:1184
std::vector< NumberType > eigenpairs_symmetric(const bool compute_eigenvectors, const std::pair< unsigned int, unsigned int > &index_limits=std::make_pair(numbers::invalid_unsigned_int, numbers::invalid_unsigned_int), const std::pair< NumberType, NumberType > &value_limits=std::make_pair(std::numeric_limits< NumberType >::quiet_NaN(), std::numeric_limits< NumberType >::quiet_NaN()))
Definition: scalapack.cc:1473
void save_serial(const std::string &filename, const std::pair< unsigned int, unsigned int > &chunk_size) const
Definition: scalapack.cc:2660
NumberType norm_general(const char type) const
Definition: scalapack.cc:2427
~ScaLAPACKMatrix() override=default
std::vector< NumberType > work
Definition: scalapack.h:940
void save(const std::string &filename, const std::pair< unsigned int, unsigned int > &chunk_size=std::make_pair(numbers::invalid_unsigned_int, numbers::invalid_unsigned_int)) const
Definition: scalapack.cc:2618
void load_parallel(const std::string &filename)
Definition: scalapack.cc:3252
NumberType l1_norm() const
Definition: scalapack.cc:2385
void compute_lu_factorization()
Definition: scalapack.cc:1273
NumberType norm_symmetric(const char type) const
Definition: scalapack.cc:2486
const char uplo
Definition: scalapack.h:957
void mult(const NumberType b, const ScaLAPACKMatrix< NumberType > &B, const NumberType c, ScaLAPACKMatrix< NumberType > &C, const bool transpose_A=false, const bool transpose_B=false) const
Definition: scalapack.cc:1067
LAPACKSupport::Property property
Definition: scalapack.h:893
void set_property(const LAPACKSupport::Property property)
Definition: scalapack.cc:304
void reinit(const size_type n_rows, const size_type n_columns, const std::shared_ptr< const Utilities::MPI::ProcessGrid > &process_grid, const size_type row_block_size=32, const size_type column_block_size=32, const LAPACKSupport::Property property=LAPACKSupport::Property::general)
Definition: scalapack.cc:217
void load_serial(const std::string &filename)
Definition: scalapack.cc:3074
size_type m() const
std::vector< NumberType > eigenpairs_symmetric_by_index_MRRR(const std::pair< unsigned int, unsigned int > &index_limits, const bool compute_eigenvectors)
Definition: scalapack.cc:1773
std::vector< int > ipiv
Definition: scalapack.h:951
NumberType reciprocal_condition_number(const NumberType a_norm) const
Definition: scalapack.cc:2323
LAPACKSupport::State state
Definition: scalapack.h:887
void TmTmult(ScaLAPACKMatrix< NumberType > &C, const ScaLAPACKMatrix< NumberType > &B, const bool adding=false) const
Definition: scalapack.cc:1226
Threads::Mutex mutex
Definition: scalapack.h:986
std::vector< NumberType > eigenpairs_symmetric_by_index(const std::pair< unsigned int, unsigned int > &index_limits, const bool compute_eigenvectors)
Definition: scalapack.cc:1430
unsigned int global_column(const unsigned int loc_column) const
Definition: scalapack.cc:515
void copy_to(FullMatrix< NumberType > &matrix) const
Definition: scalapack.cc:667
unsigned int global_row(const unsigned int loc_row) const
Definition: scalapack.cc:498
unsigned int local_m() const
void compute_cholesky_factorization()
Definition: scalapack.cc:1240
void copy_transposed(const ScaLAPACKMatrix< NumberType > &B)
Definition: scalapack.cc:981
NumberType linfty_norm() const
Definition: scalapack.cc:2399
void scale_columns(const InputVector &factors)
Definition: scalapack.cc:3506
typename AlignedVector< T >::size_type size_type
Definition: table.h:449
const TableIndices< N > & size() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
constexpr int chunk_size
Definition: cuda_size.h:35
constexpr int block_size
Definition: cuda_size.h:29
@ matrix
Contents is actually a matrix.
static const char U
@ symmetric
Matrix is symmetric.
@ general
No special properties.
SymmetricTensor< 2, dim, Number > C(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
static const unsigned int invalid_unsigned_int
Definition: types.h:213