Reference documentation for deal.II version GIT relicensing-249-g48dc7357c7 2024-03-29 12:30: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\}}\)
Loading...
Searching...
No Matches
scalapack.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2017 - 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_scalapack_h
16#define dealii_scalapack_h
17
18#include <deal.II/base/config.h>
19
20#ifdef DEAL_II_WITH_SCALAPACK
21
23# include <deal.II/base/mpi.h>
25# include <deal.II/base/mutex.h>
27
31
32# include <limits>
33# include <memory>
34
36
74template <typename NumberType>
75class ScaLAPACKMatrix : protected TransposeTable<NumberType>
76{
77public:
81 using size_type = unsigned int;
82
92 const size_type n_rows,
93 const size_type n_columns,
94 const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
95 const size_type row_block_size = 32,
98
108 const size_type size,
109 const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
110 const size_type block_size = 32,
113
127 const std::string &filename,
128 const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
129 const size_type row_block_size = 32,
130 const size_type column_block_size = 32);
131
135 ~ScaLAPACKMatrix() override = default;
136
145 void
146 reinit(
147 const size_type n_rows,
148 const size_type n_columns,
149 const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
150 const size_type row_block_size = 32,
151 const size_type column_block_size = 32,
153
161 void
162 reinit(const size_type size,
163 const std::shared_ptr<const Utilities::MPI::ProcessGrid> &process_grid,
164 const size_type block_size = 32,
167
171 void
173
178 get_property() const;
179
184 get_state() const;
185
194
204 void
206 const unsigned int rank);
207
214 void
215 copy_to(FullMatrix<NumberType> &matrix) const;
216
226 void
227 copy_to(LAPACKFullMatrix<NumberType> &matrix, const unsigned int rank) const;
228
234 void
236
261 void
263 const std::pair<unsigned int, unsigned int> &offset_A,
264 const std::pair<unsigned int, unsigned int> &offset_B,
265 const std::pair<unsigned int, unsigned int> &submatrix_size) const;
266
275 void
277
289 void
291 const NumberType a = 0.,
292 const NumberType b = 1.,
293 const bool transpose_B = false);
294
304 void
305 add(const NumberType b, const ScaLAPACKMatrix<NumberType> &B);
306
316 void
317 Tadd(const NumberType b, const ScaLAPACKMatrix<NumberType> &B);
318
339 void
340 mult(const NumberType b,
342 const NumberType c,
344 const bool transpose_A = false,
345 const bool transpose_B = false) const;
346
364 void
367 const bool adding = false) const;
368
386 void
389 const bool adding = false) const;
390
408 void
411 const bool adding = false) const;
412
431 void
434 const bool adding = false) const;
435
456 void
457 save(const std::string &filename,
458 const std::pair<unsigned int, unsigned int> &chunk_size =
459 std::make_pair(numbers::invalid_unsigned_int,
461
473 void
474 load(const std::string &filename);
475
481 void
483
489 void
491
504 void
505 invert();
506
522 std::vector<NumberType>
524 const std::pair<unsigned int, unsigned int> &index_limits,
525 const bool compute_eigenvectors);
526
535 std::vector<NumberType>
537 const std::pair<NumberType, NumberType> &value_limits,
538 const bool compute_eigenvectors);
539
556 std::vector<NumberType>
558 const std::pair<unsigned int, unsigned int> &index_limits,
559 const bool compute_eigenvectors);
560
571 std::vector<NumberType>
573 const std::pair<NumberType, NumberType> &value_limits,
574 const bool compute_eigenvectors);
575
604 std::vector<NumberType>
606 ScaLAPACKMatrix<NumberType> *VT = nullptr);
607
648 void
650
674 unsigned int
675 pseudoinverse(const NumberType ratio);
676
691 NumberType
692 reciprocal_condition_number(const NumberType a_norm) const;
693
697 NumberType
698 l1_norm() const;
699
703 NumberType
704 linfty_norm() const;
705
709 NumberType
710 frobenius_norm() const;
711
716 m() const;
717
722 n() const;
723
727 unsigned int
728 local_m() const;
729
733 unsigned int
734 local_n() const;
735
739 unsigned int
740 global_row(const unsigned int loc_row) const;
741
745 unsigned int
746 global_column(const unsigned int loc_column) const;
747
751 NumberType
752 local_el(const unsigned int loc_row, const unsigned int loc_column) const;
753
757 NumberType &
758 local_el(const unsigned int loc_row, const unsigned int loc_column);
759
770 template <class InputVector>
771 void
772 scale_columns(const InputVector &factors);
773
784 template <class InputVector>
785 void
786 scale_rows(const InputVector &factors);
787
788private:
793 NumberType
794 norm_symmetric(const char type) const;
795
800 NumberType
801 norm_general(const char type) const;
802
813 std::vector<NumberType>
815 const bool compute_eigenvectors,
816 const std::pair<unsigned int, unsigned int> &index_limits =
817 std::make_pair(numbers::invalid_unsigned_int,
819 const std::pair<NumberType, NumberType> &value_limits =
820 std::make_pair(std::numeric_limits<NumberType>::quiet_NaN(),
821 std::numeric_limits<NumberType>::quiet_NaN()));
822
842 std::vector<NumberType>
844 const bool compute_eigenvectors,
845 const std::pair<unsigned int, unsigned int> &index_limits =
846 std::make_pair(numbers::invalid_unsigned_int,
848 const std::pair<NumberType, NumberType> &value_limits =
849 std::make_pair(std::numeric_limits<NumberType>::quiet_NaN(),
850 std::numeric_limits<NumberType>::quiet_NaN()));
851
852 /*
853 * Stores the distributed matrix in @p filename
854 * using serial routines
855 */
856 void
857 save_serial(const std::string &filename,
858 const std::pair<unsigned int, unsigned int> &chunk_size) const;
859
860 /*
861 * Loads the distributed matrix from file @p filename
862 * using serial routines
863 */
864 void
865 load_serial(const std::string &filename);
866
867 /*
868 * Stores the distributed matrix in @p filename
869 * using parallel routines
870 */
871 void
872 save_parallel(const std::string &filename,
873 const std::pair<unsigned int, unsigned int> &chunk_size) const;
874
875 /*
876 * Loads the distributed matrix from file @p filename
877 * using parallel routines
878 */
879 void
880 load_parallel(const std::string &filename);
881
887
893
899 std::shared_ptr<const Utilities::MPI::ProcessGrid> grid;
900
905
910
915
920
925
930
935
939 mutable std::vector<NumberType> work;
940
944 mutable std::vector<int> iwork;
945
950 std::vector<int> ipiv;
951
956 const char uplo;
957
963
969
974 const int submatrix_row;
975
981
986};
987
988// ----------------------- inline functions ----------------------------
989
990# ifndef DOXYGEN
991
992template <typename NumberType>
993inline NumberType
994ScaLAPACKMatrix<NumberType>::local_el(const unsigned int loc_row,
995 const unsigned int loc_column) const
996{
997 return (*this)(loc_row, loc_column);
998}
999
1000
1001
1002template <typename NumberType>
1003inline NumberType &
1004ScaLAPACKMatrix<NumberType>::local_el(const unsigned int loc_row,
1005 const unsigned int loc_column)
1006{
1007 return (*this)(loc_row, loc_column);
1008}
1009
1010
1011template <typename NumberType>
1012inline unsigned int
1014{
1015 return n_rows;
1016}
1017
1018
1019
1020template <typename NumberType>
1021inline unsigned int
1023{
1024 return n_columns;
1025}
1026
1027
1028
1029template <typename NumberType>
1030unsigned int
1032{
1033 return n_local_rows;
1034}
1035
1036
1037
1038template <typename NumberType>
1039unsigned int
1041{
1042 return n_local_columns;
1043}
1044
1045
1046# endif // DOXYGEN
1047
1049
1050#endif // DEAL_II_WITH_SCALAPACK
1051
1052#endif
int descriptor[9]
Definition scalapack.h:934
std::vector< NumberType > compute_SVD(ScaLAPACKMatrix< NumberType > *U=nullptr, ScaLAPACKMatrix< NumberType > *VT=nullptr)
std::vector< NumberType > eigenpairs_symmetric_by_value(const std::pair< NumberType, NumberType > &value_limits, const bool compute_eigenvectors)
unsigned int local_n() const
NumberType frobenius_norm() const
unsigned int pseudoinverse(const NumberType ratio)
std::vector< NumberType > eigenpairs_symmetric_by_value_MRRR(const std::pair< NumberType, NumberType > &value_limits, const bool compute_eigenvectors)
const int first_process_column
Definition scalapack.h:968
void copy_from(const LAPACKFullMatrix< NumberType > &matrix, const unsigned int rank)
Definition scalapack.cc:362
void save_parallel(const std::string &filename, const std::pair< unsigned int, unsigned int > &chunk_size) const
void least_squares(ScaLAPACKMatrix< NumberType > &B, const bool transpose=false)
NumberType local_el(const unsigned int loc_row, const unsigned int loc_column) const
ScaLAPACKMatrix< NumberType > & operator=(const FullMatrix< NumberType > &)
Definition scalapack.cc:331
void Tadd(const NumberType b, const ScaLAPACKMatrix< NumberType > &B)
void mTmult(ScaLAPACKMatrix< NumberType > &C, const ScaLAPACKMatrix< NumberType > &B, const bool adding=false) const
std::vector< int > iwork
Definition scalapack.h:944
void add(const ScaLAPACKMatrix< NumberType > &B, const NumberType a=0., const NumberType b=1., const bool transpose_B=false)
Definition scalapack.cc:990
size_type n() const
LAPACKSupport::State get_state() const
Definition scalapack.cc:322
LAPACKSupport::Property get_property() const
Definition scalapack.cc:313
void Tmmult(ScaLAPACKMatrix< NumberType > &C, const ScaLAPACKMatrix< NumberType > &B, const bool adding=false) const
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()))
void scale_rows(const InputVector &factors)
const int first_process_row
Definition scalapack.h:962
NumberType & local_el(const unsigned int loc_row, const unsigned int loc_column)
std::shared_ptr< const Utilities::MPI::ProcessGrid > grid
Definition scalapack.h:899
void load(const std::string &filename)
const int submatrix_column
Definition scalapack.h:980
const int submatrix_row
Definition scalapack.h:974
void mmult(ScaLAPACKMatrix< NumberType > &C, const ScaLAPACKMatrix< NumberType > &B, const bool adding=false) const
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()))
void save_serial(const std::string &filename, const std::pair< unsigned int, unsigned int > &chunk_size) const
NumberType norm_general(const char type) const
~ScaLAPACKMatrix() override=default
std::vector< NumberType > work
Definition scalapack.h:939
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
void load_parallel(const std::string &filename)
NumberType l1_norm() const
void compute_lu_factorization()
NumberType norm_symmetric(const char type) const
const char uplo
Definition scalapack.h:956
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
LAPACKSupport::Property property
Definition scalapack.h:892
void set_property(const LAPACKSupport::Property property)
Definition scalapack.cc:303
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:216
void load_serial(const std::string &filename)
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)
std::vector< int > ipiv
Definition scalapack.h:950
NumberType reciprocal_condition_number(const NumberType a_norm) const
LAPACKSupport::State state
Definition scalapack.h:886
void TmTmult(ScaLAPACKMatrix< NumberType > &C, const ScaLAPACKMatrix< NumberType > &B, const bool adding=false) const
Threads::Mutex mutex
Definition scalapack.h:985
std::vector< NumberType > eigenpairs_symmetric_by_index(const std::pair< unsigned int, unsigned int > &index_limits, const bool compute_eigenvectors)
unsigned int global_column(const unsigned int loc_column) const
Definition scalapack.cc:514
void copy_to(FullMatrix< NumberType > &matrix) const
Definition scalapack.cc:666
unsigned int global_row(const unsigned int loc_row) const
Definition scalapack.cc:497
unsigned int local_m() const
void compute_cholesky_factorization()
void copy_transposed(const ScaLAPACKMatrix< NumberType > &B)
Definition scalapack.cc:980
NumberType linfty_norm() const
void scale_columns(const InputVector &factors)
const TableIndices< N > & size() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
@ symmetric
Matrix is symmetric.
@ general
No special properties.
static const unsigned int invalid_unsigned_int
Definition types.h:220