RcppSolveHdf5

C++ Function Reference

1 Signature

void BigDataStatMeth::RcppSolveHdf5(BigDataStatMeth::hdf5Dataset *dsA, BigDataStatMeth::hdf5Dataset *dsB, BigDataStatMeth::hdf5Dataset *dsX)

2 Description

Solves the linear system AX = B using HDF5-stored matrices.

3 Parameters

• dsA (BigDataStatMeth::hdf5Dataset *): Input HDF5 dataset containing matrix A
• dsB (BigDataStatMeth::hdf5Dataset *): Input HDF5 dataset containing matrix B
• dsX (BigDataStatMeth::hdf5Dataset *): Output HDF5 dataset for solution matrix X

4 Details

dsAInput HDF5 dataset containing matrix A dsBInput HDF5 dataset containing matrix B dsXOutput HDF5 dataset for solution matrix X This function provides an HDF5-based implementation for solving linear systems, suitable for large matrices that don’t fit in memory. It reads data in blocks and automatically handles symmetric cases.

5 Call Graph

6 Source Code

NoteImplementation

File: inst/include/hdf5Algebra/matrixEquationSolver.hpp • Lines 154-229

inline void RcppSolveHdf5(BigDataStatMeth::hdf5Dataset* dsA, BigDataStatMeth::hdf5Dataset* dsB, BigDataStatMeth::hdf5Dataset* dsX )
    {
        
        try {
            
            std::vector<hsize_t> offset = { 0, 0 },
                                 countA = { dsA->nrows(), dsA->ncols() },
                                 countB = { dsB->nrows(), dsB->ncols() },
                                 stride = { 1, 1},
                                 block = { 1, 1};
            
            std::vector<double> vdB( countB[0] * countB[1] );
            dsB->readDatasetBlock( {offset[0], offset[1]}, {countB[0], countB[1]}, stride, block, vdB.data() );
            Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>> b (vdB.data(), countB[0], countB[1]);
            
            {
                char Uchar = 'U';
                int info = 0;
                
                // Declare matrix variables
                int n = dsA->nrows();
                int nrhs = dsB->nrows();
                int lwork = std::max( 1, n );
                int lda = std::max( 1, n );
                int ldb = std::max( 1, n );
                std::vector<int> ipiv(n);
                std::vector<double> work(lwork);
                
                std::vector<double> vdA( countA[0] * countA[1] );
                dsA->readDatasetBlock( {offset[0], offset[1]}, {countA[0], countA[1]}, stride, block, vdA.data() );
                Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>> a (vdA.data(), countA[0], countA[1]);
                
                // Solve matrix equation
                if( a == a.transpose()  )
                {
                    // dsysv_( char* UPLO, int* N , int* NRHS, double* A, int* LDA, int* IPIV, double* B, int* LDB, double* WORK, int* LWORK, int* INFO);
                    dsysv_( &Uchar, &n, &nrhs, a.data(), &lda, ipiv.data(), b.data(), &ldb, work.data(), &lwork, &info);
                } else {
                    // dgesv( int N, int NRHS, double A, int LDA, int IPIV, double B, int LDB, int INFO);
                    dgesv_( &n, &nrhs, a.data(), &lda, ipiv.data(), b.data(), &ldb, &info );
                }
                
            }
            
            dsX->writeDataset(b.data());
            
        } catch( H5::FileIException& error ) { // catch failure caused by the H5File operations
            checkClose_file(dsA, dsB, dsX);
            dsA = dsB = dsX = nullptr;
            Rcpp::Rcerr<<"\nc++ exception RcppSolveHdf5 (File IException)";
            return void();
        } catch( H5::GroupIException & error ) { // catch failure caused by the DataSet operations
            checkClose_file(dsA, dsB, dsX);
            dsA = dsB = dsX = nullptr;
            Rcpp::Rcerr<<"\nc++ exception RcppSolveHdf5 (Group IException)";
            return void();
        } catch( H5::DataSetIException& error ) { // catch failure caused by the DataSet operations
            checkClose_file(dsA, dsB, dsX);
            dsA = dsB = dsX = nullptr;
            Rcpp::Rcerr<<"\nc++ exception RcppSolveHdf5 (DataSet IException)";
            return void();
        } catch(std::exception& ex) {
            checkClose_file(dsA, dsB, dsX);
            dsA = dsB = dsX = nullptr;
            Rcpp::Rcerr<<"\nc++ exception RcppSolveHdf5" << ex.what();
            return void();
        } catch (...) {
            checkClose_file(dsA, dsB, dsX);
            dsA = dsB = dsX = nullptr;
            Rcpp::Rcerr<<"\nC++ exception RcppSolveHdf5 (unknown reason)";
            return void();
        }
        
        return void();
        
    }

7 Usage Example

#include "BigDataStatMeth.hpp"

// Example usage
auto result = RcppSolveHdf5(...);