Template Numerical Library version\ main:f17d0c8
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
TNL::Solvers::IterativeSolverMonitor< Real, Index > Class Template Reference

Object for monitoring convergence of iterative solvers. More...

#include <TNL/Solvers/IterativeSolverMonitor.h>

Inheritance diagram for TNL::Solvers::IterativeSolverMonitor< Real, Index >:
Inheritance graph
[legend]
Collaboration diagram for TNL::Solvers::IterativeSolverMonitor< Real, Index >:
Collaboration graph
[legend]

Public Types

using IndexType = Index
 A type for indexing.
 
using RealType = Real
 A type of the floating-point arithmetics.
 

Public Member Functions

 IterativeSolverMonitor ()=default
 Construct with no parameters.
 
void refresh () override
 Causes that the monitor prints out the status of the solver.
 
void setIterations (const IndexType &iterations)
 Set number of the current iteration.
 
void setNodesPerIteration (const IndexType &nodes)
 Set the number of nodes of the numerical mesh or lattice.
 
void setResidue (const RealType &residue)
 Set residue of the current approximation of the solution.
 
void setStage (const std::string &stage)
 This method can be used for naming a stage of the monitored solver.
 
void setTime (const RealType &time)
 Set the time of the simulated evolution if it is time dependent.
 
void setTimeStep (const RealType &timeStep)
 Set the time step for time dependent iterative solvers.
 
void setVerbose (const IndexType &verbose)
 Set up the verbosity of the monitor.
 
- Public Member Functions inherited from TNL::Solvers::SolverMonitor
 SolverMonitor ()=default
 Basic construct with no arguments.
 
bool isStopped () const
 Checks whether the main loop was stopped.
 
void runMainLoop ()
 Starts the main loop from which the method SolverMonitor::refresh is called in given time periods.
 
void setRefreshRate (const int &refreshRate)
 Set the time interval between two consecutive calls of SolverMonitor::refresh.
 
void setTimer (Timer &timer)
 Set a timer object for the solver monitor.
 
void stopMainLoop ()
 Stops the main loop of the monitor. See runMainLoop.
 

Protected Member Functions

int getLineWidth ()
 
- Protected Member Functions inherited from TNL::Solvers::SolverMonitor
double getElapsedTime ()
 

Protected Attributes

std::atomic_bool attributes_changed { false }
 
RealType elapsed_time_before_refresh = 0
 
IndexType iterations = 0
 
IndexType iterations_before_refresh = 0
 
RealType last_mlups = 0
 
IndexType nodesPerIteration = 0
 
RealType residue = 0
 
std::atomic_bool saved { false }
 
IndexType saved_iterations = 0
 
RealType saved_residue = 0
 
std::string saved_stage
 
RealType saved_time = 0
 
RealType saved_timeStep = 0
 
std::string stage
 
RealType time = 0
 
RealType timeStep = 0
 
IndexType verbose = 2
 
- Protected Attributes inherited from TNL::Solvers::SolverMonitor
std::atomic_bool started { false }
 
std::atomic_bool stopped { false }
 
std::atomic_int timeout_milliseconds { 500 }
 
Timertimer = nullptr
 

Detailed Description

template<typename Real = double, typename Index = int>
class TNL::Solvers::IterativeSolverMonitor< Real, Index >

Object for monitoring convergence of iterative solvers.

Template Parameters
Realis a type of the floating-point arithmetics.
Indexis an indexing type.

The following example shows how to use the iterative solver monitor for monitoring convergence of linear iterative solver:

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 );
26 auto matrix_ptr = std::make_shared< MatrixType >();
27 matrix_ptr->setDimensions( size, size );
28 matrix_ptr->setRowCapacities( Vector( { 2, 3, 3, 3, 2 } ) );
29
30 auto f = [ = ] __cuda_callable__( typename MatrixType::RowView & row ) mutable
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 );
52 std::cout << *matrix_ptr << std::endl;
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;
61 std::cout << "Vector b = " << b << std::endl;
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 */
74 using IterativeSolverMonitorType = TNL::Solvers::IterativeSolverMonitor< double, int >;
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();
84 std::cout << "Vector x = " << x << std::endl;
85}
86
87int
88main( int argc, char* argv[] )
89{
90 std::cout << "Solving linear system on host: " << std::endl;
91 iterativeLinearSolverExample< TNL::Devices::Sequential >();
92
93#ifdef __CUDACC__
94 std::cout << "Solving linear system on CUDA device: " << std::endl;
95 iterativeLinearSolverExample< TNL::Devices::Cuda >();
96#endif
97}
#define __cuda_callable__
Definition Macros.h:49
Vector extends Array with algebraic operations.
Definition Vector.h:36
Implementation of sparse matrix, i.e. matrix storing only non-zero elements.
Definition SparseMatrix.h:57
Iterative solver of linear systems based on the Jacobi method.
Definition Jacobi.h:21
void setMatrix(const MatrixPointer &matrix)
Set the matrix of the linear system.
Definition LinearSolver.h:120
A RAII wrapper for launching the SolverMonitor's main loop in a separate thread.
Definition SolverMonitor.h:133
T endl(T... args)

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: 15987 RES: 0.061998
Jacobi stage: ITER: 38391 RES: 0.0019854
Jacobi stage: ITER: 54651 RES: 0.00016337
Jacobi stage: ITER: 75287 RES: 6.8639e-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: 4481 RES: 0.36909
Jacobi stage: ITER: 9113 RES: 0.17828
Jacobi stage: ITER: 13487 RES: 0.091025
Jacobi stage: ITER: 17889 RES: 0.046277
Jacobi stage: ITER: 20843 RES: 0.029406
Jacobi stage: ITER: 22867 RES: 0.021549
Jacobi stage: ITER: 27019 RES: 0.011388
Jacobi stage: ITER: 30363 RES: 0.0068136
Jacobi stage: ITER: 32649 RES: 0.0047945
Jacobi stage: ITER: 37107 RES: 0.0024182
Jacobi stage: ITER: 41563 RES: 0.0012197
Jacobi stage: ITER: 45969 RES: 0.00061971
Jacobi stage: ITER: 50447 RES: 0.0003116
Jacobi stage: ITER: 54871 RES: 0.00015794
Jacobi stage: ITER: 55735 RES: 0.00013831
Jacobi stage: ITER: 57135 RES: 0.00011155
Jacobi stage: ITER: 59075 RES: 8.2803e-05
Jacobi stage: ITER: 60795 RES: 6.3578e-05
Jacobi stage: ITER: 61291 RES: 5.8914e-05
Jacobi stage: ITER: 63599 RES: 4.1329e-05
Jacobi stage: ITER: 66387 RES: 2.6932e-05
Jacobi stage: ITER: 68423 RES: 1.9699e-05
Jacobi stage: ITER: 69113 RES: 1.7713e-05
Jacobi stage: ITER: 70263 RES: 1.4849e-05
Jacobi stage: ITER: 71015 RES: 1.323e-05
Jacobi stage: ITER: 71955 RES: 1.1451e-05
Jacobi stage: ITER: 73947 RES: 8.4325e-06
Jacobi stage: ITER: 74871 RES: 7.3168e-06
Jacobi stage: ITER: 75747 RES: 6.3956e-06
Jacobi stage: ITER: 76719 RES: 5.5086e-06
Jacobi stage: ITER: 77715 RES: 4.7272e-06
Jacobi stage: ITER: 78611 RES: 4.1194e-06
Jacobi stage: ITER: 80419 RES: 3.1205e-06
Jacobi stage: ITER: 82299 RES: 2.3378e-06
Jacobi stage: ITER: 83003 RES: 2.0982e-06
Jacobi stage: ITER: 84197 RES: 1.7461e-06
Jacobi stage: ITER: 85135 RES: 1.5123e-06
Jacobi stage: ITER: 85887 RES: 1.3473e-06
Jacobi stage: ITER: 86951 RES: 1.1441e-06
Vector x = [ 0.999999, 0.999999, 0.999998, 0.999999, 0.999999 ]

The following example shows how to employ timer (TNL::Timer) to the monitor of iterative solvers:

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 );
27 auto matrix_ptr = std::make_shared< MatrixType >();
28 matrix_ptr->setDimensions( size, size );
29 matrix_ptr->setRowCapacities( Vector( { 2, 3, 3, 3, 2 } ) );
30
31 auto f = [ = ] __cuda_callable__( typename MatrixType::RowView & row ) mutable
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 );
53 std::cout << *matrix_ptr << std::endl;
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;
62 std::cout << "Vector b = " << b << std::endl;
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 */
75 using IterativeSolverMonitorType = TNL::Solvers::IterativeSolverMonitor< double, int >;
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();
88 std::cout << "Vector x = " << x << std::endl;
89}
90
91int
92main( int argc, char* argv[] )
93{
94 std::cout << "Solving linear system on host: " << std::endl;
95 iterativeLinearSolverExample< TNL::Devices::Sequential >();
96
97#ifdef __CUDACC__
98 std::cout << "Solving linear system on CUDA device: " << std::endl;
99 iterativeLinearSolverExample< TNL::Devices::Cuda >();
100#endif
101}
Class for real time, CPU time and CPU cycles measuring.
Definition Timer.h:25
void start()
Starts timer.
Definition Timer.hpp:42

The result looks as follows:

Solving linear system on host:
Row: 0 -> 0:2.5 1:-1
Row: 1 -> 0:-1 1:2.5 2:-1
Row: 2 -> 1:-1 2:2.5 3:-1
Row: 3 -> 2:-1 3:2.5 4:-1
Row: 4 -> 3:-1 4:2.5
Vector b = [ 1.5, 0.5, 0.5, 0.5, 1.5 ]
ELA:0.001821
ELA: 0.10194 Jacobi stage: ITER: 10105 RES: 0.15303
ELA: 0.20204 Jacobi stage: ITER: 20607 RES: 0.030473
ELA: 0.30215 Jacobi stage: ITER: 32975 RES: 0.0045618
ELA: 0.40225 Jacobi stage: ITER: 53981 RES: 0.00018102
ELA: 0.50235 Jacobi stage: ITER: 75075 RES: 7.0911e-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: 0.01128
ELA: 0.1114 Jacobi stage: ITER: 2675 RES: 0.50774
ELA: 0.2115 Jacobi stage: ITER: 4419 RES: 0.37302
ELA: 0.31161 Jacobi stage: ITER: 6645 RES: 0.26126
ELA: 0.41172 Jacobi stage: ITER: 7509 RES: 0.2284
ELA: 0.51182 Jacobi stage: ITER: 8989 RES: 0.18171
ELA: 0.61192 Jacobi stage: ITER: 10035 RES: 0.15473
ELA: 0.71202 Jacobi stage: ITER: 10455 RES: 0.14505
ELA: 0.81213 Jacobi stage: ITER: 11223 RES: 0.1289
ELA: 0.91223 Jacobi stage: ITER: 12875 RES: 0.099998
ELA: 1.0123 Jacobi stage: ITER: 14695 RES: 0.075608
ELA: 1.1124 Jacobi stage: ITER: 15723 RES: 0.064564
ELA: 1.2125 Jacobi stage: ITER: 16639 RES: 0.056089
ELA: 1.3126 Jacobi stage: ITER: 17433 RES: 0.049634
ELA: 1.4127 Jacobi stage: ITER: 18371 RES: 0.042987
ELA: 1.5128 Jacobi stage: ITER: 19617 RES: 0.035489
ELA: 1.6129 Jacobi stage: ITER: 20777 RES: 0.029697
ELA: 1.713 Jacobi stage: ITER: 21827 RES: 0.025281
ELA: 1.8131 Jacobi stage: ITER: 23019 RES: 0.021051
ELA: 1.9132 Jacobi stage: ITER: 23851 RES: 0.018526
ELA: 2.0133 Jacobi stage: ITER: 26261 RES: 0.01279
ELA: 2.1134 Jacobi stage: ITER: 27637 RES: 0.010354
ELA: 2.2135 Jacobi stage: ITER: 28595 RES: 0.0089395
ELA: 2.3136 Jacobi stage: ITER: 30341 RES: 0.0068345
ELA: 2.4137 Jacobi stage: ITER: 31239 RES: 0.0059558
ELA: 2.5138 Jacobi stage: ITER: 32177 RES: 0.005155
ELA: 2.6155 Jacobi stage: ITER: 35009 RES: 0.0033367
ELA: 2.7156 Jacobi stage: ITER: 39851 RES: 0.0015865
ELA: 2.8157 Jacobi stage: ITER: 45127 RES: 0.00070549
ELA: 2.9158 Jacobi stage: ITER: 50559 RES: 0.00030629
ELA: 3.0159 Jacobi stage: ITER: 54849 RES: 0.00015842
ELA: 3.116 Jacobi stage: ITER: 59831 RES: 7.3725e-05
ELA: 3.2161 Jacobi stage: ITER: 64083 RES: 3.8368e-05
ELA: 3.3162 Jacobi stage: ITER: 68629 RES: 1.908e-05
ELA: 3.4163 Jacobi stage: ITER: 72835 RES: 1.0003e-05
ELA: 3.5164 Jacobi stage: ITER: 77689 RES: 4.7446e-06
ELA: 3.6222 Jacobi stage: ITER: 82715 RES: 2.1931e-06
ELA:Vector x = [ 3.7223 Jacobi stage: ITER:0.999999 , 87829 RES: 9.9949e-07
0.999999, 0.999998, 0.999999, 0.999999 ]

Member Function Documentation

◆ refresh()

template<typename Real , typename Index >
void TNL::Solvers::IterativeSolverMonitor< Real, Index >::refresh ( )
overridevirtual

Causes that the monitor prints out the status of the solver.

Implements TNL::Solvers::SolverMonitor.

◆ setIterations()

template<typename Real , typename Index >
void TNL::Solvers::IterativeSolverMonitor< Real, Index >::setIterations ( const IndexType & iterations)

Set number of the current iteration.

Parameters
iterationsis number of the current iteration.

◆ setNodesPerIteration()

template<typename Real , typename Index >
void TNL::Solvers::IterativeSolverMonitor< Real, Index >::setNodesPerIteration ( const IndexType & nodes)

Set the number of nodes of the numerical mesh or lattice.

This can be used to compute the number of nodes processed per one second.

Parameters
nodesis number of nodes of the numerical mesh or lattice.

◆ setResidue()

template<typename Real = double, typename Index = int>
void TNL::Solvers::IterativeSolverMonitor< Real, Index >::setResidue ( const RealType & residue)

Set residue of the current approximation of the solution.

Parameters
residueis a residue of the current approximation of the solution.

◆ setStage()

template<typename Real , typename Index >
void TNL::Solvers::IterativeSolverMonitor< Real, Index >::setStage ( const std::string & stage)

This method can be used for naming a stage of the monitored solver.

The stage name can be used to differ between various stages of iterative solvers.

Parameters
stageis name of the solver stage.

◆ setTime()

template<typename Real , typename Index >
void TNL::Solvers::IterativeSolverMonitor< Real, Index >::setTime ( const RealType & time)

Set the time of the simulated evolution if it is time dependent.

This can be used for example when solving parabolic or hyperbolic PDEs.

Parameters
timetime of the simulated evolution.

◆ setTimeStep()

template<typename Real , typename Index >
void TNL::Solvers::IterativeSolverMonitor< Real, Index >::setTimeStep ( const RealType & timeStep)

Set the time step for time dependent iterative solvers.

Parameters
timeSteptime step of the time dependent iterative solver.

◆ setVerbose()

template<typename Real , typename Index >
void TNL::Solvers::IterativeSolverMonitor< Real, Index >::setVerbose ( const IndexType & verbose)

Set up the verbosity of the monitor.

Parameters
verboseis the new value of the verbosity of the monitor.

The documentation for this class was generated from the following files: