Matrix_Sparse.cxx
1 // Copyright (C) 2001-2011 Vivien Mallet
2 // Copyright (C) 2003-2011 Marc Duruflé
3 //
4 // This file is part of the linear-algebra library Seldon,
5 // http://seldon.sourceforge.net/.
6 //
7 // Seldon is free software; you can redistribute it and/or modify it under the
8 // terms of the GNU Lesser General Public License as published by the Free
9 // Software Foundation; either version 2.1 of the License, or (at your option)
10 // any later version.
11 //
12 // Seldon is distributed in the hope that it will be useful, but WITHOUT ANY
13 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15 // more details.
16 //
17 // You should have received a copy of the GNU Lesser General Public License
18 // along with Seldon. If not, see http://www.gnu.org/licenses/.
19 
20 
21 #ifndef SELDON_FILE_MATRIX_SPARSE_CXX
22 
23 #include "Matrix_Sparse.hxx"
24 
25 #include <set>
26 
27 /*
28  Functions defined in this file:
29 
30  reads coordinate form (i, j, val) from an input stream
31  ReadCoordinateMatrix(FileStream, row_index, col_index, values, cplx)
32 
33  writes coordinate form (i, j, val) in an input stream
34  WriteCoordinateMatrix(FileStream, row_index, col_index, values, cplx)
35 
36  These functions are used by ReadText/WriteText
37 */
38 
39 namespace Seldon
40 {
41 
42 
43  /*********************
44  * MEMORY MANAGEMENT *
45  *********************/
46 
47 
49 
60  template <class T, class Prop, class Storage, class Allocator>
61  template <class Storage0, class Allocator0,
62  class Storage1, class Allocator1,
63  class Storage2, class Allocator2>
64  Matrix_Sparse<T, Prop, Storage, Allocator>::
65  Matrix_Sparse(int i, int j,
66  Vector<T, Storage0, Allocator0>& values,
67  Vector<long, Storage1, Allocator1>& ptr,
68  Vector<int, Storage2, Allocator2>& ind):
69  Matrix_Base<T, Allocator>(i, j)
70  {
71  nz_ = values.GetLength();
72 
73 #ifdef SELDON_CHECK_DIMENSIONS
74  // Checks whether vector sizes are acceptable.
75 
76  if (ind.GetLength() != nz_)
77  {
78  this->m_ = 0;
79  this->n_ = 0;
80  nz_ = 0;
81  ptr_ = NULL;
82  ind_ = NULL;
83  this->data_ = NULL;
84  throw WrongDim(string("Matrix_Sparse::Matrix_Sparse(int, int, ")
85  + "const Vector&, const Vector&, const Vector&)",
86  string("There are ") + to_str(nz_) + " values but "
87  + to_str(ind.GetLength()) + " row or column indices.");
88  }
89 
90  if (ptr.GetLength()-1 != Storage::GetFirst(i, j))
91  {
92  this->m_ = 0;
93  this->n_ = 0;
94  nz_ = 0;
95  ptr_ = NULL;
96  ind_ = NULL;
97  this->data_ = NULL;
98  throw WrongDim(string("Matrix_Sparse::Matrix_Sparse(int, int, ")
99  + "const Vector&, const Vector&, const Vector&)",
100  string("The vector of start indices contains ")
101  + to_str(ptr.GetLength()-1) + string(" row or column")
102  + string(" start indices (plus the number")
103  + " of non-zero entries) but there are "
104  + to_str(Storage::GetFirst(i, j))
105  + " rows or columns (" + to_str(i) + " by "
106  + to_str(j) + " matrix).");
107  }
108 
109  if (nz_ > 0
110  && (j == 0
111  || static_cast<long int>(nz_-1) / static_cast<long int>(j)
112  >= static_cast<long int>(i)))
113  {
114  this->m_ = 0;
115  this->n_ = 0;
116  nz_ = 0;
117  ptr_ = NULL;
118  ind_ = NULL;
119  this->data_ = NULL;
120  throw WrongDim(string("Matrix_Sparse::Matrix_Sparse(int, int, ")
121  + "const Vector&, const Vector&, const Vector&)",
122  string("There are more values (")
123  + to_str(values.GetLength())
124  + " values) than elements in the matrix ("
125  + to_str(i) + " by " + to_str(j) + ").");
126  }
127 #endif
128 
129  this->ptr_ = ptr.GetData();
130  this->ind_ = ind.GetData();
131  this->data_ = values.GetData();
132 
133  ptr.Nullify();
134  ind.Nullify();
135  values.Nullify();
136  }
137 
138 
140  template <class T, class Prop, class Storage, class Allocator>
141  Matrix_Sparse<T, Prop, Storage, Allocator>::
142  Matrix_Sparse(const Matrix_Sparse<T, Prop, Storage, Allocator>& A):
143  Matrix_Base<T, Allocator>(A)
144  {
145  this->m_ = 0;
146  this->n_ = 0;
147  this->nz_ = 0;
148  ptr_ = NULL;
149  ind_ = NULL;
150  this->Copy(A);
151  }
152 
153 
155 
158  template <class T, class Prop, class Storage, class Allocator>
159  void Matrix_Sparse<T, Prop, Storage, Allocator>::Clear()
160  {
161 #ifdef SELDON_CHECK_MEMORY
162  try
163  {
164 #endif
165 
166  if (ptr_ != NULL)
167  {
168  AllocatorLong::deallocate(ptr_, this->m_+1);
169  ptr_ = NULL;
170  }
171 
172 #ifdef SELDON_CHECK_MEMORY
173  }
174  catch (...)
175  {
176  ptr_ = NULL;
177  }
178 #endif
179 
180 #ifdef SELDON_CHECK_MEMORY
181  try
182  {
183 #endif
184 
185  if (ind_ != NULL)
186  {
187  AllocatorInt::deallocate(ind_, this->nz_);
188  ind_ = NULL;
189  }
190 
191 #ifdef SELDON_CHECK_MEMORY
192  }
193  catch (...)
194  {
195  ind_ = NULL;
196  }
197 #endif
198 
199 #ifdef SELDON_CHECK_MEMORY
200  try
201  {
202 #endif
203 
204  if (this->data_ != NULL)
205  {
206  Allocator::deallocate(this->data_, nz_);
207  this->data_ = NULL;
208  }
209 
210 #ifdef SELDON_CHECK_MEMORY
211  }
212  catch (...)
213  {
214  this->nz_ = 0;
215  this->data_ = NULL;
216  }
217 #endif
218 
219  this->m_ = 0;
220  this->n_ = 0;
221  this->nz_ = 0;
222  }
223 
224 
226 
236  template <class T, class Prop, class Storage, class Allocator>
237  template <class Storage0, class Allocator0,
238  class Storage1, class Allocator1,
239  class Storage2, class Allocator2>
240  void Matrix_Sparse<T, Prop, Storage, Allocator>::
241  SetData(int i, int j,
242  Vector<T, Storage0, Allocator0>& values,
243  Vector<long, Storage1, Allocator1>& ptr,
244  Vector<int, Storage2, Allocator2>& ind)
245  {
246  this->Clear();
247  this->m_ = i;
248  this->n_ = j;
249  this->nz_ = values.GetLength();
250 
251 #ifdef SELDON_CHECK_DIMENSIONS
252  // Checks whether vector sizes are acceptable.
253 
254  if (ind.GetM() != nz_)
255  {
256  this->m_ = 0;
257  this->n_ = 0;
258  nz_ = 0;
259  ptr_ = NULL;
260  ind_ = NULL;
261  this->data_ = NULL;
262  throw WrongDim(string("Matrix_Sparse::SetData(int, int, ")
263  + "const Vector&, const Vector&, const Vector&)",
264  string("There are ") + to_str(nz_) + " values but "
265  + to_str(ind.GetLength()) + " row or column indices.");
266  }
267 
268  if (ptr.GetM()-1 != Storage::GetFirst(i, j))
269  {
270  this->m_ = 0;
271  this->n_ = 0;
272  nz_ = 0;
273  ptr_ = NULL;
274  ind_ = NULL;
275  this->data_ = NULL;
276  throw WrongDim(string("Matrix_Sparse::SetData(int, int, ")
277  + "const Vector&, const Vector&, const Vector&)",
278  string("The vector of start indices contains ")
279  + to_str(ptr.GetLength()-1) + string(" row or column")
280  + string(" start indices (plus the number")
281  + " of non-zero entries) but there are "
282  + to_str(Storage::GetFirst(i, j))
283  + " rows or columns (" + to_str(i) + " by "
284  + to_str(j) + " matrix).");
285  }
286 
287  if (nz_ > 0
288  && (j == 0
289  || static_cast<long int>(nz_-1) / static_cast<long int>(j)
290  >= static_cast<long int>(i)))
291  {
292  this->m_ = 0;
293  this->n_ = 0;
294  nz_ = 0;
295  ptr_ = NULL;
296  ind_ = NULL;
297  this->data_ = NULL;
298  throw WrongDim(string("Matrix_Sparse::SetData(int, int, ")
299  + "const Vector&, const Vector&, const Vector&)",
300  string("There are more values (")
301  + to_str(values.GetLength())
302  + " values) than elements in the matrix ("
303  + to_str(i) + " by " + to_str(j) + ").");
304  }
305 #endif
306 
307  this->ptr_ = ptr.GetData();
308  this->ind_ = ind.GetData();
309  this->data_ = values.GetData();
310 
311  ptr.Nullify();
312  ind.Nullify();
313  values.Nullify();
314  }
315 
316 
318 
331  template <class T, class Prop, class Storage, class Allocator>
332  void Matrix_Sparse<T, Prop, Storage, Allocator>
333  ::SetData(int i, int j, long nz,
334  typename Matrix_Sparse<T, Prop, Storage, Allocator>
335  ::pointer values, long* ptr, int* ind)
336  {
337  this->Clear();
338 
339  this->m_ = i;
340  this->n_ = j;
341 
342  this->nz_ = nz;
343 
344  this->data_ = values;
345  ind_ = ind;
346  ptr_ = ptr;
347  }
348 
349 
351 
355  template <class T, class Prop, class Storage, class Allocator>
356  void Matrix_Sparse<T, Prop, Storage, Allocator>::Nullify()
357  {
358  this->data_ = NULL;
359  this->m_ = 0;
360  this->n_ = 0;
361  nz_ = 0;
362  ptr_ = NULL;
363  ind_ = NULL;
364  }
365 
366 
368 
372  template <class T, class Prop, class Storage, class Allocator>
373  void Matrix_Sparse<T, Prop, Storage, Allocator>::Reallocate(int i, int j)
374  {
375  // clearing previous entries
376  Clear();
377 
378  this->m_ = i;
379  this->n_ = j;
380 
381  // we try to allocate ptr_
382 #ifdef SELDON_CHECK_MEMORY
383  try
384  {
385 #endif
386 
387  ptr_ = reinterpret_cast<long*>( AllocatorLong::
388  allocate(Storage::GetFirst(i, j)+1, this) );
389 
390 #ifdef SELDON_CHECK_MEMORY
391  }
392  catch (...)
393  {
394  this->m_ = 0;
395  this->n_ = 0;
396  nz_ = 0;
397  ptr_ = NULL;
398  ind_ = NULL;
399  this->data_ = NULL;
400  }
401  if (ptr_ == NULL)
402  {
403  this->m_ = 0;
404  this->n_ = 0;
405  nz_ = 0;
406  ind_ = NULL;
407  this->data_ = NULL;
408  }
409  if (ptr_ == NULL && i != 0 && j != 0)
410  throw NoMemory("Matrix_Sparse::Reallocate(int, int)",
411  string("Unable to allocate ")
412  + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1) )
413  + " bytes to store " + to_str(Storage::GetFirst(i, j)+1)
414  + " row or column start indices, for a "
415  + to_str(i) + " by " + to_str(j) + " matrix.");
416 #endif
417 
418  // then filing ptr_ with 0
419  for (int k = 0; k <= Storage::GetFirst(i, j); k++)
420  ptr_[k] = 0;
421 
422  }
423 
424 
426 
431  template <class T, class Prop, class Storage, class Allocator>
432  void Matrix_Sparse<T, Prop, Storage, Allocator>::Reallocate(int i, int j, long nz)
433  {
434  // clearing previous entries
435  Clear();
436 
437  this->nz_ = nz;
438  this->m_ = i;
439  this->n_ = j;
440 
441 #ifdef SELDON_CHECK_DIMENSIONS
442  if (nz_ < 0)
443  {
444  this->m_ = 0;
445  this->n_ = 0;
446  nz_ = 0;
447  ptr_ = NULL;
448  ind_ = NULL;
449  this->data_ = NULL;
450  throw WrongDim("Matrix_Sparse::Matrix_Sparse(int, int, int)",
451  "Invalid number of non-zero elements: " + to_str(nz)
452  + ".");
453  }
454  if (nz_ > 0
455  && (j == 0
456  || static_cast<long int>(nz_-1) / static_cast<long int>(j)
457  >= static_cast<long int>(i)))
458  {
459  this->m_ = 0;
460  this->n_ = 0;
461  nz_ = 0;
462  ptr_ = NULL;
463  ind_ = NULL;
464  this->data_ = NULL;
465  throw WrongDim("Matrix_Sparse::Matrix_Sparse(int, int, int)",
466  string("There are more values (") + to_str(nz)
467  + " values) than elements in the matrix ("
468  + to_str(i) + " by " + to_str(j) + ").");
469  }
470 #endif
471 
472 #ifdef SELDON_CHECK_MEMORY
473  try
474  {
475 #endif
476 
477  ptr_ = reinterpret_cast<long*>( AllocatorLong::
478  allocate(Storage::GetFirst(i, j)+1, this) );
479 
480 #ifdef SELDON_CHECK_MEMORY
481  }
482  catch (...)
483  {
484  this->m_ = 0;
485  this->n_ = 0;
486  nz_ = 0;
487  ptr_ = NULL;
488  ind_ = NULL;
489  this->data_ = NULL;
490  }
491  if (ptr_ == NULL)
492  {
493  this->m_ = 0;
494  this->n_ = 0;
495  nz_ = 0;
496  ind_ = NULL;
497  this->data_ = NULL;
498  }
499  if (ptr_ == NULL && i != 0 && j != 0)
500  throw NoMemory("Matrix_Sparse::Matrix_Sparse(int, int, int)",
501  string("Unable to allocate ")
502  + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1) )
503  + " bytes to store " + to_str(Storage::GetFirst(i, j)+1)
504  + " row or column start indices, for a "
505  + to_str(i) + " by " + to_str(j) + " matrix.");
506 #endif
507 
508 #ifdef SELDON_CHECK_MEMORY
509  try
510  {
511 #endif
512 
513  ind_ = reinterpret_cast<int*>( AllocatorInt::
514  allocate(nz_, this) );
515 
516 #ifdef SELDON_CHECK_MEMORY
517  }
518  catch (...)
519  {
520  this->m_ = 0;
521  this->n_ = 0;
522  nz_ = 0;
523  AllocatorLong::
524  deallocate(ptr_, Storage::GetFirst(i, j)+1);
525  ptr_ = NULL;
526  ind_ = NULL;
527  this->data_ = NULL;
528  }
529  if (ind_ == NULL)
530  {
531  this->m_ = 0;
532  this->n_ = 0;
533  nz_ = 0;
534  AllocatorLong::
535  deallocate(ptr_, Storage::GetFirst(i, j)+1);
536  ptr_ = NULL;
537  this->data_ = NULL;
538  }
539  if (ind_ == NULL && i != 0 && j != 0)
540  throw NoMemory("Matrix_Sparse::Matrix_Sparse(int, int, int)",
541  string("Unable to allocate ") + to_str(sizeof(int) * nz)
542  + " bytes to store " + to_str(nz)
543  + " row or column indices, for a "
544  + to_str(i) + " by " + to_str(j) + " matrix.");
545 #endif
546 
547 #ifdef SELDON_CHECK_MEMORY
548  try
549  {
550 #endif
551 
552  this->data_ = Allocator::allocate(nz_, this);
553 
554 #ifdef SELDON_CHECK_MEMORY
555  }
556  catch (...)
557  {
558  this->m_ = 0;
559  this->n_ = 0;
560  AllocatorLong::
561  deallocate(ptr_, Storage::GetFirst(i, j)+1);
562  ptr_ = NULL;
563  AllocatorInt::deallocate(ind_, nz);
564  ind_ = NULL;
565  this->data_ = NULL;
566  }
567  if (this->data_ == NULL)
568  {
569  this->m_ = 0;
570  this->n_ = 0;
571  AllocatorLong::
572  deallocate(ptr_, Storage::GetFirst(i, j)+1);
573  ptr_ = NULL;
574  AllocatorInt::deallocate(ind_, nz);
575  ind_ = NULL;
576  }
577  if (this->data_ == NULL && i != 0 && j != 0)
578  throw NoMemory("Matrix_Sparse::Matrix_Sparse(int, int, int)",
579  string("Unable to allocate ") + to_str(sizeof(int) * nz)
580  + " bytes to store " + to_str(nz) + " values, for a "
581  + to_str(i) + " by " + to_str(j) + " matrix.");
582 #endif
583 
584  for (int k = 0; k <= Storage::GetFirst(i, j); k++)
585  ptr_[k] = 0;
586  }
587 
588 
590 
595  template <class T, class Prop, class Storage, class Allocator>
596  void Matrix_Sparse<T, Prop, Storage, Allocator>::Resize(int i, int j)
597  {
598  if (Storage::GetFirst(i, j) < Storage::GetFirst(this->m_, this->n_))
599  Resize(i, j, ptr_[Storage::GetFirst(i, j)]);
600  else
601  Resize(i, j, nz_);
602  }
603 
604 
606 
612  template <class T, class Prop, class Storage, class Allocator>
613  void Matrix_Sparse<T, Prop, Storage, Allocator>
614  ::Resize(int i, int j, long nz)
615  {
616 #ifdef SELDON_CHECK_DIMENSIONS
617  if (nz < 0)
618  {
619  this->m_ = 0;
620  this->n_ = 0;
621  nz_ = 0;
622  ptr_ = NULL;
623  ind_ = NULL;
624  this->data_ = NULL;
625  throw WrongDim("Matrix_Sparse::Resize(int, int, int)",
626  "Invalid number of non-zero elements: " + to_str(nz)
627  + ".");
628  }
629  if (nz > 0
630  && (j == 0
631  || static_cast<long int>(nz_-1) / static_cast<long int>(j)
632  >= static_cast<long int>(i)))
633  {
634  this->m_ = 0;
635  this->n_ = 0;
636  nz_ = 0;
637  ptr_ = NULL;
638  ind_ = NULL;
639  this->data_ = NULL;
640  throw WrongDim("Matrix_Sparse::Resize(int, int, int)",
641  string("There are more values (") + to_str(nz)
642  + " values) than elements in the matrix ("
643  + to_str(i) + " by " + to_str(j) + ").");
644  }
645 #endif
646 
647  if (nz != nz_)
648  {
649  // trying to resize ind_ and data_
650 #ifdef SELDON_CHECK_MEMORY
651  try
652  {
653 #endif
654  ind_
655  = reinterpret_cast<int*>( AllocatorInt::
656  reallocate(ind_, nz, this) );
657 
658 #ifdef SELDON_CHECK_MEMORY
659  }
660  catch (...)
661  {
662  this->m_ = 0;
663  this->n_ = 0;
664  nz_ = 0;
665  AllocatorLong::
666  deallocate(ptr_, Storage::GetFirst(i, j)+1);
667  ptr_ = NULL;
668  ind_ = NULL;
669  this->data_ = NULL;
670  }
671  if (ind_ == NULL)
672  {
673  this->m_ = 0;
674  this->n_ = 0;
675  nz_ = 0;
676  AllocatorLong::
677  deallocate(ptr_, Storage::GetFirst(i, j)+1);
678  ptr_ = NULL;
679  this->data_ = NULL;
680  }
681  if (ind_ == NULL && i != 0 && j != 0)
682  throw NoMemory("Matrix_Sparse::Resize(int, int, int)",
683  string("Unable to allocate ") + to_str(sizeof(int) * nz)
684  + " bytes to store " + to_str(nz)
685  + " row or column indices, for a "
686  + to_str(i) + " by " + to_str(j) + " matrix.");
687 #endif
688 
689  Vector<T, VectFull, Allocator> val;
690  val.SetData(nz_, this->data_);
691  val.Resize(nz);
692 
693  this->data_ = val.GetData();
694  nz_ = nz;
695  val.Nullify();
696  }
697 
698  if (Storage::GetFirst(this->m_, this->n_) != Storage::GetFirst(i, j))
699  {
700 #ifdef SELDON_CHECK_MEMORY
701  try
702  {
703 #endif
704  // trying to resize ptr_
705  ptr_ = reinterpret_cast<long*>( AllocatorLong::
706  reallocate(ptr_, Storage::GetFirst(i, j)+1) );
707 
708 #ifdef SELDON_CHECK_MEMORY
709  }
710  catch (...)
711  {
712  this->m_ = 0;
713  this->n_ = 0;
714  nz_ = 0;
715  ptr_ = NULL;
716  ind_ = NULL;
717  this->data_ = NULL;
718  }
719  if (ptr_ == NULL)
720  {
721  this->m_ = 0;
722  this->n_ = 0;
723  nz_ = 0;
724  ind_ = NULL;
725  this->data_ = NULL;
726  }
727  if (ptr_ == NULL && i != 0 && j != 0)
728  throw NoMemory("Matrix_Sparse::Resize(int, int)",
729  string("Unable to allocate ")
730  + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1) )
731  + " bytes to store " + to_str(Storage::GetFirst(i, j)+1)
732  + " row or column start indices, for a "
733  + to_str(i) + " by " + to_str(j) + " matrix.");
734 #endif
735 
736  // then filing last values of ptr_ with nz_
737  for (int k = Storage::GetFirst(this->m_, this->n_);
738  k <= Storage::GetFirst(i, j); k++)
739  ptr_[k] = this->nz_;
740  }
741 
742  this->m_ = i;
743  this->n_ = j;
744  }
745 
746 
748  template <class T, class Prop, class Storage, class Allocator>
749  void Matrix_Sparse<T, Prop, Storage, Allocator>::
750  Copy(const Matrix_Sparse<T, Prop, Storage, Allocator>& A)
751  {
752  this->Clear();
753  int nz = A.nz_;
754  int i = A.m_;
755  int j = A.n_;
756  this->nz_ = nz;
757  this->m_ = i;
758  this->n_ = j;
759  if ((i == 0)||(j == 0))
760  {
761  this->m_ = 0;
762  this->n_ = 0;
763  this->nz_ = 0;
764  return;
765  }
766 
767 #ifdef SELDON_CHECK_DIMENSIONS
768  if (nz_ > 0
769  && (j == 0
770  || static_cast<long int>(nz_-1) / static_cast<long int>(j)
771  >= static_cast<long int>(i)))
772  {
773  this->m_ = 0;
774  this->n_ = 0;
775  nz_ = 0;
776  ptr_ = NULL;
777  ind_ = NULL;
778  this->data_ = NULL;
779  throw WrongDim("Matrix_Sparse::Matrix_Sparse(int, int, int)",
780  string("There are more values (") + to_str(nz)
781  + " values) than elements in the matrix ("
782  + to_str(i) + " by " + to_str(j) + ").");
783  }
784 #endif
785 
786 #ifdef SELDON_CHECK_MEMORY
787  try
788  {
789 #endif
790 
791  ptr_ = reinterpret_cast<long*>( AllocatorLong::
792  allocate(Storage::GetFirst(i, j)+1) );
793 
794  AllocatorLong::memorycpy(this->ptr_, A.ptr_,
795  Storage::GetFirst(i, j) + 1);
796 #ifdef SELDON_CHECK_MEMORY
797  }
798  catch (...)
799  {
800  this->m_ = 0;
801  this->n_ = 0;
802  nz_ = 0;
803  ptr_ = NULL;
804  ind_ = NULL;
805  this->data_ = NULL;
806  }
807  if (ptr_ == NULL)
808  {
809  this->m_ = 0;
810  this->n_ = 0;
811  nz_ = 0;
812  ind_ = NULL;
813  this->data_ = NULL;
814  }
815  if (ptr_ == NULL && i != 0 && j != 0)
816  throw NoMemory("Matrix_Sparse::Matrix_Sparse(int, int, int)",
817  string("Unable to allocate ")
818  + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1) )
819  + " bytes to store " + to_str(Storage::GetFirst(i, j)+1)
820  + " row or column start indices, for a "
821  + to_str(i) + " by " + to_str(j) + " matrix.");
822 #endif
823 
824 #ifdef SELDON_CHECK_MEMORY
825  try
826  {
827 #endif
828 
829  ind_ = reinterpret_cast<int*>( AllocatorInt::
830  allocate(nz_, this) );
831  AllocatorInt::memorycpy(this->ind_, A.ind_, nz_);
832 
833 #ifdef SELDON_CHECK_MEMORY
834  }
835  catch (...)
836  {
837  this->m_ = 0;
838  this->n_ = 0;
839  nz_ = 0;
840  AllocatorLong::
841  deallocate(ptr_, Storage::GetFirst(i, j)+1);
842  ptr_ = NULL;
843  ind_ = NULL;
844  this->data_ = NULL;
845  }
846  if (ind_ == NULL)
847  {
848  this->m_ = 0;
849  this->n_ = 0;
850  nz_ = 0;
851  AllocatorLong::
852  deallocate(ptr_, Storage::GetFirst(i, j)+1);
853  ptr_ = NULL;
854  this->data_ = NULL;
855  }
856  if (ind_ == NULL && i != 0 && j != 0)
857  throw NoMemory("Matrix_Sparse::Matrix_Sparse(int, int, int)",
858  string("Unable to allocate ") + to_str(sizeof(int) * nz)
859  + " bytes to store " + to_str(nz)
860  + " row or column indices, for a "
861  + to_str(i) + " by " + to_str(j) + " matrix.");
862 #endif
863 
864 #ifdef SELDON_CHECK_MEMORY
865  try
866  {
867 #endif
868 
869  this->data_ = Allocator::allocate(nz_, this);
870  Allocator::memorycpy(this->data_, A.data_, nz_);
871 
872 #ifdef SELDON_CHECK_MEMORY
873  }
874  catch (...)
875  {
876  this->m_ = 0;
877  this->n_ = 0;
878  AllocatorLong::
879  deallocate(ptr_, Storage::GetFirst(i, j)+1);
880  ptr_ = NULL;
881  AllocatorInt::deallocate(ind_, nz);
882  ind_ = NULL;
883  this->data_ = NULL;
884  }
885  if (this->data_ == NULL)
886  {
887  this->m_ = 0;
888  this->n_ = 0;
889  AllocatorLong::
890  deallocate(ptr_, Storage::GetFirst(i, j)+1);
891  ptr_ = NULL;
892  AllocatorInt::deallocate(ind_, nz);
893  ind_ = NULL;
894  }
895  if (this->data_ == NULL && i != 0 && j != 0)
896  throw NoMemory("Matrix_Sparse::Matrix_Sparse(int, int, int)",
897  string("Unable to allocate ") + to_str(sizeof(int) * nz)
898  + " bytes to store " + to_str(nz) + " values, for a "
899  + to_str(i) + " by " + to_str(j) + " matrix.");
900 #endif
901 
902  }
903 
904 
906  template<class T, class Prop, class Storage, class Allocator>
907  size_t Matrix_Sparse<T, Prop, Storage, Allocator>::GetMemorySize() const
908  {
909  size_t taille = sizeof(*this) + this->GetPtrSize()*sizeof(long);
910  int coef = sizeof(T) + sizeof(int); // for each non-zero entry
911  taille += coef*size_t(this->nz_);
912  return taille;
913  }
914 
915 
917  template<class T, class Prop, class Storage, class Allocator>
918  void Matrix_Sparse<T, Prop, Storage, Allocator>::FillPtrInt(Vector<int>& Ptr) const
919  {
920  Ptr.Reallocate(this->m_+1);
921  for (int i = 0; i <= this->m_; i++)
922  Ptr(i) = ptr_[i];
923  }
924 
925 
926  /**********************************
927  * ELEMENT ACCESS AND AFFECTATION *
928  **********************************/
929 
930 
932 
938  template <class T, class Prop, class Storage, class Allocator>
939  const typename Matrix_Sparse<T, Prop, Storage, Allocator>::value_type
940  Matrix_Sparse<T, Prop, Storage, Allocator>::operator() (int i, int j) const
941  {
942 
943 #ifdef SELDON_CHECK_BOUNDS
944  if (i < 0 || i >= this->m_)
945  throw WrongRow("Matrix_Sparse::operator()",
946  string("Index should be in [0, ") + to_str(this->m_-1)
947  + "], but is equal to " + to_str(i) + ".");
948  if (j < 0 || j >= this->n_)
949  throw WrongCol("Matrix_Sparse::operator()",
950  string("Index should be in [0, ") + to_str(this->n_-1)
951  + "], but is equal to " + to_str(j) + ".");
952 #endif
953 
954  long k, l;
955  long a, b;
956  T zero;
957  SetComplexZero(zero);
958 
959  a = ptr_[Storage::GetFirst(i, j)];
960  b = ptr_[Storage::GetFirst(i, j) + 1];
961 
962  if (a == b)
963  return zero;
964 
965  l = Storage::GetSecond(i, j);
966 
967  for (k = a; (k < b-1) && (ind_[k] < l); k++);
968 
969  if (ind_[k] == l)
970  return this->data_[k];
971  else
972  return zero;
973  }
974 
975 
977 
985  template <class T, class Prop, class Storage, class Allocator>
986  typename Matrix_Sparse<T, Prop, Storage, Allocator>::value_type&
987  Matrix_Sparse<T, Prop, Storage, Allocator>::Val(int i, int j)
988  {
989 
990 #ifdef SELDON_CHECK_BOUNDS
991  if (i < 0 || i >= this->m_)
992  throw WrongRow("Matrix_Sparse::Val(int, int)",
993  string("Index should be in [0, ") + to_str(this->m_-1)
994  + "], but is equal to " + to_str(i) + ".");
995  if (j < 0 || j >= this->n_)
996  throw WrongCol("Matrix_Sparse::Val(int, int)",
997  string("Index should be in [0, ") + to_str(this->n_-1)
998  + "], but is equal to " + to_str(j) + ".");
999 #endif
1000 
1001  long k, l;
1002  long a, b;
1003 
1004  a = ptr_[Storage::GetFirst(i, j)];
1005  b = ptr_[Storage::GetFirst(i, j) + 1];
1006 
1007  if (a == b)
1008  throw WrongArgument("Matrix_Sparse::Val(int, int)",
1009  "No reference to element (" + to_str(i) + ", "
1010  + to_str(j)
1011  + ") can be returned: it is a zero entry.");
1012 
1013  l = Storage::GetSecond(i, j);
1014 
1015  for (k = a; (k < b-1) && (ind_[k] < l); k++);
1016 
1017  if (ind_[k] == l)
1018  return this->data_[k];
1019  else
1020  throw WrongArgument("Matrix_Sparse::Val(int, int)",
1021  "No reference to element (" + to_str(i) + ", "
1022  + to_str(j)
1023  + ") can be returned: it is a zero entry.");
1024  }
1025 
1026 
1028 
1036  template <class T, class Prop, class Storage, class Allocator>
1037  const typename Matrix_Sparse<T, Prop, Storage, Allocator>::value_type&
1038  Matrix_Sparse<T, Prop, Storage, Allocator>::Val(int i, int j) const
1039  {
1040 
1041 #ifdef SELDON_CHECK_BOUNDS
1042  if (i < 0 || i >= this->m_)
1043  throw WrongRow("Matrix_Sparse::Val(int, int)",
1044  string("Index should be in [0, ") + to_str(this->m_-1)
1045  + "], but is equal to " + to_str(i) + ".");
1046  if (j < 0 || j >= this->n_)
1047  throw WrongCol("Matrix_Sparse::Val(int, int)",
1048  string("Index should be in [0, ") + to_str(this->n_-1)
1049  + "], but is equal to " + to_str(j) + ".");
1050 #endif
1051 
1052  long k, l;
1053  long a, b;
1054 
1055  a = ptr_[Storage::GetFirst(i, j)];
1056  b = ptr_[Storage::GetFirst(i, j) + 1];
1057 
1058  if (a == b)
1059  throw WrongArgument("Matrix_Sparse::Val(int, int)",
1060  "No reference to element (" + to_str(i) + ", "
1061  + to_str(j)
1062  + ") can be returned: it is a zero entry.");
1063 
1064  l = Storage::GetSecond(i, j);
1065 
1066  for (k = a; (k < b-1) && (ind_[k] < l); k++);
1067 
1068  if (ind_[k] == l)
1069  return this->data_[k];
1070  else
1071  throw WrongArgument("Matrix_Sparse::Val(int, int)",
1072  "No reference to element (" + to_str(i) + ", "
1073  + to_str(j)
1074  + ") can be returned: it is a zero entry.");
1075  }
1076 
1077 
1079 
1086  template <class T, class Prop, class Storage, class Allocator>
1087  typename Matrix_Sparse<T, Prop, Storage, Allocator>::value_type&
1088  Matrix_Sparse<T, Prop, Storage, Allocator>::Get(int i, int j)
1089  {
1090 
1091 #ifdef SELDON_CHECK_BOUNDS
1092  if (i < 0 || i >= this->m_)
1093  throw WrongRow("Matrix_Sparse::Get(int, int)",
1094  string("Index should be in [0, ") + to_str(this->m_-1)
1095  + "], but is equal to " + to_str(i) + ".");
1096  if (j < 0 || j >= this->n_)
1097  throw WrongCol("Matrix_Sparse::Get(int, int)",
1098  string("Index should be in [0, ") + to_str(this->n_-1)
1099  + "], but is equal to " + to_str(j) + ".");
1100 #endif
1101 
1102  long k, l;
1103  long a, b;
1104 
1105  a = ptr_[Storage::GetFirst(i, j)];
1106  b = ptr_[Storage::GetFirst(i, j) + 1];
1107 
1108  if (a < b)
1109  {
1110  l = Storage::GetSecond(i, j);
1111 
1112  for (k = a; (k < b) && (ind_[k] < l); k++);
1113 
1114  if ( (k < b) && (ind_[k] == l))
1115  return this->data_[k];
1116  }
1117  else
1118  k = a;
1119 
1120  // adding a non-zero entry
1121  Resize(this->m_, this->n_, nz_+1);
1122 
1123  for (int m = Storage::GetFirst(i, j)+1;
1124  m <= Storage::GetFirst(this->m_, this->n_); m++)
1125  ptr_[m]++;
1126 
1127  for (long m = nz_-1; m >= k+1; m--)
1128  {
1129  ind_[m] = ind_[m-1];
1130  this->data_[m] = this->data_[m-1];
1131  }
1132 
1133  ind_[k] = Storage::GetSecond(i, j);
1134 
1135  // value of new non-zero entry is set to 0
1136  SetComplexZero(this->data_[k]);
1137 
1138  return this->data_[k];
1139  }
1140 
1141 
1142  /************************
1143  * CONVENIENT FUNCTIONS *
1144  ************************/
1145 
1146 
1148 
1149  template <class T, class Prop, class Storage, class Allocator>
1150  void Matrix_Sparse<T, Prop, Storage, Allocator>::Zero()
1151  {
1152  Allocator::memoryset(this->data_, char(0),
1153  this->nz_ * sizeof(value_type));
1154  }
1155 
1156 
1158 
1161  template <class T, class Prop, class Storage, class Allocator>
1162  void Matrix_Sparse<T, Prop, Storage, Allocator>::SetIdentity()
1163  {
1164  int m = this->m_;
1165  int n = this->n_;
1166  long nz = min(m, n);
1167 
1168  if (nz == 0)
1169  return;
1170 
1171  Clear();
1172 
1173  Vector<T, VectFull, Allocator> values(nz);
1174  Vector<long> ptr(Storage::GetFirst(m, n) + 1);
1175  Vector<int> ind(nz);
1176 
1177  T one; SetComplexOne(one);
1178  values.Fill(one);
1179  ind.Fill();
1180  int i;
1181  for (i = 0; i < nz + 1; i++)
1182  ptr(i) = i;
1183 
1184  for (i = nz + 1; i < ptr.GetM(); i++)
1185  ptr(i) = nz;
1186 
1187  SetData(m, n, values, ptr, ind);
1188  }
1189 
1190 
1192 
1195  template <class T, class Prop, class Storage, class Allocator>
1196  void Matrix_Sparse<T, Prop, Storage, Allocator>::Fill()
1197  {
1198  for (long i = 0; i < this->GetDataSize(); i++)
1199  SetComplexReal(i, this->data_[i]);
1200  }
1201 
1202 
1204 
1207  template <class T, class Prop, class Storage, class Allocator>
1208  template <class T0>
1209  void Matrix_Sparse<T, Prop, Storage, Allocator>::Fill(const T0& x)
1210  {
1211  T x_;
1212  SetComplexReal(x, x_);
1213  for (long i = 0; i < this->GetDataSize(); i++)
1214  this->data_[i] = x_;
1215  }
1216 
1217 
1219 
1222  template <class T, class Prop, class Storage, class Allocator>
1223  void Matrix_Sparse<T, Prop, Storage, Allocator>::FillRand()
1224  {
1225 #ifndef SELDON_WITHOUT_REINIT_RANDOM
1226  srand(time(NULL));
1227 #endif
1228  for (long i = 0; i < this->GetDataSize(); i++)
1229  SetComplexReal(rand(), this->data_[i]);
1230  }
1231 
1232 
1234 
1242  template <class T, class Prop, class Storage, class Allocator>
1243  void Matrix_Sparse<T, Prop, Storage, Allocator>
1244  ::FillRand(long Nelement)
1245  {
1246  if (this->m_ == 0 || this->n_ == 0)
1247  return;
1248 
1249  Vector<int> i(Nelement), j(Nelement);
1250  Vector<T> value(Nelement);
1251 
1252  set<pair<int, int> > skeleton;
1253  set<pair<int, int> >::iterator it;
1254 
1255 #ifndef SELDON_WITHOUT_REINIT_RANDOM
1256  srand(time(NULL));
1257 #endif
1258 
1259  // generation of triplet (i, j, value)
1260  while (static_cast<int>(skeleton.size()) != Nelement)
1261  skeleton.insert(make_pair(rand() % this->m_, rand() % this->n_));
1262 
1263  long l = 0;
1264  for (it = skeleton.begin(); it != skeleton.end(); it++)
1265  {
1266  i(l) = it->first;
1267  j(l) = it->second;
1268  SetComplexReal(rand(), value(l));
1269  l++;
1270  }
1271 
1272  // then conversion to current sparse matrix
1273  Matrix<T, Prop, Storage, Allocator>& leaf_class =
1274  static_cast<Matrix<T, Prop, Storage, Allocator>& >(*this);
1275 
1276  ConvertMatrix_from_Coordinates(i, j, value, leaf_class, 0);
1277  }
1278 
1279 
1281 
1291  template <class T, class Prop, class Storage, class Allocator>
1292  void Matrix_Sparse<T, Prop, Storage, Allocator>
1293  ::FillRand(long Nelement, const T& x)
1294  {
1295  if (this->m_ == 0 || this->n_ == 0)
1296  return;
1297 
1298  Vector<int> i(Nelement), j(Nelement);
1299  Vector<T> value(Nelement);
1300  value.Fill(x);
1301 
1302 #ifndef SELDON_WITHOUT_REINIT_RANDOM
1303  srand(time(NULL));
1304 #endif
1305 
1306  for (long l = 0; l < Nelement; l++)
1307  {
1308  i(l) = rand() % this->m_;
1309  j(l) = rand() % this->n_;
1310  }
1311 
1312  Matrix<T, Prop, Storage, Allocator>& leaf_class =
1313  static_cast<Matrix<T, Prop, Storage, Allocator>& >(*this);
1314 
1315  ConvertMatrix_from_Coordinates(i, j, value, leaf_class, 0);
1316  }
1317 
1318 
1320 
1325  template <class T, class Prop, class Storage, class Allocator>
1326  void Matrix_Sparse<T, Prop, Storage, Allocator>::Print() const
1327  {
1328  for (int i = 0; i < this->m_; i++)
1329  {
1330  for (int j = 0; j < this->n_; j++)
1331  cout << (*this)(i, j) << "\t";
1332  cout << endl;
1333  }
1334  }
1335 
1336 
1338 
1342  template <class T, class Prop, class Storage, class Allocator>
1343  void Matrix_Sparse<T, Prop, Storage, Allocator>
1344  ::Write(string FileName) const
1345  {
1346  ofstream FileStream;
1347  FileStream.open(FileName.c_str(), ofstream::binary);
1348 
1349 #ifdef SELDON_CHECK_IO
1350  // Checks if the file was opened.
1351  if (!FileStream.is_open())
1352  throw IOError("Matrix_Sparse::Write(string FileName)",
1353  string("Unable to open file \"") + FileName + "\".");
1354 #endif
1355 
1356  this->Write(FileStream);
1357 
1358  FileStream.close();
1359  }
1360 
1361 
1363 
1367  template <class T, class Prop, class Storage, class Allocator>
1368  void Matrix_Sparse<T, Prop, Storage, Allocator>
1369  ::Write(ostream& FileStream) const
1370  {
1371 #ifdef SELDON_CHECK_IO
1372  // Checks if the stream is ready.
1373  if (!FileStream.good())
1374  throw IOError("Matrix_Sparse::Write(ofstream& FileStream)",
1375  "Stream is not ready.");
1376 #endif
1377 
1378  FileStream.write(reinterpret_cast<char*>(const_cast<int*>(&this->m_)),
1379  sizeof(int));
1380  FileStream.write(reinterpret_cast<char*>(const_cast<int*>(&this->n_)),
1381  sizeof(int));
1382  FileStream.write(reinterpret_cast<char*>(const_cast<long*>(&this->nz_)),
1383  sizeof(long));
1384 
1385  FileStream.write(reinterpret_cast<char*>(this->ptr_),
1386  sizeof(long)*(Storage::GetFirst(this->m_, this->n_)+1));
1387  FileStream.write(reinterpret_cast<char*>(this->ind_),
1388  sizeof(int)*this->nz_);
1389  FileStream.write(reinterpret_cast<char*>(this->data_),
1390  sizeof(T)*this->nz_);
1391  }
1392 
1393 
1395 
1404  template <class T, class Prop, class Storage, class Allocator>
1405  void Matrix_Sparse<T, Prop, Storage, Allocator>::
1406  WriteText(string FileName, bool cplx) const
1407  {
1408  ofstream FileStream; FileStream.precision(14);
1409  FileStream.open(FileName.c_str());
1410 
1411  // changing precision
1412  FileStream.precision(cout.precision());
1413 
1414 #ifdef SELDON_CHECK_IO
1415  // Checks if the file was opened.
1416  if (!FileStream.is_open())
1417  throw IOError("Matrix_Sparse::Write(string FileName)",
1418  string("Unable to open file \"") + FileName + "\".");
1419 #endif
1420 
1421  this->WriteText(FileStream, cplx);
1422 
1423  FileStream.close();
1424  }
1425 
1426 
1428 
1437  template <class T, class Prop, class Storage, class Allocator>
1438  void Matrix_Sparse<T, Prop, Storage, Allocator>::
1439  WriteText(ostream& FileStream, bool cplx) const
1440  {
1441 
1442 #ifdef SELDON_CHECK_IO
1443  // Checks if the stream is ready.
1444  if (!FileStream.good())
1445  throw IOError("Matrix_Sparse::Write(ofstream& FileStream)",
1446  "Stream is not ready.");
1447 #endif
1448 
1449  // conversion in coordinate format (1-index convention)
1450  const Matrix<T, Prop, Storage, Allocator>& leaf_class =
1451  static_cast<const Matrix<T, Prop, Storage, Allocator>& >(*this);
1452 
1453  T zero; int index = 1;
1454  WriteCoordinateMatrix(leaf_class, FileStream, zero, index, cplx);
1455  }
1456 
1457 
1459 
1463  template <class T, class Prop, class Storage, class Allocator>
1464  void Matrix_Sparse<T, Prop, Storage, Allocator>::
1465  Read(string FileName)
1466  {
1467  ifstream FileStream;
1468  FileStream.open(FileName.c_str(), ifstream::binary);
1469 
1470 #ifdef SELDON_CHECK_IO
1471  // Checks if the file was opened.
1472  if (!FileStream.is_open())
1473  throw IOError("Matrix_Sparse::Read(string FileName)",
1474  string("Unable to open file \"") + FileName + "\".");
1475 #endif
1476 
1477  this->Read(FileStream);
1478 
1479  FileStream.close();
1480  }
1481 
1482 
1484 
1488  template <class T, class Prop, class Storage, class Allocator>
1489  void Matrix_Sparse<T, Prop, Storage, Allocator>::
1490  Read(istream& FileStream)
1491  {
1492 
1493 #ifdef SELDON_CHECK_IO
1494  // Checks if the stream is ready.
1495  if (!FileStream.good())
1496  throw IOError("Matrix_Sparse::Read(istream& FileStream)",
1497  "Stream is not ready.");
1498 #endif
1499 
1500  int m, n; long nz;
1501  FileStream.read(reinterpret_cast<char*>(&m), sizeof(int));
1502  FileStream.read(reinterpret_cast<char*>(&n), sizeof(int));
1503  FileStream.read(reinterpret_cast<char*>(&nz), sizeof(long));
1504 
1505  Reallocate(m, n, nz);
1506 
1507  FileStream.read(reinterpret_cast<char*>(ptr_),
1508  sizeof(long)*(Storage::GetFirst(m, n)+1));
1509  FileStream.read(reinterpret_cast<char*>(ind_), sizeof(int)*nz);
1510  FileStream.read(reinterpret_cast<char*>(this->data_), sizeof(T)*nz);
1511 
1512 #ifdef SELDON_CHECK_IO
1513  // Checks if data was read.
1514  if (!FileStream.good())
1515  throw IOError("Matrix_Sparse::Read(istream& FileStream)",
1516  string("Input operation failed.")
1517  + string(" The input file may have been removed")
1518  + " or may not contain enough data.");
1519 #endif
1520 
1521  }
1522 
1523 
1525 
1532  template <class T, class Prop, class Storage, class Allocator>
1533  void Matrix_Sparse<T, Prop, Storage, Allocator>::
1534  ReadText(string FileName, bool cplx)
1535  {
1536  ifstream FileStream;
1537  FileStream.open(FileName.c_str());
1538 
1539 #ifdef SELDON_CHECK_IO
1540  // Checks if the file was opened.
1541  if (!FileStream.is_open())
1542  throw IOError("Matrix_Sparse::ReadText(string FileName)",
1543  string("Unable to open file \"") + FileName + "\".");
1544 #endif
1545 
1546  this->ReadText(FileStream, cplx);
1547 
1548  FileStream.close();
1549  }
1550 
1551 
1553 
1557  template <class T, class Prop, class Storage, class Allocator>
1558  void Matrix_Sparse<T, Prop, Storage, Allocator>::
1559  ReadText(istream& FileStream, bool cplx)
1560  {
1561  Matrix<T, Prop, Storage, Allocator>& leaf_class =
1562  static_cast<Matrix<T, Prop, Storage, Allocator>& >(*this);
1563 
1564  T zero; int index = 1;
1565  ReadCoordinateMatrix(leaf_class, FileStream, zero, index, -1, cplx);
1566  }
1567 
1568 
1569 } // namespace Seldon.
1570 
1571 #define SELDON_FILE_MATRIX_SPARSE_CXX
1572 #endif
Seldon::Matrix_Base
Base class for all matrices.
Definition: Matrix_Base.hxx:143
Seldon::Matrix_Sparse::Matrix_Sparse
Matrix_Sparse()
Default constructor.
Definition: Matrix_SparseInline.cxx:53
Seldon::to_str
std::string to_str(const T &input)
Converts most types to string.
Definition: CommonInline.cxx:137
Seldon::WrongCol
Definition: Errors.hxx:138
Seldon::Vector
Definition: SeldonHeader.hxx:207
Seldon::Vector< T, VectFull, Allocator >::Nullify
void Nullify()
Clears the vector without releasing memory.
Definition: VectorInline.cxx:476
Seldon::Matrix
Definition: SeldonHeader.hxx:226
Seldon::WrongRow
Definition: Errors.hxx:126
Seldon::Vector< T, VectFull, Allocator >::Resize
void Resize(size_t i)
Changes the length of the vector, and keeps previous values.
Definition: Vector.cxx:34
Seldon::Vector< T, VectFull, Allocator >
Full vector class.
Definition: Vector.hxx:88
Seldon::Matrix_Sparse
Sparse-matrix class.
Definition: Matrix_Sparse.hxx:45
Seldon::WrongArgument
Definition: Errors.hxx:76
Seldon::WrongDim
Definition: Errors.hxx:102
Seldon
Seldon namespace.
Definition: Array.cxx:24
Seldon::IOError
Definition: Errors.hxx:150
Seldon::NoMemory
Definition: Errors.hxx:90
Seldon::Vector< T, VectFull, Allocator >::Fill
void Fill()
Fills the vector with 0, 1, 2, ...
Definition: VectorInline.cxx:736
Seldon::SetComplexReal
void SetComplexReal(int n, std::complex< T > &number)
Sets a complex number to (n, 0).
Definition: CommonInline.cxx:234
Seldon::Vector_Base::GetData
pointer GetData() const
Returns a pointer to data_ (stored data).
Definition: VectorInline.cxx:177