From 4b6e6a5334c6b4f1d037f2628788c80e46be9c24 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Tom=C3=A1=C5=A1=20Oberhuber?= <oberhuber.tomas@gmail.com>
Date: Fri, 5 Feb 2021 10:45:07 +0100
Subject: [PATCH] Added tutorial for symmetric sparse matrices.
---
.../Tutorials/Matrices/CMakeLists.txt | 5 +++
.../Matrices/SymmetricSparseMatrixExample.cpp | 35 +++++++++++++++++++
.../Matrices/SymmetricSparseMatrixExample.cu | 1 +
.../Tutorials/Matrices/tutorial_Matrices.md | 27 ++++++++++++++
src/TNL/Matrices/SparseMatrix.h | 4 ++-
5 files changed, 71 insertions(+), 1 deletion(-)
create mode 100644 Documentation/Tutorials/Matrices/SymmetricSparseMatrixExample.cpp
create mode 120000 Documentation/Tutorials/Matrices/SymmetricSparseMatrixExample.cu
diff --git a/Documentation/Tutorials/Matrices/CMakeLists.txt b/Documentation/Tutorials/Matrices/CMakeLists.txt
index 94e57ec13d..176ade630a 100644
--- a/Documentation/Tutorials/Matrices/CMakeLists.txt
+++ b/Documentation/Tutorials/Matrices/CMakeLists.txt
@@ -94,6 +94,10 @@ IF( BUILD_CUDA )
${TNL_DOCUMENTATION_OUTPUT_SNIPPETS_PATH}/SparseMatrixViewExample_setElement.out
OUTPUT SparseMatrixViewExample_setElement.out )
+ CUDA_ADD_EXECUTABLE( SymmetricSparseMatrixExample SymmetricSparseMatrixExample.cu )
+ ADD_CUSTOM_COMMAND( COMMAND SymmetricSparseMatrixExample >
+ ${TNL_DOCUMENTATION_OUTPUT_SNIPPETS_PATH}/SymmetricSparseMatrixExample.out
+ OUTPUT SymmetricSparseMatrixExample.out )
####
# THe following examples/benchmarks run for very long time
CUDA_ADD_EXECUTABLE( DenseMatrixSetup_Benchmark_cuda DenseMatrixSetup_Benchmark.cu )
@@ -128,6 +132,7 @@ ADD_CUSTOM_TARGET( TutorialsMatricesCuda ALL DEPENDS
SparseMatrixExample_forRows.out
SparseMatrixExample_rowsReduction_vectorProduct.out
SparseMatrixViewExample_setElement.out
+ SymmetricSparseMatrixExample.out
)
ELSE()
ADD_CUSTOM_TARGET( TutorialsMatrices ALL DEPENDS
diff --git a/Documentation/Tutorials/Matrices/SymmetricSparseMatrixExample.cpp b/Documentation/Tutorials/Matrices/SymmetricSparseMatrixExample.cpp
new file mode 100644
index 0000000000..27e7f5a202
--- /dev/null
+++ b/Documentation/Tutorials/Matrices/SymmetricSparseMatrixExample.cpp
@@ -0,0 +1,35 @@
+#include <iostream>
+#include <TNL/Matrices/SparseMatrix.h>
+#include <TNL/Devices/Host.h>
+
+
+template< typename Device >
+void symmetricSparseMatrixExample()
+{
+ TNL::Matrices::SparseMatrix< double, Device, int, TNL::Matrices::SymmetricMatrix > symmetricMatrix (
+ 5, // number of matrix rows
+ 5, // number of matrix columns
+ { // matrix elements definition
+ { 0, 0, 1.0 },
+ { 1, 0, 2.0 }, { 1, 1, 1.0 },
+ { 2, 0, 3.0 }, { 2, 2, 1.0 },
+ { 3, 0, 4.0 }, { 3, 3, 1.0 },
+ { 4, 0, 5.0 }, { 4, 4, 1.0 } } );
+
+ std::cout << "Symmetric sparse matrix: " << std::endl << symmetricMatrix << std::endl;
+
+ TNL::Containers::Vector< double, Device > inVector( 5, 1.0 ), outVector( 5, 0.0 );
+ symmetricMatrix.vectorProduct( inVector, outVector );
+ std::cout << "Product with vector " << inVector << " is " << outVector << std::endl << std::endl;
+}
+
+int main( int argc, char* argv[] )
+{
+ std::cout << "Creating matrix on CPU ... " << std::endl;
+ symmetricSparseMatrixExample< TNL::Devices::Host >();
+
+#ifdef HAVE_CUDA
+ std::cout << "Creating matrix on CUDA GPU ... " << std::endl;
+ symmetricSparseMatrixExample< TNL::Devices::Cuda >();
+#endif
+}
diff --git a/Documentation/Tutorials/Matrices/SymmetricSparseMatrixExample.cu b/Documentation/Tutorials/Matrices/SymmetricSparseMatrixExample.cu
new file mode 120000
index 0000000000..688efe399e
--- /dev/null
+++ b/Documentation/Tutorials/Matrices/SymmetricSparseMatrixExample.cu
@@ -0,0 +1 @@
+SymmetricSparseMatrixExample.cpp
\ No newline at end of file
diff --git a/Documentation/Tutorials/Matrices/tutorial_Matrices.md b/Documentation/Tutorials/Matrices/tutorial_Matrices.md
index b5720f5d40..5f3c01d737 100644
--- a/Documentation/Tutorials/Matrices/tutorial_Matrices.md
+++ b/Documentation/Tutorials/Matrices/tutorial_Matrices.md
@@ -630,6 +630,33 @@ would not make sense. If we pass through this test, the matrix element lies in t
\include SparseMatrixExample_forRows.out
+#### Symmetric sparse matrices
+
+For sparse [symmetric matrices](https://en.wikipedia.org/wiki/Symmetric_matrix), TNL offers a format storing only a half of the matrix elements. More precisely, ony the matrix diagonal and the elements bellow are stored in the memory. The matrix elements above the diagonal are deduced from those bellow. If such a symmetric format is used on GPU, atomic operations must be used in some matrix operations. For this reason, symmetric matrices are allowed only for when the matrix elements values are expressed with `float` and `double` type. An advantage of the symmetric formats is lower memory consumption. Since less data need to be transferred from the memory, better performance might be observed. In some cases, however, the use of atomic operations on GPU may cause performance drop. Mostly we can see approximately the same performance compared to general formats but we can profit from lower memory requirements which is appreciated especially on GPU. The following example shows how to create symmetric sparse matrix.
+
+\includelineno SymmetricSparseMatrixExample.cpp
+
+We construct matrix of the following form
+
+\f[
+\left(
+\begin{array}{ccccc}
+ 1 & \color{grey}{2} & \color{grey}{3} & \color{grey}{4} & \color{grey}{5} \\
+ 2 & 1 & & & \\
+ 3 & & 1 & & \\
+ 4 & & & 1 & \\
+ 5 & & & & 1
+\end{array}
+\right)
+\f]
+
+The elements depicted in grey color are not stored in the memory. The main difference, compared to creation of general sparse matrix, is on line 9 where we state that the matrix is symmetric by setting the matrix type to \ref TNL::Matrices::SymmetricMatrix. Next we set only the diagonal elements and those lying bellow the diagonal (lines 13-17). When we print the matrix (line 19) we can see also the symmetric part above the diagonal. Next we test product of matrix and vector (lines 21-23). The result looks as follows:
+
+\include SymmetricSparseMatrixExample.out
+
+**Warning: Assignment of symmetric sparse matrix to general sparse matrix does not give correct result, currently. Only the diagonal and the lower part of the matrix is assigned.**
+
+
### Tridiagonal matrices <a name="tridiagonal_matrices_setup"></a>
Tridiagonal matrix format serves for specific matrix pattern when the nonzero matrix elements can be placed only at the diagonal and immediately next to the diagonal. Here is an example:
diff --git a/src/TNL/Matrices/SparseMatrix.h b/src/TNL/Matrices/SparseMatrix.h
index 247ba8ef42..044c4698b2 100644
--- a/src/TNL/Matrices/SparseMatrix.h
+++ b/src/TNL/Matrices/SparseMatrix.h
@@ -903,7 +903,9 @@ class SparseMatrix : public Matrix< Real, Device, Index, RealAllocator >
/**
* \brief Assignment of any matrix type other then this and dense.
- * .
+ *
+ * **Warning: Assignment of symmetric sparse matrix to general sparse matrix does not give correct result, currently. Only the diagonal and the lower part of the matrix is assigned.**
+ *
* \param matrix is input matrix for the assignment.
* \return reference to this matrix.
*/
--
GitLab