RcppQRHdf5

C++ Function Reference

1 Signature

void BigDataStatMeth::RcppQRHdf5(BigDataStatMeth::hdf5Dataset *dsA, BigDataStatMeth::hdf5Dataset *dsQ, BigDataStatMeth::hdf5Dataset *dsR, bool bthin, Rcpp::Nullable< int > block_size, Rcpp::Nullable< int > threads)

2 Description

QR decomposition for HDF5 matrices.

3 Parameters

  • dsA (BigDataStatMeth::hdf5Dataset *): Input matrix dataset
  • dsQ (BigDataStatMeth::hdf5Dataset *): Output Q matrix dataset
  • dsR (BigDataStatMeth::hdf5Dataset *): Output R matrix dataset
  • bthin (bool): Whether to compute thin QR
  • block_size (Rcpp::Nullable< int >): Block size for processing (optional)
  • threads (Rcpp::Nullable< int >): Number of threads for parallel processing (optional)

4 Details

Computes the QR decomposition of a matrix stored in HDF5 format. Supports both full and thin decomposition with parallel processing.

5 Call Graph

Function dependencies

6 Source Code

File: inst/include/hdf5Algebra/matrixQR.hpp • Lines 156-236

inline void RcppQRHdf5( BigDataStatMeth::hdf5Dataset* dsA, 
                        BigDataStatMeth::hdf5Dataset* dsQ, 
                        BigDataStatMeth::hdf5Dataset* dsR, 
                        bool bthin,
                        Rcpp::Nullable<int> block_size = R_NilValue,
                        Rcpp::Nullable<int> threads = R_NilValue )
{
    
    
    try {

        int irank; //,
            // iblockfactor = 1;
        std::vector<hsize_t> offset = {0,0},
            count = {dsA->nrows(), dsA->ncols()},
            stride = {1,1},
            block = {1,1};
        Eigen::MatrixXd Q;
        
        
        std::vector<double> vdA( count[0] * count[1] ); 
        dsA->readDatasetBlock( {offset[0], offset[1]}, {count[0], count[1]}, stride, block, vdA.data() );
        Eigen::MatrixXd A = Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>> (vdA.data(), count[0], count[1] );
        A.transposeInPlace();
        
        
        Eigen::FullPivLU<Eigen::MatrixXd>lu_decomp(A);
        Eigen::HouseholderQR<Eigen::MatrixXd> qr(A);
        
        qr.compute(A);
        irank = lu_decomp.rank();
        
        
        //..//if (irank == count[0] + 1 || irank == count[1] + 1 )
        if (static_cast<unsigned long long>(irank) == count[0] + 1ULL ||
            static_cast<unsigned long long>(irank) == count[1] + 1ULL)
        {
            Eigen::MatrixXd R = qr.matrixQR().template triangularView<Eigen::Upper>();
            dsR->createDataset( R.rows(), R.cols(), "real" );
            dsR->writeDataset(Rcpp::wrap(R));
        } else {
            Eigen::MatrixXd R = qr.matrixQR().topLeftCorner(irank, irank).template triangularView<Eigen::Upper>();
            dsR->createDataset( R.rows(), R.cols(), "real" );
            dsR->writeDataset(Rcpp::wrap(R));
        }
        
        
        if (bthin == false)
        {
            Eigen::MatrixXd Q = qr.householderQ();
            dsQ->createDataset( Q.rows(), Q.cols(), "real" );
            dsQ->writeDataset( Rcpp::wrap(Q) );
        } else {
            
            //. 2025/01/15 error with thin calculus.// Eigen::MatrixXd Qthin = qr.householderQ() * Eigen::MatrixXd::Identity(count[0], count[1]);
            //. 2025/01/15 .// dsQ->createDataset( Qthin.rows(), Qthin.cols(), "real" );
            //. 2025/01/15 .// dsQ->writeDataset( Rcpp::wrap(Qthin));
            
            Eigen::MatrixXd Qthin = qr.householderQ() * Eigen::MatrixXd::Identity(count[1], count[0]);
            dsQ->createDataset( Qthin.rows(), Qthin.cols(), "real" );
            dsQ->writeDataset( Rcpp::wrap(Qthin));
            
        }
        
    } catch( H5::FileIException& error ) { // catch failure caused by the H5File operations
        Rcpp::Rcout<<"c++ exception RcppQRHdf5 (File IException)";
        return void();
    } catch( H5::GroupIException & error ) { // catch failure caused by the DataSet operations
        Rcpp::Rcout << "c++ exception RcppQRHdf5 (Group IException)";
        return void();
    } catch( H5::DataSetIException& error ) { // catch failure caused by the DataSet operations
        Rcpp::Rcout << "c++ exception RcppQRHdf5 (DataSet IException)";
        return void();
    } catch(std::exception& ex) {
        Rcpp::Rcout << "c++ exception RcppQRHdf5" << ex.what();
        return void();
    }
    
    return void();
    
}

7 Usage Example

#include "BigDataStatMeth.hpp"

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