Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
sparse_direct.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2001 - 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
18
22#include <deal.II/lac/vector.h>
23
24#include <cerrno>
25#include <complex>
26#include <iostream>
27#include <list>
28#include <typeinfo>
29#include <vector>
30
31
33
34
35// include UMFPACK file.
36#ifdef DEAL_II_WITH_UMFPACK
37# include <umfpack.h>
38#endif
39
40
41namespace
42{
50 template <typename SparseMatrixType>
51 unsigned int
52 parallel_grainsize(const SparseMatrixType &matrix)
53 {
54 const unsigned int avg_entries_per_row =
55 matrix.n_nonzero_elements() / matrix.m();
56 return std::max(1000 / avg_entries_per_row, 1u);
57 }
58} // namespace
59
60
62{
63 clear();
64}
65
66
67void
69{}
70
71
72#ifdef DEAL_II_WITH_UMFPACK
73
75 : n_rows(0)
76 , n_cols(0)
77 , symbolic_decomposition(nullptr)
78 , numeric_decomposition(nullptr)
79 , control(UMFPACK_CONTROL)
80{
81 umfpack_dl_defaults(control.data());
82}
83
84
85
86void
88{
89 // delete objects that haven't been deleted yet
90 if (symbolic_decomposition != nullptr)
91 {
92 umfpack_dl_free_symbolic(&symbolic_decomposition);
93 symbolic_decomposition = nullptr;
94 }
95
96 if (numeric_decomposition != nullptr)
97 {
98 umfpack_dl_free_numeric(&numeric_decomposition);
99 numeric_decomposition = nullptr;
100 }
101
102 {
103 std::vector<types::suitesparse_index> tmp;
104 tmp.swap(Ap);
105 }
106
107 {
108 std::vector<types::suitesparse_index> tmp;
109 tmp.swap(Ai);
110 }
111
112 {
113 std::vector<double> tmp;
114 tmp.swap(Ax);
115 }
116
117 {
118 std::vector<double> tmp;
119 tmp.swap(Az);
120 }
121
122 umfpack_dl_defaults(control.data());
123}
124
125
126
127template <typename number>
128void
130{
131 // do the copying around of entries so that the diagonal entry is in the
132 // right place. note that this is easy to detect: since all entries apart
133 // from the diagonal entry are sorted, we know that the diagonal entry is
134 // in the wrong place if and only if its column index is larger than the
135 // column index of the second entry in a row
136 //
137 // ignore rows with only one or no entry
139 0,
140 matrix.m(),
141 [this](const size_type row_begin, const size_type row_end) {
142 for (size_type row = row_begin; row < row_end; ++row)
143 {
144 // we may have to move some elements that are left of the diagonal
145 // but presently after the diagonal entry to the left, whereas the
146 // diagonal entry has to move to the right. we could first figure out
147 // where to move everything to, but for simplicity we just make a
148 // series of swaps instead (this is kind of a single run of
149 // bubble-sort, which gives us the desired result since the array is
150 // already "almost" sorted)
151 //
152 // in the first loop, the condition in the while-header also checks
153 // that the row has at least two entries and that the diagonal entry
154 // is really in the wrong place
155 long int cursor = Ap[row];
156 while ((cursor < Ap[row + 1] - 1) && (Ai[cursor] > Ai[cursor + 1]))
157 {
158 std::swap(Ai[cursor], Ai[cursor + 1]);
159
160 std::swap(Ax[cursor], Ax[cursor + 1]);
161 if (numbers::NumberTraits<number>::is_complex == true)
162 std::swap(Az[cursor], Az[cursor + 1]);
163
164 ++cursor;
165 }
166 }
167 },
168 parallel_grainsize(matrix));
169}
170
171
172
173template <typename number>
174void
176{
177 // same thing for SparseMatrixEZ
179 0,
180 matrix.m(),
181 [this](const size_type row_begin, const size_type row_end) {
182 for (size_type row = row_begin; row < row_end; ++row)
183 {
184 long int cursor = Ap[row];
185 while ((cursor < Ap[row + 1] - 1) && (Ai[cursor] > Ai[cursor + 1]))
186 {
187 std::swap(Ai[cursor], Ai[cursor + 1]);
188
189 std::swap(Ax[cursor], Ax[cursor + 1]);
190 if (numbers::NumberTraits<number>::is_complex == true)
191 std::swap(Az[cursor], Az[cursor + 1]);
192
193 ++cursor;
194 }
195 }
196 },
197 parallel_grainsize(matrix));
198}
199
200
201
202template <typename number>
203void
205{
206 // the case for block matrices is a bit more difficult, since all we know
207 // is that *within each block*, the diagonal of that block may come
208 // first. however, that means that there may be as many entries per row
209 // in the wrong place as there are block columns. we can do the same
210 // thing as above, but we have to do it multiple times
212 0,
213 matrix.m(),
214 [this, &matrix](const size_type row_begin, const size_type row_end) {
215 for (size_type row = row_begin; row < row_end; ++row)
216 {
217 long int cursor = Ap[row];
218 for (size_type block = 0; block < matrix.n_block_cols(); ++block)
219 {
220 // find the next out-of-order element
221 while ((cursor < Ap[row + 1] - 1) &&
222 (Ai[cursor] < Ai[cursor + 1]))
223 ++cursor;
224
225 // if there is none, then just go on
226 if (cursor == Ap[row + 1] - 1)
227 break;
228
229 // otherwise swap this entry with successive ones as long as
230 // necessary
231 long int element = cursor;
232 while ((element < Ap[row + 1] - 1) &&
233 (Ai[element] > Ai[element + 1]))
234 {
235 std::swap(Ai[element], Ai[element + 1]);
236
237 std::swap(Ax[element], Ax[element + 1]);
238 if (numbers::NumberTraits<number>::is_complex == true)
239 std::swap(Az[element], Az[element + 1]);
240
241 ++element;
242 }
243 }
244 }
245 },
246 parallel_grainsize(matrix));
247}
248
249
250
251template <class Matrix>
252void
254{
255 Assert(matrix.m() == matrix.n(), ExcNotQuadratic());
256
257 clear();
258
259 using number = typename Matrix::value_type;
260
261 n_rows = matrix.m();
262 n_cols = matrix.n();
263
264 const size_type N = matrix.m();
265
266 // copy over the data from the matrix to the data structures UMFPACK
267 // wants. note two things: first, UMFPACK wants compressed column storage
268 // whereas we always do compressed row storage; we work around this by,
269 // rather than shuffling things around, copy over the data we have, but
270 // then call the umfpack_dl_solve function with the UMFPACK_At argument,
271 // meaning that we want to solve for the transpose system
272 //
273 // second: the data we have in the sparse matrices is "almost" right
274 // already; UMFPACK wants the entries in each row (i.e. really: column)
275 // to be sorted in ascending order. we almost have that, except that we
276 // usually store the diagonal first in each row to allow for some
277 // optimizations. thus, we have to resort things a little bit, but only
278 // within each row
279 //
280 // final note: if the matrix has entries in the sparsity pattern that are
281 // actually occupied by entries that have a zero numerical value, then we
282 // keep them anyway. people are supposed to provide accurate sparsity
283 // patterns.
284 Ap.resize(N + 1);
285 Ai.resize(matrix.n_nonzero_elements());
286 Ax.resize(matrix.n_nonzero_elements());
288 Az.resize(matrix.n_nonzero_elements());
289
290 // first fill row lengths array
291 Ap[0] = 0;
292 for (size_type row = 1; row <= N; ++row)
293 Ap[row] = Ap[row - 1] + matrix.get_row_length(row - 1);
294 Assert(static_cast<size_type>(Ap.back()) == Ai.size(), ExcInternalError());
295
296 // then copy over matrix elements. note that for sparse matrices,
297 // iterators are sorted so that they traverse each row from start to end
298 // before moving on to the next row.
300 0,
301 matrix.m(),
302 [this, &matrix](const size_type row_begin, const size_type row_end) {
303 for (size_type row = row_begin; row < row_end; ++row)
304 {
305 long int index = Ap[row];
306 for (typename Matrix::const_iterator p = matrix.begin(row);
307 p != matrix.end(row);
308 ++p)
309 {
310 // write entry into the first free one for this row
311 Ai[index] = p->column();
312 Ax[index] = std::real(p->value());
313 if (numbers::NumberTraits<number>::is_complex == true)
314 Az[index] = std::imag(p->value());
315
316 // then move pointer ahead
317 ++index;
318 }
319 Assert(index == Ap[row + 1], ExcInternalError());
320 }
321 },
322 parallel_grainsize(matrix));
323
324 // make sure that the elements in each row are sorted. we have to be more
325 // careful for block sparse matrices, so ship this task out to a
326 // different function
327 sort_arrays(matrix);
328
329 int status;
331 status = umfpack_dl_symbolic(N,
332 N,
333 Ap.data(),
334 Ai.data(),
335 Ax.data(),
336 &symbolic_decomposition,
337 control.data(),
338 nullptr);
339 else
340 status = umfpack_zl_symbolic(N,
341 N,
342 Ap.data(),
343 Ai.data(),
344 Ax.data(),
345 Az.data(),
346 &symbolic_decomposition,
347 control.data(),
348 nullptr);
349 AssertThrow(status == UMFPACK_OK,
350 ExcUMFPACKError("umfpack_dl_symbolic", status));
351
353 status = umfpack_dl_numeric(Ap.data(),
354 Ai.data(),
355 Ax.data(),
356 symbolic_decomposition,
357 &numeric_decomposition,
358 control.data(),
359 nullptr);
360 else
361 status = umfpack_zl_numeric(Ap.data(),
362 Ai.data(),
363 Ax.data(),
364 Az.data(),
365 symbolic_decomposition,
366 &numeric_decomposition,
367 control.data(),
368 nullptr);
369 AssertThrow(status == UMFPACK_OK,
370 ExcUMFPACKError("umfpack_dl_numeric", status));
371
372 umfpack_dl_free_symbolic(&symbolic_decomposition);
373}
374
375
376
377void
379 const bool transpose /*=false*/) const
380{
381 // make sure that some kind of factorize() call has happened before
382 Assert(Ap.size() != 0, ExcNotInitialized());
383 Assert(Ai.size() != 0, ExcNotInitialized());
384 Assert(Ai.size() == Ax.size(), ExcNotInitialized());
385
386 Assert(Az.size() == 0,
387 ExcMessage("You have previously factored a matrix using this class "
388 "that had complex-valued entries. This then requires "
389 "applying the factored matrix to a complex-valued "
390 "vector, but you are only providing a real-valued vector "
391 "here."));
392
393 Vector<double> rhs(rhs_and_solution.size());
394 rhs = rhs_and_solution;
395
396 // solve the system. note that since UMFPACK wants compressed column
397 // storage instead of the compressed row storage format we use in
398 // deal.II's SparsityPattern classes, we solve for UMFPACK's A^T instead
399
400 // Conversely, if we solve for the transpose, we have to use UMFPACK_A
401 // instead.
402 const int status = umfpack_dl_solve(transpose ? UMFPACK_A : UMFPACK_At,
403 Ap.data(),
404 Ai.data(),
405 Ax.data(),
406 rhs_and_solution.begin(),
407 rhs.begin(),
409 control.data(),
410 nullptr);
411 AssertThrow(status == UMFPACK_OK,
412 ExcUMFPACKError("umfpack_dl_solve", status));
413}
414
415
416
417void
418SparseDirectUMFPACK::solve(Vector<std::complex<double>> &rhs_and_solution,
419 const bool transpose /*=false*/) const
420{
421# ifdef DEAL_II_WITH_COMPLEX_VALUES
422 // make sure that some kind of factorize() call has happened before
423 Assert(Ap.size() != 0, ExcNotInitialized());
424 Assert(Ai.size() != 0, ExcNotInitialized());
425 Assert(Ai.size() == Ax.size(), ExcNotInitialized());
426
427 // First see whether the matrix that was factorized was complex-valued.
428 // If so, just apply the complex factorization to the vector.
429 if (Az.size() != 0)
430 {
431 Assert(Ax.size() == Az.size(), ExcInternalError());
432
433 // It would be nice if we could just present a pointer to the
434 // first element of the complex-valued solution vector and let
435 // UMFPACK fill both the real and imaginary parts of the solution
436 // vector at that address. UMFPACK calls this the 'packed' format,
437 // and in those cases it only takes one pointer to the entire
438 // vector, rather than a pointer to the real and one pointer to
439 // the imaginary parts of the vector. The problem is that if we
440 // want to pack, then we also need to pack the matrix, and the
441 // functions above have already decided that we don't want to pack
442 // the matrix but instead deal with split format for the matrix,
443 // and then UMFPACK decides that it can't deal with a split matrix
444 // and a packed vector. We have to choose one or the other, not
445 // mix.
446 //
447 // So create four vectors, one each for the real and imaginary parts
448 // of the right hand side and solution.
449 Vector<double> rhs_re(rhs_and_solution.size());
450 Vector<double> rhs_im(rhs_and_solution.size());
451 for (unsigned int i = 0; i < rhs_and_solution.size(); ++i)
452 {
453 rhs_re(i) = std::real(rhs_and_solution(i));
454 rhs_im(i) = std::imag(rhs_and_solution(i));
455 }
456
457 Vector<double> solution_re(rhs_and_solution.size());
458 Vector<double> solution_im(rhs_and_solution.size());
459
460 // Solve the system. note that since UMFPACK wants compressed column
461 // storage instead of the compressed row storage format we use in
462 // deal.II's SparsityPattern classes, we solve for UMFPACK's A^T instead
463 //
464 // Conversely, if we solve for the transpose, we have to use UMFPACK_A
465 // instead.
466 //
467 // Note that for the complex case here, the transpose is selected using
468 // UMFPACK_Aat, not UMFPACK_At.
469 const int status = umfpack_zl_solve(transpose ? UMFPACK_A : UMFPACK_Aat,
470 Ap.data(),
471 Ai.data(),
472 Ax.data(),
473 Az.data(),
474 solution_re.data(),
475 solution_im.data(),
476 rhs_re.data(),
477 rhs_im.data(),
479 control.data(),
480 nullptr);
481 AssertThrow(status == UMFPACK_OK,
482 ExcUMFPACKError("umfpack_dl_solve", status));
483
484 // Now put things back together into the output vector
485 for (unsigned int i = 0; i < rhs_and_solution.size(); ++i)
486 rhs_and_solution(i) = {solution_re(i), solution_im(i)};
487 }
488 else
489 {
490 // We have factorized a real-valued matrix, but the rhs and solution
491 // vectors are complex-valued. UMFPACK does not natively support this
492 // case, but we can just apply the factorization to real and imaginary
493 // parts of the right hand side separately
494 const Vector<std::complex<double>> rhs = rhs_and_solution;
495
496 // Get the real part of the right hand side, solve with it, and copy the
497 // results into the result vector by just copying the real output
498 // into the complex-valued result vector (which implies setting the
499 // imaginary parts to zero):
500 Vector<double> rhs_real_or_imag(rhs_and_solution.size());
501 for (unsigned int i = 0; i < rhs.size(); ++i)
502 rhs_real_or_imag(i) = std::real(rhs(i));
503
504 solve(rhs_real_or_imag, transpose);
505
506 rhs_and_solution = rhs_real_or_imag;
507
508 // Then repeat the whole thing with the imaginary part. The copying step
509 // is more complicated now because we can only touch the imaginary
510 // component of the output vector:
511 for (unsigned int i = 0; i < rhs.size(); ++i)
512 rhs_real_or_imag(i) = std::imag(rhs(i));
513
514 solve(rhs_real_or_imag, transpose);
515
516 for (unsigned int i = 0; i < rhs.size(); ++i)
517 rhs_and_solution(i).imag(rhs_real_or_imag(i));
518 }
519
520# else
521
522 (void)rhs_and_solution;
523 (void)transpose;
524 Assert(false,
526 "This function can't be called if deal.II has been configured "
527 "with DEAL_II_WITH_COMPLEX_VALUES=FALSE."));
528# endif
529}
530
531
532void
534 const bool transpose /*=false*/) const
535{
536 // the UMFPACK functions want a contiguous array of elements, so
537 // there is no way around copying data around. thus, just copy the
538 // data into a regular vector and back
539 Vector<double> tmp(rhs_and_solution.size());
540 tmp = rhs_and_solution;
541 solve(tmp, transpose);
542 rhs_and_solution = tmp;
543}
544
545
546
547void
548SparseDirectUMFPACK::solve(BlockVector<std::complex<double>> &rhs_and_solution,
549 const bool transpose /*=false*/) const
550{
551# ifdef DEAL_II_WITH_COMPLEX_VALUES
552 // the UMFPACK functions want a contiguous array of elements, so
553 // there is no way around copying data around. thus, just copy the
554 // data into a regular vector and back
555 Vector<std::complex<double>> tmp(rhs_and_solution.size());
556 tmp = rhs_and_solution;
557 solve(tmp, transpose);
558 rhs_and_solution = tmp;
559
560# else
561 (void)rhs_and_solution;
562 (void)transpose;
563 Assert(false,
565 "This function can't be called if deal.II has been configured "
566 "with DEAL_II_WITH_COMPLEX_VALUES=FALSE."));
567# endif
568}
569
570
571
572template <class Matrix>
573void
574SparseDirectUMFPACK::solve(const Matrix & matrix,
575 Vector<double> &rhs_and_solution,
576 const bool transpose /*=false*/)
577{
578 factorize(matrix);
579 solve(rhs_and_solution, transpose);
580}
581
582
583
584template <class Matrix>
585void
586SparseDirectUMFPACK::solve(const Matrix & matrix,
587 Vector<std::complex<double>> &rhs_and_solution,
588 const bool transpose /*=false*/)
589{
590# ifdef DEAL_II_WITH_COMPLEX_VALUES
591 factorize(matrix);
592 solve(rhs_and_solution, transpose);
593
594# else
595
596 (void)matrix;
597 (void)rhs_and_solution;
598 (void)transpose;
599 Assert(false,
601 "This function can't be called if deal.II has been configured "
602 "with DEAL_II_WITH_COMPLEX_VALUES=FALSE."));
603# endif
604}
605
606
607
608template <class Matrix>
609void
610SparseDirectUMFPACK::solve(const Matrix & matrix,
611 BlockVector<double> &rhs_and_solution,
612 const bool transpose /*=false*/)
613{
614 factorize(matrix);
615 solve(rhs_and_solution, transpose);
616}
617
618
619
620template <class Matrix>
621void
622SparseDirectUMFPACK::solve(const Matrix & matrix,
623 BlockVector<std::complex<double>> &rhs_and_solution,
624 const bool transpose /*=false*/)
625{
626# ifdef DEAL_II_WITH_COMPLEX_VALUES
627 factorize(matrix);
628 solve(rhs_and_solution, transpose);
629
630# else
631
632 (void)matrix;
633 (void)rhs_and_solution;
634 (void)transpose;
635 Assert(false,
637 "This function can't be called if deal.II has been configured "
638 "with DEAL_II_WITH_COMPLEX_VALUES=FALSE."));
639# endif
640}
641
642
643#else
644
645
647 : n_rows(0)
648 , n_cols(0)
649 , symbolic_decomposition(nullptr)
650 , numeric_decomposition(nullptr)
651 , control(0)
652{}
653
654
655void
657{}
658
659
660template <class Matrix>
661void
662SparseDirectUMFPACK::factorize(const Matrix &)
663{
665 false,
667 "To call this function you need UMFPACK, but you configured deal.II "
668 "without passing the necessary switch to 'cmake'. Please consult the "
669 "installation instructions at https://dealii.org/current/readme.html"));
670}
671
672
673void
675{
677 false,
679 "To call this function you need UMFPACK, but you configured deal.II "
680 "without passing the necessary switch to 'cmake'. Please consult the "
681 "installation instructions at https://dealii.org/current/readme.html"));
682}
683
684
685
686void
687SparseDirectUMFPACK::solve(Vector<std::complex<double>> &, const bool) const
688{
690 false,
692 "To call this function you need UMFPACK, but you configured deal.II "
693 "without passing the necessary switch to 'cmake'. Please consult the "
694 "installation instructions at https://dealii.org/current/readme.html"));
695}
696
697
698
699void
701{
703 false,
705 "To call this function you need UMFPACK, but you configured deal.II "
706 "without passing the necessary switch to 'cmake'. Please consult the "
707 "installation instructions at https://dealii.org/current/readme.html"));
708}
709
710
711
712void
713SparseDirectUMFPACK::solve(BlockVector<std::complex<double>> &,
714 const bool) const
715{
717 false,
719 "To call this function you need UMFPACK, but you configured deal.II "
720 "without passing the necessary switch to 'cmake'. Please consult the "
721 "installation instructions at https://dealii.org/current/readme.html"));
722}
723
724
725
726template <class Matrix>
727void
728SparseDirectUMFPACK::solve(const Matrix &, Vector<double> &, const bool)
729{
731 false,
733 "To call this function you need UMFPACK, but you configured deal.II "
734 "without passing the necessary switch to 'cmake'. Please consult the "
735 "installation instructions at https://dealii.org/current/readme.html"));
736}
737
738
739
740template <class Matrix>
741void
742SparseDirectUMFPACK::solve(const Matrix &,
743 Vector<std::complex<double>> &,
744 const bool)
745{
747 false,
749 "To call this function you need UMFPACK, but you configured deal.II "
750 "without passing the necessary switch to 'cmake'. Please consult the "
751 "installation instructions at https://dealii.org/current/readme.html"));
752}
753
754
755
756template <class Matrix>
757void
758SparseDirectUMFPACK::solve(const Matrix &, BlockVector<double> &, const bool)
759{
761 false,
763 "To call this function you need UMFPACK, but you configured deal.II "
764 "without passing the necessary switch to 'cmake'. Please consult the "
765 "installation instructions at https://dealii.org/current/readme.html"));
766}
767
768
769
770template <class Matrix>
771void
772SparseDirectUMFPACK::solve(const Matrix &,
773 BlockVector<std::complex<double>> &,
774 const bool)
775{
777 false,
779 "To call this function you need UMFPACK, but you configured deal.II "
780 "without passing the necessary switch to 'cmake'. Please consult the "
781 "installation instructions at https://dealii.org/current/readme.html"));
782}
783
784#endif
785
786
787template <class Matrix>
788void
790{
791 this->factorize(M);
792}
793
794
795void
797{
798 dst = src;
799 this->solve(dst);
800}
801
802
803
804void
806 const BlockVector<double> &src) const
807{
808 dst = src;
809 this->solve(dst);
810}
811
812
813void
815 const Vector<double> &src) const
816{
817 dst = src;
818 this->solve(dst, /*transpose=*/true);
819}
820
821
822
823void
825 const BlockVector<double> &src) const
826{
827 dst = src;
828 this->solve(dst, /*transpose=*/true);
829}
830
833{
835 return n_rows;
836}
837
840{
842 return n_cols;
843}
844
845
846// explicit instantiations for SparseMatrixUMFPACK
847#define InstantiateUMFPACK(MatrixType) \
848 template void SparseDirectUMFPACK::factorize(const MatrixType &); \
849 template void SparseDirectUMFPACK::solve(const MatrixType &, \
850 Vector<double> &, \
851 const bool); \
852 template void SparseDirectUMFPACK::solve(const MatrixType &, \
853 Vector<std::complex<double>> &, \
854 const bool); \
855 template void SparseDirectUMFPACK::solve(const MatrixType &, \
856 BlockVector<double> &, \
857 const bool); \
858 template void SparseDirectUMFPACK::solve( \
859 const MatrixType &, BlockVector<std::complex<double>> &, const bool); \
860 template void SparseDirectUMFPACK::initialize(const MatrixType &, \
861 const AdditionalData)
862
863// Instantiate everything for real-valued matrices
870
871// Now also for complex-valued matrices
872#ifdef DEAL_II_WITH_COMPLEX_VALUES
873InstantiateUMFPACK(SparseMatrix<std::complex<double>>);
874InstantiateUMFPACK(SparseMatrix<std::complex<float>>);
877#endif
878
std::size_t size() const
std::vector< double > Az
~SparseDirectUMFPACK() override
void initialize(const SparsityPattern &sparsity_pattern)
void Tvmult(Vector< double > &dst, const Vector< double > &src) const
size_type m() const
void solve(Vector< double > &rhs_and_solution, const bool transpose=false) const
std::vector< double > Ax
void sort_arrays(const SparseMatrixEZ< number > &)
size_type n() const
void factorize(const Matrix &matrix)
std::vector< double > control
void vmult(Vector< double > &dst, const Vector< double > &src) const
std::vector< types::suitesparse_index > Ap
std::vector< types::suitesparse_index > Ai
pointer data()
size_type size() const
iterator begin()
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
static ::ExceptionBase & ExcUMFPACKError(std::string arg1, int arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNotInitialized()
static ::ExceptionBase & ExcNotQuadratic()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
@ matrix
Contents is actually a matrix.
void apply_to_subranges(const Iterator &begin, const std_cxx20::type_identity_t< Iterator > &end, const Function &f, const unsigned int grainsize)
Definition parallel.h:435
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
#define InstantiateUMFPACK(MatrixType)