omatrix_cx_tlate.h

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

Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.