CPPLapack
 All Classes Files Functions Variables Friends Pages
dgsmatrix-misc.hpp
Go to the documentation of this file.
1 //=============================================================================
2 /*! clear all the matrix data and set the sizes 0 */
3 inline void dgsmatrix::clear()
4 {CPPL_VERBOSE_REPORT;
5  m =0;
6  n =0;
7  data.clear();
8  rows.clear();
9  cols.clear();
10 }
11 
12 //=============================================================================
13 /*! change the matrix into a zero matrix */
15 {CPPL_VERBOSE_REPORT;
16  data.resize(0);
17  for(CPPL_INT i=0; i<m; i++){ rows[i].resize(0); }
18  for(CPPL_INT j=0; j<n; j++){ cols[j].resize(0); }
19  return *this;
20 }
21 
22 //=============================================================================
23 /*! change the matrix into an identity matrix */
25 {CPPL_VERBOSE_REPORT;
26 #ifdef CPPL_DEBUG
27  if(m!=n){
28  ERROR_REPORT;
29  std::cerr << "Only square matrix can be a identity matrix." << std::endl
30  << "The matrix size was " << m << "x" << n << "." << std::endl;
31  exit(1);
32  }
33 #endif//CPPL_DEBUG
34 
35  zero();
36  for(CPPL_INT i=0; i<m; i++){
37  put(i,i,1.);
38  }
39  return *this;
40 }
41 
42 //=============================================================================
43 /*! change sign(+/-) of the matrix */
44 inline void dgsmatrix::chsign()
45 {CPPL_VERBOSE_REPORT;
46  const std::vector<dcomponent>::iterator data_end =data.end();
47  for(std::vector<dcomponent>::iterator it=data.begin(); it!=data_end; it++){
48  it->v =-it->v;
49  }
50 }
51 
52 //=============================================================================
53 /*! make a deep copy of the matrix */
54 inline void dgsmatrix::copy(const dgsmatrix& mat)
55 {CPPL_VERBOSE_REPORT;
56  m =mat.m;
57  n =mat.n;
58  data =mat.data;
59  rows =mat.rows;
60  cols =mat.cols;
61 }
62 
63 //=============================================================================
64 /*! make a shallow copy of the matrix\n
65  This function is not designed to be used in project codes. */
66 inline void dgsmatrix::shallow_copy(const _dgsmatrix& mat)
67 {CPPL_VERBOSE_REPORT;
68  data.clear();
69  rows.clear();
70  cols.clear();
71 
72  m =mat.m;
73  n =mat.n;
74  data.swap(mat.data);
75  rows.swap(mat.rows);
76  cols.swap(mat.cols);
77 
78  mat.nullify();
79 }
80 
81 //=============================================================================
82 /*! resize the matrix */
83 inline dgsmatrix& dgsmatrix::resize(const CPPL_INT& _m, const CPPL_INT& _n, const CPPL_INT _c, const CPPL_INT _l)
84 {CPPL_VERBOSE_REPORT;
85 #ifdef CPPL_DEBUG
86  if( _m<0 || _n<0 || _c<0 ){
87  ERROR_REPORT;
88  std::cerr << "Matrix sizes and the length of arrays must be positive integers. " << std::endl
89  << "Your input was (" << _m << "," << _n << "," << _c << "," << _l << ")." << std::endl;
90  exit(1);
91  }
92 #endif//CPPL_DEBUG
93 
94  m =_m;
95  n =_n;
96  data.resize(0);
97  data.reserve(_c);
98  rows.resize(m);
99  for(CPPL_INT i=0; i<m; i++){
100  rows[i].resize(0);
101  rows[i].reserve(_l);
102  }
103  cols.resize(n);
104  for(CPPL_INT i=0; i<n; i++){
105  cols[i].resize(0);
106  cols[i].reserve(_l);
107  }
108 
109  return *this;
110 }
111 
112 //=============================================================================
113 /*! stretch the matrix size */
114 inline void dgsmatrix::stretch(const CPPL_INT& dm, const CPPL_INT& dn)
115 {CPPL_VERBOSE_REPORT;
116 #ifdef CPPL_DEBUG
117  if( m+dm<0 || n+dn<0 ){
118  ERROR_REPORT;
119  std::cerr << "The new matrix size must be larger than zero. " << std::endl
120  << "Your input was (" << dm << ", " << dn << ")." << std::endl;
121  exit(1);
122  }
123 #endif//CPPL_DEBUG
124 
125  //////// zero ////////
126  if(dm==0 && dn==0){ return; }
127 
128  //////// non-zero ////////
129  m +=dm;
130  n +=dn;
131 
132  //// for rows ////
133  if(dm<0){
134  //// delete components over the new size ////
135  const std::vector<dcomponent>::reverse_iterator data_rend =data.rend();
136  for(std::vector<dcomponent>::reverse_iterator it=data.rbegin(); it!=data_rend; it++){
137  if( it->i>=m ){ del( CPPL_INT(data_rend-it-1) ); }
138  }
139  //// shrink rows ////
140  for(CPPL_INT i=0; i<-dm; i++){
141  rows.pop_back();
142  }
143  }
144  else{//dm>0
145  //// expand rows ////
146  for(CPPL_INT i=0; i<dm; i++){
147  rows.push_back( std::vector<CPPL_INT>(0) );
148  }
149  }
150 
151  //// for cols ////
152  if(dn<0){
153  //// delete components over the new size ////
154  const std::vector<dcomponent>::reverse_iterator data_rend =data.rend();
155  for(std::vector<dcomponent>::reverse_iterator it=data.rbegin(); it!=data_rend; it++){
156  if( it->j>=n ){ del( CPPL_INT(data_rend-it-1) ); }
157  }
158  for(CPPL_INT j=0; j<-dn; j++){
159  cols.pop_back();
160  }
161  }
162  else{//dn>0
163  //// expand cols ////
164  for(CPPL_INT j=0; j<dn; j++){
165  cols.push_back( std::vector<CPPL_INT>(0) );
166  }
167  }
168 }
169 
170 //=============================================================================
171 /*! check if the component is listed */
172 inline bool dgsmatrix::isListed(const CPPL_INT& i, const CPPL_INT& j) const
173 {CPPL_VERBOSE_REPORT;
174 #ifdef CPPL_DEBUG
175  if( i<0 || j<0 || m<=i || n<=j ){
176  ERROR_REPORT;
177  std::cerr << "The required component is out of the matrix size." << std::endl
178  << "Your input is (" << i << "," << j << "), whereas the matrix size is " << m << "x" << n << "." << std::endl;
179  exit(1);
180  }
181 #endif//CPPL_DEBUG
182  const std::vector<CPPL_INT>::const_iterator rows_i_end =rows[i].end();
183  for(std::vector<CPPL_INT>::const_iterator p=rows[i].begin(); p!=rows_i_end; p++){
184  if(data[*p].j==j){ return 1; }
185  }
186 
187  return 0;
188 }
189 
190 
191 //=============================================================================
192 /*! return the element number of the component */
193 inline CPPL_INT dgsmatrix::number(const CPPL_INT& i, const CPPL_INT& j)
194 {CPPL_VERBOSE_REPORT;
195 #ifdef CPPL_DEBUG
196  if( i<0 || j<0 || m<=i || n<=j ){
197  ERROR_REPORT;
198  std::cerr << "The required component is out of the matrix size." << std::endl
199  << "Your input is (" << i << "," << j << "), whereas the matrix size is " << m << "x" << n << "." << std::endl;
200  exit(1);
201  }
202 #endif//CPPL_DEBUG
203  const std::vector<CPPL_INT>::iterator rows_i_end =rows[i].end();
204  for(std::vector<CPPL_INT>::iterator p=rows[i].begin(); p!=rows_i_end; p++){
205  if(data[*p].j==j){ return *p; }
206  }
207 
208  return -1;
209 }
210 
211 ///////////////////////////////////////////////////////////////////////////////
212 ///////////////////////////////////////////////////////////////////////////////
213 ///////////////////////////////////////////////////////////////////////////////
214 
215 //=============================================================================
216 /*! get row of the matrix */
217 inline _drovector dgsmatrix::row(const CPPL_INT& _m) const
218 {CPPL_VERBOSE_REPORT;
219 #ifdef CPPL_DEBUG
220  if( _m<0 || _m>m ){
221  ERROR_REPORT;
222  std::cerr << "Input row number must be between 0 and " << m << "." << std::endl
223  << "Your input was " << _m << "." << std::endl;
224  exit(1);
225  }
226 #endif//CPPL_DEBUG
227 
228  drovector vec(n);
229  vec.zero();
230 
231  const std::vector<CPPL_INT>::const_iterator rows__m_end =rows[_m].end();
232  for(std::vector<CPPL_INT>::const_iterator p=rows[_m].begin(); p!=rows__m_end; p++){
233  vec(data[*p].j) =data[*p].v;
234  }
235 
236  return _(vec);
237 }
238 
239 //=============================================================================
240 /*! get column of the matrix */
241 inline _dcovector dgsmatrix::col(const CPPL_INT& _n) const
242 {CPPL_VERBOSE_REPORT;
243 #ifdef CPPL_DEBUG
244  if( _n<0 || _n>n ){
245  ERROR_REPORT;
246  std::cerr << "Input row number must be between 0 and " << n << "." << std::endl
247  << "Your input was " << _n << "." << std::endl;
248  exit(1);
249  }
250 #endif//CPPL_DEBUG
251 
252  dcovector vec(m);
253  vec.zero();
254 
255  const std::vector<CPPL_INT>::const_iterator cols__n_end =cols[_n].end();
256  for(std::vector<CPPL_INT>::const_iterator p=cols[_n].begin(); p!=cols__n_end; p++){
257  vec(data[*p].i) =data[*p].v;
258  }
259 
260  return _(vec);
261 }
262 
263 //=============================================================================
264 /*! erase components less than DBL_MIN */
265 inline void dgsmatrix::diet(const double eps)
266 {CPPL_VERBOSE_REPORT;
267  const std::vector<dcomponent>::reverse_iterator data_rend =data.rend();
268  for(std::vector<dcomponent>::reverse_iterator it=data.rbegin(); it!=data_rend; it++){
269  if( fabs(it->v)<eps ){ del( CPPL_INT(data_rend-it-1) ); }
270  }
271 }
272 
273 ///////////////////////////////////////////////////////////////////////////////
274 ///////////////////////////////////////////////////////////////////////////////
275 ///////////////////////////////////////////////////////////////////////////////
276 
277 //=============================================================================
278 /*! health checkup */
279 inline void dgsmatrix::checkup()
280 {CPPL_VERBOSE_REPORT;
281  //////// write ////////
282  const std::vector<dcomponent>::const_iterator data_end =data.end();
283  for(std::vector<dcomponent>::const_iterator it=data.begin(); it!=data_end; it++){
284  std::cerr << "array[" << it-data.begin() << "] = (" << it->i << "," << it->j << ") = " << it->v << std::endl;
285  }
286  std::cerr << std::endl;
287 
288  for(CPPL_INT i=0; i<m; i++){
289  std::cerr << "rows[" << i << "] =" << std::flush;
290  const size_t rows_i_size =rows[i].size();
291  for(size_t k=0; k<rows_i_size; k++){
292  std::cerr << " " << rows[i][k] << std::flush;
293  }
294  std::cerr << std::endl;
295  }
296  std::cerr << std::endl;
297 
298  for(CPPL_INT j=0; j<n; j++){
299  std::cerr << "cols[" << j << "] =" << std::flush;
300  const size_t cols_j_size =cols[j].size();
301  for(size_t k=0; k<cols_j_size; k++){
302  std::cerr << " " << cols[j][k] << std::flush;
303  }
304  std::cerr << std::endl;
305  }
306 
307  //////// Elements ////////
308  for(std::vector<dcomponent>::const_iterator it=data.begin(); it!=data_end; it++){
309  //// m bound ////
310  if(it->i>=m){
311  ERROR_REPORT;
312  std::cerr << "The indx of the " << it-data.begin() << "th element is out of the matrix size." << std::endl
313  << "Its i index was " << it->i << "." << std::endl;
314  exit(1);
315  }
316 
317  //// n bound ////
318  if(it->j>=n){
319  ERROR_REPORT;
320  std::cerr << "The indx of the " << it-data.begin() << "th element is out of the matrix size." << std::endl
321  << "Its j index was " << it->j << "." << std::endl;
322  exit(1);
323  }
324 
325  //// double-listed ////
326  for(std::vector<dcomponent>::const_iterator IT=it+1; IT!=data_end; IT++){
327  if( it->i==IT->i && it->j==IT->j ){
328  ERROR_REPORT;
329  std::cerr << "The (" << it->i << ", " << it->j << ") component is double-listed at the " << it-data.begin() << "th and the" << IT-data.begin() << "the elements."<< std::endl;
330  exit(1);
331  }
332  }
333  }
334 
335  //////// ijc consistence ////////
336 
337 
338  //////// NOTE ////////
339  std::cerr << "# [NOTE]@dgsmatrix::checkup(): This sparse matrix is fine." << std::endl;
340 }
341 
342 ///////////////////////////////////////////////////////////////////////////////
343 ///////////////////////////////////////////////////////////////////////////////
344 ///////////////////////////////////////////////////////////////////////////////
345 
346 //=============================================================================
347 /*! swap two matrices */
348 inline void swap(dgsmatrix& A, dgsmatrix& B)
349 {CPPL_VERBOSE_REPORT;
350  std::swap(A.n,B.n);
351  std::swap(A.m,B.m);
352  std::swap(A.data,B.data);
353  std::swap(A.rows,B.rows);
354  std::swap(A.cols,B.cols);
355 }
356 
357 //=============================================================================
358 /*! convert user object to smart-temporary object */
359 inline _dgsmatrix _(dgsmatrix& mat)
360 {CPPL_VERBOSE_REPORT;
361  _dgsmatrix newmat;
362 
363  //////// shallow copy ////////
364  newmat.n =mat.n;
365  newmat.m =mat.m;
366  std::swap(newmat.data,mat.data);
367  std::swap(newmat.rows,mat.rows);
368  std::swap(newmat.cols,mat.cols);
369 
370  //////// nullify ////////
371  mat.m =0;
372  mat.n =0;
373 
374  return newmat;
375 }
_dgsmatrix _(dgsmatrix &mat)
void swap(dgsmatrix &A, dgsmatrix &B)
Real Double-precision General Sparse Matrix Class.
Definition: dgsmatrix.hpp:3
std::vector< std::vector< CPPL_INT > > cols
array of vector to store the entry information of component for each column
Definition: _dgsmatrix.hpp:13
void copy(const dgsmatrix &)
std::vector< dcomponent > data
matrix data
Definition: dgsmatrix.hpp:11
dgsmatrix & identity()
std::vector< std::vector< CPPL_INT > > rows
array of vector to store the entry information of component for each row
Definition: _dgsmatrix.hpp:12
CPPL_INT n
matrix column size
Definition: _dgsmatrix.hpp:10
void diet(const double=DBL_MIN)
CPPL_INT n
matrix column size
Definition: dgsmatrix.hpp:10
_drovector row(const CPPL_INT &) const
friend _dgsmatrix _(dgsmatrix &)
std::vector< std::vector< CPPL_INT > > cols
array of vector to store the entry information of component for each column
Definition: dgsmatrix.hpp:13
dcovector & zero()
_dgematrix i(const _dgbmatrix &mat)
(DO NOT USE) Smart-temporary Real Double-precision General Sparse Matrix Class
Definition: _dgsmatrix.hpp:3
Real Double-precision Row Vector Class.
Definition: drovector.hpp:3
bool isListed(const CPPL_INT &, const CPPL_INT &) const
CPPL_INT m
matrix row size
Definition: _dgsmatrix.hpp:9
drovector & zero()
(DO NOT USE) Smart-temporary Real Double-precision Row Vector Class
Definition: _drovector.hpp:3
dgsmatrix & zero()
void stretch(const CPPL_INT &, const CPPL_INT &)
dgsmatrix & del(const CPPL_INT, const CPPL_INT)
void nullify() const
std::vector< dcomponent > data
matrix data
Definition: _dgsmatrix.hpp:11
Real Double-precision Column Vector Class.
Definition: dcovector.hpp:3
dgsmatrix & put(const CPPL_INT &, const CPPL_INT &, const double &)
(DO NOT USE) Smart-temporary Real Double-precision Column Vector Class
Definition: _dcovector.hpp:3
dgsmatrix & resize(const CPPL_INT &, const CPPL_INT &, const CPPL_INT=0, const CPPL_INT=0)
std::vector< std::vector< CPPL_INT > > rows
array of vector to store the entry information of component for each row
Definition: dgsmatrix.hpp:12
_dcovector col(const CPPL_INT &) const
void shallow_copy(const _dgsmatrix &)
CPPL_INT m
matrix row size
Definition: dgsmatrix.hpp:9
void clear()
CPPL_INT number(const CPPL_INT &, const CPPL_INT &)