Warning: this documentation for the development version is under construction.
/home/vivien/public_html/.src_seldon/matrix_sparse/Matrix_ComplexSparse.cxx
00001 // Copyright (C) 2001-2009 Vivien Mallet 00002 // Copyright (C) 2003-2009 Marc Duruflé 00003 // 00004 // This file is part of the linear-algebra library Seldon, 00005 // http://seldon.sourceforge.net/. 00006 // 00007 // Seldon is free software; you can redistribute it and/or modify it under the 00008 // terms of the GNU Lesser General Public License as published by the Free 00009 // Software Foundation; either version 2.1 of the License, or (at your option) 00010 // any later version. 00011 // 00012 // Seldon is distributed in the hope that it will be useful, but WITHOUT ANY 00013 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS 00014 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for 00015 // more details. 00016 // 00017 // You should have received a copy of the GNU Lesser General Public License 00018 // along with Seldon. If not, see http://www.gnu.org/licenses/. 00019 00020 00021 #ifndef SELDON_FILE_MATRIX_COMPLEXSPARSE_CXX 00022 00023 #include "Matrix_ComplexSparse.hxx" 00024 00025 namespace Seldon 00026 { 00027 00028 00029 /**************** 00030 * CONSTRUCTORS * 00031 ****************/ 00032 00033 00035 00038 template <class T, class Prop, class Storage, class Allocator> 00039 inline Matrix_ComplexSparse<T, Prop, Storage, Allocator> 00040 ::Matrix_ComplexSparse(): Matrix_Base<T, Allocator>() 00041 { 00042 real_nz_ = 0; 00043 imag_nz_ = 0; 00044 real_ptr_ = NULL; 00045 imag_ptr_ = NULL; 00046 real_ind_ = NULL; 00047 imag_ind_ = NULL; 00048 real_data_ = NULL; 00049 imag_data_ = NULL; 00050 } 00051 00052 00054 00059 template <class T, class Prop, class Storage, class Allocator> 00060 inline Matrix_ComplexSparse<T, Prop, Storage, Allocator> 00061 ::Matrix_ComplexSparse(int i, int j): Matrix_Base<T, Allocator>() 00062 { 00063 real_nz_ = 0; 00064 imag_nz_ = 0; 00065 real_ptr_ = NULL; 00066 imag_ptr_ = NULL; 00067 real_ind_ = NULL; 00068 imag_ind_ = NULL; 00069 real_data_ = NULL; 00070 imag_data_ = NULL; 00071 00072 Reallocate(i, j); 00073 } 00074 00075 00077 00087 template <class T, class Prop, class Storage, class Allocator> 00088 inline Matrix_ComplexSparse<T, Prop, Storage, Allocator>:: 00089 Matrix_ComplexSparse(int i, int j, int real_nz, int imag_nz): 00090 Matrix_Base<T, Allocator>() 00091 { 00092 real_nz_ = 0; 00093 imag_nz_ = 0; 00094 real_ptr_ = NULL; 00095 imag_ptr_ = NULL; 00096 real_ind_ = NULL; 00097 imag_ind_ = NULL; 00098 real_data_ = NULL; 00099 imag_data_ = NULL; 00100 00101 Reallocate(i, j, real_nz, imag_nz); 00102 } 00103 00104 00106 00124 template <class T, class Prop, class Storage, class Allocator> 00125 template <class Storage0, class Allocator0, 00126 class Storage1, class Allocator1, 00127 class Storage2, class Allocator2> 00128 inline Matrix_ComplexSparse<T, Prop, Storage, Allocator>:: 00129 Matrix_ComplexSparse(int i, int j, 00130 Vector<T, Storage0, Allocator0>& real_values, 00131 Vector<int, Storage1, Allocator1>& real_ptr, 00132 Vector<int, Storage2, Allocator2>& real_ind, 00133 Vector<T, Storage0, Allocator0>& imag_values, 00134 Vector<int, Storage1, Allocator1>& imag_ptr, 00135 Vector<int, Storage2, Allocator2>& imag_ind): 00136 Matrix_Base<T, Allocator>(i, j) 00137 { 00138 00139 real_nz_ = real_values.GetLength(); 00140 imag_nz_ = imag_values.GetLength(); 00141 00142 #ifdef SELDON_CHECK_DIMENSIONS 00143 // Checks whether vector sizes are acceptable. 00144 00145 if (real_ind.GetLength() != real_nz_) 00146 { 00147 this->m_ = 0; 00148 this->n_ = 0; 00149 real_nz_ = 0; 00150 imag_nz_ = 0; 00151 real_ptr_ = NULL; 00152 imag_ptr_ = NULL; 00153 real_ind_ = NULL; 00154 imag_ind_ = NULL; 00155 this->real_data_ = NULL; 00156 this->imag_data_ = NULL; 00157 throw WrongDim(string("Matrix_ComplexSparse::") 00158 + string("Matrix_ComplexSparse(int, int, ") 00159 + string("const Vector&, const Vector&, const Vector&") 00160 + ", const Vector&, const Vector&, const Vector&)", 00161 string("There are ") + to_str(real_nz_) 00162 + " values (real part) but " 00163 + to_str(real_ind.GetLength()) 00164 + " row or column indices."); 00165 } 00166 00167 if (imag_ind.GetLength() != imag_nz_) 00168 { 00169 this->m_ = 0; 00170 this->n_ = 0; 00171 real_nz_ = 0; 00172 imag_nz_ = 0; 00173 real_ptr_ = NULL; 00174 imag_ptr_ = NULL; 00175 real_ind_ = NULL; 00176 imag_ind_ = NULL; 00177 this->real_data_ = NULL; 00178 this->imag_data_ = NULL; 00179 throw WrongDim(string("Matrix_ComplexSparse::") 00180 + string("Matrix_ComplexSparse(int, int, ") 00181 + string("const Vector&, const Vector&, const Vector&") 00182 + ", const Vector&, const Vector&, const Vector&)", 00183 string("There are ") + to_str(imag_nz_) 00184 + " values (imaginary part) but " 00185 + to_str(imag_ind.GetLength()) 00186 + " row or column indices."); 00187 } 00188 00189 if (real_ptr.GetLength()-1 != Storage::GetFirst(i, j)) 00190 { 00191 this->m_ = 0; 00192 this->n_ = 0; 00193 real_nz_ = 0; 00194 imag_nz_ = 0; 00195 real_ptr_ = NULL; 00196 imag_ptr_ = NULL; 00197 real_ind_ = NULL; 00198 imag_ind_ = NULL; 00199 this->real_data_ = NULL; 00200 this->imag_data_ = NULL; 00201 throw WrongDim(string("Matrix_ComplexSparse::") 00202 + string("Matrix_ComplexSparse(int, int, ") 00203 + string("const Vector&, const Vector&, const Vector&") 00204 + ", const Vector&, const Vector&, const Vector&)", 00205 string("The vector of start indices (real part)") 00206 + " contains " + to_str(real_ptr.GetLength()-1) 00207 + string(" row or column start indices (plus the") 00208 + " number of non-zero entries) but there are " 00209 + to_str(Storage::GetFirst(i, j)) 00210 + " rows or columns (" 00211 + to_str(i) + " by " + to_str(j) + " matrix)."); 00212 } 00213 00214 if (imag_ptr.GetLength()-1 != Storage::GetFirst(i, j)) 00215 { 00216 this->m_ = 0; 00217 this->n_ = 0; 00218 real_nz_ = 0; 00219 imag_nz_ = 0; 00220 real_ptr_ = NULL; 00221 imag_ptr_ = NULL; 00222 real_ind_ = NULL; 00223 imag_ind_ = NULL; 00224 this->real_data_ = NULL; 00225 this->imag_data_ = NULL; 00226 throw WrongDim(string("Matrix_ComplexSparse::") 00227 + string("Matrix_ComplexSparse(int, int, ") 00228 + string("const Vector&, const Vector&, const Vector&") 00229 + ", const Vector&, const Vector&, const Vector&)", 00230 string("The vector of start indices (imaginary part)") 00231 + " contains " + to_str(imag_ptr.GetLength()-1) 00232 + string(" row or column start indices (plus the") 00233 + " number of non-zero entries) but there are " 00234 + to_str(Storage::GetFirst(i, j)) 00235 + " rows or columns (" 00236 + to_str(i) + " by " + to_str(j) + " matrix)."); 00237 } 00238 00239 if ((real_nz_ > 0 00240 && (j == 0 00241 || static_cast<long int>(real_nz_-1) / static_cast<long int>(j) 00242 >= static_cast<long int>(i))) 00243 || 00244 (imag_nz_ > 0 00245 && (j == 0 00246 || static_cast<long int>(imag_nz_-1) / static_cast<long int>(j) 00247 >= static_cast<long int>(i)))) 00248 { 00249 this->m_ = 0; 00250 this->n_ = 0; 00251 real_nz_ = 0; 00252 imag_nz_ = 0; 00253 real_ptr_ = NULL; 00254 imag_ptr_ = NULL; 00255 real_ind_ = NULL; 00256 imag_ind_ = NULL; 00257 this->real_data_ = NULL; 00258 this->imag_data_ = NULL; 00259 throw WrongDim(string("Matrix_ComplexSparse::") 00260 + string("Matrix_ComplexSparse(int, int, ") 00261 + string("const Vector&, const Vector&, const Vector&") 00262 + ", const Vector&, const Vector&, const Vector&)", 00263 string("There are more values (") 00264 + to_str(real_values.GetLength()) 00265 + " values for the real part and " 00266 + to_str(real_values.GetLength()) + string(" values") 00267 + string(" for the imaginary part) than elements") 00268 + " in the matrix (" 00269 + to_str(i) + " by " + to_str(j) + ")."); 00270 } 00271 #endif 00272 00273 this->real_ptr_ = real_ptr.GetData(); 00274 this->imag_ptr_ = imag_ptr.GetData(); 00275 this->real_ind_ = real_ind.GetData(); 00276 this->imag_ind_ = imag_ind.GetData(); 00277 this->real_data_ = real_values.GetData(); 00278 this->imag_data_ = imag_values.GetData(); 00279 00280 real_ptr.Nullify(); 00281 imag_ptr.Nullify(); 00282 real_ind.Nullify(); 00283 imag_ind.Nullify(); 00284 real_values.Nullify(); 00285 imag_values.Nullify(); 00286 } 00287 00288 00290 template <class T, class Prop, class Storage, class Allocator> 00291 Matrix_ComplexSparse<T, Prop, Storage, Allocator> 00292 ::Matrix_ComplexSparse(const Matrix_ComplexSparse<T, Prop, 00293 Storage, Allocator>& A) 00294 { 00295 this->m_ = 0; 00296 this->n_ = 0; 00297 real_nz_ = 0; 00298 imag_nz_ = 0; 00299 real_ptr_ = NULL; 00300 imag_ptr_ = NULL; 00301 real_ind_ = NULL; 00302 imag_ind_ = NULL; 00303 real_data_ = NULL; 00304 imag_data_ = NULL; 00305 00306 this->Copy(A); 00307 } 00308 00309 00310 /************** 00311 * DESTRUCTOR * 00312 **************/ 00313 00314 00316 template <class T, class Prop, class Storage, class Allocator> 00317 inline Matrix_ComplexSparse<T, Prop, Storage, Allocator> 00318 ::~Matrix_ComplexSparse() 00319 { 00320 this->m_ = 0; 00321 this->n_ = 0; 00322 00323 #ifdef SELDON_CHECK_MEMORY 00324 try 00325 { 00326 #endif 00327 00328 if (real_ptr_ != NULL) 00329 { 00330 free(real_ptr_); 00331 real_ptr_ = NULL; 00332 } 00333 00334 #ifdef SELDON_CHECK_MEMORY 00335 } 00336 catch (...) 00337 { 00338 real_ptr_ = NULL; 00339 } 00340 #endif 00341 00342 #ifdef SELDON_CHECK_MEMORY 00343 try 00344 { 00345 #endif 00346 00347 if (imag_ptr_ != NULL) 00348 { 00349 free(imag_ptr_); 00350 imag_ptr_ = NULL; 00351 } 00352 00353 #ifdef SELDON_CHECK_MEMORY 00354 } 00355 catch (...) 00356 { 00357 imag_ptr_ = NULL; 00358 } 00359 #endif 00360 00361 #ifdef SELDON_CHECK_MEMORY 00362 try 00363 { 00364 #endif 00365 00366 if (real_ind_ != NULL) 00367 { 00368 free(real_ind_); 00369 real_ind_ = NULL; 00370 } 00371 00372 #ifdef SELDON_CHECK_MEMORY 00373 } 00374 catch (...) 00375 { 00376 real_ind_ = NULL; 00377 } 00378 #endif 00379 00380 #ifdef SELDON_CHECK_MEMORY 00381 try 00382 { 00383 #endif 00384 00385 if (imag_ind_ != NULL) 00386 { 00387 free(imag_ind_); 00388 imag_ind_ = NULL; 00389 } 00390 00391 #ifdef SELDON_CHECK_MEMORY 00392 } 00393 catch (...) 00394 { 00395 imag_ind_ = NULL; 00396 } 00397 #endif 00398 00399 #ifdef SELDON_CHECK_MEMORY 00400 try 00401 { 00402 #endif 00403 00404 if (this->real_data_ != NULL) 00405 { 00406 this->allocator_.deallocate(this->real_data_, real_nz_); 00407 this->real_data_ = NULL; 00408 } 00409 00410 #ifdef SELDON_CHECK_MEMORY 00411 } 00412 catch (...) 00413 { 00414 this->real_nz_ = 0; 00415 this->real_data_ = NULL; 00416 } 00417 #endif 00418 00419 #ifdef SELDON_CHECK_MEMORY 00420 try 00421 { 00422 #endif 00423 00424 if (this->imag_data_ != NULL) 00425 { 00426 this->allocator_.deallocate(this->imag_data_, imag_nz_); 00427 this->imag_data_ = NULL; 00428 } 00429 00430 #ifdef SELDON_CHECK_MEMORY 00431 } 00432 catch (...) 00433 { 00434 this->imag_nz_ = 0; 00435 this->imag_data_ = NULL; 00436 } 00437 #endif 00438 00439 this->real_nz_ = 0; 00440 this->imag_nz_ = 0; 00441 } 00442 00443 00445 00448 template <class T, class Prop, class Storage, class Allocator> 00449 inline void Matrix_ComplexSparse<T, Prop, Storage, Allocator>::Clear() 00450 { 00451 this->~Matrix_ComplexSparse(); 00452 } 00453 00454 00455 /********************* 00456 * MEMORY MANAGEMENT * 00457 *********************/ 00458 00459 00461 00478 template <class T, class Prop, class Storage, class Allocator> 00479 template <class Storage0, class Allocator0, 00480 class Storage1, class Allocator1, 00481 class Storage2, class Allocator2> 00482 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>:: 00483 SetData(int i, int j, 00484 Vector<T, Storage0, Allocator0>& real_values, 00485 Vector<int, Storage1, Allocator1>& real_ptr, 00486 Vector<int, Storage2, Allocator2>& real_ind, 00487 Vector<T, Storage0, Allocator0>& imag_values, 00488 Vector<int, Storage1, Allocator1>& imag_ptr, 00489 Vector<int, Storage2, Allocator2>& imag_ind) 00490 { 00491 this->Clear(); 00492 this->m_ = i; 00493 this->n_ = j; 00494 real_nz_ = real_values.GetLength(); 00495 imag_nz_ = imag_values.GetLength(); 00496 00497 #ifdef SELDON_CHECK_DIMENSIONS 00498 // Checks whether vector sizes are acceptable. 00499 00500 if (real_ind.GetLength() != real_nz_) 00501 { 00502 this->m_ = 0; 00503 this->n_ = 0; 00504 real_nz_ = 0; 00505 imag_nz_ = 0; 00506 real_ptr_ = NULL; 00507 imag_ptr_ = NULL; 00508 real_ind_ = NULL; 00509 imag_ind_ = NULL; 00510 this->real_data_ = NULL; 00511 this->imag_data_ = NULL; 00512 throw WrongDim(string("Matrix_ComplexSparse::SetData(int, int, ") 00513 + string("const Vector&, const Vector&, const Vector&") 00514 + ", const Vector&, const Vector&, const Vector&)", 00515 string("There are ") + to_str(real_nz_) 00516 + " values (real part) but " 00517 + to_str(real_ind.GetLength()) 00518 + " row or column indices."); 00519 } 00520 00521 if (imag_ind.GetLength() != imag_nz_) 00522 { 00523 this->m_ = 0; 00524 this->n_ = 0; 00525 real_nz_ = 0; 00526 imag_nz_ = 0; 00527 real_ptr_ = NULL; 00528 imag_ptr_ = NULL; 00529 real_ind_ = NULL; 00530 imag_ind_ = NULL; 00531 this->real_data_ = NULL; 00532 this->imag_data_ = NULL; 00533 throw WrongDim(string("Matrix_ComplexSparse::SetData(int, int, ") 00534 + string("const Vector&, const Vector&, const Vector&") 00535 + ", const Vector&, const Vector&, const Vector&)", 00536 string("There are ") + to_str(imag_nz_) 00537 + " values (imaginary part) but " 00538 + to_str(imag_ind.GetLength()) 00539 + " row or column indices."); 00540 } 00541 00542 if (real_ptr.GetLength()-1 != Storage::GetFirst(i, j)) 00543 { 00544 this->m_ = 0; 00545 this->n_ = 0; 00546 real_nz_ = 0; 00547 imag_nz_ = 0; 00548 real_ptr_ = NULL; 00549 imag_ptr_ = NULL; 00550 real_ind_ = NULL; 00551 imag_ind_ = NULL; 00552 this->real_data_ = NULL; 00553 this->imag_data_ = NULL; 00554 throw WrongDim(string("Matrix_ComplexSparse::SetData(int, int, ") 00555 + string("const Vector&, const Vector&, const Vector&") 00556 + ", const Vector&, const Vector&, const Vector&)", 00557 string("The vector of start indices (real part)") 00558 + " contains " + to_str(real_ptr.GetLength()-1) 00559 + string(" row or column start indices (plus the") 00560 + " number of non-zero entries) but there are " 00561 + to_str(Storage::GetFirst(i, j)) 00562 + " rows or columns (" 00563 + to_str(i) + " by " + to_str(j) + " matrix)."); 00564 } 00565 00566 if (imag_ptr.GetLength()-1 != Storage::GetFirst(i, j)) 00567 { 00568 this->m_ = 0; 00569 this->n_ = 0; 00570 real_nz_ = 0; 00571 imag_nz_ = 0; 00572 real_ptr_ = NULL; 00573 imag_ptr_ = NULL; 00574 real_ind_ = NULL; 00575 imag_ind_ = NULL; 00576 this->real_data_ = NULL; 00577 this->imag_data_ = NULL; 00578 throw WrongDim(string("Matrix_ComplexSparse::SetData(int, int, ") 00579 + string("const Vector&, const Vector&, const Vector&") 00580 + ", const Vector&, const Vector&, const Vector&)", 00581 string("The vector of start indices (imaginary part)") 00582 + " contains " + to_str(imag_ptr.GetLength()-1) 00583 + string(" row or column start indices (plus the") 00584 + " number of non-zero entries) but there are " 00585 + to_str(Storage::GetFirst(i, j)) 00586 + " rows or columns (" 00587 + to_str(i) + " by " + to_str(j) + " matrix)."); 00588 } 00589 00590 if ((real_nz_ > 0 00591 && (j == 0 00592 || static_cast<long int>(real_nz_-1) / static_cast<long int>(j) 00593 >= static_cast<long int>(i))) 00594 || 00595 (imag_nz_ > 0 00596 && (j == 0 00597 || static_cast<long int>(imag_nz_-1) / static_cast<long int>(j) 00598 >= static_cast<long int>(i)))) 00599 { 00600 this->m_ = 0; 00601 this->n_ = 0; 00602 real_nz_ = 0; 00603 imag_nz_ = 0; 00604 real_ptr_ = NULL; 00605 imag_ptr_ = NULL; 00606 real_ind_ = NULL; 00607 imag_ind_ = NULL; 00608 this->real_data_ = NULL; 00609 this->imag_data_ = NULL; 00610 throw WrongDim(string("Matrix_ComplexSparse::SetData(int, int, ") 00611 + string("const Vector&, const Vector&, const Vector&") 00612 + ", const Vector&, const Vector&, const Vector&)", 00613 string("There are more values (") 00614 + to_str(real_values.GetLength()) 00615 + " values for the real part and " 00616 + to_str(real_values.GetLength()) + string(" values") 00617 + string(" for the imaginary part) than elements") 00618 + " in the matrix (" 00619 + to_str(i) + " by " + to_str(j) + ")."); 00620 } 00621 #endif 00622 00623 this->real_ptr_ = real_ptr.GetData(); 00624 this->imag_ptr_ = imag_ptr.GetData(); 00625 this->real_ind_ = real_ind.GetData(); 00626 this->imag_ind_ = imag_ind.GetData(); 00627 this->real_data_ = real_values.GetData(); 00628 this->imag_data_ = imag_values.GetData(); 00629 00630 real_ptr.Nullify(); 00631 imag_ptr.Nullify(); 00632 real_ind.Nullify(); 00633 imag_ind.Nullify(); 00634 real_values.Nullify(); 00635 imag_values.Nullify(); 00636 } 00637 00638 00640 00661 template <class T, class Prop, class Storage, class Allocator> 00662 inline void Matrix_ComplexSparse<T, Prop, Storage, Allocator>:: 00663 SetData(int i, int j, int real_nz, 00664 typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 00665 ::pointer real_values, 00666 int* real_ptr, int* real_ind, int imag_nz, 00667 typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 00668 ::pointer imag_values, 00669 int* imag_ptr, int* imag_ind) 00670 { 00671 this->Clear(); 00672 00673 this->m_ = i; 00674 this->n_ = j; 00675 00676 this->real_nz_ = real_nz; 00677 this->imag_nz_ = imag_nz; 00678 00679 real_data_ = real_values; 00680 imag_data_ = imag_values; 00681 real_ind_ = real_ind; 00682 imag_ind_ = imag_ind; 00683 real_ptr_ = real_ptr; 00684 imag_ptr_ = imag_ptr; 00685 } 00686 00687 00689 00693 template <class T, class Prop, class Storage, class Allocator> 00694 inline void Matrix_ComplexSparse<T, Prop, Storage, Allocator>::Nullify() 00695 { 00696 this->data_ = NULL; 00697 this->m_ = 0; 00698 this->n_ = 0; 00699 real_nz_ = 0; 00700 real_ptr_ = NULL; 00701 real_ind_ = NULL; 00702 imag_nz_ = 0; 00703 imag_ptr_ = NULL; 00704 imag_ind_ = NULL; 00705 real_data_ = NULL; 00706 imag_data_ = NULL; 00707 } 00708 00709 00711 template <class T, class Prop, class Storage, class Allocator> 00712 inline void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 00713 ::Reallocate(int i, int j) 00714 { 00715 // previous entries are removed 00716 Clear(); 00717 00718 this->m_ = i; 00719 this->n_ = j; 00720 00721 // we try to allocate real_ptr_ and imag_ptr_ 00722 #ifdef SELDON_CHECK_MEMORY 00723 try 00724 { 00725 #endif 00726 00727 real_ptr_ = reinterpret_cast<int*>( calloc(Storage::GetFirst(i, j)+1, 00728 sizeof(int)) ); 00729 00730 imag_ptr_ = reinterpret_cast<int*>( calloc(Storage::GetFirst(i, j)+1, 00731 sizeof(int)) ); 00732 00733 #ifdef SELDON_CHECK_MEMORY 00734 } 00735 catch (...) 00736 { 00737 this->m_ = 0; 00738 this->n_ = 0; 00739 real_nz_ = 0; 00740 real_ptr_ = NULL; 00741 real_ind_ = NULL; 00742 imag_nz_ = 0; 00743 imag_ptr_ = NULL; 00744 imag_ind_ = NULL; 00745 real_data_ = NULL; 00746 imag_data_ = NULL; 00747 } 00748 if ((real_ptr_ == NULL) || (imag_ptr_ == NULL)) 00749 { 00750 this->m_ = 0; 00751 this->n_ = 0; 00752 real_nz_ = 0; 00753 real_ptr_ = NULL; 00754 real_ind_ = NULL; 00755 imag_nz_ = 0; 00756 imag_ptr_ = NULL; 00757 imag_ind_ = NULL; 00758 real_data_ = NULL; 00759 imag_data_ = NULL; 00760 } 00761 if (((real_ptr_ == NULL) || (imag_ptr_ == NULL)) && i != 0 && j != 0) 00762 throw NoMemory("Matrix_ComplexSparse::Reallocate(int, int)", 00763 string("Unable to allocate ") 00764 + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1) ) 00765 + " bytes to store " + to_str(Storage::GetFirst(i, j)+1) 00766 + " row or column start indices, for a " 00767 + to_str(i) + " by " + to_str(j) + " matrix."); 00768 #endif 00769 00770 // then filing real_ptr_ with 0 00771 for (int k = 0; k <= Storage::GetFirst(i, j); k++) 00772 { 00773 real_ptr_[k] = 0; 00774 imag_ptr_[k] = 0; 00775 } 00776 } 00777 00778 00780 template <class T, class Prop, class Storage, class Allocator> 00781 inline void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 00782 ::Reallocate(int i, int j, int real_nz, int imag_nz) 00783 { 00784 // previous entries are removed 00785 Clear(); 00786 00787 this->m_ = i; 00788 this->n_ = j; 00789 this->real_nz_ = real_nz; 00790 this->imag_nz_ = imag_nz; 00791 00792 #ifdef SELDON_CHECK_DIMENSIONS 00793 if (real_nz_ < 0 || imag_nz_ < 0) 00794 { 00795 this->m_ = 0; 00796 this->n_ = 0; 00797 real_nz_ = 0; 00798 imag_nz_ = 0; 00799 real_ptr_ = NULL; 00800 imag_ptr_ = NULL; 00801 real_ind_ = NULL; 00802 imag_ind_ = NULL; 00803 this->real_data_ = NULL; 00804 this->imag_data_ = NULL; 00805 throw WrongDim(string("Matrix_ComplexSparse::") 00806 + "Matrix_ComplexSparse(int, int, int, int)", 00807 "Invalid number of non-zero elements: " 00808 + to_str(real_nz) + " in the real part and " 00809 + to_str(imag_nz) + " in the imaginary part."); 00810 } 00811 if ((real_nz_ > 0 00812 && (j == 0 00813 || static_cast<long int>(real_nz_-1) / static_cast<long int>(j) 00814 >= static_cast<long int>(i))) 00815 || 00816 (imag_nz_ > 0 00817 && (j == 0 00818 || static_cast<long int>(imag_nz_-1) / static_cast<long int>(j) 00819 >= static_cast<long int>(i)))) 00820 { 00821 this->m_ = 0; 00822 this->n_ = 0; 00823 real_nz_ = 0; 00824 imag_nz_ = 0; 00825 real_ptr_ = NULL; 00826 imag_ptr_ = NULL; 00827 real_ind_ = NULL; 00828 imag_ind_ = NULL; 00829 this->real_data_ = NULL; 00830 this->imag_data_ = NULL; 00831 throw WrongDim(string("Matrix_ComplexSparse::") 00832 + "Matrix_ComplexSparse(int, int, int, int)", 00833 string("There are more values (") + to_str(real_nz) 00834 + " values for the real part and " + to_str(imag_nz) 00835 + string(" values for the imaginary part) than") 00836 + " elements in the matrix (" + to_str(i) + " by " 00837 + to_str(j) + ")."); 00838 } 00839 #endif 00840 00841 #ifdef SELDON_CHECK_MEMORY 00842 try 00843 { 00844 #endif 00845 00846 real_ptr_ = reinterpret_cast<int*>( calloc(Storage::GetFirst(i, j)+1, 00847 sizeof(int)) ); 00848 00849 #ifdef SELDON_CHECK_MEMORY 00850 } 00851 catch (...) 00852 { 00853 this->m_ = 0; 00854 this->n_ = 0; 00855 real_nz_ = 0; 00856 imag_nz_ = 0; 00857 real_ptr_ = NULL; 00858 imag_ptr_ = NULL; 00859 real_ind_ = NULL; 00860 imag_ind_ = NULL; 00861 this->real_data_ = NULL; 00862 this->imag_data_ = NULL; 00863 } 00864 if (real_ptr_ == NULL) 00865 { 00866 this->m_ = 0; 00867 this->n_ = 0; 00868 real_nz_ = 0; 00869 imag_nz_ = 0; 00870 imag_ptr_ = 0; 00871 real_ind_ = NULL; 00872 imag_ind_ = NULL; 00873 this->real_data_ = NULL; 00874 this->imag_data_ = NULL; 00875 } 00876 if (real_ptr_ == NULL && i != 0 && j != 0) 00877 throw NoMemory(string("Matrix_ComplexSparse::") 00878 + "Matrix_ComplexSparse(int, int, int, int)", 00879 string("Unable to allocate ") 00880 + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1)) 00881 + " bytes to store " + to_str(Storage::GetFirst(i, j)+1) 00882 + string(" row or column start indices (for the real") 00883 + " part), for a " 00884 + to_str(i) + " by " + to_str(j) + " matrix."); 00885 #endif 00886 00887 #ifdef SELDON_CHECK_MEMORY 00888 try 00889 { 00890 #endif 00891 00892 imag_ptr_ = reinterpret_cast<int*>( calloc(Storage::GetFirst(i, j)+1, 00893 sizeof(int)) ); 00894 00895 #ifdef SELDON_CHECK_MEMORY 00896 } 00897 catch (...) 00898 { 00899 this->m_ = 0; 00900 this->n_ = 0; 00901 real_nz_ = 0; 00902 imag_nz_ = 0; 00903 free(real_ptr_); 00904 real_ptr_ = NULL; 00905 imag_ptr_ = NULL; 00906 real_ind_ = NULL; 00907 imag_ind_ = NULL; 00908 this->real_data_ = NULL; 00909 this->imag_data_ = NULL; 00910 } 00911 if (imag_ptr_ == NULL) 00912 { 00913 this->m_ = 0; 00914 this->n_ = 0; 00915 real_nz_ = 0; 00916 imag_nz_ = 0; 00917 free(real_ptr_); 00918 real_ptr_ = 0; 00919 real_ind_ = NULL; 00920 imag_ind_ = NULL; 00921 this->real_data_ = NULL; 00922 this->imag_data_ = NULL; 00923 } 00924 if (imag_ptr_ == NULL && i != 0 && j != 0) 00925 throw NoMemory(string("Matrix_ComplexSparse::") 00926 + "Matrix_ComplexSparse(int, int, int, int)", 00927 string("Unable to allocate ") 00928 + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1)) 00929 + " bytes to store " + to_str(Storage::GetFirst(i, j)+1) 00930 + string(" row or column start indices (for the") 00931 + string(" imaginary part), for a ") 00932 + to_str(i) + " by " + to_str(j) + " matrix."); 00933 #endif 00934 00935 #ifdef SELDON_CHECK_MEMORY 00936 try 00937 { 00938 #endif 00939 00940 real_ind_ = reinterpret_cast<int*>( calloc(real_nz_, sizeof(int)) ); 00941 00942 #ifdef SELDON_CHECK_MEMORY 00943 } 00944 catch (...) 00945 { 00946 this->m_ = 0; 00947 this->n_ = 0; 00948 real_nz_ = 0; 00949 imag_nz_ = 0; 00950 free(real_ptr_); 00951 free(imag_ptr_); 00952 real_ptr_ = NULL; 00953 imag_ptr_ = NULL; 00954 real_ind_ = NULL; 00955 imag_ind_ = NULL; 00956 this->real_data_ = NULL; 00957 this->imag_data_ = NULL; 00958 } 00959 if (real_ind_ == NULL) 00960 { 00961 this->m_ = 0; 00962 this->n_ = 0; 00963 real_nz_ = 0; 00964 imag_nz_ = 0; 00965 free(real_ptr_); 00966 free(imag_ptr_); 00967 real_ptr_ = NULL; 00968 imag_ptr_ = NULL; 00969 this->real_data_ = NULL; 00970 this->imag_data_ = NULL; 00971 } 00972 if (real_ind_ == NULL && i != 0 && j != 0) 00973 throw NoMemory(string("Matrix_ComplexSparse::") 00974 + "Matrix_ComplexSparse(int, int, int, int)", 00975 string("Unable to allocate ") 00976 + to_str(sizeof(int) * real_nz) 00977 + " bytes to store " + to_str(real_nz) 00978 + " row or column indices (real part), for a " 00979 + to_str(i) + " by " + to_str(j) + " matrix."); 00980 #endif 00981 00982 #ifdef SELDON_CHECK_MEMORY 00983 try 00984 { 00985 #endif 00986 00987 imag_ind_ = reinterpret_cast<int*>( calloc(imag_nz_, sizeof(int)) ); 00988 00989 #ifdef SELDON_CHECK_MEMORY 00990 } 00991 catch (...) 00992 { 00993 this->m_ = 0; 00994 this->n_ = 0; 00995 real_nz_ = 0; 00996 imag_nz_ = 0; 00997 free(real_ptr_); 00998 free(imag_ptr_); 00999 real_ptr_ = NULL; 01000 imag_ptr_ = NULL; 01001 free(imag_ind_); 01002 real_ind_ = NULL; 01003 imag_ind_ = NULL; 01004 this->real_data_ = NULL; 01005 this->imag_data_ = NULL; 01006 } 01007 if (real_ind_ == NULL) 01008 { 01009 this->m_ = 0; 01010 this->n_ = 0; 01011 real_nz_ = 0; 01012 imag_nz_ = 0; 01013 free(real_ptr_); 01014 free(imag_ptr_); 01015 real_ptr_ = NULL; 01016 imag_ptr_ = NULL; 01017 free(imag_ind_); 01018 imag_ind_ = NULL; 01019 this->real_data_ = NULL; 01020 this->imag_data_ = NULL; 01021 } 01022 if (imag_ind_ == NULL && i != 0 && j != 0) 01023 throw NoMemory(string("Matrix_ComplexSparse::") 01024 + "Matrix_ComplexSparse(int, int, int, int)", 01025 string("Unable to allocate ") 01026 + to_str(sizeof(int) * imag_nz) 01027 + " bytes to store " + to_str(imag_nz) 01028 + " row or column indices (imaginary part), for a " 01029 + to_str(i) + " by " + to_str(j) + " matrix."); 01030 #endif 01031 01032 #ifdef SELDON_CHECK_MEMORY 01033 try 01034 { 01035 #endif 01036 01037 this->real_data_ = this->allocator_.allocate(real_nz_, this); 01038 01039 #ifdef SELDON_CHECK_MEMORY 01040 } 01041 catch (...) 01042 { 01043 this->m_ = 0; 01044 this->n_ = 0; 01045 free(real_ptr_); 01046 free(imag_ptr_); 01047 real_ptr_ = NULL; 01048 imag_ptr_ = NULL; 01049 free(real_ind_); 01050 free(imag_ind_); 01051 real_ind_ = NULL; 01052 imag_ind_ = NULL; 01053 this->real_data_ = NULL; 01054 this->imag_data_ = NULL; 01055 } 01056 if (real_data_ == NULL) 01057 { 01058 this->m_ = 0; 01059 this->n_ = 0; 01060 free(real_ptr_); 01061 free(imag_ptr_); 01062 real_ptr_ = NULL; 01063 imag_ptr_ = NULL; 01064 free(real_ind_); 01065 free(imag_ind_); 01066 real_ind_ = NULL; 01067 imag_ind_ = NULL; 01068 imag_data_ = NULL; 01069 } 01070 if (real_data_ == NULL && i != 0 && j != 0) 01071 throw NoMemory(string("Matrix_ComplexSparse::") 01072 + "Matrix_ComplexSparse(int, int, int, int)", 01073 string("Unable to allocate ") 01074 + to_str(sizeof(int) * real_nz) 01075 + " bytes to store " + to_str(real_nz) 01076 + " values (real part), for a " 01077 + to_str(i) + " by " + to_str(j) + " matrix."); 01078 #endif 01079 01080 #ifdef SELDON_CHECK_MEMORY 01081 try 01082 { 01083 #endif 01084 01085 this->imag_data_ = this->allocator_.allocate(imag_nz_, this); 01086 01087 #ifdef SELDON_CHECK_MEMORY 01088 } 01089 catch (...) 01090 { 01091 this->m_ = 0; 01092 this->n_ = 0; 01093 free(real_ptr_); 01094 free(imag_ptr_); 01095 real_ptr_ = NULL; 01096 imag_ptr_ = NULL; 01097 free(real_ind_); 01098 free(imag_ind_); 01099 real_ind_ = NULL; 01100 imag_ind_ = NULL; 01101 this->allocator_.deallocate(this->real_data_, real_nz_); 01102 this->real_data_ = NULL; 01103 this->imag_data_ = NULL; 01104 } 01105 if (real_data_ == NULL) 01106 { 01107 this->m_ = 0; 01108 this->n_ = 0; 01109 free(real_ptr_); 01110 free(imag_ptr_); 01111 real_ptr_ = NULL; 01112 imag_ptr_ = NULL; 01113 free(real_ind_); 01114 free(imag_ind_); 01115 real_ind_ = NULL; 01116 imag_ind_ = NULL; 01117 this->allocator_.deallocate(this->real_data_, real_nz_); 01118 real_data_ = NULL; 01119 } 01120 if (imag_data_ == NULL && i != 0 && j != 0) 01121 throw NoMemory(string("Matrix_ComplexSparse::") 01122 + "Matrix_ComplexSparse(int, int, int, int)", 01123 string("Unable to allocate ") 01124 + to_str(sizeof(int) * imag_nz) 01125 + " bytes to store " + to_str(imag_nz) 01126 + " values (imaginary part), for a " 01127 + to_str(i) + " by " + to_str(j) + " matrix."); 01128 #endif 01129 } 01130 01131 01133 01137 template <class T, class Prop, class Storage, class Allocator> 01138 void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 01139 ::Resize(int i, int j) 01140 { 01141 if (Storage::GetFirst(i, j) < Storage::GetFirst(this->m_, this->n_)) 01142 Resize(i, j, real_ptr_[Storage::GetFirst(i, j)], 01143 imag_ptr_[Storage::GetFirst(i, j)]); 01144 else 01145 Resize(i, j, real_nz_, imag_nz_); 01146 } 01147 01148 01150 01155 template <class T, class Prop, class Storage, class Allocator> 01156 void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 01157 ::Resize(int i, int j, int real_nz, int imag_nz) 01158 { 01159 #ifdef SELDON_CHECK_DIMENSIONS 01160 if (real_nz < 0 || imag_nz < 0) 01161 { 01162 this->m_ = 0; 01163 this->n_ = 0; 01164 real_nz_ = 0; 01165 imag_nz_ = 0; 01166 real_ptr_ = NULL; 01167 imag_ptr_ = NULL; 01168 real_ind_ = NULL; 01169 imag_ind_ = NULL; 01170 this->real_data_ = NULL; 01171 this->imag_data_ = NULL; 01172 throw WrongDim(string("Matrix_ComplexSparse::") 01173 + "Resize(int, int, int, int)", 01174 "Invalid number of non-zero elements: " 01175 + to_str(real_nz) + " in the real part and " 01176 + to_str(imag_nz) + " in the imaginary part."); 01177 } 01178 if ((real_nz > 0 01179 && (j == 0 01180 || static_cast<long int>(real_nz-1) / static_cast<long int>(j) 01181 >= static_cast<long int>(i))) 01182 || 01183 (imag_nz > 0 01184 && (j == 0 01185 || static_cast<long int>(imag_nz-1) / static_cast<long int>(j) 01186 >= static_cast<long int>(i)))) 01187 { 01188 this->m_ = 0; 01189 this->n_ = 0; 01190 real_nz_ = 0; 01191 imag_nz_ = 0; 01192 real_ptr_ = NULL; 01193 imag_ptr_ = NULL; 01194 real_ind_ = NULL; 01195 imag_ind_ = NULL; 01196 this->real_data_ = NULL; 01197 this->imag_data_ = NULL; 01198 throw WrongDim(string("Matrix_ComplexSparse::") 01199 + "Resize(int, int, int, int)", 01200 string("There are more values (") + to_str(real_nz) 01201 + " values for the real part and " + to_str(imag_nz) 01202 + string(" values for the imaginary part) than") 01203 + " elements in the matrix (" + to_str(i) + " by " 01204 + to_str(j) + ")."); 01205 } 01206 #endif 01207 01208 if (Storage::GetFirst(this->m_, this->n_) != Storage::GetFirst(i, j)) 01209 { 01210 #ifdef SELDON_CHECK_MEMORY 01211 try 01212 { 01213 #endif 01214 01215 real_ptr_ 01216 = reinterpret_cast<int*>( realloc(real_ptr_, 01217 (Storage::GetFirst(i, j)+1) 01218 *sizeof(int)) ); 01219 01220 #ifdef SELDON_CHECK_MEMORY 01221 } 01222 catch (...) 01223 { 01224 this->m_ = 0; 01225 this->n_ = 0; 01226 real_nz_ = 0; 01227 imag_nz_ = 0; 01228 real_ptr_ = NULL; 01229 imag_ptr_ = NULL; 01230 real_ind_ = NULL; 01231 imag_ind_ = NULL; 01232 this->real_data_ = NULL; 01233 this->imag_data_ = NULL; 01234 } 01235 if (real_ptr_ == NULL) 01236 { 01237 this->m_ = 0; 01238 this->n_ = 0; 01239 real_nz_ = 0; 01240 imag_nz_ = 0; 01241 imag_ptr_ = 0; 01242 real_ind_ = NULL; 01243 imag_ind_ = NULL; 01244 this->real_data_ = NULL; 01245 this->imag_data_ = NULL; 01246 } 01247 if (real_ptr_ == NULL && i != 0 && j != 0) 01248 throw NoMemory(string("Matrix_ComplexSparse::") 01249 + "Resize(int, int, int, int)", 01250 string("Unable to allocate ") 01251 + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1)) 01252 + " bytes to store " + to_str(Storage::GetFirst(i, j)+1) 01253 + string(" row or column start indices (for the real") 01254 + " part), for a " 01255 + to_str(i) + " by " + to_str(j) + " matrix."); 01256 #endif 01257 01258 #ifdef SELDON_CHECK_MEMORY 01259 try 01260 { 01261 #endif 01262 01263 imag_ptr_ 01264 = reinterpret_cast<int*>( realloc(imag_ptr_, 01265 (Storage::GetFirst(i, j)+1) 01266 *sizeof(int)) ); 01267 01268 #ifdef SELDON_CHECK_MEMORY 01269 } 01270 catch (...) 01271 { 01272 this->m_ = 0; 01273 this->n_ = 0; 01274 real_nz_ = 0; 01275 imag_nz_ = 0; 01276 free(real_ptr_); 01277 real_ptr_ = NULL; 01278 imag_ptr_ = NULL; 01279 real_ind_ = NULL; 01280 imag_ind_ = NULL; 01281 this->real_data_ = NULL; 01282 this->imag_data_ = NULL; 01283 } 01284 if (imag_ptr_ == NULL) 01285 { 01286 this->m_ = 0; 01287 this->n_ = 0; 01288 real_nz_ = 0; 01289 imag_nz_ = 0; 01290 free(real_ptr_); 01291 real_ptr_ = 0; 01292 real_ind_ = NULL; 01293 imag_ind_ = NULL; 01294 this->real_data_ = NULL; 01295 this->imag_data_ = NULL; 01296 } 01297 if (imag_ptr_ == NULL && i != 0 && j != 0) 01298 throw NoMemory(string("Matrix_ComplexSparse::") 01299 + "Resize(int, int, int, int)", 01300 string("Unable to allocate ") 01301 + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1)) 01302 + " bytes to store " + to_str(Storage::GetFirst(i, j)+1) 01303 + string(" row or column start indices (for the") 01304 + string(" imaginary part), for a ") 01305 + to_str(i) + " by " + to_str(j) + " matrix."); 01306 #endif 01307 } 01308 01309 if (real_nz != real_nz_) 01310 { 01311 #ifdef SELDON_CHECK_MEMORY 01312 try 01313 { 01314 #endif 01315 01316 real_ind_ 01317 = reinterpret_cast<int*>( realloc(real_ind_, 01318 real_nz*sizeof(int)) ); 01319 01320 #ifdef SELDON_CHECK_MEMORY 01321 } 01322 catch (...) 01323 { 01324 this->m_ = 0; 01325 this->n_ = 0; 01326 real_nz_ = 0; 01327 imag_nz_ = 0; 01328 free(real_ptr_); 01329 free(imag_ptr_); 01330 real_ptr_ = NULL; 01331 imag_ptr_ = NULL; 01332 real_ind_ = NULL; 01333 imag_ind_ = NULL; 01334 this->real_data_ = NULL; 01335 this->imag_data_ = NULL; 01336 } 01337 if (real_ind_ == NULL) 01338 { 01339 this->m_ = 0; 01340 this->n_ = 0; 01341 real_nz_ = 0; 01342 imag_nz_ = 0; 01343 free(real_ptr_); 01344 free(imag_ptr_); 01345 real_ptr_ = NULL; 01346 imag_ptr_ = NULL; 01347 this->real_data_ = NULL; 01348 this->imag_data_ = NULL; 01349 } 01350 if (real_ind_ == NULL && i != 0 && j != 0) 01351 throw NoMemory(string("Matrix_ComplexSparse::") 01352 + "Resize(int, int, int, int)", 01353 string("Unable to allocate ") 01354 + to_str(sizeof(int) * real_nz) 01355 + " bytes to store " + to_str(real_nz) 01356 + " row or column indices (real part), for a " 01357 + to_str(i) + " by " + to_str(j) + " matrix."); 01358 #endif 01359 } 01360 01361 if (imag_nz != imag_nz_) 01362 { 01363 #ifdef SELDON_CHECK_MEMORY 01364 try 01365 { 01366 #endif 01367 01368 imag_ind_ 01369 = reinterpret_cast<int*>( realloc(imag_ind_, 01370 imag_nz*sizeof(int)) ); 01371 01372 #ifdef SELDON_CHECK_MEMORY 01373 } 01374 catch (...) 01375 { 01376 this->m_ = 0; 01377 this->n_ = 0; 01378 real_nz_ = 0; 01379 imag_nz_ = 0; 01380 free(real_ptr_); 01381 free(imag_ptr_); 01382 real_ptr_ = NULL; 01383 imag_ptr_ = NULL; 01384 free(imag_ind_); 01385 real_ind_ = NULL; 01386 imag_ind_ = NULL; 01387 this->real_data_ = NULL; 01388 this->imag_data_ = NULL; 01389 } 01390 if (real_ind_ == NULL) 01391 { 01392 this->m_ = 0; 01393 this->n_ = 0; 01394 real_nz_ = 0; 01395 imag_nz_ = 0; 01396 free(real_ptr_); 01397 free(imag_ptr_); 01398 real_ptr_ = NULL; 01399 imag_ptr_ = NULL; 01400 free(imag_ind_); 01401 imag_ind_ = NULL; 01402 this->real_data_ = NULL; 01403 this->imag_data_ = NULL; 01404 } 01405 if (imag_ind_ == NULL && i != 0 && j != 0) 01406 throw NoMemory(string("Matrix_ComplexSparse::") 01407 + "Resize(int, int, int, int)", 01408 string("Unable to allocate ") 01409 + to_str(sizeof(int) * imag_nz) 01410 + " bytes to store " + to_str(imag_nz) 01411 + " row or column indices (imaginary part), for a " 01412 + to_str(i) + " by " + to_str(j) + " matrix."); 01413 #endif 01414 } 01415 01416 if (real_nz != real_nz_) 01417 { 01418 Vector<T, VectFull, Allocator> val; 01419 val.SetData(real_nz_, real_data_); 01420 val.Resize(real_nz); 01421 01422 real_data_ = val.GetData(); 01423 val.Nullify(); 01424 } 01425 01426 if (imag_nz != imag_nz_) 01427 { 01428 Vector<T, VectFull, Allocator> val; 01429 val.SetData(imag_nz_, imag_data_); 01430 val.Resize(imag_nz); 01431 01432 imag_data_ = val.GetData(); 01433 val.Nullify(); 01434 } 01435 01436 // then filing last values of ptr_ with nz_ 01437 for (int k = Storage::GetFirst(this->m_, this->n_); 01438 k <= Storage::GetFirst(i, j); k++) 01439 { 01440 real_ptr_[k] = real_nz_; 01441 imag_ptr_[k] = imag_nz_; 01442 } 01443 01444 this->m_ = i; 01445 this->n_ = j; 01446 real_nz_ = real_nz; 01447 imag_nz_ = imag_nz; 01448 } 01449 01450 01452 template <class T, class Prop, class Storage, class Allocator> 01453 inline void Matrix_ComplexSparse<T, Prop, Storage, Allocator>:: 01454 Copy(const Matrix_ComplexSparse<T, Prop, Storage, Allocator>& A) 01455 { 01456 this->Clear(); 01457 int i = A.m_; 01458 int j = A.n_; 01459 real_nz_ = A.real_nz_; 01460 imag_nz_ = A.imag_nz_; 01461 this->m_ = i; 01462 this->n_ = j; 01463 if ((i == 0)||(j == 0)) 01464 { 01465 this->m_ = 0; 01466 this->n_ = 0; 01467 this->real_nz_ = 0; 01468 this->imag_nz_ = 0; 01469 return; 01470 } 01471 01472 #ifdef SELDON_CHECK_DIMENSIONS 01473 if ((real_nz_ > 0 01474 && (j == 0 01475 || static_cast<long int>(real_nz_-1) / static_cast<long int>(j) 01476 >= static_cast<long int>(i))) 01477 || 01478 (imag_nz_ > 0 01479 && (j == 0 01480 || static_cast<long int>(imag_nz_-1) / static_cast<long int>(j) 01481 >= static_cast<long int>(i)))) 01482 { 01483 this->m_ = 0; 01484 this->n_ = 0; 01485 real_nz_ = 0; 01486 imag_nz_ = 0; 01487 real_ptr_ = NULL; 01488 imag_ptr_ = NULL; 01489 real_ind_ = NULL; 01490 imag_ind_ = NULL; 01491 this->real_data_ = NULL; 01492 this->imag_data_ = NULL; 01493 throw WrongDim(string("Matrix_ComplexSparse::") 01494 + "Matrix_ComplexSparse(int, int, int, int)", 01495 string("There are more values (") + to_str(real_nz_) 01496 + " values for the real part and " + to_str(imag_nz_) 01497 + string(" values for the imaginary part) than") 01498 + " elements in the matrix (" + to_str(i) + " by " 01499 + to_str(j) + ")."); 01500 } 01501 #endif 01502 01503 #ifdef SELDON_CHECK_MEMORY 01504 try 01505 { 01506 #endif 01507 01508 real_ptr_ = reinterpret_cast<int*>( calloc(Storage::GetFirst(i, j)+1, 01509 sizeof(int)) ); 01510 memcpy(this->real_ptr_, A.real_ptr_, 01511 (Storage::GetFirst(i, j) + 1) * sizeof(int)); 01512 01513 #ifdef SELDON_CHECK_MEMORY 01514 } 01515 catch (...) 01516 { 01517 this->m_ = 0; 01518 this->n_ = 0; 01519 real_nz_ = 0; 01520 imag_nz_ = 0; 01521 real_ptr_ = NULL; 01522 imag_ptr_ = NULL; 01523 real_ind_ = NULL; 01524 imag_ind_ = NULL; 01525 this->real_data_ = NULL; 01526 this->imag_data_ = NULL; 01527 } 01528 if (real_ptr_ == NULL) 01529 { 01530 this->m_ = 0; 01531 this->n_ = 0; 01532 real_nz_ = 0; 01533 imag_nz_ = 0; 01534 imag_ptr_ = 0; 01535 real_ind_ = NULL; 01536 imag_ind_ = NULL; 01537 this->real_data_ = NULL; 01538 this->imag_data_ = NULL; 01539 } 01540 if (real_ptr_ == NULL && i != 0 && j != 0) 01541 throw NoMemory(string("Matrix_ComplexSparse::") 01542 + "Matrix_ComplexSparse(int, int, int, int)", 01543 string("Unable to allocate ") 01544 + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1)) 01545 + " bytes to store " + to_str(Storage::GetFirst(i, j)+1) 01546 + string(" row or column start indices (for the real") 01547 + " part), for a " 01548 + to_str(i) + " by " + to_str(j) + " matrix."); 01549 #endif 01550 01551 #ifdef SELDON_CHECK_MEMORY 01552 try 01553 { 01554 #endif 01555 01556 imag_ptr_ = reinterpret_cast<int*>( calloc(Storage::GetFirst(i, j)+1, 01557 sizeof(int)) ); 01558 memcpy(this->imag_ptr_, A.imag_ptr_, 01559 (Storage::GetFirst(i, j) + 1) * sizeof(int)); 01560 01561 #ifdef SELDON_CHECK_MEMORY 01562 } 01563 catch (...) 01564 { 01565 this->m_ = 0; 01566 this->n_ = 0; 01567 real_nz_ = 0; 01568 imag_nz_ = 0; 01569 free(real_ptr_); 01570 real_ptr_ = NULL; 01571 imag_ptr_ = NULL; 01572 real_ind_ = NULL; 01573 imag_ind_ = NULL; 01574 this->real_data_ = NULL; 01575 this->imag_data_ = NULL; 01576 } 01577 if (imag_ptr_ == NULL) 01578 { 01579 this->m_ = 0; 01580 this->n_ = 0; 01581 real_nz_ = 0; 01582 imag_nz_ = 0; 01583 free(real_ptr_); 01584 real_ptr_ = 0; 01585 real_ind_ = NULL; 01586 imag_ind_ = NULL; 01587 this->real_data_ = NULL; 01588 this->imag_data_ = NULL; 01589 } 01590 if (imag_ptr_ == NULL && i != 0 && j != 0) 01591 throw NoMemory(string("Matrix_ComplexSparse::") 01592 + "Matrix_ComplexSparse(int, int, int, int)", 01593 string("Unable to allocate ") 01594 + to_str(sizeof(int) * (Storage::GetFirst(i, j)+1)) 01595 + " bytes to store " + to_str(Storage::GetFirst(i, j)+1) 01596 + string(" row or column start indices (for the") 01597 + string(" imaginary part), for a ") 01598 + to_str(i) + " by " + to_str(j) + " matrix."); 01599 #endif 01600 01601 #ifdef SELDON_CHECK_MEMORY 01602 try 01603 { 01604 #endif 01605 01606 real_ind_ = reinterpret_cast<int*>( calloc(real_nz_, sizeof(int)) ); 01607 memcpy(this->real_ind_, A.real_ind_, real_nz_ * sizeof(int)); 01608 01609 #ifdef SELDON_CHECK_MEMORY 01610 } 01611 catch (...) 01612 { 01613 this->m_ = 0; 01614 this->n_ = 0; 01615 real_nz_ = 0; 01616 imag_nz_ = 0; 01617 free(real_ptr_); 01618 free(imag_ptr_); 01619 real_ptr_ = NULL; 01620 imag_ptr_ = NULL; 01621 real_ind_ = NULL; 01622 imag_ind_ = NULL; 01623 this->real_data_ = NULL; 01624 this->imag_data_ = NULL; 01625 } 01626 if (real_ind_ == NULL) 01627 { 01628 this->m_ = 0; 01629 this->n_ = 0; 01630 real_nz_ = 0; 01631 imag_nz_ = 0; 01632 free(real_ptr_); 01633 free(imag_ptr_); 01634 real_ptr_ = NULL; 01635 imag_ptr_ = NULL; 01636 this->real_data_ = NULL; 01637 this->imag_data_ = NULL; 01638 } 01639 if (real_ind_ == NULL && i != 0 && j != 0) 01640 throw NoMemory(string("Matrix_ComplexSparse::") 01641 + "Matrix_ComplexSparse(int, int, int, int)", 01642 string("Unable to allocate ") 01643 + to_str(sizeof(int) * real_nz_) 01644 + " bytes to store " + to_str(real_nz_) 01645 + " row or column indices (real part), for a " 01646 + to_str(i) + " by " + to_str(j) + " matrix."); 01647 #endif 01648 01649 #ifdef SELDON_CHECK_MEMORY 01650 try 01651 { 01652 #endif 01653 01654 imag_ind_ = reinterpret_cast<int*>( calloc(imag_nz_, sizeof(int)) ); 01655 memcpy(this->imag_ind_, A.imag_ind_, imag_nz_ * sizeof(int)); 01656 01657 #ifdef SELDON_CHECK_MEMORY 01658 } 01659 catch (...) 01660 { 01661 this->m_ = 0; 01662 this->n_ = 0; 01663 real_nz_ = 0; 01664 imag_nz_ = 0; 01665 free(real_ptr_); 01666 free(imag_ptr_); 01667 real_ptr_ = NULL; 01668 imag_ptr_ = NULL; 01669 free(imag_ind_); 01670 real_ind_ = NULL; 01671 imag_ind_ = NULL; 01672 this->real_data_ = NULL; 01673 this->imag_data_ = NULL; 01674 } 01675 if (real_ind_ == NULL) 01676 { 01677 this->m_ = 0; 01678 this->n_ = 0; 01679 real_nz_ = 0; 01680 imag_nz_ = 0; 01681 free(real_ptr_); 01682 free(imag_ptr_); 01683 real_ptr_ = NULL; 01684 imag_ptr_ = NULL; 01685 free(imag_ind_); 01686 imag_ind_ = NULL; 01687 this->real_data_ = NULL; 01688 this->imag_data_ = NULL; 01689 } 01690 if (imag_ind_ == NULL && i != 0 && j != 0) 01691 throw NoMemory(string("Matrix_ComplexSparse::") 01692 + "Matrix_ComplexSparse(int, int, int, int)", 01693 string("Unable to allocate ") 01694 + to_str(sizeof(int) * imag_nz_) 01695 + " bytes to store " + to_str(imag_nz_) 01696 + " row or column indices (imaginary part), for a " 01697 + to_str(i) + " by " + to_str(j) + " matrix."); 01698 #endif 01699 01700 #ifdef SELDON_CHECK_MEMORY 01701 try 01702 { 01703 #endif 01704 01705 this->real_data_ = this->allocator_.allocate(real_nz_, this); 01706 this->allocator_.memorycpy(this->real_data_, A.real_data_, real_nz_); 01707 01708 #ifdef SELDON_CHECK_MEMORY 01709 } 01710 catch (...) 01711 { 01712 this->m_ = 0; 01713 this->n_ = 0; 01714 free(real_ptr_); 01715 free(imag_ptr_); 01716 real_ptr_ = NULL; 01717 imag_ptr_ = NULL; 01718 free(real_ind_); 01719 free(imag_ind_); 01720 real_ind_ = NULL; 01721 imag_ind_ = NULL; 01722 this->real_data_ = NULL; 01723 this->imag_data_ = NULL; 01724 } 01725 if (real_data_ == NULL) 01726 { 01727 this->m_ = 0; 01728 this->n_ = 0; 01729 free(real_ptr_); 01730 free(imag_ptr_); 01731 real_ptr_ = NULL; 01732 imag_ptr_ = NULL; 01733 free(real_ind_); 01734 free(imag_ind_); 01735 real_ind_ = NULL; 01736 imag_ind_ = NULL; 01737 imag_data_ = NULL; 01738 } 01739 if (real_data_ == NULL && i != 0 && j != 0) 01740 throw NoMemory(string("Matrix_ComplexSparse::") 01741 + "Matrix_ComplexSparse(int, int, int, int)", 01742 string("Unable to allocate ") 01743 + to_str(sizeof(int) * real_nz_) 01744 + " bytes to store " + to_str(real_nz_) 01745 + " values (real part), for a " 01746 + to_str(i) + " by " + to_str(j) + " matrix."); 01747 #endif 01748 01749 #ifdef SELDON_CHECK_MEMORY 01750 try 01751 { 01752 #endif 01753 01754 this->imag_data_ = this->allocator_.allocate(imag_nz_, this); 01755 this->allocator_.memorycpy(this->imag_data_, A.imag_data_, imag_nz_); 01756 01757 #ifdef SELDON_CHECK_MEMORY 01758 } 01759 catch (...) 01760 { 01761 this->m_ = 0; 01762 this->n_ = 0; 01763 free(real_ptr_); 01764 free(imag_ptr_); 01765 real_ptr_ = NULL; 01766 imag_ptr_ = NULL; 01767 free(real_ind_); 01768 free(imag_ind_); 01769 real_ind_ = NULL; 01770 imag_ind_ = NULL; 01771 this->allocator_.deallocate(this->real_data_, real_nz_); 01772 this->real_data_ = NULL; 01773 this->imag_data_ = NULL; 01774 } 01775 if (real_data_ == NULL) 01776 { 01777 this->m_ = 0; 01778 this->n_ = 0; 01779 free(real_ptr_); 01780 free(imag_ptr_); 01781 real_ptr_ = NULL; 01782 imag_ptr_ = NULL; 01783 free(real_ind_); 01784 free(imag_ind_); 01785 real_ind_ = NULL; 01786 imag_ind_ = NULL; 01787 this->allocator_.deallocate(this->real_data_, real_nz_); 01788 real_data_ = NULL; 01789 } 01790 if (imag_data_ == NULL && i != 0 && j != 0) 01791 throw NoMemory(string("Matrix_ComplexSparse::") 01792 + "Matrix_ComplexSparse(int, int, int, int)", 01793 string("Unable to allocate ") 01794 + to_str(sizeof(int) * imag_nz_) 01795 + " bytes to store " + to_str(imag_nz_) 01796 + " values (imaginary part), for a " 01797 + to_str(i) + " by " + to_str(j) + " matrix."); 01798 #endif 01799 01800 } 01801 01802 01803 /******************* 01804 * BASIC FUNCTIONS * 01805 *******************/ 01806 01807 01809 01815 template <class T, class Prop, class Storage, class Allocator> 01816 int Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetDataSize() const 01817 { 01818 return real_nz_ + imag_nz_; 01819 } 01820 01821 01823 01827 template <class T, class Prop, class Storage, class Allocator> 01828 int* Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetRealPtr() const 01829 { 01830 return real_ptr_; 01831 } 01832 01833 01835 01839 template <class T, class Prop, class Storage, class Allocator> 01840 int* Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetImagPtr() const 01841 { 01842 return imag_ptr_; 01843 } 01844 01845 01847 01854 template <class T, class Prop, class Storage, class Allocator> 01855 int* Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetRealInd() const 01856 { 01857 return real_ind_; 01858 } 01859 01860 01863 01870 template <class T, class Prop, class Storage, class Allocator> 01871 int* Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetImagInd() const 01872 { 01873 return imag_ind_; 01874 } 01875 01876 01878 01881 template <class T, class Prop, class Storage, class Allocator> 01882 int Matrix_ComplexSparse<T, Prop, Storage, Allocator> 01883 ::GetRealPtrSize() const 01884 { 01885 return (Storage::GetFirst(this->m_, this->n_) + 1); 01886 } 01887 01888 01890 01893 template <class T, class Prop, class Storage, class Allocator> 01894 int Matrix_ComplexSparse<T, Prop, Storage, Allocator> 01895 ::GetImagPtrSize() const 01896 { 01897 return (Storage::GetFirst(this->m_, this->n_) + 1); 01898 } 01899 01900 01903 01912 template <class T, class Prop, class Storage, class Allocator> 01913 int Matrix_ComplexSparse<T, Prop, Storage, Allocator> 01914 ::GetRealIndSize() const 01915 { 01916 return real_nz_; 01917 } 01918 01919 01922 01931 template <class T, class Prop, class Storage, class Allocator> 01932 int Matrix_ComplexSparse<T, Prop, Storage, Allocator> 01933 ::GetImagIndSize() const 01934 { 01935 return imag_nz_; 01936 } 01937 01938 01941 01950 template <class T, class Prop, class Storage, class Allocator> 01951 int Matrix_ComplexSparse<T, Prop, Storage, Allocator> 01952 ::GetRealDataSize() const 01953 { 01954 return real_nz_; 01955 } 01956 01957 01960 01969 template <class T, class Prop, class Storage, class Allocator> 01970 int Matrix_ComplexSparse<T, Prop, Storage, Allocator> 01971 ::GetImagDataSize() const 01972 { 01973 return imag_nz_; 01974 } 01975 01976 01978 01981 template <class T, class Prop, class Storage, class Allocator> 01982 T* Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetRealData() const 01983 { 01984 return real_data_; 01985 } 01986 01987 01989 01992 template <class T, class Prop, class Storage, class Allocator> 01993 T* Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetImagData() const 01994 { 01995 return imag_data_; 01996 } 01997 01998 01999 /********************************** 02000 * ELEMENT ACCESS AND AFFECTATION * 02001 **********************************/ 02002 02003 02005 02011 template <class T, class Prop, class Storage, class Allocator> 02012 inline complex<typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02013 ::value_type> 02014 Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02015 ::operator() (int i, int j) const 02016 { 02017 02018 #ifdef SELDON_CHECK_BOUNDS 02019 if (i < 0 || i >= this->m_) 02020 throw WrongRow("Matrix_ComplexSparse::operator()", 02021 string("Index should be in [0, ") + to_str(this->m_-1) 02022 + "], but is equal to " + to_str(i) + "."); 02023 if (j < 0 || j >= this->n_) 02024 throw WrongCol("Matrix_ComplexSparse::operator()", 02025 string("Index should be in [0, ") + to_str(this->n_-1) 02026 + "], but is equal to " + to_str(j) + "."); 02027 #endif 02028 02029 int real_k, imag_k, l; 02030 int real_a, real_b; 02031 int imag_a, imag_b; 02032 02033 real_a = real_ptr_[Storage::GetFirst(i, j)]; 02034 real_b = real_ptr_[Storage::GetFirst(i, j) + 1]; 02035 02036 imag_a = imag_ptr_[Storage::GetFirst(i, j)]; 02037 imag_b = imag_ptr_[Storage::GetFirst(i, j) + 1]; 02038 02039 if (real_a != real_b) 02040 { 02041 l = Storage::GetSecond(i, j); 02042 for (real_k = real_a; 02043 (real_k <real_b-1) && (real_ind_[real_k] < l); 02044 real_k++); 02045 if (imag_a != imag_b) 02046 { 02047 for (imag_k = imag_a; 02048 (imag_k < imag_b-1) && (imag_ind_[imag_k] < l); 02049 imag_k++); 02050 if (real_ind_[real_k] == l) 02051 { 02052 if (imag_ind_[imag_k] == l) 02053 return complex<T>(real_data_[real_k], imag_data_[imag_k]); 02054 else 02055 return complex<T>(real_data_[real_k], T(0)); 02056 } 02057 else 02058 if (imag_ind_[imag_k] == l) 02059 return complex<T>(T(0), imag_data_[imag_k]); 02060 else 02061 return complex<T>(T(0), T(0)); 02062 } 02063 else 02064 { 02065 if (real_ind_[real_k] == l) 02066 return complex<T>(real_data_[real_k], T(0)); 02067 else 02068 return complex<T>(T(0), T(0)); 02069 } 02070 } 02071 else 02072 { 02073 if (imag_a != imag_b) 02074 { 02075 l = Storage::GetSecond(i, j); 02076 for (imag_k = imag_a; 02077 (imag_k < imag_b-1) && (imag_ind_[imag_k] < l); 02078 imag_k++); 02079 if (imag_ind_[imag_k] == l) 02080 return complex<T>(T(0), imag_data_[imag_k]); 02081 else 02082 return complex<T>(T(0), T(0)); 02083 } 02084 else 02085 return complex<T>(T(0), T(0)); 02086 } 02087 02088 } 02089 02090 02092 02100 template <class T, class Prop, class Storage, class Allocator> 02101 inline typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02102 ::value_type& 02103 Matrix_ComplexSparse<T, Prop, Storage, Allocator>::ValReal(int i, int j) 02104 { 02105 02106 #ifdef SELDON_CHECK_BOUNDS 02107 if (i < 0 || i >= this->m_) 02108 throw WrongRow("Matrix_ComplexSparse::ValReal(int, int)", 02109 string("Index should be in [0, ") + to_str(this->m_-1) 02110 + "], but is equal to " + to_str(i) + "."); 02111 if (j < 0 || j >= this->n_) 02112 throw WrongCol("Matrix_ComplexSparse::ValReal(int, int)", 02113 string("Index should be in [0, ") + to_str(this->n_-1) 02114 + "], but is equal to " + to_str(j) + "."); 02115 #endif 02116 02117 int k, l; 02118 int a, b; 02119 02120 a = real_ptr_[Storage::GetFirst(i, j)]; 02121 b = real_ptr_[Storage::GetFirst(i, j) + 1]; 02122 02123 if (a == b) 02124 throw WrongArgument("Matrix_ComplexSparse::ValReal(int, int)", 02125 "No reference to element (" + to_str(i) + ", " 02126 + to_str(j) 02127 + ") can be returned: it is a zero entry."); 02128 02129 l = Storage::GetSecond(i, j); 02130 02131 for (k = a; (k < b-1) && (real_ind_[k] < l); k++); 02132 02133 if (real_ind_[k] == l) 02134 return this->real_data_[k]; 02135 else 02136 throw WrongArgument("Matrix_ComplexSparse::ValReal(int, int)", 02137 "No reference to element (" + to_str(i) + ", " 02138 + to_str(j) 02139 + ") can be returned: it is a zero entry."); 02140 } 02141 02142 02144 02152 template <class T, class Prop, class Storage, class Allocator> 02153 inline const typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02154 ::value_type& 02155 Matrix_ComplexSparse<T, Prop, Storage, Allocator>::ValReal(int i, int j) const 02156 { 02157 02158 #ifdef SELDON_CHECK_BOUNDS 02159 if (i < 0 || i >= this->m_) 02160 throw WrongRow("Matrix_ComplexSparse::ValReal(int, int)", 02161 string("Index should be in [0, ") + to_str(this->m_-1) 02162 + "], but is equal to " + to_str(i) + "."); 02163 if (j < 0 || j >= this->n_) 02164 throw WrongCol("Matrix_ComplexSparse::ValReal(int, int)", 02165 string("Index should be in [0, ") + to_str(this->n_-1) 02166 + "], but is equal to " + to_str(j) + "."); 02167 #endif 02168 02169 int k, l; 02170 int a, b; 02171 02172 a = real_ptr_[Storage::GetFirst(i, j)]; 02173 b = real_ptr_[Storage::GetFirst(i, j) + 1]; 02174 02175 if (a == b) 02176 throw WrongArgument("Matrix_ComplexSparse::ValReal(int, int)", 02177 "No reference to element (" + to_str(i) + ", " 02178 + to_str(j) 02179 + ") can be returned: it is a zero entry."); 02180 02181 l = Storage::GetSecond(i, j); 02182 02183 for (k = a; (k < b-1) && (real_ind_[k] < l); k++); 02184 02185 if (real_ind_[k] == l) 02186 return this->real_data_[k]; 02187 else 02188 throw WrongArgument("Matrix_ComplexSparse::ValReal(int, int)", 02189 "No reference to element (" + to_str(i) + ", " 02190 + to_str(j) 02191 + ") can be returned: it is a zero entry."); 02192 } 02193 02194 02196 02204 template <class T, class Prop, class Storage, class Allocator> 02205 inline typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02206 ::value_type& 02207 Matrix_ComplexSparse<T, Prop, Storage, Allocator>::ValImag(int i, int j) 02208 { 02209 02210 #ifdef SELDON_CHECK_BOUNDS 02211 if (i < 0 || i >= this->m_) 02212 throw WrongRow("Matrix_ComplexSparse::ValImag(int, int)", 02213 string("Index should be in [0, ") + to_str(this->m_-1) 02214 + "], but is equal to " + to_str(i) + "."); 02215 if (j < 0 || j >= this->n_) 02216 throw WrongCol("Matrix_ComplexSparse::ValImag(int, int)", 02217 string("Index should be in [0, ") + to_str(this->n_-1) 02218 + "], but is equal to " + to_str(j) + "."); 02219 #endif 02220 02221 int k, l; 02222 int a, b; 02223 02224 a = imag_ptr_[Storage::GetFirst(i, j)]; 02225 b = imag_ptr_[Storage::GetFirst(i, j) + 1]; 02226 02227 if (a == b) 02228 throw WrongArgument("Matrix_ComplexSparse::ValImag(int, int)", 02229 "No reference to element (" + to_str(i) + ", " 02230 + to_str(j) 02231 + ") can be returned: it is a zero entry."); 02232 02233 l = Storage::GetSecond(i, j); 02234 02235 for (k = a; (k < b-1) && (imag_ind_[k] < l); k++); 02236 02237 if (imag_ind_[k] == l) 02238 return this->imag_data_[k]; 02239 else 02240 throw WrongArgument("Matrix_ComplexSparse::ValImag(int, int)", 02241 "No reference to element (" + to_str(i) + ", " 02242 + to_str(j) 02243 + ") can be returned: it is a zero entry."); 02244 } 02245 02246 02248 02256 template <class T, class Prop, class Storage, class Allocator> 02257 inline const typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02258 ::value_type& 02259 Matrix_ComplexSparse<T, Prop, Storage, Allocator>::ValImag(int i, int j) const 02260 { 02261 02262 #ifdef SELDON_CHECK_BOUNDS 02263 if (i < 0 || i >= this->m_) 02264 throw WrongRow("Matrix_ComplexSparse::ValImag(int, int)", 02265 string("Index should be in [0, ") + to_str(this->m_-1) 02266 + "], but is equal to " + to_str(i) + "."); 02267 if (j < 0 || j >= this->n_) 02268 throw WrongCol("Matrix_ComplexSparse::ValImag(int, int)", 02269 string("Index should be in [0, ") + to_str(this->n_-1) 02270 + "], but is equal to " + to_str(j) + "."); 02271 #endif 02272 02273 int k, l; 02274 int a, b; 02275 02276 a = imag_ptr_[Storage::GetFirst(i, j)]; 02277 b = imag_ptr_[Storage::GetFirst(i, j) + 1]; 02278 02279 if (a == b) 02280 throw WrongArgument("Matrix_ComplexSparse::ValImag(int, int)", 02281 "No reference to element (" + to_str(i) + ", " 02282 + to_str(j) 02283 + ") can be returned: it is a zero entry."); 02284 02285 l = Storage::GetSecond(i, j); 02286 02287 for (k = a; (k < b-1) && (imag_ind_[k] < l); k++); 02288 02289 if (imag_ind_[k] == l) 02290 return this->imag_data_[k]; 02291 else 02292 throw WrongArgument("Matrix_ComplexSparse::ValImag(int, int)", 02293 "No reference to element (" + to_str(i) + ", " 02294 + to_str(j) 02295 + ") can be returned: it is a zero entry."); 02296 } 02297 02298 02300 02307 template <class T, class Prop, class Storage, class Allocator> 02308 inline typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02309 ::value_type& 02310 Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetReal(int i, int j) 02311 { 02312 02313 #ifdef SELDON_CHECK_BOUNDS 02314 if (i < 0 || i >= this->m_) 02315 throw WrongRow("Matrix_ComplexSparse::GetReal(int, int)", 02316 string("Index should be in [0, ") + to_str(this->m_-1) 02317 + "], but is equal to " + to_str(i) + "."); 02318 if (j < 0 || j >= this->n_) 02319 throw WrongCol("Matrix_ComplexSparse::GetReal(int, int)", 02320 string("Index should be in [0, ") + to_str(this->n_-1) 02321 + "], but is equal to " + to_str(j) + "."); 02322 #endif 02323 02324 int k, l; 02325 int a, b; 02326 02327 a = real_ptr_[Storage::GetFirst(i, j)]; 02328 b = real_ptr_[Storage::GetFirst(i, j) + 1]; 02329 02330 if (a < b) 02331 { 02332 l = Storage::GetSecond(i, j); 02333 02334 for (k = a; (k < b) && (real_ind_[k] < l); k++); 02335 02336 if ( (k < b) && (real_ind_[k] == l)) 02337 return this->real_data_[k]; 02338 } 02339 else 02340 k = a; 02341 02342 // adding a non-zero entry 02343 Resize(this->m_, this->n_, real_nz_+1, imag_nz_); 02344 02345 for (int m = Storage::GetFirst(i, j)+1; 02346 m <= Storage::GetFirst(this->m_, this->n_); m++) 02347 real_ptr_[m]++; 02348 02349 for (int m = real_nz_-1; m >= k+1; m--) 02350 { 02351 real_ind_[m] = real_ind_[m-1]; 02352 this->real_data_[m] = this->real_data_[m-1]; 02353 } 02354 02355 real_ind_[k] = Storage::GetSecond(i, j); 02356 02357 // value of new non-zero entry is set to 0 02358 SetComplexZero(this->real_data_[k]); 02359 02360 return this->real_data_[k]; 02361 } 02362 02363 02365 02372 template <class T, class Prop, class Storage, class Allocator> 02373 inline const typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02374 ::value_type& 02375 Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetReal(int i, int j) const 02376 { 02377 return ValReal(i, j); 02378 } 02379 02380 02382 02389 template <class T, class Prop, class Storage, class Allocator> 02390 inline typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02391 ::value_type& 02392 Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetImag(int i, int j) 02393 { 02394 02395 #ifdef SELDON_CHECK_BOUNDS 02396 if (i < 0 || i >= this->m_) 02397 throw WrongRow("Matrix_ComplexSparse::GetImag(int, int)", 02398 string("Index should be in [0, ") + to_str(this->m_-1) 02399 + "], but is equal to " + to_str(i) + "."); 02400 if (j < 0 || j >= this->n_) 02401 throw WrongCol("Matrix_ComplexSparse::GetImag(int, int)", 02402 string("Index should be in [0, ") + to_str(this->n_-1) 02403 + "], but is equal to " + to_str(j) + "."); 02404 #endif 02405 02406 int k, l; 02407 int a, b; 02408 02409 a = imag_ptr_[Storage::GetFirst(i, j)]; 02410 b = imag_ptr_[Storage::GetFirst(i, j) + 1]; 02411 02412 if (a < b) 02413 { 02414 l = Storage::GetSecond(i, j); 02415 02416 for (k = a; (k < b) && (imag_ind_[k] < l); k++); 02417 02418 if ( (k < b) && (imag_ind_[k] == l)) 02419 return this->imag_data_[k]; 02420 } 02421 else 02422 k = a; 02423 02424 // adding a non-zero entry 02425 Resize(this->m_, this->n_, real_nz_, imag_nz_+1); 02426 02427 for (int m = Storage::GetFirst(i, j)+1; 02428 m <= Storage::GetFirst(this->m_, this->n_); m++) 02429 imag_ptr_[m]++; 02430 02431 for (int m = imag_nz_-1; m >= k+1; m--) 02432 { 02433 imag_ind_[m] = imag_ind_[m-1]; 02434 this->imag_data_[m] = this->imag_data_[m-1]; 02435 } 02436 02437 imag_ind_[k] = Storage::GetSecond(i, j); 02438 02439 // value of new non-zero entry is set to 0 02440 SetComplexZero(this->imag_data_[k]); 02441 02442 return this->imag_data_[k]; 02443 } 02444 02445 02447 02454 template <class T, class Prop, class Storage, class Allocator> 02455 inline const typename Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02456 ::value_type& 02457 Matrix_ComplexSparse<T, Prop, Storage, Allocator>::GetImag(int i, int j) const 02458 { 02459 return ValImag(i, j); 02460 } 02461 02462 02464 02471 template <class T, class Prop, class Storage, class Allocator> 02472 inline void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02473 ::AddInteraction(int i, int j, const complex<T>& val) 02474 { 02475 if (real(val) != T(0)) 02476 GetReal(i, j) += real(val); 02477 02478 if (imag(val) != T(0)) 02479 GetImag(i, j) += imag(val); 02480 } 02481 02482 02484 02489 template <class T, class Prop, class Storage, class Allocator> 02490 inline void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02491 ::Set(int i, int j, const complex<T>& val) 02492 { 02493 GetReal(i, j) = real(val); 02494 GetImag(i, j) = imag(val); 02495 } 02496 02497 02499 02504 template <class T, class Prop, class Storage, class Allocator> 02505 inline Matrix_ComplexSparse<T, Prop, Storage, Allocator>& 02506 Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02507 ::operator= (const Matrix_ComplexSparse<T, Prop, Storage, Allocator>& A) 02508 { 02509 this->Copy(A); 02510 02511 return *this; 02512 } 02513 02514 02515 /************************ 02516 * CONVENIENT FUNCTIONS * 02517 ************************/ 02518 02519 02521 02522 template <class T, class Prop, class Storage, class Allocator> 02523 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>::Zero() 02524 { 02525 this->allocator_.memoryset(this->real_data_, char(0), 02526 this->real_nz_ * sizeof(value_type)); 02527 02528 this->allocator_.memoryset(this->imag_data_, char(0), 02529 this->imag_nz_ * sizeof(value_type)); 02530 } 02531 02532 02534 02537 template <class T, class Prop, class Storage, class Allocator> 02538 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>::SetIdentity() 02539 { 02540 int m = this->m_; 02541 int n = this->n_; 02542 int nz = min(m, n); 02543 02544 if (nz == 0) 02545 return; 02546 02547 Clear(); 02548 02549 Vector<T, VectFull, Allocator> real_values(nz), imag_values; 02550 Vector<int, VectFull, CallocAlloc<int> > 02551 real_ptr(Storage::GetFirst(m, n) + 1); 02552 Vector<int, VectFull, CallocAlloc<int> > real_ind(nz); 02553 Vector<int, VectFull, CallocAlloc<int> > imag_ptr(real_ptr); 02554 Vector<int, VectFull, CallocAlloc<int> > imag_ind; 02555 02556 real_values.Fill(T(1)); 02557 real_ind.Fill(); 02558 imag_ind.Zero(); 02559 int i; 02560 for (i = 0; i < nz + 1; i++) 02561 real_ptr(i) = i; 02562 for (i = nz + 1; i < real_ptr.GetLength(); i++) 02563 real_ptr(i) = nz; 02564 02565 SetData(m, n, real_values, real_ptr, real_ind, 02566 imag_values, imag_ptr, imag_ind); 02567 } 02568 02569 02571 02574 template <class T, class Prop, class Storage, class Allocator> 02575 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>::Fill() 02576 { 02577 for (int i = 0; i < this->real_nz_; i++) 02578 this->real_data_[i] = i; 02579 02580 for (int i = 0; i < this->imag_nz_; i++) 02581 this->imag_data_[i] = T(0); 02582 } 02583 02584 02586 02589 template <class T, class Prop, class Storage, class Allocator> 02590 void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02591 ::Fill(const complex<T>& x) 02592 { 02593 for (int i = 0; i < this->real_nz_; i++) 02594 this->real_data_[i] = real(x); 02595 02596 for (int i = 0; i < this->imag_nz_; i++) 02597 this->imag_data_[i] = imag(x); 02598 } 02599 02600 02602 02605 template <class T, class Prop, class Storage, class Allocator> 02606 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>::FillRand() 02607 { 02608 srand(time(NULL)); 02609 for (int i = 0; i < this->real_nz_; i++) 02610 this->real_data_[i] = rand(); 02611 02612 for (int i = 0; i < this->imag_nz_; i++) 02613 this->imag_data_[i] = rand(); 02614 } 02615 02616 02618 02623 template <class T, class Prop, class Storage, class Allocator> 02624 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>::Print() const 02625 { 02626 for (int i = 0; i < this->m_; i++) 02627 { 02628 for (int j = 0; j < this->n_; j++) 02629 cout << (*this)(i, j) << "\t"; 02630 cout << endl; 02631 } 02632 } 02633 02634 02636 02640 template <class T, class Prop, class Storage, class Allocator> 02641 void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02642 ::Write(string FileName) const 02643 { 02644 ofstream FileStream; 02645 FileStream.open(FileName.c_str()); 02646 02647 #ifdef SELDON_CHECK_IO 02648 // Checks if the file was opened. 02649 if (!FileStream.is_open()) 02650 throw IOError("Matrix_ComplexSparse::Write(string FileName)", 02651 string("Unable to open file \"") + FileName + "\"."); 02652 #endif 02653 02654 this->Write(FileStream); 02655 02656 FileStream.close(); 02657 } 02658 02659 02661 02665 template <class T, class Prop, class Storage, class Allocator> 02666 void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02667 ::Write(ostream& FileStream) const 02668 { 02669 #ifdef SELDON_CHECK_IO 02670 // Checks if the stream is ready. 02671 if (!FileStream.good()) 02672 throw IOError("Matrix_ComplexSparse::Write(ofstream& FileStream)", 02673 "Stream is not ready."); 02674 #endif 02675 02676 FileStream.write(reinterpret_cast<char*>(const_cast<int*>(&this->m_)), 02677 sizeof(int)); 02678 FileStream.write(reinterpret_cast<char*>(const_cast<int*>(&this->n_)), 02679 sizeof(int)); 02680 FileStream.write(reinterpret_cast<char*>(const_cast<int*>(&this->real_nz_)), 02681 sizeof(int)); 02682 FileStream.write(reinterpret_cast<char*>(const_cast<int*>(&this->imag_nz_)), 02683 sizeof(int)); 02684 02685 FileStream.write(reinterpret_cast<char*>(this->real_ptr_), 02686 sizeof(int)*(Storage::GetFirst(this->m_, this->n_)+1)); 02687 FileStream.write(reinterpret_cast<char*>(this->real_ind_), 02688 sizeof(int)*this->real_nz_); 02689 FileStream.write(reinterpret_cast<char*>(this->real_data_), 02690 sizeof(T)*this->real_nz_); 02691 02692 FileStream.write(reinterpret_cast<char*>(this->imag_ptr_), 02693 sizeof(int)*(Storage::GetFirst(this->m_, this->n_)+1)); 02694 FileStream.write(reinterpret_cast<char*>(this->imag_ind_), 02695 sizeof(int)*this->imag_nz_); 02696 FileStream.write(reinterpret_cast<char*>(this->imag_data_), 02697 sizeof(T)*this->imag_nz_); 02698 } 02699 02700 02702 02708 template <class T, class Prop, class Storage, class Allocator> 02709 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>:: 02710 WriteText(string FileName) const 02711 { 02712 ofstream FileStream; FileStream.precision(14); 02713 FileStream.open(FileName.c_str()); 02714 02715 #ifdef SELDON_CHECK_IO 02716 // Checks if the file was opened. 02717 if (!FileStream.is_open()) 02718 throw IOError("Matrix_ComplexSparse::Write(string FileName)", 02719 string("Unable to open file \"") + FileName + "\"."); 02720 #endif 02721 02722 this->WriteText(FileStream); 02723 02724 FileStream.close(); 02725 } 02726 02727 02729 02735 template <class T, class Prop, class Storage, class Allocator> 02736 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>:: 02737 WriteText(ostream& FileStream) const 02738 { 02739 02740 #ifdef SELDON_CHECK_IO 02741 // Checks if the stream is ready. 02742 if (!FileStream.good()) 02743 throw IOError("Matrix_ComplexSparse::Write(ofstream& FileStream)", 02744 "Stream is not ready."); 02745 #endif 02746 02747 // conversion in coordinate format (1-index convention) 02748 IVect IndRow, IndCol; Vector<complex<T> > Value; 02749 const Matrix<T, Prop, Storage, Allocator>& leaf_class = 02750 static_cast<const Matrix<T, Prop, Storage, Allocator>& >(*this); 02751 02752 ConvertMatrix_to_Coordinates(leaf_class, IndRow, IndCol, 02753 Value, 1); 02754 02755 for (int i = 0; i < IndRow.GetM(); i++) 02756 FileStream << IndRow(i) << " " << IndCol(i) << " " << Value(i) << '\n'; 02757 02758 } 02759 02760 02762 02766 template <class T, class Prop, class Storage, class Allocator> 02767 void Matrix_ComplexSparse<T, Prop, Storage, Allocator>:: 02768 Read(string FileName) 02769 { 02770 ifstream FileStream; 02771 FileStream.open(FileName.c_str()); 02772 02773 #ifdef SELDON_CHECK_IO 02774 // Checks if the file was opened. 02775 if (!FileStream.is_open()) 02776 throw IOError("Matrix_ComplexSparse::Read(string FileName)", 02777 string("Unable to open file \"") + FileName + "\"."); 02778 #endif 02779 02780 this->Read(FileStream); 02781 02782 FileStream.close(); 02783 } 02784 02785 02787 02791 template <class T, class Prop, class Storage, class Allocator> 02792 void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02793 ::Read(istream& FileStream) 02794 { 02795 02796 #ifdef SELDON_CHECK_IO 02797 // Checks if the stream is ready. 02798 if (!FileStream.good()) 02799 throw IOError("Matrix_ComplexSparse::Read(istream& FileStream)", 02800 "Stream is not ready."); 02801 #endif 02802 02803 int m, n, real_nz, imag_nz; 02804 FileStream.read(reinterpret_cast<char*>(&m), sizeof(int)); 02805 FileStream.read(reinterpret_cast<char*>(&n), sizeof(int)); 02806 FileStream.read(reinterpret_cast<char*>(&real_nz), sizeof(int)); 02807 FileStream.read(reinterpret_cast<char*>(&imag_nz), sizeof(int)); 02808 02809 Reallocate(m, n, real_nz, imag_nz); 02810 02811 FileStream.read(reinterpret_cast<char*>(real_ptr_), 02812 sizeof(int)*(Storage::GetFirst(m, n)+1)); 02813 FileStream.read(reinterpret_cast<char*>(real_ind_), sizeof(int)*real_nz); 02814 FileStream.read(reinterpret_cast<char*>(this->real_data_), sizeof(T)*real_nz); 02815 02816 FileStream.read(reinterpret_cast<char*>(imag_ptr_), 02817 sizeof(int)*(Storage::GetFirst(m, n)+1)); 02818 FileStream.read(reinterpret_cast<char*>(imag_ind_), sizeof(int)*imag_nz); 02819 FileStream.read(reinterpret_cast<char*>(this->imag_data_), sizeof(T)*imag_nz); 02820 02821 #ifdef SELDON_CHECK_IO 02822 // Checks if data was read. 02823 if (!FileStream.good()) 02824 throw IOError("Matrix_ComplexSparse::Read(istream& FileStream)", 02825 string("Input operation failed.") 02826 + string(" The input file may have been removed") 02827 + " or may not contain enough data."); 02828 #endif 02829 02830 } 02831 02832 02834 02838 template <class T, class Prop, class Storage, class Allocator> 02839 void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02840 ::ReadText(string FileName) 02841 { 02842 ifstream FileStream; 02843 FileStream.open(FileName.c_str()); 02844 02845 #ifdef SELDON_CHECK_IO 02846 // Checks if the file was opened. 02847 if (!FileStream.is_open()) 02848 throw IOError("Matrix_ComplexSparse::ReadText(string FileName)", 02849 string("Unable to open file \"") + FileName + "\"."); 02850 #endif 02851 02852 this->ReadText(FileStream); 02853 02854 FileStream.close(); 02855 } 02856 02857 02859 02863 template <class T, class Prop, class Storage, class Allocator> 02864 void Matrix_ComplexSparse<T, Prop, Storage, Allocator> 02865 ::ReadText(istream& FileStream) 02866 { 02867 Matrix<T, Prop, Storage, Allocator>& leaf_class = 02868 static_cast<Matrix<T, Prop, Storage, Allocator>& >(*this); 02869 02870 complex<T> zero; int index = 1; 02871 ReadCoordinateMatrix(leaf_class, FileStream, zero, index); 02872 } 02873 02874 02876 // MATRIX<COLCOMPLEXSPARSE> // 02878 02879 02880 /**************** 02881 * CONSTRUCTORS * 02882 ****************/ 02883 02885 02888 template <class T, class Prop, class Allocator> 02889 Matrix<T, Prop, ColComplexSparse, Allocator>::Matrix(): 02890 Matrix_ComplexSparse<T, Prop, ColComplexSparse, Allocator>() 02891 { 02892 } 02893 02894 02896 02900 template <class T, class Prop, class Allocator> 02901 Matrix<T, Prop, ColComplexSparse, Allocator>::Matrix(int i, int j): 02902 Matrix_ComplexSparse<T, Prop, ColComplexSparse, Allocator>(i, j, 0, 0) 02903 { 02904 } 02905 02906 02908 02916 template <class T, class Prop, class Allocator> 02917 Matrix<T, Prop, ColComplexSparse, Allocator>::Matrix(int i, int j, 02918 int real_nz, 02919 int imag_nz): 02920 Matrix_ComplexSparse<T, Prop, ColComplexSparse, Allocator>(i, j, 02921 real_nz, 02922 imag_nz) 02923 { 02924 } 02925 02926 02928 02946 template <class T, class Prop, class Allocator> 02947 template <class Storage0, class Allocator0, 02948 class Storage1, class Allocator1, 02949 class Storage2, class Allocator2> 02950 Matrix<T, Prop, ColComplexSparse, Allocator>:: 02951 Matrix(int i, int j, 02952 Vector<T, Storage0, Allocator0>& real_values, 02953 Vector<int, Storage1, Allocator1>& real_ptr, 02954 Vector<int, Storage2, Allocator2>& real_ind, 02955 Vector<T, Storage0, Allocator0>& imag_values, 02956 Vector<int, Storage1, Allocator1>& imag_ptr, 02957 Vector<int, Storage2, Allocator2>& imag_ind): 02958 Matrix_ComplexSparse<T, Prop, ColComplexSparse, Allocator>(i, j, 02959 real_values, 02960 real_ptr, 02961 real_ind, 02962 imag_values, 02963 imag_ptr, 02964 imag_ind) 02965 { 02966 } 02967 02968 02969 02971 // MATRIX<ROWCOMPLEXSPARSE> // 02973 02974 02975 /**************** 02976 * CONSTRUCTORS * 02977 ****************/ 02978 02980 02983 template <class T, class Prop, class Allocator> 02984 Matrix<T, Prop, RowComplexSparse, Allocator>::Matrix(): 02985 Matrix_ComplexSparse<T, Prop, RowComplexSparse, Allocator>() 02986 { 02987 } 02988 02989 02991 02995 template <class T, class Prop, class Allocator> 02996 Matrix<T, Prop, RowComplexSparse, Allocator>::Matrix(int i, int j): 02997 Matrix_ComplexSparse<T, Prop, RowComplexSparse, Allocator>(i, j, 0, 0) 02998 { 02999 } 03000 03001 03010 template <class T, class Prop, class Allocator> 03011 Matrix<T, Prop, RowComplexSparse, Allocator>::Matrix(int i, int j, 03012 int real_nz, 03013 int imag_nz): 03014 Matrix_ComplexSparse<T, Prop, RowComplexSparse, Allocator>(i, j, 03015 real_nz, 03016 imag_nz) 03017 { 03018 } 03019 03020 03022 03040 template <class T, class Prop, class Allocator> 03041 template <class Storage0, class Allocator0, 03042 class Storage1, class Allocator1, 03043 class Storage2, class Allocator2> 03044 Matrix<T, Prop, RowComplexSparse, Allocator>:: 03045 Matrix(int i, int j, 03046 Vector<T, Storage0, Allocator0>& real_values, 03047 Vector<int, Storage1, Allocator1>& real_ptr, 03048 Vector<int, Storage2, Allocator2>& real_ind, 03049 Vector<T, Storage0, Allocator0>& imag_values, 03050 Vector<int, Storage1, Allocator1>& imag_ptr, 03051 Vector<int, Storage2, Allocator2>& imag_ind): 03052 Matrix_ComplexSparse<T, Prop, RowComplexSparse, Allocator>(i, j, 03053 real_values, 03054 real_ptr, 03055 real_ind, 03056 imag_values, 03057 imag_ptr, 03058 imag_ind) 03059 { 03060 } 03061 03062 03063 } // namespace Seldon. 03064 03065 #define SELDON_FILE_MATRIX_COMPLEXSPARSE_CXX 03066 #endif
