AbstractLinAlgPack: C++ Interfaces For Vectors, Matrices And Related Linear Algebra Objects: AbstractLinAlgPack

00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
00005 //                  Copyright (2003) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #ifndef ABSTRACT_LINALG_PACK_MATRIX_WITH_OP_H
00030 #define ABSTRACT_LINALG_PACK_MATRIX_WITH_OP_H
00031 
00032 #include <iosfwd>
00033 
00034 #include "AbstractLinAlgPack_MatrixBase.hpp"
00035 #include "Teuchos_RCP.hpp"
00036 
00037 namespace AbstractLinAlgPack {
00038 
00116 class MatrixOp : public virtual MatrixBase {
00117 public:
00118 
00121 
00123   friend
00124   void Mt_S( MatrixOp* mwo_lhs, value_type alpha );
00126   friend
00127   void Mp_StM(
00128     MatrixOp* mwo_lhs, value_type alpha, const MatrixOp& M_rhs
00129     , BLAS_Cpp::Transp trans_rhs);
00131   friend
00132   void Mp_StMtP(
00133     MatrixOp* mwo_lhs, value_type alpha
00134     ,const MatrixOp& M_rhs, BLAS_Cpp::Transp M_trans
00135     ,const GenPermMatrixSlice& P_rhs, BLAS_Cpp::Transp P_rhs_trans
00136     );
00138   friend
00139   void Mp_StPtM(
00140     MatrixOp* mwo_lhs, value_type alpha
00141     ,const GenPermMatrixSlice& P_rhs, BLAS_Cpp::Transp P_rhs_trans
00142     ,const MatrixOp& M_rhs, BLAS_Cpp::Transp M_trans
00143     );
00145   friend
00146   void Mp_StPtMtP(
00147     MatrixOp* mwo_lhs, value_type alpha
00148     ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00149     ,const MatrixOp& M_rhs, BLAS_Cpp::Transp trans_rhs
00150     ,const GenPermMatrixSlice& P_rhs2, BLAS_Cpp::Transp P_rhs2_trans
00151     );
00153   friend
00154   void Vp_StMtV(
00155     VectorMutable* v_lhs, value_type alpha, const MatrixOp& M_rhs1
00156     ,BLAS_Cpp::Transp trans_rhs1, const Vector& v_rhs2, value_type beta
00157     );
00159   friend
00160   void Vp_StMtV(
00161     VectorMutable* v_lhs, value_type alpha, const MatrixOp& M_rhs1
00162     ,BLAS_Cpp::Transp trans_rhs1, const SpVectorSlice& sv_rhs2, value_type beta
00163     );
00165   friend
00166   void Vp_StPtMtV(
00167     VectorMutable* v_lhs, value_type alpha
00168     ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00169     ,const MatrixOp& M_rhs2, BLAS_Cpp::Transp M_rhs2_trans
00170     ,const Vector& v_rhs3, value_type beta
00171     );
00173   friend
00174   void Vp_StPtMtV(
00175     VectorMutable* v_lhs, value_type alpha
00176     ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00177     ,const MatrixOp& M_rhs2, BLAS_Cpp::Transp M_rhs2_trans
00178     ,const SpVectorSlice& sv_rhs3, value_type beta
00179     );
00181   friend
00182   value_type transVtMtV(
00183     const Vector& v_rhs1, const MatrixOp& M_rhs2
00184     ,BLAS_Cpp::Transp trans_rhs2, const Vector& v_rhs3
00185     );
00187   friend
00188   value_type transVtMtV(
00189     const SpVectorSlice& sv_rhs1, const MatrixOp& M_rhs2
00190     ,BLAS_Cpp::Transp trans_rhs2, const SpVectorSlice& sv_rhs3
00191     );
00193   friend
00194   void syr2k(
00195     const MatrixOp& M, BLAS_Cpp::Transp M_trans, value_type alpha
00196     ,const GenPermMatrixSlice& P1, BLAS_Cpp::Transp P1_trans
00197     ,const GenPermMatrixSlice& P2, BLAS_Cpp::Transp P2_trans
00198     ,value_type beta, MatrixSymOp* symwo_lhs
00199     );
00201   friend
00202   void Mp_StMtM(
00203     MatrixOp* mwo_lhs, value_type alpha
00204     ,const MatrixOp& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
00205     ,const MatrixOp& mwo_rhs2, BLAS_Cpp::Transp trans_rhs2
00206     ,value_type beta
00207     );
00209   friend
00210   void syrk(
00211     const MatrixOp  &mwo_rhs
00212     ,BLAS_Cpp::Transp   M_trans
00213     ,value_type         alpha
00214     ,value_type         beta
00215     ,MatrixSymOp    *sym_lhs
00216     );
00217 
00219 
00222 
00223 #ifndef DOXYGEN_COMPILE
00224 
00225   typedef Teuchos::RCP<const MatrixOp>    mat_ptr_t;
00227   typedef Teuchos::RCP<MatrixOp>          mat_mut_ptr_t;
00228 #endif
00229 
00231   enum EMatNormType {
00232     MAT_NORM_INF        
00233     ,MAT_NORM_2         
00234     ,MAT_NORM_1         
00235     ,MAT_NORM_FORB      
00236   };
00237 
00239   struct MatNorm {
00240     MatNorm(value_type _value, EMatNormType _type) : value(_value), type(_type) {}
00241     value_type    value;
00242     EMatNormType  type;
00243   };
00244 
00246   class MethodNotImplemented : public std::runtime_error
00247   {public: MethodNotImplemented(const std::string& what_arg) : std::runtime_error(what_arg) {}};
00248 
00250   class IncompatibleMatrices : public std::logic_error
00251   {public: IncompatibleMatrices(const std::string& what_arg) : std::logic_error(what_arg) {}};
00252 
00254 
00257 
00267   virtual void zero_out();
00268 
00277   virtual void Mt_S( value_type alpha );
00278 
00286   virtual MatrixOp& operator=(const MatrixOp& mwo_rhs);
00287 
00289 
00292 
00305   virtual mat_mut_ptr_t clone();
00306 
00312   virtual mat_ptr_t clone() const;
00313 
00315 
00318 
00325   virtual std::ostream& output(std::ostream& out) const;
00326 
00328 
00331 
00362   const MatNorm calc_norm(
00363     EMatNormType  requested_norm_type = MAT_NORM_1
00364     ,bool         allow_replacement   = false
00365     ) const;
00366 
00368 
00371 
00387   virtual mat_ptr_t sub_view(const Range1D& row_rng, const Range1D& col_rng) const;
00388   
00391   mat_ptr_t sub_view(
00392     const index_type& rl, const index_type& ru
00393     ,const index_type& cl, const index_type& cu
00394     ) const;
00395 
00397   
00400 
00441   virtual mat_ptr_t perm_view(
00442     const Permutation          *P_row
00443     ,const index_type          row_part[]
00444     ,int                       num_row_part
00445     ,const Permutation         *P_col
00446     ,const index_type          col_part[]
00447     ,int                       num_col_part
00448     ) const;
00449 
00479   virtual mat_ptr_t perm_view_update(
00480     const Permutation          *P_row
00481     ,const index_type          row_part[]
00482     ,int                       num_row_part
00483     ,const Permutation         *P_col
00484     ,const index_type          col_part[]
00485     ,int                       num_col_part
00486     ,const mat_ptr_t           &perm_view
00487     ) const;
00488 
00490 
00491 #ifdef TEMPLATE_FRIENDS_NOT_SUPPORTED
00492 public:
00493 #else
00494 protected:
00495 #endif
00496 
00499 
00507   virtual bool Mp_StM(
00508     MatrixOp* mwo_lhs, value_type alpha
00509     ,BLAS_Cpp::Transp trans_rhs
00510     ) const;
00511 
00519   virtual bool Mp_StM(
00520     value_type alpha,const MatrixOp& M_rhs, BLAS_Cpp::Transp trans_rhs
00521     );
00522 
00530   virtual bool Mp_StMtP(
00531     MatrixOp* mwo_lhs, value_type alpha
00532     ,BLAS_Cpp::Transp M_trans
00533     ,const GenPermMatrixSlice& P_rhs, BLAS_Cpp::Transp P_rhs_trans
00534     ) const;
00535 
00543   virtual bool Mp_StMtP(
00544     value_type alpha
00545     ,const MatrixOp& mwo_rhs, BLAS_Cpp::Transp M_trans
00546     ,const GenPermMatrixSlice& P_rhs, BLAS_Cpp::Transp P_rhs_trans
00547     );
00548 
00556   virtual bool Mp_StPtM(
00557     MatrixOp* mwo_lhs, value_type alpha
00558     ,const GenPermMatrixSlice& P_rhs, BLAS_Cpp::Transp P_rhs_trans
00559     ,BLAS_Cpp::Transp M_trans
00560     ) const;
00561 
00569   virtual bool Mp_StPtM(
00570     value_type alpha
00571     ,const GenPermMatrixSlice& P_rhs, BLAS_Cpp::Transp P_rhs_trans
00572     ,const MatrixOp& mwo_rhs, BLAS_Cpp::Transp M_trans
00573     );
00574 
00582   virtual bool Mp_StPtMtP(
00583     MatrixOp* mwo_lhs, value_type alpha
00584     ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00585     ,BLAS_Cpp::Transp M_trans
00586     ,const GenPermMatrixSlice& P_rhs2, BLAS_Cpp::Transp P_rhs2_trans
00587     ) const;
00588 
00596   virtual bool Mp_StPtMtP(
00597     value_type alpha
00598     ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00599     ,const MatrixOp& mwo_rhs, BLAS_Cpp::Transp M_trans
00600     ,const GenPermMatrixSlice& P_rhs2, BLAS_Cpp::Transp P_rhs2_trans
00601     );
00602 
00603   //    end Level-1 BLAS
00605 
00608 
00610   virtual void Vp_StMtV(
00611     VectorMutable* v_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1
00612     ,const Vector& v_rhs2, value_type beta
00613     ) const = 0;
00614 
00616   virtual void Vp_StMtV(
00617     VectorMutable* v_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1
00618     ,const SpVectorSlice& sv_rhs2, value_type beta
00619     ) const;
00620 
00622   virtual void Vp_StPtMtV(
00623     VectorMutable* v_lhs, value_type alpha
00624     ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00625     ,BLAS_Cpp::Transp M_rhs2_trans
00626     ,const Vector& v_rhs3, value_type beta
00627     ) const;
00628 
00630   virtual void Vp_StPtMtV(
00631     VectorMutable* v_lhs, value_type alpha
00632     ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00633     ,BLAS_Cpp::Transp M_rhs2_trans
00634     ,const SpVectorSlice& sv_rhs3, value_type beta
00635     ) const;
00636 
00638   virtual value_type transVtMtV(
00639     const Vector& v_rhs1, BLAS_Cpp::Transp trans_rhs2
00640     ,const Vector& v_rhs3
00641     ) const;
00642 
00644   virtual value_type transVtMtV(
00645     const SpVectorSlice& sv_rhs1, BLAS_Cpp::Transp trans_rhs2
00646     ,const SpVectorSlice& sv_rhs3
00647     ) const;
00648 
00660   virtual void syr2k(
00661     BLAS_Cpp::Transp M_trans, value_type alpha
00662     ,const GenPermMatrixSlice& P1, BLAS_Cpp::Transp P1_trans
00663     ,const GenPermMatrixSlice& P2, BLAS_Cpp::Transp P2_trans
00664     ,value_type beta, MatrixSymOp* symwo_lhs
00665     ) const;
00666 
00668 
00671 
00680   virtual bool Mp_StMtM(
00681     MatrixOp* mwo_lhs, value_type alpha
00682     ,BLAS_Cpp::Transp trans_rhs1
00683     ,const MatrixOp& mwo_rhs2, BLAS_Cpp::Transp trans_rhs2
00684     ,value_type beta
00685     ) const;
00686 
00695   virtual bool Mp_StMtM(
00696     MatrixOp* mwo_lhs, value_type alpha
00697     ,const MatrixOp& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
00698     ,BLAS_Cpp::Transp trans_rhs2
00699     ,value_type beta
00700     ) const;
00701 
00710   virtual bool Mp_StMtM(
00711     value_type alpha
00712     ,const MatrixOp& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
00713     ,const MatrixOp& mwo_rhs2,BLAS_Cpp::Transp trans_rhs2
00714     ,value_type beta
00715     );
00716   
00728   virtual bool syrk(
00729     BLAS_Cpp::Transp   M_trans
00730     ,value_type        alpha
00731     ,value_type        beta
00732     ,MatrixSymOp      *sym_lhs
00733     ) const;
00734   
00746   virtual bool syrk(
00747     const MatrixOp      &mwo_rhs
00748     ,BLAS_Cpp::Transp   M_trans
00749     ,value_type         alpha
00750     ,value_type         beta
00751     );
00752 
00753   //    end Level-3 BLAS
00755 
00756 };  // end class MatrixOp
00757 
00758 // ////////////////////////////////////////////////////////////////////////////////////////////////
00768 
00771 
00772 //
00779 void Mt_S( MatrixOp* mwo_lhs, value_type alpha );
00780 
00781 //
00796 void Mp_StM(
00797   MatrixOp* mwo_lhs, value_type alpha, const MatrixOp& M_rhs
00798   , BLAS_Cpp::Transp trans_rhs);
00799 
00806 void Mp_StMtP(
00807   MatrixOp* mwo_lhs, value_type alpha
00808   ,const MatrixOp& M_rhs, BLAS_Cpp::Transp M_trans
00809   ,const GenPermMatrixSlice& P_rhs, BLAS_Cpp::Transp P_rhs_trans
00810   );
00811 
00818 void Mp_StPtM(
00819   MatrixOp* mwo_lhs, value_type alpha
00820   ,const GenPermMatrixSlice& P_rhs, BLAS_Cpp::Transp P_rhs_trans
00821   ,const MatrixOp& M_rhs, BLAS_Cpp::Transp M_trans
00822   );
00823 
00830 void Mp_StPtMtP(
00831   MatrixOp* mwo_lhs, value_type alpha
00832   ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00833   ,const MatrixOp& M_rhs, BLAS_Cpp::Transp trans_rhs
00834   ,const GenPermMatrixSlice& P_rhs2, BLAS_Cpp::Transp P_rhs2_trans
00835   );
00836 
00837 //    end Level-1 BLAS
00839 
00842 
00844 inline void Vp_StMtV(
00845   VectorMutable* v_lhs, value_type alpha, const MatrixOp& M_rhs1
00846   ,BLAS_Cpp::Transp trans_rhs1, const Vector& v_rhs2, value_type beta = 1.0
00847   )
00848 {
00849   M_rhs1.Vp_StMtV(v_lhs,alpha,trans_rhs1,v_rhs2,beta);
00850 }
00851 
00853 inline void Vp_StMtV(
00854   VectorMutable* v_lhs, value_type alpha, const MatrixOp& M_rhs1
00855   ,BLAS_Cpp::Transp trans_rhs1, const SpVectorSlice& sv_rhs2, value_type beta = 1.0
00856   )
00857 {
00858   M_rhs1.Vp_StMtV(v_lhs,alpha,trans_rhs1,sv_rhs2,beta);
00859 }
00860 
00862 inline void Vp_StPtMtV(
00863   VectorMutable* v_lhs, value_type alpha
00864   ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00865   ,const MatrixOp& M_rhs2, BLAS_Cpp::Transp M_rhs2_trans
00866   ,const Vector& v_rhs3, value_type beta = 1.0
00867   ) 
00868 {
00869   M_rhs2.Vp_StPtMtV(v_lhs,alpha,P_rhs1,P_rhs1_trans,M_rhs2_trans,v_rhs3,beta);
00870 }
00871 
00873 inline void Vp_StPtMtV(
00874   VectorMutable* v_lhs, value_type alpha
00875   ,const GenPermMatrixSlice& P_rhs1, BLAS_Cpp::Transp P_rhs1_trans
00876   ,const MatrixOp& M_rhs2, BLAS_Cpp::Transp M_rhs2_trans
00877   ,const SpVectorSlice& sv_rhs3, value_type beta = 1.0
00878   )
00879 {
00880   M_rhs2.Vp_StPtMtV(v_lhs,alpha,P_rhs1,P_rhs1_trans,M_rhs2_trans,sv_rhs3,beta);
00881 }
00882 
00884 inline value_type transVtMtV(
00885   const Vector& v_rhs1, const MatrixOp& M_rhs2
00886   ,BLAS_Cpp::Transp trans_rhs2, const Vector& v_rhs3
00887   )
00888 {
00889   return M_rhs2.transVtMtV(v_rhs1,trans_rhs2,v_rhs3);
00890 }
00891 
00893 inline value_type transVtMtV(
00894   const SpVectorSlice& sv_rhs1, const MatrixOp& M_rhs2
00895   ,BLAS_Cpp::Transp trans_rhs2, const SpVectorSlice& sv_rhs3
00896   )
00897 {
00898   return M_rhs2.transVtMtV(sv_rhs1,trans_rhs2,sv_rhs3);
00899 }
00900 
00902 inline void syr2k(
00903   const MatrixOp& M, BLAS_Cpp::Transp M_trans, value_type alpha
00904   ,const GenPermMatrixSlice& P1, BLAS_Cpp::Transp P1_trans
00905   ,const GenPermMatrixSlice& P2, BLAS_Cpp::Transp P2_trans
00906   ,value_type beta, MatrixSymOp* symwo_lhs
00907   )
00908 {
00909   M.syr2k(M_trans,alpha,P1,P1_trans,P2,P2_trans,beta,symwo_lhs);
00910 }
00911 
00912 //    end Level-2 BLAS
00914 
00917 
00930 void Mp_StMtM(
00931   MatrixOp* mwo_lhs, value_type alpha
00932   ,const MatrixOp& mwo_rhs1, BLAS_Cpp::Transp trans_rhs1
00933   ,const MatrixOp& mwo_rhs2, BLAS_Cpp::Transp trans_rhs2
00934   ,value_type beta = 1.0
00935   );
00936 
00943 void syrk(
00944   const MatrixOp  &mwo_rhs
00945   ,BLAS_Cpp::Transp   M_trans
00946   ,value_type         alpha
00947   ,value_type         beta
00948   ,MatrixSymOp    *sym_lhs
00949   );
00950 
00951 //    end Level-3 BLAS
00953 
00954 //    end non-member functions
00956 
00957 // //////////////////////////////////////////////////
00958 // Inlined methods for MatrixOp
00959 
00960 inline
00961 MatrixOp::mat_ptr_t
00962 MatrixOp::sub_view(
00963   const index_type& rl, const index_type& ru
00964   ,const index_type& cl, const index_type& cu
00965   ) const
00966 {
00967   return this->sub_view(Range1D(rl,ru),Range1D(cl,cu));
00968 }
00969 
00970 } // end namespace AbstractLinAlgPack
00971 
00972 #endif  // ABSTRACT_LINALG_PACK_MATRIX_WITH_OP_H