29 #ifndef __VectorEigen_H
30 #define __VectorEigen_H 1
33 #include <boost/operators.hpp>
36 #include <feel/feelalg/vector.hpp>
47 struct get_real_type<double>
49 typedef double value_type;
53 struct get_real_type<std::complex<double> >
55 typedef double value_type;
72 , boost::addable<VectorEigen<T> >
73 , boost::subtractable<VectorEigen<T> >
85 typedef typename get_real_type<T>::value_type real_type;
87 typedef Eigen::Matrix<value_type,Eigen::Dynamic,1> vector_type;
89 typedef typename super1::clone_ptrtype clone_ptrtype;
92 typedef typename super1::datamap_ptrtype datamap_ptrtype;
116 clone_ptrtype
clone ()
const ;
132 const bool fast=
false );
138 const bool fast=
false );
143 void init( datamap_ptrtype
const& dm );
161 FEELPP_ASSERT ( this->
isInitialized() ).error(
"vector not initialized" );
176 FEELPP_ASSERT ( this->
isInitialized() ).error(
"vector not initialized" );
190 FEELPP_ASSERT ( this->
isInitialized() ).error(
"vector not initialized" );
205 FEELPP_ASSERT ( this->
isInitialized() ).error(
"vector not initialized" );
298 vector_type
const&
vec ()
const
325 M_vec.setConstant( v );
353 M_vec.setZero( M_vec.size() );
375 for (
int j = 0; j < n; ++j )
376 M_vec( i[j] ) += v[j];
394 const std::vector<size_type>& dof_indices )
396 FEELPP_ASSERT ( v.size() == dof_indices.size() ).error(
"invalid dof indices" );
399 this->
add ( dof_indices[i], v[i] );
409 const std::vector<size_type>& dof_indices )
411 FEELPP_ASSERT ( V.
size() == dof_indices.size() ).error(
"invalid dof indices" );
414 this->
add ( dof_indices[i], V( i ) );
425 FEELPP_ASSERT( 0 ).error(
"invalid call, not implemented yet" );
435 const std::vector<size_type>& dof_indices )
437 FEELPP_ASSERT ( V.size() == dof_indices.size() ).error(
"invalid dof indices" );
440 this->
add ( dof_indices[i], V( i ) );
448 const std::vector<size_type>& )
450 FEELPP_ASSERT( 0 ).error(
"invalid call, not implemented yet" );
460 const std::vector<size_type>& )
462 FEELPP_ASSERT( 0 ).error(
"invalid call, not implemented yet" );
473 const std::vector<size_type>& )
475 FEELPP_ASSERT( 0 ).error(
"invalid call, not implemented yet" );
484 void insert (
const ublas::vector<T>& V,
485 const std::vector<size_type>& dof_indices );
503 void printMatlab(
const std::string name=
"NULL",
bool renumber =
false )
const;
515 real_type local_min = 0;
517 real_type global_min = local_min;
521 #ifdef FEELPP_HAS_MPI
525 MPI_Allreduce ( &local_min, &global_min, 1,
526 MPI_DOUBLE, MPI_MIN, this->
comm() );
541 real_type local_max = 0;
543 real_type global_max = local_max;
546 #ifdef FEELPP_HAS_MPI
550 MPI_Allreduce ( &local_max, &global_max, 1,
551 MPI_DOUBLE, MPI_MAX, this->
comm() );
566 double local_l1 = M_vec.array().abs().sum();
568 double global_l1 = local_l1;
571 #ifdef FEELPP_HAS_MPI
575 mpi::all_reduce( this->
comm(), local_l1, global_l1, std::plus<double>() );
591 real_type local_norm2 = M_vec.squaredNorm();
592 real_type global_norm2 = local_norm2;
595 #ifdef FEELPP_HAS_MPI
599 mpi::all_reduce( this->
comm(), local_norm2, global_norm2, std::plus<real_type>() );
603 return math::sqrt( global_norm2 );
613 real_type local_norminf = M_vec.array().abs().maxCoeff();
614 real_type global_norminf = local_norminf;
617 #ifdef FEELPP_HAS_MPI
621 mpi::all_reduce( this->
comm(), local_norminf, global_norminf, mpi::maximum<real_type>() );
625 return global_norminf;
635 value_type local_sum = M_vec.array().sum();
637 value_type global_sum = local_sum;
640 #ifdef FEELPP_HAS_MPI
644 mpi::all_reduce( this->
comm(), local_sum, global_sum, std::plus<value_type>() );
656 FEELPP_DONT_INLINE this_type
sqrt()
const;
662 this_type
pow(
int n )
const;
675 M_vec.operator[]( i ) += a;
712 FEELPP_ASSERT( 0 ).error(
"invalid call, not implemented yet" );
721 const std::vector<size_type>& send_list );
735 const std::vector<size_type>& send_list )
const;
753 throw std::logic_error(
"[vetor eigen] ERROR dot function not yet implemented" );
767 void checkInvariant()
const;
780 template <
typename T>
784 return v1.
vec().array()*v2.
vec().array();
794 template <
typename T>
799 return v1->vec().array()*v2->vec().array();
807 #if !defined( FEELPP_INSTANTIATE_VECTOREIGEN )
808 extern template class VectorEigen<double>;