20 #ifndef LIBMESH_NUMERIC_VECTOR_H 21 #define LIBMESH_NUMERIC_VECTOR_H 24 #include "libmesh/libmesh_common.h" 25 #include "libmesh/enum_parallel_type.h" 26 #include "libmesh/id_types.h" 27 #include "libmesh/int_range.h" 28 #include "libmesh/reference_counted_object.h" 29 #include "libmesh/libmesh.h" 30 #include "libmesh/parallel_object.h" 31 #include "libmesh/dense_subvector.h" 32 #include "libmesh/dense_vector.h" 46 template <
typename T>
class NumericVector;
47 template <
typename T>
class DenseVector;
48 template <
typename T>
class DenseSubVector;
49 template <
typename T>
class SparseMatrix;
50 template <
typename T>
class ShellMatrix;
68 class NumericVector :
public ReferenceCountedObject<NumericVector<T>>,
105 const std::vector<numeric_index_type> & ghost,
119 virtual NumericVector<T> &
operator= (
const NumericVector<T> & v) = 0;
140 static std::unique_ptr<NumericVector<T>>
141 build(
const Parallel::Communicator &
comm,
163 #ifdef LIBMESH_ENABLE_DEPRECATED 188 virtual void close () = 0;
193 virtual void clear ();
199 virtual void zero () = 0;
207 virtual std::unique_ptr<NumericVector<T>>
zero_clone ()
const = 0;
214 virtual std::unique_ptr<NumericVector<T>>
clone ()
const = 0;
228 const bool fast =
false,
235 const bool fast =
false,
244 const std::vector<numeric_index_type> & ghost,
245 const bool fast =
false,
253 const bool fast =
false) = 0;
284 virtual T
sum()
const = 0;
383 virtual void get(
const std::vector<numeric_index_type> & index,
392 void get(
const std::vector<numeric_index_type> & index,
393 std::vector<T> & values)
const;
475 virtual void add (
const T s) = 0;
499 const std::vector<numeric_index_type> & dof_indices);
508 const std::vector<numeric_index_type> & dof_indices);
517 const std::vector<numeric_index_type> & dof_indices);
526 const std::vector<numeric_index_type> & dof_indices);
556 virtual void insert (
const T * v,
557 const std::vector<numeric_index_type> & dof_indices);
563 void insert (
const std::vector<T> & v,
564 const std::vector<numeric_index_type> & dof_indices);
571 const std::vector<numeric_index_type> & dof_indices);
578 const std::vector<numeric_index_type> & dof_indices);
585 const std::vector<numeric_index_type> & dof_indices);
590 virtual void scale (
const T factor) = 0;
595 virtual void abs() = 0;
609 virtual void localize (std::vector<T> & v_local)
const = 0;
622 const std::vector<numeric_index_type> & send_list)
const = 0;
653 virtual void localize (std::vector<T> & v_local,
654 const std::vector<numeric_index_type> & indices)
const = 0;
662 const std::vector<numeric_index_type> & send_list) = 0;
727 friend std::ostream & operator << (std::ostream & os, const NumericVector<T> & v)
738 virtual void print_matlab(
const std::string & filename =
"")
const;
751 virtual void read_matlab(
const std::string & filename);
763 const std::vector<numeric_index_type> &,
766 libmesh_not_implemented();
775 virtual std::unique_ptr<NumericVector<T>>
778 libmesh_not_implemented();
787 const std::vector<numeric_index_type> &)
789 libmesh_not_implemented();
854 template <
typename T>
867 template <
typename T>
877 libmesh_not_implemented();
883 template <
typename T>
894 libmesh_not_implemented();
900 template <
typename T>
905 const std::vector<numeric_index_type> & ,
912 libmesh_not_implemented();
918 template <
typename T>
928 template <
typename T>
933 const std::size_t num = index.size();
934 for (std::size_t i=0; i<num; i++)
936 values[i] = (*this)(index[i]);
942 template <
typename T>
945 std::vector<T> & values)
const 947 const std::size_t num = index.size();
952 this->
get(index, values.data());
957 template <
typename T>
960 const std::vector<numeric_index_type> & dof_indices)
964 this->add_vector(v.data(), dof_indices);
969 template <
typename T>
972 const std::vector<numeric_index_type> & dof_indices)
976 this->add_vector(&v(0), dof_indices);
981 template <
typename T>
984 const std::vector<numeric_index_type> & dof_indices)
988 this->insert(v.data(), dof_indices);
993 template <
typename T>
996 const std::vector<numeric_index_type> & dof_indices)
1000 this->insert(&v(0), dof_indices);
1005 template <
typename T>
1008 const std::vector<numeric_index_type> & dof_indices)
1012 this->insert(&v(0), dof_indices);
1025 os <<
"Size\tglobal = " << this->size()
1026 <<
"\t\tlocal = " << this->local_size() << std::endl;
1029 os <<
"#\tReal part\t\tImaginary part" << std::endl;
1032 << (*this)(i).
real() <<
"\t\t" 1033 << (*this)(i).
imag() << std::endl;
1038 template <
typename T>
1043 os <<
"Size\tglobal = " << this->size()
1044 <<
"\t\tlocal = " << this->local_size() << std::endl;
1046 os <<
"#\tValue" << std::endl;
1048 os << i <<
"\t" << (*this)(i) << std::endl;
1059 std::vector<Complex> v(this->size());
1063 if (this->processor_id())
1066 os <<
"Size\tglobal = " << this->size() << std::endl;
1067 os <<
"#\tReal part\t\tImaginary part" << std::endl;
1070 << v[i].real() <<
"\t\t" 1071 << v[i].imag() << std::endl;
1075 template <
typename T>
1081 std::vector<T> v(this->size());
1085 if (this->processor_id())
1088 os <<
"Size\tglobal = " << this->size() << std::endl;
1089 os <<
"#\tValue" << std::endl;
1091 os << i <<
"\t" << v[i] << std::endl;
1096 template <
typename T>
1100 std::swap(_is_closed, v._is_closed);
1102 std::swap(_type, v._type);
1105 template <
typename T>
1109 return vec.l1_norm();
1112 template <
typename T>
1116 return vec1.l1_norm_diff(vec2);
1123 #ifdef LIBMESH_DEFAULT_QUADRUPLE_PRECISION 1124 #include <boost/mpl/bool.hpp> 1126 namespace boost {
namespace multiprecision {
namespace detail {
1127 template <
typename T,
typename To>
1128 struct is_lossy_conversion<
libMesh::NumericVector<T>, To> {
1136 #endif // LIBMESH_NUMERIC_VECTOR_H virtual void insert(const T *v, const std::vector< numeric_index_type > &dof_indices)
Inserts the entries of v in *this at the locations specified by v.
virtual Real subset_l2_norm(const std::set< numeric_index_type > &indices) const
std::mutex _numeric_vector_mutex
Mutex for performing thread-safe operations.
virtual void print_global(std::ostream &os=libMesh::out) const
Prints the global contents of the vector, by default to libMesh::out.
boost::multiprecision::float128 real(const boost::multiprecision::float128 in)
virtual T el(const numeric_index_type i) const
static constexpr Real TOLERANCE
virtual void add_vector_transpose(const NumericVector< T > &v, const SparseMatrix< T > &A)=0
Computes , i.e.
virtual bool initialized() const
virtual numeric_index_type size() const =0
virtual std::unique_ptr< NumericVector< T > > zero_clone() const =0
virtual void add_vector(const T *v, const std::vector< numeric_index_type > &dof_indices)
Computes , where v is a pointer and each dof_indices[i] specifies where to add value v[i]...
virtual unsigned int size() const override final
virtual NumericVector< T > & operator=(const NumericVector< T > &v)=0
This looks like a copy assignment operator, but note that, unlike normal copy assignment operators...
Provides a uniform interface to vector storage schemes for different linear algebra libraries...
const Parallel::Communicator & comm() const
virtual std::unique_ptr< NumericVector< T > > clone() const =0
bool _is_initialized
true once init() has been called.
The libMesh namespace provides an interface to certain functionality in the library.
virtual Real subset_linfty_norm(const std::set< numeric_index_type > &indices) const
virtual T dot(const NumericVector< T > &v) const =0
virtual void init(const numeric_index_type n, const numeric_index_type n_local, const bool fast=false, const ParallelType ptype=AUTOMATIC)=0
Change the dimension of the vector to n.
virtual bool empty() const override final
virtual Real subset_l1_norm(const std::set< numeric_index_type > &indices) const
uint8_t processor_id_type
virtual void zero()=0
Set all entries to zero.
Defines a dense subvector for use in finite element computations.
virtual void create_subvector(NumericVector< T > &, const std::vector< numeric_index_type > &, bool=true) const
Fills in subvector from this vector using the indices in rows.
Real l2_norm_diff(const NumericVector< T > &other_vec) const
SolverPackage default_solver_package()
virtual void scale(const T factor)=0
Scale each element of the vector by the given factor.
virtual void localize_to_one(std::vector< T > &v_local, const processor_id_type proc_id=0) const =0
Creates a local copy of the global vector in v_local only on processor proc_id.
virtual Real l2_norm() const =0
virtual NumericVector< T > & operator+=(const NumericVector< T > &v)=0
Adds v to *this, .
virtual ~NumericVector()=default
While this class doesn't manage any memory, the derived class might and users may be deleting through...
auto l1_norm(const NumericVector< T > &vec)
dof_id_type numeric_index_type
bool _is_initialized
Flag that tells if init() has been called.
NumericVector(const Parallel::Communicator &comm_in, const ParallelType ptype=AUTOMATIC)
Dummy-Constructor.
NumericVector< T > & operator/=(const T a)
Scales the vector by 1/a, .
auto l1_norm_diff(const NumericVector< T > &vec1, const NumericVector< T > &vec2)
NumericVector< T > & operator*=(const T a)
Scales the vector by a, .
virtual bool empty() const override final
virtual void read_matlab(const std::string &filename)
Read the contents of the vector from the Matlab-script format used by PETSc.
bool compatible(const NumericVector< T > &v) 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
Prints the local contents of the vector, by default to libMesh::out.
virtual Real max() const =0
virtual Real min() const =0
virtual std::unique_ptr< NumericVector< T > > get_subvector(const std::vector< numeric_index_type > &)
Creates a view into this vector using the indices in rows.
An object whose state is distributed along a set of processors.
virtual Real l1_norm() const =0
virtual void close()=0
Calls the NumericVector's internal assembly routines, ensuring that the values are consistent across ...
ParallelType _type
Type of vector.
virtual void pointwise_mult(const NumericVector< T > &vec1, const NumericVector< T > &vec2)=0
Computes (summation not implied) i.e.
virtual bool closed() const
Real l1_norm_diff(const NumericVector< T > &other_vec) const
virtual void abs()=0
Sets for each entry in the vector.
virtual numeric_index_type first_local_index() const =0
ParallelType type() const
virtual int local_relative_compare(const NumericVector< T > &other_vector, const Real threshold=TOLERANCE) const
virtual void print_matlab(const std::string &filename="") const
Print the contents of the vector in Matlab's sparse matrix format.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual NumericVector< T > & operator-=(const NumericVector< T > &v)=0
Subtracts v from *this, .
virtual void swap(NumericVector< T > &v)
Swaps the contents of this with v.
virtual int compare(const NumericVector< T > &other_vector, const Real threshold=TOLERANCE) const
virtual numeric_index_type local_size() const =0
virtual std::size_t max_allowed_id() const =0
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
static std::unique_ptr< NumericVector< T > > build(const Parallel::Communicator &comm, SolverPackage solver_package=libMesh::default_solver_package(), ParallelType parallel_type=AUTOMATIC)
Builds a NumericVector on the processors in communicator comm using the linear solver package specifi...
virtual unsigned int size() const override final
bool initialized()
Checks that library initialization has been done.
virtual void clear()
Restores the NumericVector<T> to a pristine state.
Defines a dense vector for use in Finite Element-type computations.
boost::multiprecision::float128 imag(const boost::multiprecision::float128)
virtual void reciprocal()=0
Computes the component-wise reciprocal, .
SolverPackage
Defines an enum for various linear solver packages.
Generic shell matrix, i.e.
virtual void add(const numeric_index_type i, const T value)=0
Adds value to the vector entry specified by i.
virtual T operator()(const numeric_index_type i) const =0
virtual numeric_index_type last_local_index() const =0
virtual Real linfty_norm() const =0
void ErrorVector unsigned int
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
void set_type(ParallelType t)
Allow the user to change the ParallelType of the NumericVector under some circumstances.
const Elem & get(const ElemType type_in)
virtual void pointwise_divide(const NumericVector< T > &vec1, const NumericVector< T > &vec2)=0
Computes (summation not implied) i.e.
bool _is_closed
Flag which tracks whether the vector's values are consistent on all processors after insertion or add...
virtual void conjugate()=0
Negates the imaginary component of each entry in the vector.
ParallelType
Defines an enum for parallel data structure types.
virtual void restore_subvector(std::unique_ptr< NumericVector< T >>, const std::vector< numeric_index_type > &)
Restores a view into this vector using the indices in rows.
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.