laspack_vector.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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 // C++ includes
21 #include <algorithm> // for std::min
22 #include <limits>
23 
24 // Local Includes
26 #include "libmesh/dense_vector.h"
27 #include "libmesh/laspack_vector.h"
28 #include "libmesh/laspack_matrix.h"
29 
30 
31 #ifdef LIBMESH_HAVE_LASPACK
32 
33 namespace libMesh
34 {
35 
36 template <typename T>
38 {
39  libmesh_assert (this->closed());
40 
41  T _sum = 0;
42 
43  const numeric_index_type n = this->size();
44 
45  for (numeric_index_type i=0; i!=n; ++i)
46  _sum += (*this)(i);
47 
48  return _sum;
49 }
50 
51 
52 
53 template <typename T>
55 {
56  libmesh_assert (this->closed());
57 
58  return static_cast<Real>(l1Norm_V(const_cast<QVector*>(&_vec)));
59 }
60 
61 
62 
63 template <typename T>
65 {
66  libmesh_assert (this->closed());
67 
68  return static_cast<Real>(l2Norm_V(const_cast<QVector*>(&_vec)));
69 }
70 
71 
72 
73 template <typename T>
75 {
76  libmesh_assert (this->closed());
77 
78  return static_cast<Real>(MaxNorm_V(const_cast<QVector*>(&_vec)));
79 }
80 
81 
82 
83 template <typename T>
85 {
86  libmesh_assert (this->closed());
87 
88  this->add(1., v);
89 
90  return *this;
91 }
92 
93 
94 
95 
96 template <typename T>
98 {
99  libmesh_assert (this->closed());
100 
101  this->add(-1., v);
102 
103  return *this;
104 }
105 
106 
107 
108 template <typename T>
110 {
111  libmesh_assert_equal_to(size(), v.size());
112 
113  const numeric_index_type n = this->size();
114 
115  for (numeric_index_type i=0; i<n; i++)
116  this->set(i, (*this)(i) / v(i));
117 
118  return *this;
119 }
120 
121 
122 
123 template <typename T>
125 {
126  const numeric_index_type n = this->size();
127 
128  for (numeric_index_type i=0; i<n; i++)
129  {
130  T v = (*this)(i);
131 
132  // Don't divide by zero!
133  libmesh_assert_not_equal_to (v, T(0));
134 
135  this->set(i, 1. / v);
136  }
137 }
138 
139 
140 
141 template <typename T>
143 {
144  const numeric_index_type n = this->size();
145 
146  for (numeric_index_type i=0; i<n; i++)
147  {
148  T v = (*this)(i);
149 
150  this->set(i, libmesh_conj(v) );
151  }
152 }
153 
154 
155 template <typename T>
156 void LaspackVector<T>::add (const T v)
157 {
158  const numeric_index_type n = this->size();
159 
160  for (numeric_index_type i=0; i<n; i++)
161  this->add (i, v);
162 
163 #ifndef NDEBUG
164  this->_is_closed = false;
165 #endif
166 }
167 
168 
169 
170 
171 template <typename T>
173 {
174  this->add (1., v);
175 }
176 
177 
178 
179 template <typename T>
180 void LaspackVector<T>::add (const T a, const NumericVector<T> & v_in)
181 {
182  // Make sure the vector passed in is really a LaspackVector
183  const LaspackVector * v = cast_ptr<const LaspackVector *>(&v_in);
184 
185 #ifndef NDEBUG
186  const bool was_closed = this->_is_closed;
187 #endif
188 
189  libmesh_assert(v);
190  libmesh_assert_equal_to (this->size(), v->size());
191 
192  for (numeric_index_type i=0; i<v->size(); i++)
193  this->add (i, a*(*v)(i));
194 
195 #ifndef NDEBUG
196  this->_is_closed = was_closed;
197 #endif
198 }
199 
200 
201 
202 template <typename T>
204  const SparseMatrix<T> & mat_in)
205 {
206  // Make sure the data passed in are really in Laspack types
207  const LaspackVector<T> * vec = cast_ptr<const LaspackVector<T> *>(&vec_in);
208  const LaspackMatrix<T> * mat = cast_ptr<const LaspackMatrix<T> *>(&mat_in);
209 
210  libmesh_assert(vec);
211  libmesh_assert(mat);
212 
213  // += mat*vec
214  AddAsgn_VV (&_vec, Mul_QV(const_cast<QMatrix*>(&mat->_QMat),
215  const_cast<QVector*>(&vec->_vec)));
216 }
217 
218 
219 template <typename T>
221  const SparseMatrix<T> &)
222 {
223  libmesh_not_implemented();
224 }
225 
226 
227 
228 template <typename T>
229 void LaspackVector<T>::scale (const T factor)
230 {
231  libmesh_assert (this->initialized());
232 
233  Asgn_VV(&_vec, Mul_SV (factor, &_vec));
234 }
235 
236 template <typename T>
238 {
239  libmesh_assert (this->initialized());
240 
241  const numeric_index_type n = this->size();
242 
243  for (numeric_index_type i=0; i!=n; ++i)
244  this->set(i,std::abs((*this)(i)));
245 }
246 
247 template <typename T>
249 {
250  libmesh_assert (this->initialized());
251 
252  // Make sure the NumericVector passed in is really a LaspackVector
253  const LaspackVector<T> * v = cast_ptr<const LaspackVector<T> *>(&V);
254  libmesh_assert(v);
255 
256  return Mul_VV (const_cast<QVector*>(&(this->_vec)),
257  const_cast<QVector*>(&(v->_vec)));
258 }
259 
260 
261 
262 template <typename T>
265 {
266  libmesh_assert (this->initialized());
267  libmesh_assert (this->closed());
268 
269  V_SetAllCmp (&_vec, s);
270 
271  return *this;
272 }
273 
274 
275 
276 template <typename T>
279 {
280  // Make sure the NumericVector passed in is really a LaspackVector
281  const LaspackVector<T> * v =
282  cast_ptr<const LaspackVector<T> *>(&v_in);
283 
284  libmesh_assert(v);
285 
286  *this = *v;
287 
288  return *this;
289 }
290 
291 
292 
293 template <typename T>
296 {
297  libmesh_assert (this->initialized());
298  libmesh_assert (v.closed());
299  libmesh_assert_equal_to (this->size(), v.size());
300 
301  if (v.size() != 0)
302  Asgn_VV (const_cast<QVector*>(&_vec),
303  const_cast<QVector*>(&v._vec)
304  );
305 
306 #ifndef NDEBUG
307  this->_is_closed = true;
308 #endif
309 
310  return *this;
311 }
312 
313 
314 
315 template <typename T>
317 LaspackVector<T>::operator = (const std::vector<T> & v)
318 {
323  if (this->size() == v.size())
324  for (numeric_index_type i=0; i<v.size(); i++)
325  this->set (i, v[i]);
326 
327  else
328  libmesh_error_msg("this->size() = " << this->size() << " must be equal to v.size() = " << v.size());
329 
330  return *this;
331 }
332 
333 
334 template <typename T>
336 {
337  // Make sure the NumericVector passed in is really a LaspackVector
338  LaspackVector<T> * v_local =
339  cast_ptr<LaspackVector<T> *>(&v_local_in);
340 
341  libmesh_assert(v_local);
342 
343  *v_local = *this;
344 }
345 
346 
347 
348 template <typename T>
350  const std::vector<numeric_index_type> & libmesh_dbg_var(send_list)) const
351 {
352  // Make sure the NumericVector passed in is really a LaspackVector
353  LaspackVector<T> * v_local =
354  cast_ptr<LaspackVector<T> *>(&v_local_in);
355 
356  libmesh_assert(v_local);
357  libmesh_assert_less_equal (send_list.size(), v_local->size());
358 
359  *v_local = *this;
360 }
361 
362 
363 
364 template <typename T>
365 void LaspackVector<T>::localize (std::vector<T> & v_local,
366  const std::vector<numeric_index_type> & indices) const
367 {
368  // LaspackVectors are serial, so we can just copy values
369  v_local.resize(indices.size());
370 
371  for (numeric_index_type i=0; i<v_local.size(); i++)
372  v_local[i] = (*this)(indices[i]);
373 }
374 
375 
376 
377 template <typename T>
379  const numeric_index_type libmesh_dbg_var(last_local_idx),
380  const std::vector<numeric_index_type> & libmesh_dbg_var(send_list))
381 {
382  libmesh_assert_equal_to (first_local_idx, 0);
383  libmesh_assert_equal_to (last_local_idx+1, this->size());
384 
385  libmesh_assert_less_equal (send_list.size(), this->size());
386 
387 #ifndef NDEBUG
388  this->_is_closed = true;
389 #endif
390 }
391 
392 
393 
394 template <typename T>
395 void LaspackVector<T>::localize (std::vector<T> & v_local) const
396 
397 {
398  v_local.resize(this->size());
399 
400  for (numeric_index_type i=0; i<v_local.size(); i++)
401  v_local[i] = (*this)(i);
402 }
403 
404 
405 
406 template <typename T>
407 void LaspackVector<T>::localize_to_one (std::vector<T> & v_local,
408  const processor_id_type libmesh_dbg_var(pid)) const
409 {
410  libmesh_assert_equal_to (pid, 0);
411 
412  this->localize (v_local);
413 }
414 
415 
416 
417 template <typename T>
419  const NumericVector<T> & /*vec2*/)
420 {
421  libmesh_not_implemented();
422 }
423 
424 
425 
426 template <typename T>
428 {
429  libmesh_assert (this->initialized());
430  if (!this->size())
432 
433  Real the_max = libmesh_real((*this)(0));
434 
435  const numeric_index_type n = this->size();
436 
437  for (numeric_index_type i=1; i<n; i++)
438  the_max = std::max (the_max, libmesh_real((*this)(i)));
439 
440  return the_max;
441 }
442 
443 
444 
445 template <typename T>
447 {
448  libmesh_assert (this->initialized());
449  if (!this->size())
451 
452  Real the_min = libmesh_real((*this)(0));
453 
454  const numeric_index_type n = this->size();
455 
456  for (numeric_index_type i=1; i<n; i++)
457  the_min = std::min (the_min, libmesh_real((*this)(i)));
458 
459  return the_min;
460 }
461 
462 
463 //------------------------------------------------------------------
464 // Explicit instantiations
465 template class LaspackVector<Number>;
466 
467 } // namespace libMesh
468 
469 
470 #endif // #ifdef LIBMESH_HAVE_LASPACK
T libmesh_real(T a)
virtual bool closed() const
bool closed()
Definition: libmesh.C:279
double abs(double a)
virtual Real l1_norm() const libmesh_override
virtual numeric_index_type size() const =0
T libmesh_conj(T a)
virtual NumericVector< T > & operator/=(NumericVector< T > &v) libmesh_override
virtual T sum() const libmesh_override
virtual Real l2_norm() const libmesh_override
virtual Real max() const libmesh_override
uint8_t processor_id_type
Definition: id_types.h:99
long double max(long double a, double b)
virtual void add_vector_transpose(const NumericVector< T > &, const SparseMatrix< T > &) libmesh_override
virtual void conjugate() libmesh_override
libmesh_assert(j)
virtual void add(const numeric_index_type i, const T value) libmesh_override
virtual Real linfty_norm() const libmesh_override
virtual void abs() libmesh_override
dof_id_type numeric_index_type
Definition: id_types.h:92
virtual NumericVector< T > & operator-=(const NumericVector< T > &v) libmesh_override
PetscErrorCode Vec Mat libmesh_dbg_var(j)
virtual Real min() const libmesh_override
virtual void pointwise_mult(const NumericVector< T > &vec1, const NumericVector< T > &vec2) libmesh_override
virtual T dot(const NumericVector< T > &V) const libmesh_override
virtual void add_vector(const NumericVector< T > &, const SparseMatrix< T > &) libmesh_override
virtual numeric_index_type size() const libmesh_override
virtual void scale(const T factor) libmesh_override
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual NumericVector< T > & operator+=(const NumericVector< T > &v) libmesh_override
virtual void reciprocal() libmesh_override
bool initialized()
Definition: libmesh.C:272
virtual void localize(std::vector< T > &v_local) const libmesh_override
virtual NumericVector< T > & operator=(const T s) libmesh_override
long double min(long double a, double b)
virtual void localize_to_one(std::vector< T > &v_local, const processor_id_type proc_id=0) const libmesh_override