18#include <OMMathExceptions.H>
27 typedef std::pair<unsigned,unsigned> Index;
28 typedef std::map<Index,Matrix> Blocks;
35 size_t size()
const override {
37 for (
const auto&
block : blocks)
42 void info()
const override {
44 std::cout <<
"Empty matrix" << std::endl;
48 std::cout <<
"Symmetric block matrix" << std::endl;
49 std::cout <<
"Dimensions: " <<
nlin() <<
" x " <<
ncol() << std::endl;
50 std::cout <<
"Number of blocks: " << blocks.size() << std::endl;
51 std::cout <<
"Number of coefficients: " <<
size() << std::endl;
56 const Index& ind = find_block_indices(i,j,transposed);
62 const Matrix&
block(
const unsigned i,
const unsigned j)
const {
64 const Index& ind = find_block_indices(i,j,transposed);
67 return blocks.at(ind);
72 const Index& ind = create_block_indices(ir,jr,transposed);
80 for (
unsigned i=0; i<r.size(); ++i)
81 for (
unsigned j=i; j<r.size(); ++j)
87 const Index& ind = find_block_indices(i,j,transposed);
88 const size_t inblockindex_i = ((!transposed) ? i : j)-ranges[ind.first].start();
89 const size_t inblockindex_j = ((!transposed) ? j : i)-ranges[ind.second].start();
90 return blocks[ind](inblockindex_i,inblockindex_j);
95 const Index& ind = find_block_indices(i,j,transposed);
96 const size_t inblockindex_i = ((!transposed) ? i : j)-ranges[ind.first].start();
97 const size_t inblockindex_j = ((!transposed) ? j : i)-ranges[ind.second].start();
98 return blocks.at(ind)(inblockindex_i,inblockindex_j);
103 unsigned find_block_index(
const Range& r)
const {
return ranges.
find_index(r); }
104 unsigned find_block_index(
const unsigned i)
const {
return ranges.
find_index(i); }
106 unsigned create_block_index(
const Range& r)
try {
111 return ranges.size()-1;
114 Index create_block_indices(
const Range& ir,
const Range& jr,
bool& transposed) {
115 transposed = ir.start()>jr.start();
116 const unsigned iind = create_block_index(ir);
117 const unsigned jind = create_block_index(jr);
118 return (transposed) ?
Index({jind,iind}) :
Index({iind,jind});
121 Index find_block_indices(
const Range& ir,
const Range& jr,
bool& transposed)
const {
122 transposed = ir.start()>jr.start();
123 const unsigned iind = find_block_index(ir);
124 const unsigned jind = find_block_index(jr);
125 return (transposed) ?
Index({jind,iind}) :
Index({iind,jind});
128 Index find_block_indices(
const unsigned i,
const unsigned j,
bool& transposed)
const {
130 const unsigned iind = find_block_index(i);
131 const unsigned jind = find_block_index(j);
132 return (transposed) ?
Index({jind,iind}) :
Index({iind,jind});
virtual Dimension ncol() const
Matrix class Matrix class.
size_t size() const
Get Matrix size.
unsigned find_index(const size_t ind) const
Block symmetric matrix class Block symmetric matrix class.
void add_block(const Range &ir, const Range &jr)
const Matrix & block(const unsigned i, const unsigned j) const
Matrix & block(const unsigned i, const unsigned j)
double & operator()(const size_t i, const size_t j)
double operator()(const size_t i, const size_t j) const
void set_blocks(const Ranges &r)
size_t size() const override
void info() const override
SymmetricBlockMatrix(const size_t N)