|
AbstractLinAlgPack: C++ Interfaces For Vectors, Matrices And Related Linear Algebra Objects Version of the Day
|
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 #include "AbstractLinAlgPack_SpVectorOp.hpp" 00030 00031 namespace { 00032 // Setup some template classes to check at complile time that 00033 // the layout of SpVectorSlice::element_type is proper. 00034 template<int N, class T> 00035 class assert_compile_time { 00036 public: 00037 assert_compile_time() 00038 { 00039 // This should not compile if instantiated with a type T that 00040 // is not an integer. However, if the compiler checks this 00041 // function without instantiating it, it can not cause an error 00042 // because it does not know the type of T to see if the 00043 // conversion is legal or not. 00044 T d; 00045 static_cast<int*>(d); 00046 } 00047 }; 00048 // Template specialization for no error 00049 template<> 00050 class assert_compile_time<0,double> { 00051 public: 00052 assert_compile_time() 00053 {} 00054 }; 00055 // Validate that there is an integer stride between values 00056 assert_compile_time< 00057 ((int)sizeof(AbstractLinAlgPack::SpVectorSlice::element_type) 00058 % (int)sizeof(DenseLinAlgPack::value_type)) 00059 , double 00060 > 00061 validate_value_stride; 00062 // Validate that there is an integer stride between indexes 00063 assert_compile_time< 00064 ((int)sizeof(AbstractLinAlgPack::SpVectorSlice::element_type) 00065 % (int)sizeof(DenseLinAlgPack::index_type)) 00066 , double 00067 > 00068 validate_index_stride; 00069 } // end namespace 00070 00071 void AbstractLinAlgPack::add_elements( SpVector* sv_lhs, value_type alpha, const DVectorSlice& vs_rhs 00072 , size_type offset, bool add_zeros ) 00073 { 00074 typedef SpVector::element_type ele_t; 00075 const bool assume_sorted = !sv_lhs->nz() || ( sv_lhs->nz() && sv_lhs->is_sorted() ); 00076 DVectorSlice::const_iterator 00077 itr = vs_rhs.begin(); 00078 if(add_zeros) { 00079 for( size_type i = 1; i <= vs_rhs.dim(); ++i ) 00080 sv_lhs->add_element( ele_t( i + offset, alpha * (*itr++) ) ); 00081 } 00082 else { 00083 for( size_type i = 1; i <= vs_rhs.dim(); ++i, ++itr ) 00084 if( *itr != 0.0 ) 00085 sv_lhs->add_element( ele_t( i + offset, alpha * (*itr) ) ); 00086 } 00087 sv_lhs->assume_sorted(assume_sorted); 00088 } 00089 00090 void AbstractLinAlgPack::add_elements( SpVector* sv_lhs, value_type alpha, const SpVectorSlice& sv_rhs 00091 , size_type offset, bool add_zeros ) 00092 { 00093 typedef SpVector::element_type ele_t; 00094 const bool assume_sorted = ( !sv_lhs->nz() || ( sv_lhs->nz() && sv_lhs->is_sorted() ) ) 00095 && ( !sv_rhs.nz() || ( sv_rhs.nz() || sv_rhs.is_sorted() ) ); 00096 if(add_zeros) { 00097 for( SpVectorSlice::const_iterator itr = sv_rhs.begin(); itr != sv_rhs.end(); ++itr ) 00098 sv_lhs->add_element( ele_t( itr->index() + sv_rhs.offset() + offset, alpha * (itr->value()) ) ); 00099 } 00100 else { 00101 for( SpVectorSlice::const_iterator itr = sv_rhs.begin(); itr != sv_rhs.end(); ++itr ) 00102 if(itr->value() != 0.0 ) 00103 sv_lhs->add_element( ele_t( itr->index() + sv_rhs.offset() + offset, alpha * (itr->value()) ) ); 00104 } 00105 sv_lhs->assume_sorted(assume_sorted); 00106 }
1.7.4