dense_vector.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_DENSE_VECTOR_H
21 #define LIBMESH_DENSE_VECTOR_H
22 
23 // Local Includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/compare_types.h"
27 #include "libmesh/tensor_tools.h"
28 
29 // C++ includes
30 #include <vector>
31 
32 namespace libMesh
33 {
34 
35 // Forward Declarations
36 
37 
38 
48 // ------------------------------------------------------------
49 // DenseVector class definition
50 template<typename T>
51 class DenseVector : public DenseVectorBase<T>
52 {
53 public:
54 
58  explicit
59  DenseVector(const unsigned int n=0);
60 
64  template <typename T2>
65  DenseVector (const DenseVector<T2>& other_vector);
66 
70  template <typename T2>
71  DenseVector (const std::vector<T2>& other_vector);
72 
77 
81  virtual unsigned int size() const {
82  return libmesh_cast_int<unsigned int>(_val.size());
83  }
84 
88  virtual void zero();
89 
93  T operator() (const unsigned int i) const;
94 
98  T & operator() (const unsigned int i);
99 
103  virtual T el(const unsigned int i) const { return (*this)(i); }
104 
108  virtual T & el(const unsigned int i) { return (*this)(i); }
109 
113  template <typename T2>
114  DenseVector<T>& operator = (const DenseVector<T2>& other_vector);
115 
119  void swap(DenseVector<T>& other_vector);
120 
124  void resize (const unsigned int n);
125 
129  void scale (const T factor);
130 
134  DenseVector<T>& operator*= (const T factor);
135 
141  template <typename T2, typename T3>
142  typename boostcopy::enable_if_c<
143  ScalarTraits<T2>::value, void >::type
144  add (const T2 factor,
145  const DenseVector<T3>& vec);
146 
151  template <typename T2>
152  typename CompareTypes<T, T2>::supertype dot (const DenseVector<T2> &vec) const;
153 
158  template <typename T2>
160 
164  template <typename T2>
165  bool operator== (const DenseVector<T2> &vec) const;
166 
170  template <typename T2>
171  bool operator!= (const DenseVector<T2> &vec) const;
172 
176  template <typename T2>
178 
182  template <typename T2>
184 
190  Real min () const;
191 
197  Real max () const;
198 
203  Real l1_norm () const;
204 
210  Real l2_norm () const;
211 
217  Real linfty_norm () const;
218 
223  void get_principal_subvector (unsigned int sub_n, DenseVector<T>& dest) const;
224 
230  std::vector<T>& get_values() { return _val; }
231 
237  const std::vector<T>& get_values() const { return _val; }
238 
239 private:
240 
244  std::vector<T> _val;
245 
246 };
247 
248 
249 
250 // ------------------------------------------------------------
251 // DenseVector member functions
252 template<typename T>
253 inline
254 DenseVector<T>::DenseVector(const unsigned int n) :
255  _val (n, T(0.))
256 {
257 }
258 
259 
260 
261 template<typename T>
262 template<typename T2>
263 inline
265  DenseVectorBase<T>()
266 {
267  const std::vector<T2> &other_vals = other_vector.get_values();
268 
269  _val.clear();
270  _val.reserve(other_vals.size());
271 
272  for (unsigned int i=0; i<other_vals.size(); i++)
273  _val.push_back(other_vals[i]);
274 }
275 
276 
277 
278 template<typename T>
279 template<typename T2>
280 inline
281 DenseVector<T>::DenseVector (const std::vector<T2>& other_vector) :
282  _val(other_vector)
283 {
284 }
285 
286 
287 
288 
289 
290 template<typename T>
291 template<typename T2>
292 inline
294 {
295  // _val = other_vector._val;
296 
297  const std::vector<T2> &other_vals = other_vector.get_values();
298 
299  _val.clear();
300  _val.reserve(other_vals.size());
301 
302  for (unsigned int i=0; i<other_vals.size(); i++)
303  _val.push_back(other_vals[i]);
304 
305  return *this;
306 }
307 
308 
309 
310 template<typename T>
311 inline
313 {
314  _val.swap(other_vector._val);
315 }
316 
317 
318 
319 template<typename T>
320 inline
321 void DenseVector<T>::resize(const unsigned int n)
322 {
323  _val.resize(n);
324 
325  zero();
326 }
327 
328 
329 
330 template<typename T>
331 inline
333 {
334  std::fill (_val.begin(),
335  _val.end(),
336  T(0.));
337 }
338 
339 
340 
341 template<typename T>
342 inline
343 T DenseVector<T>::operator () (const unsigned int i) const
344 {
345  libmesh_assert_less (i, _val.size());
346 
347  return _val[i];
348 }
349 
350 
351 
352 template<typename T>
353 inline
354 T & DenseVector<T>::operator () (const unsigned int i)
355 {
356  libmesh_assert_less (i, _val.size());
357 
358  return _val[i];
359 }
360 
361 
362 
363 template<typename T>
364 inline
365 void DenseVector<T>::scale (const T factor)
366 {
367  for (unsigned int i=0; i<_val.size(); i++)
368  _val[i] *= factor;
369 }
370 
371 
372 
373 template<typename T>
374 inline
376 {
377  this->scale(factor);
378  return *this;
379 }
380 
381 
382 
383 template<typename T>
384 template<typename T2, typename T3>
385 inline
386 typename boostcopy::enable_if_c<
387  ScalarTraits<T2>::value, void >::type
388 DenseVector<T>::add (const T2 factor,
389  const DenseVector<T3>& vec)
390 {
391  libmesh_assert_equal_to (this->size(), vec.size());
392 
393  for (unsigned int i=0; i<this->size(); i++)
394  (*this)(i) += factor*vec(i);
395 }
396 
397 template<typename T>
398 template<typename T2>
399 inline
401 {
402  libmesh_assert_equal_to (this->size(), vec.size());
403 
404  typename CompareTypes<T, T2>::supertype val = 0.;
405 
406  for (unsigned int i=0; i<this->size(); i++)
407  val += (*this)(i)*libmesh_conj(vec(i));
408 
409  return val;
410 }
411 
412 template<typename T>
413 template<typename T2>
414 inline
416 {
417  libmesh_assert_equal_to (this->size(), vec.size());
418 
419  typename CompareTypes<T, T2>::supertype val = 0.;
420 
421  for (unsigned int i=0; i<this->size(); i++)
422  val += (*this)(i)*(vec(i));
423 
424  return val;
425 }
426 
427 template<typename T>
428 template<typename T2>
429 inline
431 {
432  libmesh_assert_equal_to (this->size(), vec.size());
433 
434  for (unsigned int i=0; i<this->size(); i++)
435  if ((*this)(i) != vec(i))
436  return false;
437 
438  return true;
439 }
440 
441 
442 
443 template<typename T>
444 template<typename T2>
445 inline
447 {
448  libmesh_assert_equal_to (this->size(), vec.size());
449 
450  for (unsigned int i=0; i<this->size(); i++)
451  if ((*this)(i) != vec(i))
452  return true;
453 
454  return false;
455 }
456 
457 
458 
459 template<typename T>
460 template<typename T2>
461 inline
463 {
464  libmesh_assert_equal_to (this->size(), vec.size());
465 
466  for (unsigned int i=0; i<this->size(); i++)
467  (*this)(i) += vec(i);
468 
469  return *this;
470 }
471 
472 
473 
474 template<typename T>
475 template<typename T2>
476 inline
478 {
479  libmesh_assert_equal_to (this->size(), vec.size());
480 
481  for (unsigned int i=0; i<this->size(); i++)
482  (*this)(i) -= vec(i);
483 
484  return *this;
485 }
486 
487 
488 
489 template<typename T>
490 inline
492 {
493  libmesh_assert (this->size());
494  Real my_min = libmesh_real((*this)(0));
495 
496  for (unsigned int i=1; i!=this->size(); i++)
497  {
498  Real current = libmesh_real((*this)(i));
499  my_min = (my_min < current? my_min : current);
500  }
501  return my_min;
502 }
503 
504 
505 
506 template<typename T>
507 inline
509 {
510  libmesh_assert (this->size());
511  Real my_max = libmesh_real((*this)(0));
512 
513  for (unsigned int i=1; i!=this->size(); i++)
514  {
515  Real current = libmesh_real((*this)(i));
516  my_max = (my_max > current? my_max : current);
517  }
518  return my_max;
519 }
520 
521 
522 
523 template<typename T>
524 inline
526 {
527  Real my_norm = 0.;
528  for (unsigned int i=0; i!=this->size(); i++)
529  {
530  my_norm += std::abs((*this)(i));
531  }
532  return my_norm;
533 }
534 
535 
536 
537 template<typename T>
538 inline
540 {
541  Real my_norm = 0.;
542  for (unsigned int i=0; i!=this->size(); i++)
543  {
544  my_norm += TensorTools::norm_sq((*this)(i));
545  }
546  return sqrt(my_norm);
547 }
548 
549 
550 
551 template<typename T>
552 inline
554 {
555  if (!this->size())
556  return 0.;
557  Real my_norm = TensorTools::norm_sq((*this)(0));
558 
559  for (unsigned int i=1; i!=this->size(); i++)
560  {
561  Real current = TensorTools::norm_sq((*this)(i));
562  my_norm = (my_norm > current? my_norm : current);
563  }
564  return sqrt(my_norm);
565 }
566 
567 template<typename T>
568 inline
570  DenseVector<T>& dest) const
571 {
572  libmesh_assert_less_equal ( sub_n, this->size() );
573 
574  dest.resize(sub_n);
575  for(unsigned int i=0; i<sub_n; i++)
576  dest(i) = (*this)(i);
577 }
578 
579 } // namespace libMesh
580 
581 #endif // LIBMESH_DENSE_VECTOR_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:04 UTC

Hosted By:
SourceForge.net Logo