Reference documentation for deal.II version Git 2276b94619 2020-02-28 10:11:55 +0100
\(\newcommand{\vcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\vcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
tensor_product_matrix.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2019 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_tensor_product_matrix_h
17 #define dealii_tensor_product_matrix_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/array_view.h>
23 #include <deal.II/base/thread_management.h>
24 
25 #include <deal.II/lac/lapack_full_matrix.h>
26 
27 #include <deal.II/matrix_free/tensor_product_kernels.h>
28 
29 DEAL_II_NAMESPACE_OPEN
30 
31 // Forward declarations
32 #ifndef DOXYGEN
33 template <typename>
34 class Vector;
35 template <typename>
36 class FullMatrix;
37 #endif
38 
75 template <int dim, typename Number, int n_rows_1d = -1>
77 {
78 public:
83  using value_type = Number;
84 
89  static constexpr int n_rows_1d_static = n_rows_1d;
90 
96  unsigned int
97  m() const;
98 
104  unsigned int
105  n() const;
106 
113  void
114  vmult(const ArrayView<Number> &dst, const ArrayView<const Number> &src) const;
115 
122  void
123  apply_inverse(const ArrayView<Number> & dst,
124  const ArrayView<const Number> &src) const;
125 
126 protected:
131 
135  std::array<Table<2, Number>, dim> mass_matrix;
136 
140  std::array<Table<2, Number>, dim> derivative_matrix;
141 
146  std::array<AlignedVector<Number>, dim> eigenvalues;
147 
152  std::array<Table<2, Number>, dim> eigenvectors;
153 
154 private:
159 
164 };
165 
166 
167 
241 template <int dim, typename Number, int n_rows_1d = -1>
243  : public TensorProductMatrixSymmetricSumBase<dim, Number, n_rows_1d>
244 {
245 public:
250 
258  const std::array<Table<2, Number>, dim> &mass_matrix,
259  const std::array<Table<2, Number>, dim> &derivative_matrix);
260 
268  const std::array<FullMatrix<Number>, dim> &mass_matrix,
269  const std::array<FullMatrix<Number>, dim> &derivative_matrix);
270 
277 
289  void
290  reinit(const std::array<Table<2, Number>, dim> &mass_matrix,
291  const std::array<Table<2, Number>, dim> &derivative_matrix);
292 
298  void
299  reinit(const std::array<FullMatrix<Number>, dim> &mass_matrix,
300  const std::array<FullMatrix<Number>, dim> &derivative_matrix);
301 
307  void
308  reinit(const Table<2, Number> &mass_matrix,
309  const Table<2, Number> &derivative_matrix);
310 
311 private:
320  template <typename MatrixArray>
321  void
322  reinit_impl(MatrixArray &&mass_matrix, MatrixArray &&derivative_matrix);
323 };
324 
325 
326 
335 template <int dim, typename Number, int n_rows_1d>
336 class TensorProductMatrixSymmetricSum<dim, VectorizedArray<Number>, n_rows_1d>
338  VectorizedArray<Number>,
339  n_rows_1d>
340 {
341 public:
346 
354  const std::array<Table<2, VectorizedArray<Number>>, dim> &mass_matrix,
355  const std::array<Table<2, VectorizedArray<Number>>, dim>
357 
367 
379  void
380  reinit(const std::array<Table<2, VectorizedArray<Number>>, dim> &mass_matrix,
381  const std::array<Table<2, VectorizedArray<Number>>, dim>
383 
389  void
390  reinit(const Table<2, VectorizedArray<Number>> &mass_matrix,
392 
393 private:
402  template <typename MatrixArray>
403  void
404  reinit_impl(MatrixArray &&mass_matrix, MatrixArray &&derivative_matrix);
405 };
406 
407 
408 /*----------------------- Inline functions ----------------------------------*/
409 
410 #ifndef DOXYGEN
411 
412 namespace internal
413 {
414  namespace TensorProductMatrix
415  {
424  template <typename Number>
425  void
426  spectral_assembly(const Number * mass_matrix,
427  const Number * derivative_matrix,
428  const unsigned int n_rows,
429  const unsigned int n_cols,
430  Number * eigenvalues,
431  Number * eigenvectors)
432  {
433  Assert(n_rows == n_cols, ExcNotImplemented());
434 
435  auto &&transpose_fill_nm = [](Number * out,
436  const Number * in,
437  const unsigned int n,
438  const unsigned int m) {
439  for (unsigned int mm = 0; mm < m; ++mm)
440  for (unsigned int nn = 0; nn < n; ++nn)
441  out[mm + nn * m] = *(in++);
442  };
443 
444  std::vector<::Vector<Number>> eigenvecs(n_rows);
445  LAPACKFullMatrix<Number> mass_copy(n_rows, n_cols);
446  LAPACKFullMatrix<Number> deriv_copy(n_rows, n_cols);
447 
448  transpose_fill_nm(&(mass_copy(0, 0)), mass_matrix, n_rows, n_cols);
449  transpose_fill_nm(&(deriv_copy(0, 0)), derivative_matrix, n_rows, n_cols);
450 
451  deriv_copy.compute_generalized_eigenvalues_symmetric(mass_copy,
452  eigenvecs);
453  AssertDimension(eigenvecs.size(), n_rows);
454  for (unsigned int i = 0; i < n_rows; ++i)
455  for (unsigned int j = 0; j < n_cols; ++j, ++eigenvectors)
456  *eigenvectors = eigenvecs[j][i];
457 
458  for (unsigned int i = 0; i < n_rows; ++i, ++eigenvalues)
459  *eigenvalues = deriv_copy.eigenvalue(i).real();
460  }
461  } // namespace TensorProductMatrix
462 } // namespace internal
463 
464 
465 template <int dim, typename Number, int n_rows_1d>
466 inline unsigned int
468 {
469  unsigned int m = mass_matrix[0].n_rows();
470  for (unsigned int d = 1; d < dim; ++d)
471  m *= mass_matrix[d].n_rows();
472  return m;
473 }
474 
475 
476 
477 template <int dim, typename Number, int n_rows_1d>
478 inline unsigned int
480 {
481  unsigned int n = mass_matrix[0].n_cols();
482  for (unsigned int d = 1; d < dim; ++d)
483  n *= mass_matrix[d].n_cols();
484  return n;
485 }
486 
487 
488 
489 template <int dim, typename Number, int n_rows_1d>
490 inline void
492  const ArrayView<Number> & dst_view,
493  const ArrayView<const Number> &src_view) const
494 {
495  AssertDimension(dst_view.size(), this->m());
496  AssertDimension(src_view.size(), this->n());
497  std::lock_guard<std::mutex> lock(this->mutex);
498  const unsigned int n = Utilities::fixed_power<dim>(
499  n_rows_1d > 0 ? n_rows_1d : eigenvalues[0].size());
500  tmp_array.resize_fast(n * 2);
501  constexpr int kernel_size = n_rows_1d > 0 ? n_rows_1d : 0;
503  dim,
504  kernel_size,
505  kernel_size,
506  Number>
507  eval(AlignedVector<Number>{},
510  mass_matrix[0].n_rows(),
511  mass_matrix[0].n_rows());
512  Number * t = tmp_array.begin();
513  const Number *src = src_view.begin();
514  Number * dst = dst_view.data();
515 
516  if (dim == 1)
517  {
518  const Number *A = &derivative_matrix[0](0, 0);
519  eval.template apply<0, false, false>(A, src, dst);
520  }
521 
522  else if (dim == 2)
523  {
524  const Number *A0 = &derivative_matrix[0](0, 0);
525  const Number *M0 = &mass_matrix[0](0, 0);
526  const Number *A1 = &derivative_matrix[1](0, 0);
527  const Number *M1 = &mass_matrix[1](0, 0);
528  eval.template apply<0, false, false>(M0, src, t);
529  eval.template apply<1, false, false>(A1, t, dst);
530  eval.template apply<0, false, false>(A0, src, t);
531  eval.template apply<1, false, true>(M1, t, dst);
532  }
533 
534  else if (dim == 3)
535  {
536  const Number *A0 = &derivative_matrix[0](0, 0);
537  const Number *M0 = &mass_matrix[0](0, 0);
538  const Number *A1 = &derivative_matrix[1](0, 0);
539  const Number *M1 = &mass_matrix[1](0, 0);
540  const Number *A2 = &derivative_matrix[2](0, 0);
541  const Number *M2 = &mass_matrix[2](0, 0);
542  eval.template apply<0, false, false>(M0, src, t + n);
543  eval.template apply<1, false, false>(M1, t + n, t);
544  eval.template apply<2, false, false>(A2, t, dst);
545  eval.template apply<1, false, false>(A1, t + n, t);
546  eval.template apply<0, false, false>(A0, src, t + n);
547  eval.template apply<1, false, true>(M1, t + n, t);
548  eval.template apply<2, false, true>(M2, t, dst);
549  }
550 
551  else
552  AssertThrow(false, ExcNotImplemented());
553 }
554 
555 
556 
557 template <int dim, typename Number, int n_rows_1d>
558 inline void
560  const ArrayView<Number> & dst_view,
561  const ArrayView<const Number> &src_view) const
562 {
563  AssertDimension(dst_view.size(), this->n());
564  AssertDimension(src_view.size(), this->m());
565  std::lock_guard<std::mutex> lock(this->mutex);
566  const unsigned int n = n_rows_1d > 0 ? n_rows_1d : eigenvalues[0].size();
567  tmp_array.resize_fast(Utilities::fixed_power<dim>(n));
568  constexpr int kernel_size = n_rows_1d > 0 ? n_rows_1d : 0;
570  dim,
571  kernel_size,
572  kernel_size,
573  Number>
574  eval(AlignedVector<Number>(),
577  mass_matrix[0].n_rows(),
578  mass_matrix[0].n_rows());
579  Number * t = tmp_array.begin();
580  const Number *src = src_view.data();
581  Number * dst = dst_view.data();
582 
583  // NOTE: dof_to_quad has to be interpreted as 'dof to eigenvalue index'
584  // --> apply<.,true,.> (S,src,dst) calculates dst = S^T * src,
585  // --> apply<.,false,.> (S,src,dst) calculates dst = S * src,
586  // while the eigenvectors are stored column-wise in S, i.e.
587  // rows correspond to dofs whereas columns to eigenvalue indices!
588  if (dim == 1)
589  {
590  const Number *S = &eigenvectors[0](0, 0);
591  eval.template apply<0, true, false>(S, src, t);
592  for (unsigned int i = 0; i < n; ++i)
593  t[i] /= eigenvalues[0][i];
594  eval.template apply<0, false, false>(S, t, dst);
595  }
596 
597  else if (dim == 2)
598  {
599  const Number *S0 = &(eigenvectors[0](0, 0));
600  const Number *S1 = &(eigenvectors[1](0, 0));
601  eval.template apply<0, true, false>(S0, src, t);
602  eval.template apply<1, true, false>(S1, t, dst);
603  for (unsigned int i1 = 0, c = 0; i1 < n; ++i1)
604  for (unsigned int i0 = 0; i0 < n; ++i0, ++c)
605  dst[c] /= (eigenvalues[1][i1] + eigenvalues[0][i0]);
606  eval.template apply<0, false, false>(S0, dst, t);
607  eval.template apply<1, false, false>(S1, t, dst);
608  }
609 
610  else if (dim == 3)
611  {
612  const Number *S0 = &eigenvectors[0](0, 0);
613  const Number *S1 = &eigenvectors[1](0, 0);
614  const Number *S2 = &eigenvectors[2](0, 0);
615  eval.template apply<0, true, false>(S0, src, t);
616  eval.template apply<1, true, false>(S1, t, dst);
617  eval.template apply<2, true, false>(S2, dst, t);
618  for (unsigned int i2 = 0, c = 0; i2 < n; ++i2)
619  for (unsigned int i1 = 0; i1 < n; ++i1)
620  for (unsigned int i0 = 0; i0 < n; ++i0, ++c)
621  t[c] /=
622  (eigenvalues[2][i2] + eigenvalues[1][i1] + eigenvalues[0][i0]);
623  eval.template apply<0, false, false>(S0, t, dst);
624  eval.template apply<1, false, false>(S1, dst, t);
625  eval.template apply<2, false, false>(S2, t, dst);
626  }
627 
628  else
629  Assert(false, ExcNotImplemented());
630 }
631 
632 
633 //---------------------- TensorProductMatrixSymmetricSum ----------------------
634 
635 template <int dim, typename Number, int n_rows_1d>
638  const std::array<Table<2, Number>, dim> &mass_matrix,
639  const std::array<Table<2, Number>, dim> &derivative_matrix)
640 {
642 }
643 
644 
645 
646 template <int dim, typename Number, int n_rows_1d>
649  const std::array<FullMatrix<Number>, dim> &mass_matrix,
650  const std::array<FullMatrix<Number>, dim> &derivative_matrix)
651 {
653 }
654 
655 
656 
657 template <int dim, typename Number, int n_rows_1d>
661 {
662  reinit(mass_matrix, derivative_matrix);
663 }
664 
665 
666 
667 template <int dim, typename Number, int n_rows_1d>
668 template <typename MatrixArray>
669 inline void
671  MatrixArray &&mass_matrices_,
672  MatrixArray &&derivative_matrices_)
673 {
674  auto &&mass_matrices = std::forward<MatrixArray>(mass_matrices_);
675  auto &&derivative_matrices = std::forward<MatrixArray>(derivative_matrices_);
676  this->mass_matrix = mass_matrices;
677  this->derivative_matrix = derivative_matrices;
678 
679  for (int dir = 0; dir < dim; ++dir)
680  {
681  Assert(n_rows_1d == -1 ||
682  (n_rows_1d > 0 && static_cast<unsigned int>(n_rows_1d) ==
683  mass_matrices[dir].n_rows()),
684  ExcDimensionMismatch(n_rows_1d, mass_matrices[dir].n_rows()));
685  AssertDimension(mass_matrices[dir].n_rows(), mass_matrices[dir].n_cols());
686  AssertDimension(mass_matrices[dir].n_rows(),
687  derivative_matrices[dir].n_rows());
688  AssertDimension(mass_matrices[dir].n_rows(),
689  derivative_matrices[dir].n_cols());
690 
691  this->eigenvectors[dir].reinit(mass_matrices[dir].n_cols(),
692  mass_matrices[dir].n_rows());
693  this->eigenvalues[dir].resize(mass_matrices[dir].n_cols());
694  internal::TensorProductMatrix ::spectral_assembly<Number>(
695  &(mass_matrices[dir](0, 0)),
696  &(derivative_matrices[dir](0, 0)),
697  mass_matrices[dir].n_rows(),
698  mass_matrices[dir].n_cols(),
699  this->eigenvalues[dir].begin(),
700  &(this->eigenvectors[dir](0, 0)));
701  }
702 }
703 
704 
705 
706 template <int dim, typename Number, int n_rows_1d>
707 inline void
709  const std::array<Table<2, Number>, dim> &mass_matrix,
710  const std::array<Table<2, Number>, dim> &derivative_matrix)
711 {
712  reinit_impl(mass_matrix, derivative_matrix);
713 }
714 
715 
716 
717 template <int dim, typename Number, int n_rows_1d>
718 inline void
720  const std::array<FullMatrix<Number>, dim> &mass_matrix,
721  const std::array<FullMatrix<Number>, dim> &derivative_matrix)
722 {
723  std::array<Table<2, Number>, dim> mass_copy;
724  std::array<Table<2, Number>, dim> deriv_copy;
725 
726  std::transform(mass_matrix.cbegin(),
727  mass_matrix.cend(),
728  mass_copy.begin(),
729  [](const FullMatrix<Number> &m) -> Table<2, Number> {
730  return m;
731  });
732  std::transform(derivative_matrix.cbegin(),
733  derivative_matrix.cend(),
734  deriv_copy.begin(),
735  [](const FullMatrix<Number> &m) -> Table<2, Number> {
736  return m;
737  });
738 
739  reinit_impl(std::move(mass_copy), std::move(deriv_copy));
740 }
741 
742 
743 
744 template <int dim, typename Number, int n_rows_1d>
745 inline void
747  const Table<2, Number> &mass_matrix,
748  const Table<2, Number> &derivative_matrix)
749 {
750  std::array<Table<2, Number>, dim> mass_matrices;
751  std::array<Table<2, Number>, dim> derivative_matrices;
752 
753  std::fill(mass_matrices.begin(), mass_matrices.end(), mass_matrix);
754  std::fill(derivative_matrices.begin(),
755  derivative_matrices.end(),
757 
758  reinit_impl(std::move(mass_matrices), std::move(derivative_matrices));
759 }
760 
761 
762 
763 //------------- vectorized spec.: TensorProductMatrixSymmetricSum -------------
764 
765 template <int dim, typename Number, int n_rows_1d>
768  n_rows_1d>::
769  TensorProductMatrixSymmetricSum(
770  const std::array<Table<2, VectorizedArray<Number>>, dim> &mass_matrix,
771  const std::array<Table<2, VectorizedArray<Number>>, dim> &derivative_matrix)
772 {
773  reinit(mass_matrix, derivative_matrix);
774 }
775 
776 
777 
778 template <int dim, typename Number, int n_rows_1d>
781  n_rows_1d>::
782  TensorProductMatrixSymmetricSum(
783  const Table<2, VectorizedArray<Number>> &mass_matrix,
784  const Table<2, VectorizedArray<Number>> &derivative_matrix)
785 {
786  reinit(mass_matrix, derivative_matrix);
787 }
788 
789 
790 
791 template <int dim, typename Number, int n_rows_1d>
792 template <typename MatrixArray>
793 inline void
795  reinit_impl(MatrixArray &&mass_matrices_, MatrixArray &&derivative_matrices_)
796 {
797  auto &&mass_matrix = std::forward<MatrixArray>(mass_matrices_);
798  auto &&derivative_matrix = std::forward<MatrixArray>(derivative_matrices_);
799  this->mass_matrix = mass_matrix;
800  this->derivative_matrix = derivative_matrix;
801 
802  constexpr unsigned int macro_size = VectorizedArray<Number>::n_array_elements;
803  std::size_t n_rows_max = (n_rows_1d > 0) ? n_rows_1d : 0;
804  if (n_rows_1d == -1)
805  for (unsigned int d = 0; d < dim; ++d)
806  n_rows_max = std::max(n_rows_max, mass_matrix[d].n_rows());
807  const std::size_t nm_flat_size_max = n_rows_max * n_rows_max * macro_size;
808  const std::size_t n_flat_size_max = n_rows_max * macro_size;
809 
810  std::vector<Number> mass_matrix_flat;
811  std::vector<Number> deriv_matrix_flat;
812  std::vector<Number> eigenvalues_flat;
813  std::vector<Number> eigenvectors_flat;
814  mass_matrix_flat.resize(nm_flat_size_max);
815  deriv_matrix_flat.resize(nm_flat_size_max);
816  eigenvalues_flat.resize(n_flat_size_max);
817  eigenvectors_flat.resize(nm_flat_size_max);
818  std::array<unsigned int, macro_size> offsets_nm;
819  std::array<unsigned int, macro_size> offsets_n;
820  for (int dir = 0; dir < dim; ++dir)
821  {
822  Assert(n_rows_1d == -1 ||
823  (n_rows_1d > 0 && static_cast<unsigned int>(n_rows_1d) ==
824  mass_matrix[dir].n_rows()),
825  ExcDimensionMismatch(n_rows_1d, mass_matrix[dir].n_rows()));
826  AssertDimension(mass_matrix[dir].n_rows(), mass_matrix[dir].n_cols());
827  AssertDimension(mass_matrix[dir].n_rows(),
828  derivative_matrix[dir].n_rows());
829  AssertDimension(mass_matrix[dir].n_rows(),
830  derivative_matrix[dir].n_cols());
831 
832  const unsigned int n_rows = mass_matrix[dir].n_rows();
833  const unsigned int n_cols = mass_matrix[dir].n_cols();
834  const unsigned int nm = n_rows * n_cols;
835  for (unsigned int vv = 0; vv < macro_size; ++vv)
836  offsets_nm[vv] = nm * vv;
837 
838  vectorized_transpose_and_store(false,
839  nm,
840  &(mass_matrix[dir](0, 0)),
841  offsets_nm.cbegin(),
842  mass_matrix_flat.data());
843  vectorized_transpose_and_store(false,
844  nm,
845  &(derivative_matrix[dir](0, 0)),
846  offsets_nm.cbegin(),
847  deriv_matrix_flat.data());
848 
849  const Number *mass_cbegin = mass_matrix_flat.data();
850  const Number *deriv_cbegin = deriv_matrix_flat.data();
851  Number * eigenvec_begin = eigenvectors_flat.data();
852  Number * eigenval_begin = eigenvalues_flat.data();
853  for (unsigned int lane = 0; lane < macro_size; ++lane)
854  internal::TensorProductMatrix ::spectral_assembly<Number>(
855  mass_cbegin + nm * lane,
856  deriv_cbegin + nm * lane,
857  n_rows,
858  n_cols,
859  eigenval_begin + n_rows * lane,
860  eigenvec_begin + nm * lane);
861 
862  this->eigenvalues[dir].resize(n_rows);
863  this->eigenvectors[dir].reinit(n_rows, n_cols);
864  for (unsigned int vv = 0; vv < macro_size; ++vv)
865  offsets_n[vv] = n_rows * vv;
866  vectorized_load_and_transpose(n_rows,
867  eigenvalues_flat.data(),
868  offsets_n.cbegin(),
869  this->eigenvalues[dir].begin());
870  vectorized_load_and_transpose(nm,
871  eigenvectors_flat.data(),
872  offsets_nm.cbegin(),
873  &(this->eigenvectors[dir](0, 0)));
874  }
875 }
876 
877 
878 
879 template <int dim, typename Number, int n_rows_1d>
880 inline void
882  reinit(
883  const std::array<Table<2, VectorizedArray<Number>>, dim> &mass_matrix,
884  const std::array<Table<2, VectorizedArray<Number>>, dim> &derivative_matrix)
885 {
886  reinit_impl(mass_matrix, derivative_matrix);
887 }
888 
889 
890 
891 template <int dim, typename Number, int n_rows_1d>
892 inline void
894  reinit(const Table<2, VectorizedArray<Number>> &mass_matrix,
895  const Table<2, VectorizedArray<Number>> &derivative_matrix)
896 {
897  std::array<Table<2, VectorizedArray<Number>>, dim> mass_matrices;
898  std::array<Table<2, VectorizedArray<Number>>, dim> derivative_matrices;
899 
900  std::fill(mass_matrices.begin(), mass_matrices.end(), mass_matrix);
901  std::fill(derivative_matrices.begin(),
902  derivative_matrices.end(),
904 
905  reinit_impl(std::move(mass_matrices), std::move(derivative_matrices));
906 }
907 
908 
909 
910 #endif
911 
912 DEAL_II_NAMESPACE_CLOSE
913 
914 #endif
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1571
std::array< AlignedVector< Number >, dim > eigenvalues
std::complex< number > eigenvalue(const size_type i) const
void vmult(const ArrayView< Number > &dst, const ArrayView< const Number > &src) const
#define AssertThrow(cond, exc)
Definition: exceptions.h:1523
LinearAlgebra::distributed::Vector< Number > Vector
std::size_t size() const
Definition: array_view.h:484
void apply_inverse(const ArrayView< Number > &dst, const ArrayView< const Number > &src) const
#define Assert(cond, exc)
Definition: exceptions.h:1411
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void reinit(const std::array< Table< 2, Number >, dim > &mass_matrix, const std::array< Table< 2, Number >, dim > &derivative_matrix)
value_type * data() const noexcept
Definition: array_view.h:461
void reinit_impl(MatrixArray &&mass_matrix, MatrixArray &&derivative_matrix)
std::array< Table< 2, Number >, dim > eigenvectors
void resize_fast(const size_type size)
iterator begin()
void compute_generalized_eigenvalues_symmetric(LAPACKFullMatrix< number > &B, const number lower_bound, const number upper_bound, const number abs_accuracy, Vector< number > &eigenvalues, std::vector< Vector< number >> &eigenvectors, const types::blas_int itype=1)
static ::ExceptionBase & ExcNotImplemented()
std::array< Table< 2, Number >, dim > derivative_matrix
iterator begin() const
Definition: array_view.h:493
std::array< Table< 2, Number >, dim > mass_matrix