Linear solvers
Table of Contents
Introduction
Solvers of linear systems are one of the most important algorithms in scientific computations. TNL offers the following iterative methods:
- Stationary methods
- Krylov subspace methods
- Conjugate gradient method, CG (TNL::Solvers::Linear::CG)
- Biconjugate gradient stabilized method, BICGStab (TNL::Solvers::Linear::BICGStab)
- Biconjugate gradient stabilized method, BICGStab(l) (TNL::Solvers::Linear::BICGStabL)
- Transpose-free quasi-minimal residual method, TFQMR (TNL::Solvers::Linear::TFQMR)
- Generalized minimal residual method, GMRES (TNL::Solvers::Linear::GMRES) with various methods of orthogonalization:
- Classical Gramm-Schmidt, CGS
- Classical Gramm-Schmidt with reorthogonalization, CGSR
- Modified Gramm-Schmidt, MGS
- Modified Gramm-Schmidt with reorthogonalization, MGSR
- Compact WY form of the Householder reflections, CWY
The iterative solvers (not the stationary solvers like TNL::Solvers::Linear::Jacobi and TNL::Solvers::Linear::SOR) can be combined with the following preconditioners:
- Diagonal or Jacobi (TNL::Solvers::Linear::Preconditioners::Diagonal)
- ILU (Incomplete LU) - CPU only currently
Iterative solvers of linear systems
Basic setup
All iterative solvers for linear systems can be found in the namespace TNL::Solvers::Linear. The following example shows the use the iterative solvers:
1#include <iostream>
2#include <memory>
3#include <TNL/Matrices/SparseMatrix.h>
4#include <TNL/Devices/Host.h>
5#include <TNL/Devices/Cuda.h>
6#include <TNL/Solvers/Linear/TFQMR.h>
7
8template< typename Device >
9void
10iterativeLinearSolverExample()
11{
12 /***
13 * Set the following matrix (dots represent zero matrix elements):
14 *
15 * / 2.5 -1 . . . \
16 * | -1 2.5 -1 . . |
17 * | . -1 2.5 -1. . |
18 * | . . -1 2.5 -1 |
19 * \ . . . -1 2.5 /
20 */
23 const int size( 5 );
25 matrix_ptr->setDimensions( size, size );
26 matrix_ptr->setRowCapacities( Vector( { 2, 3, 3, 3, 2 } ) );
27
29 {
30 const int rowIdx = row.getRowIndex();
31 if( rowIdx == 0 ) {
32 row.setElement( 0, rowIdx, 2.5 ); // diagonal element
33 row.setElement( 1, rowIdx + 1, -1 ); // element above the diagonal
34 }
35 else if( rowIdx == size - 1 ) {
36 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
37 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
38 }
39 else {
40 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
41 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
42 row.setElement( 2, rowIdx + 1, -1.0 ); // element above the diagonal
43 }
44 };
45
46 /***
47 * Set the matrix elements.
48 */
49 matrix_ptr->forAllRows( f );
51
52 /***
53 * Set the right-hand side vector.
54 */
55 Vector x( size, 1.0 );
56 Vector b( size );
57 matrix_ptr->vectorProduct( x, b );
58 x = 0.0;
60
61 /***
62 * Solve the linear system.
63 */
65 LinearSolver solver;
66 solver.setMatrix( matrix_ptr );
67 solver.setConvergenceResidue( 1.0e-6 );
68 solver.solve( b, x );
70}
71
72int
73main( int argc, char* argv[] )
74{
76 iterativeLinearSolverExample< TNL::Devices::Sequential >();
77
78#ifdef __CUDACC__
80 iterativeLinearSolverExample< TNL::Devices::Cuda >();
81#endif
82}
Implementation of sparse matrix, i.e. matrix storing only non-zero elements.
Definition SparseMatrix.h:55
virtual void setMatrix(const MatrixPointer &matrix)
Set the matrix of the linear system.
Definition LinearSolver.h:120
Iterative solver of linear systems based on the Transpose-free quasi-minimal residual (TFQMR) method.
Definition TFQMR.h:21
T endl(T... args)
In this example we solve a linear system \( A \vec x = \vec b \) where
\[A = \left( \begin{array}{cccc} 2.5 & -1 & & & \\ -1 & 2.5 & -1 & & \\ & -1 & 2.5 & -1 & \\ & & -1 & 2.5 & -1 \\ & & & -1 & 2.5 \\ \end{array} \right) \]
The right-hand side vector \(\vec b \) is set to \(( 1.5, 0.5, 0.5, 0.5, 1.5 )^T \) so that the exact solution is \( \vec x = ( 1, 1, 1, 1, 1 )^T \). The elements of the matrix \( A \) are set using the method TNL::Matrices::SparseMatrix::forAllRows. In this example, we use the sparse matrix but any other matrix type can be used as well (see the namespace TNL::Matrices). Next we set the solution vector \( \vec x = ( 1, 1, 1, 1, 1 )^T \) and multiply it with matrix \( A \) to get the right-hand side vector \( \vec b \). Finally, we reset the vector \( \vec x \) to zero.
To solve the linear system in the example, we use the TFQMR solver. Other solvers can be used as well (see the namespace TNL::Solvers::Linear). The solver needs only one template parameter which is the matrix type. Next we create an instance of the solver and set the matrix of the linear system. Note that the matrix is passed to the solver as a std::shared_ptr. Then we set the stopping criterion for the iterative method in terms of the relative residue, i.e. \( \lVert \vec b - A \vec x \rVert / \lVert b \rVert \). The solver is executed by calling the TNL::Solvers::Linear::LinearSolver::solve method which accepts the right-hand side vector \( \vec b \) and the solution vector \( \vec x \).
The result looks as follows:
Solving linear system on host:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
Vector x = [ 1, 1, 1, 1, 1 ]
Solving linear system on CUDA device:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
Vector x = [ 1, 1, 1, 1, 1 ]
Setup with a solver monitor
Solution of large linear systems may take a lot of time. In such situations, it is useful to be able to monitor the convergence of the solver or the solver status in general. For this purpose, TNL provides a solver monitor which can show various metrics in real time, such as current number of iterations, current residue of the approximate solution, etc. The solver monitor in TNL runs in a separate thread and it refreshes the status of the solver with a configurable refresh rate (once per 500 ms by default). The use of the solver monitor is demonstrated in the following example.
1#include <iostream>
2#include <memory>
3#include <chrono>
4#include <thread>
5#include <TNL/Matrices/SparseMatrix.h>
6#include <TNL/Devices/Sequential.h>
7#include <TNL/Devices/Cuda.h>
8#include <TNL/Solvers/Linear/Jacobi.h>
9
10template< typename Device >
11void
12iterativeLinearSolverExample()
13{
14 /***
15 * Set the following matrix (dots represent zero matrix elements):
16 *
17 * / 2.5 -1 . . . \
18 * | -1 2.5 -1 . . |
19 * | . -1 2.5 -1. . |
20 * | . . -1 2.5 -1 |
21 * \ . . . -1 2.5 /
22 */
25 const int size( 5 );
27 matrix_ptr->setDimensions( size, size );
28 matrix_ptr->setRowCapacities( Vector( { 2, 3, 3, 3, 2 } ) );
29
31 {
32 const int rowIdx = row.getRowIndex();
33 if( rowIdx == 0 ) {
34 row.setElement( 0, rowIdx, 2.5 ); // diagonal element
35 row.setElement( 1, rowIdx + 1, -1 ); // element above the diagonal
36 }
37 else if( rowIdx == size - 1 ) {
38 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
39 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
40 }
41 else {
42 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
43 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
44 row.setElement( 2, rowIdx + 1, -1.0 ); // element above the diagonal
45 }
46 };
47
48 /***
49 * Set the matrix elements.
50 */
51 matrix_ptr->forAllRows( f );
53
54 /***
55 * Set the right-hand side vector.
56 */
57 Vector x( size, 1.0 );
58 Vector b( size );
59 matrix_ptr->vectorProduct( x, b );
60 x = 0.0;
62
63 /***
64 * Setup solver of the linear system.
65 */
67 LinearSolver solver;
68 solver.setMatrix( matrix_ptr );
69 solver.setOmega( 0.0005 );
70
71 /***
72 * Setup monitor of the iterative solver.
73 */
75 IterativeSolverMonitorType monitor;
76 TNL::Solvers::SolverMonitorThread monitorThread( monitor );
77 monitor.setRefreshRate( 10 ); // refresh rate in milliseconds
78 monitor.setVerbose( 1 );
79 monitor.setStage( "Jacobi stage:" );
80 solver.setSolverMonitor( monitor );
81 solver.setConvergenceResidue( 1.0e-6 );
82 solver.solve( b, x );
83 monitor.stopMainLoop();
85}
86
87int
88main( int argc, char* argv[] )
89{
91 iterativeLinearSolverExample< TNL::Devices::Sequential >();
92
93#ifdef __CUDACC__
95 iterativeLinearSolverExample< TNL::Devices::Cuda >();
96#endif
97}
Iterative solver of linear systems based on the Jacobi method.
Definition Jacobi.h:21
A RAII wrapper for launching the SolverMonitor's main loop in a separate thread.
Definition SolverMonitor.h:133
First, we set up the same linear system as in the previous example, we create an instance of the Jacobi solver and we pass the matrix of the linear system to the solver. Then, we set the relaxation parameter \( \omega \) of the Jacobi solver to 0.0005. The reason is to artificially slow down the convergence, because we want to see some iterations in this example. Next, we create an instance of the solver monitor and a special thread for the monitor (an instance of the TNL::Solvers::SolverMonitorThread class). We use the following methods to configure the solver monitor:
- TNL::Solvers::SolverMonitor::setRefreshRate sets the refresh rate of the monitor to 10 milliseconds.
- TNL::Solvers::IterativeSolverMonitor::setVerbose sets verbosity of the monitor to 1.
- TNL::Solvers::IterativeSolverMonitor::setStage sets a name of the solver stage. The monitor stages serve for distinguishing between different stages of more complex solvers (for example when the linear system solver is embedded into a time-dependent PDE solver).
Next, we call TNL::Solvers::IterativeSolver::setSolverMonitor to connect the solver with the monitor and we set the convergence criterion based on the relative residue. Finally, we start the solver using the TNL::Solvers::Linear::Jacobi::solve method and when the solver finishes, we stop the monitor using TNL::Solvers::SolverMonitor::stopMainLoop.
The result looks as follows:
Solving linear system on host:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
Jacobi stage: ITER: 23807 RES: 0.018652
Jacobi stage: ITER: 47923 RES: 0.00045917
Jacobi stage: ITER: 68059 RES: 2.0832e-05
Vector x = [ 0.999999, 0.999999, 0.999998, 0.999999, 0.999999 ]
Solving linear system on CUDA device:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
Jacobi stage: ITER: 1383 RES: 0.67008
Jacobi stage: ITER: 1871 RES: 0.59876
Jacobi stage: ITER: 2219 RES: 0.55604
Jacobi stage: ITER: 2551 RES: 0.52013
Jacobi stage: ITER: 2883 RES: 0.48804
Jacobi stage: ITER: 3217 RES: 0.45874
Jacobi stage: ITER: 3551 RES: 0.4324
Jacobi stage: ITER: 3885 RES: 0.40795
Jacobi stage: ITER: 4213 RES: 0.38592
Jacobi stage: ITER: 4551 RES: 0.36498
Jacobi stage: ITER: 4887 RES: 0.34546
Jacobi stage: ITER: 5183 RES: 0.32933
Jacobi stage: ITER: 5427 RES: 0.3167
Jacobi stage: ITER: 5763 RES: 0.30022
Jacobi stage: ITER: 6091 RES: 0.28506
Jacobi stage: ITER: 6411 RES: 0.27109
Jacobi stage: ITER: 6747 RES: 0.25722
Jacobi stage: ITER: 7023 RES: 0.24639
Jacobi stage: ITER: 7357 RES: 0.23385
Jacobi stage: ITER: 7647 RES: 0.22364
Jacobi stage: ITER: 7993 RES: 0.21191
Jacobi stage: ITER: 8341 RES: 0.20082
Jacobi stage: ITER: 8627 RES: 0.19221
Jacobi stage: ITER: 8947 RES: 0.18295
Jacobi stage: ITER: 9151 RES: 0.17729
Jacobi stage: ITER: 9267 RES: 0.17415
Jacobi stage: ITER: 9415 RES: 0.17023
Jacobi stage: ITER: 9759 RES: 0.16144
Jacobi stage: ITER: 10027 RES: 0.15492
Jacobi stage: ITER: 10291 RES: 0.14876
Jacobi stage: ITER: 10569 RES: 0.14249
Jacobi stage: ITER: 10805 RES: 0.13741
Jacobi stage: ITER: 11115 RES: 0.13105
Jacobi stage: ITER: 11435 RES: 0.12476
Jacobi stage: ITER: 11739 RES: 0.11907
Jacobi stage: ITER: 12051 RES: 0.11349
Jacobi stage: ITER: 12355 RES: 0.10832
Jacobi stage: ITER: 12671 RES: 0.10318
Jacobi stage: ITER: 13003 RES: 0.098051
Jacobi stage: ITER: 13341 RES: 0.093061
Jacobi stage: ITER: 13673 RES: 0.088433
Jacobi stage: ITER: 14011 RES: 0.083985
Jacobi stage: ITER: 14349 RES: 0.079711
Jacobi stage: ITER: 14691 RES: 0.075654
Jacobi stage: ITER: 15015 RES: 0.071981
Jacobi stage: ITER: 15335 RES: 0.068529
Jacobi stage: ITER: 15655 RES: 0.065242
Jacobi stage: ITER: 15971 RES: 0.06215
Jacobi stage: ITER: 16309 RES: 0.058988
Jacobi stage: ITER: 16651 RES: 0.055986
Jacobi stage: ITER: 16991 RES: 0.053137
Jacobi stage: ITER: 17311 RES: 0.050589
Jacobi stage: ITER: 17633 RES: 0.048133
Jacobi stage: ITER: 17973 RES: 0.045683
Jacobi stage: ITER: 18277 RES: 0.043599
Jacobi stage: ITER: 18611 RES: 0.041432
Jacobi stage: ITER: 18947 RES: 0.039348
Jacobi stage: ITER: 19285 RES: 0.037345
Jacobi stage: ITER: 19625 RES: 0.035445
Jacobi stage: ITER: 19947 RES: 0.033745
Jacobi stage: ITER: 20281 RES: 0.032048
Jacobi stage: ITER: 20589 RES: 0.030567
Jacobi stage: ITER: 20899 RES: 0.029154
Jacobi stage: ITER: 21237 RES: 0.027671
Jacobi stage: ITER: 21571 RES: 0.026295
Jacobi stage: ITER: 21905 RES: 0.024972
Jacobi stage: ITER: 22247 RES: 0.023702
Jacobi stage: ITER: 22575 RES: 0.022537
Jacobi stage: ITER: 22907 RES: 0.021417
Jacobi stage: ITER: 23215 RES: 0.020427
Jacobi stage: ITER: 23539 RES: 0.019435
Jacobi stage: ITER: 23881 RES: 0.018435
Jacobi stage: ITER: 24209 RES: 0.017529
Jacobi stage: ITER: 24539 RES: 0.016668
Jacobi stage: ITER: 24879 RES: 0.01582
Jacobi stage: ITER: 25219 RES: 0.015015
Jacobi stage: ITER: 25523 RES: 0.01433
Jacobi stage: ITER: 25839 RES: 0.013651
Jacobi stage: ITER: 26173 RES: 0.012964
Jacobi stage: ITER: 26511 RES: 0.012312
Jacobi stage: ITER: 26851 RES: 0.011686
Jacobi stage: ITER: 27187 RES: 0.011098
Jacobi stage: ITER: 27527 RES: 0.010533
Jacobi stage: ITER: 27855 RES: 0.010016
Jacobi stage: ITER: 28159 RES: 0.0095587
Jacobi stage: ITER: 28475 RES: 0.0091058
Jacobi stage: ITER: 28795 RES: 0.0086691
Jacobi stage: ITER: 29115 RES: 0.0082533
Jacobi stage: ITER: 29439 RES: 0.0078526
Jacobi stage: ITER: 29777 RES: 0.007453
Jacobi stage: ITER: 30079 RES: 0.0071174
Jacobi stage: ITER: 30419 RES: 0.0067552
Jacobi stage: ITER: 30735 RES: 0.0064352
Jacobi stage: ITER: 31017 RES: 0.0061605
Jacobi stage: ITER: 31329 RES: 0.0058722
Jacobi stage: ITER: 31659 RES: 0.0055837
Jacobi stage: ITER: 31993 RES: 0.0053028
Jacobi stage: ITER: 32329 RES: 0.0050361
Jacobi stage: ITER: 32667 RES: 0.0047828
Jacobi stage: ITER: 33003 RES: 0.0045422
Jacobi stage: ITER: 33339 RES: 0.0043137
Jacobi stage: ITER: 33675 RES: 0.0040967
Jacobi stage: ITER: 34011 RES: 0.0038907
Jacobi stage: ITER: 34347 RES: 0.0036949
Jacobi stage: ITER: 34683 RES: 0.0035091
Jacobi stage: ITER: 35019 RES: 0.0033326
Jacobi stage: ITER: 35355 RES: 0.0031649
Jacobi stage: ITER: 35691 RES: 0.0030057
Jacobi stage: ITER: 36025 RES: 0.0028546
Jacobi stage: ITER: 36363 RES: 0.002711
Jacobi stage: ITER: 36697 RES: 0.0025746
Jacobi stage: ITER: 37035 RES: 0.0024451
Jacobi stage: ITER: 37369 RES: 0.0023221
Jacobi stage: ITER: 37707 RES: 0.0022053
Jacobi stage: ITER: 38041 RES: 0.0020944
Jacobi stage: ITER: 38377 RES: 0.001989
Jacobi stage: ITER: 38719 RES: 0.0018878
Jacobi stage: ITER: 39055 RES: 0.0017929
Jacobi stage: ITER: 39389 RES: 0.0017027
Jacobi stage: ITER: 39723 RES: 0.001618
Jacobi stage: ITER: 40059 RES: 0.0015366
Jacobi stage: ITER: 40395 RES: 0.0014593
Jacobi stage: ITER: 40731 RES: 0.0013859
Jacobi stage: ITER: 41065 RES: 0.0013162
Jacobi stage: ITER: 41401 RES: 0.00125
Jacobi stage: ITER: 41735 RES: 0.0011879
Jacobi stage: ITER: 42071 RES: 0.0011281
Jacobi stage: ITER: 42407 RES: 0.0010714
Jacobi stage: ITER: 42741 RES: 0.0010175
Jacobi stage: ITER: 43077 RES: 0.0009663
Jacobi stage: ITER: 43411 RES: 0.00091826
Jacobi stage: ITER: 43747 RES: 0.00087207
Jacobi stage: ITER: 44083 RES: 0.0008282
Jacobi stage: ITER: 44419 RES: 0.00078654
Jacobi stage: ITER: 44753 RES: 0.00074698
Jacobi stage: ITER: 45089 RES: 0.0007094
Jacobi stage: ITER: 45423 RES: 0.00067413
Jacobi stage: ITER: 45759 RES: 0.00064023
Jacobi stage: ITER: 46095 RES: 0.00060802
Jacobi stage: ITER: 46429 RES: 0.00057744
Jacobi stage: ITER: 46765 RES: 0.00054839
Jacobi stage: ITER: 47005 RES: 0.00052854
Jacobi stage: ITER: 47495 RES: 0.00049037
Jacobi stage: ITER: 47831 RES: 0.00046571
Jacobi stage: ITER: 48167 RES: 0.00044228
Jacobi stage: ITER: 48503 RES: 0.00042003
Jacobi stage: ITER: 48837 RES: 0.00039891
Jacobi stage: ITER: 49173 RES: 0.00037884
Jacobi stage: ITER: 49507 RES: 0.00036001
Jacobi stage: ITER: 49843 RES: 0.0003419
Jacobi stage: ITER: 50179 RES: 0.0003247
Jacobi stage: ITER: 50515 RES: 0.00030837
Jacobi stage: ITER: 50849 RES: 0.00029286
Jacobi stage: ITER: 51185 RES: 0.00027812
Jacobi stage: ITER: 51519 RES: 0.0002643
Jacobi stage: ITER: 51855 RES: 0.000251
Jacobi stage: ITER: 52191 RES: 0.00023838
Jacobi stage: ITER: 52527 RES: 0.00022639
Jacobi stage: ITER: 52861 RES: 0.000215
Jacobi stage: ITER: 53197 RES: 0.00020418
Jacobi stage: ITER: 53531 RES: 0.00019403
Jacobi stage: ITER: 53867 RES: 0.00018427
Jacobi stage: ITER: 54203 RES: 0.000175
Jacobi stage: ITER: 54539 RES: 0.0001662
Jacobi stage: ITER: 54873 RES: 0.00015784
Jacobi stage: ITER: 55209 RES: 0.0001499
Jacobi stage: ITER: 55543 RES: 0.00014245
Jacobi stage: ITER: 55879 RES: 0.00013528
Jacobi stage: ITER: 56215 RES: 0.00012848
Jacobi stage: ITER: 56549 RES: 0.00012202
Jacobi stage: ITER: 56885 RES: 0.00011588
Jacobi stage: ITER: 57221 RES: 0.00011005
Jacobi stage: ITER: 57555 RES: 0.00010458
Jacobi stage: ITER: 57891 RES: 9.9318e-05
Jacobi stage: ITER: 58227 RES: 9.4322e-05
Jacobi stage: ITER: 58563 RES: 8.9578e-05
Jacobi stage: ITER: 58897 RES: 8.5072e-05
Jacobi stage: ITER: 59231 RES: 8.0842e-05
Jacobi stage: ITER: 59567 RES: 7.6776e-05
Jacobi stage: ITER: 59903 RES: 7.2914e-05
Jacobi stage: ITER: 60239 RES: 6.9246e-05
Jacobi stage: ITER: 60573 RES: 6.5763e-05
Jacobi stage: ITER: 60909 RES: 6.2455e-05
Jacobi stage: ITER: 61243 RES: 5.935e-05
Jacobi stage: ITER: 61579 RES: 5.6365e-05
Jacobi stage: ITER: 61915 RES: 5.3529e-05
Jacobi stage: ITER: 62251 RES: 5.0837e-05
Jacobi stage: ITER: 62587 RES: 4.828e-05
Jacobi stage: ITER: 62921 RES: 4.5851e-05
Jacobi stage: ITER: 63255 RES: 4.3572e-05
Jacobi stage: ITER: 63591 RES: 4.138e-05
Jacobi stage: ITER: 63927 RES: 3.9298e-05
Jacobi stage: ITER: 64263 RES: 3.7322e-05
Jacobi stage: ITER: 64597 RES: 3.5444e-05
Jacobi stage: ITER: 64933 RES: 3.3662e-05
Jacobi stage: ITER: 65269 RES: 3.1968e-05
Jacobi stage: ITER: 65603 RES: 3.0379e-05
Jacobi stage: ITER: 65939 RES: 2.8851e-05
Jacobi stage: ITER: 66275 RES: 2.74e-05
Jacobi stage: ITER: 66611 RES: 2.6021e-05
Jacobi stage: ITER: 66945 RES: 2.4712e-05
Jacobi stage: ITER: 67281 RES: 2.3469e-05
Jacobi stage: ITER: 67615 RES: 2.2303e-05
Jacobi stage: ITER: 67951 RES: 2.1181e-05
Jacobi stage: ITER: 68287 RES: 2.0115e-05
Jacobi stage: ITER: 68621 RES: 1.9104e-05
Jacobi stage: ITER: 68957 RES: 1.8143e-05
Jacobi stage: ITER: 69291 RES: 1.7241e-05
Jacobi stage: ITER: 69627 RES: 1.6373e-05
Jacobi stage: ITER: 69963 RES: 1.555e-05
Jacobi stage: ITER: 70299 RES: 1.4768e-05
Jacobi stage: ITER: 70633 RES: 1.4025e-05
Jacobi stage: ITER: 70969 RES: 1.3319e-05
Jacobi stage: ITER: 71303 RES: 1.2657e-05
Jacobi stage: ITER: 71639 RES: 1.202e-05
Jacobi stage: ITER: 71975 RES: 1.1416e-05
Jacobi stage: ITER: 72309 RES: 1.0842e-05
Jacobi stage: ITER: 72645 RES: 1.0296e-05
Jacobi stage: ITER: 72979 RES: 9.7843e-06
Jacobi stage: ITER: 73315 RES: 9.2922e-06
Jacobi stage: ITER: 73651 RES: 8.8248e-06
Jacobi stage: ITER: 73987 RES: 8.3809e-06
Jacobi stage: ITER: 74321 RES: 7.9593e-06
Jacobi stage: ITER: 74657 RES: 7.559e-06
Jacobi stage: ITER: 74991 RES: 7.1831e-06
Jacobi stage: ITER: 75327 RES: 6.8218e-06
Jacobi stage: ITER: 75663 RES: 6.4787e-06
Jacobi stage: ITER: 75999 RES: 6.1528e-06
Jacobi stage: ITER: 76333 RES: 5.8433e-06
Jacobi stage: ITER: 76669 RES: 5.5494e-06
Jacobi stage: ITER: 77003 RES: 5.2735e-06
Jacobi stage: ITER: 77339 RES: 5.0082e-06
Jacobi stage: ITER: 77675 RES: 4.7563e-06
Jacobi stage: ITER: 78009 RES: 4.5171e-06
Jacobi stage: ITER: 78343 RES: 4.2925e-06
Jacobi stage: ITER: 78679 RES: 4.0766e-06
Jacobi stage: ITER: 79015 RES: 3.8715e-06
Jacobi stage: ITER: 79351 RES: 3.6768e-06
Jacobi stage: ITER: 79685 RES: 3.4918e-06
Jacobi stage: ITER: 80019 RES: 3.3162e-06
Jacobi stage: ITER: 80355 RES: 3.1513e-06
Jacobi stage: ITER: 80691 RES: 2.9928e-06
Jacobi stage: ITER: 81027 RES: 2.8423e-06
Jacobi stage: ITER: 81363 RES: 2.6993e-06
Jacobi stage: ITER: 81697 RES: 2.5635e-06
Jacobi stage: ITER: 82031 RES: 2.4361e-06
Jacobi stage: ITER: 82367 RES: 2.3135e-06
Jacobi stage: ITER: 82703 RES: 2.1972e-06
Jacobi stage: ITER: 83039 RES: 2.0866e-06
Jacobi stage: ITER: 83373 RES: 1.9817e-06
Jacobi stage: ITER: 83709 RES: 1.882e-06
Jacobi stage: ITER: 84043 RES: 1.7884e-06
Jacobi stage: ITER: 84379 RES: 1.6985e-06
Jacobi stage: ITER: 84715 RES: 1.613e-06
Jacobi stage: ITER: 85049 RES: 1.5319e-06
Jacobi stage: ITER: 85385 RES: 1.4548e-06
Jacobi stage: ITER: 85719 RES: 1.3825e-06
Jacobi stage: ITER: 86055 RES: 1.313e-06
Jacobi stage: ITER: 86391 RES: 1.2469e-06
Jacobi stage: ITER: 86727 RES: 1.1842e-06
Jacobi stage: ITER: 87063 RES: 1.1246e-06
Jacobi stage: ITER: 87397 RES: 1.0681e-06
Jacobi stage: ITER: 87733 RES: 1.0143e-06
Vector x = [ 0.999999, 0.999999, 0.999998, 0.999999, 0.999999 ]
The monitoring of the solver can be improved by time elapsed since the beginning of the computation as demonstrated in the following example:
1#include <iostream>
2#include <memory>
3#include <chrono>
4#include <thread>
5#include <TNL/Timer.h>
6#include <TNL/Matrices/SparseMatrix.h>
7#include <TNL/Devices/Sequential.h>
8#include <TNL/Devices/Cuda.h>
9#include <TNL/Solvers/Linear/Jacobi.h>
10
11template< typename Device >
12void
13iterativeLinearSolverExample()
14{
15 /***
16 * Set the following matrix (dots represent zero matrix elements):
17 *
18 * / 2.5 -1 . . . \
19 * | -1 2.5 -1 . . |
20 * | . -1 2.5 -1. . |
21 * | . . -1 2.5 -1 |
22 * \ . . . -1 2.5 /
23 */
26 const int size( 5 );
28 matrix_ptr->setDimensions( size, size );
29 matrix_ptr->setRowCapacities( Vector( { 2, 3, 3, 3, 2 } ) );
30
32 {
33 const int rowIdx = row.getRowIndex();
34 if( rowIdx == 0 ) {
35 row.setElement( 0, rowIdx, 2.5 ); // diagonal element
36 row.setElement( 1, rowIdx + 1, -1 ); // element above the diagonal
37 }
38 else if( rowIdx == size - 1 ) {
39 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
40 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
41 }
42 else {
43 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
44 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
45 row.setElement( 2, rowIdx + 1, -1.0 ); // element above the diagonal
46 }
47 };
48
49 /***
50 * Set the matrix elements.
51 */
52 matrix_ptr->forAllRows( f );
54
55 /***
56 * Set the right-hand side vector.
57 */
58 Vector x( size, 1.0 );
59 Vector b( size );
60 matrix_ptr->vectorProduct( x, b );
61 x = 0.0;
63
64 /***
65 * Setup solver of the linear system.
66 */
68 LinearSolver solver;
69 solver.setMatrix( matrix_ptr );
70 solver.setOmega( 0.0005 );
71
72 /***
73 * Setup monitor of the iterative solver.
74 */
76 IterativeSolverMonitorType monitor;
77 TNL::Solvers::SolverMonitorThread mmonitorThread( monitor );
78 monitor.setRefreshRate( 10 ); // refresh rate in milliseconds
79 monitor.setVerbose( 1 );
80 monitor.setStage( "Jacobi stage:" );
81 TNL::Timer timer;
82 monitor.setTimer( timer );
83 timer.start();
84 solver.setSolverMonitor( monitor );
85 solver.setConvergenceResidue( 1.0e-6 );
86 solver.solve( b, x );
87 monitor.stopMainLoop();
89}
90
91int
92main( int argc, char* argv[] )
93{
95 iterativeLinearSolverExample< TNL::Devices::Sequential >();
96
97#ifdef __CUDACC__
99 iterativeLinearSolverExample< TNL::Devices::Cuda >();
100#endif
101}
The only changes are around the lines where we create an instance of TNL::Timer, connect it with the monitor using TNL::Solvers::SolverMonitor::setTimer and start the timer with TNL::Timer::start.
The result looks as follows:
Solving linear system on host:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
ELA: 0
ELA: 0.10009 Jacobi stage: ITER: 25707 RES: 0.013931
ELA: 0.20188 Jacobi stage: ITER: 51015 RES: 0.00028557
ELA: 0.30199 Jacobi stage: ITER: 78113 RES: 4.4455e-06
Vector x = [ 0.999999, 0.999999, 0.999998, 0.999999, 0.999999 ]
Solving linear system on CUDA device:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
ELA:1.7664e-
ELA: 0.10046 Jacobi stage: ITER: 499 RES: 0.85129
ELA: 0.20057 Jacobi stage: ITER: 955 RES: 0.74746
ELA: 0.30068 Jacobi stage: ITER: 1287 RES: 0.68606
ELA: 0.40079 Jacobi stage: ITER: 1619 RES: 0.63368
ELA: 0.50089 Jacobi stage: ITER: 1951 RES: 0.58841
ELA: 0.60099 Jacobi stage: ITER: 2287 RES: 0.54833
ELA: 0.70114 Jacobi stage: ITER: 2619 RES: 0.51327
ELA: 0.80124 Jacobi stage: ITER: 2955 RES: 0.48151
ELA: 0.90134 Jacobi stage: ITER: 3291 RES: 0.45282
ELA: 1.0014 Jacobi stage: ITER: 3589 RES: 0.42938
ELA: 1.1024 Jacobi stage: ITER: 3921 RES: 0.40545
ELA: 1.2034 Jacobi stage: ITER: 4289 RES: 0.38105
ELA: 1.3036 Jacobi stage: ITER: 4655 RES: 0.35879
ELA: 1.4037 Jacobi stage: ITER: 4993 RES: 0.33947
ELA: 1.5038 Jacobi stage: ITER: 5319 RES: 0.32222
ELA: 1.6039 Jacobi stage: ITER: 5631 RES: 0.30657
ELA: 1.704 Jacobi stage: ITER: 5975 RES: 0.29032
ELA: 1.8041 Jacobi stage: ITER: 6325 RES: 0.27469
ELA: 1.9042 Jacobi stage: ITER: 6653 RES: 0.26094
ELA: 2.0043 Jacobi stage: ITER: 6991 RES: 0.24762
ELA: 2.1044 Jacobi stage: ITER: 7327 RES: 0.23502
ELA: 2.2045 Jacobi stage: ITER: 7595 RES: 0.22545
ELA: 2.3046 Jacobi stage: ITER: 7697 RES: 0.22184
ELA: 2.4047 Jacobi stage: ITER: 7813 RES: 0.2179
ELA: 2.5048 Jacobi stage: ITER: 7897 RES: 0.21508
ELA: 2.6049 Jacobi stage: ITER: 8099 RES: 0.20853
ELA: 2.705 Jacobi stage: ITER: 8243 RES: 0.20394
ELA: 2.8094 Jacobi stage: ITER: 8409 RES: 0.19872
ELA: 2.9095 Jacobi stage: ITER: 8669 RES: 0.19091
ELA: 3.0096 Jacobi stage: ITER: 8793 RES: 0.18729
ELA: 3.1097 Jacobi stage: ITER: 8925 RES: 0.18352
ELA: 3.2098 Jacobi stage: ITER: 9121 RES: 0.17806
ELA: 3.3099 Jacobi stage: ITER: 9375 RES: 0.17128
ELA: 3.41 Jacobi stage: ITER: 9525 RES: 0.16731
ELA: 3.5101 Jacobi stage: ITER: 9747 RES: 0.16174
ELA: 3.6102 Jacobi stage: ITER: 10007 RES: 0.1554
ELA: 3.7104 Jacobi stage: ITER: 10205 RES: 0.15069
ELA: 3.8105 Jacobi stage: ITER: 10507 RES: 0.1439
ELA: 3.9136 Jacobi stage: ITER: 10837 RES: 0.13673
ELA: 4.0137 Jacobi stage: ITER: 11173 RES: 0.12985
ELA: 4.1138 Jacobi stage: ITER: 11505 RES: 0.12339
ELA: 4.2139 Jacobi stage: ITER: 11845 RES: 0.11711
ELA: 4.314 Jacobi stage: ITER: 12183 RES: 0.11122
ELA: 4.4154 Jacobi stage: ITER: 12515 RES: 0.10568
ELA: 4.5155 Jacobi stage: ITER: 12835 RES: 0.10061
ELA: 4.6156 Jacobi stage: ITER: 13157 RES: 0.095729
ELA: 4.7167 Jacobi stage: ITER: 13463 RES: 0.091361
ELA: 4.8168 Jacobi stage: ITER: 13791 RES: 0.086871
ELA: 4.9169 Jacobi stage: ITER: 14129 RES: 0.082451
ELA: 5.0173 Jacobi stage: ITER: 14469 RES: 0.078255
ELA: 5.1174 Jacobi stage: ITER: 14793 RES: 0.074455
ELA: 5.2175 Jacobi stage: ITER: 15127 RES: 0.070753
ELA: 5.3177 Jacobi stage: ITER: 15439 RES: 0.067442
ELA: 5.4179 Jacobi stage: ITER: 15777 RES: 0.064011
ELA: 5.518 Jacobi stage: ITER: 16075 RES: 0.061165
ELA: 5.6182 Jacobi stage: ITER: 16405 RES: 0.058124
ELA: 5.7183 Jacobi stage: ITER: 16741 RES: 0.055201
ELA: 5.8184 Jacobi stage: ITER: 17079 RES: 0.052424
ELA: 5.9184 Jacobi stage: ITER: 17407 RES: 0.049848
ELA: 6.0188 Jacobi stage: ITER: 17739 RES: 0.04737
ELA: 6.119 Jacobi stage: ITER: 18077 RES: 0.044959
ELA: 6.2191 Jacobi stage: ITER: 18381 RES: 0.042908
ELA: 6.3193 Jacobi stage: ITER: 18687 RES: 0.040951
ELA: 6.4195 Jacobi stage: ITER: 19019 RES: 0.038915
ELA: 6.5196 Jacobi stage: ITER: 19351 RES: 0.03698
ELA: 6.6197 Jacobi stage: ITER: 19687 RES: 0.03512
ELA: 6.7198 Jacobi stage: ITER: 20027 RES: 0.033333
ELA: 6.8199 Jacobi stage: ITER: 20349 RES: 0.031715
ELA: 6.92 Jacobi stage: ITER: 20671 RES: 0.030193
ELA: 7.0205 Jacobi stage: ITER: 20987 RES: 0.028763
ELA: 7.1205 Jacobi stage: ITER: 21307 RES: 0.027383
ELA: 7.2215 Jacobi stage: ITER: 21647 RES: 0.02599
ELA: 7.3216 Jacobi stage: ITER: 21971 RES: 0.024728
ELA: 7.4217 Jacobi stage: ITER: 22299 RES: 0.023513
ELA: 7.5218 Jacobi stage: ITER: 22637 RES: 0.022317
ELA: 7.6219 Jacobi stage: ITER: 22959 RES: 0.021246
ELA: 7.722 Jacobi stage: ITER: 23271 RES: 0.020252
ELA: 7.8221 Jacobi stage: ITER: 23587 RES: 0.019293
ELA: 7.9222 Jacobi stage: ITER: 23919 RES: 0.018333
ELA: 8.0224 Jacobi stage: ITER: 24257 RES: 0.017401
ELA: 8.1225 Jacobi stage: ITER: 24595 RES: 0.016525
ELA: 8.2226 Jacobi stage: ITER: 24935 RES: 0.015684
ELA: 8.323 Jacobi stage: ITER: 25271 RES: 0.014895
ELA: 8.4232 Jacobi stage: ITER: 25583 RES: 0.014198
ELA: 8.5234 Jacobi stage: ITER: 25883 RES: 0.013559
ELA: 8.6235 Jacobi stage: ITER: 26207 RES: 0.012901
ELA: 8.7244 Jacobi stage: ITER: 26515 RES: 0.012305
ELA: 8.8264 Jacobi stage: ITER: 26855 RES: 0.011679
ELA: 8.9265 Jacobi stage: ITER: 27171 RES: 0.011125
ELA: 9.0266 Jacobi stage: ITER: 27509 RES: 0.010559
ELA: 9.1267 Jacobi stage: ITER: 27847 RES: 0.010028
ELA: 9.2268 Jacobi stage: ITER: 28179 RES: 0.0095294
ELA: 9.3269 Jacobi stage: ITER: 28523 RES: 0.0090389
ELA: 9.427 Jacobi stage: ITER: 28851 RES: 0.0085948
ELA: 9.5271 Jacobi stage: ITER: 29155 RES: 0.0082027
ELA: 9.6272 Jacobi stage: ITER: 29491 RES: 0.0077901
ELA: 9.7273 Jacobi stage: ITER: 29827 RES: 0.0073983
ELA: 9.8274 Jacobi stage: ITER: 30161 RES: 0.0070261
ELA: 9.9284 Jacobi stage: ITER: 30499 RES: 0.0066727
ELA: 10.029 Jacobi stage: ITER: 30835 RES: 0.0063371
ELA: 10.129 Jacobi stage: ITER: 31171 RES: 0.0060183
ELA: 10.229 Jacobi stage: ITER: 31507 RES: 0.0057156
ELA: 10.33 Jacobi stage: ITER: 31843 RES: 0.0054281
ELA: 10.43 Jacobi stage: ITER: 32179 RES: 0.005155
ELA: 10.53 Jacobi stage: ITER: 32515 RES: 0.0048957
ELA: 10.63 Jacobi stage: ITER: 32849 RES: 0.0046495
ELA: 10.73 Jacobi stage: ITER: 33185 RES: 0.0044156
ELA: 10.831 Jacobi stage: ITER: 33521 RES: 0.0041935
ELA: 10.931 Jacobi stage: ITER: 33855 RES: 0.003985
ELA: 11.032 Jacobi stage: ITER: 34197 RES: 0.0037799
ELA: 11.133 Jacobi stage: ITER: 34531 RES: 0.003592
ELA: 11.233 Jacobi stage: ITER: 34869 RES: 0.0034092
ELA: 11.333 Jacobi stage: ITER: 35203 RES: 0.0032397
ELA: 11.433 Jacobi stage: ITER: 35541 RES: 0.0030749
ELA: 11.534 Jacobi stage: ITER: 35875 RES: 0.002922
ELA: 11.634 Jacobi stage: ITER: 36215 RES: 0.0027733
ELA: 11.735 Jacobi stage: ITER: 36555 RES: 0.0026322
ELA: 11.835 Jacobi stage: ITER: 36893 RES: 0.0024982
ELA: 11.936 Jacobi stage: ITER: 37229 RES: 0.0023726
ELA: 12.036 Jacobi stage: ITER: 37563 RES: 0.0022546
ELA: 12.136 Jacobi stage: ITER: 37899 RES: 0.0021412
ELA: 12.236 Jacobi stage: ITER: 38235 RES: 0.0020335
ELA: 12.336 Jacobi stage: ITER: 38571 RES: 0.0019312
ELA: 12.437 Jacobi stage: ITER: 38907 RES: 0.0018341
ELA: 12.537 Jacobi stage: ITER: 39243 RES: 0.0017418
ELA: 12.637 Jacobi stage: ITER: 39579 RES: 0.0016542
ELA: 12.738 Jacobi stage: ITER: 39915 RES: 0.001571
ELA: 12.838 Jacobi stage: ITER: 40251 RES: 0.001492
ELA: 12.938 Jacobi stage: ITER: 40585 RES: 0.0014169
ELA: 13.038 Jacobi stage: ITER: 40919 RES: 0.0013465
ELA: 13.138 Jacobi stage: ITER: 41255 RES: 0.0012788
ELA: 13.238 Jacobi stage: ITER: 41591 RES: 0.0012144
ELA: 13.338 Jacobi stage: ITER: 41927 RES: 0.0011533
ELA: 13.438 Jacobi stage: ITER: 42261 RES: 0.0010953
ELA: 13.538 Jacobi stage: ITER: 42597 RES: 0.0010402
ELA: 13.638 Jacobi stage: ITER: 42931 RES: 0.00098852
ELA: 13.739 Jacobi stage: ITER: 43267 RES: 0.00093879
ELA: 13.839 Jacobi stage: ITER: 43603 RES: 0.00089157
ELA: 13.939 Jacobi stage: ITER: 43937 RES: 0.00084672
ELA: 14.039 Jacobi stage: ITER: 44273 RES: 0.00080413
ELA: 14.143 Jacobi stage: ITER: 44629 RES: 0.00076134
ELA: 14.244 Jacobi stage: ITER: 45005 RES: 0.00071862
ELA: 14.356 Jacobi stage: ITER: 45383 RES: 0.00067829
ELA: 14.456 Jacobi stage: ITER: 45719 RES: 0.00064417
ELA: 14.556 Jacobi stage: ITER: 46053 RES: 0.00061177
ELA: 14.657 Jacobi stage: ITER: 46389 RES: 0.000581
ELA: 14.757 Jacobi stage: ITER: 46723 RES: 0.00055211
ELA: 14.857 Jacobi stage: ITER: 47059 RES: 0.00052434
ELA: 14.957 Jacobi stage: ITER: 47395 RES: 0.00049796
ELA: 15.057 Jacobi stage: ITER: 47731 RES: 0.00047292
ELA: 15.157 Jacobi stage: ITER: 48065 RES: 0.00044913
ELA: 15.257 Jacobi stage: ITER: 48401 RES: 0.00042654
ELA: 15.357 Jacobi stage: ITER: 48735 RES: 0.00040533
ELA: 15.457 Jacobi stage: ITER: 49071 RES: 0.00038494
ELA: 15.557 Jacobi stage: ITER: 49407 RES: 0.00036558
ELA: 15.658 Jacobi stage: ITER: 49743 RES: 0.00034719
ELA: 15.758 Jacobi stage: ITER: 50077 RES: 0.00032973
ELA: 15.858 Jacobi stage: ITER: 50413 RES: 0.00031314
ELA: 15.958 Jacobi stage: ITER: 50747 RES: 0.00029757
ELA: 16.058 Jacobi stage: ITER: 51083 RES: 0.0002826
ELA: 16.158 Jacobi stage: ITER: 51419 RES: 0.00026839
ELA: 16.258 Jacobi stage: ITER: 51753 RES: 0.00025489
ELA: 16.358 Jacobi stage: ITER: 52089 RES: 0.00024207
ELA: 16.458 Jacobi stage: ITER: 52425 RES: 0.00022989
ELA: 16.558 Jacobi stage: ITER: 52759 RES: 0.00021846
ELA: 16.659 Jacobi stage: ITER: 53095 RES: 0.00020747
ELA: 16.759 Jacobi stage: ITER: 53431 RES: 0.00019704
ELA: 16.859 Jacobi stage: ITER: 53765 RES: 0.00018713
ELA: 16.959 Jacobi stage: ITER: 54101 RES: 0.00017771
ELA: 17.059 Jacobi stage: ITER: 54435 RES: 0.00016888
ELA: 17.159 Jacobi stage: ITER: 54771 RES: 0.00016038
ELA: 17.259 Jacobi stage: ITER: 55107 RES: 0.00015232
ELA: 17.359 Jacobi stage: ITER: 55443 RES: 0.00014465
ELA: 17.459 Jacobi stage: ITER: 55779 RES: 0.00013738
ELA: 17.56 Jacobi stage: ITER: 56113 RES: 0.00013047
ELA: 17.66 Jacobi stage: ITER: 56447 RES: 0.00012398
ELA: 17.76 Jacobi stage: ITER: 56783 RES: 0.00011774
ELA: 17.86 Jacobi stage: ITER: 57119 RES: 0.00011182
ELA: 17.96 Jacobi stage: ITER: 57455 RES: 0.0001062
ELA: 18.06 Jacobi stage: ITER: 57789 RES: 0.00010086
ELA: 18.16 Jacobi stage: ITER: 58125 RES: 9.5782e-05
ELA: 18.26 Jacobi stage: ITER: 58459 RES: 9.102e-05
ELA: 18.36 Jacobi stage: ITER: 58795 RES: 8.6442e-05
ELA: 18.46 Jacobi stage: ITER: 59131 RES: 8.2094e-05
ELA: 18.561 Jacobi stage: ITER: 59467 RES: 7.7964e-05
ELA: 18.661 Jacobi stage: ITER: 59801 RES: 7.4042e-05
ELA: 18.761 Jacobi stage: ITER: 60137 RES: 7.0318e-05
ELA: 18.861 Jacobi stage: ITER: 60471 RES: 6.6822e-05
ELA: 18.961 Jacobi stage: ITER: 60807 RES: 6.3461e-05
ELA: 19.061 Jacobi stage: ITER: 61143 RES: 6.0269e-05
ELA: 19.161 Jacobi stage: ITER: 61479 RES: 5.7237e-05
ELA: 19.261 Jacobi stage: ITER: 61813 RES: 5.4358e-05
ELA: 19.361 Jacobi stage: ITER: 62149 RES: 5.1624e-05
ELA: 19.462 Jacobi stage: ITER: 62483 RES: 4.9057e-05
ELA: 19.562 Jacobi stage: ITER: 62819 RES: 4.659e-05
ELA: 19.662 Jacobi stage: ITER: 63155 RES: 4.4246e-05
ELA: 19.762 Jacobi stage: ITER: 63491 RES: 4.202e-05
ELA: 19.862 Jacobi stage: ITER: 63827 RES: 3.9907e-05
ELA: 19.962 Jacobi stage: ITER: 64161 RES: 3.7899e-05
ELA: 20.062 Jacobi stage: ITER: 64495 RES: 3.6015e-05
ELA: 20.162 Jacobi stage: ITER: 64831 RES: 3.4204e-05
ELA: 20.262 Jacobi stage: ITER: 65167 RES: 3.2483e-05
ELA: 20.362 Jacobi stage: ITER: 65503 RES: 3.0849e-05
ELA: 20.463 Jacobi stage: ITER: 65837 RES: 2.9297e-05
ELA: 20.563 Jacobi stage: ITER: 66173 RES: 2.7824e-05
ELA: 20.663 Jacobi stage: ITER: 66507 RES: 2.644e-05
ELA: 20.763 Jacobi stage: ITER: 66843 RES: 2.511e-05
ELA: 20.863 Jacobi stage: ITER: 67179 RES: 2.3847e-05
ELA: 20.963 Jacobi stage: ITER: 67515 RES: 2.2648e-05
ELA: 21.063 Jacobi stage: ITER: 67849 RES: 2.1509e-05
ELA: 21.163 Jacobi stage: ITER: 68185 RES: 2.0427e-05
ELA: 21.263 Jacobi stage: ITER: 68519 RES: 1.9411e-05
ELA: 21.363 Jacobi stage: ITER: 68855 RES: 1.8435e-05
ELA: 21.464 Jacobi stage: ITER: 69191 RES: 1.7507e-05
ELA: 21.564 Jacobi stage: ITER: 69525 RES: 1.6627e-05
ELA: 21.664 Jacobi stage: ITER: 69861 RES: 1.579e-05
ELA: 21.764 Jacobi stage: ITER: 70195 RES: 1.5005e-05
ELA: 21.864 Jacobi stage: ITER: 70531 RES: 1.4251e-05
ELA: 21.964 Jacobi stage: ITER: 70867 RES: 1.3534e-05
ELA: 22.064 Jacobi stage: ITER: 71203 RES: 1.2853e-05
ELA: 22.164 Jacobi stage: ITER: 71537 RES: 1.2207e-05
ELA: 22.264 Jacobi stage: ITER: 71873 RES: 1.1593e-05
ELA: 22.365 Jacobi stage: ITER: 72207 RES: 1.1016e-05
ELA: 22.465 Jacobi stage: ITER: 72543 RES: 1.0462e-05
ELA: 22.565 Jacobi stage: ITER: 72879 RES: 9.9358e-06
ELA: 22.665 Jacobi stage: ITER: 73213 RES: 9.436e-06
ELA: 22.765 Jacobi stage: ITER: 73549 RES: 8.9614e-06
ELA: 22.865 Jacobi stage: ITER: 73885 RES: 8.5106e-06
ELA: 22.965 Jacobi stage: ITER: 74219 RES: 8.0875e-06
ELA: 23.065 Jacobi stage: ITER: 74555 RES: 7.6807e-06
ELA: 23.165 Jacobi stage: ITER: 74891 RES: 7.2943e-06
ELA: 23.265 Jacobi stage: ITER: 75225 RES: 6.9274e-06
ELA: 23.366 Jacobi stage: ITER: 75561 RES: 6.579e-06
ELA: 23.466 Jacobi stage: ITER: 75895 RES: 6.2519e-06
ELA: 23.566 Jacobi stage: ITER: 76231 RES: 5.9374e-06
ELA: 23.666 Jacobi stage: ITER: 76567 RES: 5.6387e-06
ELA: 23.766 Jacobi stage: ITER: 76901 RES: 5.3551e-06
ELA: 23.866 Jacobi stage: ITER: 77235 RES: 5.0889e-06
ELA: 23.966 Jacobi stage: ITER: 77571 RES: 4.8329e-06
ELA: 24.066 Jacobi stage: ITER: 77907 RES: 4.5898e-06
ELA: 24.166 Jacobi stage: ITER: 78243 RES: 4.3589e-06
ELA: 24.267 Jacobi stage: ITER: 78577 RES: 4.1397e-06
ELA: 24.367 Jacobi stage: ITER: 78913 RES: 3.9314e-06
ELA: 24.467 Jacobi stage: ITER: 79247 RES: 3.736e-06
ELA: 24.567 Jacobi stage: ITER: 79583 RES: 3.5481e-06
ELA: 24.667 Jacobi stage: ITER: 79919 RES: 3.3696e-06
ELA: 24.767 Jacobi stage: ITER: 80255 RES: 3.2001e-06
ELA: 24.867 Jacobi stage: ITER: 80589 RES: 3.0391e-06
ELA: 24.967 Jacobi stage: ITER: 80925 RES: 2.8862e-06
ELA: 25.067 Jacobi stage: ITER: 81259 RES: 2.7428e-06
ELA: 25.167 Jacobi stage: ITER: 81595 RES: 2.6048e-06
ELA: 25.268 Jacobi stage: ITER: 81931 RES: 2.4738e-06
ELA: 25.368 Jacobi stage: ITER: 82265 RES: 2.3493e-06
ELA: 25.468 Jacobi stage: ITER: 82601 RES: 2.2312e-06
ELA: 25.568 Jacobi stage: ITER: 82935 RES: 2.1202e-06
ELA: 25.668 Jacobi stage: ITER: 83271 RES: 2.0136e-06
ELA: 25.768 Jacobi stage: ITER: 83607 RES: 1.9123e-06
ELA: 25.868 Jacobi stage: ITER: 83943 RES: 1.8161e-06
ELA: 25.968 Jacobi stage: ITER: 84279 RES: 1.7248e-06
ELA: 26.068 Jacobi stage: ITER: 84613 RES: 1.638e-06
ELA: 26.168 Jacobi stage: ITER: 84949 RES: 1.5556e-06
ELA: 26.269 Jacobi stage: ITER: 85283 RES: 1.4783e-06
ELA: 26.369 Jacobi stage: ITER: 85619 RES: 1.4039e-06
ELA: 26.469 Jacobi stage: ITER: 86075 RES: 1.3089e-06
ELA: 26.569 Jacobi stage: ITER: 86583 RES: 1.2107e-06
ELA: 26.669 Jacobi stage: ITER: 87087 RES: 1.1205e-06
ELA: 26.769 Jacobi stage: ITER: 87591 RES: 1.037e-06
Vector x = [ 0.999999, 0.999999, 0.999998, 0.999999, 0.999999 ]
Setup with preconditioner
Preconditioners of iterative solvers can significantly improve the performance of the solver. In the case of the linear systems, they are used mainly with the Krylov subspace methods. Preconditioners cannot be used with the starionary methods (TNL::Solvers::Linear::Jacobi and TNL::Solvers::Linear::SOR). The following example shows how to setup an iterative solver of linear systems with preconditioning.
1#include <iostream>
2#include <memory>
3#include <TNL/Matrices/SparseMatrix.h>
4#include <TNL/Devices/Host.h>
5#include <TNL/Devices/Cuda.h>
6#include <TNL/Solvers/Linear/TFQMR.h>
7#include <TNL/Solvers/Linear/Preconditioners/Diagonal.h>
8
9template< typename Device >
10void
11iterativeLinearSolverExample()
12{
13 /***
14 * Set the following matrix (dots represent zero matrix elements):
15 *
16 * / 2.5 -1 . . . \
17 * | -1 2.5 -1 . . |
18 * | . -1 2.5 -1. . |
19 * | . . -1 2.5 -1 |
20 * \ . . . -1 2.5 /
21 */
24 const int size( 5 );
26 matrix_ptr->setDimensions( size, size );
27 matrix_ptr->setRowCapacities( Vector( { 2, 3, 3, 3, 2 } ) );
28
30 {
31 const int rowIdx = row.getRowIndex();
32 if( rowIdx == 0 ) {
33 row.setElement( 0, rowIdx, 2.5 ); // diagonal element
34 row.setElement( 1, rowIdx + 1, -1 ); // element above the diagonal
35 }
36 else if( rowIdx == size - 1 ) {
37 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
38 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
39 }
40 else {
41 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
42 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
43 row.setElement( 2, rowIdx + 1, -1.0 ); // element above the diagonal
44 }
45 };
46
47 /***
48 * Set the matrix elements.
49 */
50 matrix_ptr->forAllRows( f );
52
53 /***
54 * Set the right-hand side vector.
55 */
56 Vector x( size, 1.0 );
57 Vector b( size );
58 matrix_ptr->vectorProduct( x, b );
59 x = 0.0;
61
62 /***
63 * Solve the linear system using diagonal (Jacobi) preconditioner.
64 */
68 preconditioner_ptr->update( matrix_ptr );
69 LinearSolver solver;
70 solver.setMatrix( matrix_ptr );
71 solver.setPreconditioner( preconditioner_ptr );
72 solver.setConvergenceResidue( 1.0e-6 );
73 solver.solve( b, x );
75}
76
77int
78main( int argc, char* argv[] )
79{
81 iterativeLinearSolverExample< TNL::Devices::Sequential >();
82
83#ifdef __CUDACC__
85 iterativeLinearSolverExample< TNL::Devices::Cuda >();
86#endif
87}
Diagonal (Jacobi) preconditioner for iterative solvers of linear systems.
Definition Diagonal.h:21
In this example, we solve the same problem as in all other examples in this section. The only differences concerning the preconditioner happen in the solver setup. Similarly to the matrix of the linear system, the preconditioner needs to be passed to the solver as a std::shared_ptr. When the preconditioner object is created, we have to initialize it using the update method, which has to be called everytime the matrix of the linear system changes. This is important, for example, when solving time-dependent PDEs, but it does not happen in this example. Finally, we need to connect the solver with the preconditioner using the setPreconditioner method.
The result looks as follows:
Solving linear system on host:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
Vector x = [ 1, 1, 1, 1, 1 ]
Solving linear system on CUDA device:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
Vector x = [ 1, 1, 1, 1, 1 ]
Choosing the solver and preconditioner type at runtime
When developing a numerical solver, one often has to search for a combination of various methods and algorithms that fit given requirements the best. To make this easier, TNL provides the functions TNL::Solvers::getLinearSolver and TNL::Solvers::getPreconditioner for selecting the linear solver and preconditioner at runtime. The following example shows how to use these functions:
1#include <iostream>
2#include <memory>
3#include <TNL/Matrices/SparseMatrix.h>
4#include <TNL/Devices/Host.h>
5#include <TNL/Devices/Cuda.h>
6#include <TNL/Solvers/LinearSolverTypeResolver.h>
7
8template< typename Device >
9void
10iterativeLinearSolverExample()
11{
12 /***
13 * Set the following matrix (dots represent zero matrix elements):
14 *
15 * / 2.5 -1 . . . \
16 * | -1 2.5 -1 . . |
17 * | . -1 2.5 -1. . |
18 * | . . -1 2.5 -1 |
19 * \ . . . -1 2.5 /
20 */
23 const int size( 5 );
25 matrix_ptr->setDimensions( size, size );
26 matrix_ptr->setRowCapacities( Vector( { 2, 3, 3, 3, 2 } ) );
27
29 {
30 const int rowIdx = row.getRowIndex();
31 if( rowIdx == 0 ) {
32 row.setElement( 0, rowIdx, 2.5 ); // diagonal element
33 row.setElement( 1, rowIdx + 1, -1 ); // element above the diagonal
34 }
35 else if( rowIdx == size - 1 ) {
36 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
37 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
38 }
39 else {
40 row.setElement( 0, rowIdx - 1, -1.0 ); // element below the diagonal
41 row.setElement( 1, rowIdx, 2.5 ); // diagonal element
42 row.setElement( 2, rowIdx + 1, -1.0 ); // element above the diagonal
43 }
44 };
45
46 /***
47 * Set the matrix elements.
48 */
49 matrix_ptr->forAllRows( f );
51
52 /***
53 * Set the right-hand side vector.
54 */
55 Vector x( size, 1.0 );
56 Vector b( size );
57 matrix_ptr->vectorProduct( x, b );
58 x = 0.0;
60
61 /***
62 * Solve the linear system using diagonal (Jacobi) preconditioner.
63 */
66 preconditioner_ptr->update( matrix_ptr );
67 solver_ptr->setMatrix( matrix_ptr );
68 solver_ptr->setPreconditioner( preconditioner_ptr );
69 solver_ptr->setConvergenceResidue( 1.0e-6 );
70 solver_ptr->solve( b, x );
72}
73
74int
75main( int argc, char* argv[] )
76{
78 iterativeLinearSolverExample< TNL::Devices::Sequential >();
79
80#ifdef __CUDACC__
82 iterativeLinearSolverExample< TNL::Devices::Cuda >();
83#endif
84}
std::shared_ptr< Linear::Preconditioners::Preconditioner< MatrixType > > getPreconditioner(const std::string &name)
Function returning shared pointer with linear preconditioner given by its name in a form of a string.
Definition LinearSolverTypeResolver.h:144
std::shared_ptr< Linear::LinearSolver< MatrixType > > getLinearSolver(const std::string &name)
Function returning shared pointer with linear solver given by its name in a form of a string.
Definition LinearSolverTypeResolver.h:87
We still stay with the same problem and the only changes are in the solver setup. We first use TNL::Solvers::getLinearSolver to get a shared pointer holding the solver and then TNL::Solvers::getPreconditioner to get a shared pointer holding the preconditioner. The rest of the code is the same as in the previous examples with the only difference that we work with the pointer solver_ptr instead of the direct instance solver of the solver type.
The result looks as follows:
Solving linear system on host:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
Vector x = [ 1, 1, 1, 1, 1 ]
Solving linear system on CUDA device:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
Vector x = [ 1, 1, 1, 1, 1 ]
Generated on for Template Numerical Library by
1.14.0