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

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).

Get Object-oriented Scientific Computing
Lib at SourceForge.net. Fast, secure and Free Open Source software
downloads.