Logo  0.95.0-final
Finite Element Embedded Library and Language in C++
Feel++ Feel++ on Github Feel++ community
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
matrixepetra.hpp
Go to the documentation of this file.
1 /* -*- mode: c++; coding: utf-8; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; show-trailing-whitespace: t -*- vim:fenc=utf-8:ft=tcl:et:sw=4:ts=4:sts=4
2 
3  This file is part of the Feel library
4 
5  Author(s):
6  Christophe Prud'homme <christophe.prudhomme@feelpp.org>
7  Klaus Sapelza <klaus.sapelza@epfl.ch>
8  Date: 2006-09-14
9 
10  Copyright (C) 2006,2007 EPFL
11  Copyright (C) 2006,2007 Universit� Joseph Fourier (Grenoble I)
12 
13  This library is free software; you can redistribute it and/or
14  modify it under the terms of the GNU Lesser General Public
15  License as published by the Free Software Foundation; either
16  version 3.0 of the License, or (at your option) any later version.
17 
18  This library is distributed in the hope that it will be useful,
19  but WITHOUT ANY WARRANTY; without even the implied warranty of
20  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
21  Lesser General Public License for more details.
22 
23  You should have received a copy of the GNU Lesser General Public
24  License along with this library; if not, write to the Free Software
25  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
26 */
34 #ifndef __MatrixEpetra_H
35 #define __MatrixEpetra_H 1
36 
37 #include <feel/feelconfig.h>
38 
39 
41 
44 
45 #if defined( FEELPP_HAS_TRILINOS_EPETRA )
46 #undef PACKAGE_BUGREPORT
47 #undef PACKAGE_NAME
48 #undef PACKAGE_STRING
49 #undef PACKAGE_TARNAME
50 #undef PACKAGE_VERSION
51 #if defined(FEELPP_HAS_MPI)
52 #include <Epetra_MpiComm.h>
53 #else
54 #include <Epetra_SerialComm.h>
55 #endif /* FEELPP_HAS_MPI */
56 #include <Epetra_Map.h>
57 #include <EpetraExt_MatrixMatrix.h>
58 #include <EpetraExt_RowMatrixOut.h>
59 #include <Epetra_FECrsMatrix.h>
60 #include <Epetra_Vector.h>
61 #undef PACKAGE_BUGREPORT
62 #undef PACKAGE_NAME
63 #undef PACKAGE_STRING
64 #undef PACKAGE_TARNAME
65 #undef PACKAGE_VERSION
66 
67 namespace Feel
68 {
69 template<typename T> class VectorEpetra;
70 
84 class MatrixEpetra : public MatrixSparse<double>
85 {
86  typedef MatrixSparse<double> super;
87 
88 
89 public:
94  static const bool is_row_major = true;
95  typedef super::value_type value_type;
96  typedef super::real_type real_type;
97  typedef std::vector<std::set<size_type> > pattern_type;
98 
99  typedef super::graph_type graph_type;
100  typedef super::graph_ptrtype graph_ptrtype;
101 
102  typedef Vector<value_type> vector_type;
103  typedef boost::shared_ptr<Vector<value_type> > vector_ptrtype;
104 
105  typedef VectorEpetra<value_type> epetra_vector_type;
106  typedef boost::shared_ptr<epetra_vector_type> epetra_vector_ptrtype;
107 
108  typedef MatrixEpetra epetra_sparse_matrix_type;
109 
111 
112 
116 
125  MatrixEpetra ( Epetra_FECrsMatrix const& m );
126 
130  MatrixEpetra ( MatrixEpetra const& T );
131 
135  MatrixEpetra ( Epetra_Map const& emap, int nnz = 50 );
136 
137  MatrixEpetra ( Epetra_Vector const& x );
138 
139  MatrixEpetra ( Epetra_Map const& rowmap, Epetra_Map const& colmap );
140 
141  MatrixEpetra( Epetra_Map const& row_emap, Epetra_Map const& col_emap,
142  Epetra_Map const& dom_map, Epetra_Map const& range_map );
143 
144 
145  MatrixEpetra ( const size_type m,
146  const size_type n,
147  const size_type m_l,
148  const size_type n_l,
149  const size_type nnz,
150  const size_type /*noz*/ );
151 
158  ~MatrixEpetra();
159 
161 
165 
166  MatrixEpetra & operator = ( MatrixSparse<value_type> const& M )
167  {
168  return *this;
169  }
181  value_type operator () ( const size_type i,
182  const size_type j ) const;
183 
184 
190  MatrixEpetra& operator=( const MatrixEpetra& T )
191  {
192  if ( &T != this )
193  {
194  M_emap = T.getRowMap();
195  M_col_emap = T.getColMap();
196  M_dom_map = T.mat().DomainMap();
197  M_range_map = T.mat().RangeMap();
198  *M_mat = T.mat();
199 
200  this->setInitialized( T.isInitialized() );
201  }
202 
203  return *this;
204  }
205 
207 
211 
216  size_type size1 () const;
217 
222  size_type size2 () const;
223 
228  size_type rowStart () const;
229 
234  size_type rowStop () const;
235 
236 
240  Epetra_Map getRowMap() const
241  {
242  return M_emap;
243  }
244 
248  Epetra_Map getColMap() const
249  {
250  return M_col_emap;
251  }
252 
256  Epetra_Map getDomainMap() const
257  {
258  return M_mat->DomainMap();
259  }
260 
264  Epetra_Map getRangeMap() const
265  {
266  return M_mat->RangeMap();
267  }
269 
270 
271  std::vector<size_type> bcIndices()
272  {
273  return M_bc_index;
274  }
275 
276  boost::shared_ptr<Epetra_FECrsMatrix> matrix()
277  {
278  return M_mat;
279  }
280 
284 
285 
287 
291 
300  void init ( const size_type m,
301  const size_type n,
302  const size_type m_l,
303  const size_type n_l,
304  const size_type nnz=10,
305  const size_type noz=10 );
306 
310  void init ( const size_type m,
311  const size_type n,
312  const size_type m_l,
313  const size_type n_l,
314  graph_ptrtype const& graph );
315 
316 
317 
321  void init ();
322 
326  // void resize( size_type m, size_type n, bool /*preserve*/ )
327  // {
328  // this->setInitialized( false );
329  // init( m, n, m, n );
330  // }
331 
332  // void fill( pattern_type const& patt )
333  // {}
334 
341  void clear ();
342 
343 
348  void zero ();
349 
355  void zero ( size_type /*start1*/, size_type /*stop1*/, size_type /*start2*/, size_type /*stop2*/ )
356  {
357  }
358 
362  void diagonal ( Vector<double>& dest ) const;
363 
364  void setDiagonal( Epetra_Vector const& x );
365 
382  void close () const;
383 
384 
389  bool closed() const;
390 
391  bool EpetraIndicesAreLocal() const;
392 
393  bool HaveColMap() const;
394 
400  void transpose( MatrixSparse<value_type>& Mt ) const;
401 
405  real_type energy ( vector_type const& v1, vector_type const& v2, bool tranpose = false ) const;
406 
414  real_type l1Norm () const;
415 
426  real_type linftyNorm () const;
427 
434  void set ( const size_type i,
435  const size_type j,
436  const value_type& value );
437 
441  void multiplyMatrix ( const MatrixEpetra& A, const MatrixEpetra& B );
442 
443  template<typename T>
444  void multiply ( bool trans, const Vector<T>& v, Vector<T>& r ) const
445  {
446  epetra_vector_type const& ev( dynamic_cast<epetra_vector_type const&>( v ) );
447  epetra_vector_type& er( dynamic_cast<epetra_vector_type&>( r ) );
448  M_mat->Multiply( trans, ev.vec(), er.vec() );
449  }
450 
451  void multiply ( bool trans, const Epetra_MultiVector &X, Epetra_MultiVector &Y ) const
452  {
453  M_mat->Multiply( trans, X, Y );
454  }
455 
456 
465  void add ( const size_type i,
466  const size_type j,
467  const value_type& value );
468 
475  void addMatrix ( double coeff, const MatrixEpetra& _m )
476  {
477  EpetraExt::MatrixMatrix::Add( _m.mat(), false, coeff, ( *this ).mat(), 1. );
478  }
479 
480  void addMatrix ( const ublas::matrix<value_type> &/*dm*/,
481  const std::vector<size_type> &/*rows*/,
482  const std::vector<size_type> &/*cols*/ )
483  {
484  }
485 
494  void addMatrix ( int* rows, int nrows,
495  int* cols, int ncols,
496  value_type* data );
497 
502  void addMatrix ( const ublas::matrix<value_type> &dm,
503  const std::vector<size_type> &dof_indices )
504  {
505  this->addMatrix ( dm, dof_indices, dof_indices );
506  }
507 
508 
520  void addMatrix ( const value_type a, super &X );
521 
522 
526  void scale ( double const scalar );
527 
528 
532  Epetra_FECrsMatrix& mat ()
533  {
534  //FEELPP_ASSERT (M_mat != NULL).error("null epetra matrix");
535  return *M_mat;
536  }
537 
541  Epetra_FECrsMatrix const& mat () const
542  {
543  //FEELPP_ASSERT (M_mat != NULL).error("null epetra matrix");
544  return *M_mat;
545  }
546 
553  void printMatlab( const std::string name="NULL" ) const;
554  void printKonsole() const;
555 
556 
563 
571  void zeroRows( std::vector<int> const& rows, std::vector<value_type> const& values, Vector<value_type>& rhs, Context const& on_context );
572 
573 
577  void updateBlockMat( boost::shared_ptr<MatrixSparse<value_type> > m, std::vector<size_type> start_i, std::vector<size_type> start_j );
578 
579 protected:
580 
581 private:
582 
583  mpi::communicator M_comm;
584 
600  MatrixEpetra();
601 
606  Epetra_Map M_emap;
607  Epetra_Map M_col_emap;
608  Epetra_Map M_dom_map;
609  Epetra_Map M_range_map;
610 
611  mutable boost::shared_ptr<Epetra_FECrsMatrix> M_mat;
612 
613  std::vector<size_type> M_bc_index;
614 };
615 
616 
617 
618 inline
619 MatrixEpetra::MatrixEpetra()
620  :
621  super(),
622 
623 #ifdef FEELPP_HAS_MPI
624  M_emap( Epetra_Map( -1, 0, 0, Epetra_MpiComm( M_comm ) ) ),
625  M_col_emap( Epetra_Map( -1, 0, 0, Epetra_MpiComm( M_comm ) ) ),
626  M_dom_map( Epetra_Map( -1, 0, 0, Epetra_MpiComm( M_comm ) ) ),
627  M_range_map( Epetra_Map( -1, 0, 0, Epetra_MpiComm( M_comm ) ) ),
628  M_mat( new Epetra_FECrsMatrix( Copy, M_emap, 0 ) )
629 #else
630  M_emap( Epetra_Map( -1, 0, 0, Epetra_SerialComm ) ),
631  M_col_emap( Epetra_Map( -1, 0, 0, Epetra_SerialComm ) ),
632  M_dom_map( Epetra_Map( -1, 0, 0, Epetra_SerialComm ) ),
633  M_range_map( Epetra_Map( -1, 0, 0, Epetra_SerialComm ) ),
634  M_mat( new Epetra_FECrsMatrix( Copy, M_emap, 0 ) )
635 
636 #endif
637 {}
638 
639 inline
640 MatrixEpetra::MatrixEpetra( Epetra_Map const& emap, int nnz )
641  :
642  super(),
643  M_emap( emap ),
644  M_col_emap( emap ),
645  M_dom_map( emap ),
646  M_range_map( emap ),
647  M_mat( new Epetra_FECrsMatrix( Copy, emap, 0 ) )
648 {
649  //this->setInitialized( true );
650 }
651 
652 
653 inline
654 MatrixEpetra::MatrixEpetra( Epetra_Map const& row_emap, Epetra_Map const& col_emap )
655  :
656  super(),
657  M_emap( row_emap ),
658  M_col_emap( col_emap ),
659  M_dom_map( col_emap ),
660  M_range_map( row_emap ),
661  M_mat( new Epetra_FECrsMatrix( Copy, row_emap, col_emap, 0 ) )
662 {
663  //this->setInitialized( true );
664 
665 }
666 
667 inline
668 MatrixEpetra::MatrixEpetra( Epetra_Map const& row_emap, Epetra_Map const& col_emap,
669  Epetra_Map const& dom_map, Epetra_Map const& range_map )
670  :
671  super(),
672  M_emap( row_emap ),
673  M_col_emap( col_emap ),
674  M_dom_map( dom_map ),
675  M_range_map( range_map ),
676  M_mat( new Epetra_FECrsMatrix( Copy, M_emap, M_col_emap, 0 ) )
677 {
678  //this->setInitialized( true );
679 
680 }
681 
682 
683 inline
684 MatrixEpetra::MatrixEpetra( Epetra_FECrsMatrix const& M )
685  :
686  super(),
687  M_emap( M.RowMap() ),
688  M_col_emap( M.ColMap() ),
689  M_dom_map( M.DomainMap() ),
690  M_range_map( M.RangeMap() ),
691  M_mat( new Epetra_FECrsMatrix( Copy, M_emap, M_col_emap, 0 ) )
692 {
693 }
694 
695 
696 inline
697 MatrixEpetra::MatrixEpetra( MatrixEpetra const& T )
698  :
699  super(),
700  M_emap( T.M_emap ),
701  M_col_emap( T.M_col_emap ),
702  M_dom_map( T.M_dom_map ),
703  M_range_map( T.M_range_map ),
704  M_mat( new Epetra_FECrsMatrix( T.mat() ) )
705 
706 {
707 }
708 
709 
710 inline
711 MatrixEpetra::MatrixEpetra( const size_type m,
712  const size_type n,
713  const size_type m_l,
714  const size_type /*n_l*/,
715  const size_type nnz,
716  const size_type /*noz*/ )
717  :
718  super(),
719  M_emap( Epetra_Map( m, m_l, 0, Epetra_MpiComm( M_comm ) ) ),
720  M_col_emap( Epetra_Map( n, n, 0, Epetra_MpiComm( M_comm ) ) ),
721  M_dom_map( M_emap ),
722  M_range_map( M_emap ),
723  M_mat( new Epetra_FECrsMatrix( Copy, M_emap, M_col_emap, nnz ) )
724 {
725 }
726 
727 
728 inline
729 MatrixEpetra::~MatrixEpetra()
730 {
731 
732 }
733 
734 
735 inline
736 void MatrixEpetra::init ( const size_type m,
737  const size_type n,
738  const size_type /*m_l*/,
739  const size_type /*n_l*/,
740  const size_type /*nnz*/,
741  const size_type /*noz*/ )
742 {
743  if ( ( m==0 ) || ( n==0 ) )
744  return;
745 
746  {
747  // Clear initialized matrices
748  if ( this->isInitialized() )
749  this->clear();
750 
751  this->setInitialized( true );
752  }
753 }
754 
755 
756 
757 
758 
759 inline
760 void MatrixEpetra::init ()
761 {
762  if ( this->isInitialized() )
763  this->clear();
764 
765  this->setInitialized( true );
766 }
767 
768 inline
769 void MatrixEpetra::zero ()
770 {
771  FEELPP_ASSERT ( this->isInitialized() ).error( "epetra matrix not properly initialized" ) ;
772 
773  M_bc_index.resize( 0 );
774  M_mat->PutScalar( 0.0 );
775 
776 }
777 
778 
779 inline
780 void MatrixEpetra::clear ()
781 {
782  if ( this->isInitialized() )
783  {
784  M_mat = boost::shared_ptr<Epetra_FECrsMatrix>( new Epetra_FECrsMatrix( Copy, M_emap, M_col_emap, 50 ) );
785 
786  this->setInitialized( false );
787  }
788 }
789 
790 
791 
792 
793 inline
794 void MatrixEpetra::close () const
795 {
796  int ierr=0;
797 
798  if ( !this->closed() )
799  ierr = M_mat->GlobalAssemble( M_dom_map, M_range_map, true );
800 
801  if ( ierr != 0 )
802  {
803  DVLOG(2) << "ERRCODE GlobalAssemble: " << ierr << "\n";
804  }
805 }
806 
807 
808 inline
809 void MatrixEpetra::setDiagonal ( Epetra_Vector const& x )
810 {
811  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );;
812 
813  const Epetra_Map& rowMap( M_mat->RowMatrixRowMap() );
814  const Epetra_Map& colMap( M_mat->RowMatrixColMap() );
815 
816  size_type L = x.MyLength();
817 
818  value_type zero_value = static_cast<value_type>( 0.0 );
819 
820  for ( size_type i=0; i< L; i++ )
821  {
822  int i_val = static_cast<int>( rowMap.GID( i ) );
823  int j_val = static_cast<int>( colMap.GID( i ) );
824 
825  M_mat->InsertGlobalValues( 1, &i_val, 1, &j_val, &zero_value );
826  }
827 
828  this->close();
829 
830  M_mat->ReplaceDiagonalValues( x );
831 }
832 
833 
834 
835 inline
836 size_type MatrixEpetra::size1 () const
837 {
838  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );;
839 
840  DVLOG(2) << "Size in size1(): " << M_mat->NumGlobalRows() << "\n";
841 
842  int epetra_m = M_mat->NumGlobalRows();
843 
844  return static_cast<size_type>( epetra_m );
845 }
846 
847 
848 
849 inline
850 size_type MatrixEpetra::size2 () const
851 {
852  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );;
853 
854  int epetra_n = M_mat->NumGlobalCols();
855 
856  return static_cast<size_type>( epetra_n );
857 }
858 
859 
860 
861 inline
862 size_type MatrixEpetra::rowStart () const
863 {
864  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );;
865 
866  int start = M_emap.MinMyGID();
867 
868  return static_cast<size_type>( start );
869 }
870 
871 
872 
873 
874 inline
875 size_type MatrixEpetra::rowStop () const
876 {
877  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );;
878 
879  int stop = M_emap.MaxMyGID();
880 
881  return static_cast<size_type>( stop );
882 }
883 
884 
885 inline
886 void MatrixEpetra::set ( const size_type i,
887  const size_type j,
888  const value_type& value )
889 {
890  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );;
891 
892  int i_val = static_cast<int>( i );
893  int j_val = static_cast<int>( j );
894 
895  int ierr = 0;
896 
897  value_type epetra_value = static_cast<value_type>( value );
898 
899  ierr = M_mat->ReplaceGlobalValues( 1, &i_val, 1, &j_val, &epetra_value );
900 
901  if ( ierr )
902  {
903  ierr = M_mat->InsertGlobalValues( 1, &i_val, 1, &j_val, &epetra_value );
904 
905  if ( ierr != 0 )
906  {
907  DVLOG(2) << "ERRORCODE InsertGlobalValues: " << ierr << " in M(" << i_val << "," << j_val << ") for value "<< epetra_value << "." << "\n";
908  }
909  }
910 
911  //FEELPP_ASSERT( ierr == 0 )( ierr ).warn ( "invalid MatrixEpetra::set operation" );
912 }
913 
914 
915 
916 
917 inline
918 bool MatrixEpetra::closed() const
919 {
920  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );
921 
922  bool filled;
923  filled = M_mat->Filled();
924 
925  return static_cast<bool>( filled );
926 }
927 
928 
929 inline
930 bool MatrixEpetra::EpetraIndicesAreLocal() const
931 {
932  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );
933 
934  bool IsLocal;
935  IsLocal = M_mat->IndicesAreLocal();
936 
937  return static_cast<bool>( IsLocal );
938 }
939 
940 
941 inline
942 bool MatrixEpetra::HaveColMap() const
943 {
944  FEELPP_ASSERT ( this->isInitialized() ).error( "MatrixEpetra<> not properly initialized" );
945 
946  bool HaveMap;
947  HaveMap = M_mat->HaveColMap();
948 
949  return static_cast<bool>( HaveMap );
950 }
951 
952 inline
953 void
954 MatrixEpetra::addMatrix ( const value_type coeff, super &X )
955 {
956  FEELPP_ASSERT ( this->isInitialized() ).error( "epetra matrix not initialized" );
957 
958  epetra_sparse_matrix_type* A_ptr = dynamic_cast< epetra_sparse_matrix_type*>( &X );
959 
960  this->addMatrix( coeff, *A_ptr );
961 }
962 
963 inline
964 void
965 MatrixEpetra::scale ( double scalar )
966 {
967  FEELPP_ASSERT ( this->isInitialized() ).error( "epetra matrix not initialized" );
968  this->close();
969  M_mat->Scale( scalar );
970 }
971 
972 inline
973 MatrixEpetra::real_type
974 MatrixEpetra::l1Norm() const
975 {
976  FEELPP_ASSERT ( this->isInitialized() ).error( "epetra matrix not initialized" );
977  this->close();
978 
979  real_type NormOne = M_mat->NormOne();
980 
981  return static_cast<real_type>( NormOne );
982 }
983 
984 inline
985 MatrixEpetra::real_type
986 MatrixEpetra::linftyNorm() const
987 {
988  FEELPP_ASSERT ( this->isInitialized() ).error( "epetra matrix not initialized" );
989  this->close();
990 
991  real_type NormInf = M_mat->NormInf();
992 
993  return static_cast<real_type>( NormInf );
994 }
995 
996 
997 // This is for a GLOBAL matrix!
998 
999 inline
1000 MatrixEpetra::value_type //doesnt work in parallel yet
1001 MatrixEpetra::operator () ( const size_type i,
1002  const size_type j ) const
1003 {
1004  FEELPP_ASSERT ( this->isInitialized() ).error( "epetra matrix not initialized" );
1005 
1006 
1007  int i_val=static_cast<int>( i ),
1008  j_val=static_cast<int>( j );
1009 
1010  int NumEntries;
1011  double* Values;
1012  int* Indices;
1013 
1014 
1015  // int ierr = M_mat->ExtractMyRowView( i_val, NumEntries, Values, Indices);
1016  M_mat->ExtractMyRowView( i_val, NumEntries, Values, Indices );
1017 
1018  for ( int k = 0; k < NumEntries; ++k )
1019  {
1020  if ( Indices[k] == j_val )
1021  return static_cast<double> ( Values[k] );
1022  }
1023 
1024  /*
1025  int *epetra_cols;
1026  double *epetra_row;
1027 
1028  int
1029  ncols=0,
1030  i_val=static_cast<int>(i),
1031  j_val=static_cast<int>(j);
1032 
1033  // the matrix must NOT be closed
1034  // this->close();
1035 
1036 
1037  M_mat->ExtractGlobalRowView(i_val, ncols, epetra_row, epetra_cols);
1038 
1039  // Perform a binary search to find the contiguous index in
1040  // petsc_cols (resp. petsc_row) corresponding to global index j_val
1041  std::pair<const int*, const int*> p = std::equal_range (&epetra_cols[0], &epetra_cols[0] + ncols, j_val);
1042 
1043  std::cout << "Epetra row: " << *epetra_row << "\n";
1044  std::cout << "Epetra col: " << *epetra_cols << "\n";
1045 
1046 
1047  // Found an entry for j_val
1048  if (p.first != p.second)
1049  {
1050  // The entry in the contiguous row corresponding
1051  // to the j_val column of interest
1052  const int j = std::distance (const_cast<int*>(&epetra_cols[0]),
1053  const_cast<int*>(p.first));
1054 
1055  assert (j < ncols);
1056  assert (epetra_cols[j] == j_val);
1057 
1058  return static_cast<double> (epetra_row[j]);
1059 
1060  //return value;
1061  }
1062 
1063  std::cout << "Element " << epetra_row[j] << "\n";*/
1064 
1065  // Otherwise the entry is not in the sparse matrix,
1066  // i.e. it is 0.
1067  return 0.;
1068 
1069  //return M_mat[i_val][j_val]; //no good: locally indexed!
1070 }
1071 
1072 inline
1073 void
1074 MatrixEpetra::updateBlockMat( boost::shared_ptr<MatrixSparse<value_type> > m, std::vector<size_type> start_i, std::vector<size_type> start_j )
1075 {
1076  LOG(ERROR) << "Invalid call to updateBlockMat, not yet implemented\n";
1077 }
1078 
1079 
1080 inline
1081 NdebugStream&
1082 operator<<( NdebugStream& __os, MatrixEpetra const& /*__n*/ )
1083 {
1084  return __os;
1085 }
1086 
1087 } // Feel
1088 
1089 #endif /* FEELPP_HAS_TRILINOS_EPETRA */
1090 #endif /* __MatrixEpetra_H */

Generated on Sun Oct 20 2013 08:25:01 for Feel++ by doxygen 1.8.4