Added test and fix of MultidiagionalMatrix::setElements. (9ab855e2) · Commits · TNL / tnl-dev

Documentation/Examples/Matrices/MultidiagonalMatrix/MultidiagonalMatrixExample_Constructor_init_list_2.cpp

+12 −29

Original line number	Diff line number	Diff line
		@@ -12,40 +12,23 @@ void laplaceOperatorMatrix()
		const int matrixSize = gridSize * gridSize;
		TNL::Matrices::MultidiagonalMatrix< double, Device > matrix(
		matrixSize, { - gridSize, -1, 0, 1, gridSize }, {
		{ 1.0 },
		{ 1.0 },
		{ 1.0 },
		{ 1.0 },
		{ 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ -1.0, -1.0, 4.0, -1.0, -1.0 },
		{ -1.0, -1.0, 4.0, -1.0, -1.0 },
		{ 1.0 },
		{ 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ -1.0, -1.0, 4.0, -1.0, -1.0 },
		{ -1.0, -1.0, 4.0, -1.0, -1.0 },
		{ 1.0 },
		{ 1.0 },
		{ 1.0 },
		{ 1.0 },
		{ 1.0 }
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 }
		} );
		auto matrixView = matrix.getView();
		auto f = [=] __cuda_callable__ ( int i, int j ) mutable {
		const int elementIdx = i * gridSize + j;
		auto row = matrixView.getRow( elementIdx );
		if( i == 0 \|\| j == 0 \|\| i == gridSize - 1 \|\| j == gridSize - 1 )
		row.setElement( 0, 1.0 );
		else
		{
		row.setElement( 0, -1.0 );
		row.setElement( 1, -1.0 );
		row.setElement( 2, 4.0 );
		row.setElement( 3, -1.0 );
		row.setElement( 4, -1.0 );
		}
		};
		TNL::Algorithms::ParallelFor2D< Device >::exec( 0, 0, gridSize, gridSize, f );

		std::cout << "Laplace operator matrix: " << std::endl << matrix << std::endl;
		}

src/TNL/Matrices/MultidiagonalMatrix.hpp

+4 −8

Original line number	Diff line number	Diff line
		@@ -208,18 +208,14 @@ setElements( const std::initializer_list< std::initializer_list< ListReal > >& d
		{
		if( std::is_same< DeviceType, Devices::Host >::value )
		{
		this->getValues() = 0.0;
		auto row_it = data.begin();
		for( size_t rowIdx = 0; rowIdx < data.size(); rowIdx++ )
		{
		auto data_it = row_it->begin();
		for( IndexType i = 0; i < this->diagonalsShifts.getSize(); i++ )
		{
		const auto columnIdx = this->diagonalsShifts[ i ] + rowIdx;
		if( columnIdx >= 0 && columnIdx < this->getColumns() && data_it != row_it->end() )
		this->getRow( rowIdx ).setElement( i, *data_it++ );
		else
		this->getRow( rowIdx ).setElement( i, 0 );
		}
		IndexType i = 0;
		while( data_it != row_it->end() )
		this->getRow( rowIdx ).setElement( i++, *data_it++ );
		row_it ++;
		}
		}

src/TNL/Matrices/MultidiagonalMatrixView.hpp

+3 −3

Original line number	Diff line number	Diff line
		@@ -320,9 +320,9 @@ getElement( const IndexType row, const IndexType column ) const
		TNL_ASSERT_GE( column, 0, "" );
		TNL_ASSERT_LT( column, this->getColumns(), "" );

		for( IndexType i = 0; i < hostDiagonalsShifts.getSize(); i++ )
		if( row + hostDiagonalsShifts[ i ] == column )
		return this->values.getElement( this->getElementIndex( row, i ) );
		for( IndexType localIdx = 0; localIdx < hostDiagonalsShifts.getSize(); localIdx++ )
		if( row + hostDiagonalsShifts[ localIdx ] == column )
		return this->values.getElement( this->indexer.getGlobalIndex( row, localIdx ) );
		return 0.0;
		}

src/UnitTests/Matrices/MultidiagonalMatrixTest.h

+69 −1

Original line number	Diff line number	Diff line
		@@ -89,6 +89,68 @@ void test_SetLike()
		EXPECT_EQ( m1.getColumns(), m2.getColumns() );
		}

		template< typename Matrix >
		void test_SetElements()
		{
		using RealType = typename Matrix::RealType;
		using DeviceType = typename Matrix::DeviceType;
		using IndexType = typename Matrix::IndexType;

		const int gridSize( 4 );
		const int matrixSize( gridSize * gridSize );
		Matrix matrix( matrixSize, matrixSize, { - gridSize, -1, 0, 1, gridSize } );
		matrix.setElements( {
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ -1.0, -1.0, 4.0, -1.0, -1.0 },
		{ -1.0, -1.0, 4.0, -1.0, -1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ -1.0, -1.0, 4.0, -1.0, -1.0 },
		{ -1.0, -1.0, 4.0, -1.0, -1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 },
		{ 0.0, 0.0, 1.0 }
		} );

		std::cerr << matrix << std::endl;
		for( int i = 0; i < gridSize; i++ )
		for( int j = 0; j < gridSize; j++ )
		{
		const int elementIdx = i * gridSize + j;
		if( i == 0 \|\| j == 0 \|\| i == gridSize - 1 \|\| j == gridSize - 1 ) // check matrix elements corresponding to boundary grid nodes
		{
		for( int k = 0; k < matrixSize; k++ )
		{
		if( elementIdx == k )
		EXPECT_EQ( matrix.getElement( elementIdx, k ), 1.0 );
		else
		EXPECT_EQ( matrix.getElement( elementIdx, k ), 0.0 );
		}
		}
		else // check matrix elements corresponding to inner grid nodes
		{
		for( int k = 0; k < matrixSize; k++ )
		{
		if( k == elementIdx - gridSize \|\|
		k == elementIdx - 1 \|\|
		k == elementIdx + 1 \|\|
		k == elementIdx + gridSize )
		EXPECT_EQ( matrix.getElement( elementIdx, k ), -1.0 );
		else if( k == elementIdx )
		EXPECT_EQ( matrix.getElement( elementIdx, k ), 4.0 );
		else
		EXPECT_EQ( matrix.getElement( elementIdx, k ), 0.0 );
		}
		}
		}
		}

		template< typename Matrix >
		void test_GetNonemptyRowsCount()
		{
		@@ -1381,6 +1443,13 @@ TYPED_TEST( MatrixTest, setLikeTest )
		test_SetLike< MatrixType, MatrixType >();
		}

		TYPED_TEST( MatrixTest, setElements )
		{
		using MatrixType = typename TestFixture::MatrixType;

		test_SetElements< MatrixType >();
		}

		TYPED_TEST( MatrixTest, getNonemptyRowsCountTest )
		{
		using MatrixType = typename TestFixture::MatrixType;
		@@ -1388,7 +1457,6 @@ TYPED_TEST( MatrixTest, getNonemptyRowsCountTest )
		test_GetNonemptyRowsCount< MatrixType >();
		}


		TYPED_TEST( MatrixTest, getCompressedRowLengthTest )
		{
		using MatrixType = typename TestFixture::MatrixType;