Thyra_SerialVectorBase.hpp

// @HEADER
// ***********************************************************************
// 
//    Thyra: Interfaces and Support for Abstract Numerical Algorithms
//                 Copyright (2004) Sandia Corporation
// 
// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
// license for use of this work by or on behalf of the U.S. Government.
// 
// This library is free software; you can redistribute it and/or modify
// it under the terms of the GNU Lesser General Public License as
// published by the Free Software Foundation; either version 2.1 of the
// License, or (at your option) any later version.
//  
// This library is distributed in the hope that it will be useful, but
// WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.
//  
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
// USA
// Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
// 
// ***********************************************************************
// @HEADER

#ifndef THYRA_VECTOR_SERIAL_BASE_HPP
#define THYRA_VECTOR_SERIAL_BASE_HPP

#include "Thyra_SerialVectorBaseDecl.hpp"
#include "Thyra_apply_op_helper.hpp"
#include "Teuchos_Workspace.hpp"
#include "Teuchos_TestForException.hpp"

namespace Thyra {

template<class Scalar>
SerialVectorBase<Scalar>::SerialVectorBase()
  :in_applyOp_(false)
{}

// Virtual methods with default implementations

template<class Scalar>
void SerialVectorBase<Scalar>::getData( const Scalar** values, Index* stride ) const
{
  const_cast<SerialVectorBase<Scalar>*>(this)->getData(const_cast<Scalar**>(values),stride);
}

template<class Scalar>
void SerialVectorBase<Scalar>::freeData( const Scalar** values ) const
{
  const_cast<SerialVectorBase<Scalar>*>(this)->commitData(const_cast<Scalar**>(values));
}

// Overridden from VectorBase

template<class Scalar>
void SerialVectorBase<Scalar>::applyOp(
  const RTOpPack::RTOpT<Scalar>   &op
  ,const int                      num_vecs
  ,const VectorBase<Scalar>*      vecs[]
  ,const int                      num_targ_vecs
  ,VectorBase<Scalar>*            targ_vecs[]
  ,RTOpPack::ReductTarget         *reduct_obj
  ,const Index                    first_ele
  ,const Index                    sub_dim
  ,const Index                    global_offset
  ) const
{
#ifdef _DEBUG
  TEST_FOR_EXCEPTION(
    in_applyOp_, std::logic_error
    ,"SerialVectorBase::applyOp(...): Error, something is not right here!" );
  Thyra::apply_op_validate_input(
    "SerialVectorBase<>::applyOp(...)",*this->space()
    ,op,num_vecs,vecs,num_targ_vecs,targ_vecs,reduct_obj,first_ele,sub_dim,global_offset
    );
#endif
  in_applyOp_ = true;
  Thyra::apply_op_serial(
    *(this->space())
    ,op,num_vecs,vecs,num_targ_vecs,targ_vecs,reduct_obj
    ,first_ele,sub_dim,global_offset
    );
  in_applyOp_ = false;
}

template<class Scalar>
void SerialVectorBase<Scalar>::getSubVector( const Range1D& rng_in, RTOpPack::SubVectorT<Scalar>* sub_vec ) const
{
  const Index   this_dim = this->space()->dim(); // ToDo: Cache this!
  const Range1D rng      = RangePack::full_range(rng_in,1,this_dim);
#ifdef _DEBUG
  TEST_FOR_EXCEPTION(
    rng.ubound() > this_dim, std::out_of_range
    ,"SerialVectorBase<Scalar>::getSubVector(...) : Error, "
    "rng = ["<<rng.lbound()<<","<<rng.ubound()<<"] "
    "is not in the range [1,this->dim()] = [1," << this_dim << "]!" );
#endif
  const Scalar *values = NULL; Index stride = 0;
  this->getData(&values,&stride);
  sub_vec->initialize(
    rng.lbound()-1                   // globalOffset
    ,rng.size()                      // subDim
    ,values+stride*(rng.lbound()-1)  // values
    ,stride                          // stride
    );
}

template<class Scalar>
void SerialVectorBase<Scalar>::freeSubVector( RTOpPack::SubVectorT<Scalar>* sub_vec ) const
{
  sub_vec->set_uninitialized();  // Nothing to deallocate!
}

template<class Scalar>
void SerialVectorBase<Scalar>::getSubVector( const Range1D& rng_in, RTOpPack::MutableSubVectorT<Scalar>* sub_vec )
{
  const Index   this_dim = this->space()->dim(); // ToDo: Cache this!
  const Range1D rng      = RangePack::full_range(rng_in,1,this_dim);
#ifdef _DEBUG
  TEST_FOR_EXCEPTION(
    rng.ubound() > this_dim, std::out_of_range
    ,"SerialVectorBase<Scalar>::getSubVector(...) : Error, "
    "rng = ["<<rng.lbound()<<","<<rng.ubound()<<"] "
    "is not in the range [1,this->dim()] = [1," << this_dim << "]!" );
#endif
  Scalar *values = NULL; Index stride = 0;
  this->getData(&values,&stride);
  sub_vec->initialize(
    rng.lbound()-1                   // globalOffset
    ,rng.size()                      // subDim
    ,values+stride*(rng.lbound()-1)  // values
    ,stride                          // stride
    );
}

template<class Scalar>
void SerialVectorBase<Scalar>::commitSubVector( RTOpPack::MutableSubVectorT<Scalar>* sub_vec )
{
  sub_vec->set_uninitialized();  // Nothing to deallocate!
}

template<class Scalar>
void SerialVectorBase<Scalar>::setSubVector( const RTOpPack::SparseSubVectorT<Scalar>& sub_vec )
{
  VectorBase<Scalar>::setSubVector(sub_vec);  // This implementation is okay
}

} // end namespace Thyra

#endif // THYRA_VECTOR_SERIAL_BASE_HPP

---