Reference documentation for deal.II version GIT 29f9da0a34 2023-12-07 10:00:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
tridiagonal_matrix.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2005 - 2023 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_tridiagonal_matrix_h
17 #define dealii_tridiagonal_matrix_h
18 
19 #include <deal.II/base/config.h>
20 
22 
24 
25 #include <iomanip>
26 #include <vector>
27 
29 
30 // forward declarations
31 #ifndef DOXYGEN
32 template <typename number>
33 class Vector;
34 #endif
35 
52 template <typename number>
54 {
55 public:
57 
62 
69  TridiagonalMatrix(size_type n = 0, bool symmetric = false);
70 
75  void
76  reinit(size_type n, bool symmetric = false);
77 
78 
82 
88  size_type
89  m() const;
90 
95  size_type
96  n() const;
97 
103  bool
104  all_zero() const;
105 
109 
114  number
116 
126  number &
128 
132 
143  void
145  const Vector<number> &v,
146  const bool adding = false) const;
147 
154  void
155  vmult_add(Vector<number> &w, const Vector<number> &v) const;
156 
166  void
168  const Vector<number> &v,
169  const bool adding = false) const;
170 
178  void
179  Tvmult_add(Vector<number> &w, const Vector<number> &v) const;
180 
186  number
187  matrix_scalar_product(const Vector<number> &u, const Vector<number> &v) const;
188 
198  number
199  matrix_norm_square(const Vector<number> &v) const;
200 
204 
211  void
216  number
217  eigenvalue(const size_type i) const;
221 
225  template <class OutputStream>
226  void
227  print(OutputStream &s,
228  const unsigned int width = 5,
229  const unsigned int precision = 2) const;
232 private:
236  std::vector<number> diagonal;
246  std::vector<number> left;
252  std::vector<number> right;
253 
259 
267 };
268 
271 //---------------------------------------------------------------------------
272 #ifndef DOXYGEN
273 
274 template <typename number>
277 {
278  return diagonal.size();
279 }
280 
281 
282 
283 template <typename number>
286 {
287  return diagonal.size();
288 }
289 
290 
291 template <typename number>
292 inline number
294 {
295  AssertIndexRange(i, n());
296  AssertIndexRange(j, n());
297  Assert(i <= j + 1, ExcIndexRange(i, j - 1, j + 2));
298  Assert(j <= i + 1, ExcIndexRange(j, i - 1, i + 2));
299 
300  if (j == i)
301  return diagonal[i];
302  if (j == i - 1)
303  {
304  if (is_symmetric)
305  return right[i - 1];
306  else
307  return left[i];
308  }
309 
310  if (j == i + 1)
311  return right[i];
312 
313  Assert(false, ExcInternalError());
314  return 0;
315 }
316 
317 
318 template <typename number>
319 inline number &
321 {
322  AssertIndexRange(i, n());
323  AssertIndexRange(j, n());
324  Assert(i <= j + 1, ExcIndexRange(i, j - 1, j + 2));
325  Assert(j <= i + 1, ExcIndexRange(j, i - 1, i + 2));
326 
327  if (j == i)
328  return diagonal[i];
329  if (j == i - 1)
330  {
331  if (is_symmetric)
332  return right[i - 1];
333  else
334  return left[i];
335  }
336 
337  if (j == i + 1)
338  return right[i];
339 
340  Assert(false, ExcInternalError());
341  return diagonal[0];
342 }
343 
344 
345 template <typename number>
346 template <class OutputStream>
347 void
348 TridiagonalMatrix<number>::print(OutputStream &s,
349  const unsigned int width,
350  const unsigned int) const
351 {
352  for (size_type i = 0; i < n(); ++i)
353  {
354  if (i > 0)
355  s << std::setw(width) << (*this)(i, i - 1);
356  else
357  s << std::setw(width) << "";
358 
359  s << ' ' << (*this)(i, i) << ' ';
360 
361  if (i < n() - 1)
362  s << std::setw(width) << (*this)(i, i + 1);
363 
364  s << std::endl;
365  }
366 }
367 
368 
369 #endif // DOXYGEN
370 
372 
373 #endif
number operator()(size_type i, size_type j) const
LAPACKSupport::State state
void Tvmult_add(Vector< number > &w, const Vector< number > &v) const
void vmult(Vector< number > &w, const Vector< number > &v, const bool adding=false) const
void print(OutputStream &s, const unsigned int width=5, const unsigned int precision=2) const
size_type n() const
number & operator()(size_type i, size_type j)
number matrix_norm_square(const Vector< number > &v) const
std::vector< number > diagonal
void reinit(size_type n, bool symmetric=false)
void vmult_add(Vector< number > &w, const Vector< number > &v) const
std::vector< number > right
types::global_dof_index size_type
void Tvmult(Vector< number > &w, const Vector< number > &v, const bool adding=false) const
size_type m() const
number eigenvalue(const size_type i) const
TridiagonalMatrix(size_type n=0, bool symmetric=false)
number matrix_scalar_product(const Vector< number > &u, const Vector< number > &v) const
std::vector< number > left
Definition: vector.h:110
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1631
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
static ::ExceptionBase & ExcIndexRange(std::size_t arg1, std::size_t arg2, std::size_t arg3)
@ symmetric
Matrix is symmetric.
@ diagonal
Matrix is diagonal.
types::global_dof_index size_type
Definition: cuda_kernels.h:45
Tensor< 2, dim, Number > w(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
unsigned int global_dof_index
Definition: types.h:82