|
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 <assert.h> 00030 00031 #include <stdexcept> 00032 00033 #include "AbstractLinAlgPack_VectorMutableSubView.hpp" 00034 #include "Teuchos_TestForException.hpp" 00035 #include "Teuchos_Workspace.hpp" 00036 #include "Teuchos_dyn_cast.hpp" 00037 00038 namespace AbstractLinAlgPack { 00039 00040 VectorSubView::VectorSubView( const vec_ptr_t& vec, const Range1D& rng ) 00041 { 00042 initialize(vec,rng); 00043 } 00044 00045 void VectorSubView::initialize( const vec_ptr_t& vec, const Range1D& rng ) 00046 { 00047 namespace rcp = MemMngPack; 00048 #ifdef TEUCHOS_DEBUG 00049 TEST_FOR_EXCEPTION( 00050 vec.get() == NULL, std::invalid_argument 00051 ,"VectorSubView::initialize(...) : Error!" ); 00052 #endif 00053 typedef VectorSpace::space_ptr_t space_ptr_t; 00054 space_.initialize( Teuchos::rcp(&vec->space(),false), rng ); 00055 full_vec_ = vec; 00056 this->has_changed(); 00057 } 00058 00059 void VectorSubView::set_uninitialized() 00060 { 00061 space_.set_uninitialized(); 00062 full_vec_ = Teuchos::null; 00063 this->has_changed(); 00064 } 00065 00066 // Overridden from Vector 00067 00068 const VectorSpace& VectorSubView::space() const 00069 { 00070 return space_; 00071 } 00072 00073 index_type VectorSubView::dim() const 00074 { 00075 return space_.dim(); 00076 } 00077 00078 void VectorSubView::apply_op( 00079 const RTOpPack::RTOp& op 00080 ,const size_t num_vecs, const Vector* vecs[] 00081 ,const size_t num_targ_vecs, VectorMutable* targ_vecs[] 00082 ,RTOpPack::ReductTarget* reduct_obj 00083 ,const index_type first_ele_in, const index_type sub_dim_in, const index_type global_offset_in 00084 ) const 00085 { 00086 using Teuchos::dyn_cast; 00087 using Teuchos::Workspace; 00088 Teuchos::WorkspaceStore* wss = Teuchos::get_default_workspace_store().get(); 00089 // If these are in-core vectors then just let a default implementation 00090 // take care of this. 00091 if( this->space().is_in_core() ) { 00092 this->apply_op_serial( 00093 op,num_vecs,vecs,num_targ_vecs,targ_vecs,reduct_obj 00094 ,first_ele_in,sub_dim_in,global_offset_in 00095 ); 00096 return; 00097 } 00098 // These are not in-core vectors so ... 00099 int k; 00100 const index_type this_dim = this->dim(); 00101 #ifdef TEUCHOS_DEBUG 00102 TEST_FOR_EXCEPTION( 00103 sub_dim_in < 0 00104 || !(1 <= first_ele_in && first_ele_in <= this_dim) 00105 || ( sub_dim_in > 0 && (sub_dim_in - (first_ele_in - 1) > this_dim) ) 00106 , std::logic_error 00107 ,"VectorSubView::apply_op(...): Error, first_ele_in = " 00108 << first_ele_in << ", global_offset_in = " << global_offset_in 00109 << ", sub_dim_in = " << sub_dim_in << " and this->dim() = this_dim = " 00110 << this_dim << " are not compatible." ); 00111 #endif 00112 const index_type this_offset = space_impl().rng().lbound() - 1; 00113 const index_type 00114 this_sub_dim = ( sub_dim_in 00115 ? sub_dim_in 00116 : space_impl().rng().size() - (first_ele_in - 1) 00117 ); 00118 Workspace<const Vector*> vecs_full(wss,num_vecs); 00119 for( k = 0; k < num_vecs; ++k ) 00120 vecs_full[k] = dyn_cast<const VectorSubView>(*vecs[k]).full_vec().get(); 00121 Workspace<VectorMutable*> targ_vecs_full(wss,num_targ_vecs); 00122 for( k = 0; k < num_targ_vecs; ++k ) 00123 targ_vecs_full[k] = dyn_cast<VectorMutableSubView>(*targ_vecs[k]).full_vec().get(); 00124 AbstractLinAlgPack::apply_op( 00125 op 00126 ,num_vecs, num_vecs ? &vecs_full[0] : NULL 00127 ,num_targ_vecs, num_targ_vecs ? &targ_vecs_full[0] : NULL 00128 ,reduct_obj 00129 ,this_offset + first_ele_in // first_ele 00130 ,this_sub_dim // sub_dim 00131 ,global_offset_in // global_dim 00132 ); 00133 } 00134 00135 value_type VectorSubView::get_ele(index_type i) const 00136 { 00137 space_.validate_range(Range1D(i,i)); 00138 return full_vec_->get_ele( space_.rng().lbound() + i - 1 ); 00139 } 00140 00141 Vector::vec_ptr_t 00142 VectorSubView::sub_view( const Range1D& rng_in ) const 00143 { 00144 namespace rcp = MemMngPack; 00145 const index_type this_dim = this->dim(); 00146 const Range1D rng = RangePack::full_range(rng_in,1,this_dim); 00147 space_.validate_range(rng); 00148 const index_type this_offset = space_.rng().lbound() - 1; 00149 return Teuchos::rcp( 00150 new VectorSubView( 00151 full_vec_ 00152 ,Range1D( 00153 this_offset + rng.lbound() 00154 ,this_offset + rng.ubound() ) 00155 ) ); 00156 } 00157 00158 void VectorSubView::get_sub_vector( const Range1D& rng_in, RTOpPack::SubVector* sub_vec ) const 00159 { 00160 #ifdef TEUCHOS_DEBUG 00161 TEST_FOR_EXCEPTION( !sub_vec, std::logic_error, "VectorSubView::get_sub_vector(...): Error!" ); 00162 #endif 00163 const index_type this_dim = this->dim(); 00164 const Range1D rng = RangePack::full_range(rng_in,1,this_dim); 00165 space_.validate_range(rng); 00166 const index_type this_offset = space_.rng().lbound() - 1; 00167 full_vec_->get_sub_vector( rng + this_offset, sub_vec ); 00168 sub_vec->setGlobalOffset( sub_vec->globalOffset() - this_offset ); 00169 } 00170 00171 void VectorSubView::free_sub_vector( RTOpPack::SubVector* sub_vec ) const 00172 { 00173 #ifdef TEUCHOS_DEBUG 00174 TEST_FOR_EXCEPTION( !sub_vec, std::logic_error, "VectorSubView::free_sub_vector(...): Error!" ); 00175 #endif 00176 const index_type this_offset = space_.rng().lbound() - 1; 00177 sub_vec->setGlobalOffset( sub_vec->globalOffset() + this_offset ); 00178 full_vec_->free_sub_vector( sub_vec ); 00179 } 00180 00181 } // end namespace AbstractLinAlgPack
1.7.4