Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
diagonal_matrix.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_diagonal_matrix_h
17 #define dealii_diagonal_matrix_h
18 
19 
20 #include <deal.II/lac/vector.h>
21 #include <deal.II/lac/vector_operation.h>
22 
23 DEAL_II_NAMESPACE_OPEN
24 
47 template <typename VectorType = Vector<double>>
49 {
50 public:
51  using value_type = typename VectorType::value_type;
52  using size_type = typename VectorType::size_type;
53 
57  DiagonalMatrix() = default;
58 
63  void
64  reinit(const VectorType &vec);
65 
72  void
74 
79  VectorType &
80  get_vector();
81 
85  void
86  clear();
87 
91  const VectorType &
92  get_vector() const;
93 
98  size_type
99  m() const;
100 
105  size_type
106  n() const;
107 
117  value_type
118  operator()(const size_type i, const size_type j) const;
119 
129  value_type &
130  operator()(const size_type i, const size_type j);
131 
143  template <typename number2>
144  void
145  add(const size_type row,
146  const size_type n_cols,
147  const size_type *col_indices,
148  const number2 * values,
149  const bool elide_zero_values = true,
150  const bool col_indices_are_sorted = false);
151 
158  void
159  add(const size_type i, const size_type j, const value_type value);
160 
164  void
165  vmult(VectorType &dst, const VectorType &src) const;
166 
172  void
173  Tvmult(VectorType &dst, const VectorType &src) const;
174 
180  void
181  vmult_add(VectorType &dst, const VectorType &src) const;
182 
188  void
189  Tvmult_add(VectorType &dst, const VectorType &src) const;
190 
198  void
199  initialize_dof_vector(VectorType &dst) const;
200 
204  std::size_t
205  memory_consumption() const;
206 
207 private:
211  VectorType diagonal;
212 };
213 
214 /* ---------------------------------- Inline functions ------------------- */
215 
216 #ifndef DOXYGEN
217 
218 template <typename VectorType>
219 void
221 {
222  diagonal.reinit(0);
223 }
224 
225 
226 
227 template <typename VectorType>
228 std::size_t
230 {
231  return diagonal.memory_consumption();
232 }
233 
234 
235 
236 template <typename VectorType>
237 void
238 DiagonalMatrix<VectorType>::reinit(const VectorType &vec)
239 {
240  diagonal = vec;
241 }
242 
243 
244 
245 template <typename VectorType>
246 void
248 {
249  dst.reinit(diagonal);
250 }
251 
252 
253 
254 template <typename VectorType>
255 void
257 {
258  diagonal.compress(operation);
259 }
260 
261 
262 
263 template <typename VectorType>
264 VectorType &
266 {
267  return diagonal;
268 }
269 
270 
271 
272 template <typename VectorType>
273 const VectorType &
275 {
276  return diagonal;
277 }
278 
279 
280 
281 template <typename VectorType>
282 typename VectorType::size_type
284 {
285  return diagonal.size();
286 }
287 
288 
289 
290 template <typename VectorType>
291 typename VectorType::size_type
293 {
294  return diagonal.size();
295 }
296 
297 
298 
299 template <typename VectorType>
300 typename VectorType::value_type
302  const size_type j) const
303 {
304  Assert(i == j, ExcIndexRange(j, i, i + 1));
305  (void)j;
306  return diagonal(i);
307 }
308 
309 
310 
311 template <typename VectorType>
312 typename VectorType::value_type &
313 DiagonalMatrix<VectorType>::operator()(const size_type i, const size_type j)
314 {
315  Assert(i == j, ExcIndexRange(j, i, i + 1));
316  (void)j;
317  return diagonal(i);
318 }
319 
320 
321 
322 template <typename VectorType>
323 template <typename number2>
324 void
325 DiagonalMatrix<VectorType>::add(const size_type row,
326  const size_type n_cols,
327  const size_type *col_indices,
328  const number2 * values,
329  const bool,
330  const bool)
331 {
332  for (size_type i = 0; i < n_cols; ++i)
333  if (col_indices[i] == row)
334  diagonal(row) += values[i];
335 }
336 
337 
338 
339 template <typename VectorType>
340 void
341 DiagonalMatrix<VectorType>::add(const size_type i,
342  const size_type j,
343  const value_type value)
344 {
345  if (i == j)
346  diagonal(i) += value;
347 }
348 
349 
350 
351 template <typename VectorType>
352 void
353 DiagonalMatrix<VectorType>::vmult(VectorType &dst, const VectorType &src) const
354 {
355  dst = src;
356  dst.scale(diagonal);
357 }
358 
359 
360 
361 template <typename VectorType>
362 void
363 DiagonalMatrix<VectorType>::Tvmult(VectorType &dst, const VectorType &src) const
364 {
365  vmult(dst, src);
366 }
367 
368 
369 
370 template <typename VectorType>
371 void
373  const VectorType &src) const
374 {
375  VectorType tmp(src);
376  tmp.scale(diagonal);
377  dst += tmp;
378 }
379 
380 
381 
382 template <typename VectorType>
383 void
385  const VectorType &src) const
386 {
387  vmult_add(dst, src);
388 }
389 
390 
391 #endif
392 
393 DEAL_II_NAMESPACE_CLOSE
394 
395 #endif
std::size_t memory_consumption() const
VectorType & get_vector()
void Tvmult(VectorType &dst, const VectorType &src) const
void compress(VectorOperation::values operation)
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
value_type operator()(const size_type i, const size_type j) const
void add(const size_type row, const size_type n_cols, const size_type *col_indices, const number2 *values, const bool elide_zero_values=true, const bool col_indices_are_sorted=false)
Matrix is diagonal.
size_type n() const
void vmult(VectorType &dst, const VectorType &src) const
#define Assert(cond, exc)
Definition: exceptions.h:1407
size_type m() const
VectorType diagonal
DiagonalMatrix()=default
void initialize_dof_vector(VectorType &dst) const
void reinit(const VectorType &vec)
void vmult_add(VectorType &dst, const VectorType &src) const
void Tvmult_add(VectorType &dst, const VectorType &src) const