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: 24453 RES: 0.016884
Jacobi stage: ITER: 49249 RES: 0.00037444
Jacobi stage: ITER: 74247 RES: 8.0478e-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 ]
Jacobi stage: ITER: 2227 RES: 0.55512
Jacobi stage: ITER: 3139 RES: 0.46546
Jacobi stage: ITER: 3823 RES: 0.41245
Jacobi stage: ITER: 5209 RES: 0.32785
Jacobi stage: ITER: 5835 RES: 0.29682
Jacobi stage: ITER: 6091 RES: 0.28506
Jacobi stage: ITER: 6361 RES: 0.27314
Jacobi stage: ITER: 6627 RES: 0.26208
Jacobi stage: ITER: 6917 RES: 0.25041
Jacobi stage: ITER: 7243 RES: 0.2381
Jacobi stage: ITER: 7741 RES: 0.22034
Jacobi stage: ITER: 8463 RES: 0.19713
Jacobi stage: ITER: 9101 RES: 0.17861
Jacobi stage: ITER: 9503 RES: 0.16793
Jacobi stage: ITER: 10095 RES: 0.15331
Jacobi stage: ITER: 10443 RES: 0.14532
Jacobi stage: ITER: 11607 RES: 0.12151
Jacobi stage: ITER: 11971 RES: 0.1149
Jacobi stage: ITER: 12159 RES: 0.11163
Jacobi stage: ITER: 12903 RES: 0.099569
Jacobi stage: ITER: 13183 RES: 0.095377
Jacobi stage: ITER: 13407 RES: 0.09215
Jacobi stage: ITER: 13701 RES: 0.088054
Jacobi stage: ITER: 14337 RES: 0.079858
Jacobi stage: ITER: 14575 RES: 0.077015
Jacobi stage: ITER: 14735 RES: 0.075145
Jacobi stage: ITER: 14911 RES: 0.07314
Jacobi stage: ITER: 15083 RES: 0.071233
Jacobi stage: ITER: 15283 RES: 0.069078
Jacobi stage: ITER: 15449 RES: 0.067318
Jacobi stage: ITER: 15763 RES: 0.064168
Jacobi stage: ITER: 15987 RES: 0.061998
Jacobi stage: ITER: 16215 RES: 0.059864
Jacobi stage: ITER: 16431 RES: 0.05791
Jacobi stage: ITER: 16665 RES: 0.055849
Jacobi stage: ITER: 16955 RES: 0.053432
Jacobi stage: ITER: 17227 RES: 0.051246
Jacobi stage: ITER: 17553 RES: 0.048728
Jacobi stage: ITER: 17795 RES: 0.046964
Jacobi stage: ITER: 18121 RES: 0.044657
Jacobi stage: ITER: 18459 RES: 0.04241
Jacobi stage: ITER: 18725 RES: 0.0407
Jacobi stage: ITER: 19047 RES: 0.038748
Jacobi stage: ITER: 19361 RES: 0.036912
Jacobi stage: ITER: 19687 RES: 0.03512
Jacobi stage: ITER: 20023 RES: 0.033353
Jacobi stage: ITER: 20359 RES: 0.031676
Jacobi stage: ITER: 20683 RES: 0.030138
Jacobi stage: ITER: 21019 RES: 0.028622
Jacobi stage: ITER: 21355 RES: 0.027182
Jacobi stage: ITER: 21801 RES: 0.025375
Jacobi stage: ITER: 22127 RES: 0.024143
Jacobi stage: ITER: 22493 RES: 0.022816
Jacobi stage: ITER: 22827 RES: 0.021681
Jacobi stage: ITER: 23185 RES: 0.020515
Jacobi stage: ITER: 23515 RES: 0.019507
Jacobi stage: ITER: 23849 RES: 0.018526
Jacobi stage: ITER: 24185 RES: 0.017594
Jacobi stage: ITER: 24523 RES: 0.016709
Jacobi stage: ITER: 24837 RES: 0.015917
Jacobi stage: ITER: 25145 RES: 0.015182
Jacobi stage: ITER: 25483 RES: 0.014418
Jacobi stage: ITER: 25825 RES: 0.013676
Jacobi stage: ITER: 26169 RES: 0.012972
Jacobi stage: ITER: 26487 RES: 0.012358
Jacobi stage: ITER: 26823 RES: 0.011736
Jacobi stage: ITER: 27155 RES: 0.011153
Jacobi stage: ITER: 27493 RES: 0.010585
Jacobi stage: ITER: 27803 RES: 0.010096
Jacobi stage: ITER: 28127 RES: 0.0096058
Jacobi stage: ITER: 28463 RES: 0.0091226
Jacobi stage: ITER: 28751 RES: 0.0087279
Jacobi stage: ITER: 29011 RES: 0.0083862
Jacobi stage: ITER: 29335 RES: 0.007979
Jacobi stage: ITER: 29615 RES: 0.0076431
Jacobi stage: ITER: 29837 RES: 0.0073846
Jacobi stage: ITER: 30119 RES: 0.0070738
Jacobi stage: ITER: 30375 RES: 0.006801
Jacobi stage: ITER: 30637 RES: 0.0065307
Jacobi stage: ITER: 30967 RES: 0.0062099
Jacobi stage: ITER: 31297 RES: 0.0059011
Jacobi stage: ITER: 31635 RES: 0.0056043
Jacobi stage: ITER: 31971 RES: 0.0053224
Jacobi stage: ITER: 32303 RES: 0.0050578
Jacobi stage: ITER: 32635 RES: 0.0048063
Jacobi stage: ITER: 32931 RES: 0.0045927
Jacobi stage: ITER: 33219 RES: 0.004394
Jacobi stage: ITER: 33551 RES: 0.0041755
Jacobi stage: ITER: 33831 RES: 0.0039997
Jacobi stage: ITER: 34043 RES: 0.0038716
Jacobi stage: ITER: 34369 RES: 0.0036814
Jacobi stage: ITER: 34693 RES: 0.0035026
Jacobi stage: ITER: 35031 RES: 0.0033264
Jacobi stage: ITER: 35313 RES: 0.0031845
Jacobi stage: ITER: 35611 RES: 0.0030429
Jacobi stage: ITER: 35943 RES: 0.0028916
Jacobi stage: ITER: 36275 RES: 0.0027479
Jacobi stage: ITER: 36607 RES: 0.0026112
Jacobi stage: ITER: 36935 RES: 0.0024829
Jacobi stage: ITER: 37263 RES: 0.002361
Jacobi stage: ITER: 37593 RES: 0.0022436
Jacobi stage: ITER: 37823 RES: 0.0021664
Jacobi stage: ITER: 38035 RES: 0.002097
Jacobi stage: ITER: 38367 RES: 0.0019927
Jacobi stage: ITER: 38699 RES: 0.0018936
Jacobi stage: ITER: 39033 RES: 0.0017984
Jacobi stage: ITER: 39367 RES: 0.001709
Jacobi stage: ITER: 39699 RES: 0.001624
Jacobi stage: ITER: 40031 RES: 0.0015433
Jacobi stage: ITER: 40359 RES: 0.0014674
Jacobi stage: ITER: 40627 RES: 0.0014083
Jacobi stage: ITER: 40959 RES: 0.0013382
Jacobi stage: ITER: 41293 RES: 0.0012709
Jacobi stage: ITER: 41625 RES: 0.0012077
Jacobi stage: ITER: 41957 RES: 0.0011477
Jacobi stage: ITER: 42291 RES: 0.0010906
Jacobi stage: ITER: 42623 RES: 0.0010364
Jacobi stage: ITER: 42915 RES: 0.00099095
Jacobi stage: ITER: 43215 RES: 0.00094632
Jacobi stage: ITER: 43523 RES: 0.00090259
Jacobi stage: ITER: 43819 RES: 0.00086248
Jacobi stage: ITER: 44143 RES: 0.0008206
Jacobi stage: ITER: 44473 RES: 0.00077981
Jacobi stage: ITER: 44785 RES: 0.00074332
Jacobi stage: ITER: 45085 RES: 0.00070984
Jacobi stage: ITER: 45403 RES: 0.00067621
Jacobi stage: ITER: 45717 RES: 0.00064417
Jacobi stage: ITER: 46041 RES: 0.0006129
Jacobi stage: ITER: 46351 RES: 0.00058458
Jacobi stage: ITER: 46649 RES: 0.00055825
Jacobi stage: ITER: 46995 RES: 0.00052952
Jacobi stage: ITER: 47347 RES: 0.00050165
Jacobi stage: ITER: 47699 RES: 0.00047525
Jacobi stage: ITER: 48031 RES: 0.00045162
Jacobi stage: ITER: 48363 RES: 0.00042917
Jacobi stage: ITER: 48695 RES: 0.00040783
Jacobi stage: ITER: 49029 RES: 0.00038731
Jacobi stage: ITER: 49361 RES: 0.00036806
Jacobi stage: ITER: 49695 RES: 0.00034976
Jacobi stage: ITER: 50031 RES: 0.00033217
Jacobi stage: ITER: 50365 RES: 0.00031546
Jacobi stage: ITER: 50697 RES: 0.00029977
Jacobi stage: ITER: 51029 RES: 0.00028487
Jacobi stage: ITER: 51363 RES: 0.00027071
Jacobi stage: ITER: 51695 RES: 0.00025725
Jacobi stage: ITER: 52027 RES: 0.00024446
Jacobi stage: ITER: 52359 RES: 0.0002323
Jacobi stage: ITER: 52703 RES: 0.00022035
Jacobi stage: ITER: 53047 RES: 0.00020901
Jacobi stage: ITER: 53389 RES: 0.00019825
Jacobi stage: ITER: 53727 RES: 0.00018828
Jacobi stage: ITER: 54063 RES: 0.00017881
Jacobi stage: ITER: 54395 RES: 0.00016992
Jacobi stage: ITER: 54729 RES: 0.00016137
Jacobi stage: ITER: 55061 RES: 0.00015335
Jacobi stage: ITER: 55393 RES: 0.00014572
Jacobi stage: ITER: 55727 RES: 0.00013848
Jacobi stage: ITER: 56059 RES: 0.00013159
Jacobi stage: ITER: 56391 RES: 0.00012505
Jacobi stage: ITER: 56727 RES: 0.00011876
Jacobi stage: ITER: 57059 RES: 0.00011286
Jacobi stage: ITER: 57391 RES: 0.00010725
Jacobi stage: ITER: 57723 RES: 0.00010191
Jacobi stage: ITER: 58055 RES: 9.6847e-05
Jacobi stage: ITER: 58387 RES: 9.2032e-05
Jacobi stage: ITER: 58721 RES: 8.7403e-05
Jacobi stage: ITER: 59053 RES: 8.3057e-05
Jacobi stage: ITER: 59387 RES: 7.8928e-05
Jacobi stage: ITER: 59719 RES: 7.5004e-05
Jacobi stage: ITER: 60051 RES: 7.1275e-05
Jacobi stage: ITER: 60383 RES: 6.7731e-05
Jacobi stage: ITER: 60717 RES: 6.4324e-05
Jacobi stage: ITER: 61049 RES: 6.1126e-05
Jacobi stage: ITER: 61383 RES: 5.8087e-05
Jacobi stage: ITER: 61715 RES: 5.5199e-05
Jacobi stage: ITER: 62047 RES: 5.2455e-05
Jacobi stage: ITER: 62379 RES: 4.9847e-05
Jacobi stage: ITER: 62711 RES: 4.7369e-05
Jacobi stage: ITER: 63045 RES: 4.4986e-05
Jacobi stage: ITER: 63377 RES: 4.275e-05
Jacobi stage: ITER: 63711 RES: 4.0624e-05
Jacobi stage: ITER: 64043 RES: 3.8604e-05
Jacobi stage: ITER: 64375 RES: 3.6685e-05
Jacobi stage: ITER: 64707 RES: 3.4861e-05
Jacobi stage: ITER: 65039 RES: 3.3128e-05
Jacobi stage: ITER: 65371 RES: 3.1481e-05
Jacobi stage: ITER: 65705 RES: 2.9897e-05
Jacobi stage: ITER: 66037 RES: 2.8411e-05
Jacobi stage: ITER: 66369 RES: 2.6999e-05
Jacobi stage: ITER: 66703 RES: 2.5656e-05
Jacobi stage: ITER: 67035 RES: 2.4381e-05
Jacobi stage: ITER: 67367 RES: 2.3169e-05
Jacobi stage: ITER: 67699 RES: 2.2017e-05
Jacobi stage: ITER: 68031 RES: 2.0922e-05
Jacobi stage: ITER: 68365 RES: 1.987e-05
Jacobi stage: ITER: 68697 RES: 1.8882e-05
Jacobi stage: ITER: 69031 RES: 1.7943e-05
Jacobi stage: ITER: 69363 RES: 1.7051e-05
Jacobi stage: ITER: 69695 RES: 1.6203e-05
Jacobi stage: ITER: 70027 RES: 1.5398e-05
Jacobi stage: ITER: 70359 RES: 1.4632e-05
Jacobi stage: ITER: 70693 RES: 1.3896e-05
Jacobi stage: ITER: 71027 RES: 1.3205e-05
Jacobi stage: ITER: 71359 RES: 1.2549e-05
Jacobi stage: ITER: 71691 RES: 1.1925e-05
Jacobi stage: ITER: 72023 RES: 1.1332e-05
Jacobi stage: ITER: 72355 RES: 1.0769e-05
Jacobi stage: ITER: 72687 RES: 1.0233e-05
Jacobi stage: ITER: 73021 RES: 9.7184e-06
Jacobi stage: ITER: 73353 RES: 9.2353e-06
Jacobi stage: ITER: 73687 RES: 8.7761e-06
Jacobi stage: ITER: 74019 RES: 8.3398e-06
Jacobi stage: ITER: 74351 RES: 7.9252e-06
Jacobi stage: ITER: 74683 RES: 7.5311e-06
Jacobi stage: ITER: 75017 RES: 7.1523e-06
Jacobi stage: ITER: 75349 RES: 6.7967e-06
Jacobi stage: ITER: 75681 RES: 6.4588e-06
Jacobi stage: ITER: 76015 RES: 6.1377e-06
Jacobi stage: ITER: 76347 RES: 5.8325e-06
Jacobi stage: ITER: 76679 RES: 5.5426e-06
Jacobi stage: ITER: 77011 RES: 5.267e-06
Jacobi stage: ITER: 77343 RES: 5.0051e-06
Jacobi stage: ITER: 77675 RES: 4.7563e-06
Jacobi stage: ITER: 78009 RES: 4.5171e-06
Jacobi stage: ITER: 78341 RES: 4.2925e-06
Jacobi stage: ITER: 78675 RES: 4.0791e-06
Jacobi stage: ITER: 79007 RES: 3.8763e-06
Jacobi stage: ITER: 79339 RES: 3.6836e-06
Jacobi stage: ITER: 79671 RES: 3.5004e-06
Jacobi stage: ITER: 80003 RES: 3.3264e-06
Jacobi stage: ITER: 80335 RES: 3.161e-06
Jacobi stage: ITER: 80669 RES: 3.002e-06
Jacobi stage: ITER: 81055 RES: 2.8301e-06
Jacobi stage: ITER: 81577 RES: 2.6112e-06
Jacobi stage: ITER: 82077 RES: 2.4182e-06
Jacobi stage: ITER: 82577 RES: 2.2394e-06
Jacobi stage: ITER: 83075 RES: 2.0751e-06
Jacobi stage: ITER: 83573 RES: 1.9217e-06
Jacobi stage: ITER: 84073 RES: 1.7797e-06
Jacobi stage: ITER: 84571 RES: 1.6491e-06
Jacobi stage: ITER: 85071 RES: 1.5272e-06
Jacobi stage: ITER: 85571 RES: 1.4143e-06
Jacobi stage: ITER: 86071 RES: 1.3097e-06
Jacobi stage: ITER: 86571 RES: 1.2129e-06
Jacobi stage: ITER: 87071 RES: 1.1233e-06
Jacobi stage: ITER: 87569 RES: 1.0402e-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:6.778e-0
ELA: 0.10009 Jacobi stage: ITER: 25677 RES: 0.013991
ELA: 0.2002 Jacobi stage: ITER: 52169 RES: 0.00023911
ELA: 0.3003 Jacobi stage: ITER: 79019 RES: 3.8691e-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.702e-0
ELA: 0.10012 Jacobi stage: ITER: 2211 RES: 0.5565
ELA: 0.20024 Jacobi stage: ITER: 4383 RES: 0.37526
ELA: 0.30034 Jacobi stage: ITER: 6127 RES: 0.28345
ELA: 0.40045 Jacobi stage: ITER: 7395 RES: 0.23255
ELA: 0.50055 Jacobi stage: ITER: 8325 RES: 0.20131
ELA: 0.60065 Jacobi stage: ITER: 9315 RES: 0.17287
ELA: 0.70075 Jacobi stage: ITER: 10161 RES: 0.15171
ELA: 0.80104 Jacobi stage: ITER: 10535 RES: 0.14328
ELA: 0.90114 Jacobi stage: ITER: 11035 RES: 0.13268
ELA: 1.0012 Jacobi stage: ITER: 11461 RES: 0.12423
ELA: 1.102 Jacobi stage: ITER: 11991 RES: 0.11455
ELA: 1.2021 Jacobi stage: ITER: 12327 RES: 0.10878
ELA: 1.303 Jacobi stage: ITER: 12643 RES: 0.10363
ELA: 1.4031 Jacobi stage: ITER: 12915 RES: 0.099386
ELA: 1.5032 Jacobi stage: ITER: 13437 RES: 0.091698
ELA: 1.6033 Jacobi stage: ITER: 14919 RES: 0.073051
ELA: 1.7034 Jacobi stage: ITER: 15139 RES: 0.070623
ELA: 1.8038 Jacobi stage: ITER: 15587 RES: 0.065927
ELA: 1.906 Jacobi stage: ITER: 16203 RES: 0.059975
ELA: 2.0062 Jacobi stage: ITER: 16379 RES: 0.058375
ELA: 2.1063 Jacobi stage: ITER: 16555 RES: 0.056818
ELA: 2.211 Jacobi stage: ITER: 16985 RES: 0.05317
ELA: 2.3111 Jacobi stage: ITER: 17167 RES: 0.05172
ELA: 2.4112 Jacobi stage: ITER: 17335 RES: 0.050402
ELA: 2.5113 Jacobi stage: ITER: 18003 RES: 0.045487
ELA: 2.6114 Jacobi stage: ITER: 18341 RES: 0.043173
ELA: 2.7115 Jacobi stage: ITER: 18775 RES: 0.040401
ELA: 2.8116 Jacobi stage: ITER: 19143 RES: 0.038181
ELA: 2.9117 Jacobi stage: ITER: 19309 RES: 0.037208
ELA: 3.0118 Jacobi stage: ITER: 19611 RES: 0.035532
ELA: 3.1119 Jacobi stage: ITER: 20163 RES: 0.032644
ELA: 3.212 Jacobi stage: ITER: 20511 RES: 0.030945
ELA: 3.3121 Jacobi stage: ITER: 20873 RES: 0.029262
ELA: 3.4122 Jacobi stage: ITER: 21195 RES: 0.027858
ELA: 3.5123 Jacobi stage: ITER: 21875 RES: 0.025095
ELA: 3.6125 Jacobi stage: ITER: 22259 RES: 0.023658
ELA: 3.7126 Jacobi stage: ITER: 22639 RES: 0.022317
ELA: 3.8127 Jacobi stage: ITER: 22991 RES: 0.021142
ELA: 3.9128 Jacobi stage: ITER: 23359 RES: 0.01998
ELA: 4.0129 Jacobi stage: ITER: 23683 RES: 0.01901
ELA: 4.113 Jacobi stage: ITER: 24017 RES: 0.018054
ELA: 4.213 Jacobi stage: ITER: 24375 RES: 0.017093
ELA: 4.3132 Jacobi stage: ITER: 24701 RES: 0.016253
ELA: 4.4133 Jacobi stage: ITER: 25011 RES: 0.015502
ELA: 4.5134 Jacobi stage: ITER: 25335 RES: 0.01475
ELA: 4.6135 Jacobi stage: ITER: 25671 RES: 0.014008
ELA: 4.7136 Jacobi stage: ITER: 26005 RES: 0.013303
ELA: 4.8137 Jacobi stage: ITER: 26337 RES: 0.012642
ELA: 4.9138 Jacobi stage: ITER: 26669 RES: 0.012013
ELA: 5.0139 Jacobi stage: ITER: 27007 RES: 0.011409
ELA: 5.114 Jacobi stage: ITER: 27375 RES: 0.010782
ELA: 5.2141 Jacobi stage: ITER: 27695 RES: 0.010265
ELA: 5.3142 Jacobi stage: ITER: 28057 RES: 0.0097067
ELA: 5.4143 Jacobi stage: ITER: 28391 RES: 0.0092241
ELA: 5.516 Jacobi stage: ITER: 28743 RES: 0.0087386
ELA: 5.617 Jacobi stage: ITER: 29075 RES: 0.0083041
ELA: 5.7172 Jacobi stage: ITER: 29411 RES: 0.0078864
ELA: 5.8173 Jacobi stage: ITER: 29747 RES: 0.0074897
ELA: 5.918 Jacobi stage: ITER: 30085 RES: 0.0071086
ELA: 6.022 Jacobi stage: ITER: 30413 RES: 0.0067594
ELA: 6.1221 Jacobi stage: ITER: 30739 RES: 0.0064312
ELA: 6.2222 Jacobi stage: ITER: 31075 RES: 0.0061077
ELA: 6.3241 Jacobi stage: ITER: 31357 RES: 0.005847
ELA: 6.4242 Jacobi stage: ITER: 31549 RES: 0.0056771
ELA: 6.5243 Jacobi stage: ITER: 31793 RES: 0.0054682
ELA: 6.6244 Jacobi stage: ITER: 32131 RES: 0.0051932
ELA: 6.7245 Jacobi stage: ITER: 32463 RES: 0.004935
ELA: 6.8246 Jacobi stage: ITER: 32801 RES: 0.0046839
ELA: 6.9247 Jacobi stage: ITER: 33111 RES: 0.0044675
ELA: 7.0248 Jacobi stage: ITER: 33435 RES: 0.0042506
ELA: 7.1249 Jacobi stage: ITER: 33767 RES: 0.0040392
ELA: 7.228 Jacobi stage: ITER: 34107 RES: 0.0038337
ELA: 7.3281 Jacobi stage: ITER: 34459 RES: 0.0036319
ELA: 7.4283 Jacobi stage: ITER: 34789 RES: 0.0034514
ELA: 7.5284 Jacobi stage: ITER: 35159 RES: 0.0032617
ELA: 7.629 Jacobi stage: ITER: 35557 RES: 0.0030673
ELA: 7.7292 Jacobi stage: ITER: 35883 RES: 0.0029184
ELA: 7.8293 Jacobi stage: ITER: 36221 RES: 0.0027699
ELA: 7.9294 Jacobi stage: ITER: 36591 RES: 0.0026177
ELA: 8.0295 Jacobi stage: ITER: 36923 RES: 0.0024875
ELA: 8.1296 Jacobi stage: ITER: 37247 RES: 0.0023668
ELA: 8.2297 Jacobi stage: ITER: 37585 RES: 0.0022463
ELA: 8.3298 Jacobi stage: ITER: 37921 RES: 0.0021333
ELA: 8.4299 Jacobi stage: ITER: 38255 RES: 0.0020273
ELA: 8.53 Jacobi stage: ITER: 38587 RES: 0.0019265
ELA: 8.6301 Jacobi stage: ITER: 38875 RES: 0.0018431
ELA: 8.7302 Jacobi stage: ITER: 39195 RES: 0.0017547
ELA: 8.8303 Jacobi stage: ITER: 39521 RES: 0.0016685
ELA: 8.9304 Jacobi stage: ITER: 39889 RES: 0.0015768
ELA: 9.0305 Jacobi stage: ITER: 40267 RES: 0.0014883
ELA: 9.1306 Jacobi stage: ITER: 40599 RES: 0.0014143
ELA: 9.2307 Jacobi stage: ITER: 40927 RES: 0.0013448
ELA: 9.3308 Jacobi stage: ITER: 41265 RES: 0.0012764
ELA: 9.4309 Jacobi stage: ITER: 41599 RES: 0.0012129
ELA: 9.531 Jacobi stage: ITER: 41899 RES: 0.0011583
ELA: 9.6311 Jacobi stage: ITER: 42239 RES: 0.0010994
ELA: 9.7312 Jacobi stage: ITER: 42571 RES: 0.0010447
ELA: 9.8313 Jacobi stage: ITER: 42903 RES: 0.00099278
ELA: 9.9314 Jacobi stage: ITER: 43239 RES: 0.00094284
ELA: 10.031 Jacobi stage: ITER: 43573 RES: 0.00089541
ELA: 10.132 Jacobi stage: ITER: 43907 RES: 0.0008509
ELA: 10.232 Jacobi stage: ITER: 44303 RES: 0.00080068
ELA: 10.332 Jacobi stage: ITER: 44683 RES: 0.00075529
ELA: 10.432 Jacobi stage: ITER: 45015 RES: 0.00071773
ELA: 10.532 Jacobi stage: ITER: 45349 RES: 0.00068163
ELA: 10.632 Jacobi stage: ITER: 45681 RES: 0.00064774
ELA: 10.732 Jacobi stage: ITER: 46015 RES: 0.00061554
ELA: 10.832 Jacobi stage: ITER: 46347 RES: 0.00058494
ELA: 10.932 Jacobi stage: ITER: 46679 RES: 0.00055585
ELA: 11.032 Jacobi stage: ITER: 47019 RES: 0.00052757
ELA: 11.133 Jacobi stage: ITER: 47379 RES: 0.00049919
ELA: 11.233 Jacobi stage: ITER: 47711 RES: 0.00047437
ELA: 11.333 Jacobi stage: ITER: 48043 RES: 0.00045079
ELA: 11.433 Jacobi stage: ITER: 48377 RES: 0.00042811
ELA: 11.533 Jacobi stage: ITER: 48709 RES: 0.00040683
ELA: 11.634 Jacobi stage: ITER: 49041 RES: 0.0003866
ELA: 11.734 Jacobi stage: ITER: 49375 RES: 0.00036738
ELA: 11.834 Jacobi stage: ITER: 49725 RES: 0.00034804
ELA: 11.934 Jacobi stage: ITER: 50083 RES: 0.00032952
ELA: 12.034 Jacobi stage: ITER: 50423 RES: 0.00031276
ELA: 12.134 Jacobi stage: ITER: 50773 RES: 0.00029629
ELA: 12.234 Jacobi stage: ITER: 51171 RES: 0.00027881
ELA: 12.334 Jacobi stage: ITER: 51507 RES: 0.00026478
ELA: 12.434 Jacobi stage: ITER: 51855 RES: 0.000251
ELA: 12.534 Jacobi stage: ITER: 52193 RES: 0.00023823
ELA: 12.635 Jacobi stage: ITER: 52521 RES: 0.00022653
ELA: 12.735 Jacobi stage: ITER: 52859 RES: 0.00021513
ELA: 12.835 Jacobi stage: ITER: 53199 RES: 0.00020418
ELA: 12.935 Jacobi stage: ITER: 53543 RES: 0.00019368
ELA: 13.035 Jacobi stage: ITER: 53897 RES: 0.00018337
ELA: 13.136 Jacobi stage: ITER: 54169 RES: 0.00017587
ELA: 13.236 Jacobi stage: ITER: 54415 RES: 0.0001694
ELA: 13.336 Jacobi stage: ITER: 54691 RES: 0.00016237
ELA: 13.436 Jacobi stage: ITER: 55023 RES: 0.00015429
ELA: 13.536 Jacobi stage: ITER: 55355 RES: 0.00014662
ELA: 13.636 Jacobi stage: ITER: 55687 RES: 0.00013933
ELA: 13.736 Jacobi stage: ITER: 56021 RES: 0.00013232
ELA: 13.836 Jacobi stage: ITER: 56353 RES: 0.00012575
ELA: 13.936 Jacobi stage: ITER: 56685 RES: 0.00011949
ELA: 14.036 Jacobi stage: ITER: 57023 RES: 0.00011348
ELA: 14.137 Jacobi stage: ITER: 57355 RES: 0.00010784
ELA: 14.237 Jacobi stage: ITER: 57687 RES: 0.00010248
ELA: 14.337 Jacobi stage: ITER: 58019 RES: 9.7384e-05
ELA: 14.437 Jacobi stage: ITER: 58351 RES: 9.2543e-05
ELA: 14.537 Jacobi stage: ITER: 58685 RES: 8.7888e-05
ELA: 14.637 Jacobi stage: ITER: 59017 RES: 8.3518e-05
ELA: 14.737 Jacobi stage: ITER: 59349 RES: 7.9366e-05
ELA: 14.837 Jacobi stage: ITER: 59695 RES: 7.5281e-05
ELA: 14.937 Jacobi stage: ITER: 60035 RES: 7.145e-05
ELA: 15.037 Jacobi stage: ITER: 60383 RES: 6.7731e-05
ELA: 15.138 Jacobi stage: ITER: 60719 RES: 6.4324e-05
ELA: 15.238 Jacobi stage: ITER: 61047 RES: 6.1164e-05
ELA: 15.338 Jacobi stage: ITER: 61379 RES: 5.8123e-05
ELA: 15.438 Jacobi stage: ITER: 61711 RES: 5.5233e-05
ELA: 15.538 Jacobi stage: ITER: 62045 RES: 5.2455e-05
ELA: 15.638 Jacobi stage: ITER: 62377 RES: 4.9847e-05
ELA: 15.738 Jacobi stage: ITER: 62709 RES: 4.7369e-05
ELA: 15.839 Jacobi stage: ITER: 63043 RES: 4.5014e-05
ELA: 15.939 Jacobi stage: ITER: 63375 RES: 4.2776e-05
ELA: 16.039 Jacobi stage: ITER: 63701 RES: 4.0674e-05
ELA: 16.139 Jacobi stage: ITER: 64033 RES: 3.8652e-05
ELA: 16.239 Jacobi stage: ITER: 64365 RES: 3.673e-05
ELA: 16.339 Jacobi stage: ITER: 64699 RES: 3.4904e-05
ELA: 16.439 Jacobi stage: ITER: 65031 RES: 3.3169e-05
ELA: 16.539 Jacobi stage: ITER: 65363 RES: 3.152e-05
ELA: 16.639 Jacobi stage: ITER: 65695 RES: 2.9953e-05
ELA: 16.74 Jacobi stage: ITER: 66027 RES: 2.8463e-05
ELA: 16.84 Jacobi stage: ITER: 66359 RES: 2.7048e-05
ELA: 16.94 Jacobi stage: ITER: 66691 RES: 2.5704e-05
ELA: 17.04 Jacobi stage: ITER: 67025 RES: 2.4411e-05
ELA: 17.14 Jacobi stage: ITER: 67357 RES: 2.3197e-05
ELA: 17.24 Jacobi stage: ITER: 67689 RES: 2.2044e-05
ELA: 17.34 Jacobi stage: ITER: 68023 RES: 2.0948e-05
ELA: 17.44 Jacobi stage: ITER: 68355 RES: 1.9906e-05
ELA: 17.54 Jacobi stage: ITER: 68687 RES: 1.8917e-05
ELA: 17.64 Jacobi stage: ITER: 69019 RES: 1.7976e-05
ELA: 17.741 Jacobi stage: ITER: 69353 RES: 1.7072e-05
ELA: 17.841 Jacobi stage: ITER: 69685 RES: 1.6223e-05
ELA: 17.941 Jacobi stage: ITER: 70017 RES: 1.5417e-05
ELA: 18.041 Jacobi stage: ITER: 70349 RES: 1.465e-05
ELA: 18.141 Jacobi stage: ITER: 70683 RES: 1.3922e-05
ELA: 18.241 Jacobi stage: ITER: 71015 RES: 1.323e-05
ELA: 18.341 Jacobi stage: ITER: 71347 RES: 1.2572e-05
ELA: 18.441 Jacobi stage: ITER: 71679 RES: 1.1947e-05
ELA: 18.541 Jacobi stage: ITER: 72011 RES: 1.1353e-05
ELA: 18.642 Jacobi stage: ITER: 72345 RES: 1.0782e-05
ELA: 18.742 Jacobi stage: ITER: 72677 RES: 1.0246e-05
ELA: 18.842 Jacobi stage: ITER: 73011 RES: 9.7364e-06
ELA: 18.942 Jacobi stage: ITER: 73343 RES: 9.2523e-06
ELA: 19.042 Jacobi stage: ITER: 73675 RES: 8.7923e-06
ELA: 19.142 Jacobi stage: ITER: 74007 RES: 8.3552e-06
ELA: 19.242 Jacobi stage: ITER: 74339 RES: 7.9398e-06
ELA: 19.342 Jacobi stage: ITER: 74671 RES: 7.545e-06
ELA: 19.442 Jacobi stage: ITER: 75005 RES: 7.1655e-06
ELA: 19.542 Jacobi stage: ITER: 75337 RES: 6.8093e-06
ELA: 19.643 Jacobi stage: ITER: 75671 RES: 6.4707e-06
ELA: 19.743 Jacobi stage: ITER: 76003 RES: 6.149e-06
ELA: 19.843 Jacobi stage: ITER: 76335 RES: 5.8433e-06
ELA: 19.943 Jacobi stage: ITER: 76667 RES: 5.5528e-06
ELA: 20.043 Jacobi stage: ITER: 76999 RES: 5.2735e-06
ELA: 20.143 Jacobi stage: ITER: 77333 RES: 5.0113e-06
ELA: 20.243 Jacobi stage: ITER: 77667 RES: 4.7621e-06
ELA: 20.343 Jacobi stage: ITER: 77999 RES: 4.5254e-06
ELA: 20.443 Jacobi stage: ITER: 78331 RES: 4.3004e-06
ELA: 20.544 Jacobi stage: ITER: 78663 RES: 4.0866e-06
ELA: 20.644 Jacobi stage: ITER: 78995 RES: 3.8834e-06
ELA: 20.744 Jacobi stage: ITER: 79329 RES: 3.6881e-06
ELA: 20.844 Jacobi stage: ITER: 79661 RES: 3.5047e-06
ELA: 20.944 Jacobi stage: ITER: 79995 RES: 3.3305e-06
ELA: 21.044 Jacobi stage: ITER: 80327 RES: 3.1649e-06
ELA: 21.144 Jacobi stage: ITER: 80659 RES: 3.0075e-06
ELA: 21.244 Jacobi stage: ITER: 80991 RES: 2.858e-06
ELA: 21.344 Jacobi stage: ITER: 81323 RES: 2.7159e-06
ELA: 21.444 Jacobi stage: ITER: 81657 RES: 2.5793e-06
ELA: 21.545 Jacobi stage: ITER: 81989 RES: 2.4511e-06
ELA: 21.645 Jacobi stage: ITER: 82323 RES: 2.3292e-06
ELA: 21.745 Jacobi stage: ITER: 82655 RES: 2.2134e-06
ELA: 21.845 Jacobi stage: ITER: 82987 RES: 2.1034e-06
ELA: 21.945 Jacobi stage: ITER: 83319 RES: 1.9988e-06
ELA: 22.045 Jacobi stage: ITER: 83651 RES: 1.8994e-06
ELA: 22.145 Jacobi stage: ITER: 83985 RES: 1.8039e-06
ELA: 22.245 Jacobi stage: ITER: 84317 RES: 1.7142e-06
ELA: 22.345 Jacobi stage: ITER: 84649 RES: 1.629e-06
ELA: 22.446 Jacobi stage: ITER: 84983 RES: 1.548e-06
ELA: 22.546 Jacobi stage: ITER: 85315 RES: 1.471e-06
ELA: 22.646 Jacobi stage: ITER: 85647 RES: 1.3979e-06
ELA: 22.746 Jacobi stage: ITER: 85979 RES: 1.3284e-06
ELA: 22.846 Jacobi stage: ITER: 86311 RES: 1.2623e-06
ELA: 22.946 Jacobi stage: ITER: 86645 RES: 1.1988e-06
ELA: 23.046 Jacobi stage: ITER: 86977 RES: 1.1392e-06
ELA: 23.146 Jacobi stage: ITER: 87309 RES: 1.0826e-06
ELA: 23.246 Jacobi stage: ITER: 87643 RES: 1.0288e-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.15.0