libMesh::PetscVector< T > Class Template Reference

#include <petsc_vector.h>

Inheritance diagram for libMesh::PetscVector< T >:

Public Member Functions

 PetscVector (const Parallel::Communicator &comm, const ParallelType type=AUTOMATIC)
 
 PetscVector (const Parallel::Communicator &comm, const numeric_index_type n, const ParallelType type=AUTOMATIC)
 
 PetscVector (const Parallel::Communicator &comm, const numeric_index_type n, const numeric_index_type n_local, const ParallelType type=AUTOMATIC)
 
 PetscVector (const Parallel::Communicator &comm, const numeric_index_type N, const numeric_index_type n_local, const std::vector< numeric_index_type > &ghost, const ParallelType type=AUTOMATIC)
 
 PetscVector (Vec v, const Parallel::Communicator &comm LIBMESH_CAN_DEFAULT_TO_COMMWORLD)
 
 ~PetscVector ()
 
void close ()
 
void clear ()
 
void zero ()
 
virtual AutoPtr< NumericVector
< T > > 
zero_clone () const
 
AutoPtr< NumericVector< T > > clone () const
 
void init (const numeric_index_type N, const numeric_index_type n_local, const bool fast=false, const ParallelType type=AUTOMATIC)
 
void init (const numeric_index_type N, const bool fast=false, const ParallelType type=AUTOMATIC)
 
virtual void init (const numeric_index_type, const numeric_index_type, const std::vector< numeric_index_type > &, const bool=false, const ParallelType=AUTOMATIC)
 
virtual void init (const NumericVector< T > &other, const bool fast=false)
 
NumericVector< T > & operator= (const T s)
 
NumericVector< T > & operator= (const NumericVector< T > &V)
 
PetscVector< T > & operator= (const PetscVector< T > &V)
 
NumericVector< T > & operator= (const std::vector< T > &v)
 
Real min () const
 
Real max () const
 
sum () const
 
Real l1_norm () const
 
Real l2_norm () const
 
Real linfty_norm () const
 
numeric_index_type size () const
 
numeric_index_type local_size () const
 
numeric_index_type first_local_index () const
 
numeric_index_type last_local_index () const
 
numeric_index_type map_global_to_local_index (const numeric_index_type i) const
 
operator() (const numeric_index_type i) const
 
virtual void get (const std::vector< numeric_index_type > &index, std::vector< T > &values) const
 
NumericVector< T > & operator+= (const NumericVector< T > &V)
 
NumericVector< T > & operator-= (const NumericVector< T > &V)
 
virtual void reciprocal ()
 
virtual void conjugate ()
 
void set (const numeric_index_type i, const T value)
 
void add (const numeric_index_type i, const T value)
 
void add (const T s)
 
void add (const NumericVector< T > &V)
 
void add (const T a, const NumericVector< T > &v)
 
void add_vector (const std::vector< T > &v, const std::vector< numeric_index_type > &dof_indices)
 
void add_vector (const NumericVector< T > &V, const std::vector< numeric_index_type > &dof_indices)
 
void add_vector (const NumericVector< T > &V, const SparseMatrix< T > &A)
 
void add_vector (const DenseVector< T > &V, const std::vector< numeric_index_type > &dof_indices)
 
void add_vector_transpose (const NumericVector< T > &V, const SparseMatrix< T > &A_trans)
 
void add_vector_conjugate_transpose (const NumericVector< T > &V, const SparseMatrix< T > &A_trans)
 
virtual void insert (const std::vector< T > &v, const std::vector< numeric_index_type > &dof_indices)
 
virtual void insert (const NumericVector< T > &V, const std::vector< numeric_index_type > &dof_indices)
 
virtual void insert (const DenseVector< T > &V, const std::vector< numeric_index_type > &dof_indices)
 
virtual void insert (const DenseSubVector< T > &V, const std::vector< numeric_index_type > &dof_indices)
 
void scale (const T factor)
 
virtual NumericVector< T > & operator/= (NumericVector< T > &v)
 
virtual void abs ()
 
virtual T dot (const NumericVector< T > &) const
 
virtual T indefinite_dot (const NumericVector< T > &) const
 
void localize (std::vector< T > &v_local) const
 
void localize (NumericVector< T > &v_local) const
 
void localize (NumericVector< T > &v_local, const std::vector< numeric_index_type > &send_list) const
 
void localize (const numeric_index_type first_local_idx, const numeric_index_type last_local_idx, const std::vector< numeric_index_type > &send_list)
 
void localize_to_one (std::vector< T > &v_local, const processor_id_type proc_id=0) const
 
virtual void pointwise_mult (const NumericVector< T > &vec1, const NumericVector< T > &vec2)
 
void print_matlab (const std::string name="NULL") const
 
virtual void create_subvector (NumericVector< T > &subvector, const std::vector< numeric_index_type > &rows) const
 
virtual void swap (NumericVector< T > &v)
 
Vec vec ()
 
template<>
void localize_to_one (std::vector< Real > &v_local, const processor_id_type pid) const
 
template<>
void localize_to_one (std::vector< Complex > &v_local, const processor_id_type pid) const
 
virtual bool initialized () const
 
ParallelType type () const
 
ParallelType & type ()
 
virtual bool closed () const
 
virtual Real subset_l1_norm (const std::set< numeric_index_type > &indices) const
 
virtual Real subset_l2_norm (const std::set< numeric_index_type > &indices) const
 
virtual Real subset_linfty_norm (const std::set< numeric_index_type > &indices) const
 
virtual T el (const numeric_index_type i) const
 
NumericVector< T > & operator*= (const T a)
 
NumericVector< T > & operator/= (const T a)
 
void add_vector (const NumericVector< T > &v, const ShellMatrix< T > &a)
 
virtual int compare (const NumericVector< T > &other_vector, const Real threshold=TOLERANCE) const
 
virtual int local_relative_compare (const NumericVector< T > &other_vector, const Real threshold=TOLERANCE) const
 
virtual int global_relative_compare (const NumericVector< T > &other_vector, const Real threshold=TOLERANCE) const
 
virtual void print (std::ostream &os=libMesh::out) const
 
template<>
void print (std::ostream &os) const
 
virtual void print_global (std::ostream &os=libMesh::out) const
 
template<>
void print_global (std::ostream &os) const
 
const Parallel::Communicatorcomm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static AutoPtr< NumericVector
< T > > 
build (const Parallel::Communicator &comm, const SolverPackage solver_package=libMesh::default_solver_package())
 
static AutoPtr< NumericVector
< T > > 
build (const SolverPackage solver_package=libMesh::default_solver_package())
 
static std::string get_info ()
 
static void print_info (std::ostream &out=libMesh::out)
 
static unsigned int n_objects ()
 
static void enable_print_counter_info ()
 
static void disable_print_counter_info ()
 

Protected Types

typedef std::map< std::string,
std::pair< unsigned int,
unsigned int > > 
Counts
 

Protected Member Functions

void increment_constructor_count (const std::string &name)
 
void increment_destructor_count (const std::string &name)
 

Protected Attributes

bool _is_closed
 
bool _is_initialized
 
ParallelType _type
 
const Parallel::Communicator_communicator
 

Static Protected Attributes

static Counts _counts
 
static Threads::atomic
< unsigned int > 
_n_objects
 
static Threads::spin_mutex _mutex
 
static bool _enable_print_counter = true
 

Private Types

typedef std::map
< numeric_index_type,
numeric_index_type
GlobalToLocalMap
 

Private Member Functions

void _get_array (void) const
 
void _restore_array (void) const
 

Private Attributes

Vec _vec
 
bool _array_is_present
 
numeric_index_type _local_size
 
Vec _local_form
 
PetscScalar * _values
 
GlobalToLocalMap _global_to_local_map
 
bool _destroy_vec_on_exit
 

Detailed Description

template<typename T>
class libMesh::PetscVector< T >

Petsc vector. Provides a nice interface to the Petsc C-based data structures for parallel vectors.

Author
Benjamin S. Kirk, 2002

Definition at line 61 of file petsc_vector.h.

Member Typedef Documentation

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts
protectedinherited

Data structure to log the information. The log is identified by the class name.

Definition at line 113 of file reference_counter.h.

template<typename T>
typedef std::map<numeric_index_type,numeric_index_type> libMesh::PetscVector< T >::GlobalToLocalMap
private

Type for map that maps global to local ghost cells.

Definition at line 596 of file petsc_vector.h.

Constructor & Destructor Documentation

template<typename T >
libMesh::PetscVector< T >::PetscVector ( const Parallel::Communicator comm,
const ParallelType  type = AUTOMATIC 
)
inlineexplicit

Dummy-Constructor. Dimension=0

Definition at line 618 of file petsc_vector.h.

References libMesh::NumericVector< T >::_type.

619  : NumericVector<T>(comm, ptype),
620  _array_is_present(false),
621  _local_form(NULL),
622  _values(NULL),
625 {
626  this->_type = ptype;
627 }
template<typename T >
libMesh::PetscVector< T >::PetscVector ( const Parallel::Communicator comm,
const numeric_index_type  n,
const ParallelType  type = AUTOMATIC 
)
inlineexplicit

Constructor. Set dimension to n and initialize all elements with zero.

Definition at line 633 of file petsc_vector.h.

References libMesh::PetscVector< T >::init().

636  : NumericVector<T>(comm, ptype),
637  _array_is_present(false),
638  _local_form(NULL),
639  _values(NULL),
642 {
643  this->init(n, n, false, ptype);
644 }
template<typename T >
libMesh::PetscVector< T >::PetscVector ( const Parallel::Communicator comm,
const numeric_index_type  n,
const numeric_index_type  n_local,
const ParallelType  type = AUTOMATIC 
)
inline

Constructor. Set local dimension to n_local, the global dimension to n, and initialize all elements with zero.

Definition at line 650 of file petsc_vector.h.

References libMesh::PetscVector< T >::init().

654  : NumericVector<T>(comm, ptype),
655  _array_is_present(false),
656  _local_form(NULL),
657  _values(NULL),
660 {
661  this->init(n, n_local, false, ptype);
662 }
template<typename T >
libMesh::PetscVector< T >::PetscVector ( const Parallel::Communicator comm,
const numeric_index_type  N,
const numeric_index_type  n_local,
const std::vector< numeric_index_type > &  ghost,
const ParallelType  type = AUTOMATIC 
)
inline

Constructor. Set local dimension to n_local, the global dimension to n, but additionally reserve memory for the indices specified by the ghost argument.

Definition at line 668 of file petsc_vector.h.

References libMesh::PetscVector< T >::init().

673  : NumericVector<T>(comm, ptype),
674  _array_is_present(false),
675  _local_form(NULL),
676  _values(NULL),
679 {
680  this->init(n, n_local, ghost, false, ptype);
681 }
template<typename T>
libMesh::PetscVector< T >::PetscVector ( Vec  v,
const Parallel::Communicator &comm  LIBMESH_CAN_DEFAULT_TO_COMMWORLD 
)

Constructor. Creates a PetscVector assuming you already have a valid PETSc Vec object. In this case, v is NOT destroyed by the PetscVector constructor when this object goes out of scope. This allows ownership of v to remain with the original creator, and to simply provide additional functionality with the PetscVector.

template<typename T >
libMesh::PetscVector< T >::~PetscVector ( )
inline

Destructor, deallocates memory. Made virtual to allow for derived classes to behave properly.

Definition at line 773 of file petsc_vector.h.

774 {
775  this->clear ();
776 }

Member Function Documentation

template<typename T >
void libMesh::PetscVector< T >::_get_array ( void  ) const
inlineprivate

Queries the array (and the local form if the vector is ghosted) from Petsc.

Definition at line 1308 of file petsc_vector.h.

References libMeshEnums::GHOSTED, libMesh::initialized(), and libMesh::libmesh_assert().

1309 {
1310  libmesh_assert (this->initialized());
1311  if(!_array_is_present)
1312  {
1313  PetscErrorCode ierr=0;
1314  if(this->type() != GHOSTED)
1315  {
1316  ierr = VecGetArray(_vec, &_values);
1317  LIBMESH_CHKERRABORT(ierr);
1318  }
1319  else
1320  {
1321  ierr = VecGhostGetLocalForm (_vec,&_local_form);
1322  LIBMESH_CHKERRABORT(ierr);
1323  ierr = VecGetArray(_local_form, &_values);
1324  LIBMESH_CHKERRABORT(ierr);
1325 #ifndef NDEBUG
1326  PetscInt my_local_size = 0;
1327  ierr = VecGetLocalSize(_local_form, &my_local_size);
1328  LIBMESH_CHKERRABORT(ierr);
1329  _local_size = static_cast<numeric_index_type>(my_local_size);
1330 #endif
1331  }
1332  _array_is_present = true;
1333  }
1334 }
template<typename T >
void libMesh::PetscVector< T >::_restore_array ( void  ) const
inlineprivate

Restores the array (and the local form if the vector is ghosted) to Petsc.

Definition at line 1340 of file petsc_vector.h.

References libMeshEnums::GHOSTED, libMesh::initialized(), and libMesh::libmesh_assert().

Referenced by libMesh::PetscVector< T >::add(), libMesh::PetscVector< T >::create_subvector(), libMesh::PetscVector< T >::init(), and libMesh::PetscVector< T >::operator=().

1341 {
1342  libmesh_assert (this->initialized());
1343  if(_array_is_present)
1344  {
1345  PetscErrorCode ierr=0;
1346  if(this->type() != GHOSTED)
1347  {
1348  ierr = VecRestoreArray (_vec, &_values);
1349  LIBMESH_CHKERRABORT(ierr);
1350  _values = NULL;
1351  }
1352  else
1353  {
1354  ierr = VecRestoreArray (_local_form, &_values);
1355  LIBMESH_CHKERRABORT(ierr);
1356  _values = NULL;
1357  ierr = VecGhostRestoreLocalForm (_vec,&_local_form);
1358  LIBMESH_CHKERRABORT(ierr);
1359  _local_form = NULL;
1360 #ifndef NDEBUG
1361  _local_size = 0;
1362 #endif
1363  }
1364  _array_is_present = false;
1365  }
1366 }
template<typename T >
void libMesh::PetscVector< T >::abs ( )
virtual

v = abs(v)... that is, each entry in v is replaced by its absolute value.

Implements libMesh::NumericVector< T >.

Definition at line 565 of file petsc_vector.C.

References libMeshEnums::GHOSTED.

566 {
567  this->_restore_array();
568 
569  PetscErrorCode ierr = 0;
570 
571  if(this->type() != GHOSTED)
572  {
573  ierr = VecAbs(_vec);
574  LIBMESH_CHKERRABORT(ierr);
575  }
576  else
577  {
578  Vec loc_vec;
579  ierr = VecGhostGetLocalForm (_vec,&loc_vec);
580  LIBMESH_CHKERRABORT(ierr);
581 
582  ierr = VecAbs(loc_vec);
583  LIBMESH_CHKERRABORT(ierr);
584 
585  ierr = VecGhostRestoreLocalForm (_vec,&loc_vec);
586  LIBMESH_CHKERRABORT(ierr);
587  }
588 }
template<typename T >
void libMesh::PetscVector< T >::add ( const numeric_index_type  i,
const T  value 
)
virtual

v(i) += value

Implements libMesh::NumericVector< T >.

Definition at line 191 of file petsc_vector.C.

References libMesh::ierr.

192 {
193  this->_restore_array();
194  libmesh_assert_less (i, size());
195 
197  PetscInt i_val = static_cast<PetscInt>(i);
198  PetscScalar petsc_value = static_cast<PetscScalar>(value);
199 
200  ierr = VecSetValues (_vec, 1, &i_val, &petsc_value, ADD_VALUES);
201  LIBMESH_CHKERRABORT(ierr);
202 
203  this->_is_closed = false;
204 }
template<typename T >
void libMesh::PetscVector< T >::add ( const T  s)
virtual

$U(0-LIBMESH_DIM)+=s$. Addition of s to all components. Note that s is a scalar and not a vector.

Implements libMesh::NumericVector< T >.

Definition at line 336 of file petsc_vector.C.

References libMeshEnums::GHOSTED.

337 {
338  this->_restore_array();
339 
341  PetscScalar* values;
342  const PetscScalar v = static_cast<PetscScalar>(v_in);
343 
344  if(this->type() != GHOSTED)
345  {
346  const PetscInt n = static_cast<PetscInt>(this->local_size());
347  const PetscInt fli = static_cast<PetscInt>(this->first_local_index());
348 
349  for (PetscInt i=0; i<n; i++)
350  {
351  ierr = VecGetArray (_vec, &values);
352  LIBMESH_CHKERRABORT(ierr);
353 
354  PetscInt ig = fli + i;
355 
356  PetscScalar value = (values[i] + v);
357 
358  ierr = VecRestoreArray (_vec, &values);
359  LIBMESH_CHKERRABORT(ierr);
360 
361  ierr = VecSetValues (_vec, 1, &ig, &value, INSERT_VALUES);
362  LIBMESH_CHKERRABORT(ierr);
363  }
364  }
365  else
366  {
367  /* Vectors that include ghost values require a special
368  handling. */
369  Vec loc_vec;
370  ierr = VecGhostGetLocalForm (_vec,&loc_vec);
371  LIBMESH_CHKERRABORT(ierr);
372 
373  PetscInt n=0;
374  ierr = VecGetSize(loc_vec, &n);
375  LIBMESH_CHKERRABORT(ierr);
376 
377  for (PetscInt i=0; i<n; i++)
378  {
379  ierr = VecGetArray (loc_vec, &values);
380  LIBMESH_CHKERRABORT(ierr);
381 
382  PetscScalar value = (values[i] + v);
383 
384  ierr = VecRestoreArray (loc_vec, &values);
385  LIBMESH_CHKERRABORT(ierr);
386 
387  ierr = VecSetValues (loc_vec, 1, &i, &value, INSERT_VALUES);
388  LIBMESH_CHKERRABORT(ierr);
389  }
390 
391  ierr = VecGhostRestoreLocalForm (_vec,&loc_vec);
392  LIBMESH_CHKERRABORT(ierr);
393  }
394 
395  this->_is_closed = false;
396 }
template<typename T >
void libMesh::PetscVector< T >::add ( const NumericVector< T > &  V)
virtual

$ U+=V $ . Simple vector addition, equal to the operator +=.

Implements libMesh::NumericVector< T >.

Definition at line 401 of file petsc_vector.C.

402 {
403  this->add (1., v);
404 }
template<typename T >
void libMesh::PetscVector< T >::add ( const T  a,
const NumericVector< T > &  v 
)
virtual

$ U+=a*V $ . Simple vector addition, equal to the operator +=.

Implements libMesh::NumericVector< T >.

Definition at line 409 of file petsc_vector.C.

References libMesh::PetscVector< T >::_restore_array(), libMesh::PetscVector< T >::_vec, libMeshEnums::GHOSTED, and libMesh::PetscVector< T >::size().

410 {
411  this->_restore_array();
412 
413  PetscErrorCode ierr = 0;
414  PetscScalar a = static_cast<PetscScalar>(a_in);
415 
416  // Make sure the NumericVector passed in is really a PetscVector
417  const PetscVector<T>* v = libmesh_cast_ptr<const PetscVector<T>*>(&v_in);
418  v->_restore_array();
419 
420  libmesh_assert_equal_to (this->size(), v->size());
421 
422  if(this->type() != GHOSTED)
423  {
424 #if PETSC_VERSION_LESS_THAN(2,3,0)
425  // 2.2.x & earlier style
426  ierr = VecAXPY(&a, v->_vec, _vec);
427  LIBMESH_CHKERRABORT(ierr);
428 #else
429  // 2.3.x & later style
430  ierr = VecAXPY(_vec, a, v->_vec);
431  LIBMESH_CHKERRABORT(ierr);
432 #endif
433  }
434  else
435  {
436  Vec loc_vec;
437  Vec v_loc_vec;
438  ierr = VecGhostGetLocalForm (_vec,&loc_vec);
439  LIBMESH_CHKERRABORT(ierr);
440  ierr = VecGhostGetLocalForm (v->_vec,&v_loc_vec);
441  LIBMESH_CHKERRABORT(ierr);
442 
443 #if PETSC_VERSION_LESS_THAN(2,3,0)
444  // 2.2.x & earlier style
445  ierr = VecAXPY(&a, v_loc_vec, loc_vec);
446  LIBMESH_CHKERRABORT(ierr);
447 #else
448  // 2.3.x & later style
449  ierr = VecAXPY(loc_vec, a, v_loc_vec);
450  LIBMESH_CHKERRABORT(ierr);
451 #endif
452 
453  ierr = VecGhostRestoreLocalForm (v->_vec,&v_loc_vec);
454  LIBMESH_CHKERRABORT(ierr);
455  ierr = VecGhostRestoreLocalForm (_vec,&loc_vec);
456  LIBMESH_CHKERRABORT(ierr);
457  }
458 }
template<typename T >
void libMesh::PetscVector< T >::add_vector ( const std::vector< T > &  v,
const std::vector< numeric_index_type > &  dof_indices 
)
virtual

$ U+=v $ where v is a std::vector<T> and you want to specify WHERE to add it

Implements libMesh::NumericVector< T >.

Definition at line 209 of file petsc_vector.C.

References libMesh::ierr.

211 {
212  // If we aren't adding anything just return
213  if(v.empty() || dof_indices.empty())
214  return;
215 
216  this->_restore_array();
217  libmesh_assert_equal_to (v.size(), dof_indices.size());
218 
220  const PetscInt * i_val = reinterpret_cast<const PetscInt*>(&dof_indices[0]);
221  const PetscScalar * petsc_value = static_cast<const PetscScalar*>(&v[0]);
222 
223  ierr = VecSetValues (_vec, v.size(), i_val, petsc_value, ADD_VALUES);
224  LIBMESH_CHKERRABORT(ierr);
225 
226  this->_is_closed = false;
227 }
template<typename T >
void libMesh::PetscVector< T >::add_vector ( const NumericVector< T > &  V,
const std::vector< numeric_index_type > &  dof_indices 
)
virtual

$ U+=V $ where U and V are type NumericVector<T> and you want to specify WHERE to add the NumericVector<T> V

Implements libMesh::NumericVector< T >.

Definition at line 232 of file petsc_vector.C.

References libMesh::NumericVector< T >::size().

234 {
235  libmesh_assert_equal_to (V.size(), dof_indices.size());
236 
237  for (unsigned int i=0; i<V.size(); i++)
238  this->add (dof_indices[i], V(i));
239 }
template<typename T >
void libMesh::PetscVector< T >::add_vector ( const NumericVector< T > &  V,
const SparseMatrix< T > &  A 
)
virtual

$U+=A*V$, add the product of a SparseMatrix A and a NumericVector V to this, where this=U.

Implements libMesh::NumericVector< T >.

Definition at line 244 of file petsc_vector.C.

References libMesh::PetscVector< T >::_vec, A, and libMesh::ierr.

246 {
247  this->_restore_array();
248  // Make sure the data passed in are really of Petsc types
249  const PetscVector<T>* V = libmesh_cast_ptr<const PetscVector<T>*>(&V_in);
250  const PetscMatrix<T>* A = libmesh_cast_ptr<const PetscMatrix<T>*>(&A_in);
251 
253 
254  A->close();
255 
256  // The const_cast<> is not elegant, but it is required since PETSc
257  // is not const-correct.
258  ierr = MatMultAdd(const_cast<PetscMatrix<T>*>(A)->mat(), V->_vec, _vec, _vec);
259  LIBMESH_CHKERRABORT(ierr);
260 }
template<typename T >
void libMesh::PetscVector< T >::add_vector ( const DenseVector< T > &  V,
const std::vector< numeric_index_type > &  dof_indices 
)
virtual

$U+=V $ where U and V are type DenseVector<T> and you want to specify WHERE to add the DenseVector<T> V

Implements libMesh::NumericVector< T >.

Definition at line 265 of file petsc_vector.C.

References libMesh::DenseVector< T >::size().

267 {
268  libmesh_assert_equal_to (V.size(), dof_indices.size());
269 
270  for (unsigned int i=0; i<V.size(); i++)
271  this->add (dof_indices[i], V(i));
272 }
template<typename T>
void libMesh::NumericVector< T >::add_vector ( const NumericVector< T > &  v,
const ShellMatrix< T > &  a 
)
inherited

$U+=A*V$, add the product of a ShellMatrix A and a NumericVector V to this, where this=U.

Definition at line 380 of file numeric_vector.C.

References libMesh::ShellMatrix< T >::vector_mult_add().

382 {
383  a.vector_mult_add(*this,v);
384 }
template<typename T >
void libMesh::PetscVector< T >::add_vector_conjugate_transpose ( const NumericVector< T > &  V,
const SparseMatrix< T > &  A_trans 
)

$U+=A^H*V$, add the product of the conjugate-transpose of SparseMatrix A_trans and a NumericVector V to this, where this=U.

Definition at line 297 of file petsc_vector.C.

References libMesh::out.

299 {
300 
301  libMesh::out << "MatMultHermitianTranspose was introduced in PETSc 3.1.0,"
302  << "No one has made it backwards compatible with older "
303  << "versions of PETSc so far." << std::endl;
304  libmesh_error();
305 }
template<typename T >
void libMesh::PetscVector< T >::add_vector_transpose ( const NumericVector< T > &  V,
const SparseMatrix< T > &  A_trans 
)
virtual

$U+=A^T*V$, add the product of the transpose of SparseMatrix A_trans and a NumericVector V to this, where this=U.

Implements libMesh::NumericVector< T >.

Definition at line 276 of file petsc_vector.C.

References libMesh::PetscVector< T >::_vec, A, and libMesh::ierr.

278 {
279  this->_restore_array();
280  // Make sure the data passed in are really of Petsc types
281  const PetscVector<T>* V = libmesh_cast_ptr<const PetscVector<T>*>(&V_in);
282  const PetscMatrix<T>* A = libmesh_cast_ptr<const PetscMatrix<T>*>(&A_in);
283 
285 
286  A->close();
287 
288  // The const_cast<> is not elegant, but it is required since PETSc
289  // is not const-correct.
290  ierr = MatMultTransposeAdd(const_cast<PetscMatrix<T>*>(A)->mat(), V->_vec, _vec, _vec);
291  LIBMESH_CHKERRABORT(ierr);
292 }
template<typename T >
AutoPtr< NumericVector< T > > libMesh::NumericVector< T >::build ( const Parallel::Communicator comm,
const SolverPackage  solver_package = libMesh::default_solver_package() 
)
staticinherited

Builds a NumericVector on the processors in communicator comm using the linear solver package specified by solver_package

Definition at line 46 of file numeric_vector.C.

References libMeshEnums::AUTOMATIC, libMesh::EIGEN_SOLVERS, libMesh::LASPACK_SOLVERS, libMeshEnums::PETSC_SOLVERS, and libMesh::TRILINOS_SOLVERS.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::System::add_vector(), libMesh::NumericVector< T >::build(), libMesh::EquationSystems::build_solution_vector(), libMesh::System::calculate_norm(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::for(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::EquationSystems::get_solution(), and libMesh::System::project_vector().

47 {
48  // Build the appropriate vector
49  switch (solver_package)
50  {
51 
52 
53 #ifdef LIBMESH_HAVE_LASPACK
54  case LASPACK_SOLVERS:
55  {
56  AutoPtr<NumericVector<T> > ap(new LaspackVector<T>(comm, AUTOMATIC));
57  return ap;
58  }
59 #endif
60 
61 
62 #ifdef LIBMESH_HAVE_PETSC
63  case PETSC_SOLVERS:
64  {
65  AutoPtr<NumericVector<T> > ap(new PetscVector<T>(comm, AUTOMATIC));
66  return ap;
67  }
68 #endif
69 
70 
71 #ifdef LIBMESH_HAVE_TRILINOS
72  case TRILINOS_SOLVERS:
73  {
74  AutoPtr<NumericVector<T> > ap(new EpetraVector<T>(comm, AUTOMATIC));
75  return ap;
76  }
77 #endif
78 
79 
80 #ifdef LIBMESH_HAVE_EIGEN
81  case EIGEN_SOLVERS:
82  {
83  AutoPtr<NumericVector<T> > ap(new EigenSparseVector<T>(comm, AUTOMATIC));
84  return ap;
85  }
86 #endif
87 
88 
89  default:
90  AutoPtr<NumericVector<T> > ap(new DistributedVector<T>(comm, AUTOMATIC));
91  return ap;
92 
93  }
94 
95  AutoPtr<NumericVector<T> > ap(NULL);
96  return ap;
97 }
template<typename T >
AutoPtr< NumericVector< T > > libMesh::NumericVector< T >::build ( const SolverPackage  solver_package = libMesh::default_solver_package())
staticinherited

Builds a NumericVector on the processors in communicator CommWorld using the linear solver package specified by solver_package. Deprecated.

Definition at line 103 of file numeric_vector.C.

References libMesh::NumericVector< T >::build(), and libMesh::CommWorld.

104 {
105  libmesh_deprecated();
106  return NumericVector<T>::build(CommWorld, solver_package);
107 }
template<typename T >
void libMesh::PetscVector< T >::clear ( )
inlinevirtual
Returns
the PetscVector<T> to a pristine state.

Reimplemented from libMesh::NumericVector< T >.

Definition at line 988 of file petsc_vector.h.

References libMesh::libMeshPrivateData::_is_initialized, and libMesh::initialized().

989 {
990  if (this->initialized())
991  this->_restore_array();
992 
993  if ((this->initialized()) && (this->_destroy_vec_on_exit))
994  {
996 
997  ierr = LibMeshVecDestroy(&_vec);
998  LIBMESH_CHKERRABORT(ierr);
999  }
1000 
1001  this->_is_closed = this->_is_initialized = false;
1002 
1003  _global_to_local_map.clear();
1004 }
template<typename T >
AutoPtr< NumericVector< T > > libMesh::PetscVector< T >::clone ( ) const
inlinevirtual

Creates a copy of this vector and returns it in an AutoPtr.

Implements libMesh::NumericVector< T >.

Definition at line 1071 of file petsc_vector.h.

References libMesh::comm.

1072 {
1073  AutoPtr<NumericVector<T> > cloned_vector
1074  (new PetscVector<T>(this->comm(), this->type()));
1075 
1076  cloned_vector->init(*this, true);
1077 
1078  *cloned_vector = *this;
1079 
1080  return cloned_vector;
1081 }
template<typename T >
void libMesh::PetscVector< T >::close ( )
inlinevirtual

Call the assemble functions

Implements libMesh::NumericVector< T >.

Definition at line 962 of file petsc_vector.h.

References libMeshEnums::GHOSTED.

Referenced by libMesh::__libmesh_petsc_diff_solver_monitor(), libMesh::SlepcEigenSolver< T >::get_eigenpair(), libMesh::PetscVector< T >::localize(), and libMesh::PetscLinearSolver< T >::solve().

963 {
964  this->_restore_array();
965 
967 
968  ierr = VecAssemblyBegin(_vec);
969  LIBMESH_CHKERRABORT(ierr);
970  ierr = VecAssemblyEnd(_vec);
971  LIBMESH_CHKERRABORT(ierr);
972 
973  if(this->type() == GHOSTED)
974  {
975  ierr = VecGhostUpdateBegin(_vec,INSERT_VALUES,SCATTER_FORWARD);
976  LIBMESH_CHKERRABORT(ierr);
977  ierr = VecGhostUpdateEnd(_vec,INSERT_VALUES,SCATTER_FORWARD);
978  LIBMESH_CHKERRABORT(ierr);
979  }
980 
981  this->_is_closed = true;
982 }
template<typename T>
virtual bool libMesh::NumericVector< T >::closed ( ) const
inlinevirtualinherited
const Parallel::Communicator& libMesh::ParallelObject::comm ( ) const
inlineinherited
Returns
a reference to the Parallel::Communicator object used by this mesh.

Definition at line 86 of file parallel_object.h.

References libMesh::ParallelObject::_communicator.

Referenced by libMesh::__libmesh_petsc_diff_solver_jacobian(), libMesh::__libmesh_petsc_diff_solver_monitor(), libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::__libmesh_petsc_snes_jacobian(), libMesh::__libmesh_petsc_snes_residual(), libMesh::MeshRefinement::_coarsen_elements(), libMesh::ExactSolution::_compute_error(), libMesh::MetisPartitioner::_do_partition(), libMesh::ParmetisPartitioner::_do_repartition(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_get_diagonal(), libMesh::SlepcEigenSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult(), libMesh::PetscLinearSolver< T >::_petsc_shell_matrix_mult_add(), libMesh::EquationSystems::_read_impl(), libMesh::MeshRefinement::_refine_elements(), libMesh::ParallelMesh::add_elem(), libMesh::ImplicitSystem::add_matrix(), libMesh::ParallelMesh::add_node(), libMesh::System::add_vector(), libMesh::UnstructuredMesh::all_second_order(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::FEMSystem::assemble_qoi(), libMesh::MeshCommunication::assign_global_indices(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::DofMap::attach_matrix(), libMesh::MeshTools::bounding_box(), libMesh::System::calculate_norm(), libMesh::MeshRefinement::coarsen_elements(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Problem_Interface::computeF(), libMesh::Problem_Interface::computeJacobian(), libMesh::Problem_Interface::computePreconditioner(), libMesh::MeshTools::correct_node_proc_ids(), libMesh::MeshCommunication::delete_remote_elements(), libMesh::DofMap::distribute_dofs(), DMlibMeshFunction(), DMlibMeshJacobian(), DMLibMeshSetSystem(), DMVariableBounds_libMesh(), libMesh::MeshRefinement::eliminate_unrefined_patches(), libMesh::WeightedPatchRecoveryErrorEstimator::estimate_error(), libMesh::PatchRecoveryErrorEstimator::estimate_error(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::AdjointRefinementEstimator::estimate_error(), libMesh::MeshRefinement::flag_elements_by_elem_fraction(), libMesh::MeshRefinement::flag_elements_by_error_fraction(), libMesh::MeshRefinement::flag_elements_by_nelem_target(), libMesh::for(), libMesh::CondensedEigenSystem::get_eigenpair(), libMesh::ImplicitSystem::get_linear_solver(), libMesh::LocationMap< T >::init(), libMesh::PetscDiffSolver::init(), libMesh::TimeSolver::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::EigenSystem::init_data(), libMesh::EigenSystem::init_matrices(), libMesh::ParmetisPartitioner::initialize(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::ParallelMesh::libmesh_assert_valid_parallel_flags(), libMesh::MeshTools::libmesh_assert_valid_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshRefinement::limit_level_mismatch_at_edge(), libMesh::MeshRefinement::limit_level_mismatch_at_node(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshCommunication::make_elems_parallel_consistent(), libMesh::MeshRefinement::make_flags_parallel_consistent(), libMesh::MeshCommunication::make_node_ids_parallel_consistent(), libMesh::MeshCommunication::make_node_proc_ids_parallel_consistent(), libMesh::MeshCommunication::make_nodes_parallel_consistent(), libMesh::MeshRefinement::make_refinement_compatible(), libMesh::FEMSystem::mesh_position_set(), libMesh::MeshSerializer::MeshSerializer(), libMesh::ParallelMesh::n_active_elem(), libMesh::MeshTools::n_active_levels(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::CondensedEigenSystem::n_global_non_condensed_dofs(), libMesh::MeshTools::n_levels(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::MeshTools::n_p_levels(), libMesh::ParallelMesh::parallel_max_elem_id(), libMesh::ParallelMesh::parallel_max_node_id(), libMesh::ParallelMesh::parallel_n_elem(), libMesh::ParallelMesh::parallel_n_nodes(), libMesh::Partitioner::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshBase::prepare_for_use(), libMesh::System::project_vector(), libMesh::Nemesis_IO::read(), libMesh::XdrIO::read(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::MeshBase::recalculate_n_partitions(), libMesh::MeshRefinement::refine_and_coarsen_elements(), libMesh::MeshRefinement::refine_elements(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::MeshBase::subdomain_ids(), libMesh::BoundaryInfo::sync(), libMesh::Parallel::sync_element_data_by_parent_id(), libMesh::MeshRefinement::test_level_one(), libMesh::MeshRefinement::test_unflagged(), libMesh::MeshTools::total_weight(), libMesh::CheckpointIO::write(), libMesh::XdrIO::write(), libMesh::UnstructuredMesh::write(), libMesh::LegacyXdrIO::write_mesh(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), and libMesh::DivaIO::write_stream().

87  { return _communicator; }
template<typename T>
int libMesh::NumericVector< T >::compare ( const NumericVector< T > &  other_vector,
const Real  threshold = TOLERANCE 
) const
virtualinherited
Returns
-1 when this is equivalent to other_vector, up to the given threshold. When differences occur, the return value contains the first index i where the difference (a[i]-b[i]) exceeded the threshold. When no threshold is given, the libMesh TOLERANCE is used.

Definition at line 112 of file numeric_vector.C.

References std::abs(), libMesh::comm, libMesh::NumericVector< T >::first_local_index(), libMesh::NumericVector< T >::initialized(), libMesh::initialized(), libMesh::NumericVector< T >::last_local_index(), libMesh::libmesh_assert(), and std::max().

114 {
115  libmesh_assert (this->initialized());
116  libmesh_assert (other_vector.initialized());
117  libmesh_assert_equal_to (this->first_local_index(), other_vector.first_local_index());
118  libmesh_assert_equal_to (this->last_local_index(), other_vector.last_local_index());
119 
120  int first_different_i = std::numeric_limits<int>::max();
122 
123  do
124  {
125  if ( std::abs( (*this)(i) - other_vector(i) ) > threshold )
126  first_different_i = i;
127  else
128  i++;
129  }
130  while (first_different_i==std::numeric_limits<int>::max()
131  && i<last_local_index());
132 
133  // Find the correct first differing index in parallel
134  this->comm().min(first_different_i);
135 
136  if (first_different_i == std::numeric_limits<int>::max())
137  return -1;
138 
139  return first_different_i;
140 }
template<typename T >
void libMesh::PetscVector< T >::conjugate ( )
virtual

Replace each entry v_i = real(v_i) + imag(v_i) of this vector by its complex conjugate, real(v_i) - imag(v_i)

Implements libMesh::NumericVector< T >.

Definition at line 179 of file petsc_vector.C.

References libMesh::ierr.

180 {
181  PetscErrorCode ierr = 0;
182 
183  // We just call the PETSc VecConjugate
184  ierr = VecConjugate(_vec);
185  LIBMESH_CHKERRABORT(ierr);
186 }
template<typename T >
void libMesh::PetscVector< T >::create_subvector ( NumericVector< T > &  subvector,
const std::vector< numeric_index_type > &  rows 
) const
virtual

Creates a "subvector" from this vector using the rows indices of the "rows" array.

Reimplemented from libMesh::NumericVector< T >.

Definition at line 1391 of file petsc_vector.C.

References libMesh::NumericVector< T >::_is_initialized, libMesh::PetscVector< T >::_restore_array(), libMesh::PetscVector< T >::_vec, libMesh::comm, libMesh::MeshTools::Generation::Private::idx(), libMesh::NumericVector< T >::initialized(), libMesh::Utility::iota(), and PETSC_USE_POINTER.

1393 {
1394  this->_restore_array();
1395 
1396  // PETSc data structures
1397  IS parent_is, subvector_is;
1398  VecScatter scatter;
1399  PetscErrorCode ierr = 0;
1400 
1401  // Make sure the passed in subvector is really a PetscVector
1402  PetscVector<T>* petsc_subvector = libmesh_cast_ptr<PetscVector<T>*>(&subvector);
1403 
1404  // If the petsc_subvector is already initialized, we assume that the
1405  // user has already allocated the *correct* amount of space for it.
1406  // If not, we use the appropriate PETSc routines to initialize it.
1407  if (!petsc_subvector->initialized())
1408  {
1409  // Initialize the petsc_subvector to have enough space to hold
1410  // the entries which will be scattered into it. Note: such an
1411  // init() function (where we let PETSc decide the number of local
1412  // entries) is not currently offered by the PetscVector
1413  // class. Should we differentiate here between sequential and
1414  // parallel vector creation based on this->n_processors() ?
1415  ierr = VecCreateMPI(this->comm().get(),
1416  PETSC_DECIDE, // n_local
1417  rows.size(), // n_global
1418  &(petsc_subvector->_vec)); LIBMESH_CHKERRABORT(ierr);
1419 
1420  ierr = VecSetFromOptions (petsc_subvector->_vec); LIBMESH_CHKERRABORT(ierr);
1421 
1422  // Mark the subvector as initialized
1423  petsc_subvector->_is_initialized = true;
1424  }
1425  else
1426  {
1427  petsc_subvector->_restore_array();
1428  }
1429 
1430  // Use iota to fill an array with entries [0,1,2,3,4,...rows.size()]
1431  std::vector<PetscInt> idx(rows.size());
1432  Utility::iota (idx.begin(), idx.end(), 0);
1433 
1434  // Construct index sets
1435  ierr = ISCreateLibMesh(this->comm().get(),
1436  rows.size(),
1437  (PetscInt*) &rows[0],
1439  &parent_is); LIBMESH_CHKERRABORT(ierr);
1440 
1441  ierr = ISCreateLibMesh(this->comm().get(),
1442  rows.size(),
1443  (PetscInt*) &idx[0],
1445  &subvector_is); LIBMESH_CHKERRABORT(ierr);
1446 
1447  // Construct the scatter object
1448  ierr = VecScatterCreate(this->_vec,
1449  parent_is,
1450  petsc_subvector->_vec,
1451  subvector_is,
1452  &scatter); LIBMESH_CHKERRABORT(ierr);
1453 
1454  // Actually perform the scatter
1455 #if PETSC_VERSION_LESS_THAN(2,3,3)
1456  ierr = VecScatterBegin(this->_vec,
1457  petsc_subvector->_vec,
1458  INSERT_VALUES,
1459  SCATTER_FORWARD,
1460  scatter); LIBMESH_CHKERRABORT(ierr);
1461 
1462  ierr = VecScatterEnd(this->_vec,
1463  petsc_subvector->_vec,
1464  INSERT_VALUES,
1465  SCATTER_FORWARD,
1466  scatter); LIBMESH_CHKERRABORT(ierr);
1467 #else
1468  // API argument order change in PETSc 2.3.3
1469  ierr = VecScatterBegin(scatter,
1470  this->_vec,
1471  petsc_subvector->_vec,
1472  INSERT_VALUES,
1473  SCATTER_FORWARD); LIBMESH_CHKERRABORT(ierr);
1474 
1475  ierr = VecScatterEnd(scatter,
1476  this->_vec,
1477  petsc_subvector->_vec,
1478  INSERT_VALUES,
1479  SCATTER_FORWARD); LIBMESH_CHKERRABORT(ierr);
1480 #endif
1481 
1482  // Clean up
1483  ierr = LibMeshISDestroy(&parent_is); LIBMESH_CHKERRABORT(ierr);
1484  ierr = LibMeshISDestroy(&subvector_is); LIBMESH_CHKERRABORT(ierr);
1485  ierr = LibMeshVecScatterDestroy(&scatter); LIBMESH_CHKERRABORT(ierr);
1486 
1487 }
void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited

Definition at line 106 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

107 {
108  _enable_print_counter = false;
109  return;
110 }
template<typename T >
T libMesh::PetscVector< T >::dot ( const NumericVector< T > &  V) const
virtual

Computes the dot product, p = U.V. Use complex-conjugate of V in the complex-valued case.

Implements libMesh::NumericVector< T >.

Definition at line 591 of file petsc_vector.C.

References libMesh::PetscVector< T >::_vec.

592 {
593  this->_restore_array();
594 
595  // Error flag
596  PetscErrorCode ierr = 0;
597 
598  // Return value
599  PetscScalar value=0.;
600 
601  // Make sure the NumericVector passed in is really a PetscVector
602  const PetscVector<T>* v = libmesh_cast_ptr<const PetscVector<T>*>(&V);
603 
604  // 2.3.x (at least) style. Untested for previous versions.
605  ierr = VecDot(this->_vec, v->_vec, &value);
606  LIBMESH_CHKERRABORT(ierr);
607 
608  return static_cast<T>(value);
609 }
template<typename T>
virtual T libMesh::NumericVector< T >::el ( const numeric_index_type  i) const
inlinevirtualinherited
Returns
the element U(i)

Definition at line 342 of file numeric_vector.h.

342 { return (*this)(i); }
void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

Methods to enable/disable the reference counter output from print_info()

Definition at line 100 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

101 {
102  _enable_print_counter = true;
103  return;
104 }
template<typename T >
numeric_index_type libMesh::PetscVector< T >::first_local_index ( ) const
inlinevirtual
Returns
the index of the first vector element actually stored on this processor

Implements libMesh::NumericVector< T >.

Definition at line 1124 of file petsc_vector.h.

References libMesh::initialized(), and libMesh::libmesh_assert().

1125 {
1126  libmesh_assert (this->initialized());
1127 
1128  PetscErrorCode ierr=0;
1129  PetscInt petsc_first=0, petsc_last=0;
1130 
1131  ierr = VecGetOwnershipRange (_vec, &petsc_first, &petsc_last);
1132  LIBMESH_CHKERRABORT(ierr);
1133 
1134  return static_cast<numeric_index_type>(petsc_first);
1135 }
template<typename T >
void libMesh::PetscVector< T >::get ( const std::vector< numeric_index_type > &  index,
std::vector< T > &  values 
) const
inlinevirtual

Access multiple components at once. Overloaded method that should be faster (probably much faster) than calling operator() individually for each index.

Reimplemented from libMesh::NumericVector< T >.

Definition at line 1228 of file petsc_vector.h.

References libMeshEnums::GHOSTED.

1229 {
1230  this->_get_array();
1231 
1232  const std::size_t num = index.size();
1233  values.resize(num);
1234 
1235  for(std::size_t i=0; i<num; i++)
1236  {
1237  const numeric_index_type local_index = this->map_global_to_local_index(index[i]);
1238 #ifndef NDEBUG
1239  if(this->type() == GHOSTED)
1240  {
1241  libmesh_assert_less (local_index, _local_size);
1242  }
1243 #endif
1244  values[i] = static_cast<T>(_values[local_index]);
1245  }
1246 }
std::string libMesh::ReferenceCounter::get_info ( )
staticinherited

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

References libMesh::ReferenceCounter::_counts, and libMesh::Quality::name().

Referenced by libMesh::ReferenceCounter::print_info().

48 {
49 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
50 
51  std::ostringstream oss;
52 
53  oss << '\n'
54  << " ---------------------------------------------------------------------------- \n"
55  << "| Reference count information |\n"
56  << " ---------------------------------------------------------------------------- \n";
57 
58  for (Counts::iterator it = _counts.begin();
59  it != _counts.end(); ++it)
60  {
61  const std::string name(it->first);
62  const unsigned int creations = it->second.first;
63  const unsigned int destructions = it->second.second;
64 
65  oss << "| " << name << " reference count information:\n"
66  << "| Creations: " << creations << '\n'
67  << "| Destructions: " << destructions << '\n';
68  }
69 
70  oss << " ---------------------------------------------------------------------------- \n";
71 
72  return oss.str();
73 
74 #else
75 
76  return "";
77 
78 #endif
79 }
template<typename T>
int libMesh::NumericVector< T >::global_relative_compare ( const NumericVector< T > &  other_vector,
const Real  threshold = TOLERANCE 
) const
virtualinherited
Returns
-1 when this is equivalent to other_vector, up to the given local relative threshold. When differences occur, the return value contains the first index where the difference (a[i]-b[i])/max_j(a[j],b[j]) exceeded the threshold. When no threshold is given, the libMesh TOLERANCE is used.

Definition at line 177 of file numeric_vector.C.

References std::abs(), libMesh::comm, libMesh::NumericVector< T >::first_local_index(), libMesh::NumericVector< T >::initialized(), libMesh::initialized(), libMesh::NumericVector< T >::last_local_index(), libMesh::libmesh_assert(), libMesh::NumericVector< T >::linfty_norm(), std::max(), and libMesh::Real.

179 {
180  libmesh_assert (this->initialized());
181  libmesh_assert (other_vector.initialized());
182  libmesh_assert_equal_to (this->first_local_index(), other_vector.first_local_index());
183  libmesh_assert_equal_to (this->last_local_index(), other_vector.last_local_index());
184 
185  int first_different_i = std::numeric_limits<int>::max();
187 
188  const Real my_norm = this->linfty_norm();
189  const Real other_norm = other_vector.linfty_norm();
190  const Real abs_threshold = std::max(my_norm, other_norm) * threshold;
191 
192  do
193  {
194  if ( std::abs( (*this)(i) - other_vector(i) ) > abs_threshold )
195  first_different_i = i;
196  else
197  i++;
198  }
199  while (first_different_i==std::numeric_limits<int>::max()
200  && i<last_local_index());
201 
202  // Find the correct first differing index in parallel
203  this->comm().min(first_different_i);
204 
205  if (first_different_i == std::numeric_limits<int>::max())
206  return -1;
207 
208  return first_different_i;
209 }
void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name)
inlineprotectedinherited

Increments the construction counter. Should be called in the constructor of any derived class that will be reference counted.

Definition at line 163 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject().

164 {
165  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
166  std::pair<unsigned int, unsigned int>& p = _counts[name];
167 
168  p.first++;
169 }
void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name)
inlineprotectedinherited

Increments the destruction counter. Should be called in the destructor of any derived class that will be reference counted.

Definition at line 176 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

177 {
178  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
179  std::pair<unsigned int, unsigned int>& p = _counts[name];
180 
181  p.second++;
182 }
template<typename T >
T libMesh::PetscVector< T >::indefinite_dot ( const NumericVector< T > &  V) const
virtual

Computes the dot product, p = U.V. Do not use complex-conjugate of V in the complex-valued case.

Definition at line 612 of file petsc_vector.C.

References libMesh::PetscVector< T >::_vec.

613 {
614  this->_restore_array();
615 
616  // Error flag
617  PetscErrorCode ierr = 0;
618 
619  // Return value
620  PetscScalar value=0.;
621 
622  // Make sure the NumericVector passed in is really a PetscVector
623  const PetscVector<T>* v = libmesh_cast_ptr<const PetscVector<T>*>(&V);
624 
625  // 2.3.x (at least) style. Untested for previous versions.
626  ierr = VecTDot(this->_vec, v->_vec, &value);
627  LIBMESH_CHKERRABORT(ierr);
628 
629  return static_cast<T>(value);
630 }
template<typename T >
void libMesh::PetscVector< T >::init ( const numeric_index_type  N,
const numeric_index_type  n_local,
const bool  fast = false,
const ParallelType  type = AUTOMATIC 
)
inlinevirtual

Change the dimension of the vector to N. The reserved memory for this vector remains unchanged if possible, to make things faster, but this may waste some memory, so take this in the back of your head. However, if N==0 all memory is freed, i.e. if you want to resize the vector and release the memory not needed, you have to first call init(0) and then init(N). This cited behaviour is analogous to that of the STL containers.

On fast==false, the vector is filled by zeros.

Implements libMesh::NumericVector< T >.

Definition at line 782 of file petsc_vector.h.

References libMesh::libMeshPrivateData::_is_initialized, libMeshEnums::AUTOMATIC, libMesh::CHKERRABORT(), libMesh::comm, libMesh::initialized(), libMesh::libmesh_assert(), libMesh::n_local, libMeshEnums::PARALLEL, libMeshEnums::SERIAL, and libMesh::zero.

Referenced by libMesh::PetscVector< T >::localize(), and libMesh::PetscVector< T >::PetscVector().

786 {
788  PetscInt petsc_n=static_cast<PetscInt>(n);
789  PetscInt petsc_n_local=static_cast<PetscInt>(n_local);
790 
791 
792  // Clear initialized vectors
793  if (this->initialized())
794  this->clear();
795 
796  if (ptype == AUTOMATIC)
797  {
798  if (n == n_local)
799  this->_type = SERIAL;
800  else
801  this->_type = PARALLEL;
802  }
803  else
804  this->_type = ptype;
805 
806  libmesh_assert ((this->_type==SERIAL && n==n_local) ||
807  this->_type==PARALLEL);
808 
809  // create a sequential vector if on only 1 processor
810  if (this->_type == SERIAL)
811  {
812  ierr = VecCreateSeq (PETSC_COMM_SELF, petsc_n, &_vec);
813  CHKERRABORT(PETSC_COMM_SELF,ierr);
814 
815  ierr = VecSetFromOptions (_vec);
816  CHKERRABORT(PETSC_COMM_SELF,ierr);
817  }
818  // otherwise create an MPI-enabled vector
819  else if (this->_type == PARALLEL)
820  {
821 #ifdef LIBMESH_HAVE_MPI
822  libmesh_assert_less_equal (n_local, n);
823  ierr = VecCreateMPI (this->comm().get(), petsc_n_local, petsc_n,
824  &_vec);
825  LIBMESH_CHKERRABORT(ierr);
826 #else
827  libmesh_assert_equal_to (n_local, n);
828  ierr = VecCreateSeq (PETSC_COMM_SELF, petsc_n, &_vec);
829  CHKERRABORT(PETSC_COMM_SELF,ierr);
830 #endif
831 
832  ierr = VecSetFromOptions (_vec);
833  LIBMESH_CHKERRABORT(ierr);
834  }
835  else
836  libmesh_error();
837 
838  this->_is_initialized = true;
839  this->_is_closed = true;
840 
841 
842  if (fast == false)
843  this->zero ();
844 }
template<typename T >
void libMesh::PetscVector< T >::init ( const numeric_index_type  N,
const bool  fast = false,
const ParallelType  type = AUTOMATIC 
)
inlinevirtual

call init with n_local = N,

Implements libMesh::NumericVector< T >.

Definition at line 850 of file petsc_vector.h.

References libMesh::TriangleWrapper::init().

853 {
854  this->init(n,n,fast,ptype);
855 }
template<typename T>
virtual void libMesh::PetscVector< T >::init ( const numeric_index_type  ,
const numeric_index_type  ,
const std::vector< numeric_index_type > &  ,
const bool  = false,
const ParallelType  = AUTOMATIC 
)
virtual

Create a vector that holds tha local indices plus those specified in the ghost argument.

Implements libMesh::NumericVector< T >.

template<typename T >
void libMesh::PetscVector< T >::init ( const NumericVector< T > &  other,
const bool  fast = false 
)
inlinevirtual

Creates a vector that has the same dimension and storage type as other, including ghost dofs.

Implements libMesh::NumericVector< T >.

Definition at line 921 of file petsc_vector.h.

References libMesh::PetscVector< T >::_global_to_local_map, libMesh::libMeshPrivateData::_is_initialized, libMesh::PetscVector< T >::_restore_array(), libMesh::NumericVector< T >::_type, libMesh::PetscVector< T >::_vec, libMesh::NumericVector< T >::initialized(), libMesh::initialized(), libMesh::PetscVector< T >::size(), and libMesh::zero.

923 {
924  // Clear initialized vectors
925  if (this->initialized())
926  this->clear();
927 
928  const PetscVector<T>& v = libmesh_cast_ref<const PetscVector<T>&>(other);
929 
930  // Other vector should restore array.
931  if(v.initialized())
932  {
933  v._restore_array();
934  }
935 
936  this->_global_to_local_map = v._global_to_local_map;
937 
938  // Even if we're initializeing sizes based on an uninitialized or
939  // unclosed vector, *this* vector is being initialized now and is
940  // initially closed.
941  this->_is_closed = true; // v._is_closed;
942  this->_is_initialized = true; // v._is_initialized;
943 
944  this->_type = v._type;
945 
946  if (v.size() != 0)
947  {
948  PetscErrorCode ierr = 0;
949 
950  ierr = VecDuplicate (v._vec, &this->_vec);
951  LIBMESH_CHKERRABORT(ierr);
952  }
953 
954  if (fast == false)
955  this->zero ();
956 }
template<typename T>
virtual bool libMesh::NumericVector< T >::initialized ( ) const
inlinevirtualinherited
template<typename T >
void libMesh::PetscVector< T >::insert ( const std::vector< T > &  v,
const std::vector< numeric_index_type > &  dof_indices 
)
virtual

$ U=v $ where v is a std::vector<T> and you want to specify WHERE to insert it

Implements libMesh::NumericVector< T >.

Definition at line 463 of file petsc_vector.C.

465 {
466  libmesh_assert_equal_to (v.size(), dof_indices.size());
467 
468  for (unsigned int i=0; i<v.size(); i++)
469  this->set (dof_indices[i], v[i]);
470 }
template<typename T >
void libMesh::PetscVector< T >::insert ( const NumericVector< T > &  V,
const std::vector< numeric_index_type > &  dof_indices 
)
virtual

$U=V$, where U and V are type NumericVector<T> and you want to specify WHERE to insert the NumericVector<T> V

Implements libMesh::NumericVector< T >.

Definition at line 475 of file petsc_vector.C.

References libMesh::NumericVector< T >::size().

477 {
478  libmesh_assert_equal_to (V.size(), dof_indices.size());
479 
480  for (unsigned int i=0; i<V.size(); i++)
481  this->set (dof_indices[i], V(i));
482 }
template<typename T >
void libMesh::PetscVector< T >::insert ( const DenseVector< T > &  V,
const std::vector< numeric_index_type > &  dof_indices 
)
virtual

$ U=V $ where V is type DenseVector<T> and you want to specify WHERE to insert it

Implements libMesh::NumericVector< T >.

Definition at line 487 of file petsc_vector.C.

References libMesh::DenseVector< T >::size().

489 {
490  libmesh_assert_equal_to (V.size(), dof_indices.size());
491 
492  for (unsigned int i=0; i<V.size(); i++)
493  this->set (dof_indices[i], V(i));
494 }
template<typename T >
void libMesh::PetscVector< T >::insert ( const DenseSubVector< T > &  V,
const std::vector< numeric_index_type > &  dof_indices 
)
virtual

$ U=V $ where V is type DenseSubVector<T> and you want to specify WHERE to insert it

Implements libMesh::NumericVector< T >.

Definition at line 499 of file petsc_vector.C.

References libMesh::DenseVectorBase< T >::size().

501 {
502  libmesh_assert_equal_to (V.size(), dof_indices.size());
503 
504  for (unsigned int i=0; i<V.size(); i++)
505  this->set (dof_indices[i], V(i));
506 }
template<typename T >
Real libMesh::PetscVector< T >::l1_norm ( ) const
virtual
Returns
the $l_1$-norm of the vector, i.e. the sum of the absolute values.

Implements libMesh::NumericVector< T >.

Definition at line 68 of file petsc_vector.C.

References libMesh::closed(), libMesh::ierr, libMesh::libmesh_assert(), and libMesh::Real.

69 {
70  this->_restore_array();
71  libmesh_assert(this->closed());
72 
74  PetscReal value=0.;
75 
76  ierr = VecNorm (_vec, NORM_1, &value);
77  LIBMESH_CHKERRABORT(ierr);
78 
79  return static_cast<Real>(value);
80 }
template<typename T >
Real libMesh::PetscVector< T >::l2_norm ( ) const
virtual
Returns
the $l_2$-norm of the vector, i.e. the square root of the sum of the squares of the elements.

Implements libMesh::NumericVector< T >.

Definition at line 85 of file petsc_vector.C.

References libMesh::closed(), libMesh::ierr, libMesh::libmesh_assert(), and libMesh::Real.

86 {
87  this->_restore_array();
88  libmesh_assert(this->closed());
89 
91  PetscReal value=0.;
92 
93  ierr = VecNorm (_vec, NORM_2, &value);
94  LIBMESH_CHKERRABORT(ierr);
95 
96  return static_cast<Real>(value);
97 }
template<typename T >
numeric_index_type libMesh::PetscVector< T >::last_local_index ( ) const
inlinevirtual
Returns
the index of the last vector element actually stored on this processor

Implements libMesh::NumericVector< T >.

Definition at line 1141 of file petsc_vector.h.

References libMesh::initialized(), and libMesh::libmesh_assert().

1142 {
1143  libmesh_assert (this->initialized());
1144 
1145  PetscErrorCode ierr=0;
1146  PetscInt petsc_first=0, petsc_last=0;
1147 
1148  ierr = VecGetOwnershipRange (_vec, &petsc_first, &petsc_last);
1149  LIBMESH_CHKERRABORT(ierr);
1150 
1151  return static_cast<numeric_index_type>(petsc_last);
1152 }
template<typename T >
Real libMesh::PetscVector< T >::linfty_norm ( ) const
virtual
Returns
the maximum absolute value of the elements of this vector, which is the $l_\infty$-norm of a vector.

Implements libMesh::NumericVector< T >.

Definition at line 103 of file petsc_vector.C.

References libMesh::closed(), libMesh::ierr, libMesh::libmesh_assert(), and libMesh::Real.

104 {
105  this->_restore_array();
106  libmesh_assert(this->closed());
107 
109  PetscReal value=0.;
110 
111  ierr = VecNorm (_vec, NORM_INFINITY, &value);
112  LIBMESH_CHKERRABORT(ierr);
113 
114  return static_cast<Real>(value);
115 }
template<typename T>
int libMesh::NumericVector< T >::local_relative_compare ( const NumericVector< T > &  other_vector,
const Real  threshold = TOLERANCE 
) const
virtualinherited
Returns
-1 when this is equivalent to other_vector, up to the given local relative threshold. When differences occur, the return value contains the first index where the difference (a[i]-b[i])/max(a[i],b[i]) exceeded the threshold. When no threshold is given, the libMesh TOLERANCE is used.

Definition at line 144 of file numeric_vector.C.

References std::abs(), libMesh::comm, libMesh::NumericVector< T >::first_local_index(), libMesh::NumericVector< T >::initialized(), libMesh::initialized(), libMesh::NumericVector< T >::last_local_index(), libMesh::libmesh_assert(), and std::max().

146 {
147  libmesh_assert (this->initialized());
148  libmesh_assert (other_vector.initialized());
149  libmesh_assert_equal_to (this->first_local_index(), other_vector.first_local_index());
150  libmesh_assert_equal_to (this->last_local_index(), other_vector.last_local_index());
151 
152  int first_different_i = std::numeric_limits<int>::max();
154 
155  do
156  {
157  if ( std::abs( (*this)(i) - other_vector(i) ) > threshold *
158  std::max(std::abs((*this)(i)), std::abs(other_vector(i))))
159  first_different_i = i;
160  else
161  i++;
162  }
163  while (first_different_i==std::numeric_limits<int>::max()
164  && i<last_local_index());
165 
166  // Find the correct first differing index in parallel
167  this->comm().min(first_different_i);
168 
169  if (first_different_i == std::numeric_limits<int>::max())
170  return -1;
171 
172  return first_different_i;
173 }
template<typename T >
numeric_index_type libMesh::PetscVector< T >::local_size ( ) const
inlinevirtual
Returns
the local size of the vector (index_stop-index_start)

Implements libMesh::NumericVector< T >.

Definition at line 1107 of file petsc_vector.h.

References libMesh::initialized(), and libMesh::libmesh_assert().

Referenced by libMesh::PetscVector< T >::operator=().

1108 {
1109  libmesh_assert (this->initialized());
1110 
1111  PetscErrorCode ierr=0;
1112  PetscInt petsc_size=0;
1113 
1114  ierr = VecGetLocalSize(_vec, &petsc_size);
1115  LIBMESH_CHKERRABORT(ierr);
1116 
1117  return static_cast<numeric_index_type>(petsc_size);
1118 }
template<typename T >
void libMesh::PetscVector< T >::localize ( std::vector< T > &  v_local) const
virtual

Creates a copy of the global vector in the local vector v_local.

Implements libMesh::NumericVector< T >.

Definition at line 1060 of file petsc_vector.C.

References libMesh::comm.

Referenced by libMesh::PetscVector< T >::localize().

1061 {
1062  this->_restore_array();
1063 
1064  // This function must be run on all processors at once
1065  parallel_object_only();
1066 
1067  PetscErrorCode ierr=0;
1068  const PetscInt n = this->size();
1069  const PetscInt nl = this->local_size();
1070  PetscScalar *values;
1071 
1072  v_local.clear();
1073  v_local.resize(n, 0.);
1074 
1075  ierr = VecGetArray (_vec, &values);
1076  LIBMESH_CHKERRABORT(ierr);
1077 
1079 
1080  for (PetscInt i=0; i<nl; i++)
1081  v_local[i+ioff] = static_cast<T>(values[i]);
1082 
1083  ierr = VecRestoreArray (_vec, &values);
1084  LIBMESH_CHKERRABORT(ierr);
1085 
1086  this->comm().sum(v_local);
1087 }
template<typename T >
void libMesh::PetscVector< T >::localize ( NumericVector< T > &  v_local) const
virtual

Same, but fills a NumericVector<T> instead of a std::vector.

Implements libMesh::NumericVector< T >.

Definition at line 817 of file petsc_vector.C.

References libMesh::PetscVector< T >::_vec, libMesh::PetscVector< T >::close(), libMesh::comm, libMeshEnums::GHOSTED, libMesh::MeshTools::Generation::Private::idx(), libMesh::Utility::iota(), libMesh::libmesh_assert(), PETSC_USE_POINTER, libMesh::PetscVector< T >::size(), and libMesh::NumericVector< T >::type().

818 {
819  this->_restore_array();
820 
821  // Make sure the NumericVector passed in is really a PetscVector
822  PetscVector<T>* v_local = libmesh_cast_ptr<PetscVector<T>*>(&v_local_in);
823 
824  libmesh_assert(v_local);
825  libmesh_assert_equal_to (v_local->size(), this->size());
826 
827  PetscErrorCode ierr = 0;
828  const PetscInt n = this->size();
829 
830  IS is;
831  VecScatter scatter;
832 
833  // Create idx, idx[i] = i;
834  std::vector<PetscInt> idx(n); Utility::iota (idx.begin(), idx.end(), 0);
835 
836  // Create the index set & scatter object
837  ierr = ISCreateLibMesh(this->comm().get(), n, &idx[0], PETSC_USE_POINTER, &is);
838  LIBMESH_CHKERRABORT(ierr);
839 
840  ierr = VecScatterCreate(_vec, is,
841  v_local->_vec, is,
842  &scatter);
843  LIBMESH_CHKERRABORT(ierr);
844 
845  // Perform the scatter
846 #if PETSC_VERSION_LESS_THAN(2,3,3)
847 
848  ierr = VecScatterBegin(_vec, v_local->_vec, INSERT_VALUES,
849  SCATTER_FORWARD, scatter);
850  LIBMESH_CHKERRABORT(ierr);
851 
852  ierr = VecScatterEnd (_vec, v_local->_vec, INSERT_VALUES,
853  SCATTER_FORWARD, scatter);
854  LIBMESH_CHKERRABORT(ierr);
855 #else
856  // API argument order change in PETSc 2.3.3
857  ierr = VecScatterBegin(scatter, _vec, v_local->_vec,
858  INSERT_VALUES, SCATTER_FORWARD);
859  LIBMESH_CHKERRABORT(ierr);
860 
861  ierr = VecScatterEnd (scatter, _vec, v_local->_vec,
862  INSERT_VALUES, SCATTER_FORWARD);
863  LIBMESH_CHKERRABORT(ierr);
864 #endif
865 
866  // Clean up
867  ierr = LibMeshISDestroy (&is);
868  LIBMESH_CHKERRABORT(ierr);
869 
870  ierr = LibMeshVecScatterDestroy(&scatter);
871  LIBMESH_CHKERRABORT(ierr);
872 
873  // Make sure ghost dofs are up to date
874  if (v_local->type() == GHOSTED)
875  v_local->close();
876 }
template<typename T >
void libMesh::PetscVector< T >::localize ( NumericVector< T > &  v_local,
const std::vector< numeric_index_type > &  send_list 
) const
virtual

Creates a local vector v_local containing only information relevant to this processor, as defined by the send_list.

Implements libMesh::NumericVector< T >.

Definition at line 881 of file petsc_vector.C.

References libMesh::PetscVector< T >::_vec, libMesh::PetscVector< T >::close(), libMesh::comm, libMeshEnums::GHOSTED, libMesh::MeshTools::Generation::Private::idx(), libMesh::libmesh_assert(), libMeshEnums::PARALLEL, PETSC_USE_POINTER, libMesh::PetscVector< T >::size(), and libMesh::NumericVector< T >::type().

883 {
884  // FIXME: Workaround for a strange bug at large-scale.
885  // If we have ghosting, PETSc lets us just copy the solution, and
886  // doing so avoids a segfault?
887  if (v_local_in.type() == GHOSTED &&
888  this->type() == PARALLEL)
889  {
890  v_local_in = *this;
891  return;
892  }
893 
894  // Normal code path begins here
895 
896  this->_restore_array();
897 
898  // Make sure the NumericVector passed in is really a PetscVector
899  PetscVector<T>* v_local = libmesh_cast_ptr<PetscVector<T>*>(&v_local_in);
900 
901  libmesh_assert(v_local);
902  libmesh_assert_equal_to (v_local->size(), this->size());
903  libmesh_assert_less_equal (send_list.size(), v_local->size());
904 
906  const numeric_index_type n_sl = send_list.size();
907 
908  IS is;
909  VecScatter scatter;
910 
911  std::vector<PetscInt> idx(n_sl + this->local_size());
912 
913  for (numeric_index_type i=0; i<n_sl; i++)
914  idx[i] = static_cast<PetscInt>(send_list[i]);
915  for (numeric_index_type i = 0; i != this->local_size(); ++i)
916  idx[n_sl+i] = i + this->first_local_index();
917 
918  // Create the index set & scatter object
919  if (idx.empty())
920  ierr = ISCreateLibMesh(this->comm().get(),
921  n_sl+this->local_size(), PETSC_NULL, PETSC_USE_POINTER, &is);
922  else
923  ierr = ISCreateLibMesh(this->comm().get(),
924  n_sl+this->local_size(), &idx[0], PETSC_USE_POINTER, &is);
925  LIBMESH_CHKERRABORT(ierr);
926 
927  ierr = VecScatterCreate(_vec, is,
928  v_local->_vec, is,
929  &scatter);
930  LIBMESH_CHKERRABORT(ierr);
931 
932 
933  // Perform the scatter
934 #if PETSC_VERSION_LESS_THAN(2,3,3)
935 
936  ierr = VecScatterBegin(_vec, v_local->_vec, INSERT_VALUES,
937  SCATTER_FORWARD, scatter);
938  LIBMESH_CHKERRABORT(ierr);
939 
940  ierr = VecScatterEnd (_vec, v_local->_vec, INSERT_VALUES,
941  SCATTER_FORWARD, scatter);
942  LIBMESH_CHKERRABORT(ierr);
943 
944 #else
945 
946  // API argument order change in PETSc 2.3.3
947  ierr = VecScatterBegin(scatter, _vec, v_local->_vec,
948  INSERT_VALUES, SCATTER_FORWARD);
949  LIBMESH_CHKERRABORT(ierr);
950 
951  ierr = VecScatterEnd (scatter, _vec, v_local->_vec,
952  INSERT_VALUES, SCATTER_FORWARD);
953  LIBMESH_CHKERRABORT(ierr);
954 
955 #endif
956 
957 
958  // Clean up
959  ierr = LibMeshISDestroy (&is);
960  LIBMESH_CHKERRABORT(ierr);
961 
962  ierr = LibMeshVecScatterDestroy(&scatter);
963  LIBMESH_CHKERRABORT(ierr);
964 
965  // Make sure ghost dofs are up to date
966  if (v_local->type() == GHOSTED)
967  v_local->close();
968 }
template<typename T >
void libMesh::PetscVector< T >::localize ( const numeric_index_type  first_local_idx,
const numeric_index_type  last_local_idx,
const std::vector< numeric_index_type > &  send_list 
)
virtual

Updates a local vector with selected values from neighboring processors, as defined by send_list.

Implements libMesh::NumericVector< T >.

Definition at line 972 of file petsc_vector.C.

References libMesh::PetscVector< T >::_vec, libMesh::PetscVector< T >::close(), libMesh::comm, libMesh::MeshTools::Generation::Private::idx(), libMesh::PetscVector< T >::init(), libMesh::Utility::iota(), libMesh::PetscVector< T >::localize(), libMesh::n_processors(), libMeshEnums::PARALLEL, and PETSC_USE_POINTER.

975 {
976  this->_restore_array();
977 
978  libmesh_assert_less_equal (send_list.size(), this->size());
979  libmesh_assert_less_equal (last_local_idx+1, this->size());
980 
981  const numeric_index_type my_size = this->size();
982  const numeric_index_type my_local_size = (last_local_idx - first_local_idx + 1);
984 
985  // Don't bother for serial cases
986 // if ((first_local_idx == 0) &&
987 // (my_local_size == my_size))
988  // But we do need to stay in sync for degenerate cases
989  if (this->n_processors() == 1)
990  return;
991 
992 
993  // Build a parallel vector, initialize it with the local
994  // parts of (*this)
995  PetscVector<T> parallel_vec(this->comm(), PARALLEL);
996 
997  parallel_vec.init (my_size, my_local_size, true, PARALLEL);
998 
999 
1000  // Copy part of *this into the parallel_vec
1001  {
1002  IS is;
1003  VecScatter scatter;
1004 
1005  // Create idx, idx[i] = i+first_local_idx;
1006  std::vector<PetscInt> idx(my_local_size);
1007  Utility::iota (idx.begin(), idx.end(), first_local_idx);
1008 
1009  // Create the index set & scatter object
1010  ierr = ISCreateLibMesh(this->comm().get(), my_local_size,
1011  my_local_size ? &idx[0] : NULL, PETSC_USE_POINTER, &is);
1012  LIBMESH_CHKERRABORT(ierr);
1013 
1014  ierr = VecScatterCreate(_vec, is,
1015  parallel_vec._vec, is,
1016  &scatter);
1017  LIBMESH_CHKERRABORT(ierr);
1018 
1019  // Perform the scatter
1020 #if PETSC_VERSION_LESS_THAN(2,3,3)
1021 
1022  ierr = VecScatterBegin(_vec, parallel_vec._vec, INSERT_VALUES,
1023  SCATTER_FORWARD, scatter);
1024  LIBMESH_CHKERRABORT(ierr);
1025 
1026  ierr = VecScatterEnd (_vec, parallel_vec._vec, INSERT_VALUES,
1027  SCATTER_FORWARD, scatter);
1028  LIBMESH_CHKERRABORT(ierr);
1029 
1030 #else
1031 
1032  // API argument order change in PETSc 2.3.3
1033  ierr = VecScatterBegin(scatter, _vec, parallel_vec._vec,
1034  INSERT_VALUES, SCATTER_FORWARD);
1035  LIBMESH_CHKERRABORT(ierr);
1036 
1037  ierr = VecScatterEnd (scatter, _vec, parallel_vec._vec,
1038  INSERT_VALUES, SCATTER_FORWARD);
1039  LIBMESH_CHKERRABORT(ierr);
1040 
1041 #endif
1042 
1043  // Clean up
1044  ierr = LibMeshISDestroy (&is);
1045  LIBMESH_CHKERRABORT(ierr);
1046 
1047  ierr = LibMeshVecScatterDestroy(&scatter);
1048  LIBMESH_CHKERRABORT(ierr);
1049  }
1050 
1051  // localize like normal
1052  parallel_vec.close();
1053  parallel_vec.localize (*this, send_list);
1054  this->close();
1055 }
template<typename T>
void libMesh::PetscVector< T >::localize_to_one ( std::vector< T > &  v_local,
const processor_id_type  proc_id = 0 
) const
virtual

Creates a local copy of the global vector in v_local only on processor proc_id. By default the data is sent to processor 0. This method is useful for outputting data from one processor.

Implements libMesh::NumericVector< T >.

template<>
void libMesh::PetscVector< Real >::localize_to_one ( std::vector< Real > &  v_local,
const processor_id_type  pid 
) const

Definition at line 1095 of file petsc_vector.C.

References libMesh::comm, libMesh::n_processors(), and libMesh::processor_id().

1097 {
1098  this->_restore_array();
1099 
1100  PetscErrorCode ierr=0;
1101  const PetscInt n = size();
1102  const PetscInt nl = local_size();
1103  PetscScalar *values;
1104 
1105 
1106  // only one processor
1107  if (n_processors() == 1)
1108  {
1109  v_local.resize(n);
1110 
1111  ierr = VecGetArray (_vec, &values);
1112  LIBMESH_CHKERRABORT(ierr);
1113 
1114  for (PetscInt i=0; i<n; i++)
1115  v_local[i] = static_cast<Real>(values[i]);
1116 
1117  ierr = VecRestoreArray (_vec, &values);
1118  LIBMESH_CHKERRABORT(ierr);
1119  }
1120 
1121  // otherwise multiple processors
1122  else
1123  {
1124  if(pid == 0) // optimized version for localizing to 0
1125  {
1126  Vec vout;
1127  VecScatter ctx;
1128 
1129  ierr = VecScatterCreateToZero(_vec, &ctx, &vout);
1130  LIBMESH_CHKERRABORT(ierr);
1131 
1132  ierr = VecScatterBegin(ctx, _vec, vout, INSERT_VALUES, SCATTER_FORWARD);
1133  LIBMESH_CHKERRABORT(ierr);
1134  ierr = VecScatterEnd(ctx, _vec, vout, INSERT_VALUES, SCATTER_FORWARD);
1135  LIBMESH_CHKERRABORT(ierr);
1136 
1137  if(processor_id() == 0)
1138  {
1139  v_local.resize(n);
1140 
1141  ierr = VecGetArray (vout, &values);
1142  LIBMESH_CHKERRABORT(ierr);
1143 
1144  for (PetscInt i=0; i<n; i++)
1145  v_local[i] = static_cast<Real>(values[i]);
1146 
1147  ierr = VecRestoreArray (vout, &values);
1148  LIBMESH_CHKERRABORT(ierr);
1149  }
1150 
1151  ierr = LibMeshVecScatterDestroy(&ctx);
1152  LIBMESH_CHKERRABORT(ierr);
1153  ierr = LibMeshVecDestroy(&vout);
1154  LIBMESH_CHKERRABORT(ierr);
1155 
1156  }
1157  else
1158  {
1159  v_local.resize(n);
1160 
1161  numeric_index_type ioff = this->first_local_index();
1162  std::vector<Real> local_values (n, 0.);
1163 
1164  {
1165  ierr = VecGetArray (_vec, &values);
1166  LIBMESH_CHKERRABORT(ierr);
1167 
1168  for (PetscInt i=0; i<nl; i++)
1169  local_values[i+ioff] = static_cast<Real>(values[i]);
1170 
1171  ierr = VecRestoreArray (_vec, &values);
1172  LIBMESH_CHKERRABORT(ierr);
1173  }
1174 
1175 
1176  MPI_Reduce (&local_values[0], &v_local[0], n, MPI_REAL, MPI_SUM,
1177  pid, this->comm().get());
1178  }
1179  }
1180 }
template<>
void libMesh::PetscVector< Complex >::localize_to_one ( std::vector< Complex > &  v_local,
const processor_id_type  pid 
) const

Definition at line 1189 of file petsc_vector.C.

References libMesh::comm.

1191 {
1192  this->_restore_array();
1193 
1194  PetscErrorCode ierr=0;
1195  const PetscInt n = size();
1196  const PetscInt nl = local_size();
1197  PetscScalar *values;
1198 
1199 
1200  v_local.resize(n);
1201 
1202 
1203  for (PetscInt i=0; i<n; i++)
1204  v_local[i] = 0.;
1205 
1206  // only one processor
1207  if (n == nl)
1208  {
1209  ierr = VecGetArray (_vec, &values);
1210  LIBMESH_CHKERRABORT(ierr);
1211 
1212  for (PetscInt i=0; i<n; i++)
1213  v_local[i] = static_cast<Complex>(values[i]);
1214 
1215  ierr = VecRestoreArray (_vec, &values);
1216  LIBMESH_CHKERRABORT(ierr);
1217  }
1218 
1219  // otherwise multiple processors
1220  else
1221  {
1222  numeric_index_type ioff = this->first_local_index();
1223 
1224  /* in here the local values are stored, acting as send buffer for MPI
1225  * initialize to zero, since we collect using MPI_SUM
1226  */
1227  std::vector<Real> real_local_values(n, 0.);
1228  std::vector<Real> imag_local_values(n, 0.);
1229 
1230  {
1231  ierr = VecGetArray (_vec, &values);
1232  LIBMESH_CHKERRABORT(ierr);
1233 
1234  // provide my local share to the real and imag buffers
1235  for (PetscInt i=0; i<nl; i++)
1236  {
1237  real_local_values[i+ioff] = static_cast<Complex>(values[i]).real();
1238  imag_local_values[i+ioff] = static_cast<Complex>(values[i]).imag();
1239  }
1240 
1241  ierr = VecRestoreArray (_vec, &values);
1242  LIBMESH_CHKERRABORT(ierr);
1243  }
1244 
1245  /* have buffers of the real and imaginary part of v_local.
1246  * Once MPI_Reduce() collected all the real and imaginary
1247  * parts in these std::vector<Real>, the values can be
1248  * copied to v_local
1249  */
1250  std::vector<Real> real_v_local(n);
1251  std::vector<Real> imag_v_local(n);
1252 
1253  // collect entries from other proc's in real_v_local, imag_v_local
1254  MPI_Reduce (&real_local_values[0], &real_v_local[0], n,
1255  MPI_REAL, MPI_SUM,
1256  pid, this->comm().get());
1257 
1258  MPI_Reduce (&imag_local_values[0], &imag_v_local[0], n,
1259  MPI_REAL, MPI_SUM,
1260  pid, this->comm().get());
1261 
1262  // copy real_v_local and imag_v_local to v_local
1263  for (PetscInt i=0; i<n; i++)
1264  v_local[i] = Complex(real_v_local[i], imag_v_local[i]);
1265  }
1266 }
template<typename T >
numeric_index_type libMesh::PetscVector< T >::map_global_to_local_index ( const numeric_index_type  i) const
inline

Maps the global index i to the corresponding global index. If the index is not a ghost cell, this is done by subtraction the number of the first local index. If it is a ghost cell, it has to be looked up in the map.

Definition at line 1158 of file petsc_vector.h.

References end, libMesh::err, libMesh::initialized(), and libMesh::libmesh_assert().

1159 {
1160  libmesh_assert (this->initialized());
1161 
1162  PetscErrorCode ierr=0;
1163  PetscInt petsc_first=0, petsc_last=0;
1164  ierr = VecGetOwnershipRange (_vec, &petsc_first, &petsc_last);
1165  LIBMESH_CHKERRABORT(ierr);
1166  const numeric_index_type first = static_cast<numeric_index_type>(petsc_first);
1167  const numeric_index_type last = static_cast<numeric_index_type>(petsc_last);
1168 
1169  if((i>=first) && (i<last))
1170  {
1171  return i-first;
1172  }
1173 
1174  GlobalToLocalMap::const_iterator it = _global_to_local_map.find(i);
1175 #ifndef NDEBUG
1176  const GlobalToLocalMap::const_iterator end = _global_to_local_map.end();
1177  if (it == end)
1178  {
1179  std::ostringstream error_message;
1180  error_message << "No index " << i << " in ghosted vector.\n"
1181  << "Vector contains [" << first << ',' << last << ")\n";
1182  GlobalToLocalMap::const_iterator b = _global_to_local_map.begin();
1183  if (b == end)
1184  {
1185  error_message << "And empty ghost array.\n";
1186  }
1187  else
1188  {
1189  error_message << "And ghost array {" << b->first;
1190  for (b++; b != end; b++)
1191  error_message << ',' << b->first;
1192  error_message << "}\n";
1193  }
1194 
1195  libMesh::err << error_message.str();
1196 
1197  libmesh_error();
1198  }
1199  libmesh_assert (it != _global_to_local_map.end());
1200 #endif
1201  return it->second+last-first;
1202 }
template<typename T >
Real libMesh::PetscVector< T >::max ( ) const
inlinevirtual
Returns
the maximum element in the vector. In case of complex numbers, this returns the maximum Real part.

Implements libMesh::NumericVector< T >.

Definition at line 1271 of file petsc_vector.h.

References libMesh::Real.

1272 {
1273  this->_restore_array();
1274 
1275  PetscErrorCode ierr=0;
1276  PetscInt index=0;
1277  PetscReal returnval=0.;
1278 
1279  ierr = VecMax (_vec, &index, &returnval);
1280  LIBMESH_CHKERRABORT(ierr);
1281 
1282  // this return value is correct: VecMax returns a PetscReal
1283  return static_cast<Real>(returnval);
1284 }
template<typename T >
Real libMesh::PetscVector< T >::min ( ) const
inlinevirtual
Returns
the minimum element in the vector. In case of complex numbers, this returns the minimum Real part.

Implements libMesh::NumericVector< T >.

Definition at line 1252 of file petsc_vector.h.

References libMesh::Real.

1253 {
1254  this->_restore_array();
1255 
1256  PetscErrorCode ierr=0;
1257  PetscInt index=0;
1258  PetscReal returnval=0.;
1259 
1260  ierr = VecMin (_vec, &index, &returnval);
1261  LIBMESH_CHKERRABORT(ierr);
1262 
1263  // this return value is correct: VecMin returns a PetscReal
1264  return static_cast<Real>(returnval);
1265 }
static unsigned int libMesh::ReferenceCounter::n_objects ( )
inlinestaticinherited

Prints the number of outstanding (created, but not yet destroyed) objects.

Definition at line 79 of file reference_counter.h.

References libMesh::ReferenceCounter::_n_objects.

80  { return _n_objects; }
processor_id_type libMesh::ParallelObject::n_processors ( ) const
inlineinherited
Returns
the number of processors in the group.

Definition at line 92 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and libMesh::Parallel::Communicator::size().

Referenced by libMesh::ParmetisPartitioner::_do_repartition(), libMesh::ParallelMesh::add_elem(), libMesh::ParallelMesh::add_node(), libMesh::LaplaceMeshSmoother::allgather_graph(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::ParallelMesh::assign_unique_ids(), libMesh::AztecLinearSolver< T >::AztecLinearSolver(), libMesh::ParallelMesh::clear(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::UnstructuredMesh::create_pid_mesh(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::EnsightIO::EnsightIO(), libMesh::MeshBase::get_info(), libMesh::EquationSystems::init(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), libMesh::Nemesis_IO_Helper::initialize(), libMesh::MeshTools::libmesh_assert_valid_dof_ids(), libMesh::MeshTools::libmesh_assert_valid_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_procids< Node >(), libMesh::MeshTools::libmesh_assert_valid_refinement_flags(), libMesh::MeshBase::n_active_elem_on_proc(), libMesh::MeshBase::n_elem_on_proc(), libMesh::MeshBase::n_nodes_on_proc(), libMesh::Partitioner::partition(), libMesh::MeshBase::partition(), libMesh::Partitioner::partition_unpartitioned_elements(), libMesh::PetscLinearSolver< T >::PetscLinearSolver(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::MeshTools::processor_bounding_box(), libMesh::System::project_vector(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::UnstructuredMesh::read(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::System::read_serialized_vector(), libMesh::Partitioner::repartition(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::BoundaryInfo::sync(), libMesh::ParallelMesh::update_parallel_id_counts(), libMesh::CheckpointIO::write(), libMesh::GMVIO::write_binary(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::System::write_parallel_data(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bcs(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::XdrIO::write_serialized_nodes(), and libMesh::XdrIO::write_serialized_nodesets().

93  { return libmesh_cast_int<processor_id_type>(_communicator.size()); }
template<typename T >
T libMesh::PetscVector< T >::operator() ( const numeric_index_type  i) const
inlinevirtual

Access components, returns U(i).

Implements libMesh::NumericVector< T >.

Definition at line 1208 of file petsc_vector.h.

References libMeshEnums::GHOSTED.

1209 {
1210  this->_get_array();
1211 
1212  const numeric_index_type local_index = this->map_global_to_local_index(i);
1213 
1214 #ifndef NDEBUG
1215  if(this->type() == GHOSTED)
1216  {
1217  libmesh_assert_less (local_index, _local_size);
1218  }
1219 #endif
1220 
1221  return static_cast<T>(_values[local_index]);
1222 }
template<typename T>
NumericVector<T>& libMesh::NumericVector< T >::operator*= ( const T  a)
inlineinherited

Multiplication operator. Equivalent to U.scale(a)

Definition at line 368 of file numeric_vector.h.

368 { this->scale(a); return *this; }
template<typename T >
NumericVector< T > & libMesh::PetscVector< T >::operator+= ( const NumericVector< T > &  V)
virtual

Addition operator. Fast equivalent to U.add(1, V).

Implements libMesh::NumericVector< T >.

Definition at line 122 of file petsc_vector.C.

References libMesh::closed(), and libMesh::libmesh_assert().

123 {
124  this->_restore_array();
125  libmesh_assert(this->closed());
126 
127  this->add(1., v);
128 
129  return *this;
130 }
template<typename T >
NumericVector< T > & libMesh::PetscVector< T >::operator-= ( const NumericVector< T > &  V)
virtual

Subtraction operator. Fast equivalent to U.add(-1, V).

Implements libMesh::NumericVector< T >.

Definition at line 136 of file petsc_vector.C.

References libMesh::closed(), and libMesh::libmesh_assert().

137 {
138  this->_restore_array();
139  libmesh_assert(this->closed());
140 
141  this->add(-1., v);
142 
143  return *this;
144 }
template<typename T>
NumericVector<T>& libMesh::NumericVector< T >::operator/= ( const T  a)
inlineinherited

Division operator. Equivalent to U.scale(1./a)

Definition at line 374 of file numeric_vector.h.

374 { this->scale(1./a); return *this; }
template<typename T >
NumericVector< T > & libMesh::PetscVector< T >::operator/= ( NumericVector< T > &  v)
virtual

Pointwise Division operator. ie divide every entry in this vector by the entry in v

Implements libMesh::NumericVector< T >.

Definition at line 552 of file petsc_vector.C.

References libMesh::PetscVector< T >::_vec.

553 {
554  PetscErrorCode ierr = 0;
555 
556  const PetscVector<T>* v_vec = libmesh_cast_ptr<const PetscVector<T>*>(&v);
557 
558  ierr = VecPointwiseDivide(_vec, _vec, v_vec->_vec);
559  LIBMESH_CHKERRABORT(ierr);
560 
561  return *this;
562 }
template<typename T >
NumericVector< T > & libMesh::PetscVector< T >::operator= ( const T  s)
virtual

Change the dimension to that of the vector V. The same applies as for the other init function.

The elements of V are not copied, i.e. this function is the same as calling init(V.size(),fast). $U(0-N) = s$: fill all components.

Implements libMesh::NumericVector< T >.

Definition at line 635 of file petsc_vector.C.

References libMesh::closed(), libMeshEnums::GHOSTED, and libMesh::libmesh_assert().

636 {
637  this->_restore_array();
638  libmesh_assert(this->closed());
639 
640  PetscErrorCode ierr = 0;
641  PetscScalar s = static_cast<PetscScalar>(s_in);
642 
643  if (this->size() != 0)
644  {
645  if(this->type() != GHOSTED)
646  {
647 #if PETSC_VERSION_LESS_THAN(2,3,0)
648  // 2.2.x & earlier style
649  ierr = VecSet(&s, _vec);
650  LIBMESH_CHKERRABORT(ierr);
651 #else
652  // 2.3.x & later style
653  ierr = VecSet(_vec, s);
654  LIBMESH_CHKERRABORT(ierr);
655 #endif
656  }
657  else
658  {
659  Vec loc_vec;
660  ierr = VecGhostGetLocalForm (_vec,&loc_vec);
661  LIBMESH_CHKERRABORT(ierr);
662 
663 #if PETSC_VERSION_LESS_THAN(2,3,0)
664  // 2.2.x & earlier style
665  ierr = VecSet(&s, loc_vec);
666  LIBMESH_CHKERRABORT(ierr);
667 #else
668  // 2.3.x & later style
669  ierr = VecSet(loc_vec, s);
670  LIBMESH_CHKERRABORT(ierr);
671 #endif
672 
673  ierr = VecGhostRestoreLocalForm (_vec,&loc_vec);
674  LIBMESH_CHKERRABORT(ierr);
675  }
676  }
677 
678  return *this;
679 }
template<typename T >
NumericVector< T > & libMesh::PetscVector< T >::operator= ( const NumericVector< T > &  V)
virtual

$U = V$: copy all components.

Implements libMesh::NumericVector< T >.

Definition at line 685 of file petsc_vector.C.

686 {
687  // Make sure the NumericVector passed in is really a PetscVector
688  const PetscVector<T>* v = libmesh_cast_ptr<const PetscVector<T>*>(&v_in);
689 
690  *this = *v;
691 
692  return *this;
693 }
template<typename T >
PetscVector< T > & libMesh::PetscVector< T >::operator= ( const PetscVector< T > &  V)

$U = V$: copy all components.

Definition at line 699 of file petsc_vector.C.

References libMesh::PetscVector< T >::_global_to_local_map, libMesh::PetscVector< T >::_restore_array(), libMesh::NumericVector< T >::_type, libMesh::PetscVector< T >::_vec, libMesh::NumericVector< T >::closed(), libMeshEnums::GHOSTED, libMesh::libmesh_assert(), libMesh::PetscVector< T >::local_size(), libMeshEnums::PARALLEL, libMeshEnums::SERIAL, libMesh::PetscVector< T >::size(), and libMesh::NumericVector< T >::type().

700 {
701  this->_restore_array();
702  v._restore_array();
703 
704  libmesh_assert_equal_to (this->size(), v.size());
705  libmesh_assert_equal_to (this->local_size(), v.local_size());
706  libmesh_assert (v.closed());
707 
708  PetscErrorCode ierr = 0;
709 
710  if (((this->type()==PARALLEL) && (v.type()==GHOSTED)) ||
711  ((this->type()==GHOSTED) && (v.type()==PARALLEL)) ||
712  ((this->type()==GHOSTED) && (v.type()==SERIAL)) ||
713  ((this->type()==SERIAL) && (v.type()==GHOSTED)))
714  {
715  /* Allow assignment of a ghosted to a parallel vector since this
716  causes no difficulty. See discussion in libmesh-devel of
717  June 24, 2010. */
718  ierr = VecCopy (v._vec, this->_vec);
719  LIBMESH_CHKERRABORT(ierr);
720  }
721  else
722  {
723  /* In all other cases, we assert that both vectors are of equal
724  type. */
725  libmesh_assert_equal_to (this->_type, v._type);
726  libmesh_assert (this->_global_to_local_map == v._global_to_local_map);
727 
728  if (v.size() != 0)
729  {
730  if(this->type() != GHOSTED)
731  {
732  ierr = VecCopy (v._vec, this->_vec);
733  LIBMESH_CHKERRABORT(ierr);
734  }
735  else
736  {
737  Vec loc_vec;
738  Vec v_loc_vec;
739  ierr = VecGhostGetLocalForm (_vec,&loc_vec);
740  LIBMESH_CHKERRABORT(ierr);
741  ierr = VecGhostGetLocalForm (v._vec,&v_loc_vec);
742  LIBMESH_CHKERRABORT(ierr);
743 
744  ierr = VecCopy (v_loc_vec, loc_vec);
745  LIBMESH_CHKERRABORT(ierr);
746 
747  ierr = VecGhostRestoreLocalForm (v._vec,&v_loc_vec);
748  LIBMESH_CHKERRABORT(ierr);
749  ierr = VecGhostRestoreLocalForm (_vec,&loc_vec);
750  LIBMESH_CHKERRABORT(ierr);
751  }
752  }
753  }
754 
755  close();
756 
757  return *this;
758 }
template<typename T >
NumericVector< T > & libMesh::PetscVector< T >::operator= ( const std::vector< T > &  v)
virtual

$U = V$: copy all components.

Case 1: The vector is the same size of The global vector. Only add the local components.

Case 2: The vector is the same size as our local piece. Insert directly to the local piece.

Implements libMesh::NumericVector< T >.

Definition at line 764 of file petsc_vector.C.

References libMeshEnums::GHOSTED.

765 {
766  this->_restore_array();
767 
768  const numeric_index_type nl = this->local_size();
769  const numeric_index_type ioff = this->first_local_index();
771  PetscScalar* values;
772 
777  if (this->size() == v.size())
778  {
779  ierr = VecGetArray (_vec, &values);
780  LIBMESH_CHKERRABORT(ierr);
781 
782  for (numeric_index_type i=0; i<nl; i++)
783  values[i] = static_cast<PetscScalar>(v[i+ioff]);
784 
785  ierr = VecRestoreArray (_vec, &values);
786  LIBMESH_CHKERRABORT(ierr);
787  }
788 
793  else
794  {
795  libmesh_assert_equal_to (this->local_size(), v.size());
796 
797  ierr = VecGetArray (_vec, &values);
798  LIBMESH_CHKERRABORT(ierr);
799 
800  for (numeric_index_type i=0; i<nl; i++)
801  values[i] = static_cast<PetscScalar>(v[i]);
802 
803  ierr = VecRestoreArray (_vec, &values);
804  LIBMESH_CHKERRABORT(ierr);
805  }
806 
807  // Make sure ghost dofs are up to date
808  if (this->type() == GHOSTED)
809  this->close();
810 
811  return *this;
812 }
template<typename T >
void libMesh::PetscVector< T >::pointwise_mult ( const NumericVector< T > &  vec1,
const NumericVector< T > &  vec2 
)
virtual

Computes the pointwise (i.e. component-wise) product of vec1 and vec2 and stores the result in *this.

Implements libMesh::NumericVector< T >.

Definition at line 1273 of file petsc_vector.C.

References libMeshEnums::GHOSTED, and libMesh::out.

1275 {
1276  this->_restore_array();
1277 
1278  PetscErrorCode ierr = 0;
1279 
1280  // Convert arguments to PetscVector*.
1281  const PetscVector<T>* vec1_petsc = libmesh_cast_ptr<const PetscVector<T>*>(&vec1);
1282  const PetscVector<T>* vec2_petsc = libmesh_cast_ptr<const PetscVector<T>*>(&vec2);
1283 
1284  // Call PETSc function.
1285 
1286 #if PETSC_VERSION_LESS_THAN(2,3,1)
1287 
1288  libMesh::out << "This method has been developed with PETSc 2.3.1. "
1289  << "No one has made it backwards compatible with older "
1290  << "versions of PETSc so far; however, it might work "
1291  << "without any change with some older version." << std::endl;
1292  libmesh_error();
1293 
1294 #else
1295 
1296  if(this->type() != GHOSTED)
1297  {
1298  ierr = VecPointwiseMult(this->vec(),
1299  const_cast<PetscVector<T>*>(vec1_petsc)->vec(),
1300  const_cast<PetscVector<T>*>(vec2_petsc)->vec());
1301  LIBMESH_CHKERRABORT(ierr);
1302  }
1303  else
1304  {
1305  Vec loc_vec;
1306  Vec v1_loc_vec;
1307  Vec v2_loc_vec;
1308  ierr = VecGhostGetLocalForm (_vec,&loc_vec);
1309  LIBMESH_CHKERRABORT(ierr);
1310  ierr = VecGhostGetLocalForm (const_cast<PetscVector<T>*>(vec1_petsc)->vec(),&v1_loc_vec);
1311  LIBMESH_CHKERRABORT(ierr);
1312  ierr = VecGhostGetLocalForm (const_cast<PetscVector<T>*>(vec2_petsc)->vec(),&v2_loc_vec);
1313  LIBMESH_CHKERRABORT(ierr);
1314 
1315  ierr = VecPointwiseMult(loc_vec,v1_loc_vec,v2_loc_vec);
1316  LIBMESH_CHKERRABORT(ierr);
1317 
1318  ierr = VecGhostRestoreLocalForm (const_cast<PetscVector<T>*>(vec1_petsc)->vec(),&v1_loc_vec);
1319  LIBMESH_CHKERRABORT(ierr);
1320  ierr = VecGhostRestoreLocalForm (const_cast<PetscVector<T>*>(vec2_petsc)->vec(),&v2_loc_vec);
1321  LIBMESH_CHKERRABORT(ierr);
1322  ierr = VecGhostRestoreLocalForm (_vec,&loc_vec);
1323  LIBMESH_CHKERRABORT(ierr);
1324  }
1325 
1326 #endif
1327 
1328 }
template<typename T >
void libMesh::NumericVector< T >::print ( std::ostream &  os = libMesh::out) const
inlinevirtualinherited

Prints the local contents of the vector, by default to libMesh::out

Definition at line 832 of file numeric_vector.h.

References libMesh::initialized(), and libMesh::libmesh_assert().

833 {
834  libmesh_assert (this->initialized());
835  os << "Size\tglobal = " << this->size()
836  << "\t\tlocal = " << this->local_size() << std::endl;
837 
838  os << "#\tValue" << std::endl;
839  for (numeric_index_type i=this->first_local_index(); i<this->last_local_index(); i++)
840  os << i << "\t" << (*this)(i) << std::endl;
841 }
template<>
void libMesh::NumericVector< Complex >::print ( std::ostream &  os) const
inlineinherited

Definition at line 814 of file numeric_vector.h.

References libMesh::initialized(), and libMesh::libmesh_assert().

815 {
816  libmesh_assert (this->initialized());
817  os << "Size\tglobal = " << this->size()
818  << "\t\tlocal = " << this->local_size() << std::endl;
819 
820  // std::complex<>::operator<<() is defined, but use this form
821  os << "#\tReal part\t\tImaginary part" << std::endl;
822  for (numeric_index_type i=this->first_local_index(); i<this->last_local_index(); i++)
823  os << i << "\t"
824  << (*this)(i).real() << "\t\t"
825  << (*this)(i).imag() << std::endl;
826 }
template<typename T >
void libMesh::NumericVector< T >::print_global ( std::ostream &  os = libMesh::out) const
inlinevirtualinherited

Prints the global contents of the vector, by default to libMesh::out

Definition at line 869 of file numeric_vector.h.

References libMesh::initialized(), libMesh::libmesh_assert(), and libMesh::processor_id().

870 {
871  libmesh_assert (this->initialized());
872 
873  std::vector<T> v(this->size());
874  this->localize(v);
875 
876  // Right now we only want one copy of the output
877  if (this->processor_id())
878  return;
879 
880  os << "Size\tglobal = " << this->size() << std::endl;
881  os << "#\tValue" << std::endl;
882  for (numeric_index_type i=0; i!=v.size(); i++)
883  os << i << "\t" << v[i] << std::endl;
884 }
template<>
void libMesh::NumericVector< Complex >::print_global ( std::ostream &  os) const
inlineinherited

Definition at line 847 of file numeric_vector.h.

References libMesh::initialized(), libMesh::libmesh_assert(), and libMesh::processor_id().

848 {
849  libmesh_assert (this->initialized());
850 
851  std::vector<Complex> v(this->size());
852  this->localize(v);
853 
854  // Right now we only want one copy of the output
855  if (this->processor_id())
856  return;
857 
858  os << "Size\tglobal = " << this->size() << std::endl;
859  os << "#\tReal part\t\tImaginary part" << std::endl;
860  for (numeric_index_type i=0; i!=v.size(); i++)
861  os << i << "\t"
862  << v[i].real() << "\t\t"
863  << v[i].imag() << std::endl;
864 }
void libMesh::ReferenceCounter::print_info ( std::ostream &  out = libMesh::out)
staticinherited

Prints the reference information, by default to libMesh::out.

Definition at line 88 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter, and libMesh::ReferenceCounter::get_info().

89 {
91 }
template<typename T >
void libMesh::PetscVector< T >::print_matlab ( const std::string  name = "NULL") const
virtual

Print the contents of the vector in Matlab format. Optionally prints the matrix to the file named name. If name is not specified it is dumped to the screen.

Create an ASCII file containing the matrix if a filename was provided.

Otherwise the matrix will be dumped to the screen.

Destroy the viewer.

Reimplemented from libMesh::NumericVector< T >.

Definition at line 1333 of file petsc_vector.C.

References libMesh::closed(), libMesh::comm, and libMesh::libmesh_assert().

1334 {
1335  this->_restore_array();
1336  libmesh_assert (this->closed());
1337 
1338  PetscErrorCode ierr=0;
1339  PetscViewer petsc_viewer;
1340 
1341 
1342  ierr = PetscViewerCreate (this->comm().get(),
1343  &petsc_viewer);
1344  LIBMESH_CHKERRABORT(ierr);
1345 
1350  if (name != "NULL")
1351  {
1352  ierr = PetscViewerASCIIOpen( this->comm().get(),
1353  name.c_str(),
1354  &petsc_viewer);
1355  LIBMESH_CHKERRABORT(ierr);
1356 
1357  ierr = PetscViewerSetFormat (petsc_viewer,
1358  PETSC_VIEWER_ASCII_MATLAB);
1359  LIBMESH_CHKERRABORT(ierr);
1360 
1361  ierr = VecView (_vec, petsc_viewer);
1362  LIBMESH_CHKERRABORT(ierr);
1363  }
1364 
1368  else
1369  {
1370  ierr = PetscViewerSetFormat (PETSC_VIEWER_STDOUT_WORLD,
1371  PETSC_VIEWER_ASCII_MATLAB);
1372  LIBMESH_CHKERRABORT(ierr);
1373 
1374  ierr = VecView (_vec, PETSC_VIEWER_STDOUT_WORLD);
1375  LIBMESH_CHKERRABORT(ierr);
1376  }
1377 
1378 
1382  ierr = LibMeshPetscViewerDestroy (&petsc_viewer);
1383  LIBMESH_CHKERRABORT(ierr);
1384 }
processor_id_type libMesh::ParallelObject::processor_id ( ) const
inlineinherited
Returns
the rank of this processor in the group.

Definition at line 98 of file parallel_object.h.

References libMesh::ParallelObject::_communicator, and libMesh::Parallel::Communicator::rank().

Referenced by libMesh::MetisPartitioner::_do_partition(), libMesh::EquationSystems::_read_impl(), libMesh::SerialMesh::active_local_elements_begin(), libMesh::ParallelMesh::active_local_elements_begin(), libMesh::SerialMesh::active_local_elements_end(), libMesh::ParallelMesh::active_local_elements_end(), libMesh::SerialMesh::active_local_subdomain_elements_begin(), libMesh::ParallelMesh::active_local_subdomain_elements_begin(), libMesh::SerialMesh::active_local_subdomain_elements_end(), libMesh::ParallelMesh::active_local_subdomain_elements_end(), libMesh::SerialMesh::active_not_local_elements_begin(), libMesh::ParallelMesh::active_not_local_elements_begin(), libMesh::SerialMesh::active_not_local_elements_end(), libMesh::ParallelMesh::active_not_local_elements_end(), libMesh::ParallelMesh::add_elem(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::ParallelMesh::add_node(), libMesh::UnstructuredMesh::all_second_order(), libMesh::ParmetisPartitioner::assign_partitioning(), libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::Nemesis_IO_Helper::build_element_and_node_maps(), libMesh::ParmetisPartitioner::build_graph(), libMesh::InfElemBuilder::build_inf_elem(), libMesh::DofMap::build_sparsity(), libMesh::ParallelMesh::clear(), libMesh::ExodusII_IO_Helper::close(), libMesh::Nemesis_IO_Helper::compute_border_node_ids(), libMesh::Nemesis_IO_Helper::compute_communication_map_parameters(), libMesh::Nemesis_IO_Helper::compute_internal_and_border_elems_and_internal_nodes(), libMesh::Nemesis_IO_Helper::compute_node_communication_maps(), libMesh::Nemesis_IO_Helper::compute_num_global_elem_blocks(), libMesh::Nemesis_IO_Helper::compute_num_global_nodesets(), libMesh::Nemesis_IO_Helper::compute_num_global_sidesets(), libMesh::Nemesis_IO_Helper::construct_nemesis_filename(), libMesh::ExodusII_IO_Helper::create(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DofMap::end_dof(), libMesh::DofMap::end_old_dof(), libMesh::EnsightIO::EnsightIO(), libMesh::UnstructuredMesh::find_neighbors(), libMesh::DofMap::first_dof(), libMesh::DofMap::first_old_dof(), libMesh::Nemesis_IO_Helper::get_cmap_params(), libMesh::Nemesis_IO_Helper::get_eb_info_global(), libMesh::Nemesis_IO_Helper::get_elem_cmap(), libMesh::Nemesis_IO_Helper::get_elem_map(), libMesh::MeshBase::get_info(), libMesh::Nemesis_IO_Helper::get_init_global(), libMesh::Nemesis_IO_Helper::get_init_info(), libMesh::Nemesis_IO_Helper::get_loadbal_param(), libMesh::Nemesis_IO_Helper::get_node_cmap(), libMesh::Nemesis_IO_Helper::get_node_map(), libMesh::Nemesis_IO_Helper::get_ns_param_global(), libMesh::Nemesis_IO_Helper::get_ss_param_global(), libMesh::MeshFunction::gradient(), libMesh::MeshFunction::hessian(), libMesh::SystemSubsetBySubdomain::init(), libMesh::ParmetisPartitioner::initialize(), libMesh::ExodusII_IO_Helper::initialize(), libMesh::ExodusII_IO_Helper::initialize_discontinuous(), libMesh::ExodusII_IO_Helper::initialize_element_variables(), libMesh::ExodusII_IO_Helper::initialize_global_variables(), libMesh::ExodusII_IO_Helper::initialize_nodal_variables(), libMesh::SparsityPattern::Build::join(), libMesh::DofMap::last_dof(), libMesh::MeshTools::libmesh_assert_valid_procids< Elem >(), libMesh::MeshTools::libmesh_assert_valid_procids< Node >(), libMesh::SerialMesh::local_elements_begin(), libMesh::ParallelMesh::local_elements_begin(), libMesh::SerialMesh::local_elements_end(), libMesh::ParallelMesh::local_elements_end(), libMesh::SerialMesh::local_level_elements_begin(), libMesh::ParallelMesh::local_level_elements_begin(), libMesh::SerialMesh::local_level_elements_end(), libMesh::ParallelMesh::local_level_elements_end(), libMesh::SerialMesh::local_nodes_begin(), libMesh::ParallelMesh::local_nodes_begin(), libMesh::SerialMesh::local_nodes_end(), libMesh::ParallelMesh::local_nodes_end(), libMesh::SerialMesh::local_not_level_elements_begin(), libMesh::ParallelMesh::local_not_level_elements_begin(), libMesh::SerialMesh::local_not_level_elements_end(), libMesh::ParallelMesh::local_not_level_elements_end(), libMesh::MeshRefinement::make_coarsening_compatible(), libMesh::MeshBase::n_active_local_elem(), libMesh::BoundaryInfo::n_boundary_conds(), libMesh::BoundaryInfo::n_edge_conds(), libMesh::DofMap::n_local_dofs(), libMesh::System::n_local_dofs(), libMesh::MeshBase::n_local_elem(), libMesh::MeshBase::n_local_nodes(), libMesh::BoundaryInfo::n_nodeset_conds(), libMesh::SerialMesh::not_local_elements_begin(), libMesh::ParallelMesh::not_local_elements_begin(), libMesh::SerialMesh::not_local_elements_end(), libMesh::ParallelMesh::not_local_elements_end(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::SparsityPattern::Build::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::MeshFunction::operator()(), libMesh::ParallelMesh::ParallelMesh(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::System::point_value(), libMesh::System::project_vector(), libMesh::Nemesis_IO_Helper::put_cmap_params(), libMesh::Nemesis_IO_Helper::put_elem_cmap(), libMesh::Nemesis_IO_Helper::put_elem_map(), libMesh::Nemesis_IO_Helper::put_loadbal_param(), libMesh::Nemesis_IO_Helper::put_node_cmap(), libMesh::Nemesis_IO_Helper::put_node_map(), libMesh::Nemesis_IO::read(), libMesh::CheckpointIO::read(), libMesh::XdrIO::read(), libMesh::UnstructuredMesh::read(), libMesh::CheckpointIO::read_connectivity(), libMesh::ExodusII_IO_Helper::read_elem_num_map(), libMesh::System::read_header(), libMesh::System::read_legacy_data(), libMesh::ExodusII_IO_Helper::read_node_num_map(), libMesh::System::read_parallel_data(), libMesh::System::read_SCALAR_dofs(), libMesh::XdrIO::read_serialized_bc_names(), libMesh::XdrIO::read_serialized_bcs(), libMesh::System::read_serialized_blocked_dof_objects(), libMesh::XdrIO::read_serialized_connectivity(), libMesh::System::read_serialized_data(), libMesh::XdrIO::read_serialized_nodes(), libMesh::XdrIO::read_serialized_nodesets(), libMesh::XdrIO::read_serialized_subdomain_names(), libMesh::System::read_serialized_vector(), libMesh::System::read_serialized_vectors(), libMesh::MeshData::read_xdr(), libMesh::Partitioner::set_node_processor_ids(), libMesh::DofMap::set_nonlocal_dof_objects(), libMesh::LaplaceMeshSmoother::smooth(), libMesh::BoundaryInfo::sync(), libMesh::MeshTools::total_weight(), libMesh::ParallelMesh::update_parallel_id_counts(), libMesh::MeshTools::weight(), libMesh::ExodusII_IO::write(), libMesh::CheckpointIO::write(), libMesh::XdrIO::write(), libMesh::UnstructuredMesh::write(), libMesh::EquationSystems::write(), libMesh::GMVIO::write_discontinuous_gmv(), libMesh::ExodusII_IO::write_element_data(), libMesh::ExodusII_IO_Helper::write_element_values(), libMesh::ExodusII_IO_Helper::write_elements(), libMesh::ExodusII_IO_Helper::write_elements_discontinuous(), libMesh::ExodusII_IO::write_global_data(), libMesh::ExodusII_IO_Helper::write_global_values(), libMesh::System::write_header(), libMesh::ExodusII_IO::write_information_records(), libMesh::ExodusII_IO_Helper::write_information_records(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates(), libMesh::ExodusII_IO_Helper::write_nodal_coordinates_discontinuous(), libMesh::UCDIO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data(), libMesh::ExodusII_IO::write_nodal_data_discontinuous(), libMesh::ExodusII_IO_Helper::write_nodal_values(), libMesh::ExodusII_IO_Helper::write_nodesets(), libMesh::Nemesis_IO_Helper::write_nodesets(), libMesh::System::write_parallel_data(), libMesh::System::write_SCALAR_dofs(), libMesh::XdrIO::write_serialized_bc_names(), libMesh::XdrIO::write_serialized_bcs(), libMesh::System::write_serialized_blocked_dof_objects(), libMesh::XdrIO::write_serialized_connectivity(), libMesh::System::write_serialized_data(), libMesh::XdrIO::write_serialized_nodes(), libMesh::XdrIO::write_serialized_nodesets(), libMesh::XdrIO::write_serialized_subdomain_names(), libMesh::System::write_serialized_vector(), libMesh::System::write_serialized_vectors(), libMesh::ExodusII_IO_Helper::write_sidesets(), libMesh::Nemesis_IO_Helper::write_sidesets(), libMesh::ExodusII_IO::write_timestep(), and libMesh::ExodusII_IO_Helper::write_timestep().

99  { return libmesh_cast_int<processor_id_type>(_communicator.rank()); }
template<typename T >
void libMesh::PetscVector< T >::reciprocal ( )
virtual

Replace each entry v_i of this vector by its reciprocal, 1/v_i.

Implements libMesh::NumericVector< T >.

Definition at line 167 of file petsc_vector.C.

References libMesh::ierr.

168 {
169  PetscErrorCode ierr = 0;
170 
171  // VecReciprocal has been in PETSc since at least 2.3.3 days
172  ierr = VecReciprocal(_vec);
173  LIBMESH_CHKERRABORT(ierr);
174 }
template<typename T >
void libMesh::PetscVector< T >::scale ( const T  factor)
virtual

Scale each element of the vector by the given factor.

Implements libMesh::NumericVector< T >.

Definition at line 511 of file petsc_vector.C.

References libMeshEnums::GHOSTED.

512 {
513  this->_restore_array();
514 
515  PetscErrorCode ierr = 0;
516  PetscScalar factor = static_cast<PetscScalar>(factor_in);
517 
518  if(this->type() != GHOSTED)
519  {
520 #if PETSC_VERSION_LESS_THAN(2,3,0)
521  // 2.2.x & earlier style
522  ierr = VecScale(&factor, _vec);
523  LIBMESH_CHKERRABORT(ierr);
524 #else
525  // 2.3.x & later style
526  ierr = VecScale(_vec, factor);
527  LIBMESH_CHKERRABORT(ierr);
528 #endif
529  }
530  else
531  {
532  Vec loc_vec;
533  ierr = VecGhostGetLocalForm (_vec,&loc_vec);
534  LIBMESH_CHKERRABORT(ierr);
535 
536 #if PETSC_VERSION_LESS_THAN(2,3,0)
537  // 2.2.x & earlier style
538  ierr = VecScale(&factor, loc_vec);
539  LIBMESH_CHKERRABORT(ierr);
540 #else
541  // 2.3.x & later style
542  ierr = VecScale(loc_vec, factor);
543  LIBMESH_CHKERRABORT(ierr);
544 #endif
545 
546  ierr = VecGhostRestoreLocalForm (_vec,&loc_vec);
547  LIBMESH_CHKERRABORT(ierr);
548  }
549 }
template<typename T >
void libMesh::PetscVector< T >::set ( const numeric_index_type  i,
const T  value 
)
virtual

v(i) = value

Implements libMesh::NumericVector< T >.

Definition at line 149 of file petsc_vector.C.

References libMesh::ierr.

150 {
151  this->_restore_array();
152  libmesh_assert_less (i, size());
153 
155  PetscInt i_val = static_cast<PetscInt>(i);
156  PetscScalar petsc_value = static_cast<PetscScalar>(value);
157 
158  ierr = VecSetValues (_vec, 1, &i_val, &petsc_value, INSERT_VALUES);
159  LIBMESH_CHKERRABORT(ierr);
160 
161  this->_is_closed = false;
162 }
template<typename T >
numeric_index_type libMesh::PetscVector< T >::size ( ) const
inlinevirtual
Returns
dimension of the vector. This function was formerly called n(), but was renamed to get the PetscVector<T> class closer to the C++ standard library's std::vector container.

Implements libMesh::NumericVector< T >.

Definition at line 1087 of file petsc_vector.h.

References libMesh::initialized(), and libMesh::libmesh_assert().

Referenced by libMesh::PetscVector< T >::add(), libMesh::PetscVector< T >::init(), libMesh::PetscVector< T >::localize(), and libMesh::PetscVector< T >::operator=().

1088 {
1089  libmesh_assert (this->initialized());
1090 
1091  PetscErrorCode ierr=0;
1092  PetscInt petsc_size=0;
1093 
1094  if (!this->initialized())
1095  return 0;
1096 
1097  ierr = VecGetSize(_vec, &petsc_size);
1098  LIBMESH_CHKERRABORT(ierr);
1099 
1100  return static_cast<numeric_index_type>(petsc_size);
1101 }
template<class T >
Real libMesh::NumericVector< T >::subset_l1_norm ( const std::set< numeric_index_type > &  indices) const
virtualinherited
Returns
the $l_1$-norm of the vector, i.e. the sum of the absolute values for the specified entries in the vector.

Note that the indices must necessary live on this processor.

Definition at line 320 of file numeric_vector.C.

References std::abs(), libMesh::comm, and libMesh::Real.

Referenced by libMesh::System::discrete_var_norm().

321 {
322  const NumericVector<T> & v = *this;
323 
324  std::set<numeric_index_type>::const_iterator it = indices.begin();
325  const std::set<numeric_index_type>::const_iterator it_end = indices.end();
326 
327  Real norm = 0;
328 
329  for(; it!=it_end; ++it)
330  norm += std::abs(v(*it));
331 
332  this->comm().sum(norm);
333 
334  return norm;
335 }
template<class T >
Real libMesh::NumericVector< T >::subset_l2_norm ( const std::set< numeric_index_type > &  indices) const
virtualinherited
Returns
the $l_2$-norm of the vector, i.e. the square root of the sum of the squares of the elements for the specified entries in the vector.

Note that the indices must necessary live on this processor.

Definition at line 338 of file numeric_vector.C.

References libMesh::comm, libMesh::TensorTools::norm_sq(), and libMesh::Real.

Referenced by libMesh::System::discrete_var_norm().

339 {
340  const NumericVector<T> & v = *this;
341 
342  std::set<numeric_index_type>::const_iterator it = indices.begin();
343  const std::set<numeric_index_type>::const_iterator it_end = indices.end();
344 
345  Real norm = 0;
346 
347  for(; it!=it_end; ++it)
348  norm += TensorTools::norm_sq(v(*it));
349 
350  this->comm().sum(norm);
351 
352  return std::sqrt(norm);
353 }
template<class T >
Real libMesh::NumericVector< T >::subset_linfty_norm ( const std::set< numeric_index_type > &  indices) const
virtualinherited
Returns
the maximum absolute value of the specified entries of this vector, which is the $l_\infty$-norm of a vector.

Note that the indices must necessary live on this processor.

Definition at line 356 of file numeric_vector.C.

References std::abs(), libMesh::comm, and libMesh::Real.

Referenced by libMesh::System::discrete_var_norm().

357 {
358  const NumericVector<T> & v = *this;
359 
360  std::set<numeric_index_type>::const_iterator it = indices.begin();
361  const std::set<numeric_index_type>::const_iterator it_end = indices.end();
362 
363  Real norm = 0;
364 
365  for(; it!=it_end; ++it)
366  {
367  Real value = std::abs(v(*it));
368  if(value > norm)
369  norm = value;
370  }
371 
372  this->comm().max(norm);
373 
374  return norm;
375 }
template<typename T >
T libMesh::PetscVector< T >::sum ( ) const
virtual
Returns
the sum of values in a vector

Implements libMesh::NumericVector< T >.

Definition at line 52 of file petsc_vector.C.

References libMesh::closed(), libMesh::ierr, and libMesh::libmesh_assert().

53 {
54  this->_restore_array();
55  libmesh_assert(this->closed());
56 
58  PetscScalar value=0.;
59 
60  ierr = VecSum (_vec, &value);
61  LIBMESH_CHKERRABORT(ierr);
62 
63  return static_cast<T>(value);
64 }
template<typename T >
void libMesh::PetscVector< T >::swap ( NumericVector< T > &  v)
inlinevirtual

Swaps the raw PETSc vector context pointers.

Reimplemented from libMesh::NumericVector< T >.

Definition at line 1290 of file petsc_vector.h.

References libMesh::PetscVector< T >::_array_is_present, libMesh::PetscVector< T >::_destroy_vec_on_exit, libMesh::PetscVector< T >::_global_to_local_map, libMesh::PetscVector< T >::_local_form, libMesh::PetscVector< T >::_values, libMesh::PetscVector< T >::_vec, swap(), and libMesh::NumericVector< T >::swap().

Referenced by libMesh::__libmesh_petsc_diff_solver_residual(), libMesh::__libmesh_petsc_snes_jacobian(), libMesh::__libmesh_petsc_snes_residual(), DMlibMeshFunction(), and DMlibMeshJacobian().

1291 {
1292  NumericVector<T>::swap(other);
1293 
1294  PetscVector<T>& v = libmesh_cast_ref<PetscVector<T>&>(other);
1295 
1296  std::swap(_vec, v._vec);
1297  std::swap(_destroy_vec_on_exit, v._destroy_vec_on_exit);
1298  std::swap(_global_to_local_map, v._global_to_local_map);
1299  std::swap(_array_is_present, v._array_is_present);
1300  std::swap(_local_form, v._local_form);
1301  std::swap(_values, v._values);
1302 }
template<typename T>
ParallelType& libMesh::NumericVector< T >::type ( )
inlineinherited
Returns
the type (SERIAL, PARALLEL, GHOSTED) of the vector.

Definition at line 139 of file numeric_vector.h.

139 { return _type; }
template<typename T>
Vec libMesh::PetscVector< T >::vec ( )
inline

Returns the raw PETSc vector context pointer. Note this is generally not required in user-level code. Just don't do anything crazy like calling LibMeshVecDestroy()!

Definition at line 540 of file petsc_vector.h.

References libMesh::PetscVector< T >::_vec, and libMesh::libmesh_assert().

Referenced by libMesh::PetscPreconditioner< T >::apply(), libMesh::PetscNonlinearSolver< T >::build_mat_null_space(), DMCreateGlobalVector_libMesh(), libMesh::PetscMatrix< T >::get_diagonal(), libMesh::SlepcEigenSolver< T >::get_eigenpair(), libMesh::PetscDiffSolver::solve(), and libMesh::PetscLinearSolver< T >::solve().

540 { libmesh_assert (_vec); return _vec; }
template<typename T >
void libMesh::PetscVector< T >::zero ( )
inlinevirtual

Set all entries to zero. Equivalent to v = 0, but more obvious and faster.

Implements libMesh::NumericVector< T >.

Definition at line 1010 of file petsc_vector.h.

References libMesh::closed(), libMeshEnums::GHOSTED, and libMesh::libmesh_assert().

1011 {
1012  libmesh_assert(this->closed());
1013 
1014  this->_restore_array();
1015 
1016  PetscErrorCode ierr=0;
1017 
1018  PetscScalar z=0.;
1019 
1020  if(this->type() != GHOSTED)
1021  {
1022 #if PETSC_VERSION_LESS_THAN(2,3,0)
1023  // 2.2.x & earlier style
1024  ierr = VecSet (&z, _vec);
1025  LIBMESH_CHKERRABORT(ierr);
1026 #else
1027  // 2.3.x & newer
1028  ierr = VecSet (_vec, z);
1029  LIBMESH_CHKERRABORT(ierr);
1030 #endif
1031  }
1032  else
1033  {
1034  /* Vectors that include ghost values require a special
1035  handling. */
1036  Vec loc_vec;
1037  ierr = VecGhostGetLocalForm (_vec,&loc_vec);
1038  LIBMESH_CHKERRABORT(ierr);
1039 #if PETSC_VERSION_LESS_THAN(2,3,0)
1040  // 2.2.x & earlier style
1041  ierr = VecSet (&z, loc_vec);
1042  LIBMESH_CHKERRABORT(ierr);
1043 #else
1044  // 2.3.x & newer
1045  ierr = VecSet (loc_vec, z);
1046  LIBMESH_CHKERRABORT(ierr);
1047 #endif
1048  ierr = VecGhostRestoreLocalForm (_vec,&loc_vec);
1049  LIBMESH_CHKERRABORT(ierr);
1050  }
1051 }
template<typename T >
AutoPtr< NumericVector< T > > libMesh::PetscVector< T >::zero_clone ( ) const
inlinevirtual

Creates a vector which has the same type, size and partitioning as this vector, but whose data is all zero. Returns it in an AutoPtr.

Implements libMesh::NumericVector< T >.

Definition at line 1057 of file petsc_vector.h.

References libMesh::comm.

1058 {
1059  AutoPtr<NumericVector<T> > cloned_vector
1060  (new PetscVector<T>(this->comm(), this->type()));
1061 
1062  cloned_vector->init(*this);
1063 
1064  return cloned_vector;
1065 }

Member Data Documentation

template<typename T>
bool libMesh::PetscVector< T >::_array_is_present
mutableprivate

If true, the actual Petsc array of the values of the vector is currently accessible. That means that the members _local_form and _values are valid.

Definition at line 557 of file petsc_vector.h.

Referenced by libMesh::PetscVector< T >::swap().

ReferenceCounter::Counts libMesh::ReferenceCounter::_counts
staticprotectedinherited
template<typename T>
bool libMesh::PetscVector< T >::_destroy_vec_on_exit
private

This boolean value should only be set to false for the constructor which takes a PETSc Vec object.

Definition at line 608 of file petsc_vector.h.

Referenced by libMesh::PetscVector< T >::swap().

bool libMesh::ReferenceCounter::_enable_print_counter = true
staticprotectedinherited

Flag to control whether reference count information is printed when print_info is called.

Definition at line 137 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::disable_print_counter_info(), libMesh::ReferenceCounter::enable_print_counter_info(), and libMesh::ReferenceCounter::print_info().

template<typename T>
GlobalToLocalMap libMesh::PetscVector< T >::_global_to_local_map
private

Map that maps global to local ghost cells (will be empty if not in ghost cell mode).

Definition at line 602 of file petsc_vector.h.

Referenced by libMesh::PetscVector< T >::init(), libMesh::PetscVector< T >::operator=(), and libMesh::PetscVector< T >::swap().

template<typename T>
bool libMesh::NumericVector< T >::_is_closed
protectedinherited

Flag to see if the Numeric assemble routines have been called yet

Definition at line 657 of file numeric_vector.h.

Referenced by libMesh::NumericVector< Number >::closed(), libMesh::DistributedVector< T >::localize(), and libMesh::DistributedVector< T >::operator=().

template<typename T>
bool libMesh::NumericVector< T >::_is_initialized
protectedinherited
template<typename T>
Vec libMesh::PetscVector< T >::_local_form
mutableprivate

Petsc vector datatype to hold the local form of a ghosted vector. The contents of this field are only valid if the vector is ghosted and _array_is_present is true.

Definition at line 573 of file petsc_vector.h.

Referenced by libMesh::PetscVector< T >::swap().

template<typename T>
numeric_index_type libMesh::PetscVector< T >::_local_size
mutableprivate

Size of the local form, for being used in assertations. The contents of this field are only valid if the vector is ghosted and _array_is_present is true.

Definition at line 565 of file petsc_vector.h.

Threads::spin_mutex libMesh::ReferenceCounter::_mutex
staticprotectedinherited

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 131 of file reference_counter.h.

Threads::atomic< unsigned int > libMesh::ReferenceCounter::_n_objects
staticprotectedinherited

The number of objects. Print the reference count information when the number returns to 0.

Definition at line 126 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::n_objects(), libMesh::ReferenceCounter::ReferenceCounter(), and libMesh::ReferenceCounter::~ReferenceCounter().

template<typename T>
PetscScalar* libMesh::PetscVector< T >::_values
mutableprivate

Pointer to the actual Petsc array of the values of the vector. This pointer is only valid if _array_is_present is true.

Definition at line 579 of file petsc_vector.h.

Referenced by libMesh::PetscVector< T >::swap().


The documentation for this class was generated from the following files:

Site Created By: libMesh Developers
Last modified: February 07 2014 16:58:01 UTC

Hosted By:
SourceForge.net Logo