Dune Core Modules (2.4.1)

dynmatrixev.hh
Go to the documentation of this file.
1// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=4 sw=2 sts=2:
3#ifndef DUNE_DYNMATRIXEIGENVALUES_HH
4#define DUNE_DYNMATRIXEIGENVALUES_HH
5
6#include <memory>
7
8#include <dune/common/std/memory.hh>
9
10#include "dynmatrix.hh"
11
20namespace Dune {
21
22 namespace DynamicMatrixHelp {
23
24 // defined in fmatrixev_ext.cpp
25 extern void eigenValuesNonsymLapackCall(
26 const char* jobvl, const char* jobvr, const long
27 int* n, double* a, const long int* lda, double* wr, double* wi, double* vl,
28 const long int* ldvl, double* vr, const long int* ldvr, double* work,
29 const long int* lwork, const long int* info);
30
38 template <typename K, class C>
39 static void eigenValuesNonSym(const DynamicMatrix<K>& matrix,
40 DynamicVector<C>& eigenValues)
41 {
42 {
43 const long int N = matrix.rows();
44 const char jobvl = 'n';
45 const char jobvr = 'n';
46
47
48 // matrix to put into dgeev
49 std::unique_ptr<double[]> matrixVector = Std::make_unique<double[]>(N*N);
50
51 // copy matrix
52 int row = 0;
53 for(int i=0; i<N; ++i)
54 {
55 for(int j=0; j<N; ++j, ++row)
56 {
57 matrixVector[ row ] = matrix[ i ][ j ];
58 }
59 }
60
61 // working memory
62 std::unique_ptr<double[]> eigenR = Std::make_unique<double[]>(N);
63 std::unique_ptr<double[]> eigenI = Std::make_unique<double[]>(N);
64 std::unique_ptr<double[]> work = Std::make_unique<double[]>(3*N);
65
66 // return value information
67 long int info = 0;
68 long int lwork = 3*N;
69
70 // call LAPACK routine (see fmatrixev_ext.cc)
71 eigenValuesNonsymLapackCall(&jobvl, &jobvr, &N, &matrixVector[0], &N,
72 &eigenR[0], &eigenI[0], 0, &N, 0, &N, &work[0],
73 &lwork, &info);
74
75 if( info != 0 )
76 {
77 std::cerr << "For matrix " << matrix << " eigenvalue calculation failed! " << std::endl;
78 DUNE_THROW(InvalidStateException,"eigenValues: Eigenvalue calculation failed!");
79 }
80
81 eigenValues.resize(N);
82 for (int i=0; i<N; ++i)
83 eigenValues[i] = std::complex<double>(eigenR[i], eigenI[i]);
84 }
85 }
86
87 }
88
89}
91#endif
size_type rows() const
number of rows
Definition: densematrix.hh:695
Construct a matrix with a dynamic size.
Definition: dynmatrix.hh:56
Construct a vector with a dynamic size.
Definition: dynvector.hh:55
Default exception if a function was called while the object is not in a valid state for that function...
Definition: exceptions.hh:306
This file implements a dense matrix with dynamic numbers of rows and columns.
static void eigenValuesNonSym(const DynamicMatrix< K > &matrix, DynamicVector< C > &eigenValues)
calculates the eigenvalues of a symmetric field matrix
Definition: dynmatrixev.hh:39
#define DUNE_THROW(E, m)
Definition: exceptions.hh:243
Dune namespace.
Definition: alignment.hh:10
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 24, 23:30, 2024)