AbstractLinAlgPack: C++ Interfaces For Vectors, Matrices And Related Linear Algebra Objects  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
AbstractLinAlgPack_apply_op_helper.cpp
1 /*
2 // @HEADER
3 // ***********************************************************************
4 //
5 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
6 // Copyright (2003) Sandia Corporation
7 //
8 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
9 // license for use of this work by or on behalf of the U.S. Government.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov)
39 //
40 // ***********************************************************************
41 // @HEADER
42 */
43 
44 // /////////////////////////////////////////////////////////////////////////////
45 // apply_op_helper.cpp
46 
47 #include <assert.h>
48 
49 #include "AbstractLinAlgPack_apply_op_helper.hpp"
50 #include "AbstractLinAlgPack_VectorSpace.hpp"
51 #include "AbstractLinAlgPack_VectorMutable.hpp"
52 #include "Teuchos_Workspace.hpp"
53 #include "Teuchos_Assert.hpp"
54 
55 void AbstractLinAlgPack::apply_op_validate_input(
56  const char func_name[]
57  ,const RTOpPack::RTOp& op
58  ,const size_t num_vecs, const Vector* vecs[]
59  ,const size_t num_targ_vecs, VectorMutable* targ_vecs[]
60  ,RTOpPack::ReductTarget *reduct_obj
61  ,const index_type first_ele_in, const index_type sub_dim_in, const index_type global_offset_in
62  )
63 {
64  const VectorSpace
65  &space = (num_vecs ? vecs[0]->space() : targ_vecs[0]->space() );
66  const index_type
67  dim = space.dim();
69  global_offset_in < 0, std::logic_error
70  ,func_name << " : Error! global_offset_in = "
71  <<global_offset_in<<" is not valid" );
73  first_ele_in > dim, std::logic_error
74  ,func_name << " : Error! first_ele_in = "
75  <<first_ele_in<<" is not compatible with space.dim() = " << dim );
77  sub_dim_in < 0 || (sub_dim_in > 0 && sub_dim_in > dim-(first_ele_in-1)), std::logic_error
78  ,func_name << " : Error! first_ele_in = "
79  <<first_ele_in<<" and sub_dim_in = "<<sub_dim_in
80  <<" is not compatible with space.dim() = " << dim );
81  {for(size_t k = 0; k < num_vecs; ++k) {
82  const bool is_compatible = space.is_compatible(vecs[k]->space());
84  !is_compatible, VectorSpace::IncompatibleVectorSpaces
85  ,func_name << " : Error! vecs["<<k<<"]->space() of type \'"
86  << typeName(vecs[k]->space()) << "\' "
87  << " with dimension vecs["<<k<<"].dim() = " << vecs[k]->dim()
88  << " is not compatible with space of type \'"
89  << typeName(space) << " with dimmension space.dim() = " << dim );
90  }}
91  {for(size_t k = 0; k < num_targ_vecs; ++k) {
92  const bool is_compatible = space.is_compatible(targ_vecs[k]->space());
94  !is_compatible, VectorSpace::IncompatibleVectorSpaces
95  ,func_name << " : Error! targ_vecs["<<k<<"]->space() of type \'"
96  << typeName(targ_vecs[k]->space()) << "\' "
97  << " with dimension targ_vecs["<<k<<"].dim() = " << targ_vecs[k]->dim()
98  << " is not compatible with space of type \'"
99  << typeName(space) << " with dimmension space.dim() = " << dim );
100  }}
101 }
102 
103 void AbstractLinAlgPack::apply_op_serial(
104  const RTOpPack::RTOp& op
105  ,const size_t num_vecs, const Vector* vecs[]
106  ,const size_t num_targ_vecs, VectorMutable* targ_vecs[]
107  ,RTOpPack::ReductTarget *reduct_obj
108  ,const index_type first_ele_in, const index_type sub_dim_in, const index_type global_offset_in
109  )
110 {
111  using Teuchos::Workspace;
113 
114  // Dimension of global sub-vector
115  const VectorSpace
116  &space = ( num_vecs ? vecs[0]->space() : targ_vecs[0]->space() );
117  const index_type
118  full_dim = space.dim(),
119  global_sub_dim = sub_dim_in ? sub_dim_in : full_dim - (first_ele_in-1);
120  const Range1D
121  global_sub_rng = Range1D(first_ele_in,(first_ele_in-1)+global_sub_dim);
122 
123  //
124  // Get explicit views of the vector elements
125  //
126 
127  Workspace<RTOpPack::ConstSubVectorView<value_type> > local_vecs(wss,num_vecs);
128  Workspace<RTOpPack::SubVectorView<value_type> > local_targ_vecs(wss,num_targ_vecs);
129  size_t k;
130  for(k = 0; k < num_vecs; ++k) {
131  RTOpPack::SubVector _v;
132  vecs[k]->get_sub_vector( global_sub_rng, &_v );
133  (local_vecs[k] = _v).setGlobalOffset( _v.globalOffset() + global_offset_in );
134  }
135  for(k = 0; k < num_targ_vecs; ++k) {
136  RTOpPack::MutableSubVector _v;
137  targ_vecs[k]->get_sub_vector( global_sub_rng, &_v );
138  (local_targ_vecs[k] = _v).setGlobalOffset( _v.globalOffset() + global_offset_in );
139  }
140 
141  //
142  // Apply the reduction/transformation operator on all elements all at once!
143  //
144 
145  op.apply_op( local_vecs(), local_targ_vecs(), Teuchos::ptr(reduct_obj) );
146 
147  //
148  // Free (and commit) the explicit views of the vector elements
149  // which should also inform the vectors that they have
150  // changed.
151  //
152 
153  for(k = 0; k < num_vecs; ++k) {
154  RTOpPack::ConstSubVectorView<value_type> &v = local_vecs[k];
155  v.setGlobalOffset( v.globalOffset() - global_offset_in );
156  RTOpPack::SubVector _v = v;
157  vecs[k]->free_sub_vector(&_v);
158  }
159  for(k = 0; k < num_targ_vecs; ++k) {
160  RTOpPack::SubVectorView<value_type> &v = local_targ_vecs[k];
161  v.setGlobalOffset( v.globalOffset() - global_offset_in );
162  RTOpPack::MutableSubVector _v = v;
163  targ_vecs[k]->commit_sub_vector(&_v);
164  }
165 
166 }
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Ordinal globalOffset() const
void setGlobalOffset(Ordinal globalOffset_in)
TEUCHOSCORE_LIB_DLL_EXPORT Teuchos::RCP< WorkspaceStore > get_default_workspace_store()
std::string typeName(const T &t)