00001 /* 00002 ------------------------------------------------------------------- 00003 00004 Copyright (C) 2006, 2007, 2008, Andrew W. Steiner 00005 00006 This file is part of O2scl. 00007 00008 O2scl is free software; you can redistribute it and/or modify 00009 it under the terms of the GNU General Public License as published by 00010 the Free Software Foundation; either version 3 of the License, or 00011 (at your option) any later version. 00012 00013 O2scl is distributed in the hope that it will be useful, 00014 but WITHOUT ANY WARRANTY; without even the implied warranty of 00015 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00016 GNU General Public License for more details. 00017 00018 You should have received a copy of the GNU General Public License 00019 along with O2scl. If not, see <http://www.gnu.org/licenses/>. 00020 00021 ------------------------------------------------------------------- 00022 */ 00023 #ifndef O2SCL_OMATRIX_CX_TLATE_H 00024 #define O2SCL_OMATRIX_CX_TLATE_H 00025 00026 /** \file omatrix_cx_tlate.h 00027 \brief File for definitions of complex matrices 00028 */ 00029 00030 #include <iostream> 00031 #include <cstdlib> 00032 #include <string> 00033 #include <fstream> 00034 #include <sstream> 00035 #include <o2scl/err_hnd.h> 00036 #include <gsl/gsl_matrix.h> 00037 #include <gsl/gsl_complex.h> 00038 #include <o2scl/ovector_tlate.h> 00039 #include <o2scl/ovector_cx_tlate.h> 00040 00041 #ifndef DOXYGENP 00042 namespace o2scl { 00043 #endif 00044 00045 /** 00046 \brief A matrix view of double-precision numbers 00047 */ 00048 template<class data_t, class parent_t, class block_t, class complex_t> 00049 class omatrix_cx_view_tlate : 00050 public parent_t { 00051 00052 public: 00053 00054 /// \name Copy constructors 00055 //@{ 00056 /// So2scllow copy constructor - create a new view of the same matrix 00057 omatrix_cx_view_tlate(const omatrix_cx_view_tlate &v) { 00058 parent_t::block=NULL; 00059 parent_t::data=v.data; 00060 parent_t::size1=v.size1; 00061 parent_t::size2=v.size2; 00062 parent_t::tda=v.tda; 00063 parent_t::owner=0; 00064 } 00065 00066 /// So2scllow copy constructor - create a new view of the same matrix 00067 omatrix_cx_view_tlate& operator=(const omatrix_cx_view_tlate &v) { 00068 parent_t::block=NULL; 00069 parent_t::data=v.data; 00070 parent_t::size1=v.size1; 00071 parent_t::size2=v.size2; 00072 parent_t::tda=v.tda; 00073 parent_t::owner=0; 00074 00075 return *this; 00076 } 00077 //@} 00078 00079 ~omatrix_cx_view_tlate() {}; 00080 00081 /// \name Get and set methods 00082 //@{ 00083 /** 00084 \brief Array-like indexing 00085 */ 00086 complex_t *operator[](size_t i) { 00087 #if O2SCL_NO_RANGE_CHECK 00088 #else 00089 if (i>=parent_t::size1) { 00090 set_err((((std::string)"Array index ")+itos(i)+" out of bounds" 00091 +" in omatrix_cx_view_tlate::operator[]. Size: "+ 00092 itos(parent_t::size1)+ 00093 " (index should be less than size).").c_str(),gsl_index); 00094 return (complex_t *)(parent_t::data); 00095 } 00096 #endif 00097 return (complex_t *)(parent_t::data+i*parent_t::tda*2); 00098 } 00099 00100 /** 00101 \brief Array-like indexing 00102 */ 00103 const complex_t *operator[](size_t i) const { 00104 #if O2SCL_NO_RANGE_CHECK 00105 #else 00106 if (i>=parent_t::size1) { 00107 set_err((((std::string)"Array index ")+itos(i)+" out of bounds" 00108 +" in omatrix_cx_view_tlate::operator[] const. Size: "+ 00109 itos(parent_t::size1)+ 00110 " (index should be less than size).").c_str(),gsl_index); 00111 return (const complex_t *)(parent_t::data); 00112 } 00113 #endif 00114 return (const complex_t *)(parent_t::data+i*parent_t::tda*2); 00115 } 00116 00117 /** 00118 \brief Array-like indexing 00119 */ 00120 complex_t &operator()(size_t i, size_t j) { 00121 #if O2SCL_NO_RANGE_CHECK 00122 #else 00123 if (i>=parent_t::size1 || j>=parent_t::size2) { 00124 set_err((((std::string)"Array indices ")+itos(i)+ 00125 ", "+itos(j)+" out of bounds" 00126 +" in omatrix_cx_view_tlate::operator(). Sizes: "+ 00127 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00128 " (indices should be less than sizes).").c_str(), 00129 gsl_index); 00130 return *(complex_t *)(parent_t::data); 00131 } 00132 #endif 00133 return *(complex_t *)(parent_t::data+(i*parent_t::tda+j)*2); 00134 } 00135 00136 /** 00137 \brief Array-like indexing 00138 */ 00139 const complex_t &operator()(size_t i, size_t j) const { 00140 #if O2SCL_NO_RANGE_CHECK 00141 #else 00142 if (i>=parent_t::size1 || j>=parent_t::size2) { 00143 set_err((((std::string)"Array indices ")+itos(i)+ 00144 ", "+itos(j)+" out of bounds" 00145 +" in omatrix_cx_view_tlate::operator() const. Sizes: "+ 00146 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00147 " (indices should be less than sizes).").c_str(), 00148 gsl_index); 00149 return *(const complex_t *)(parent_t::data); 00150 } 00151 #endif 00152 return *(const complex_t *)(parent_t::data+(i*parent_t::tda+j)*2); 00153 } 00154 00155 /** \brief Get (with optional range-checking) */ 00156 complex_t get(size_t i, size_t j) const { 00157 #if O2SCL_NO_RANGE_CHECK 00158 #else 00159 if (i>=parent_t::size1 || j>=parent_t::size2) { 00160 set_err((((std::string)"Array indices ")+itos(i)+ 00161 ", "+itos(j)+" out of bounds" 00162 +" in omatrix_cx_view_tlate::get(). Sizes: "+ 00163 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00164 " (indices should be less than sizes).").c_str(), 00165 gsl_index); 00166 complex_t zero={{0,0}}; 00167 return zero; 00168 } 00169 #endif 00170 return *(complex_t *)(parent_t::data+(i*parent_t::tda+j)*2); 00171 } 00172 00173 /** \brief Get STL-like complex number (with optional range-checking) */ 00174 std::complex<data_t> get_stl(size_t i, size_t j) const { 00175 #if O2SCL_NO_RANGE_CHECK 00176 #else 00177 if (i>=parent_t::size1 || j>=parent_t::size2) { 00178 set_err((((std::string)"Array indices ")+itos(i)+ 00179 ", "+itos(j)+" out of bounds" 00180 +" in omatrix_cx_view_tlate::get_stl(). Sizes: "+ 00181 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00182 " (indices should be less than sizes).").c_str(), 00183 gsl_index); 00184 std::complex<data_t> zero=0; 00185 return zero; 00186 } 00187 #endif 00188 return *(std::complex<data_t> *)(parent_t::data+(i*parent_t::tda+j)*2); 00189 } 00190 00191 /** \brief Get real part (with optional range-checking) */ 00192 data_t real(size_t i, size_t j) const { 00193 #if O2SCL_NO_RANGE_CHECK 00194 #else 00195 if (i>=parent_t::size1 || j>=parent_t::size2) { 00196 set_err((((std::string)"Array indices ")+itos(i)+ 00197 ", "+itos(j)+" out of bounds" 00198 +" in omatrix_cx_view_tlate::real(). Sizes: "+ 00199 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00200 " (indices should be less than sizes).").c_str(), 00201 gsl_index); 00202 return 0; 00203 } 00204 #endif 00205 return *(parent_t::data+(i*parent_t::tda+j)*2); 00206 } 00207 00208 /** \brief Get imaginary part(with optional range-checking) */ 00209 data_t imag(size_t i, size_t j) const { 00210 #if O2SCL_NO_RANGE_CHECK 00211 #else 00212 if (i>=parent_t::size1 || j>=parent_t::size2) { 00213 set_err((((std::string)"Array indices ")+itos(i)+ 00214 ", "+itos(j)+" out of bounds" 00215 +" in omatrix_cx_view_tlate::imag(). Sizes: "+ 00216 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00217 " (indices should be less than sizes).").c_str(), 00218 gsl_index); 00219 return 0; 00220 } 00221 #endif 00222 return *(parent_t::data+(i*parent_t::tda+j)*2+1); 00223 } 00224 00225 /** \brief Get pointer (with optional range-checking) */ 00226 complex_t *get_ptr(size_t i, size_t j) { 00227 #if O2SCL_NO_RANGE_CHECK 00228 #else 00229 if (i>=parent_t::size1 || j>=parent_t::size2) { 00230 set_err((((std::string)"Array indices ")+itos(i)+ 00231 ", "+itos(j)+" out of bounds" 00232 +" in omatrix_cx_view_tlate::get_ptr(). Sizes: "+ 00233 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00234 " (indices should be less than sizes).").c_str(), 00235 gsl_index); 00236 return (complex_t *)(parent_t::data); 00237 } 00238 #endif 00239 return (complex_t *)(parent_t::data+(i*parent_t::tda+j)*2); 00240 } 00241 00242 /** \brief Get pointer (with optional range-checking) */ 00243 const complex_t *get_const_ptr(size_t i, size_t j) const { 00244 #if O2SCL_NO_RANGE_CHECK 00245 #else 00246 if (i>=parent_t::size1 || j>=parent_t::size2) { 00247 set_err((((std::string)"Array indices ")+itos(i)+ 00248 ", "+itos(j)+" out of bounds" 00249 +" in omatrix_cx_view_tlate::get_const_ptr(). Sizes: "+ 00250 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00251 " (indices should be less than sizes).").c_str(), 00252 gsl_index); 00253 return (const complex_t *)(parent_t::data); 00254 } 00255 #endif 00256 return (const complex_t *)(parent_t::data+(i*parent_t::tda+j)*2); 00257 } 00258 00259 /** \brief Set (with optional range-checking) */ 00260 int set(size_t i, size_t j, complex_t &val) { 00261 #if O2SCL_NO_RANGE_CHECK 00262 #else 00263 if (i>=parent_t::size1 || j>=parent_t::size2) { 00264 set_err_ret((((std::string)"Array indices ")+itos(i)+ 00265 ", "+itos(j)+" out of bounds" 00266 +" in omatrix_cx_view_tlate::set(). Sizes: "+ 00267 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00268 " (indices should be less than sizes).").c_str(), 00269 gsl_index); 00270 } 00271 #endif 00272 parent_t::data[(i*parent_t::tda+j)*2]=GSL_REAL(val); 00273 parent_t::data[(i*parent_t::tda+j)*2+1]=GSL_IMAG(val); 00274 return 0; 00275 } 00276 00277 /** \brief Set (with optional range-checking) */ 00278 int set(size_t i, size_t j, data_t vr, data_t vi) { 00279 #if O2SCL_NO_RANGE_CHECK 00280 #else 00281 if (i>=parent_t::size1 || j>=parent_t::size2) { 00282 set_err_ret((((std::string)"Array indices ")+itos(i)+ 00283 ", "+itos(j)+" out of bounds" 00284 +" in omatrix_cx_view_tlate::set(). Sizes: "+ 00285 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00286 " (indices should be less than sizes).").c_str(), 00287 gsl_index); 00288 } 00289 #endif 00290 parent_t::data[(i*parent_t::tda+j)*2]=vr; 00291 parent_t::data[(i*parent_t::tda+j)*2+1]=vi; 00292 return 0; 00293 } 00294 00295 /** \brief Set (with optional range-checking) */ 00296 int set_real(size_t i, size_t j, data_t vr) { 00297 #if O2SCL_NO_RANGE_CHECK 00298 #else 00299 if (i>=parent_t::size1 || j>=parent_t::size2) { 00300 set_err_ret((((std::string)"Array indices ")+itos(i)+ 00301 ", "+itos(j)+" out of bounds" 00302 +" in omatrix_cx_view_tlate::set_real(). Sizes: "+ 00303 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00304 " (indices should be less than sizes).").c_str(), 00305 gsl_index); 00306 } 00307 #endif 00308 parent_t::data[(i*parent_t::tda+j)*2]=vr; 00309 return 0; 00310 } 00311 00312 /** \brief Set (with optional range-checking) */ 00313 int set_imag(size_t i, size_t j, data_t vi) { 00314 #if O2SCL_NO_RANGE_CHECK 00315 #else 00316 if (i>=parent_t::size1 || j>=parent_t::size2) { 00317 set_err_ret((((std::string)"Array indices ")+itos(i)+ 00318 ", "+itos(j)+" out of bounds" 00319 +" in omatrix_cx_view_tlate::set_imag(). Sizes: "+ 00320 itos(parent_t::size1)+","+itos(parent_t::size2)+ 00321 " (indices should be less than sizes).").c_str(), 00322 gsl_index); 00323 } 00324 #endif 00325 parent_t::data[(i*parent_t::tda+j)*2+1]=vi; 00326 return 0; 00327 } 00328 00329 /** \brief Set all*/ 00330 int set_all(complex_t &val) { 00331 for(size_t i=0;i<parent_t::size1;i++) { 00332 for(size_t j=0;j<parent_t::size2;j++) { 00333 parent_t::data[(i*parent_t::tda+j)*2]=GSL_REAL(val); 00334 parent_t::data[(i*parent_t::tda+j)*2+1]=GSL_IMAG(val); 00335 } 00336 } 00337 return 0; 00338 } 00339 00340 /** 00341 \brief Method to return number of rows 00342 00343 If no memory has been allocated, this will quietly 00344 return zero. 00345 */ 00346 size_t rows() const { 00347 return parent_t::size1; 00348 } 00349 00350 /** 00351 \brief Method to return number of columns 00352 00353 If no memory has been allocated, this will quietly 00354 return zero. 00355 */ 00356 size_t cols() const { 00357 return parent_t::size2; 00358 } 00359 00360 /** 00361 \brief Method to return matrix tda 00362 00363 If no memory has been allocated, this will quietly 00364 return zero. 00365 */ 00366 size_t tda() const { 00367 return parent_t::tda; 00368 } 00369 //@} 00370 00371 /// \name Other methods 00372 //@{ 00373 /// Return true if this object owns the data it refers to 00374 bool is_owner() const { 00375 if (parent_t::owner==1) return true; 00376 return false; 00377 } 00378 //@} 00379 00380 #ifndef DOXYGENP 00381 00382 protected: 00383 00384 /** \brief Empty constructor provided for use by 00385 omatrix_cx_tlate(const omatrix_cx_tlate &v) 00386 */ 00387 omatrix_cx_view_tlate() {}; 00388 00389 #endif 00390 00391 }; 00392 00393 /** 00394 \brief A matrix of double-precision numbers 00395 00396 */ 00397 template<class data_t, class parent_t, class block_t, class complex_t> 00398 class omatrix_cx_tlate : 00399 public omatrix_cx_view_tlate<data_t,parent_t,block_t,complex_t> { 00400 public: 00401 00402 /// \name Standard constructor 00403 //@{ 00404 /** \brief Create an omatrix of size \c n with owner as 'true' */ 00405 omatrix_cx_tlate(size_t r=0, size_t c=0) { 00406 00407 parent_t::data=0; 00408 00409 // This must be set to 1 even if n=0 so that future 00410 // calls to operator= work properly 00411 parent_t::owner=1; 00412 00413 if (r>0 && c>0) { 00414 parent_t::block=(block_t *)malloc(sizeof(block_t)); 00415 if (parent_t::block) { 00416 parent_t::block->data=(data_t *)malloc(2*r*c*sizeof(data_t)); 00417 if (parent_t::block->data) { 00418 parent_t::block->size=r*c; 00419 parent_t::data=parent_t::block->data; 00420 parent_t::size1=r; 00421 parent_t::size2=c; 00422 parent_t::tda=c; 00423 } else { 00424 std::free(parent_t::block); 00425 set_err("No memory for data in omatrix_cx_tlate constructor", 00426 gsl_enomem); 00427 } 00428 } else { 00429 set_err("No memory for block in omatrix_cx_tlate contructor", 00430 gsl_enomem); 00431 } 00432 } else { 00433 parent_t::size1=0; 00434 parent_t::size2=0; 00435 parent_t::tda=0; 00436 } 00437 } 00438 //@} 00439 00440 /// \name Copy constructors 00441 //@{ 00442 /// Deep copy constructor, allocate new space and make a copy 00443 omatrix_cx_tlate(const omatrix_cx_tlate &v) : 00444 omatrix_cx_view_tlate<data_t,parent_t,block_t,complex_t>() { 00445 size_t n=v.size1; 00446 size_t n2=v.size2; 00447 if (n>0 && n2>0) { 00448 parent_t::block=(block_t *)malloc(sizeof(block_t)); 00449 if (parent_t::block) { 00450 parent_t::block->data=(data_t *)malloc(2*n*n2*sizeof(data_t)); 00451 if (parent_t::block->data) { 00452 parent_t::block->size=n*n2; 00453 parent_t::data=parent_t::block->data; 00454 parent_t::size1=n; 00455 parent_t::size2=n2; 00456 parent_t::tda=n2; 00457 parent_t::owner=1; 00458 for(size_t i=0;i<n;i++) { 00459 for(size_t j=0;j<n2;j++) { 00460 *(parent_t::data+i*parent_t::tda+j)=*(v.data+i*v.tda+j); 00461 *(parent_t::data+i*parent_t::tda+j+1)=*(v.data+i*v.tda+j+1); 00462 } 00463 } 00464 } else { 00465 std::free(parent_t::block); 00466 set_err("No memory for data in omatrix_cx_tlate constructor", 00467 gsl_enomem); 00468 } 00469 } else { 00470 set_err("No memory for block in omatrix_cx_tlate contructor", 00471 gsl_enomem); 00472 } 00473 } else { 00474 parent_t::size1=0; 00475 parent_t::size2=0; 00476 parent_t::tda=0; 00477 } 00478 } 00479 00480 /// Deep copy constructor, allocate new space and make a copy 00481 omatrix_cx_tlate 00482 (const omatrix_cx_view_tlate<data_t,parent_t,block_t,complex_t> &v) : 00483 omatrix_cx_view_tlate<data_t,parent_t,block_t,complex_t>() { 00484 size_t r=v.rows(); 00485 size_t c=v.cols(); 00486 if (r>0 && c>0) { 00487 parent_t::block=(block_t *)malloc(sizeof(block_t)); 00488 if (parent_t::block) { 00489 parent_t::block->data=(data_t *)malloc(2*r*c*sizeof(data_t)); 00490 if (parent_t::block->data) { 00491 parent_t::block->size=r*c; 00492 parent_t::data=parent_t::block->data; 00493 parent_t::size1=v.size1; 00494 parent_t::size2=v.size2; 00495 parent_t::tda=v.tda; 00496 parent_t::owner=1; 00497 for(size_t i=0;i<r;i++) { 00498 for(size_t j=0;i<c;j++) { 00499 *(parent_t::data+i*parent_t::tda+j)=*(v.data+i*v.tda+j); 00500 *(parent_t::data+i*parent_t::tda+j+1)= 00501 *(v.data+i*v.tda+j+1); 00502 } 00503 } 00504 } else { 00505 std::free(parent_t::block); 00506 set_err("No memory for data in omatrix_cx_tlate constructor", 00507 gsl_enomem); 00508 } 00509 } else { 00510 set_err("No memory for block in omatrix_cx_tlate contructor", 00511 gsl_enomem); 00512 } 00513 } else { 00514 parent_t::size1=0; 00515 parent_t::size2=0; 00516 parent_t::tda=0; 00517 } 00518 } 00519 //@} 00520 00521 ~omatrix_cx_tlate() { 00522 if (parent_t::size1>0) { 00523 if (parent_t::owner==1) { 00524 if (parent_t::block->size>0) { 00525 std::free(parent_t::block->data); 00526 } 00527 std::free(parent_t::block); 00528 parent_t::size1=0; 00529 parent_t::size2=0; 00530 } 00531 } 00532 } 00533 00534 /// \name Memory allocation 00535 //@{ 00536 /** 00537 \brief Allocate memory for size \c n after freeing any memory 00538 presently in use 00539 */ 00540 int allocate(size_t nrows, size_t ncols) { 00541 if (parent_t::size1>0 || parent_t::size2>0) free(); 00542 00543 if (nrows>0 && ncols>0) { 00544 parent_t::block=(block_t *)malloc(sizeof(block_t)); 00545 if (parent_t::block) { 00546 parent_t::block->data=(data_t *) 00547 malloc(2*nrows*ncols*sizeof(data_t)); 00548 if (parent_t::block->data) { 00549 parent_t::block->size=nrows*ncols; 00550 parent_t::data=parent_t::block->data; 00551 parent_t::size1=nrows; 00552 parent_t::size2=ncols; 00553 parent_t::tda=ncols; 00554 parent_t::owner=1; 00555 } else { 00556 std::free(parent_t::block); 00557 set_err_ret 00558 ("No memory for data in omatrix_cx_tlate::allocate()", 00559 gsl_enomem); 00560 } 00561 } else { 00562 set_err_ret 00563 ("No memory for block in omatrix_cx_tlate::allocate()", 00564 gsl_enomem); 00565 } 00566 } else { 00567 set_err_ret("Zero size in omatrix::allocate()",gsl_einval); 00568 } 00569 return 0; 00570 } 00571 00572 /** 00573 \brief Free the memory 00574 00575 This function will safely do nothing if used without first 00576 allocating memory or if called multiple times in succession. 00577 */ 00578 int free() { 00579 if (parent_t::size1>0) { 00580 if (parent_t::owner==1) { 00581 if (parent_t::block->size>0) { 00582 std::free(parent_t::block->data); 00583 } 00584 std::free(parent_t::block); 00585 } 00586 parent_t::size1=0; 00587 parent_t::size2=0; 00588 parent_t::tda=0; 00589 } 00590 return 0; 00591 } 00592 //@} 00593 00594 /// \name Other methods 00595 //@{ 00596 /// \brief Compute the transpose 00597 omatrix_cx_tlate<data_t,parent_t,block_t,complex_t> transpose() { 00598 omatrix_cx_tlate<data_t,parent_t,block_t,complex_t> 00599 result(parent_t::size2,parent_t::size1); 00600 for(size_t i=0;i<parent_t::size1;i++) { 00601 for(size_t j=0;j<parent_t::size2;j++) { 00602 result[j][i]=(*this)[i][j]; 00603 } 00604 } 00605 } 00606 00607 /// Compute the conjugate transpose 00608 omatrix_cx_tlate<data_t,parent_t,block_t,complex_t> htranspose() { 00609 omatrix_cx_tlate<data_t,parent_t,block_t,complex_t> 00610 result(parent_t::size2,parent_t::size1); 00611 for(size_t i=0;i<parent_t::size1;i++) { 00612 for(size_t j=0;j<parent_t::size2;j++) { 00613 result[j][i]=gsl_complex_conjugate((*this)[i][j]); 00614 } 00615 } 00616 } 00617 //@} 00618 00619 }; 00620 00621 /** \brief Create a vector from a row of a matrix 00622 */ 00623 template<class data_t, class mparent_t, class vparent_t, class block_t, 00624 class complex_t> class omatrix_cx_row_tlate : 00625 public ovector_cx_view_tlate<data_t,vparent_t,block_t,complex_t> { 00626 public: 00627 /** \brief Create a vector from row \c i of matrix \c m */ 00628 omatrix_cx_row_tlate 00629 (omatrix_cx_view_tlate<data_t,mparent_t,block_t,complex_t> &m, 00630 size_t i) { 00631 if (i<m.size1) { 00632 vparent_t::size=m.size2; 00633 vparent_t::stride=1; 00634 vparent_t::data=m.data+m.tda()*i; 00635 vparent_t::owner=0; 00636 vparent_t::block=NULL; 00637 } 00638 } 00639 }; 00640 00641 /** \brief Create a vector from a row of a matrix 00642 */ 00643 template<class data_t, class mparent_t, class vparent_t, class block_t, 00644 class complex_t> class omatrix_cx_const_row_tlate : 00645 public ovector_cx_view_tlate<data_t,vparent_t,block_t,complex_t> { 00646 public: 00647 /** \brief Create a vector from row \c i of matrix \c m */ 00648 omatrix_cx_const_row_tlate 00649 (const omatrix_cx_view_tlate<data_t,mparent_t,block_t,complex_t> &m, 00650 size_t i) { 00651 if (i<m.size1) { 00652 vparent_t::size=m.size2; 00653 vparent_t::stride=1; 00654 vparent_t::data=m.data+m.tda*i; 00655 vparent_t::owner=0; 00656 vparent_t::block=NULL; 00657 } 00658 } 00659 }; 00660 00661 /** \brief Create a vector from a column of a matrix 00662 */ 00663 template<class data_t, class mparent_t, class vparent_t, class block_t, 00664 class complex_t> class omatrix_cx_col_tlate : 00665 public ovector_cx_view_tlate<data_t,vparent_t,block_t,complex_t> { 00666 public: 00667 /** \brief Create a vector from col \c i of matrix \c m */ 00668 omatrix_cx_col_tlate 00669 (omatrix_cx_view_tlate<data_t,mparent_t,block_t,complex_t> &m, 00670 size_t i) { 00671 if (i<m.size2) { 00672 vparent_t::size=m.size1; 00673 vparent_t::stride=m.tda(); 00674 vparent_t::data=m.data+2*i; 00675 vparent_t::owner=0; 00676 vparent_t::block=NULL; 00677 } 00678 } 00679 }; 00680 00681 /** \brief Create a vector from a column of a matrix 00682 */ 00683 template<class data_t, class mparent_t, class vparent_t, class block_t, 00684 class complex_t> class omatrix_cx_const_col_tlate : 00685 public ovector_cx_view_tlate<data_t,vparent_t,block_t,complex_t> { 00686 public: 00687 /** \brief Create a vector from col \c i of matrix \c m */ 00688 omatrix_cx_const_col_tlate 00689 (omatrix_cx_view_tlate<data_t,mparent_t,block_t,complex_t> &m, 00690 size_t i) { 00691 if (i<m.size2) { 00692 vparent_t::size=m.size1; 00693 vparent_t::stride=m.tda(); 00694 vparent_t::data=m.data+2*i; 00695 vparent_t::owner=0; 00696 vparent_t::block=NULL; 00697 } 00698 } 00699 }; 00700 00701 /// omatrix_cx typedef 00702 typedef omatrix_cx_tlate<double,gsl_matrix_complex,gsl_block_complex, 00703 gsl_complex> omatrix_cx; 00704 /// omatrix_cx_view typedef 00705 typedef omatrix_cx_view_tlate<double,gsl_matrix_complex,gsl_block_complex, 00706 gsl_complex> omatrix_cx_view; 00707 /// omatrix_cx_row typedef 00708 typedef omatrix_cx_row_tlate<double,gsl_matrix_complex,gsl_vector_complex, 00709 gsl_block_complex,gsl_complex> omatrix_cx_row; 00710 /// omatrix_cx_col typedef 00711 typedef omatrix_cx_col_tlate<double,gsl_matrix_complex,gsl_vector_complex, 00712 gsl_block_complex,gsl_complex> omatrix_cx_col; 00713 /// omatrix_cx_const_row typedef 00714 typedef omatrix_cx_const_row_tlate<double,gsl_matrix_complex, 00715 gsl_vector_complex,gsl_block_complex,gsl_complex> omatrix_cx_const_row; 00716 /// omatrix_cx_const_col typedef 00717 typedef omatrix_cx_const_col_tlate<double,gsl_matrix_complex, 00718 gsl_vector_complex,gsl_block_complex,gsl_complex> omatrix_cx_const_col; 00719 00720 #ifndef DOXYGENP 00721 } 00722 #endif 00723 00724 #endif 00725
Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.
Project hosting provided by
,
O2scl Sourceforge Project Page