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
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
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 <TNL/Matrices/SparseMatrix.h>
4#include <TNL/Devices/Sequential.h>
5#include <TNL/Devices/Cuda.h>
6#include <TNL/Solvers/Linear/Jacobi.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 * Setup solver of the linear system.
63 */
65 LinearSolver solver;
66 solver.setMatrix( matrix_ptr );
67 solver.setOmega( 0.0005 );
68
69 /***
70 * Setup monitor of the iterative solver.
71 */
73 IterativeSolverMonitorType monitor;
74 TNL::Solvers::SolverMonitorThread monitorThread( monitor );
75 monitor.setRefreshRate( 10 ); // refresh rate in milliseconds
76 monitor.setVerbose( 1 );
77 monitor.setStage( "Jacobi stage:" );
78 solver.setSolverMonitor( monitor );
79 solver.setConvergenceResidue( 1.0e-6 );
80 solver.solve( b, x );
81 monitor.stopMainLoop();
83}
84
85int
86main( int argc, char* argv[] )
87{
89 iterativeLinearSolverExample< TNL::Devices::Sequential >();
90
91#ifdef __CUDACC__
93 iterativeLinearSolverExample< TNL::Devices::Cuda >();
94#endif
95}
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:142
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: 48115 RES: 0.00044583
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: 2391 RES: 0.53691
Jacobi stage: ITER: 4335 RES: 0.37826
Jacobi stage: ITER: 7687 RES: 0.22226
Jacobi stage: ITER: 11409 RES: 0.12523
Jacobi stage: ITER: 15399 RES: 0.067858
Jacobi stage: ITER: 17579 RES: 0.048548
Jacobi stage: ITER: 20783 RES: 0.029678
Jacobi stage: ITER: 22735 RES: 0.02199
Jacobi stage: ITER: 25455 RES: 0.01448
Jacobi stage: ITER: 28467 RES: 0.009117
Jacobi stage: ITER: 29487 RES: 0.0077949
Jacobi stage: ITER: 30059 RES: 0.0071393
Jacobi stage: ITER: 30607 RES: 0.0065629
Jacobi stage: ITER: 30965 RES: 0.0062099
Jacobi stage: ITER: 31351 RES: 0.0058542
Jacobi stage: ITER: 31581 RES: 0.0056492
Jacobi stage: ITER: 31965 RES: 0.0053257
Jacobi stage: ITER: 32349 RES: 0.0050206
Jacobi stage: ITER: 32599 RES: 0.004833
Jacobi stage: ITER: 32971 RES: 0.0045646
Jacobi stage: ITER: 33307 RES: 0.004335
Jacobi stage: ITER: 33587 RES: 0.0041525
Jacobi stage: ITER: 33971 RES: 0.0039146
Jacobi stage: ITER: 34191 RES: 0.0037846
Jacobi stage: ITER: 34563 RES: 0.0035744
Jacobi stage: ITER: 34907 RES: 0.0033904
Jacobi stage: ITER: 35147 RES: 0.0032677
Jacobi stage: ITER: 35519 RES: 0.0030862
Jacobi stage: ITER: 35749 RES: 0.0029782
Jacobi stage: ITER: 36123 RES: 0.0028128
Jacobi stage: ITER: 36495 RES: 0.0026566
Jacobi stage: ITER: 36751 RES: 0.0025541
Jacobi stage: ITER: 37123 RES: 0.0024123
Jacobi stage: ITER: 37511 RES: 0.0022727
Jacobi stage: ITER: 37727 RES: 0.0021985
Jacobi stage: ITER: 38099 RES: 0.0020764
Jacobi stage: ITER: 38399 RES: 0.0019829
Jacobi stage: ITER: 38659 RES: 0.0019053
Jacobi stage: ITER: 38899 RES: 0.0018363
Jacobi stage: ITER: 39035 RES: 0.0017984
Jacobi stage: ITER: 39271 RES: 0.0017343
Jacobi stage: ITER: 39491 RES: 0.0016767
Jacobi stage: ITER: 39633 RES: 0.00164
Jacobi stage: ITER: 39875 RES: 0.0015807
Jacobi stage: ITER: 40087 RES: 0.00153
Jacobi stage: ITER: 40267 RES: 0.0014883
Jacobi stage: ITER: 40507 RES: 0.0014344
Jacobi stage: ITER: 40655 RES: 0.0014022
Jacobi stage: ITER: 40903 RES: 0.0013498
Jacobi stage: ITER: 41203 RES: 0.001289
Jacobi stage: ITER: 41441 RES: 0.0012424
Jacobi stage: ITER: 41833 RES: 0.0011698
Jacobi stage: ITER: 42059 RES: 0.0011302
Jacobi stage: ITER: 42427 RES: 0.0010681
Jacobi stage: ITER: 42801 RES: 0.0010081
Jacobi stage: ITER: 43023 RES: 0.00097465
Jacobi stage: ITER: 43393 RES: 0.00092052
Jacobi stage: ITER: 43691 RES: 0.0008796
Jacobi stage: ITER: 44003 RES: 0.00083844
Jacobi stage: ITER: 44387 RES: 0.00079042
Jacobi stage: ITER: 44615 RES: 0.00076322
Jacobi stage: ITER: 44973 RES: 0.00072216
Jacobi stage: ITER: 45347 RES: 0.00068205
Jacobi stage: ITER: 45573 RES: 0.00065858
Jacobi stage: ITER: 45967 RES: 0.00062009
Jacobi stage: ITER: 46203 RES: 0.00059802
Jacobi stage: ITER: 46559 RES: 0.0005662
Jacobi stage: ITER: 46941 RES: 0.00053377
Jacobi stage: ITER: 47167 RES: 0.00051571
Jacobi stage: ITER: 47543 RES: 0.00048677
Jacobi stage: ITER: 47941 RES: 0.00045776
Jacobi stage: ITER: 48175 RES: 0.00044174
Jacobi stage: ITER: 48561 RES: 0.00041618
Jacobi stage: ITER: 48787 RES: 0.00040211
Jacobi stage: ITER: 49171 RES: 0.00037907
Jacobi stage: ITER: 49527 RES: 0.0003589
Jacobi stage: ITER: 49739 RES: 0.0003474
Jacobi stage: ITER: 50089 RES: 0.00032912
Jacobi stage: ITER: 50319 RES: 0.00031779
Jacobi stage: ITER: 50689 RES: 0.00030014
Jacobi stage: ITER: 51071 RES: 0.00028312
Jacobi stage: ITER: 51299 RES: 0.00027338
Jacobi stage: ITER: 51689 RES: 0.00025741
Jacobi stage: ITER: 52083 RES: 0.00024236
Jacobi stage: ITER: 52315 RES: 0.00023388
Jacobi stage: ITER: 52689 RES: 0.00022076
Jacobi stage: ITER: 52907 RES: 0.00021355
Jacobi stage: ITER: 53285 RES: 0.00020144
Jacobi stage: ITER: 53633 RES: 0.00019096
Jacobi stage: ITER: 53883 RES: 0.00018382
Jacobi stage: ITER: 54251 RES: 0.00017372
Jacobi stage: ITER: 54469 RES: 0.00016795
Jacobi stage: ITER: 54843 RES: 0.00015862
Jacobi stage: ITER: 55209 RES: 0.0001499
Jacobi stage: ITER: 55431 RES: 0.00014492
Jacobi stage: ITER: 55815 RES: 0.00013662
Jacobi stage: ITER: 56159 RES: 0.00012959
Jacobi stage: ITER: 56415 RES: 0.00012459
Jacobi stage: ITER: 56791 RES: 0.0001176
Jacobi stage: ITER: 57019 RES: 0.00011355
Jacobi stage: ITER: 57405 RES: 0.00010698
Jacobi stage: ITER: 57753 RES: 0.00010141
Jacobi stage: ITER: 57983 RES: 9.7924e-05
Jacobi stage: ITER: 58343 RES: 9.2656e-05
Jacobi stage: ITER: 58553 RES: 8.9688e-05
Jacobi stage: ITER: 58927 RES: 8.4707e-05
Jacobi stage: ITER: 59291 RES: 8.0101e-05
Jacobi stage: ITER: 59533 RES: 7.7154e-05
Jacobi stage: ITER: 59919 RES: 7.2735e-05
Jacobi stage: ITER: 60225 RES: 6.9374e-05
Jacobi stage: ITER: 60529 RES: 6.6209e-05
Jacobi stage: ITER: 60923 RES: 6.234e-05
Jacobi stage: ITER: 61145 RES: 6.0232e-05
Jacobi stage: ITER: 61525 RES: 5.6817e-05
Jacobi stage: ITER: 61887 RES: 5.376e-05
Jacobi stage: ITER: 62139 RES: 5.1719e-05
Jacobi stage: ITER: 62515 RES: 4.8817e-05
Jacobi stage: ITER: 62767 RES: 4.6963e-05
Jacobi stage: ITER: 63115 RES: 4.4519e-05
Jacobi stage: ITER: 63475 RES: 4.2124e-05
Jacobi stage: ITER: 63699 RES: 4.0699e-05
Jacobi stage: ITER: 64081 RES: 3.8368e-05
Jacobi stage: ITER: 64467 RES: 3.617e-05
Jacobi stage: ITER: 64719 RES: 3.4797e-05
Jacobi stage: ITER: 65097 RES: 3.2824e-05
Jacobi stage: ITER: 65401 RES: 3.1327e-05
Jacobi stage: ITER: 65709 RES: 2.9879e-05
Jacobi stage: ITER: 66095 RES: 2.8168e-05
Jacobi stage: ITER: 66323 RES: 2.7198e-05
Jacobi stage: ITER: 66699 RES: 2.5672e-05
Jacobi stage: ITER: 67039 RES: 2.4366e-05
Jacobi stage: ITER: 67283 RES: 2.3469e-05
Jacobi stage: ITER: 67633 RES: 2.2234e-05
Jacobi stage: ITER: 67857 RES: 2.1482e-05
Jacobi stage: ITER: 68207 RES: 2.0364e-05
Jacobi stage: ITER: 68587 RES: 1.9209e-05
Jacobi stage: ITER: 68835 RES: 1.8491e-05
Jacobi stage: ITER: 69215 RES: 1.7443e-05
Jacobi stage: ITER: 69541 RES: 1.6586e-05
Jacobi stage: ITER: 69795 RES: 1.5956e-05
Jacobi stage: ITER: 70151 RES: 1.5107e-05
Jacobi stage: ITER: 70367 RES: 1.4614e-05
Jacobi stage: ITER: 70733 RES: 1.3811e-05
Jacobi stage: ITER: 71001 RES: 1.3254e-05
Jacobi stage: ITER: 71319 RES: 1.2626e-05
Jacobi stage: ITER: 71703 RES: 1.1903e-05
Jacobi stage: ITER: 71921 RES: 1.1507e-05
Jacobi stage: ITER: 72289 RES: 1.0875e-05
Jacobi stage: ITER: 72649 RES: 1.029e-05
Jacobi stage: ITER: 72891 RES: 9.9175e-06
Jacobi stage: ITER: 73259 RES: 9.3725e-06
Jacobi stage: ITER: 73625 RES: 8.8574e-06
Jacobi stage: ITER: 73845 RES: 8.5631e-06
Jacobi stage: ITER: 74219 RES: 8.0875e-06
Jacobi stage: ITER: 74595 RES: 7.6336e-06
Jacobi stage: ITER: 74819 RES: 7.3754e-06
Jacobi stage: ITER: 75191 RES: 6.9658e-06
Jacobi stage: ITER: 75491 RES: 6.6521e-06
Jacobi stage: ITER: 75795 RES: 6.3486e-06
Jacobi stage: ITER: 76171 RES: 5.9924e-06
Jacobi stage: ITER: 76391 RES: 5.7933e-06
Jacobi stage: ITER: 76759 RES: 5.4749e-06
Jacobi stage: ITER: 77129 RES: 5.1708e-06
Jacobi stage: ITER: 77371 RES: 4.9837e-06
Jacobi stage: ITER: 77707 RES: 4.733e-06
Jacobi stage: ITER: 78079 RES: 4.4701e-06
Jacobi stage: ITER: 78299 RES: 4.3216e-06
Jacobi stage: ITER: 78667 RES: 4.0841e-06
Jacobi stage: ITER: 78995 RES: 3.8834e-06
Jacobi stage: ITER: 79259 RES: 3.7291e-06
Jacobi stage: ITER: 79625 RES: 3.5242e-06
Jacobi stage: ITER: 79845 RES: 3.4071e-06
Jacobi stage: ITER: 80203 RES: 3.2257e-06
Jacobi stage: ITER: 80565 RES: 3.0503e-06
Jacobi stage: ITER: 80779 RES: 2.9526e-06
Jacobi stage: ITER: 81143 RES: 2.7921e-06
Jacobi stage: ITER: 81427 RES: 2.6729e-06
Jacobi stage: ITER: 81733 RES: 2.5494e-06
Jacobi stage: ITER: 82091 RES: 2.4137e-06
Jacobi stage: ITER: 82311 RES: 2.3335e-06
Jacobi stage: ITER: 82677 RES: 2.2053e-06
Jacobi stage: ITER: 82949 RES: 2.115e-06
Jacobi stage: ITER: 83253 RES: 2.0185e-06
Jacobi stage: ITER: 83613 RES: 1.9099e-06
Jacobi stage: ITER: 83821 RES: 1.8499e-06
Jacobi stage: ITER: 84183 RES: 1.7504e-06
Jacobi stage: ITER: 84481 RES: 1.6715e-06
Jacobi stage: ITER: 84777 RES: 1.5973e-06
Jacobi stage: ITER: 85147 RES: 1.5095e-06
Jacobi stage: ITER: 85371 RES: 1.4584e-06
Jacobi stage: ITER: 85729 RES: 1.38e-06
Jacobi stage: ITER: 86089 RES: 1.3057e-06
Jacobi stage: ITER: 86305 RES: 1.2631e-06
Jacobi stage: ITER: 86673 RES: 1.1937e-06
Jacobi stage: ITER: 86891 RES: 1.1547e-06
Jacobi stage: ITER: 87263 RES: 1.0906e-06
Jacobi stage: ITER: 87633 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 <TNL/Timer.h>
4#include <TNL/Matrices/SparseMatrix.h>
5#include <TNL/Devices/Sequential.h>
6#include <TNL/Devices/Cuda.h>
7#include <TNL/Solvers/Linear/Jacobi.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 * Setup solver of the linear system.
64 */
66 LinearSolver solver;
67 solver.setMatrix( matrix_ptr );
68 solver.setOmega( 0.0005 );
69
70 /***
71 * Setup monitor of the iterative solver.
72 */
74 IterativeSolverMonitorType monitor;
75 TNL::Solvers::SolverMonitorThread mmonitorThread( monitor );
76 monitor.setRefreshRate( 10 ); // refresh rate in milliseconds
77 monitor.setVerbose( 1 );
78 monitor.setStage( "Jacobi stage:" );
79 TNL::Timer timer;
80 monitor.setTimer( timer );
81 timer.start();
82 solver.setSolverMonitor( monitor );
83 solver.setConvergenceResidue( 1.0e-6 );
84 solver.solve( b, x );
85 monitor.stopMainLoop();
87}
88
89int
90main( int argc, char* argv[] )
91{
93 iterativeLinearSolverExample< TNL::Devices::Sequential >();
94
95#ifdef __CUDACC__
97 iterativeLinearSolverExample< TNL::Devices::Cuda >();
98#endif
99}
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:3.4931e-
ELA: 0.10014 Jacobi stage: ITER: 48787 RES: 0.00040211
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.8561e-
ELA: 0.10012 Jacobi stage: ITER: 567 RES: 0.83409
ELA: 0.20023 Jacobi stage: ITER: 1153 RES: 0.70926
ELA: 0.30031 Jacobi stage: ITER: 1491 RES: 0.65294
ELA: 0.40045 Jacobi stage: ITER: 1931 RES: 0.59097
ELA: 0.50059 Jacobi stage: ITER: 2249 RES: 0.55239
ELA: 0.60067 Jacobi stage: ITER: 2543 RES: 0.52095
ELA: 0.70076 Jacobi stage: ITER: 2931 RES: 0.48367
ELA: 0.80085 Jacobi stage: ITER: 3207 RES: 0.45974
ELA: 0.90097 Jacobi stage: ITER: 3559 RES: 0.43179
ELA: 1.0011 Jacobi stage: ITER: 3939 RES: 0.40435
ELA: 1.1012 Jacobi stage: ITER: 4169 RES: 0.38878
ELA: 1.2013 Jacobi stage: ITER: 4543 RES: 0.36546
ELA: 1.3014 Jacobi stage: ITER: 4837 RES: 0.34817
ELA: 1.4015 Jacobi stage: ITER: 5145 RES: 0.33124
ELA: 1.5016 Jacobi stage: ITER: 5523 RES: 0.31189
ELA: 1.6017 Jacobi stage: ITER: 5739 RES: 0.30136
ELA: 1.7018 Jacobi stage: ITER: 6109 RES: 0.28417
ELA: 1.8018 Jacobi stage: ITER: 6423 RES: 0.27058
ELA: 1.9019 Jacobi stage: ITER: 6707 RES: 0.25883
ELA: 2.002 Jacobi stage: ITER: 7093 RES: 0.24364
ELA: 2.1021 Jacobi stage: ITER: 7455 RES: 0.23039
ELA: 2.2022 Jacobi stage: ITER: 7731 RES: 0.22075
ELA: 2.3023 Jacobi stage: ITER: 8117 RES: 0.20789
ELA: 2.4024 Jacobi stage: ITER: 8351 RES: 0.20057
ELA: 2.5025 Jacobi stage: ITER: 8707 RES: 0.18985
ELA: 2.6026 Jacobi stage: ITER: 9091 RES: 0.17894
ELA: 2.7027 Jacobi stage: ITER: 9315 RES: 0.17287
ELA: 2.8028 Jacobi stage: ITER: 9561 RES: 0.16639
ELA: 2.9029 Jacobi stage: ITER: 9715 RES: 0.16254
ELA: 3.003 Jacobi stage: ITER: 9927 RES: 0.15732
ELA: 3.1032 Jacobi stage: ITER: 10159 RES: 0.15181
ELA: 3.2033 Jacobi stage: ITER: 10295 RES: 0.14867
ELA: 3.3033 Jacobi stage: ITER: 10539 RES: 0.14319
ELA: 3.4034 Jacobi stage: ITER: 10789 RES: 0.13775
ELA: 3.5035 Jacobi stage: ITER: 10937 RES: 0.13465
ELA: 3.6036 Jacobi stage: ITER: 11171 RES: 0.12993
ELA: 3.7037 Jacobi stage: ITER: 11357 RES: 0.12623
ELA: 3.8038 Jacobi stage: ITER: 11565 RES: 0.12226
ELA: 3.9039 Jacobi stage: ITER: 11861 RES: 0.11682
ELA: 4.004 Jacobi stage: ITER: 12077 RES: 0.11301
ELA: 4.1041 Jacobi stage: ITER: 12459 RES: 0.1066
ELA: 4.2042 Jacobi stage: ITER: 12777 RES: 0.10148
ELA: 4.3043 Jacobi stage: ITER: 13061 RES: 0.097151
ELA: 4.4044 Jacobi stage: ITER: 13447 RES: 0.091586
ELA: 4.5044 Jacobi stage: ITER: 13669 RES: 0.088488
ELA: 4.6045 Jacobi stage: ITER: 14049 RES: 0.08347
ELA: 4.7046 Jacobi stage: ITER: 14427 RES: 0.078785
ELA: 4.8047 Jacobi stage: ITER: 14653 RES: 0.076074
ELA: 4.9048 Jacobi stage: ITER: 15033 RES: 0.07176
ELA: 5.0049 Jacobi stage: ITER: 15299 RES: 0.068909
ELA: 5.105 Jacobi stage: ITER: 15631 RES: 0.065482
ELA: 5.2051 Jacobi stage: ITER: 16011 RES: 0.06177
ELA: 5.3052 Jacobi stage: ITER: 16237 RES: 0.059644
ELA: 5.4053 Jacobi stage: ITER: 16631 RES: 0.056158
ELA: 5.5054 Jacobi stage: ITER: 16975 RES: 0.053268
ELA: 5.6055 Jacobi stage: ITER: 17231 RES: 0.051214
ELA: 5.7056 Jacobi stage: ITER: 17615 RES: 0.048281
ELA: 5.8057 Jacobi stage: ITER: 17837 RES: 0.046648
ELA: 5.9058 Jacobi stage: ITER: 18215 RES: 0.04403
ELA: 6.0059 Jacobi stage: ITER: 18603 RES: 0.041483
ELA: 6.106 Jacobi stage: ITER: 18851 RES: 0.039932
ELA: 6.2061 Jacobi stage: ITER: 19235 RES: 0.037645
ELA: 6.3062 Jacobi stage: ITER: 19563 RES: 0.035795
ELA: 6.4063 Jacobi stage: ITER: 19843 RES: 0.034288
ELA: 6.5064 Jacobi stage: ITER: 20215 RES: 0.032384
ELA: 6.6065 Jacobi stage: ITER: 20435 RES: 0.031308
ELA: 6.7066 Jacobi stage: ITER: 20791 RES: 0.029642
ELA: 6.8067 Jacobi stage: ITER: 21117 RES: 0.028186
ELA: 6.9068 Jacobi stage: ITER: 21387 RES: 0.027049
ELA: 7.0068 Jacobi stage: ITER: 21773 RES: 0.025484
ELA: 7.1069 Jacobi stage: ITER: 22003 RES: 0.024607
ELA: 7.207 Jacobi stage: ITER: 22391 RES: 0.023183
ELA: 7.3071 Jacobi stage: ITER: 22775 RES: 0.021855
ELA: 7.4072 Jacobi stage: ITER: 23023 RES: 0.021038
ELA: 7.5074 Jacobi stage: ITER: 23391 RES: 0.019882
ELA: 7.6074 Jacobi stage: ITER: 23685 RES: 0.018999
ELA: 7.7075 Jacobi stage: ITER: 23999 RES: 0.018109
ELA: 7.8076 Jacobi stage: ITER: 24381 RES: 0.017072
ELA: 7.9077 Jacobi stage: ITER: 24603 RES: 0.016505
ELA: 8.0078 Jacobi stage: ITER: 24979 RES: 0.015579
ELA: 8.1079 Jacobi stage: ITER: 25281 RES: 0.014868
ELA: 8.208 Jacobi stage: ITER: 25561 RES: 0.014242
ELA: 8.3081 Jacobi stage: ITER: 25937 RES: 0.013443
ELA: 8.4082 Jacobi stage: ITER: 26159 RES: 0.012996
ELA: 8.5083 Jacobi stage: ITER: 26539 RES: 0.012259
ELA: 8.6084 Jacobi stage: ITER: 26919 RES: 0.011564
ELA: 8.7085 Jacobi stage: ITER: 27141 RES: 0.011173
ELA: 8.8085 Jacobi stage: ITER: 27517 RES: 0.010546
ELA: 8.9086 Jacobi stage: ITER: 27757 RES: 0.010164
ELA: 9.0087 Jacobi stage: ITER: 28135 RES: 0.009594
ELA: 9.1088 Jacobi stage: ITER: 28511 RES: 0.0090556
ELA: 9.2089 Jacobi stage: ITER: 28725 RES: 0.0087601
ELA: 9.309 Jacobi stage: ITER: 29075 RES: 0.0083041
ELA: 9.4091 Jacobi stage: ITER: 29347 RES: 0.0079643
ELA: 9.5092 Jacobi stage: ITER: 29653 RES: 0.0075963
ELA: 9.6093 Jacobi stage: ITER: 30033 RES: 0.0071656
ELA: 9.7094 Jacobi stage: ITER: 30259 RES: 0.0069233
ELA: 9.8095 Jacobi stage: ITER: 30639 RES: 0.0065307
ELA: 9.9095 Jacobi stage: ITER: 31025 RES: 0.0061529
ELA: 10.01 Jacobi stage: ITER: 31255 RES: 0.0059412
ELA: 10.11 Jacobi stage: ITER: 31651 RES: 0.0055905
ELA: 10.21 Jacobi stage: ITER: 31891 RES: 0.0053882
ELA: 10.31 Jacobi stage: ITER: 32259 RES: 0.0050921
ELA: 10.41 Jacobi stage: ITER: 32639 RES: 0.0048034
ELA: 10.51 Jacobi stage: ITER: 32867 RES: 0.0046381
ELA: 10.61 Jacobi stage: ITER: 33247 RES: 0.0043751
ELA: 10.71 Jacobi stage: ITER: 33607 RES: 0.0041397
ELA: 10.811 Jacobi stage: ITER: 33857 RES: 0.0039826
ELA: 10.911 Jacobi stage: ITER: 34219 RES: 0.0037683
ELA: 11.011 Jacobi stage: ITER: 34449 RES: 0.0036364
ELA: 11.111 Jacobi stage: ITER: 34823 RES: 0.0034344
ELA: 11.211 Jacobi stage: ITER: 35211 RES: 0.0032337
ELA: 11.311 Jacobi stage: ITER: 35455 RES: 0.0031167
ELA: 11.411 Jacobi stage: ITER: 35847 RES: 0.0029346
ELA: 11.511 Jacobi stage: ITER: 36233 RES: 0.0027648
ELA: 11.612 Jacobi stage: ITER: 36457 RES: 0.0026713
ELA: 11.712 Jacobi stage: ITER: 36847 RES: 0.0025167
ELA: 11.812 Jacobi stage: ITER: 37095 RES: 0.0024227
ELA: 11.912 Jacobi stage: ITER: 37455 RES: 0.0022923
ELA: 12.012 Jacobi stage: ITER: 37833 RES: 0.0021624
ELA: 12.112 Jacobi stage: ITER: 38049 RES: 0.0020918
ELA: 12.212 Jacobi stage: ITER: 38395 RES: 0.0019841
ELA: 12.312 Jacobi stage: ITER: 38711 RES: 0.0018901
ELA: 12.412 Jacobi stage: ITER: 38993 RES: 0.0018095
ELA: 12.513 Jacobi stage: ITER: 39353 RES: 0.0017121
ELA: 12.613 Jacobi stage: ITER: 39651 RES: 0.001636
ELA: 12.713 Jacobi stage: ITER: 39987 RES: 0.0015537
ELA: 12.813 Jacobi stage: ITER: 40351 RES: 0.0014692
ELA: 12.913 Jacobi stage: ITER: 40571 RES: 0.0014204
ELA: 13.013 Jacobi stage: ITER: 40927 RES: 0.0013448
ELA: 13.113 Jacobi stage: ITER: 41163 RES: 0.001297
ELA: 13.213 Jacobi stage: ITER: 41509 RES: 0.0012294
ELA: 13.313 Jacobi stage: ITER: 41871 RES: 0.0011633
ELA: 13.413 Jacobi stage: ITER: 42087 RES: 0.0011253
ELA: 13.513 Jacobi stage: ITER: 42477 RES: 0.0010596
ELA: 13.614 Jacobi stage: ITER: 42767 RES: 0.0010137
ELA: 13.714 Jacobi stage: ITER: 43063 RES: 0.00096868
ELA: 13.814 Jacobi stage: ITER: 43441 RES: 0.00091375
ELA: 13.914 Jacobi stage: ITER: 43753 RES: 0.000871
ELA: 14.014 Jacobi stage: ITER: 44051 RES: 0.00083228
ELA: 14.114 Jacobi stage: ITER: 44419 RES: 0.00078654
ELA: 14.214 Jacobi stage: ITER: 44649 RES: 0.00075901
ELA: 14.314 Jacobi stage: ITER: 45019 RES: 0.00071729
ELA: 14.414 Jacobi stage: ITER: 45391 RES: 0.00067746
ELA: 14.514 Jacobi stage: ITER: 45623 RES: 0.00065374
ELA: 14.614 Jacobi stage: ITER: 45991 RES: 0.00061781
ELA: 14.714 Jacobi stage: ITER: 46353 RES: 0.00058422
ELA: 14.815 Jacobi stage: ITER: 46583 RES: 0.00056411
ELA: 14.915 Jacobi stage: ITER: 46959 RES: 0.00053245
ELA: 15.015 Jacobi stage: ITER: 47265 RES: 0.00050785
ELA: 15.115 Jacobi stage: ITER: 47553 RES: 0.00048588
ELA: 15.215 Jacobi stage: ITER: 47927 RES: 0.00045889
ELA: 15.315 Jacobi stage: ITER: 48257 RES: 0.00043608
ELA: 15.415 Jacobi stage: ITER: 48539 RES: 0.00041772
ELA: 15.515 Jacobi stage: ITER: 48925 RES: 0.00039355
ELA: 15.615 Jacobi stage: ITER: 49145 RES: 0.00038047
ELA: 15.715 Jacobi stage: ITER: 49513 RES: 0.00035956
ELA: 15.816 Jacobi stage: ITER: 49885 RES: 0.00033959
ELA: 15.916 Jacobi stage: ITER: 50103 RES: 0.00032851
ELA: 16.016 Jacobi stage: ITER: 50475 RES: 0.00031027
ELA: 16.116 Jacobi stage: ITER: 50775 RES: 0.00029629
ELA: 16.216 Jacobi stage: ITER: 51051 RES: 0.000284
ELA: 16.316 Jacobi stage: ITER: 51415 RES: 0.00026855
ELA: 16.416 Jacobi stage: ITER: 51631 RES: 0.00025979
ELA: 16.516 Jacobi stage: ITER: 51993 RES: 0.00024566
ELA: 16.616 Jacobi stage: ITER: 52359 RES: 0.0002323
ELA: 16.716 Jacobi stage: ITER: 52587 RES: 0.00022431
ELA: 16.816 Jacobi stage: ITER: 52951 RES: 0.00021211
ELA: 16.917 Jacobi stage: ITER: 53171 RES: 0.00020506
ELA: 17.017 Jacobi stage: ITER: 53533 RES: 0.00019391
ELA: 17.117 Jacobi stage: ITER: 53899 RES: 0.00018337
ELA: 17.217 Jacobi stage: ITER: 54115 RES: 0.00017739
ELA: 17.317 Jacobi stage: ITER: 54471 RES: 0.00016795
ELA: 17.417 Jacobi stage: ITER: 54695 RES: 0.00016227
ELA: 17.517 Jacobi stage: ITER: 55049 RES: 0.00015363
ELA: 17.617 Jacobi stage: ITER: 55423 RES: 0.0001451
ELA: 17.717 Jacobi stage: ITER: 55639 RES: 0.00014036
ELA: 17.817 Jacobi stage: ITER: 56013 RES: 0.00013249
ELA: 17.917 Jacobi stage: ITER: 56321 RES: 0.00012636
ELA: 18.018 Jacobi stage: ITER: 56599 RES: 0.00012112
ELA: 18.118 Jacobi stage: ITER: 56963 RES: 0.00011453
ELA: 18.218 Jacobi stage: ITER: 57179 RES: 0.0001108
ELA: 18.318 Jacobi stage: ITER: 57543 RES: 0.00010477
ELA: 18.418 Jacobi stage: ITER: 57855 RES: 9.9869e-05
ELA: 18.518 Jacobi stage: ITER: 58135 RES: 9.5664e-05
ELA: 18.618 Jacobi stage: ITER: 58503 RES: 9.0407e-05
ELA: 18.718 Jacobi stage: ITER: 58723 RES: 8.7403e-05
ELA: 18.818 Jacobi stage: ITER: 59255 RES: 8.0545e-05
ELA: 18.918 Jacobi stage: ITER: 59787 RES: 7.4225e-05
ELA: 19.018 Jacobi stage: ITER: 60171 RES: 6.9973e-05
ELA: 19.119 Jacobi stage: ITER: 60731 RES: 6.4206e-05
ELA: 19.219 Jacobi stage: ITER: 61079 RES: 6.0864e-05
ELA: 19.319 Jacobi stage: ITER: 61635 RES: 5.5882e-05
ELA: 19.419 Jacobi stage: ITER: 62207 RES: 5.1182e-05
ELA: 19.519 Jacobi stage: ITER: 62545 RES: 4.8577e-05
ELA: 19.619 Jacobi stage: ITER: 63125 RES: 4.4437e-05
ELA: 19.719 Jacobi stage: ITER: 63463 RES: 4.2202e-05
ELA: 19.819 Jacobi stage: ITER: 64009 RES: 3.8795e-05
ELA: 19.919 Jacobi stage: ITER: 64565 RES: 3.5619e-05
ELA: 20.02 Jacobi stage: ITER: 64889 RES: 3.389e-05
ELA: 20.12 Jacobi stage: ITER: 65443 RES: 3.1135e-05
ELA: 20.22 Jacobi stage: ITER: 65789 RES: 2.9514e-05
ELA: 20.32 Jacobi stage: ITER: 66341 RES: 2.7115e-05
ELA: 20.42 Jacobi stage: ITER: 66903 RES: 2.488e-05
ELA: 20.52 Jacobi stage: ITER: 67251 RES: 2.3585e-05
ELA: 20.62 Jacobi stage: ITER: 67815 RES: 2.1628e-05
ELA: 20.72 Jacobi stage: ITER: 68179 RES: 2.0452e-05
ELA: 20.82 Jacobi stage: ITER: 68701 RES: 1.887e-05
ELA: 20.92 Jacobi stage: ITER: 69241 RES: 1.7368e-05
ELA: 21.02 Jacobi stage: ITER: 69567 RES: 1.6525e-05
ELA: 21.121 Jacobi stage: ITER: 70117 RES: 1.5182e-05
ELA: 21.221 Jacobi stage: ITER: 70477 RES: 1.4365e-05
ELA: 21.321 Jacobi stage: ITER: 71015 RES: 1.323e-05
ELA: 21.421 Jacobi stage: ITER: 71563 RES: 1.2162e-05
ELA: 21.521 Jacobi stage: ITER: 71883 RES: 1.1578e-05
ELA: 21.621 Jacobi stage: ITER: 72439 RES: 1.0631e-05
ELA: 21.721 Jacobi stage: ITER: 72787 RES: 1.0077e-05
ELA: 21.821 Jacobi stage: ITER: 73343 RES: 9.2523e-06
ELA: 21.921 Jacobi stage: ITER: 73889 RES: 8.5054e-06
ELA: 22.021 Jacobi stage: ITER: 74217 RES: 8.0875e-06
ELA: 22.122 Jacobi stage: ITER: 74777 RES: 7.4209e-06
ELA: 22.222 Jacobi stage: ITER: 75099 RES: 7.065e-06
ELA: 22.322 Jacobi stage: ITER: 75665 RES: 6.4747e-06
ELA: 22.422 Jacobi stage: ITER: 76211 RES: 5.9557e-06
ELA: 22.522 Jacobi stage: ITER: 76573 RES: 5.6318e-06
ELA: 22.622 Jacobi stage: ITER: 77131 RES: 5.1708e-06
ELA: 22.722 Jacobi stage: ITER: 77459 RES: 4.9168e-06
ELA: 22.822 Jacobi stage: ITER: 78021 RES: 4.5087e-06
ELA: 22.922 Jacobi stage: ITER: 78577 RES: 4.1397e-06
ELA: 23.023 Jacobi stage: ITER: 78919 RES: 3.929e-06
ELA: 23.123 Jacobi stage: ITER: 79483 RES: 3.603e-06
ELA: 23.223 Jacobi stage: ITER: 79815 RES: 3.4238e-06
ELA: 23.323 Jacobi stage: ITER: 80387 RES: 3.1359e-06
ELA: 23.423 Jacobi stage: ITER: 80849 RES: 2.9201e-06
ELA: 23.523 Jacobi stage: ITER: 81279 RES: 2.7343e-06
ELA: 23.623 Jacobi stage: ITER: 81833 RES: 2.5105e-06
ELA: 23.723 Jacobi stage: ITER: 82171 RES: 2.3842e-06
ELA: 23.823 Jacobi stage: ITER: 82721 RES: 2.1904e-06
ELA: 23.923 Jacobi stage: ITER: 83131 RES: 2.0574e-06
ELA: 24.023 Jacobi stage: ITER: 83611 RES: 1.9111e-06
ELA: 24.124 Jacobi stage: ITER: 84153 RES: 1.7579e-06
ELA: 24.224 Jacobi stage: ITER: 84527 RES: 1.6603e-06
ELA: 24.324 Jacobi stage: ITER: 85087 RES: 1.5234e-06
ELA: 24.424 Jacobi stage: ITER: 85651 RES: 1.397e-06
ELA: 24.524 Jacobi stage: ITER: 85993 RES: 1.3251e-06
ELA: 24.624 Jacobi stage: ITER: 86563 RES: 1.2144e-06
ELA: 24.724 Jacobi stage: ITER: 86893 RES: 1.154e-06
ELA: 24.824 Jacobi stage: ITER: 87455 RES: 1.0589e-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.16.1