MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
vector.hpp
Go to the documentation of this file.
1// Copyright (c) 2010-2024, Lawrence Livermore National Security, LLC. Produced
2// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3// LICENSE and NOTICE for details. LLNL-CODE-806117.
4//
5// This file is part of the MFEM library. For more information and source code
6// availability visit https://mfem.org.
7//
8// MFEM is free software; you can redistribute it and/or modify it under the
9// terms of the BSD-3 license. We welcome feedback and contributions, see file
10// CONTRIBUTING.md for details.
11
12#ifndef MFEM_VECTOR
13#define MFEM_VECTOR
14
15#include "../general/array.hpp"
16#ifdef MFEM_USE_ADIOS2
18#endif
21#include "../general/device.hpp"
22
23#include <cmath>
24#include <cstdlib>
25#include <iostream>
26#include <limits>
27#if defined(_MSC_VER) && (_MSC_VER < 1800)
28#include <float.h>
29#define isfinite _finite
30#endif
31
32#ifdef MFEM_USE_MPI
33#include <mpi.h>
34#endif
35
36namespace mfem
37{
38
39/** Count the number of entries in an array of doubles for which isfinite
40 is false, i.e. the entry is a NaN or +/-Inf. */
41inline int CheckFinite(const real_t *v, const int n);
42
43/// Define a shortcut for std::numeric_limits<double>::infinity()
44#ifndef __CYGWIN__
46{
47 return std::numeric_limits<real_t>::infinity();
48}
49#else
50// On Cygwin math.h defines a function 'infinity()' which will conflict with the
51// above definition if we have 'using namespace mfem;' and try to use something
52// like 'double a = infinity();'. This 'infinity()' function is non-standard and
53// is defined by the Newlib C standard library implementation used by Cygwin,
54// see https://en.wikipedia.org/wiki/Newlib, http://www.sourceware.org/newlib.
55using ::infinity;
56#endif
57
58/// Generate a random `real_t` number in the interval [0,1) using rand().
60{
61 constexpr real_t max = (real_t)(RAND_MAX) + 1_r;
62#if defined(MFEM_USE_SINGLE)
63 // Note: For RAND_MAX = 2^31-1, float(RAND_MAX) = 2^31, and
64 // max = float(RAND_MAX)+1.0f is equal to 2^31 too; this is the actual value
65 // that we want for max. However, this rounding behavior means that for a
66 // range of values of rand() close to RAND_MAX, the expression
67 // float(rand())/max will give 1.0f.
68 //
69 // Therefore, to ensure we return a number less than 1, we take the minimum
70 // of float(rand())/max and (1 - 2^(-24)), where the latter number is the
71 // largest float less than 1.
72 return std::fmin(real_t(rand())/max, 0.9999999403953552_r);
73#else
74 return real_t(rand())/max;
75#endif
76}
77
78/// Vector data type.
79class Vector
80{
81protected:
82
84 int size;
85
86public:
87
88 /** Default constructor for Vector. Sets size = 0, and calls Memory::Reset on
89 data through Memory<double>'s default constructor. */
90 Vector(): size(0) { }
91
92 /// Copy constructor. Allocates a new data array and copies the data.
93 Vector(const Vector &);
94
95 /// Move constructor. "Steals" data from its argument.
96 Vector(Vector&& v);
97
98 /// @brief Creates vector of size s.
99 /// @warning Entries are not initialized to zero!
100 explicit Vector(int s);
101
102 /// Creates a vector referencing an array of doubles, owned by someone else.
103 /** The pointer @a data_ can be NULL. The data array can be replaced later
104 with SetData(). */
105 Vector(real_t *data_, int size_)
106 { data.Wrap(data_, size_, false); size = size_; }
107
108 /** @brief Create a Vector referencing a sub-vector of the Vector @a base
109 starting at the given offset, @a base_offset, and size @a size_. */
110 Vector(Vector &base, int base_offset, int size_)
111 : data(base.data, base_offset, size_), size(size_) { }
112
113 /// Create a Vector of size @a size_ using MemoryType @a mt.
114 Vector(int size_, MemoryType mt)
115 : data(size_, mt), size(size_) { }
116
117 /** @brief Create a Vector of size @a size_ using host MemoryType @a h_mt and
118 device MemoryType @a d_mt. */
119 Vector(int size_, MemoryType h_mt, MemoryType d_mt)
120 : data(size_, h_mt, d_mt), size(size_) { }
121
122 /// Create a vector using a braced initializer list
123 template <int N, typename T = real_t>
124 explicit Vector(const T (&values)[N]) : Vector(N)
125 { std::copy(values, values + N, GetData()); }
126
127 /// Enable execution of Vector operations using the mfem::Device.
128 /** The default is to use Backend::CPU (serial execution on each MPI rank),
129 regardless of the mfem::Device configuration.
130
131 When appropriate, MFEM functions and class methods will enable the use
132 of the mfem::Device for their Vector parameters.
133
134 Some derived classes, e.g. GridFunction, enable the use of the
135 mfem::Device by default. */
136 virtual void UseDevice(bool use_dev) const { data.UseDevice(use_dev); }
137
138 /// Return the device flag of the Memory object used by the Vector
139 virtual bool UseDevice() const { return data.UseDevice(); }
140
141 /// Reads a vector from multiple files
142 void Load(std::istream ** in, int np, int * dim);
143
144 /// Load a vector from an input stream.
145 void Load(std::istream &in, int Size);
146
147 /// Load a vector from an input stream, reading the size from the stream.
148 void Load(std::istream &in) { int s; in >> s; Load(in, s); }
149
150 /// @brief Resize the vector to size @a s.
151 /** If the new size is less than or equal to Capacity() then the internal
152 data array remains the same. Otherwise, the old array is deleted, if
153 owned, and a new array of size @a s is allocated without copying the
154 previous content of the Vector.
155 @warning In the second case above (new size greater than current one),
156 the vector will allocate new data array, even if it did not own the
157 original data! Also, new entries are not initialized! */
158 void SetSize(int s);
159
160 /// Resize the vector to size @a s using MemoryType @a mt.
161 void SetSize(int s, MemoryType mt);
162
163 /// Resize the vector to size @a s using the MemoryType of @a v.
164 void SetSize(int s, const Vector &v) { SetSize(s, v.GetMemory().GetMemoryType()); }
165
166 /// Set the Vector data.
167 /// @warning This method should be called only when OwnsData() is false.
168 void SetData(real_t *d) { data.Wrap(d, data.Capacity(), false); }
169
170 /// Set the Vector data and size.
171 /** The Vector does not assume ownership of the new data. The new size is
172 also used as the new Capacity().
173 @warning This method should be called only when OwnsData() is false.
174 @sa NewDataAndSize(). */
175 void SetDataAndSize(real_t *d, int s) { data.Wrap(d, s, false); size = s; }
176
177 /// Set the Vector data and size, deleting the old data, if owned.
178 /** The Vector does not assume ownership of the new data. The new size is
179 also used as the new Capacity().
180 @sa SetDataAndSize(). */
181 void NewDataAndSize(real_t *d, int s)
182 {
183 data.Delete();
184 SetDataAndSize(d, s);
185 }
186
187 /// Reset the Vector to use the given external Memory @a mem and size @a s.
188 /** If @a own_mem is false, the Vector will not own any of the pointers of
189 @a mem.
190
191 Note that when @a own_mem is true, the @a mem object can be destroyed
192 immediately by the caller but `mem.Delete()` should NOT be called since
193 the Vector object takes ownership of all pointers owned by @a mem.
194
195 @sa NewDataAndSize(). */
196 inline void NewMemoryAndSize(const Memory<real_t> &mem, int s, bool own_mem);
197
198 /// Reset the Vector to be a reference to a sub-vector of @a base.
199 inline void MakeRef(Vector &base, int offset, int size);
200
201 /** @brief Reset the Vector to be a reference to a sub-vector of @a base
202 without changing its current size. */
203 inline void MakeRef(Vector &base, int offset);
204
205 /// Set the Vector data (host pointer) ownership flag.
206 void MakeDataOwner() const { data.SetHostPtrOwner(true); }
207
208 /// Destroy a vector
209 void Destroy();
210
211 /** @brief Delete the device pointer, if owned. If @a copy_to_host is true
212 and the data is valid only on device, move it to host before deleting.
213 Invalidates the device memory. */
214 void DeleteDevice(bool copy_to_host = true)
215 { data.DeleteDevice(copy_to_host); }
216
217 /// Returns the size of the vector.
218 inline int Size() const { return size; }
219
220 /// Return the size of the currently allocated data array.
221 /** It is always true that Capacity() >= Size(). */
222 inline int Capacity() const { return data.Capacity(); }
223
224 /// Return a pointer to the beginning of the Vector data.
225 /** @warning This method should be used with caution as it gives write access
226 to the data of const-qualified Vector%s. */
227 inline real_t *GetData() const
228 { return const_cast<real_t*>((const real_t*)data); }
229
230 /// Conversion to `double *`. Deprecated.
231 MFEM_DEPRECATED inline operator real_t *() { return data; }
232
233 /// Conversion to `const double *`. Deprecated.
234 MFEM_DEPRECATED inline operator const real_t *() const { return data; }
235
236 /// STL-like begin.
237 inline real_t *begin() { return data; }
238
239 /// STL-like end.
240 inline real_t *end() { return data + size; }
241
242 /// STL-like begin (const version).
243 inline const real_t *begin() const { return data; }
244
245 /// STL-like end (const version).
246 inline const real_t *end() const { return data + size; }
247
248 /// Return a reference to the Memory object used by the Vector.
250
251 /** @brief Return a reference to the Memory object used by the Vector, const
252 version. */
253 const Memory<real_t> &GetMemory() const { return data; }
254
255 /// Update the memory location of the vector to match @a v.
256 void SyncMemory(const Vector &v) const { GetMemory().Sync(v.GetMemory()); }
257
258 /// Update the alias memory location of the vector to match @a v.
259 void SyncAliasMemory(const Vector &v) const
260 { GetMemory().SyncAlias(v.GetMemory(),Size()); }
261
262 /// Read the Vector data (host pointer) ownership flag.
263 inline bool OwnsData() const { return data.OwnsHostPtr(); }
264
265 /// Changes the ownership of the data; after the call the Vector is empty
266 inline void StealData(real_t **p)
267 { *p = data; data.Reset(); size = 0; }
268
269 /// Changes the ownership of the data; after the call the Vector is empty
270 inline real_t *StealData() { real_t *p; StealData(&p); return p; }
271
272 /// Access Vector entries. Index i = 0 .. size-1.
273 real_t &Elem(int i);
274
275 /// Read only access to Vector entries. Index i = 0 .. size-1.
276 const real_t &Elem(int i) const;
277
278 /// Access Vector entries using () for 0-based indexing.
279 /** @note If MFEM_DEBUG is enabled, bounds checking is performed. */
280 inline real_t &operator()(int i);
281
282 /// Read only access to Vector entries using () for 0-based indexing.
283 /** @note If MFEM_DEBUG is enabled, bounds checking is performed. */
284 inline const real_t &operator()(int i) const;
285
286 /// Access Vector entries using [] for 0-based indexing.
287 /** @note If MFEM_DEBUG is enabled, bounds checking is performed. */
288 inline real_t &operator[](int i) { return (*this)(i); }
289
290 /// Read only access to Vector entries using [] for 0-based indexing.
291 /** @note If MFEM_DEBUG is enabled, bounds checking is performed. */
292 inline const real_t &operator[](int i) const { return (*this)(i); }
293
294 /// Dot product with a `double *` array.
295 real_t operator*(const real_t *) const;
296
297 /// Return the inner-product.
298 real_t operator*(const Vector &v) const;
299
300 /// Copy Size() entries from @a v.
301 Vector &operator=(const real_t *v);
302
303 /// Copy assignment.
304 /** @note Defining this method overwrites the implicitly defined copy
305 assignment operator. */
306 Vector &operator=(const Vector &v);
307
308 /// Move assignment
309 Vector &operator=(Vector&& v);
310
311 /// Redefine '=' for vector = constant.
312 Vector &operator=(real_t value);
313
315
316 /// Component-wise scaling: (*this)(i) *= v(i)
317 Vector &operator*=(const Vector &v);
318
320
321 /// Component-wise division: (*this)(i) /= v(i)
322 Vector &operator/=(const Vector &v);
323
325
326 Vector &operator-=(const Vector &v);
327
329
330 Vector &operator+=(const Vector &v);
331
332 /// (*this) += a * Va
333 Vector &Add(const real_t a, const Vector &Va);
334
335 /// (*this) = a * x
336 Vector &Set(const real_t a, const Vector &x);
337
338 void SetVector(const Vector &v, int offset);
339
340 void AddSubVector(const Vector &v, int offset);
341
342 /// (*this) = -(*this)
343 void Neg();
344
345 /// (*this)(i) = 1.0 / (*this)(i)
346 void Reciprocal();
347
348 /// Swap the contents of two Vectors
349 inline void Swap(Vector &other);
350
351 /// Set v = v1 + v2.
352 friend void add(const Vector &v1, const Vector &v2, Vector &v);
353
354 /// Set v = v1 + alpha * v2.
355 friend void add(const Vector &v1, real_t alpha, const Vector &v2, Vector &v);
356
357 /// z = a * (x + y)
358 friend void add(const real_t a, const Vector &x, const Vector &y, Vector &z);
359
360 /// z = a * x + b * y
361 friend void add(const real_t a, const Vector &x,
362 const real_t b, const Vector &y, Vector &z);
363
364 /// Set v = v1 - v2.
365 friend void subtract(const Vector &v1, const Vector &v2, Vector &v);
366
367 /// z = a * (x - y)
368 friend void subtract(const real_t a, const Vector &x,
369 const Vector &y, Vector &z);
370
371 /// Computes cross product of this vector with another 3D vector.
372 /// vout = this x vin.
373 void cross3D(const Vector &vin, Vector &vout) const;
374
375 /// v = median(v,lo,hi) entrywise. Implementation assumes lo <= hi.
376 void median(const Vector &lo, const Vector &hi);
377
378 /// Extract entries listed in @a dofs to the output Vector @a elemvect.
379 /** Negative dof values cause the -dof-1 position in @a elemvect to receive
380 the -val in from this Vector. */
381 void GetSubVector(const Array<int> &dofs, Vector &elemvect) const;
382
383 /// Extract entries listed in @a dofs to the output array @a elem_data.
384 /** Negative dof values cause the -dof-1 position in @a elem_data to receive
385 the -val in from this Vector. */
386 void GetSubVector(const Array<int> &dofs, real_t *elem_data) const;
387
388 /// Set the entries listed in @a dofs to the given @a value.
389 /** Negative dof values cause the -dof-1 position in this Vector to receive
390 the -value. */
391 void SetSubVector(const Array<int> &dofs, const real_t value);
392
393 /** @brief Set the entries listed in @a dofs to the values given in the @a
394 elemvect Vector. Negative dof values cause the -dof-1 position in this
395 Vector to receive the -val from @a elemvect. */
396 void SetSubVector(const Array<int> &dofs, const Vector &elemvect);
397
398 /** @brief Set the entries listed in @a dofs to the values given the @a ,
399 elem_data array. Negative dof values cause the -dof-1 position in this
400 Vector to receive the -val from @a elem_data. */
401 void SetSubVector(const Array<int> &dofs, real_t *elem_data);
402
403 /** @brief Add elements of the @a elemvect Vector to the entries listed in @a
404 dofs. Negative dof values cause the -dof-1 position in this Vector to add
405 the -val from @a elemvect. */
406 void AddElementVector(const Array<int> & dofs, const Vector & elemvect);
407
408 /** @brief Add elements of the @a elem_data array to the entries listed in @a
409 dofs. Negative dof values cause the -dof-1 position in this Vector to add
410 the -val from @a elem_data. */
411 void AddElementVector(const Array<int> & dofs, real_t *elem_data);
412
413 /** @brief Add @a times the elements of the @a elemvect Vector to the entries
414 listed in @a dofs. Negative dof values cause the -dof-1 position in this
415 Vector to add the -a*val from @a elemvect. */
416 void AddElementVector(const Array<int> & dofs, const real_t a,
417 const Vector & elemvect);
418
419 /// Set all vector entries NOT in the @a dofs Array to the given @a val.
420 void SetSubVectorComplement(const Array<int> &dofs, const real_t val);
421
422 /// Prints vector to stream out.
423 void Print(std::ostream &out = mfem::out, int width = 8) const;
424
425#ifdef MFEM_USE_ADIOS2
426 /// Prints vector to stream out.
427 /// @param out adios2stream output
428 /// @param variable_name variable name associated with current Vector
429 void Print(adios2stream & out, const std::string& variable_name) const;
430#endif
431
432 /// Prints vector to stream out in HYPRE_Vector format.
433 void Print_HYPRE(std::ostream &out) const;
434
435 /// Print the Vector size and hash of its data.
436 /** This is a compact text representation of the Vector contents that can be
437 used to compare vectors from different runs without the need to save the
438 whole vector. */
439 void PrintHash(std::ostream &out) const;
440
441 /// Set random values in the vector.
442 void Randomize(int seed = 0);
443 /// Returns the l2 norm of the vector.
444 real_t Norml2() const;
445 /// Returns the l_infinity norm of the vector.
446 real_t Normlinf() const;
447 /// Returns the l_1 norm of the vector.
448 real_t Norml1() const;
449 /// Returns the l_p norm of the vector.
450 real_t Normlp(real_t p) const;
451 /// Returns the maximal element of the vector.
452 real_t Max() const;
453 /// Returns the minimal element of the vector.
454 real_t Min() const;
455 /// Return the sum of the vector entries
456 real_t Sum() const;
457 /// Compute the square of the Euclidean distance to another vector.
458 inline real_t DistanceSquaredTo(const real_t *p) const;
459 /// Compute the square of the Euclidean distance to another vector.
460 inline real_t DistanceSquaredTo(const Vector &p) const;
461 /// Compute the Euclidean distance to another vector.
462 inline real_t DistanceTo(const real_t *p) const;
463 /// Compute the Euclidean distance to another vector.
464 inline real_t DistanceTo(const Vector &p) const;
465
466 /** @brief Count the number of entries in the Vector for which isfinite
467 is false, i.e. the entry is a NaN or +/-Inf. */
468 int CheckFinite() const { return mfem::CheckFinite(HostRead(), size); }
469
470 /// Destroys vector.
471 virtual ~Vector();
472
473 /// Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), on_dev).
474 virtual const real_t *Read(bool on_dev = true) const
475 { return mfem::Read(data, size, on_dev); }
476
477 /// Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), false).
478 virtual const real_t *HostRead() const
479 { return mfem::Read(data, size, false); }
480
481 /// Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), on_dev).
482 virtual real_t *Write(bool on_dev = true)
483 { return mfem::Write(data, size, on_dev); }
484
485 /// Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), false).
486 virtual real_t *HostWrite()
487 { return mfem::Write(data, size, false); }
488
489 /// Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), on_dev).
490 virtual real_t *ReadWrite(bool on_dev = true)
491 { return mfem::ReadWrite(data, size, on_dev); }
492
493 /// Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), false).
495 { return mfem::ReadWrite(data, size, false); }
496
497};
498
499// Inline methods
500
501template <typename T>
502inline T ZeroSubnormal(T val)
503{
504 return (std::fpclassify(val) == FP_SUBNORMAL) ? 0.0 : val;
505}
506
507inline bool IsFinite(const real_t &val)
508{
509 // isfinite didn't appear in a standard until C99, and later C++11. It wasn't
510 // standard in C89 or C++98. PGI as of 14.7 still defines it as a macro.
511#ifdef isfinite
512 return isfinite(val);
513#else
514 return std::isfinite(val);
515#endif
516}
517
518inline int CheckFinite(const real_t *v, const int n)
519{
520 int bad = 0;
521 for (int i = 0; i < n; i++)
522 {
523 if (!IsFinite(v[i])) { bad++; }
524 }
525 return bad;
526}
527
528inline Vector::Vector(int s)
529{
530 MFEM_ASSERT(s>=0,"Unexpected negative size.");
531 size = s;
532 if (s > 0)
533 {
534 data.New(s);
535 }
536}
537
538inline void Vector::SetSize(int s)
539{
540 if (s == size)
541 {
542 return;
543 }
544 if (s <= data.Capacity())
545 {
546 size = s;
547 return;
548 }
549 // preserve a valid MemoryType and device flag
550 const MemoryType mt = data.GetMemoryType();
551 const bool use_dev = data.UseDevice();
552 data.Delete();
553 size = s;
554 data.New(s, mt);
555 data.UseDevice(use_dev);
556}
557
558inline void Vector::SetSize(int s, MemoryType mt)
559{
560 if (mt == data.GetMemoryType())
561 {
562 if (s == size)
563 {
564 return;
565 }
566 if (s <= data.Capacity())
567 {
568 size = s;
569 return;
570 }
571 }
572 const bool use_dev = data.UseDevice();
573 data.Delete();
574 if (s > 0)
575 {
576 data.New(s, mt);
577 size = s;
578 }
579 else
580 {
581 data.Reset();
582 size = 0;
583 }
584 data.UseDevice(use_dev);
585}
586
587inline void Vector::NewMemoryAndSize(const Memory<real_t> &mem, int s,
588 bool own_mem)
589{
590 data.Delete();
591 size = s;
592 if (own_mem)
593 {
594 data = mem;
595 }
596 else
597 {
598 data.MakeAlias(mem, 0, s);
599 }
600}
601
602inline void Vector::MakeRef(Vector &base, int offset, int s)
603{
604 data.Delete();
605 size = s;
606 data.MakeAlias(base.GetMemory(), offset, s);
607}
608
609inline void Vector::MakeRef(Vector &base, int offset)
610{
611 data.Delete();
612 data.MakeAlias(base.GetMemory(), offset, size);
613}
614
615inline void Vector::Destroy()
616{
617 const bool use_dev = data.UseDevice();
618 data.Delete();
619 size = 0;
620 data.Reset();
621 data.UseDevice(use_dev);
622}
623
625{
626 MFEM_ASSERT(data && i >= 0 && i < size,
627 "index [" << i << "] is out of range [0," << size << ")");
628
629 return data[i];
630}
631
632inline const real_t &Vector::operator()(int i) const
633{
634 MFEM_ASSERT(data && i >= 0 && i < size,
635 "index [" << i << "] is out of range [0," << size << ")");
636
637 return data[i];
638}
639
640inline void Vector::Swap(Vector &other)
641{
642 mfem::Swap(data, other.data);
643 mfem::Swap(size, other.size);
644}
645
646/// Specialization of the template function Swap<> for class Vector
647template<> inline void Swap<Vector>(Vector &a, Vector &b)
648{
649 a.Swap(b);
650}
651
653{
654 data.Delete();
655}
656
657inline real_t DistanceSquared(const real_t *x, const real_t *y, const int n)
658{
659 real_t d = 0.0;
660
661 for (int i = 0; i < n; i++)
662 {
663 d += (x[i]-y[i])*(x[i]-y[i]);
664 }
665
666 return d;
667}
668
669inline real_t Distance(const real_t *x, const real_t *y, const int n)
670{
671 return std::sqrt(DistanceSquared(x, y, n));
672}
673
674inline real_t Distance(const Vector &x, const Vector &y)
675{
676 return x.DistanceTo(y);
677}
678
680{
681 return DistanceSquared(data, p, size);
682}
683
685{
686 MFEM_ASSERT(p.Size() == Size(), "Incompatible vector sizes.");
687 return DistanceSquared(data, p.data, size);
688}
689
690inline real_t Vector::DistanceTo(const real_t *p) const
691{
692 return Distance(data, p, size);
693}
694
695inline real_t Vector::DistanceTo(const Vector &p) const
696{
697 MFEM_ASSERT(p.Size() == Size(), "Incompatible vector sizes.");
698 return Distance(data, p.data, size);
699}
700
701/// Returns the inner product of x and y
702/** In parallel this computes the inner product of the local vectors,
703 producing different results on each MPI rank.
704*/
705inline real_t InnerProduct(const Vector &x, const Vector &y)
706{
707 return x * y;
708}
709
710#ifdef MFEM_USE_MPI
711/// Returns the inner product of x and y in parallel
712/** In parallel this computes the inner product of the global vectors,
713 producing identical results on each MPI rank.
714*/
715inline real_t InnerProduct(MPI_Comm comm, const Vector &x, const Vector &y)
716{
717 real_t loc_prod = x * y;
718 real_t glb_prod;
719 MPI_Allreduce(&loc_prod, &glb_prod, 1, MFEM_MPI_REAL_T, MPI_SUM, comm);
720 return glb_prod;
721}
722#endif
723
724} // namespace mfem
725
726#endif
Class used by MFEM to store pointers to host and/or device memory.
void SetHostPtrOwner(bool own) const
Set/clear the ownership flag for the host pointer. Ownership indicates whether the pointer will be de...
int Capacity() const
Return the size of the allocated memory.
void MakeAlias(const Memory &base, int offset, int size)
Create a memory object that points inside the memory object base.
bool UseDevice() const
Read the internal device flag.
void SyncAlias(const Memory &base, int alias_size) const
Update the alias Memory *this to match the memory location (all valid locations) of its base Memory,...
bool OwnsHostPtr() const
Return true if the host pointer is owned. Ownership indicates whether the pointer will be deleted by ...
MemoryType GetMemoryType() const
Return a MemoryType that is currently valid. If both the host and the device pointers are currently v...
void Sync(const Memory &other) const
Copy the host/device pointer validity flags from other to *this.
void DeleteDevice(bool copy_to_host=true)
Delete the device pointer, if owned. If copy_to_host is true and the data is valid only on device,...
void Reset()
Reset the memory to be empty, ensuring that Delete() will be a no-op.
void Wrap(T *ptr, int size, bool own)
Wrap an externally allocated host pointer, ptr with the current host memory type returned by MemoryMa...
void Delete()
Delete the owned pointers and reset the Memory object.
void New(int size)
Allocate host memory for size entries with the current host memory type returned by MemoryManager::Ge...
Vector data type.
Definition vector.hpp:80
void Randomize(int seed=0)
Set random values in the vector.
Definition vector.cpp:816
void PrintHash(std::ostream &out) const
Print the Vector size and hash of its data.
Definition vector.cpp:808
real_t & Elem(int i)
Access Vector entries. Index i = 0 .. size-1.
Definition vector.cpp:95
Vector(int size_, MemoryType mt)
Create a Vector of size size_ using MemoryType mt.
Definition vector.hpp:114
void StealData(real_t **p)
Changes the ownership of the data; after the call the Vector is empty.
Definition vector.hpp:266
void SetSize(int s, const Vector &v)
Resize the vector to size s using the MemoryType of v.
Definition vector.hpp:164
void SetVector(const Vector &v, int offset)
Definition vector.cpp:274
virtual const real_t * HostRead() const
Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), false).
Definition vector.hpp:478
virtual const real_t * Read(bool on_dev=true) const
Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), on_dev).
Definition vector.hpp:474
const real_t & operator[](int i) const
Read only access to Vector entries using [] for 0-based indexing.
Definition vector.hpp:292
void median(const Vector &lo, const Vector &hi)
v = median(v,lo,hi) entrywise. Implementation assumes lo <= hi.
Definition vector.cpp:554
void Load(std::istream &in)
Load a vector from an input stream, reading the size from the stream.
Definition vector.hpp:148
void Neg()
(*this) = -(*this)
Definition vector.cpp:300
real_t * begin()
STL-like begin.
Definition vector.hpp:237
void Print(std::ostream &out=mfem::out, int width=8) const
Prints vector to stream out.
Definition vector.cpp:755
virtual real_t * ReadWrite(bool on_dev=true)
Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), on_dev).
Definition vector.hpp:490
void SetDataAndSize(real_t *d, int s)
Set the Vector data and size.
Definition vector.hpp:175
Vector(Vector &base, int base_offset, int size_)
Create a Vector referencing a sub-vector of the Vector base starting at the given offset,...
Definition vector.hpp:110
real_t Normlinf() const
Returns the l_infinity norm of the vector.
Definition vector.cpp:850
void MakeDataOwner() const
Set the Vector data (host pointer) ownership flag.
Definition vector.hpp:206
real_t Norml1() const
Returns the l_1 norm of the vector.
Definition vector.cpp:861
bool OwnsData() const
Read the Vector data (host pointer) ownership flag.
Definition vector.hpp:263
void SetSubVector(const Array< int > &dofs, const real_t value)
Set the entries listed in dofs to the given value.
Definition vector.cpp:604
friend void subtract(const Vector &v1, const Vector &v2, Vector &v)
Set v = v1 - v2.
Definition vector.cpp:472
Memory< real_t > & GetMemory()
Return a reference to the Memory object used by the Vector.
Definition vector.hpp:249
void AddElementVector(const Array< int > &dofs, const Vector &elemvect)
Add elements of the elemvect Vector to the entries listed in dofs. Negative dof values cause the -dof...
Definition vector.cpp:670
real_t operator*(const real_t *) const
Dot product with a double * array.
Definition vector.cpp:105
Memory< real_t > data
Definition vector.hpp:83
real_t & operator[](int i)
Access Vector entries using [] for 0-based indexing.
Definition vector.hpp:288
int CheckFinite() const
Count the number of entries in the Vector for which isfinite is false, i.e. the entry is a NaN or +/-...
Definition vector.hpp:468
void AddSubVector(const Vector &v, int offset)
Definition vector.cpp:287
const real_t * end() const
STL-like end (const version).
Definition vector.hpp:246
void Swap(Vector &other)
Swap the contents of two Vectors.
Definition vector.hpp:640
Vector & operator*=(real_t c)
Definition vector.cpp:162
void SyncAliasMemory(const Vector &v) const
Update the alias memory location of the vector to match v.
Definition vector.hpp:259
real_t Norml2() const
Returns the l2 norm of the vector.
Definition vector.cpp:832
void SyncMemory(const Vector &v) const
Update the memory location of the vector to match v.
Definition vector.hpp:256
void Load(std::istream **in, int np, int *dim)
Reads a vector from multiple files.
Definition vector.cpp:49
real_t DistanceSquaredTo(const real_t *p) const
Compute the square of the Euclidean distance to another vector.
Definition vector.hpp:679
Vector & Set(const real_t a, const Vector &x)
(*this) = a * x
Definition vector.cpp:262
real_t * end()
STL-like end.
Definition vector.hpp:240
const Memory< real_t > & GetMemory() const
Return a reference to the Memory object used by the Vector, const version.
Definition vector.hpp:253
void NewMemoryAndSize(const Memory< real_t > &mem, int s, bool own_mem)
Reset the Vector to use the given external Memory mem and size s.
Definition vector.hpp:587
real_t Max() const
Returns the maximal element of the vector.
Definition vector.cpp:922
virtual bool UseDevice() const
Return the device flag of the Memory object used by the Vector.
Definition vector.hpp:139
Vector(int size_, MemoryType h_mt, MemoryType d_mt)
Create a Vector of size size_ using host MemoryType h_mt and device MemoryType d_mt.
Definition vector.hpp:119
void Destroy()
Destroy a vector.
Definition vector.hpp:615
Vector(const T(&values)[N])
Create a vector using a braced initializer list.
Definition vector.hpp:124
int Size() const
Returns the size of the vector.
Definition vector.hpp:218
virtual void UseDevice(bool use_dev) const
Enable execution of Vector operations using the mfem::Device.
Definition vector.hpp:136
real_t Sum() const
Return the sum of the vector entries.
Definition vector.cpp:1283
void Print_HYPRE(std::ostream &out) const
Prints vector to stream out in HYPRE_Vector format.
Definition vector.cpp:789
void SetSize(int s)
Resize the vector to size s.
Definition vector.hpp:538
void Reciprocal()
(*this)(i) = 1.0 / (*this)(i)
Definition vector.cpp:308
virtual real_t * HostWrite()
Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), false).
Definition vector.hpp:486
void NewDataAndSize(real_t *d, int s)
Set the Vector data and size, deleting the old data, if owned.
Definition vector.hpp:181
real_t Normlp(real_t p) const
Returns the l_p norm of the vector.
Definition vector.cpp:872
Vector & operator-=(real_t c)
Definition vector.cpp:205
const real_t * begin() const
STL-like begin (const version).
Definition vector.hpp:243
int Capacity() const
Return the size of the currently allocated data array.
Definition vector.hpp:222
real_t * GetData() const
Return a pointer to the beginning of the Vector data.
Definition vector.hpp:227
Vector & operator=(const real_t *v)
Copy Size() entries from v.
Definition vector.cpp:118
void SetSubVectorComplement(const Array< int > &dofs, const real_t val)
Set all vector entries NOT in the dofs Array to the given val.
Definition vector.cpp:739
void SetData(real_t *d)
Definition vector.hpp:168
virtual real_t * HostReadWrite()
Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), false).
Definition vector.hpp:494
virtual ~Vector()
Destroys vector.
Definition vector.hpp:652
real_t Min() const
Returns the minimal element of the vector.
Definition vector.cpp:1212
void GetSubVector(const Array< int > &dofs, Vector &elemvect) const
Extract entries listed in dofs to the output Vector elemvect.
Definition vector.cpp:578
Vector & operator+=(real_t c)
Definition vector.cpp:226
virtual real_t * Write(bool on_dev=true)
Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), on_dev).
Definition vector.hpp:482
real_t * StealData()
Changes the ownership of the data; after the call the Vector is empty.
Definition vector.hpp:270
Vector & Add(const real_t a, const Vector &Va)
(*this) += a * Va
Definition vector.cpp:247
real_t & operator()(int i)
Access Vector entries using () for 0-based indexing.
Definition vector.hpp:624
void DeleteDevice(bool copy_to_host=true)
Delete the device pointer, if owned. If copy_to_host is true and the data is valid only on device,...
Definition vector.hpp:214
Vector(real_t *data_, int size_)
Creates a vector referencing an array of doubles, owned by someone else.
Definition vector.hpp:105
void cross3D(const Vector &vin, Vector &vout) const
Definition vector.cpp:541
void MakeRef(Vector &base, int offset, int size)
Reset the Vector to be a reference to a sub-vector of base.
Definition vector.hpp:602
Vector & operator/=(real_t c)
Definition vector.cpp:183
friend void add(const Vector &v1, const Vector &v2, Vector &v)
Set v = v1 + v2.
Definition vector.cpp:316
real_t DistanceTo(const real_t *p) const
Compute the Euclidean distance to another vector.
Definition vector.hpp:690
const real_t alpha
Definition ex15.cpp:369
int dim
Definition ex24.cpp:53
real_t b
Definition lissajous.cpp:42
real_t a
Definition lissajous.cpp:41
const T * Read(const Memory< T > &mem, int size, bool on_dev=true)
Get a pointer for read access to mem with the mfem::Device's DeviceMemoryClass, if on_dev = true,...
Definition device.hpp:320
void Swap(Array< T > &, Array< T > &)
Definition array.hpp:648
real_t infinity()
Define a shortcut for std::numeric_limits<double>::infinity()
Definition vector.hpp:45
T * Write(Memory< T > &mem, int size, bool on_dev=true)
Get a pointer for write access to mem with the mfem::Device's DeviceMemoryClass, if on_dev = true,...
Definition device.hpp:337
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition globals.hpp:66
real_t rand_real()
Generate a random real_t number in the interval [0,1) using rand().
Definition vector.hpp:59
bool IsFinite(const real_t &val)
Definition vector.hpp:507
real_t Distance(const real_t *x, const real_t *y, const int n)
Definition vector.hpp:669
real_t InnerProduct(HypreParVector *x, HypreParVector *y)
Definition hypre.cpp:439
real_t DistanceSquared(const real_t *x, const real_t *y, const int n)
Definition vector.hpp:657
T * ReadWrite(Memory< T > &mem, int size, bool on_dev=true)
Get a pointer for read+write access to mem with the mfem::Device's DeviceMemoryClass,...
Definition device.hpp:354
int CheckFinite(const real_t *v, const int n)
Definition vector.hpp:518
T ZeroSubnormal(T val)
Definition vector.hpp:502
float real_t
Definition config.hpp:43
MemoryType
Memory types supported by MFEM.
void Swap< Vector >(Vector &a, Vector &b)
Specialization of the template function Swap<> for class Vector.
Definition vector.hpp:647
real_t p(const Vector &x, real_t t)
RefCoord s[3]