RcppbdEigen_hdf5_Block

C++ Function Reference

1 Signature

void BigDataStatMeth::RcppbdEigen_hdf5_Block(T *dsA, BigDataStatMeth::hdf5Dataset *dsd, BigDataStatMeth::hdf5Dataset *dsu, int k, const std::string &which, int ncv, bool bcenter, bool bscale, double tol, int max_iter, bool compute_vectors, Rcpp::Nullable< int > threads=R_NilValue)

2 Description

Block-wise eigendecomposition for large HDF5 matrices using Spectra.

3 Parameters

  • dsA (T *): Input dataset (square matrix)
  • dsd (BigDataStatMeth::hdf5Dataset *): Output dataset for eigenvalues (1 x k)
  • dsu (BigDataStatMeth::hdf5Dataset *): Output dataset for eigenvectors (n x k), may be nullptr
  • k (int): Number of eigenvalues (must be < n, caller ensures this)

4 Returns

Type: class T

5 Details

dsAInput dataset (square matrix) dsdOutput dataset for eigenvalues (1 x k) dsuOutput dataset for eigenvectors (n x k), may be nullptr kNumber of eigenvalues (must be < n, caller ensures this)

6 Call Graph

Function dependencies

7 Source Code

File: inst/include/hdf5Algebra/matrixEigenDecomposition.hpp • Lines 189-265

inline void RcppbdEigen_hdf5_Block(T* dsA,
                                        BigDataStatMeth::hdf5Dataset* dsd,
                                        BigDataStatMeth::hdf5Dataset* dsu,
                                        int k, const std::string& which, int ncv,
                                        bool bcenter, bool bscale,
                                        double tol, int max_iter,
                                        bool compute_vectors,
                                        Rcpp::Nullable<int> threads = R_NilValue)
    {
        static_assert(std::is_same<T*, BigDataStatMeth::hdf5Dataset*>::value ||
                      std::is_same<T*, BigDataStatMeth::hdf5DatasetInternal*>::value,
                      "Error - type not allowed");
        try {
            std::vector<hsize_t> stride = {1, 1}, block = {1, 1};
            
            hsize_t n = dsA->nrows();
            if (n != dsA->ncols())
                throw std::runtime_error("Matrix must be square for eigendecomposition");
            
            std::tie(k, ncv) = validateSpectraParams((int)n, k, ncv);
            
            std::vector<double> vdA(n * n);
            dsA->readDatasetBlock({0, 0}, {n, n}, stride, block, vdA.data());
            Eigen::MatrixXd X = Eigen::Map<Eigen::Matrix<double, Eigen::Dynamic,
                                                          Eigen::Dynamic, Eigen::RowMajor>>(
                vdA.data(), n, n);
            if (bcenter || bscale)
                X = RcppNormalize_Data(X, bcenter, bscale, false);
            
            eigdecomp reteig = RcppbdEigen_spectra(X, k, which, ncv, false, false, tol, max_iter);
            if (!reteig.bconv)
                throw std::runtime_error("Eigendecomposition failed to converge");
            
            int comp = (int)dsA->getCompressionLevel();
            
            dsd->inheritCompressionLevel(comp);
            dsd->createDataset(1, reteig.eigenvalues_real.size(), "real");
            dsd->writeDataset(Rcpp::wrap(reteig.eigenvalues_real));
            
            if (compute_vectors && reteig.bcomputevectors && dsu != nullptr) {
                dsu->inheritCompressionLevel(comp);
                dsu->createDataset(reteig.eigenvectors_real.rows(),
                                   reteig.eigenvectors_real.cols(), "real");
                dsu->writeDataset(Rcpp::wrap(reteig.eigenvectors_real));
            }
            
            if (!reteig.is_symmetric && reteig.eigenvalues_imag.cwiseAbs().maxCoeff() > 1e-14) {
                BigDataStatMeth::hdf5Dataset* dsd_imag = new BigDataStatMeth::hdf5Dataset(
                    dsA->getFullPath(), "EIGEN/" + dsA->getDatasetName(), "values_imag", true);
                dsd_imag->inheritCompressionLevel(comp);
                dsd_imag->createDataset(1, reteig.eigenvalues_imag.size(), "real");
                dsd_imag->writeDataset(Rcpp::wrap(reteig.eigenvalues_imag));
                delete dsd_imag;
                
                if (compute_vectors && reteig.bcomputevectors && dsu != nullptr &&
                    reteig.eigenvectors_imag.cwiseAbs().maxCoeff() > 1e-14) {
                    BigDataStatMeth::hdf5Dataset* dsu_imag = new BigDataStatMeth::hdf5Dataset(
                        dsA->getFullPath(), "EIGEN/" + dsA->getDatasetName(), "vectors_imag", true);
                    dsu_imag->inheritCompressionLevel(comp);
                    dsu_imag->createDataset(reteig.eigenvectors_imag.rows(),
                                            reteig.eigenvectors_imag.cols(), "real");
                    dsu_imag->writeDataset(Rcpp::wrap(reteig.eigenvectors_imag));
                    delete dsu_imag;
                }
            }
        } catch(H5::FileIException& error) {
            throw std::runtime_error("c++ exception RcppbdEigen_hdf5_Block (File IException)");
        } catch(H5::DataSetIException& error) {
            throw std::runtime_error("c++ exception RcppbdEigen_hdf5_Block (DataSet IException)");
        } catch(H5::GroupIException& error) {
            throw std::runtime_error("c++ exception RcppbdEigen_hdf5_Block (Group IException)");
        } catch(std::exception &ex) {
            throw std::runtime_error(std::string("C++ exception RcppbdEigen_hdf5_Block: ") + ex.what());
        } catch (...) {
            throw std::runtime_error("C++ exception RcppbdEigen_hdf5_Block (unknown reason)");
        }
    }

8 Usage Example

#include "BigDataStatMeth.hpp"

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