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:57
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: 20865 RES: 0.029298
Jacobi stage: ITER: 45357 RES: 0.00068079
Jacobi stage: ITER: 70495 RES: 1.433e-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: 3899 RES: 0.40712
Jacobi stage: ITER: 8815 RES: 0.18671
Jacobi stage: ITER: 9803 RES: 0.16036
Jacobi stage: ITER: 10267 RES: 0.14931
Jacobi stage: ITER: 10639 RES: 0.141
Jacobi stage: ITER: 10967 RES: 0.13407
Jacobi stage: ITER: 11299 RES: 0.1274
Jacobi stage: ITER: 11625 RES: 0.12114
Jacobi stage: ITER: 11939 RES: 0.11547
Jacobi stage: ITER: 12267 RES: 0.10979
Jacobi stage: ITER: 12597 RES: 0.10433
Jacobi stage: ITER: 12929 RES: 0.099142
Jacobi stage: ITER: 13237 RES: 0.09456
Jacobi stage: ITER: 13563 RES: 0.089968
Jacobi stage: ITER: 13895 RES: 0.085495
Jacobi stage: ITER: 14227 RES: 0.081244
Jacobi stage: ITER: 14555 RES: 0.077251
Jacobi stage: ITER: 14859 RES: 0.073727
Jacobi stage: ITER: 15187 RES: 0.070104
Jacobi stage: ITER: 15511 RES: 0.066701
Jacobi stage: ITER: 15831 RES: 0.063501
Jacobi stage: ITER: 16165 RES: 0.060307
Jacobi stage: ITER: 16495 RES: 0.057344
Jacobi stage: ITER: 16815 RES: 0.054593
Jacobi stage: ITER: 17133 RES: 0.051975
Jacobi stage: ITER: 17451 RES: 0.049512
Jacobi stage: ITER: 17779 RES: 0.04708
Jacobi stage: ITER: 18107 RES: 0.044766
Jacobi stage: ITER: 18427 RES: 0.042619
Jacobi stage: ITER: 18755 RES: 0.040525
Jacobi stage: ITER: 19093 RES: 0.038463
Jacobi stage: ITER: 19407 RES: 0.036663
Jacobi stage: ITER: 19735 RES: 0.034862
Jacobi stage: ITER: 20043 RES: 0.033251
Jacobi stage: ITER: 20373 RES: 0.031598
Jacobi stage: ITER: 20671 RES: 0.030193
Jacobi stage: ITER: 20991 RES: 0.028745
Jacobi stage: ITER: 21307 RES: 0.027383
Jacobi stage: ITER: 21637 RES: 0.026022
Jacobi stage: ITER: 21963 RES: 0.024759
Jacobi stage: ITER: 22295 RES: 0.023528
Jacobi stage: ITER: 22627 RES: 0.022358
Jacobi stage: ITER: 22947 RES: 0.021285
Jacobi stage: ITER: 23273 RES: 0.02024
Jacobi stage: ITER: 23595 RES: 0.019269
Jacobi stage: ITER: 23923 RES: 0.018322
Jacobi stage: ITER: 24231 RES: 0.017476
Jacobi stage: ITER: 24511 RES: 0.01674
Jacobi stage: ITER: 24835 RES: 0.015927
Jacobi stage: ITER: 25149 RES: 0.015173
Jacobi stage: ITER: 25479 RES: 0.014427
Jacobi stage: ITER: 25807 RES: 0.013718
Jacobi stage: ITER: 26135 RES: 0.013044
Jacobi stage: ITER: 26467 RES: 0.012396
Jacobi stage: ITER: 26797 RES: 0.011779
Jacobi stage: ITER: 27129 RES: 0.011194
Jacobi stage: ITER: 27461 RES: 0.010637
Jacobi stage: ITER: 27759 RES: 0.010164
Jacobi stage: ITER: 28033 RES: 0.0097425
Jacobi stage: ITER: 28355 RES: 0.0092752
Jacobi stage: ITER: 28673 RES: 0.0088304
Jacobi stage: ITER: 28999 RES: 0.0084017
Jacobi stage: ITER: 29333 RES: 0.007979
Jacobi stage: ITER: 29659 RES: 0.0075917
Jacobi stage: ITER: 29979 RES: 0.0072275
Jacobi stage: ITER: 30299 RES: 0.0068809
Jacobi stage: ITER: 30603 RES: 0.006567
Jacobi stage: ITER: 30931 RES: 0.0062443
Jacobi stage: ITER: 31259 RES: 0.0059375
Jacobi stage: ITER: 31585 RES: 0.0056458
Jacobi stage: ITER: 31913 RES: 0.0053684
Jacobi stage: ITER: 32237 RES: 0.0051078
Jacobi stage: ITER: 32563 RES: 0.0048598
Jacobi stage: ITER: 32891 RES: 0.004621
Jacobi stage: ITER: 33219 RES: 0.004394
Jacobi stage: ITER: 33545 RES: 0.0041781
Jacobi stage: ITER: 33873 RES: 0.0039728
Jacobi stage: ITER: 34199 RES: 0.0037799
Jacobi stage: ITER: 34527 RES: 0.0035942
Jacobi stage: ITER: 34855 RES: 0.0034176
Jacobi stage: ITER: 35181 RES: 0.0032497
Jacobi stage: ITER: 35509 RES: 0.00309
Jacobi stage: ITER: 35835 RES: 0.00294
Jacobi stage: ITER: 36163 RES: 0.0027955
Jacobi stage: ITER: 36491 RES: 0.0026582
Jacobi stage: ITER: 36817 RES: 0.0025276
Jacobi stage: ITER: 37145 RES: 0.0024034
Jacobi stage: ITER: 37471 RES: 0.0022867
Jacobi stage: ITER: 37799 RES: 0.0021744
Jacobi stage: ITER: 38127 RES: 0.0020675
Jacobi stage: ITER: 38453 RES: 0.0019659
Jacobi stage: ITER: 38779 RES: 0.0018705
Jacobi stage: ITER: 39107 RES: 0.0017786
Jacobi stage: ITER: 39435 RES: 0.0016912
Jacobi stage: ITER: 39763 RES: 0.0016081
Jacobi stage: ITER: 40089 RES: 0.0015291
Jacobi stage: ITER: 40415 RES: 0.0014549
Jacobi stage: ITER: 40743 RES: 0.0013834
Jacobi stage: ITER: 41071 RES: 0.0013154
Jacobi stage: ITER: 41399 RES: 0.0012508
Jacobi stage: ITER: 41725 RES: 0.0011893
Jacobi stage: ITER: 42051 RES: 0.0011316
Jacobi stage: ITER: 42379 RES: 0.001076
Jacobi stage: ITER: 42707 RES: 0.0010231
Jacobi stage: ITER: 43035 RES: 0.00097285
Jacobi stage: ITER: 43361 RES: 0.00092505
Jacobi stage: ITER: 43689 RES: 0.0008796
Jacobi stage: ITER: 44015 RES: 0.0008369
Jacobi stage: ITER: 44343 RES: 0.00079578
Jacobi stage: ITER: 44671 RES: 0.00075668
Jacobi stage: ITER: 44997 RES: 0.0007195
Jacobi stage: ITER: 45323 RES: 0.00068415
Jacobi stage: ITER: 45651 RES: 0.00065093
Jacobi stage: ITER: 45979 RES: 0.00061895
Jacobi stage: ITER: 46315 RES: 0.00058782
Jacobi stage: ITER: 46643 RES: 0.00055894
Jacobi stage: ITER: 46971 RES: 0.00053147
Jacobi stage: ITER: 47297 RES: 0.00050536
Jacobi stage: ITER: 47625 RES: 0.00048053
Jacobi stage: ITER: 47951 RES: 0.0004572
Jacobi stage: ITER: 48279 RES: 0.00043474
Jacobi stage: ITER: 48607 RES: 0.00041338
Jacobi stage: ITER: 48933 RES: 0.00039307
Jacobi stage: ITER: 49261 RES: 0.00037375
Jacobi stage: ITER: 49587 RES: 0.00035561
Jacobi stage: ITER: 49915 RES: 0.00033814
Jacobi stage: ITER: 50243 RES: 0.00032152
Jacobi stage: ITER: 50569 RES: 0.00030573
Jacobi stage: ITER: 50895 RES: 0.00029088
Jacobi stage: ITER: 51223 RES: 0.00027659
Jacobi stage: ITER: 51551 RES: 0.000263
Jacobi stage: ITER: 51877 RES: 0.00025008
Jacobi stage: ITER: 52205 RES: 0.00023779
Jacobi stage: ITER: 52531 RES: 0.00022625
Jacobi stage: ITER: 52859 RES: 0.00021513
Jacobi stage: ITER: 53187 RES: 0.00020456
Jacobi stage: ITER: 53515 RES: 0.00019451
Jacobi stage: ITER: 53841 RES: 0.00018495
Jacobi stage: ITER: 54167 RES: 0.00017597
Jacobi stage: ITER: 54495 RES: 0.00016733
Jacobi stage: ITER: 54823 RES: 0.00015911
Jacobi stage: ITER: 55149 RES: 0.00015129
Jacobi stage: ITER: 55477 RES: 0.00014386
Jacobi stage: ITER: 55803 RES: 0.00013687
Jacobi stage: ITER: 56131 RES: 0.00013015
Jacobi stage: ITER: 56459 RES: 0.00012375
Jacobi stage: ITER: 56785 RES: 0.00011767
Jacobi stage: ITER: 57111 RES: 0.00011196
Jacobi stage: ITER: 57439 RES: 0.00010646
Jacobi stage: ITER: 57767 RES: 0.00010123
Jacobi stage: ITER: 58093 RES: 9.6254e-05
Jacobi stage: ITER: 58421 RES: 9.1525e-05
Jacobi stage: ITER: 58747 RES: 8.7081e-05
Jacobi stage: ITER: 59075 RES: 8.2803e-05
Jacobi stage: ITER: 59403 RES: 7.8734e-05
Jacobi stage: ITER: 59729 RES: 7.4866e-05
Jacobi stage: ITER: 60055 RES: 7.1231e-05
Jacobi stage: ITER: 60383 RES: 6.7731e-05
Jacobi stage: ITER: 60711 RES: 6.4404e-05
Jacobi stage: ITER: 61037 RES: 6.1239e-05
Jacobi stage: ITER: 61365 RES: 5.823e-05
Jacobi stage: ITER: 61691 RES: 5.5403e-05
Jacobi stage: ITER: 62019 RES: 5.2681e-05
Jacobi stage: ITER: 62347 RES: 5.0093e-05
Jacobi stage: ITER: 62673 RES: 4.7632e-05
Jacobi stage: ITER: 63001 RES: 4.5291e-05
Jacobi stage: ITER: 63327 RES: 4.3092e-05
Jacobi stage: ITER: 63655 RES: 4.0975e-05
Jacobi stage: ITER: 63983 RES: 3.8962e-05
Jacobi stage: ITER: 64309 RES: 3.7048e-05
Jacobi stage: ITER: 64637 RES: 3.5227e-05
Jacobi stage: ITER: 64963 RES: 3.3517e-05
Jacobi stage: ITER: 65291 RES: 3.187e-05
Jacobi stage: ITER: 65619 RES: 3.0304e-05
Jacobi stage: ITER: 65947 RES: 2.8815e-05
Jacobi stage: ITER: 66273 RES: 2.74e-05
Jacobi stage: ITER: 66601 RES: 2.6053e-05
Jacobi stage: ITER: 66927 RES: 2.4789e-05
Jacobi stage: ITER: 67255 RES: 2.3571e-05
Jacobi stage: ITER: 67583 RES: 2.2412e-05
Jacobi stage: ITER: 67909 RES: 2.1311e-05
Jacobi stage: ITER: 68235 RES: 2.0277e-05
Jacobi stage: ITER: 68563 RES: 1.928e-05
Jacobi stage: ITER: 68891 RES: 1.8333e-05
Jacobi stage: ITER: 69219 RES: 1.7432e-05
Jacobi stage: ITER: 69639 RES: 1.6343e-05
Jacobi stage: ITER: 69967 RES: 1.554e-05
Jacobi stage: ITER: 70295 RES: 1.4777e-05
Jacobi stage: ITER: 70621 RES: 1.4051e-05
Jacobi stage: ITER: 70949 RES: 1.336e-05
Jacobi stage: ITER: 71275 RES: 1.2712e-05
Jacobi stage: ITER: 71603 RES: 1.2087e-05
Jacobi stage: ITER: 71931 RES: 1.1493e-05
Jacobi stage: ITER: 72257 RES: 1.0929e-05
Jacobi stage: ITER: 72583 RES: 1.0398e-05
Jacobi stage: ITER: 72911 RES: 9.8871e-06
Jacobi stage: ITER: 73239 RES: 9.4013e-06
Jacobi stage: ITER: 73567 RES: 8.9394e-06
Jacobi stage: ITER: 73893 RES: 8.5002e-06
Jacobi stage: ITER: 74219 RES: 8.0875e-06
Jacobi stage: ITER: 74547 RES: 7.6901e-06
Jacobi stage: ITER: 74875 RES: 7.3123e-06
Jacobi stage: ITER: 75203 RES: 6.953e-06
Jacobi stage: ITER: 75529 RES: 6.6114e-06
Jacobi stage: ITER: 75855 RES: 6.2904e-06
Jacobi stage: ITER: 76183 RES: 5.9813e-06
Jacobi stage: ITER: 76511 RES: 5.6875e-06
Jacobi stage: ITER: 76837 RES: 5.408e-06
Jacobi stage: ITER: 77165 RES: 5.1423e-06
Jacobi stage: ITER: 77491 RES: 4.8926e-06
Jacobi stage: ITER: 77819 RES: 4.6523e-06
Jacobi stage: ITER: 78147 RES: 4.4237e-06
Jacobi stage: ITER: 78473 RES: 4.2063e-06
Jacobi stage: ITER: 78801 RES: 3.9997e-06
Jacobi stage: ITER: 79127 RES: 3.8055e-06
Jacobi stage: ITER: 79455 RES: 3.6185e-06
Jacobi stage: ITER: 79783 RES: 3.4407e-06
Jacobi stage: ITER: 80109 RES: 3.2717e-06
Jacobi stage: ITER: 80435 RES: 3.1128e-06
Jacobi stage: ITER: 80763 RES: 2.9599e-06
Jacobi stage: ITER: 81091 RES: 2.8145e-06
Jacobi stage: ITER: 81419 RES: 2.6762e-06
Jacobi stage: ITER: 81745 RES: 2.5447e-06
Jacobi stage: ITER: 82071 RES: 2.4211e-06
Jacobi stage: ITER: 82399 RES: 2.3022e-06
Jacobi stage: ITER: 82727 RES: 2.1891e-06
Jacobi stage: ITER: 83053 RES: 2.0815e-06
Jacobi stage: ITER: 83381 RES: 1.9792e-06
Jacobi stage: ITER: 83707 RES: 1.8831e-06
Jacobi stage: ITER: 84035 RES: 1.7906e-06
Jacobi stage: ITER: 84363 RES: 1.7026e-06
Jacobi stage: ITER: 84689 RES: 1.619e-06
Jacobi stage: ITER: 85017 RES: 1.5394e-06
Jacobi stage: ITER: 85343 RES: 1.4647e-06
Jacobi stage: ITER: 85671 RES: 1.3927e-06
Jacobi stage: ITER: 85999 RES: 1.3243e-06
Jacobi stage: ITER: 86327 RES: 1.2592e-06
Jacobi stage: ITER: 86653 RES: 1.1974e-06
Jacobi stage: ITER: 86979 RES: 1.1392e-06
Jacobi stage: ITER: 87307 RES: 1.0833e-06
Jacobi stage: ITER: 87635 RES: 1.03e-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:1.2412e-
ELA: 0.1001 Jacobi stage: ITER: 21439 RES: 0.026834
ELA: 0.20017 Jacobi stage: ITER: 44269 RES: 0.00080463
ELA: 0.30025 Jacobi stage: ITER: 61951 RES: 5.3234e-05
ELA: 0.40035 Jacobi stage: ITER: 82947 RES: 2.1163e-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:2.2224e-
ELA: 0.10011 Jacobi stage: ITER: 431 RES: 0.86917
ELA: 0.20021 Jacobi stage: ITER: 899 RES: 0.75887
ELA: 0.30032 Jacobi stage: ITER: 1265 RES: 0.68948
ELA: 0.40042 Jacobi stage: ITER: 1591 RES: 0.6378
ELA: 0.50189 Jacobi stage: ITER: 1925 RES: 0.59148
ELA: 0.60198 Jacobi stage: ITER: 2251 RES: 0.55239
ELA: 0.70209 Jacobi stage: ITER: 2559 RES: 0.51931
ELA: 0.8022 Jacobi stage: ITER: 2885 RES: 0.48767
ELA: 0.90269 Jacobi stage: ITER: 3217 RES: 0.45874
ELA: 1.0028 Jacobi stage: ITER: 3547 RES: 0.4327
ELA: 1.1029 Jacobi stage: ITER: 3857 RES: 0.40991
ELA: 1.203 Jacobi stage: ITER: 4173 RES: 0.38852
ELA: 1.3031 Jacobi stage: ITER: 4505 RES: 0.36763
ELA: 1.4032 Jacobi stage: ITER: 4835 RES: 0.3484
ELA: 1.5033 Jacobi stage: ITER: 5163 RES: 0.33039
ELA: 1.6034 Jacobi stage: ITER: 5465 RES: 0.31468
ELA: 1.7035 Jacobi stage: ITER: 5791 RES: 0.29889
ELA: 1.8036 Jacobi stage: ITER: 6109 RES: 0.28417
ELA: 1.9037 Jacobi stage: ITER: 6431 RES: 0.27024
ELA: 2.0038 Jacobi stage: ITER: 6759 RES: 0.25674
ELA: 2.1039 Jacobi stage: ITER: 7087 RES: 0.24395
ELA: 2.204 Jacobi stage: ITER: 7411 RES: 0.23197
ELA: 2.3042 Jacobi stage: ITER: 7727 RES: 0.22088
ELA: 2.4043 Jacobi stage: ITER: 8041 RES: 0.21035
ELA: 2.5044 Jacobi stage: ITER: 8365 RES: 0.20007
ELA: 2.6045 Jacobi stage: ITER: 8691 RES: 0.19032
ELA: 2.7046 Jacobi stage: ITER: 9007 RES: 0.18127
ELA: 2.8047 Jacobi stage: ITER: 9335 RES: 0.17234
ELA: 2.9049 Jacobi stage: ITER: 9655 RES: 0.16405
ELA: 3.0053 Jacobi stage: ITER: 9967 RES: 0.15636
ELA: 3.1053 Jacobi stage: ITER: 10295 RES: 0.14867
ELA: 3.2055 Jacobi stage: ITER: 10603 RES: 0.14179
ELA: 3.3056 Jacobi stage: ITER: 10933 RES: 0.13473
ELA: 3.4057 Jacobi stage: ITER: 11229 RES: 0.12874
ELA: 3.5058 Jacobi stage: ITER: 11545 RES: 0.12263
ELA: 3.6059 Jacobi stage: ITER: 11857 RES: 0.11689
ELA: 3.7059 Jacobi stage: ITER: 12183 RES: 0.11122
ELA: 3.806 Jacobi stage: ITER: 12511 RES: 0.10575
ELA: 3.9061 Jacobi stage: ITER: 12839 RES: 0.10055
ELA: 4.0062 Jacobi stage: ITER: 13171 RES: 0.095553
ELA: 4.1063 Jacobi stage: ITER: 13487 RES: 0.091025
ELA: 4.2064 Jacobi stage: ITER: 13815 RES: 0.086552
ELA: 4.3065 Jacobi stage: ITER: 14131 RES: 0.082451
ELA: 4.4066 Jacobi stage: ITER: 14459 RES: 0.078399
ELA: 4.5067 Jacobi stage: ITER: 14767 RES: 0.074776
ELA: 4.6069 Jacobi stage: ITER: 15047 RES: 0.071628
ELA: 4.707 Jacobi stage: ITER: 15365 RES: 0.068193
ELA: 4.8071 Jacobi stage: ITER: 15679 RES: 0.065001
ELA: 4.9072 Jacobi stage: ITER: 16005 RES: 0.061808
ELA: 5.0073 Jacobi stage: ITER: 16329 RES: 0.058807
ELA: 5.1074 Jacobi stage: ITER: 16657 RES: 0.055917
ELA: 5.2075 Jacobi stage: ITER: 16987 RES: 0.05317
ELA: 5.3076 Jacobi stage: ITER: 17319 RES: 0.050526
ELA: 5.4077 Jacobi stage: ITER: 17651 RES: 0.048014
ELA: 5.5078 Jacobi stage: ITER: 17981 RES: 0.045627
ELA: 5.6079 Jacobi stage: ITER: 18279 RES: 0.043599
ELA: 5.708 Jacobi stage: ITER: 18555 RES: 0.041789
ELA: 5.8081 Jacobi stage: ITER: 18863 RES: 0.039858
ELA: 5.9082 Jacobi stage: ITER: 19181 RES: 0.037947
ELA: 6.0083 Jacobi stage: ITER: 19505 RES: 0.036104
ELA: 6.1088 Jacobi stage: ITER: 19837 RES: 0.034309
ELA: 6.209 Jacobi stage: ITER: 20163 RES: 0.032644
ELA: 6.3091 Jacobi stage: ITER: 20479 RES: 0.031097
ELA: 6.4092 Jacobi stage: ITER: 20801 RES: 0.029587
ELA: 6.5093 Jacobi stage: ITER: 21101 RES: 0.028255
ELA: 6.6094 Jacobi stage: ITER: 21429 RES: 0.026867
ELA: 6.7095 Jacobi stage: ITER: 21755 RES: 0.025562
ELA: 6.8096 Jacobi stage: ITER: 22083 RES: 0.024306
ELA: 6.9097 Jacobi stage: ITER: 22411 RES: 0.023112
ELA: 7.0097 Jacobi stage: ITER: 22739 RES: 0.021977
ELA: 7.1098 Jacobi stage: ITER: 23067 RES: 0.020897
ELA: 7.2099 Jacobi stage: ITER: 23395 RES: 0.01987
ELA: 7.31 Jacobi stage: ITER: 23723 RES: 0.018894
ELA: 7.4101 Jacobi stage: ITER: 24049 RES: 0.017965
ELA: 7.5102 Jacobi stage: ITER: 24375 RES: 0.017093
ELA: 7.6103 Jacobi stage: ITER: 24703 RES: 0.016253
ELA: 7.7104 Jacobi stage: ITER: 25031 RES: 0.015455
ELA: 7.8104 Jacobi stage: ITER: 25359 RES: 0.014695
ELA: 7.9105 Jacobi stage: ITER: 25685 RES: 0.013973
ELA: 8.0106 Jacobi stage: ITER: 26011 RES: 0.013295
ELA: 8.1107 Jacobi stage: ITER: 26339 RES: 0.012642
ELA: 8.2108 Jacobi stage: ITER: 26667 RES: 0.012021
ELA: 8.3109 Jacobi stage: ITER: 26995 RES: 0.01143
ELA: 8.411 Jacobi stage: ITER: 27321 RES: 0.010868
ELA: 8.511 Jacobi stage: ITER: 27647 RES: 0.010341
ELA: 8.6111 Jacobi stage: ITER: 27975 RES: 0.0098327
ELA: 8.7112 Jacobi stage: ITER: 28303 RES: 0.0093496
ELA: 8.8113 Jacobi stage: ITER: 28631 RES: 0.0088902
ELA: 8.9116 Jacobi stage: ITER: 28957 RES: 0.0084534
ELA: 9.0116 Jacobi stage: ITER: 29285 RES: 0.0080381
ELA: 9.1117 Jacobi stage: ITER: 29611 RES: 0.0076478
ELA: 9.2118 Jacobi stage: ITER: 29939 RES: 0.0072721
ELA: 9.3119 Jacobi stage: ITER: 30267 RES: 0.0069148
ELA: 9.412 Jacobi stage: ITER: 30593 RES: 0.006575
ELA: 9.5121 Jacobi stage: ITER: 30921 RES: 0.006252
ELA: 9.6122 Jacobi stage: ITER: 31247 RES: 0.0059485
ELA: 9.7123 Jacobi stage: ITER: 31575 RES: 0.0056562
ELA: 9.8124 Jacobi stage: ITER: 31903 RES: 0.0053783
ELA: 9.9125 Jacobi stage: ITER: 32229 RES: 0.005114
ELA: 10.013 Jacobi stage: ITER: 32557 RES: 0.0048628
ELA: 10.113 Jacobi stage: ITER: 32883 RES: 0.0046267
ELA: 10.213 Jacobi stage: ITER: 33211 RES: 0.0043994
ELA: 10.313 Jacobi stage: ITER: 33539 RES: 0.0041832
ELA: 10.413 Jacobi stage: ITER: 33867 RES: 0.0039777
ELA: 10.513 Jacobi stage: ITER: 34193 RES: 0.0037822
ELA: 10.613 Jacobi stage: ITER: 34519 RES: 0.0035986
ELA: 10.713 Jacobi stage: ITER: 34847 RES: 0.0034218
ELA: 10.813 Jacobi stage: ITER: 35175 RES: 0.0032537
ELA: 10.913 Jacobi stage: ITER: 35503 RES: 0.0030938
ELA: 11.013 Jacobi stage: ITER: 35829 RES: 0.0029418
ELA: 11.113 Jacobi stage: ITER: 36155 RES: 0.002799
ELA: 11.213 Jacobi stage: ITER: 36483 RES: 0.0026615
ELA: 11.314 Jacobi stage: ITER: 36811 RES: 0.0025307
ELA: 11.414 Jacobi stage: ITER: 37137 RES: 0.0024063
ELA: 11.514 Jacobi stage: ITER: 37465 RES: 0.0022881
ELA: 11.614 Jacobi stage: ITER: 37791 RES: 0.002177
ELA: 11.714 Jacobi stage: ITER: 38119 RES: 0.0020701
ELA: 11.814 Jacobi stage: ITER: 38447 RES: 0.0019684
ELA: 11.914 Jacobi stage: ITER: 38773 RES: 0.0018716
ELA: 12.014 Jacobi stage: ITER: 39099 RES: 0.0017808
ELA: 12.114 Jacobi stage: ITER: 39427 RES: 0.0016933
ELA: 12.214 Jacobi stage: ITER: 39755 RES: 0.0016101
ELA: 12.314 Jacobi stage: ITER: 40083 RES: 0.001531
ELA: 12.414 Jacobi stage: ITER: 40409 RES: 0.0014558
ELA: 12.515 Jacobi stage: ITER: 40735 RES: 0.0013851
ELA: 12.615 Jacobi stage: ITER: 41063 RES: 0.001317
ELA: 12.715 Jacobi stage: ITER: 41391 RES: 0.0012523
ELA: 12.815 Jacobi stage: ITER: 41717 RES: 0.0011908
ELA: 12.915 Jacobi stage: ITER: 42045 RES: 0.0011323
ELA: 13.015 Jacobi stage: ITER: 42371 RES: 0.0010773
ELA: 13.115 Jacobi stage: ITER: 42699 RES: 0.0010244
ELA: 13.215 Jacobi stage: ITER: 43027 RES: 0.00097405
ELA: 13.315 Jacobi stage: ITER: 43353 RES: 0.00092619
ELA: 13.415 Jacobi stage: ITER: 43681 RES: 0.00088068
ELA: 13.515 Jacobi stage: ITER: 44007 RES: 0.00083793
ELA: 13.615 Jacobi stage: ITER: 44335 RES: 0.00079676
ELA: 13.715 Jacobi stage: ITER: 44663 RES: 0.00075761
ELA: 13.816 Jacobi stage: ITER: 44989 RES: 0.00072039
ELA: 13.916 Jacobi stage: ITER: 45317 RES: 0.00068499
ELA: 14.016 Jacobi stage: ITER: 45643 RES: 0.00065173
ELA: 14.116 Jacobi stage: ITER: 45971 RES: 0.00061971
ELA: 14.216 Jacobi stage: ITER: 46299 RES: 0.00058926
ELA: 14.316 Jacobi stage: ITER: 46625 RES: 0.00056031
ELA: 14.416 Jacobi stage: ITER: 46953 RES: 0.00053278
ELA: 14.516 Jacobi stage: ITER: 47279 RES: 0.00050692
ELA: 14.616 Jacobi stage: ITER: 47607 RES: 0.00048201
ELA: 14.716 Jacobi stage: ITER: 47935 RES: 0.00045833
ELA: 14.816 Jacobi stage: ITER: 48261 RES: 0.00043581
ELA: 14.916 Jacobi stage: ITER: 48587 RES: 0.00041465
ELA: 15.016 Jacobi stage: ITER: 48915 RES: 0.00039428
ELA: 15.117 Jacobi stage: ITER: 49243 RES: 0.0003749
ELA: 15.217 Jacobi stage: ITER: 49571 RES: 0.00035648
ELA: 15.317 Jacobi stage: ITER: 49897 RES: 0.00033897
ELA: 15.417 Jacobi stage: ITER: 50223 RES: 0.00032251
ELA: 15.517 Jacobi stage: ITER: 50551 RES: 0.00030667
ELA: 15.617 Jacobi stage: ITER: 50879 RES: 0.0002916
ELA: 15.717 Jacobi stage: ITER: 51205 RES: 0.00027727
ELA: 15.817 Jacobi stage: ITER: 51533 RES: 0.00026365
ELA: 15.917 Jacobi stage: ITER: 51859 RES: 0.00025085
ELA: 16.017 Jacobi stage: ITER: 52187 RES: 0.00023852
ELA: 16.117 Jacobi stage: ITER: 52515 RES: 0.0002268
ELA: 16.217 Jacobi stage: ITER: 52841 RES: 0.00021566
ELA: 16.317 Jacobi stage: ITER: 53167 RES: 0.00020519
ELA: 16.418 Jacobi stage: ITER: 53495 RES: 0.00019511
ELA: 16.518 Jacobi stage: ITER: 53823 RES: 0.00018552
ELA: 16.618 Jacobi stage: ITER: 54151 RES: 0.00017641
ELA: 16.718 Jacobi stage: ITER: 54477 RES: 0.00016774
ELA: 16.818 Jacobi stage: ITER: 54805 RES: 0.0001595
ELA: 16.918 Jacobi stage: ITER: 55131 RES: 0.00015175
ELA: 17.018 Jacobi stage: ITER: 55459 RES: 0.0001443
ELA: 17.118 Jacobi stage: ITER: 55787 RES: 0.00013721
ELA: 17.218 Jacobi stage: ITER: 56115 RES: 0.00013047
ELA: 17.318 Jacobi stage: ITER: 56441 RES: 0.00012406
ELA: 17.418 Jacobi stage: ITER: 56769 RES: 0.00011796
ELA: 17.518 Jacobi stage: ITER: 57095 RES: 0.00011223
ELA: 17.619 Jacobi stage: ITER: 57423 RES: 0.00010672
ELA: 17.719 Jacobi stage: ITER: 57751 RES: 0.00010148
ELA: 17.819 Jacobi stage: ITER: 58077 RES: 9.6491e-05
ELA: 17.919 Jacobi stage: ITER: 58405 RES: 9.175e-05
ELA: 18.019 Jacobi stage: ITER: 58731 RES: 8.7296e-05
ELA: 18.119 Jacobi stage: ITER: 59059 RES: 8.3006e-05
ELA: 18.219 Jacobi stage: ITER: 59387 RES: 7.8928e-05
ELA: 18.319 Jacobi stage: ITER: 59713 RES: 7.505e-05
ELA: 18.45 Jacobi stage: ITER: 60143 RES: 7.0275e-05
ELA: 18.551 Jacobi stage: ITER: 60471 RES: 6.6822e-05
ELA: 18.651 Jacobi stage: ITER: 60797 RES: 6.3539e-05
ELA: 18.751 Jacobi stage: ITER: 61123 RES: 6.0454e-05
ELA: 18.851 Jacobi stage: ITER: 61451 RES: 5.7484e-05
ELA: 18.951 Jacobi stage: ITER: 61779 RES: 5.4659e-05
ELA: 19.051 Jacobi stage: ITER: 62105 RES: 5.1974e-05
ELA: 19.151 Jacobi stage: ITER: 62433 RES: 4.942e-05
ELA: 19.251 Jacobi stage: ITER: 62759 RES: 4.7021e-05
ELA: 19.351 Jacobi stage: ITER: 63087 RES: 4.4711e-05
ELA: 19.451 Jacobi stage: ITER: 63415 RES: 4.2514e-05
ELA: 19.551 Jacobi stage: ITER: 63741 RES: 4.0425e-05
ELA: 19.651 Jacobi stage: ITER: 64069 RES: 3.8439e-05
ELA: 19.751 Jacobi stage: ITER: 64395 RES: 3.6573e-05
ELA: 19.852 Jacobi stage: ITER: 64723 RES: 3.4776e-05
ELA: 19.952 Jacobi stage: ITER: 65051 RES: 3.3067e-05
ELA: 20.052 Jacobi stage: ITER: 65377 RES: 3.1442e-05
ELA: 20.152 Jacobi stage: ITER: 65705 RES: 2.9897e-05
ELA: 20.252 Jacobi stage: ITER: 66031 RES: 2.8446e-05
ELA: 20.352 Jacobi stage: ITER: 66359 RES: 2.7048e-05
ELA: 20.452 Jacobi stage: ITER: 66687 RES: 2.5719e-05
ELA: 20.552 Jacobi stage: ITER: 67013 RES: 2.4456e-05
ELA: 20.652 Jacobi stage: ITER: 67341 RES: 2.3254e-05
ELA: 20.752 Jacobi stage: ITER: 67667 RES: 2.2125e-05
ELA: 20.852 Jacobi stage: ITER: 67995 RES: 2.1038e-05
ELA: 20.952 Jacobi stage: ITER: 68323 RES: 2.0004e-05
ELA: 21.052 Jacobi stage: ITER: 68649 RES: 1.9022e-05
ELA: 21.153 Jacobi stage: ITER: 68975 RES: 1.8098e-05
ELA: 21.253 Jacobi stage: ITER: 69303 RES: 1.7209e-05
ELA: 21.353 Jacobi stage: ITER: 69631 RES: 1.6363e-05
ELA: 21.453 Jacobi stage: ITER: 69959 RES: 1.5559e-05
ELA: 21.553 Jacobi stage: ITER: 70285 RES: 1.4795e-05
ELA: 21.653 Jacobi stage: ITER: 70611 RES: 1.4077e-05
ELA: 21.753 Jacobi stage: ITER: 70939 RES: 1.3385e-05
ELA: 21.853 Jacobi stage: ITER: 71267 RES: 1.2727e-05
ELA: 21.953 Jacobi stage: ITER: 71593 RES: 1.2102e-05
ELA: 22.053 Jacobi stage: ITER: 71921 RES: 1.1507e-05
ELA: 22.153 Jacobi stage: ITER: 72247 RES: 1.0949e-05
ELA: 22.253 Jacobi stage: ITER: 72575 RES: 1.0411e-05
ELA: 22.354 Jacobi stage: ITER: 72903 RES: 9.8992e-06
ELA: 22.454 Jacobi stage: ITER: 73229 RES: 9.4129e-06
ELA: 22.554 Jacobi stage: ITER: 73557 RES: 8.9504e-06
ELA: 22.654 Jacobi stage: ITER: 73883 RES: 8.5158e-06
ELA: 22.754 Jacobi stage: ITER: 74211 RES: 8.0974e-06
ELA: 22.854 Jacobi stage: ITER: 74539 RES: 7.6996e-06
ELA: 22.954 Jacobi stage: ITER: 74865 RES: 7.3213e-06
ELA: 23.054 Jacobi stage: ITER: 75191 RES: 6.9658e-06
ELA: 23.154 Jacobi stage: ITER: 75519 RES: 6.6236e-06
ELA: 23.254 Jacobi stage: ITER: 75847 RES: 6.2981e-06
ELA: 23.354 Jacobi stage: ITER: 76175 RES: 5.9887e-06
ELA: 23.455 Jacobi stage: ITER: 76503 RES: 5.6944e-06
ELA: 23.555 Jacobi stage: ITER: 76829 RES: 5.4147e-06
ELA: 23.655 Jacobi stage: ITER: 77155 RES: 5.1518e-06
ELA: 23.755 Jacobi stage: ITER: 77483 RES: 4.8987e-06
ELA: 23.855 Jacobi stage: ITER: 77811 RES: 4.658e-06
ELA: 23.955 Jacobi stage: ITER: 78137 RES: 4.4291e-06
ELA: 24.055 Jacobi stage: ITER: 78531 RES: 4.1703e-06
ELA: 24.155 Jacobi stage: ITER: 79027 RES: 3.8644e-06
ELA: 24.255 Jacobi stage: ITER: 79519 RES: 3.5831e-06
ELA: 24.355 Jacobi stage: ITER: 80011 RES: 3.3223e-06
ELA: 24.455 Jacobi stage: ITER: 80503 RES: 3.0805e-06
ELA: 24.556 Jacobi stage: ITER: 80995 RES: 2.8563e-06
ELA: 24.656 Jacobi stage: ITER: 81487 RES: 2.6484e-06
ELA: 24.756 Jacobi stage: ITER: 81979 RES: 2.4556e-06
ELA: 24.856 Jacobi stage: ITER: 82471 RES: 2.2769e-06
ELA: 24.956 Jacobi stage: ITER: 82963 RES: 2.1111e-06
ELA: 25.056 Jacobi stage: ITER: 83455 RES: 1.9575e-06
ELA: 25.156 Jacobi stage: ITER: 83947 RES: 1.815e-06
ELA: 25.256 Jacobi stage: ITER: 84439 RES: 1.6829e-06
ELA: 25.356 Jacobi stage: ITER: 84931 RES: 1.5604e-06
ELA: 25.457 Jacobi stage: ITER: 85423 RES: 1.4468e-06
ELA: 25.557 Jacobi stage: ITER: 85915 RES: 1.3415e-06
ELA: 25.657 Jacobi stage: ITER: 86425 RES: 1.24e-06
ELA: 25.757 Jacobi stage: ITER: 86917 RES: 1.1498e-06
ELA: 25.857 Jacobi stage: ITER: 87409 RES: 1.0661e-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::LinearSolver< MatrixType > > getLinearSolver(std::string name)
Function returning shared pointer with linear solver given by its name in a form of a string.
Definition LinearSolverTypeResolver.h:84
std::shared_ptr< Linear::Preconditioners::Preconditioner< MatrixType > > getPreconditioner(std::string name)
Function returning shared pointer with linear preconditioner given by its name in a form of a string.
Definition LinearSolverTypeResolver.h:138
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 Tue Mar 25 2025 19:47:25 for Template Numerical Library by
1.13.2