diff --git a/src/core/arrays/tnlArray.h b/src/core/arrays/tnlArray.h index d9e0c6f33e4e003885017453630e4f2c1eea6f40..b51184959845a9d7c356dd52114d18a139d74d4c 100644 --- a/src/core/arrays/tnlArray.h +++ b/src/core/arrays/tnlArray.h @@ -30,7 +30,7 @@ class tnlSharedArray; template< typename Element, typename Device = tnlHost, typename Index = int > -class tnlArray : public tnlObject +class tnlArray : public virtual tnlObject { public: diff --git a/src/core/tnlCuda.h b/src/core/tnlCuda.h index 48069908c0186a5a290593447e0740016227fa3d..89c67df3d3a4f01cb9174c16be08a4df9212ea75 100644 --- a/src/core/tnlCuda.h +++ b/src/core/tnlCuda.h @@ -19,6 +19,7 @@ #define TNLCUDA_H_ #include <iostream> +#include <unistd.h> #include <core/tnlDevice.h> #include <core/tnlString.h> #include <core/tnlAssert.h> @@ -43,6 +44,8 @@ class tnlCuda static bool checkDevice( const char* file_name, int line ); + static size_t getFreeMemory(); + protected: static int maxGridSize, maxBlockSize; diff --git a/src/core/tnlHost.h b/src/core/tnlHost.h index 8adcd3ed6a83275515bc7caa3b56f919358c7af1..16ad2fa9329a855956a089683b2c269bbc3f98bb 100644 --- a/src/core/tnlHost.h +++ b/src/core/tnlHost.h @@ -18,6 +18,7 @@ #ifndef TNLHOST_H_ #define TNLHOST_H_ +#include <unistd.h> #include <core/tnlDevice.h> #include <core/tnlString.h> @@ -28,6 +29,8 @@ class tnlHost static tnlString getDeviceType(); static tnlDeviceEnum getDevice(); + + static size_t getFreeMemory(); }; #endif /* TNLHOST_H_ */ diff --git a/src/implementation/core/tnlCuda.cpp b/src/implementation/core/tnlCuda.cpp index 36856fe5c25581ac1d85f4d0f34dcd9366f0cb99..a06f45ed37b16d222d0637a511d9082145b8235a 100644 --- a/src/implementation/core/tnlCuda.cpp +++ b/src/implementation/core/tnlCuda.cpp @@ -421,3 +421,8 @@ bool tnlCuda::checkDevice( const char* file_name, int line ) #endif } +size_t tnlCuda::getFreeMemory() +{ + +} + diff --git a/src/implementation/core/tnlHost.cpp b/src/implementation/core/tnlHost.cpp index ea839b96213339f663f54a459734a83fa47663a3..837323be293497fe36fa648ea767d05be2d0572d 100644 --- a/src/implementation/core/tnlHost.cpp +++ b/src/implementation/core/tnlHost.cpp @@ -20,14 +20,21 @@ #include <core/tnlHost.h> -tnlString tnlHost :: getDeviceType() +tnlString tnlHost::getDeviceType() { return tnlString( "tnlHost" ); }; -tnlDeviceEnum tnlHost :: getDevice() +tnlDeviceEnum tnlHost::getDevice() { return tnlHostDevice; }; +size_t tnlHost::getFreeMemory() +{ + long pages = sysconf(_SC_PHYS_PAGES); + long page_size = sysconf(_SC_PAGE_SIZE); + return pages * page_size; +}; + #endif /* TNLHOST_H_ */ diff --git a/src/implementation/matrices/CMakeLists.txt b/src/implementation/matrices/CMakeLists.txt index e31fb4bb42b9f0df17b5cb972fbe8c715f2e5e1f..b4a222734d75bddd400df5184eae38f024b9808c 100755 --- a/src/implementation/matrices/CMakeLists.txt +++ b/src/implementation/matrices/CMakeLists.txt @@ -2,6 +2,7 @@ SET( headers tnlMatrix_impl.h tnlDenseMatrix_impl.h tnlTridiagonalMatrix_impl.h tnlMultidiagonalMatrix_impl.h + tnlSparseMatrix_impl.h tnlEllpackMatrix_impl.h tnlSlicedEllpackMatrix_impl.h tnlChunkedEllpackMatrix_impl.h diff --git a/src/implementation/matrices/tnlCSRMatrix_impl.h b/src/implementation/matrices/tnlCSRMatrix_impl.h index 76d5e963a1dc29fc0bd1518de2ae52012d41504d..580eb7a8e28dcee0b1e7bd89c89950bb3c8ee32e 100644 --- a/src/implementation/matrices/tnlCSRMatrix_impl.h +++ b/src/implementation/matrices/tnlCSRMatrix_impl.h @@ -26,8 +26,6 @@ template< typename Real, typename Device, typename Index > tnlCSRMatrix< Real, Device, Index >::tnlCSRMatrix() -: rows( 0 ), - columns( 0 ) { }; @@ -55,14 +53,10 @@ template< typename Real, typename Device, typename Index > bool tnlCSRMatrix< Real, Device, Index >::setDimensions( const IndexType rows, - const IndexType columns ) + const IndexType columns ) { - tnlAssert( rows > 0 && columns > 0, - cerr << "rows = " << rows - << " columns = " << columns << endl ); - this->rows = rows; - this->columns = columns; - if( ! this->rowPointers.setSize( this->rows + 1 ) ) + if( ! tnlSparseMatrix< Real, Device, Index >::setDimensions( rows, columns ) || + ! this->rowPointers.setSize( this->rows + 1 ) ) return false; this->rowPointers.setValue( 0 ); return true; @@ -71,8 +65,7 @@ bool tnlCSRMatrix< Real, Device, Index >::setDimensions( const IndexType rows, template< typename Real, typename Device, typename Index > - template< typename Vector > -bool tnlCSRMatrix< Real, Device, Index >::setRowLengths( const Vector& rowLengths ) +bool tnlCSRMatrix< Real, Device, Index >::setRowLengths( const RowLengthsVector& rowLengths ) { /**** * Compute the rows pointers. The last one is @@ -103,50 +96,21 @@ template< typename Real, typename Index2 > bool tnlCSRMatrix< Real, Device, Index >::setLike( const tnlCSRMatrix< Real2, Device2, Index2 >& matrix ) { - if( ! this->setDimensions( matrix.getRows(), matrix.getColumns() ) || - ! this->values.setLike( matrix.values ) || - ! this->columnIndexes.setLike( matrix.columnIndexes ) || + if( ! tnlSparseMatrix< Real, Device, Index >::setLike( matrix ) || ! this->rowPointers.setLike( matrix.rowPointers ) ) return false; return true; } -template< typename Real, - typename Device, - typename Index > -Index tnlCSRMatrix< Real, Device, Index >::getNumberOfAllocatedElements() const -{ - return this->values.getSize(); -} - template< typename Real, typename Device, typename Index > void tnlCSRMatrix< Real, Device, Index >::reset() { - this->columns = 0; - this->rows = 0; - this->values.reset(); - this->columnIndexes.reset(); + tnlSparseMatrix< Real, Device, Index >::reset(); this->rowPointers.reset(); } -template< typename Real, - typename Device, - typename Index > -Index tnlCSRMatrix< Real, Device, Index >::getRows() const -{ - return this->rows; -} - -template< typename Real, - typename Device, - typename Index > -Index tnlCSRMatrix< Real, Device, Index >::getColumns() const -{ - return this->columns; -} - template< typename Real, typename Device, typename Index > @@ -197,7 +161,7 @@ Real tnlCSRMatrix< Real, Device, Index >::getElement( const IndexType row, const IndexType rowEnd = this->rowPointers[ row + 1 ]; while( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] < column ) elementPtr++; - if( this->columnIndexes[ elementPtr ] == column ) + if( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] == column ) return this->values[ elementPtr ]; return 0.0; } @@ -347,11 +311,9 @@ template< typename Real, typename Index > bool tnlCSRMatrix< Real, Device, Index >::save( tnlFile& file ) const { - if( ! file.write( &this->rows ) ) return false; - if( ! file.write( &this->columns ) ) return false; - if( ! this->values.save( file ) ) return false; - if( ! this->columnIndexes.save( file ) ) return false; - if( ! this->rowPointers.save( file ) ) return false; + if( ! tnlSparseMatrix< Real, Device, Index >::save( file ) || + ! this->rowPointers.save( file ) ) + return false; return true; } @@ -360,11 +322,9 @@ template< typename Real, typename Index > bool tnlCSRMatrix< Real, Device, Index >::load( tnlFile& file ) { - if( ! file.read( &this->rows ) ) return false; - if( ! file.read( &this->columns ) ) return false; - if( ! this->values.load( file ) ) return false; - if( ! this->columnIndexes.load( file ) ) return false; - if( ! this->rowPointers.load( file ) ) return false; + if( ! tnlSparseMatrix< Real, Device, Index >::load( file ) || + ! this->rowPointers.load( file ) ) + return false; return true; } diff --git a/src/implementation/matrices/tnlChunkedEllpackMatrix_impl.h b/src/implementation/matrices/tnlChunkedEllpackMatrix_impl.h index c44679a0587aa57ad91bfd095158b3124a96acb2..46ab5408e39fc4952ae05c0e2e08ab9f92502e09 100644 --- a/src/implementation/matrices/tnlChunkedEllpackMatrix_impl.h +++ b/src/implementation/matrices/tnlChunkedEllpackMatrix_impl.h @@ -27,13 +27,11 @@ template< typename Real, typename Device, typename Index > tnlChunkedEllpackMatrix< Real, Device, Index >::tnlChunkedEllpackMatrix() -: rows( 0 ), - columns( 0 ), - chunksInSlice( 256 ), +: chunksInSlice( 256 ), desiredChunkSize( 16 ) { - values.setName( "tnlChunkedEllpackMatrix::values" ); - columnIndexes.setName( "tnlChunkedEllpackMatrix::columnIndexes" ); + this->values.setName( "tnlChunkedEllpackMatrix::values" ); + this->columnIndexes.setName( "tnlChunkedEllpackMatrix::columnIndexes" ); chunksToRowsMapping.setName( "tnlChunkedEllpackMatrix::chunksToRowsMapping" ); slicesToRowsMapping.setName( "tnlChunkedEllpackMatrix::slicesToRowsMapping" ); rowPointers.setName( "tnlChunkedEllpackMatrix::rowPointers" ); @@ -69,8 +67,8 @@ bool tnlChunkedEllpackMatrix< Real, Device, Index >::setDimensions( const IndexT tnlAssert( rows > 0 && columns > 0, cerr << "rows = " << rows << " columns = " << columns << endl ); - this->rows = rows; - this->columns = columns; + if( ! tnlSparseMatrix< Real, Device, Index >::setDimensions( rows, columns ) ) + return false; /**** * Allocate slice info array. Note that there cannot be @@ -87,8 +85,7 @@ bool tnlChunkedEllpackMatrix< Real, Device, Index >::setDimensions( const IndexT template< typename Real, typename Device, typename Index > - template< typename Vector > -bool tnlChunkedEllpackMatrix< Real, Device, Index >::setRowLengths( const Vector& rowLengths ) +bool tnlChunkedEllpackMatrix< Real, Device, Index >::setRowLengths( const RowLengthsVector& rowLengths ) { /**** * Iterate over rows and allocate slices so that each slice has @@ -137,9 +134,14 @@ bool tnlChunkedEllpackMatrix< Real, Device, Index >::setRowLengths( const Vector RealType rowRatio( 0.0 ); if( allocatedElementsInSlice != 0 ) rowRatio = ( RealType ) rowLengths[ i ] / ( RealType ) allocatedElementsInSlice; - freeChunks -= this->chunksToRowsMapping[ i ] = ceil( freeChunks * rowRatio ); + const IndexType addedChunks = ceil( freeChunks * rowRatio ); + freeChunks -= addedChunks; + this->chunksToRowsMapping[ i ] += addedChunks; + tnlAssert( chunksToRowsMapping[ i ] > 0, + cerr << " chunksToRowsMapping[ i ] = " << chunksToRowsMapping[ i ] << endl ); } tnlAssert( freeChunks >= 0, ); + //cout << "freeChunks = " << freeChunks << endl; } /**** @@ -150,6 +152,9 @@ bool tnlChunkedEllpackMatrix< Real, Device, Index >::setRowLengths( const Vector maxChunkInSlice = Max( maxChunkInSlice, ceil( ( RealType ) rowLengths[ i ] / ( RealType ) this->chunksToRowsMapping[ i ] ) ); + tnlAssert( maxChunkInSlice > 0, + cerr << " maxChunkInSlice = " << maxChunkInSlice << endl ); + /**** * Set-up the slice info. @@ -165,8 +170,14 @@ bool tnlChunkedEllpackMatrix< Real, Device, Index >::setRowLengths( const Vector sliceIndex++; for( IndexType i = sliceBegin; i < sliceEnd; i++ ) + { this->rowPointers[ i + 1 ] = this->rowPointers[ i ] + maxChunkInSlice*chunksToRowsMapping[ i ]; + tnlAssert( this->rowPointers[ i ] >= 0, + cerr << "this->rowPointers[ i ] = " << this->rowPointers[ i ] ); + tnlAssert( this->rowPointers[ i + 1 ] >= 0, + cerr << "this->rowPointers[ i + 1 ] = " << this->rowPointers[ i + 1 ] ); + } /**** * Finish the chunks to rows mapping by computing the prefix sum. @@ -177,18 +188,12 @@ bool tnlChunkedEllpackMatrix< Real, Device, Index >::setRowLengths( const Vector * Proceed to the next row */ sliceSize = 0; + allocatedElementsInSlice = 0; if( row < this->rows ) continue; else break; } - /**** - * Allocate values and column indexes - */ - if( ! this->values.setSize( elementsToAllocation ) || - ! this->columnIndexes.setSize( elementsToAllocation ) ) - return false; - this->columnIndexes.setValue( this->columns ); - return true; + return tnlSparseMatrix< Real, Device, Index >::allocateMatrixElements( elementsToAllocation ); } template< typename Real, @@ -201,52 +206,23 @@ bool tnlChunkedEllpackMatrix< Real, Device, Index >::setLike( const tnlChunkedEl { this->chunksInSlice = matrix.chunksInSlice; this->desiredChunkSize = matrix.desiredChunkSize; - if( ! this->setDimensions( matrix.getRows(), matrix.getColumns() ) || - ! this->values.setLike( matrix.values ) || - ! this->columnIndexes.setLike( matrix.columnIndexes ) || + if( ! tnlSparseMatrix< Real, Device, Index >::setLike( matrix ) || ! this->chunksToRowsMapping.setLike( matrix.chunksToRowsMapping ) || ! this->slicesToRowsMapping.setLike( matrix.slicesToRowsMapping ) ) return false; return true; } -template< typename Real, - typename Device, - typename Index > -Index tnlChunkedEllpackMatrix< Real, Device, Index >::getNumberOfAllocatedElements() const -{ - return this->values.getSize(); -} - template< typename Real, typename Device, typename Index > void tnlChunkedEllpackMatrix< Real, Device, Index >::reset() { - this->columns = 0; - this->rows = 0; - this->values.reset(); - this->columnIndexes.reset(); + tnlSparseMatrix< Real, Device, Index >::reset(); this->chunksToRowsMapping.reset(); this->slicesToRowsMapping.reset(); } -template< typename Real, - typename Device, - typename Index > -Index tnlChunkedEllpackMatrix< Real, Device, Index >::getRows() const -{ - return this->rows; -} - -template< typename Real, - typename Device, - typename Index > -Index tnlChunkedEllpackMatrix< Real, Device, Index >::getColumns() const -{ - return this->columns; -} - template< typename Real, typename Device, typename Index > @@ -330,9 +306,11 @@ Real tnlChunkedEllpackMatrix< Real, Device, Index >::getElement( const IndexType const IndexType& chunkSize = slices.getElement( sliceIndex ).chunkSize; IndexType elementPtr = rowPointers[ row ]; const IndexType rowEnd = rowPointers[ row + 1 ]; + tnlAssert( rowEnd <= this->columnIndexes.getSize(), + cerr << "rowEnd = " << rowEnd << " this->columnIndexes.getSize() = " << this->columnIndexes.getSize() ); while( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] < column ) elementPtr++; - if( this->columnIndexes[ elementPtr ] == column ) + if( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] == column ) return this->values[ elementPtr ]; return 0.0; } @@ -358,6 +336,11 @@ bool tnlChunkedEllpackMatrix< Real, Device, Index >::addElement( const IndexType IndexType elementPtr = rowPointers[ row ]; const IndexType rowEnd = rowPointers[ row + 1 ]; + tnlAssert( elementPtr >= 0, + cerr << "elementPtr = " << elementPtr ); + tnlAssert( rowEnd <= this->columnIndexes.getSize(), + cerr << "rowEnd = " << rowEnd << " this->columnIndexes.getSize() = " << this->columnIndexes.getSize() ); + while( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] < column ) elementPtr++; if( elementPtr == rowEnd ) return false; @@ -496,14 +479,12 @@ template< typename Real, typename Index > bool tnlChunkedEllpackMatrix< Real, Device, Index >::save( tnlFile& file ) const { - if( ! file.write( &this->rows ) ) return false; - if( ! file.write( &this->columns ) ) return false; - if( ! this->values.save( file ) ) return false; - if( ! this->columnIndexes.save( file ) ) return false; - if( ! this->chunksToRowsMapping.save( file ) ) return false; - if( ! this->slicesToRowsMapping.save( file ) ) return false; - if( ! this->rowPointers.save( file ) ) return false; - if( ! this->slices.save( file ) ) return false; + if( ! tnlSparseMatrix< Real, Device, Index >::save( file ) || + ! this->chunksToRowsMapping.save( file ) || + ! this->slicesToRowsMapping.save( file ) || + ! this->rowPointers.save( file ) || + ! this->slices.save( file ) ) + return false; return true; } @@ -512,14 +493,12 @@ template< typename Real, typename Index > bool tnlChunkedEllpackMatrix< Real, Device, Index >::load( tnlFile& file ) { - if( ! file.read( &this->rows ) ) return false; - if( ! file.read( &this->columns ) ) return false; - if( ! this->values.load( file ) ) return false; - if( ! this->columnIndexes.load( file ) ) return false; - if( ! this->chunksToRowsMapping.load( file ) ) return false; - if( ! this->slicesToRowsMapping.load( file ) ) return false; - if( ! this->rowPointers.load( file ) ) return false; - if( ! this->slices.load( file ) ) return false; + if( ! tnlSparseMatrix< Real, Device, Index >::load( file ) || + ! this->chunksToRowsMapping.load( file ) || + ! this->slicesToRowsMapping.load( file ) || + ! this->rowPointers.load( file ) || + ! this->slices.load( file ) ) + return false; return true; } diff --git a/src/implementation/matrices/tnlDenseMatrix_impl.h b/src/implementation/matrices/tnlDenseMatrix_impl.h index 5983c6f00b7be79bfa7a71e0ba9b80ea03db2f69..c610061bab14a466c46f9617e0a03ba0bd81e024 100644 --- a/src/implementation/matrices/tnlDenseMatrix_impl.h +++ b/src/implementation/matrices/tnlDenseMatrix_impl.h @@ -53,15 +53,28 @@ template< typename Real, bool tnlDenseMatrix< Real, Device, Index >::setDimensions( const IndexType rows, const IndexType columns ) { - return tnlMultiArray< 2, Real, Device, Index >::setDimensions( columns, rows ); + if( ! tnlMatrix< Real, Device, Index >::setDimensions( rows, columns ) || + ! tnlMultiArray< 2, Real, Device, Index >::setDimensions( rows, columns ) ) + return false; tnlMultiArray< 2, Real, Device, Index >::setValue( 0.0 ); + return true; } template< typename Real, typename Device, typename Index > - template< typename Vector > -bool tnlDenseMatrix< Real, Device, Index >::setRowLengths( const Vector& rowLengths ) + template< typename Real2, + typename Device2, + typename Index2 > +bool tnlDenseMatrix< Real, Device, Index >::setLike( const tnlDenseMatrix< Real2, Device2, Index2 >& matrix ) +{ + return this->setDimensions( matrix.getRows(), matrix.getColumns() ); +} + +template< typename Real, + typename Device, + typename Index > +bool tnlDenseMatrix< Real, Device, Index >::setRowLengths( const RowLengthsVector& rowLengths ) { return true; } @@ -69,7 +82,7 @@ bool tnlDenseMatrix< Real, Device, Index >::setRowLengths( const Vector& rowLeng template< typename Real, typename Device, typename Index > -Index tnlDenseMatrix< Real, Device, Index >::getNumberOfAllocatedElements() const +Index tnlDenseMatrix< Real, Device, Index >::getNumberOfMatrixElements() const { return this->getRows() * this->getColumns(); } @@ -77,17 +90,30 @@ Index tnlDenseMatrix< Real, Device, Index >::getNumberOfAllocatedElements() cons template< typename Real, typename Device, typename Index > -Index tnlDenseMatrix< Real, Device, Index >::getRows() const +void tnlDenseMatrix< Real, Device, Index >::reset() { - return this->getDimensions().x(); + tnlMatrix< Real, Device, Index >::reset(); + tnlMultiArray< 2, Real, Device, Index >::reset(); } template< typename Real, typename Device, typename Index > -Index tnlDenseMatrix< Real, Device, Index >::getColumns() const +bool tnlDenseMatrix< Real, Device, Index >::setElement( const IndexType row, + const IndexType column, + const RealType& value ) { - return this->getDimensions().y(); + tnlMultiArray< 2, Real, Device, Index >::setElement( row, column, value ); + return true; +} + +template< typename Real, + typename Device, + typename Index > +Real tnlDenseMatrix< Real, Device, Index >::getElement( const IndexType row, + const IndexType column ) const +{ + return tnlMultiArray< 2, Real, Device, Index >::getElement( row, column ); } template< typename Real, @@ -293,6 +319,44 @@ void tnlDenseMatrix< Real, Device, Index >::performSORIteration( const Vector& b x[ row ] += omega / this->operator()( row, row )( b[ row ] - sum ); } +template< typename Real, + typename Device, + typename Index > +bool tnlDenseMatrix< Real, Device, Index >::save( const tnlString& fileName ) const +{ + return tnlObject::save( fileName ); +} + +template< typename Real, + typename Device, + typename Index > +bool tnlDenseMatrix< Real, Device, Index >::load( const tnlString& fileName ) +{ + return tnlObject::load( fileName ); +} + +template< typename Real, + typename Device, + typename Index > +bool tnlDenseMatrix< Real, Device, Index >::save( tnlFile& file ) const +{ + if( ! tnlMatrix< Real, Device, Index >::save( file ) || + ! tnlMultiArray< 2, Real, Device, Index >::save( file ) ) + return false; + return true; +} + +template< typename Real, + typename Device, + typename Index > +bool tnlDenseMatrix< Real, Device, Index >::load( tnlFile& file ) +{ + if( ! tnlMatrix< Real, Device, Index >::load( file ) || + ! tnlMultiArray< 2, Real, Device, Index >::load( file ) ) + return false; + return true; +} + template< typename Real, typename Device, typename Index > diff --git a/src/implementation/matrices/tnlEllpackMatrix_impl.h b/src/implementation/matrices/tnlEllpackMatrix_impl.h index 61a3fa6f74eb296d238ca40dc24ed99806f3333b..611d0d83ecdfbf67bd55f2c25ff98f21e6a40ed1 100644 --- a/src/implementation/matrices/tnlEllpackMatrix_impl.h +++ b/src/implementation/matrices/tnlEllpackMatrix_impl.h @@ -26,9 +26,7 @@ template< typename Real, typename Device, typename Index > tnlEllpackMatrix< Real, Device, Index > :: tnlEllpackMatrix() -: rows( 0 ), - columns( 0 ), - rowLengths( 0 ) +: rowLengths( 0 ) { }; @@ -55,8 +53,8 @@ tnlString tnlEllpackMatrix< Real, Device, Index >::getTypeVirtual() const template< typename Real, typename Device, typename Index > -bool tnlEllpackMatrix< Real, Device, Index > :: setDimensions( const IndexType rows, - const IndexType columns ) +bool tnlEllpackMatrix< Real, Device, Index >::setDimensions( const IndexType rows, + const IndexType columns ) { tnlAssert( rows > 0 && columns > 0, cerr << "rows = " << rows @@ -71,8 +69,7 @@ bool tnlEllpackMatrix< Real, Device, Index > :: setDimensions( const IndexType r template< typename Real, typename Device, typename Index > - template< typename Vector > -bool tnlEllpackMatrix< Real, Device, Index > :: setRowLengths( const Vector& rowLengths ) +bool tnlEllpackMatrix< Real, Device, Index >::setRowLengths( const RowLengthsVector& rowLengths ) { this->rowLengths = 0; for( IndexType i = 0; i < rowLengths.getSize(); i++ ) @@ -87,7 +84,7 @@ bool tnlEllpackMatrix< Real, Device, Index > :: setRowLengths( const Vector& row template< typename Real, typename Device, typename Index > -bool tnlEllpackMatrix< Real, Device, Index > :: setConstantRowLengths( const IndexType& rowLengths ) +bool tnlEllpackMatrix< Real, Device, Index >::setConstantRowLengths( const IndexType& rowLengths ) { tnlAssert( rowLengths > 0, cerr << " rowLengths = " << rowLengths ); @@ -103,49 +100,21 @@ template< typename Real, template< typename Real2, typename Device2, typename Index2 > -bool tnlEllpackMatrix< Real, Device, Index > :: setLike( const tnlEllpackMatrix< Real2, Device2, Index2 >& matrix ) +bool tnlEllpackMatrix< Real, Device, Index >::setLike( const tnlEllpackMatrix< Real2, Device2, Index2 >& matrix ) { - this->rowLengths = 0; - this->setDimensions( matrix.getRows(), matrix.getColumns() ); - if( ! this->setConstantRowLengths( matrix.rowLengths ) ) + if( ! tnlSparseMatrix< Real, Device, Index >::setLike( matrix ) ) return false; + this->rowLengths = matrix.rowLengths; return true; } -template< typename Real, - typename Device, - typename Index > -Index tnlEllpackMatrix< Real, Device, Index > :: getNumberOfAllocatedElements() const -{ - return this->values.getSize(); -} - template< typename Real, typename Device, typename Index > void tnlEllpackMatrix< Real, Device, Index > :: reset() { - this->rows = 0; - this->columns = 0; + tnlSparseMatrix< Real, Device, Index >::reset(); this->rowLengths = 0; - this->values.reset(); - this->columnIndexes.reset(); -} - -template< typename Real, - typename Device, - typename Index > -Index tnlEllpackMatrix< Real, Device, Index >::getRows() const -{ - return this->rows; -} - -template< typename Real, - typename Device, - typename Index > -Index tnlEllpackMatrix< Real, Device, Index >::getColumns() const -{ - return this->columns; } template< typename Real, @@ -194,12 +163,12 @@ template< typename Real, Real tnlEllpackMatrix< Real, Device, Index >::getElement( const IndexType row, const IndexType column ) const { - IndexType i( row * this->rowLengths ); - const IndexType rowEnd( i + this->rowLengths ); - while( i < rowEnd && this->columnIndexes[ i ] < column ) i++; - if( i == rowEnd || this->columnIndexes[ i ] != column ) - return 0.0; - return this->values.getElement( i ); + IndexType elementPtr( row * this->rowLengths ); + const IndexType rowEnd( elementPtr + this->rowLengths ); + while( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] < column ) elementPtr++; + if( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] == column ) + return this->values.getElement( elementPtr ); + return 0.0; } template< typename Real, @@ -339,11 +308,8 @@ template< typename Real, typename Index > bool tnlEllpackMatrix< Real, Device, Index >::save( tnlFile& file ) const { - if( ! file.write( &this->rows ) ) return false; - if( ! file.write( &this->columns ) ) return false; + if( ! tnlSparseMatrix< Real, Device, Index >::save( file) ) return false; if( ! file.write( &this->rowLengths ) ) return false; - if( ! this->values.save( file ) ) return false; - if( ! this->columnIndexes.save( file ) ) return false; return true; } @@ -352,11 +318,8 @@ template< typename Real, typename Index > bool tnlEllpackMatrix< Real, Device, Index >::load( tnlFile& file ) { - if( ! file.read( &this->rows ) ) return false; - if( ! file.read( &this->columns ) ) return false; + if( ! tnlSparseMatrix< Real, Device, Index >::load( file) ) return false; if( ! file.read( &this->rowLengths ) ) return false; - if( ! this->values.load( file ) ) return false; - if( ! this->columnIndexes.load( file ) ) return false; return true; } @@ -401,15 +364,8 @@ template< typename Real, typename Index > bool tnlEllpackMatrix< Real, Device, Index >::allocateElements() { - if( ! this->values.setSize( this->rows * this->rowLengths ) || - ! this->columnIndexes.setSize( this->rows * this->rowLengths ) ) + if( ! tnlSparseMatrix< Real, Device, Index >::allocateMatrixElements( this->rows * this->rowLengths ) ) return false; - - /**** - * Setting a column index to this->columns means that the - * index is undefined. - */ - this->columnIndexes.setValue( this->columns ); return true; } diff --git a/src/implementation/matrices/tnlMatrixReader_impl.h b/src/implementation/matrices/tnlMatrixReader_impl.h index ecfc08867390255c83f759040d87b95f0f20b697..a03a38676400dd525f686ec66828dd0ff35f3b5b 100644 --- a/src/implementation/matrices/tnlMatrixReader_impl.h +++ b/src/implementation/matrices/tnlMatrixReader_impl.h @@ -78,7 +78,7 @@ bool tnlMatrixReader< Matrix >::verifyMtxFile( std::istream& file, dimensionsLine = true; continue; } - IndexType row, column; + IndexType row( 1 ), column( 1 ); RealType value; if( ! parseMtxLineWithElement( line, row, column, value ) ) return false; @@ -95,17 +95,52 @@ bool tnlMatrixReader< Matrix >::verifyMtxFile( std::istream& file, if( symmetricMatrix && row != column ) processedElements++; if( verbose ) - cout << " Verifying the matrix elements ... " << processedElements << " / " << matrix.getNumberOfAllocatedElements() << " \r" << flush; + cout << " Verifying the matrix elements ... " << processedElements << " / " << matrix.getNumberOfMatrixElements() << " \r" << flush; } file.clear(); long int fileSize = file.tellg(); if( verbose ) - cout << " Verifying the matrix elements ... " << processedElements << " / " << matrix.getNumberOfAllocatedElements() + cout << " Verifying the matrix elements ... " << processedElements << " / " << matrix.getNumberOfMatrixElements() << " -> " << timer.GetTime() << " sec. i.e. " << fileSize / ( timer.GetTime() * ( 1 << 20 )) << "MB/s." << endl; return true; } +template< typename Matrix > +bool tnlMatrixReader< Matrix >::findLineByElement( std::istream& file, + const IndexType& row, + const IndexType& column, + tnlString& line, + IndexType& lineNumber ) +{ + file.clear(); + file.seekg( 0, ios::beg ); + bool symmetricMatrix( false ); + bool dimensionsLine( false ); + lineNumber = 0; + tnlTimerRT timer; + IndexType currentRow, currentColumn; + RealType value; + while( line.getLine( file ) ) + { + lineNumber++; + if( line[ 0 ] == '%' ) continue; + if( ! dimensionsLine ) + { + dimensionsLine = true; + continue; + } + IndexType currentRow( 1 ), currentColumn( 1 ); + RealType value; + if( ! parseMtxLineWithElement( line, currentRow, currentColumn, value ) ) + return false; + if( ( currentRow == row + 1 && currentColumn == column + 1 ) || + ( symmetricMatrix && currentRow == column + 1 && currentColumn == row + 1 ) ) + return true; + } + return false; +} + template< typename Matrix > bool tnlMatrixReader< Matrix >::checkMtxHeader( const tnlString& header, bool& symmetric ) @@ -221,7 +256,7 @@ bool tnlMatrixReader< Matrix >::computeRowLengthsFromMtxFile( std::istream& file dimensionsLine = true; continue; } - IndexType row, column; + IndexType row( 1 ), column( 1 ); RealType value; if( ! parseMtxLineWithElement( line, row, column, value ) ) return false; @@ -261,7 +296,7 @@ bool tnlMatrixReader< Matrix >::readMatrixElementsFromMtxFile( std::istream& fil dimensionsLine = true; continue; } - IndexType row, column; + IndexType row( 1 ), column( 1 ); RealType value; if( ! parseMtxLineWithElement( line, row, column, value ) ) return false; @@ -273,12 +308,12 @@ bool tnlMatrixReader< Matrix >::readMatrixElementsFromMtxFile( std::istream& fil processedElements++; } if( verbose ) - cout << " Reading the matrix elements ... " << processedElements << " / " << matrix.getNumberOfAllocatedElements() << " \r" << flush; + cout << " Reading the matrix elements ... " << processedElements << " / " << matrix.getNumberOfMatrixElements() << " \r" << flush; } file.clear(); long int fileSize = file.tellg(); if( verbose ) - cout << " Reading the matrix elements ... " << processedElements << " / " << matrix.getNumberOfAllocatedElements() + cout << " Reading the matrix elements ... " << processedElements << " / " << matrix.getNumberOfMatrixElements() << " -> " << timer.GetTime() << " sec. i.e. " << fileSize / ( timer.GetTime() * ( 1 << 20 )) << "MB/s." << endl; return true; diff --git a/src/implementation/matrices/tnlMatrix_impl.h b/src/implementation/matrices/tnlMatrix_impl.h index 346eb0fce8a57f07b17bbf51dd66fe17f4ec67fb..0b8fa5f6ad745399bbc8364198a4de4b34581a0c 100644 --- a/src/implementation/matrices/tnlMatrix_impl.h +++ b/src/implementation/matrices/tnlMatrix_impl.h @@ -18,7 +18,88 @@ #ifndef TNLMATRIX_IMPL_H_ #define TNLMATRIX_IMPL_H_ +#include <matrices/tnlMatrix.h> +template< typename Real, + typename Device, + typename Index > +tnlMatrix< Real, Device, Index >::tnlMatrix() +: rows( 0 ), + columns( 0 ) +{ +} +template< typename Real, + typename Device, + typename Index > + bool tnlMatrix< Real, Device, Index >::setDimensions( const IndexType rows, + const IndexType columns ) +{ + tnlAssert( rows > 0 && columns > 0, + cerr << " rows = " << rows << " columns = " << columns ); + this->rows = rows; + this->columns = columns; + return true; +} + +template< typename Real, + typename Device, + typename Index > + template< typename Real2, + typename Device2, + typename Index2 > +bool tnlMatrix< Real, Device, Index >::setLike( const tnlMatrix< Real2, Device2, Index2 >& matrix ) +{ + return setDimensions( matrix.getRows(), matrix.getColumns() ); +} + +template< typename Real, + typename Device, + typename Index > +Index tnlMatrix< Real, Device, Index >::getRows() const +{ + return this->rows; +} + +template< typename Real, + typename Device, + typename Index > +Index tnlMatrix< Real, Device, Index >::getColumns() const +{ + return this->columns; +} + +template< typename Real, + typename Device, + typename Index > +void tnlMatrix< Real, Device, Index >::reset() +{ + this->rows = 0; + this->columns = 0; +} + +template< typename Real, + typename Device, + typename Index > +bool tnlMatrix< Real, Device, Index >::save( tnlFile& file ) const +{ + if( ! tnlObject::save( file ) || + ! file.write( &this->rows ) || + ! file.write( &this->columns ) ) + return false; + return true; +} + +template< typename Real, + typename Device, + typename Index > +bool tnlMatrix< Real, Device, Index >::load( tnlFile& file ) +{ + if( ! tnlObject::load( file ) || + ! file.read( &this->rows ) || + ! file.read( &this->columns ) ) + return false; + return true; +} #endif /* TNLMATRIX_IMPL_H_ */ diff --git a/src/implementation/matrices/tnlMultidiagonalMatrix_impl.h b/src/implementation/matrices/tnlMultidiagonalMatrix_impl.h index 78735c15425050ccaa8c14d2470e073f43e65dd6..512f12ff5b55869f02d76e20eeb920de51c983bf 100644 --- a/src/implementation/matrices/tnlMultidiagonalMatrix_impl.h +++ b/src/implementation/matrices/tnlMultidiagonalMatrix_impl.h @@ -26,8 +26,6 @@ template< typename Real, typename Device, typename Index > tnlMultidiagonalMatrix< Real, Device, Index > :: tnlMultidiagonalMatrix() -: rows( 0 ), - columns( 0 ) { }; @@ -54,14 +52,14 @@ tnlString tnlMultidiagonalMatrix< Real, Device, Index >::getTypeVirtual() const template< typename Real, typename Device, typename Index > -bool tnlMultidiagonalMatrix< Real, Device, Index > :: setDimensions( const IndexType rows, - const IndexType columns ) +bool tnlMultidiagonalMatrix< Real, Device, Index >::setDimensions( const IndexType rows, + const IndexType columns ) { tnlAssert( rows > 0 && columns > 0, cerr << "rows = " << rows << " columns = " << columns << endl ); - this->rows = rows; - this->columns = columns; + if( ! tnlMatrix< Real, Device, Index >::setDimensions( rows, columns ) ) + return false; if( this->diagonalsShift.getSize() != 0 ) { if( ! this->values.setSize( Min( this->rows, this->columns ) * this->diagonalsShift.getSize() ) ) @@ -71,6 +69,17 @@ bool tnlMultidiagonalMatrix< Real, Device, Index > :: setDimensions( const Index return true; } +template< typename Real, + typename Device, + typename Index > +bool tnlMultidiagonalMatrix< Real, Device, Index >::setRowLengths( const RowLengthsVector& rowLengths ) +{ + /**** + * TODO: implement some check here similar to the one in the tridiagonal matrix + */ + return true; +} + template< typename Real, typename Device, typename Index > @@ -116,7 +125,7 @@ bool tnlMultidiagonalMatrix< Real, Device, Index > :: setLike( const tnlMultidia template< typename Real, typename Device, typename Index > -Index tnlMultidiagonalMatrix< Real, Device, Index > :: getNumberOfAllocatedElements() const +Index tnlMultidiagonalMatrix< Real, Device, Index > :: getNumberOfMatrixElements() const { return this->values.getSize(); } @@ -131,22 +140,6 @@ void tnlMultidiagonalMatrix< Real, Device, Index > :: reset() this->values.reset(); } -template< typename Real, - typename Device, - typename Index > -Index tnlMultidiagonalMatrix< Real, Device, Index >::getRows() const -{ - return this->rows; -} - -template< typename Real, - typename Device, - typename Index > -Index tnlMultidiagonalMatrix< Real, Device, Index >::getColumns() const -{ - return this->columns; -} - template< typename Real, typename Device, typename Index > @@ -318,8 +311,7 @@ template< typename Real, typename Index > bool tnlMultidiagonalMatrix< Real, Device, Index >::save( tnlFile& file ) const { - if( ! file.write( &this->rows ) ) return false; - if( ! file.write( &this->columns ) ) return false; + if( ! tnlMatrix< Real, Device, Index >::save( file ) ) return false; if( ! this->values.save( file ) ) return false; if( ! this->diagonalsShift.save( file ) ) return false; return true; @@ -330,8 +322,7 @@ template< typename Real, typename Index > bool tnlMultidiagonalMatrix< Real, Device, Index >::load( tnlFile& file ) { - if( ! file.read( &this->rows ) ) return false; - if( ! file.read( &this->columns ) ) return false; + if( ! tnlMatrix< Real, Device, Index >::load( file ) ) return false; if( ! this->values.load( file ) ) return false; if( ! this->diagonalsShift.load( file ) ) return false; return true; @@ -364,7 +355,7 @@ void tnlMultidiagonalMatrix< Real, Device, Index >::print( ostream& str ) const for( IndexType i = 0; i < this->diagonalsShift.getSize(); i++ ) { const IndexType column = row + diagonalsShift[ i ]; - if( column >=0 && columns < this-columns ) + if( column >=0 && column < this->columns ) str << " Col:" << column << "->" << this->operator()( row, column ) << "\t"; } str << endl; diff --git a/src/implementation/matrices/tnlSlicedEllpackMatrix_impl.h b/src/implementation/matrices/tnlSlicedEllpackMatrix_impl.h index 55e2fc7e4969b55d2da0b9e2baca0079f111fe26..66553809b4a91754b7ddced900df891edeb1c368 100644 --- a/src/implementation/matrices/tnlSlicedEllpackMatrix_impl.h +++ b/src/implementation/matrices/tnlSlicedEllpackMatrix_impl.h @@ -27,8 +27,6 @@ template< typename Real, typename Index, int SliceSize > tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::tnlSlicedEllpackMatrix() -: rows( 0 ), - columns( 0 ) { }; @@ -64,17 +62,14 @@ bool tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::setDimensions( co tnlAssert( rows > 0 && columns > 0, cerr << "rows = " << rows << " columns = " << columns << endl ); - this->rows = rows; - this->columns = columns; - return true; + return tnlSparseMatrix< Real, Device, Index >::setDimensions( rows, columns ); } template< typename Real, typename Device, typename Index, int SliceSize > - template< typename Vector > -bool tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::setRowLengths( const Vector& rowLengths ) +bool tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::setRowLengths( const RowLengthsVector& rowLengths ) { const IndexType slices = roundUpDivision( this->rows, SliceSize ); if( ! this->sliceRowLengths.setSize( slices ) || @@ -107,14 +102,7 @@ bool tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::setRowLengths( co this->slicePointers.setElement( slices, 0 ); this->slicePointers.computeExclusivePrefixSum(); - /**** - * Allocate values and column indexes - */ - if( ! this->values.setSize( this->slicePointers[ slices ] ) || - ! this->columnIndexes.setSize( this->slicePointers[ slices ] ) ) - return false; - this->columnIndexes.setValue( this->columns ); - return true; + return this->allocateMatrixElements( this->slicePointers[ slices ] ); } template< typename Real, @@ -126,56 +114,24 @@ template< typename Real, typename Index2 > bool tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::setLike( const tnlSlicedEllpackMatrix< Real2, Device2, Index2, SliceSize >& matrix ) { - if( ! this->setDimensions( matrix.getRows(), matrix.getColumns() ) || - ! this->values.setLike( matrix.values ) || - ! this->columnIndexes.setLike( matrix.columnIndexes ) || + if( !tnlSparseMatrix< Real, Device, Index >::setLike( matrix ) || ! this->slicePointers.setLike( matrix.slicePointers ) || ! this->sliceRowLengths.setLike( matrix.sliceRowLengths ) ) return false; return true; } -template< typename Real, - typename Device, - typename Index, - int SliceSize > -Index tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::getNumberOfAllocatedElements() const -{ - return this->values.getSize(); -} - template< typename Real, typename Device, typename Index, int SliceSize > void tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::reset() { - this->columns = 0; - this->rows = 0; - this->values.reset(); - this->columnIndexes.reset(); + tnlSparseMatrix< Real, Device, Index >::reset(); this->slicePointers.reset(); this->sliceRowLengths.reset(); } -template< typename Real, - typename Device, - typename Index, - int SliceSize > -Index tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::getRows() const -{ - return this->rows; -} - -template< typename Real, - typename Device, - typename Index, - int SliceSize > -Index tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::getColumns() const -{ - return this->columns; -} - template< typename Real, typename Device, typename Index, @@ -233,7 +189,7 @@ Real tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::getElement( const const IndexType rowEnd = elementPtr + rowLength; while( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] < column ) elementPtr++; - if( this->columnIndexes[ elementPtr ] == column ) + if( elementPtr < rowEnd && this->columnIndexes[ elementPtr ] == column ) return this->values[ elementPtr ]; return 0.0; } @@ -391,12 +347,9 @@ template< typename Real, int SliceSize > bool tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::save( tnlFile& file ) const { - if( ! file.write( &this->rows ) ) return false; - if( ! file.write( &this->columns ) ) return false; - if( ! this->values.save( file ) ) return false; - if( ! this->columnIndexes.save( file ) ) return false; - if( ! this->slicePointers.save( file ) ) return false; - if( ! this->sliceRowLengths.save( file ) ) return false; + if( ! tnlSparseMatrix< Real, Device, Index >::save( file ) || + ! this->slicePointers.save( file ) || + ! this->sliceRowLengths.save( file ) ) return true; } @@ -406,12 +359,9 @@ template< typename Real, int SliceSize > bool tnlSlicedEllpackMatrix< Real, Device, Index, SliceSize >::load( tnlFile& file ) { - if( ! file.read( &this->rows ) ) return false; - if( ! file.read( &this->columns ) ) return false; - if( ! this->values.load( file ) ) return false; - if( ! this->columnIndexes.load( file ) ) return false; - if( ! this->slicePointers.load( file ) ) return false; - if( ! this->sliceRowLengths.load( file ) ) return false; + if( ! tnlSparseMatrix< Real, Device, Index >::load( file ) || + ! this->slicePointers.load( file ) || + ! this->sliceRowLengths.load( file ) ) return true; } diff --git a/src/implementation/matrices/tnlSparseMatrix_impl.h b/src/implementation/matrices/tnlSparseMatrix_impl.h new file mode 100644 index 0000000000000000000000000000000000000000..eb71680da0a9b202994960d154069e9785c6adfc --- /dev/null +++ b/src/implementation/matrices/tnlSparseMatrix_impl.h @@ -0,0 +1,102 @@ +/*************************************************************************** + tnlSparseMatrix_impl.h - description + ------------------- + begin : Dec 21, 2013 + copyright : (C) 2013 by Tomas Oberhuber + email : tomas.oberhuber@fjfi.cvut.cz + ***************************************************************************/ + +/*************************************************************************** + * * + * This program is free software; you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation; either version 2 of the License, or * + * (at your option) any later version. * + * * + ***************************************************************************/ + +#ifndef TNLSPARSEMATRIX_IMPL_H_ +#define TNLSPARSEMATRIX_IMPL_H_ + +template< typename Real, + typename Device, + typename Index > +tnlSparseMatrix< Real, Device, Index >::tnlSparseMatrix() +{ +} + +template< typename Real, + typename Device, + typename Index > + template< typename Real2, + typename Device2, + typename Index2 > +bool tnlSparseMatrix< Real, Device, Index >::setLike( const tnlSparseMatrix< Real2, Device2, Index2 >& matrix ) +{ + if( ! tnlMatrix< Real, Device, Index >::setLike( matrix ) || + ! this->allocateMatrixElements( matrix.getNumberOfMatrixElements() ) ) + return false; + return true; +} + +template< typename Real, + typename Device, + typename Index > +Index tnlSparseMatrix< Real, Device, Index >::getNumberOfMatrixElements() const +{ + return this->values.getSize(); +} + +template< typename Real, + typename Device, + typename Index > +void tnlSparseMatrix< Real, Device, Index >::reset() +{ + tnlMatrix< Real, Device, Index >::reset(); + this->values.reset(); + this->columnIndexes.reset(); +} + +template< typename Real, + typename Device, + typename Index > +bool tnlSparseMatrix< Real, Device, Index >::save( tnlFile& file ) const +{ + if( ! tnlMatrix< Real, Device, Index >::save( file ) || + ! this->values.save( file ) || + ! this->columnIndexes.save( file ) ) + return false; + return true; +} + +template< typename Real, + typename Device, + typename Index > +bool tnlSparseMatrix< Real, Device, Index >::load( tnlFile& file ) +{ + if( ! tnlMatrix< Real, Device, Index >::load( file ) || + ! this->values.load( file ) || + ! this->columnIndexes.load( file ) ) + return false; + return true; +} + +template< typename Real, + typename Device, + typename Index > +bool tnlSparseMatrix< Real, Device, Index >::allocateMatrixElements( const IndexType& numberOfMatrixElements ) +{ + if( ! this->values.setSize( numberOfMatrixElements ) || + ! this->columnIndexes.setSize( numberOfMatrixElements ) ) + return false; + + /**** + * Setting a column index to this->columns means that the + * index is undefined. + */ + this->columnIndexes.setValue( this->columns ); + return true; +} + + +#endif /* TNLSPARSEMATRIX_IMPL_H_ */ diff --git a/src/implementation/matrices/tnlTridiagonalMatrix_impl.h b/src/implementation/matrices/tnlTridiagonalMatrix_impl.h index 704e35988b229f5c9ed63f81b609319a4734f573..bd051a3a83c295dd2afb1d598175e4e8d3e05ce1 100644 --- a/src/implementation/matrices/tnlTridiagonalMatrix_impl.h +++ b/src/implementation/matrices/tnlTridiagonalMatrix_impl.h @@ -25,7 +25,6 @@ template< typename Real, typename Device, typename Index > tnlTridiagonalMatrix< Real, Device, Index >::tnlTridiagonalMatrix() -: rows( 0 ) { } @@ -51,55 +50,63 @@ tnlString tnlTridiagonalMatrix< Real, Device, Index >::getTypeVirtual() const template< typename Real, typename Device, typename Index > -bool tnlTridiagonalMatrix< Real, Device, Index >::setDimensions( const IndexType rows ) +bool tnlTridiagonalMatrix< Real, Device, Index >::setDimensions( const IndexType rows, + const IndexType columns ) { - if( ! values.setSize( 3*rows - 2 ) ) + if( ! tnlMatrix< Real, Device, Index >::setDimensions( rows, columns ) ) + return false; + if( ! values.setSize( 3*Min( rows, columns ) ) ) return false; - this->rows = rows; - this->columns = rows; this->values.setValue( 0.0 ); } template< typename Real, typename Device, typename Index > - template< typename Real2, typename Device2, typename Index2 > -bool tnlTridiagonalMatrix< Real, Device, Index >::setLike( const tnlTridiagonalMatrix< Real2, Device2, Index2 >& m ) +bool tnlTridiagonalMatrix< Real, Device, Index >::setRowLengths( const RowLengthsVector& rowLengths ) { - return this->setDimensions( m.getRows() ); -} - -template< typename Real, - typename Device, - typename Index > -Index tnlTridiagonalMatrix< Real, Device, Index >::getNumberOfAllocatedElements() const -{ - return 3 * rows - 2; + if( rowLengths[ 0 ] > 2 ) + return false; + const IndexType diagonalLength = Min( this->getRows(), this->getColumns() ); + for( Index i = 1; i < diagonalLength-1; i++ ) + if( rowLengths[ i ] > 3 ) + return false; + if( this->getRows() > this->getColumns() ) + if( rowLengths[ this->getRows()-1 ] > 1 ) + return false; + if( this->getRows() == this->getColumns() ) + if( rowLengths[ this->getRows()-1 ] > 2 ) + return false; + if( this->getRows() < this->getColumns() ) + if( rowLengths[ this->getRows()-1 ] > 3 ) + return false; + return true; } template< typename Real, typename Device, typename Index > -void tnlTridiagonalMatrix< Real, Device, Index >::reset() + template< typename Real2, typename Device2, typename Index2 > +bool tnlTridiagonalMatrix< Real, Device, Index >::setLike( const tnlTridiagonalMatrix< Real2, Device2, Index2 >& m ) { - this->rows = 0; - this->values.reset(); + return this->setDimensions( m.getRows(), m.getColumns() ); } template< typename Real, typename Device, typename Index > -Index tnlTridiagonalMatrix< Real, Device, Index >::getRows() const +Index tnlTridiagonalMatrix< Real, Device, Index >::getNumberOfMatrixElements() const { - return this->rows; + return 3 * Min( this->getRows(), this->getColumns() ); } template< typename Real, typename Device, typename Index > -Index tnlTridiagonalMatrix< Real, Device, Index >::getColumns() const +void tnlTridiagonalMatrix< Real, Device, Index >::reset() { - return this->rows; + tnlMatrix< Real, Device, Index >::reset(); + this->values.reset(); } template< typename Real, @@ -131,11 +138,12 @@ void tnlTridiagonalMatrix< Real, Device, Index >::setValue( const RealType& v ) template< typename Real, typename Device, typename Index > -void tnlTridiagonalMatrix< Real, Device, Index >::setElement( const IndexType row, +bool tnlTridiagonalMatrix< Real, Device, Index >::setElement( const IndexType row, const IndexType column, const RealType& value ) { this->values.setElement( this->getElementIndex( row, column ), value ); + return true; } template< typename Real, @@ -332,7 +340,7 @@ template< typename Real, typename Index > bool tnlTridiagonalMatrix< Real, Device, Index >::save( tnlFile& file ) const { - if( ! file.write( &this->rows ) || + if( ! tnlMatrix< Real, Device, Index >::save( file ) || ! this->values.save( file ) ) { cerr << "Unable to save the tridiagonal matrix " << this->getName() << "." << endl; @@ -346,7 +354,7 @@ template< typename Real, typename Index > bool tnlTridiagonalMatrix< Real, Device, Index >::load( tnlFile& file ) { - if( ! file.read( &this->rows ) || + if( ! tnlMatrix< Real, Device, Index >::load( file ) || ! this->values.load( file ) ) { cerr << "Unable to save the tridiagonal matrix " << this->getName() << "." << endl; diff --git a/src/legacy/matrices/tnlCusparseCSRMatrix.h b/src/legacy/matrices/tnlCusparseCSRMatrix.h index f40d59ef298ea1874a9d327c140d439e1f0df699..80b1b4d82fc095fdde8f893bb968b2fe5f47fe2e 100644 --- a/src/legacy/matrices/tnlCusparseCSRMatrix.h +++ b/src/legacy/matrices/tnlCusparseCSRMatrix.h @@ -26,7 +26,6 @@ #include <core/tnlAssert.h> #include <core/mfuncs.h> #include <matrices/tnlCSRMatrix.h> -#include <matrices/tnlRgCSRMatrix.h> #include <debug/tnlDebug.h> #ifdef HAVE_CUSPARSE diff --git a/src/matrices/CMakeLists.txt b/src/matrices/CMakeLists.txt index 0cdd3036542b1ebef0a5b90d29a51f25cfb2399b..00caba1fb8f49b6b566934af9eae8648618a46de 100755 --- a/src/matrices/CMakeLists.txt +++ b/src/matrices/CMakeLists.txt @@ -2,6 +2,7 @@ SET( headers tnlMatrix.h tnlDenseMatrix.h tnlTridiagonalMatrix.h tnlMultidiagonalMatrix.h + tnlSparseMatrix.h tnlEllpackMatrix.h tnlSlicedEllpackMatrix.h tnlChunkedEllpackMatrix.h diff --git a/src/matrices/tnlCSRMatrix.h b/src/matrices/tnlCSRMatrix.h index 3382173a0dd7f3497e432da540cc77c4dcc8f76c..11fdffd9d1aac831f7b8637b9287800354bd0774 100644 --- a/src/matrices/tnlCSRMatrix.h +++ b/src/matrices/tnlCSRMatrix.h @@ -18,16 +18,18 @@ #ifndef TNLCSRMATRIX_H_ #define TNLCSRMATRIX_H_ +#include <matrices/tnlSparseMatrix.h> #include <core/vectors/tnlVector.h> template< typename Real, typename Device = tnlHost, typename Index = int > -class tnlCSRMatrix : public tnlObject +class tnlCSRMatrix : public tnlSparseMatrix< Real, Device, Index > { public: typedef Real RealType; typedef Device DeviceType; typedef Index IndexType; + typedef typename tnlSparseMatrix< RealType, DeviceType, IndexType >:: RowLengthsVector RowLengthsVector; tnlCSRMatrix(); @@ -38,20 +40,13 @@ class tnlCSRMatrix : public tnlObject bool setDimensions( const IndexType rows, const IndexType columns ); - template< typename Vector > - bool setRowLengths( const Vector& rowLengths ); + bool setRowLengths( const RowLengthsVector& rowLengths ); template< typename Real2, typename Device2, typename Index2 > bool setLike( const tnlCSRMatrix< Real2, Device2, Index2 >& matrix ); - IndexType getNumberOfAllocatedElements() const; - void reset(); - IndexType getRows() const; - - IndexType getColumns() const; - template< typename Real2, typename Device2, typename Index2 > bool operator == ( const tnlCSRMatrix< Real2, Device2, Index2 >& matrix ) const; @@ -110,11 +105,7 @@ class tnlCSRMatrix : public tnlObject protected: - IndexType rows, columns; - - tnlVector< Real, Device, Index > values; - - tnlVector< Index, Device, Index > columnIndexes, rowPointers; + tnlVector< Index, Device, Index > rowPointers; }; diff --git a/src/matrices/tnlChunkedEllpackMatrix.h b/src/matrices/tnlChunkedEllpackMatrix.h index c906b3560b64d2360d6f868368bd8dc401a35111..e78b958166107a9f3e01bfd7aa8b857396337f31 100644 --- a/src/matrices/tnlChunkedEllpackMatrix.h +++ b/src/matrices/tnlChunkedEllpackMatrix.h @@ -18,16 +18,18 @@ #ifndef TNLCHUNKEDELLPACKMATRIX_H_ #define TNLCHUNKEDELLPACKMATRIX_H_ +#include <matrices/tnlSparseMatrix.h> #include <core/vectors/tnlVector.h> template< typename Real, typename Device = tnlHost, typename Index = int > -class tnlChunkedEllpackMatrix : public tnlObject +class tnlChunkedEllpackMatrix : public tnlSparseMatrix< Real, Device, Index > { public: typedef Real RealType; typedef Device DeviceType; typedef Index IndexType; + typedef typename tnlSparseMatrix< RealType, DeviceType, IndexType >:: RowLengthsVector RowLengthsVector; tnlChunkedEllpackMatrix(); @@ -38,20 +40,13 @@ class tnlChunkedEllpackMatrix : public tnlObject bool setDimensions( const IndexType rows, const IndexType columns ); - template< typename Vector > - bool setRowLengths( const Vector& rowLengths ); + bool setRowLengths( const RowLengthsVector& rowLengths ); template< typename Real2, typename Device2, typename Index2 > bool setLike( const tnlChunkedEllpackMatrix< Real2, Device2, Index2 >& matrix ); - IndexType getNumberOfAllocatedElements() const; - void reset(); - IndexType getRows() const; - - IndexType getColumns() const; - template< typename Real2, typename Device2, typename Index2 > bool operator == ( const tnlChunkedEllpackMatrix< Real2, Device2, Index2 >& matrix ) const; @@ -129,13 +124,9 @@ class tnlChunkedEllpackMatrix : public tnlObject { return tnlString( "tnlChunkedEllpackSliceInfo" ); }; }; - IndexType rows, columns; - IndexType chunksInSlice, desiredChunkSize; - tnlVector< Real, Device, Index > values; - - tnlVector< Index, Device, Index > columnIndexes, chunksToRowsMapping, slicesToRowsMapping, rowPointers; + tnlVector< Index, Device, Index > chunksToRowsMapping, slicesToRowsMapping, rowPointers; tnlArray< tnlChunkedEllpackSliceInfo, Device, Index > slices; diff --git a/src/matrices/tnlDenseMatrix.h b/src/matrices/tnlDenseMatrix.h index fbdc7cc755777cd5f864f797014c393ebdbc5c28..0c22c13011387aecd42ec124ea54e2d70407fc9b 100644 --- a/src/matrices/tnlDenseMatrix.h +++ b/src/matrices/tnlDenseMatrix.h @@ -19,18 +19,21 @@ #define TNLDENSEMATRIX_H_ #include <core/tnlHost.h> +#include <matrices/tnlMatrix.h> #include <core/arrays/tnlMultiArray.h> template< typename Real = double, typename Device = tnlHost, typename Index = int > -class tnlDenseMatrix : public tnlMultiArray< 2, Real, Device, Index > +class tnlDenseMatrix : public tnlMatrix< Real, Device, Index >, + public tnlMultiArray< 2, Real, Device, Index > { public: typedef Real RealType; typedef Device DeviceType; typedef Index IndexType; + typedef typename tnlMatrix< Real, Device, Index >::RowLengthsVector RowLengthsVector; tnlDenseMatrix(); @@ -41,17 +44,24 @@ class tnlDenseMatrix : public tnlMultiArray< 2, Real, Device, Index > bool setDimensions( const IndexType rows, const IndexType columns ); + template< typename Real2, typename Device2, typename Index2 > + bool setLike( const tnlDenseMatrix< Real2, Device2, Index2 >& matrix ); + /**** * This method is only for the compatibility with the sparse matrices. */ - template< typename Vector > - bool setRowLengths( const Vector& rowLengths ); + bool setRowLengths( const RowLengthsVector& rowLengths ); - IndexType getNumberOfAllocatedElements() const; + IndexType getNumberOfMatrixElements() const; - IndexType getRows() const; + void reset(); + + bool setElement( const IndexType row, + const IndexType column, + const RealType& value ); - IndexType getColumns() const; + Real getElement( const IndexType row, + const IndexType column ) const; bool addElement( const IndexType row, const IndexType column, @@ -98,6 +108,14 @@ class tnlDenseMatrix : public tnlMultiArray< 2, Real, Device, Index > Vector& x, const RealType& omega = 1.0 ) const; + bool save( const tnlString& fileName ) const; + + bool load( const tnlString& fileName ); + + bool save( tnlFile& file ) const; + + bool load( tnlFile& file ); + void print( ostream& str ) const; }; diff --git a/src/matrices/tnlEllpackMatrix.h b/src/matrices/tnlEllpackMatrix.h index ae848ad19e57972a9382d5617f1c86d417cc22ff..194b3971e6881716f75e28102dcbc339501fed56 100644 --- a/src/matrices/tnlEllpackMatrix.h +++ b/src/matrices/tnlEllpackMatrix.h @@ -18,16 +18,18 @@ #ifndef TNLELLPACKMATRIX_H_ #define TNLELLPACKMATRIX_H_ +#include <matrices/tnlSparseMatrix.h> #include <core/vectors/tnlVector.h> template< typename Real, typename Device = tnlHost, typename Index = int > -class tnlEllpackMatrix : public tnlObject +class tnlEllpackMatrix : public tnlSparseMatrix< Real, Device, Index > { public: typedef Real RealType; typedef Device DeviceType; typedef Index IndexType; + typedef typename tnlSparseMatrix< RealType, DeviceType, IndexType >:: RowLengthsVector RowLengthsVector; tnlEllpackMatrix(); @@ -38,22 +40,15 @@ class tnlEllpackMatrix : public tnlObject bool setDimensions( const IndexType rows, const IndexType columns ); - template< typename Vector > - bool setRowLengths( const Vector& rowLengths ); + bool setRowLengths( const RowLengthsVector& rowLengths ); bool setConstantRowLengths( const IndexType& rowLengths ); template< typename Real2, typename Device2, typename Index2 > bool setLike( const tnlEllpackMatrix< Real2, Device2, Index2 >& matrix ); - IndexType getNumberOfAllocatedElements() const; - void reset(); - IndexType getRows() const; - - IndexType getColumns() const; - template< typename Real2, typename Device2, typename Index2 > bool operator == ( const tnlEllpackMatrix< Real2, Device2, Index2 >& matrix ) const; @@ -110,12 +105,7 @@ class tnlEllpackMatrix : public tnlObject bool allocateElements(); - IndexType rows, columns, rowLengths; - - tnlVector< Real, Device, Index > values; - - tnlVector< Index, Device, Index > columnIndexes; - + IndexType rowLengths; }; #include <implementation/matrices/tnlEllpackMatrix_impl.h> diff --git a/src/matrices/tnlMatrix.h b/src/matrices/tnlMatrix.h index 32a0369eb890887c21d53e771ffe0d6885a7cbf0..498c18289c03500d48134cd53efc9e953476986e 100644 --- a/src/matrices/tnlMatrix.h +++ b/src/matrices/tnlMatrix.h @@ -20,49 +20,60 @@ #include <core/tnlObject.h> #include <core/tnlHost.h> +#include <core/vectors/tnlVector.h> template< typename Real = double, typename Device = tnlHost, typename Index = int > -class tnlMatrix : public tnlObject +class tnlMatrix : public virtual tnlObject { public: typedef Real RealType; typedef Device DeviceType; typedef Index IndexType; + typedef tnlVector< IndexType, DeviceType, IndexType > RowLengthsVector; + + tnlMatrix(); virtual bool setDimensions( const IndexType rows, const IndexType columns ); + virtual bool setRowLengths( const RowLengthsVector& rowLengths ) = 0; + + template< typename Real2, typename Device2, typename Index2 > + bool setLike( const tnlMatrix< Real2, Device2, Index2 >& matrix ); + + virtual IndexType getNumberOfMatrixElements() const = 0; + + void reset(); + IndexType getRows() const; IndexType getColumns() const; virtual bool setElement( const IndexType row, const IndexType column, - const RealType& value ); + const RealType& value ) = 0; virtual bool addElement( const IndexType row, const IndexType column, const RealType& value, - const RealType& thisElementMultiplicator = 1.0 ); + const RealType& thisElementMultiplicator = 1.0 ) = 0; + virtual Real getElement( const IndexType row, + const IndexType column ) const = 0; virtual bool save( tnlFile& file ) const; virtual bool load( tnlFile& file ); - virtual bool save( const tnlString& fileName ) const; - - virtual bool load( const tnlString& fileName ); - protected: IndexType rows, columns; }; -#include <implementation/matrices/tnlMatrix.h> +#include <implementation/matrices/tnlMatrix_impl.h> #endif /* TNLMATRIX_H_ */ diff --git a/src/matrices/tnlMatrixReader.h b/src/matrices/tnlMatrixReader.h index 3d3b8d401f4a7fd6544cd27c4e54996c3088d138..5f2ee48cc546476b6aa2484b88bf0598bf7a8a02 100644 --- a/src/matrices/tnlMatrixReader.h +++ b/src/matrices/tnlMatrixReader.h @@ -37,6 +37,12 @@ class tnlMatrixReader static bool verifyMtxFile( std::istream& file, const Matrix& matrix, bool verbose = false ); + + static bool findLineByElement( std::istream& file, + const IndexType& row, + const IndexType& column, + tnlString& line, + IndexType& lineNumber ); protected: static bool checkMtxHeader( const tnlString& header, diff --git a/src/matrices/tnlMultidiagonalMatrix.h b/src/matrices/tnlMultidiagonalMatrix.h index 09d91bc96bffbee219304668f1a3abb9d8fe5e40..f853096173392af8e90574531470433be4c16b95 100644 --- a/src/matrices/tnlMultidiagonalMatrix.h +++ b/src/matrices/tnlMultidiagonalMatrix.h @@ -18,16 +18,18 @@ #ifndef TNLMULTIDIAGONALMATRIX_H_ #define TNLMULTIDIAGONALMATRIX_H_ +#include <matrices/tnlMatrix.h> #include <core/vectors/tnlVector.h> template< typename Real, typename Device = tnlHost, typename Index = int > -class tnlMultidiagonalMatrix : public tnlObject +class tnlMultidiagonalMatrix : public tnlMatrix< Real, Device, Index > { public: typedef Real RealType; typedef Device DeviceType; typedef Index IndexType; + typedef typename tnlMatrix< Real, Device, Index >::RowLengthsVector RowLengthsVector; tnlMultidiagonalMatrix(); @@ -38,6 +40,8 @@ class tnlMultidiagonalMatrix : public tnlObject bool setDimensions( const IndexType rows, const IndexType columns ); + bool setRowLengths( const RowLengthsVector& rowLengths ); + bool setDiagonals( const IndexType diagonalsNumber, const IndexType* diagonalsShift ); @@ -46,14 +50,10 @@ class tnlMultidiagonalMatrix : public tnlObject template< typename Real2, typename Device2, typename Index2 > bool setLike( const tnlMultidiagonalMatrix< Real2, Device2, Index2 >& matrix ); - IndexType getNumberOfAllocatedElements() const; + IndexType getNumberOfMatrixElements() const; void reset(); - IndexType getRows() const; - - IndexType getColumns() const; - template< typename Real2, typename Device2, typename Index2 > bool operator == ( const tnlMultidiagonalMatrix< Real2, Device2, Index2 >& matrix ) const; @@ -109,8 +109,6 @@ class tnlMultidiagonalMatrix : public tnlObject const IndexType column, IndexType& index ) const; - IndexType rows, columns; - tnlVector< Real, Device, Index > values; tnlVector< Index, Device, Index > diagonalsShift; diff --git a/src/matrices/tnlSlicedEllpackMatrix.h b/src/matrices/tnlSlicedEllpackMatrix.h index 49de206492fe0e11e666784b2ef8675d1df2e086..17ad1b95086372f46745d47491ce86fb30df60ae 100644 --- a/src/matrices/tnlSlicedEllpackMatrix.h +++ b/src/matrices/tnlSlicedEllpackMatrix.h @@ -18,16 +18,21 @@ #ifndef TNLSLICEDELLPACKMATRIX_H_ #define TNLSLICEDELLPACKMATRIX_H_ +#include <matrices/tnlSparseMatrix.h> #include <core/vectors/tnlVector.h> -template< typename Real, typename Device = tnlHost, typename Index = int, int SliceSize = 32 > -class tnlSlicedEllpackMatrix : public tnlObject +template< typename Real = double, + typename Device = tnlHost, + typename Index = int, + int SliceSize = 32 > +class tnlSlicedEllpackMatrix : public tnlSparseMatrix< Real, Device, Index > { public: typedef Real RealType; typedef Device DeviceType; typedef Index IndexType; + typedef typename tnlSparseMatrix< RealType, DeviceType, IndexType >:: RowLengthsVector RowLengthsVector; tnlSlicedEllpackMatrix(); @@ -38,20 +43,13 @@ class tnlSlicedEllpackMatrix : public tnlObject bool setDimensions( const IndexType rows, const IndexType columns ); - template< typename Vector > - bool setRowLengths( const Vector& rowLengths ); + bool setRowLengths( const RowLengthsVector& rowLengths ); template< typename Real2, typename Device2, typename Index2 > bool setLike( const tnlSlicedEllpackMatrix< Real2, Device2, Index2, SliceSize >& matrix ); - IndexType getNumberOfAllocatedElements() const; - void reset(); - IndexType getRows() const; - - IndexType getColumns() const; - template< typename Real2, typename Device2, typename Index2 > bool operator == ( const tnlSlicedEllpackMatrix< Real2, Device2, Index2 >& matrix ) const; @@ -106,11 +104,7 @@ class tnlSlicedEllpackMatrix : public tnlObject protected: - IndexType rows, columns; - - tnlVector< Real, Device, Index > values; - - tnlVector< Index, Device, Index > columnIndexes, slicePointers, sliceRowLengths; + tnlVector< Index, Device, Index > slicePointers, sliceRowLengths; }; diff --git a/src/matrices/tnlSparseMatrix.h b/src/matrices/tnlSparseMatrix.h new file mode 100644 index 0000000000000000000000000000000000000000..b9ccc26e0c1e6f96796ade6dc59c0704e63f1718 --- /dev/null +++ b/src/matrices/tnlSparseMatrix.h @@ -0,0 +1,59 @@ +/*************************************************************************** + tnlSparseMatrix.h - description + ------------------- + begin : Dec 21, 2013 + copyright : (C) 2013 by Tomas Oberhuber + email : tomas.oberhuber@fjfi.cvut.cz + ***************************************************************************/ + +/*************************************************************************** + * * + * This program is free software; you can redistribute it and/or modify * + * it under the terms of the GNU General Public License as published by * + * the Free Software Foundation; either version 2 of the License, or * + * (at your option) any later version. * + * * + ***************************************************************************/ + +#ifndef TNLSPARSEMATRIX_H_ +#define TNLSPARSEMATRIX_H_ + +#include <matrices/tnlMatrix.h> + +template< typename Real, + typename Device, + typename Index > +class tnlSparseMatrix : public tnlMatrix< Real, Device, Index > +{ + public: + + typedef Real RealType; + typedef Device DeviceType; + typedef Index IndexType; + typedef typename tnlMatrix< RealType, DeviceType, IndexType >:: RowLengthsVector RowLengthsVector; + + tnlSparseMatrix(); + + template< typename Real2, typename Device2, typename Index2 > + bool setLike( const tnlSparseMatrix< Real2, Device2, Index2 >& matrix ); + + IndexType getNumberOfMatrixElements() const; + + void reset(); + + bool save( tnlFile& file ) const; + + bool load( tnlFile& file ); + + protected: + + bool allocateMatrixElements( const IndexType& numberOfMatrixElements ); + + tnlVector< Index, Device, Index > columnIndexes; + + tnlVector< Real, Device, Index > values; +}; + +#include <implementation/matrices/tnlSparseMatrix_impl.h> + +#endif /* TNLSPARSEMATRIX_H_ */ diff --git a/src/matrices/tnlTridiagonalMatrix.h b/src/matrices/tnlTridiagonalMatrix.h index 1b6bb1472808f6025428289ec9d6c28fb300e8f4..5d7328d229ec644feb4cb1b3ea9cdd08392236be 100644 --- a/src/matrices/tnlTridiagonalMatrix.h +++ b/src/matrices/tnlTridiagonalMatrix.h @@ -18,20 +18,20 @@ #ifndef TNLTRIDIAGONALMATRIX_H_ #define TNLTRIDIAGONALMATRIX_H_ -#include <core/tnlObject.h> -#include <core/tnlHost.h> +#include <matrices/tnlMatrix.h> #include <core/vectors/tnlVector.h> template< typename Real = double, typename Device = tnlHost, typename Index = int > -class tnlTridiagonalMatrix : public tnlObject +class tnlTridiagonalMatrix : public tnlMatrix< Real, Device, Index > { public: typedef Real RealType; typedef Device DeviceType; typedef Index IndexType; + typedef typename tnlMatrix< Real, Device, Index >::RowLengthsVector RowLengthsVector; tnlTridiagonalMatrix(); @@ -39,19 +39,18 @@ class tnlTridiagonalMatrix : public tnlObject tnlString getTypeVirtual() const; - bool setDimensions( const IndexType rows ); + bool setDimensions( const IndexType rows, + const IndexType columns ); + + bool setRowLengths( const RowLengthsVector& rowLengths ); template< typename Real2, typename Device2, typename Index2 > bool setLike( const tnlTridiagonalMatrix< Real2, Device2, Index2 >& m ); - IndexType getNumberOfAllocatedElements() const; + IndexType getNumberOfMatrixElements() const; void reset(); - IndexType getRows() const; - - IndexType getColumns() const; - template< typename Real2, typename Device2, typename Index2 > bool operator == ( const tnlTridiagonalMatrix< Real2, Device2, Index2 >& matrix ) const; @@ -60,7 +59,7 @@ class tnlTridiagonalMatrix : public tnlObject void setValue( const RealType& v ); - void setElement( const IndexType row, + bool setElement( const IndexType row, const IndexType column, const RealType& value ); @@ -127,8 +126,6 @@ class tnlTridiagonalMatrix : public tnlObject IndexType getElementIndex( const IndexType row, const IndexType column ) const; - IndexType rows, columns; - tnlVector< RealType, DeviceType, IndexType > values; }; diff --git a/tests/benchmarks/CMakeLists.txt b/tests/benchmarks/CMakeLists.txt index 462ad789e9cc5e0b807389d3593e4265dad4b7af..19e8a854bbcf4f74cb60a8dffa288a8c4977a6a4 100755 --- a/tests/benchmarks/CMakeLists.txt +++ b/tests/benchmarks/CMakeLists.txt @@ -9,13 +9,13 @@ SET( tnlSpmvBenchmark_headers sparse-matrix-benchmark.h tnlSpmvBenchmarkRgCSRMatrix.h ) -IF( BUILD_CUDA ) - #CUDA_ADD_EXECUTABLE( tnl-sparse-matrix-benchmark${debugExt} sparse-matrix-benchmark.cu ) - #SET_TARGET_PROPERTIES( tnl-sparse-matrix-benchmark${debugExt} PROPERTIES CUDA_COMPILE_FLAGS "${CXX_OPTIMIZE_FLAGS}" ) -ELSE() - #ADD_EXECUTABLE( tnl-sparse-matrix-benchmark${debugExt} sparse-matrix-benchmark.cpp ) - #SET_TARGET_PROPERTIES( tnl-sparse-matrix-benchmark${debugExt} PROPERTIES COMPILE_FLAGS "${CXX_OPTIMIZE_FLAGS}" ) -ENDIF() +#IF( BUILD_CUDA ) +# CUDA_ADD_EXECUTABLE( tnl-sparse-matrix-benchmark${debugExt} sparse-matrix-benchmark.cu ) +# SET_TARGET_PROPERTIES( tnl-sparse-matrix-benchmark${debugExt} PROPERTIES CUDA_COMPILE_FLAGS "${CXX_OPTIMIZE_FLAGS}" ) +#ELSE() +# ADD_EXECUTABLE( tnl-sparse-matrix-benchmark${debugExt} sparse-matrix-benchmark.cpp ) +# SET_TARGET_PROPERTIES( tnl-sparse-matrix-benchmark${debugExt} PROPERTIES COMPILE_FLAGS "${CXX_OPTIMIZE_FLAGS}" ) +#ENDIF() #TARGET_LINK_LIBRARIES( tnl-sparse-matrix-benchmark${debugExt} tnl${debugExt}-${tnlVersion} # ${CUSPARSE_LIBRARY} ) diff --git a/tests/benchmarks/matrix-solvers-benchmark.h b/tests/benchmarks/matrix-solvers-benchmark.h index 213e7ef01ff0946b93b8bf04b117ff409377dae6..d8f0d9ee7fee139540001cc9b335b92bb98f4ad7 100644 --- a/tests/benchmarks/matrix-solvers-benchmark.h +++ b/tests/benchmarks/matrix-solvers-benchmark.h @@ -72,7 +72,7 @@ bool benchmarkSolver( const tnlParameterContainer& parameters, const RealType& maxResidue = parameters. GetParameter< double >( "max-residue" ); const IndexType& size = matrix. getRows(); - const IndexType nonZeros = matrix. getNumberOfAllocatedElements(); + const IndexType nonZeros = matrix. getNumberOfMatrixElements(); //const IndexType maxIterations = size * ( ( double ) size * size / ( double ) nonZeros ); const IndexType maxIterations = size; cout << "Setting max. number of iterations to " << maxIterations << endl; @@ -277,7 +277,7 @@ bool benchmarkMatrix( const tnlParameterContainer& parameters ) return false; } matrixStatsFile << " <td> " << csrMatrix. getRows() << " </td> " << endl - << " <td> " << csrMatrix. getNumberOfAllocatedElements() << " </td> " << endl; + << " <td> " << csrMatrix. getNumberOfMatrixElements() << " </td> " << endl; matrixStatsFile. close(); } diff --git a/tests/benchmarks/sparse-matrix-benchmark.h b/tests/benchmarks/sparse-matrix-benchmark.h index 6b5743233f7666f95a44089238c38d3024bb1c8b..104ff9a294094b8699dff29fe48612b758cdcecb 100644 --- a/tests/benchmarks/sparse-matrix-benchmark.h +++ b/tests/benchmarks/sparse-matrix-benchmark.h @@ -23,12 +23,11 @@ #include <iomanip> #include <config/tnlConfigDescription.h> #include <config/tnlParameterContainer.h> -#include <legacy/matrices/tnlFullMatrix.h> -#include <legacy/matrices/tnlFastCSRMatrix.h> -#include <legacy/matrices/tnlFastRgCSRMatrix.h> -#include <legacy/matrices/tnlFastRgCSRMatrixCUDA.h> -#include <legacy/matrices/tnlEllpackMatrix.h> -#include <legacy/matrices/tnlEllpackMatrixCUDA.h> +#include <matrices/tnlDenseMatrix.h> +#include <matrices/tnlEllpackMatrix.h> +#include <matrices/tnlSlicedEllpackMatrix.h> +#include <matrices/tnlChunkedEllpackMatrix.h> +#include <matrices/tnlCSRMatrix.h> #include <core/mfuncs.h> #include "tnlSpmvBenchmarkCSRMatrix.h" #include "tnlSpmvBenchmarkCusparseCSRMatrix.h" @@ -102,38 +101,25 @@ void benchmarkRgCSRFormat( const tnlCSRMatrix< Real, tnlHost, int >& csrMatrix, } } -template< class Real > -bool benchmarkMatrix( const tnlString& inputFile, - const tnlString& inputMtxFile, - const tnlString& pdfFile, - const tnlString& logFileName, - bool formatTest, - int maxIterations, - int verbose ) +template< typename RealType > +bool benchmarkMatrix( const tnlParameterContainer& parameters ) { /**** * Read the CSR matrix ... */ - tnlCSRMatrix< Real > csrMatrix; - csrMatrix.setName( "csr-matrix" ); - tnlString inputMtxSortedFile( inputMtxFile ); - inputMtxSortedFile += tnlString( ".sort" ); - tnlFile binaryFile; - if( ! binaryFile. open( inputFile, tnlReadMode ) ) - { - cerr << "I am not able to open the file " << inputFile << "." << endl; - return 1; - } - if( verbose ) - cout << "Reading the CSR matrix ... " << flush; - if( ! csrMatrix. load( binaryFile ) ) + typedef tnlCSRMatrix< RealType, tnlHost, int > CsrMatrix; + CsrMatrix csrMatrix; + + const tnlString& inputFileName = parameters.GetParameter< tnlString >( "input-file" ); + fstream inputFile; + inputFile.open( inputFileName.getString(), ios::in ); + if( ! inputFile ) { - cerr << "Unable to restore the CSR matrix." << endl; + cerr << "I am not able to open the file " << inputFileName << "." << endl; return false; } - if( verbose ) - cout << " OK." << endl; - binaryFile. close(); + if( ! tnlMatrixReader< CsrMatrix >::readMtxFile( inputFile, csrMatrix ) ) + return false; /**** * Check the number of the non-zero elements @@ -436,57 +422,14 @@ int main( int argc, char* argv[] ) conf_desc. PrintUsage( argv[ 0 ] ); return 1; } - tnlString inputFile = parameters. GetParameter< tnlString >( "input-file" ); - tnlString inputMtxFile = parameters. GetParameter< tnlString >( "input-mtx-file" ); - tnlString pdfFile = parameters. GetParameter< tnlString >( "pdf-file" ); - tnlString logFileName = parameters. GetParameter< tnlString >( "log-file" ); - double stop_time = parameters. GetParameter< double >( "stop-time" ); - bool formatTest = parameters. GetParameter< bool >( "format-test" ); - int maxIterations = parameters. GetParameter< int >( "max-iterations" ); - int verbose = parameters. GetParameter< int >( "verbose"); - - - tnlFile binaryFile; - if( ! binaryFile. open( inputFile, tnlReadMode ) ) - { - cerr << "I am not able to open the file " << inputFile << "." << endl; - return 1; - } - tnlString object_type; - if( ! getObjectType( binaryFile, object_type ) ) - { - cerr << "Unknown object ... SKIPPING!" << endl; - return EXIT_FAILURE; - } - if( verbose ) - cout << object_type << " detected ... " << endl; - binaryFile. close(); - - if( object_type == "tnlCSRMatrix< float, tnlHost >") - benchmarkMatrix< float >( inputFile, - inputMtxFile, - pdfFile, - logFileName, - formatTest, - maxIterations, - verbose ); - - if( object_type == "tnlCSRMatrix< double, tnlHost >" ) - { - benchmarkMatrix< double >( inputFile, - inputMtxFile, - pdfFile, - logFileName, - formatTest, - maxIterations, - verbose ); - } - //binaryFile. close(); - - - + const tnlString& precision = parameters.GetParameter< tnlString >( "precision" ); + if( precision == "float" ) + if( ! benchmarkMatrix< float >( parameters ) ) + return EXIT_FAILURE; + if( precision == "double" ) + if( ! benchmarkMatrix< double >( parameters ) ) + return EXIT_FAILURE; return EXIT_SUCCESS; } - #endif /* SPARSEMATRIXBENCHMARK_H_ */ diff --git a/tests/benchmarks/tnlSpmvBenchmarkAdaptiveRgCSRMatrix.h b/tests/benchmarks/tnlSpmvBenchmarkAdaptiveRgCSRMatrix.h index 81e17151b80f09fd25c05ab213aaf5071e413eec..5d51a2723101db419640f18f9ccbc0888c26f6cb 100644 --- a/tests/benchmarks/tnlSpmvBenchmarkAdaptiveRgCSRMatrix.h +++ b/tests/benchmarks/tnlSpmvBenchmarkAdaptiveRgCSRMatrix.h @@ -19,7 +19,6 @@ #define TNLSPMVBENCHMARKADAPTIVERGCSRMATRIX_H_ #include "tnlSpmvBenchmark.h" -#include <matrices/tnlAdaptiveRgCSRMatrix.h> #include <core/tnlAssert.h> template< typename Real, typename Device, typename Index> diff --git a/tests/benchmarks/tnlSpmvBenchmarkRgCSRMatrix.h b/tests/benchmarks/tnlSpmvBenchmarkRgCSRMatrix.h index 67be6da56215329462677d92eb4d71c1cb3c1a80..5209fededda3b0bd79bf46a5b30c421318d1e5a1 100644 --- a/tests/benchmarks/tnlSpmvBenchmarkRgCSRMatrix.h +++ b/tests/benchmarks/tnlSpmvBenchmarkRgCSRMatrix.h @@ -19,7 +19,6 @@ #define TNLSPMVBENCHMARKRGCSRMATRIX_H_ #include "tnlSpmvBenchmark.h" -#include <matrices/tnlRgCSRMatrix.h> template< typename Real, typename Device, typename Index> class tnlSpmvBenchmarkRgCSRMatrix : public tnlSpmvBenchmark< Real, Device, Index, tnlRgCSRMatrix > diff --git a/tests/data/tnl-test-matrix-formats.cfg.desc b/tests/data/tnl-test-matrix-formats.cfg.desc index e166bc61dbc41f9b024d925e6ff21c6025dd2e47..bb1ceab8081c9d0bd0e5d640b8b31dccaf937b72 100644 --- a/tests/data/tnl-test-matrix-formats.cfg.desc +++ b/tests/data/tnl-test-matrix-formats.cfg.desc @@ -1,7 +1,9 @@ group Arguments { - string input-file [Input file name.]; - string matrix-format [Matrix format. It can be: dense, ellpack, sliced-ellpack, chunked-ellpack and csr.]; - bool verbose [Verbose mode.]; + string input-file [Input file name.]; + string matrix-format [Matrix format. It can be: dense, ellpack, sliced-ellpack, chunked-ellpack and csr.]; + bool hard-test( no ) [Comparison against the dense matrix.]; + bool multiplication-test( no ) [Matrix-vector multiplication test.]; + bool verbose [Verbose mode.]; },[Arguments of the matrix format test.]; diff --git a/tests/long-time-unit-tests/matrix-formats-test.h b/tests/long-time-unit-tests/matrix-formats-test.h index 0d7591f130853391b32d9abe18181c02a6372403..5bdf13ac6768d198a21d5c70631fd9d7a778ac91 100644 --- a/tests/long-time-unit-tests/matrix-formats-test.h +++ b/tests/long-time-unit-tests/matrix-formats-test.h @@ -40,6 +40,10 @@ template< typename Matrix > bool testMatrix( const tnlParameterContainer& parameters ) { Matrix matrix; + typedef typename Matrix::RealType RealType; + typedef typename Matrix::DeviceType DeviceType; + typedef typename Matrix::IndexType IndexType; + const tnlString& fileName = parameters.GetParameter< tnlString >( "input-file" ); bool verbose = parameters.GetParameter< bool >( "verbose" ); fstream file; @@ -50,16 +54,65 @@ bool testMatrix( const tnlParameterContainer& parameters ) return false; } if( ! tnlMatrixReader< Matrix >::readMtxFile( file, matrix, verbose ) ) - { - file.close(); return false; - } if( ! tnlMatrixReader< Matrix >::verifyMtxFile( file, matrix, verbose ) ) - { - file.close(); return false; + if( parameters.GetParameter< bool >( "hard-test" ) ) + { + typedef tnlDenseMatrix< RealType, DeviceType, IndexType > DenseMatrix; + DenseMatrix denseMatrix; + if( ! tnlMatrixReader< DenseMatrix >::readMtxFile( file, denseMatrix, verbose ) ) + return false; + if( ! tnlMatrixReader< DenseMatrix >::verifyMtxFile( file, denseMatrix, verbose ) ) + return false; + //matrix.print( cout ); + //denseMatrix.print( cout ); + for( IndexType i = 0; i < matrix.getRows(); i++ ) + { + for( IndexType j = 0; j < matrix.getColumns(); j++ ) + if( matrix.getElement( i, j ) != denseMatrix.getElement( i, j ) ) + { + cerr << "The matrices differ at position " << i << ", " << j << "." << endl + << " The values are " << matrix.getElement( i, j ) << " (sparse) and " + << denseMatrix.getElement( i, j ) << " (dense)." << endl; + tnlString line; + IndexType lineNumber; + if( tnlMatrixReader< Matrix >::findLineByElement( file, i, j, line, lineNumber ) ) + cerr << "The mtx file says ( line " << lineNumber << " ): " << line << endl; + else + cerr << "The element is missing in the file. Should be zero therefore." << endl; + return false; + } + if( verbose ) + cout << " Comparing the sparse matrix with the dense matrix ... " << i << " / " << matrix.getRows() << " \r" << flush; + } + if( verbose ) + cout << " Comparing the sparse matrix with the dense matrix ... OK. " << endl; + } + if( parameters.GetParameter< bool >( "multiplication-test" ) ) + { + tnlVector< RealType, DeviceType, IndexType > x, b; + x.setSize( matrix.getColumns() ); + b.setSize( matrix.getRows() ); + for( IndexType i = 0; i < x.getSize(); i++ ) + { + x.setValue( 0 ); + x.setElement( i, 1.0 ); + matrix.vectorProduct( x, b ); + for( IndexType j = 0; j < b.getSize(); j++ ) + if( b.getElement( j ) != matrix.getElement( j, i ) ) + { + cerr << "The matrix-vector multiplication gives wrong result at positions " + << j << ", " << i << ". The result is " << b.getElement( j ) << " and it should be " + << matrix.getElement( j, i ) << "." << endl; + return false; + } + if( verbose ) + cerr << " Testing the matrix-vector multiplication ... " << i << " / " << matrix.getRows() << " \r" << flush; + } + if( verbose ) + cerr << " Testing the matrix-vector multiplication ... OK. " << endl; } - file.close(); return true; } diff --git a/tests/long-time-unit-tests/tnl-run-matrix-formats-test b/tests/long-time-unit-tests/tnl-run-matrix-formats-test index 563ee83be733abd2703dc547c1bc9954341b74d5..74c74700d33298f75f04db119736a3862964ae02 100644 --- a/tests/long-time-unit-tests/tnl-run-matrix-formats-test +++ b/tests/long-time-unit-tests/tnl-run-matrix-formats-test @@ -18,7 +18,7 @@ test_matrix() # $1 - matrix name FILESIZE=$(stat -c%s "$1") echo "Testing $1 ( `echo \"$FILESIZE / 1024 / 1024 \" | bc` MB)." - $TNL_MATRIX_TEST --input-file $1 --matrix-format csr --verbose $VERBOSE + $TNL_MATRIX_TEST --input-file $1 --matrix-format chunked-ellpack --hard-test yes --verbose $VERBOSE echo "-------------------------------------------------------------------------------------" } diff --git a/tests/unit-tests/matrices/tnlTridiagonalMatrixTester.h b/tests/unit-tests/matrices/tnlTridiagonalMatrixTester.h index 13d3854fd9576ae0dac01d99031a10a7cdb5e565..ae232530df292979c4b15a8c0b2975e470dfa2e8 100644 --- a/tests/unit-tests/matrices/tnlTridiagonalMatrixTester.h +++ b/tests/unit-tests/matrices/tnlTridiagonalMatrixTester.h @@ -62,7 +62,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase void setDimensionsTest() { MatrixType m; - m.setDimensions( 10 ); + m.setDimensions( 10, 10 ); CPPUNIT_ASSERT( m.getRows() == 10 ); CPPUNIT_ASSERT( m.getColumns() == 10 ); } @@ -70,7 +70,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase void setLikeTest() { MatrixType m1, m2; - m1.setDimensions( 10 ); + m1.setDimensions( 10, 10 ); m2.setLike( m1 ); CPPUNIT_ASSERT( m1.getRows() == m2.getRows() ); } @@ -79,7 +79,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase { const int size( 10 ); MatrixType m; - m.setDimensions( size ); + m.setDimensions( size, size ); m.setValue( 1.0 ); for( int i = 0; i < size; i++ ) for( int j = 0; j < size; j++ ) @@ -92,7 +92,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase void setElementTest() { MatrixType m; - m.setDimensions( 10 ); + m.setDimensions( 10, 10 ); for( int i = 0; i < 10; i++ ) m.setElement( i, i, i ); for( int i = 0; i < 10; i++ ) @@ -107,7 +107,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase void addElementTest() { MatrixType m; - m.setDimensions( 10 ); + m.setDimensions( 10, 10 ); for( int i = 0; i < 10; i++ ) m.setElement( i, i, i ); for( int i = 0; i < 10; i++ ) @@ -129,7 +129,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase void setRowTest() { MatrixType m; - m.setDimensions( 10 ); + m.setDimensions( 10, 10 ); m.setValue( 0.0 ); for( int i = 0; i < 10; i++ ) m.setElement( i, i, i ); @@ -165,7 +165,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase v.setSize( size ); w.setSize( size ); MatrixType m; - m.setDimensions( size ); + m.setDimensions( size, size ); for( int i = 0; i < size; i++ ) { v.setElement( i, i ); @@ -181,7 +181,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase { const int size = 10; MatrixType m; - m.setDimensions( 10); + m.setDimensions( 10, 10 ); for( int i = 0; i < size; i++ ) for( int j = 0; j < size; j++ ) if( abs( i - j ) <= 1 ) @@ -213,7 +213,7 @@ class tnlTridiagonalMatrixTester : public CppUnit :: TestCase { const int size = 10; MatrixType m; - m.setDimensions( 10 ); + m.setDimensions( 10, 10 ); for( int i = 0; i < size; i++ ) for( int j = 0; j < size; j++ ) if( abs( i - j ) <= 1 )