As we planned, to keep the focus more on a working framework initially, I started with a direct logarithm of dense matrices that Eigen3 unsupported module provides (available for Eigen3.1.0 and later) instead of the rational approximation. This week, I have implemented the framework for dense matrix linear operators using Gaussian trace samples () and that shows that for a small matrix, say,
we can approximate the by estimating up to a certain precision with a large number of estimates. For example, , and for log-det estimates with Gaussian samples, we get an expected estimate , .
The framework has been developed in a modular way, with all independent base classes and implementations first and then the dependent ones. I added a unit-test for (almost) every component that has been added to ensure that they behave as we might expect them to do. Currently the framework has –
- CLinearOperator<T> base class, CDenseMatrixOperator<T> implementation of this, a unit-test.
- CJobResult base class, CScalarResult<T> and CVectorResult<T> implementation of this, a unit-test for scalar.
- CJobResultAggregator base class, CStoreScalarAggregator<T> implementation of this and a unit-test.
- CIndependentJob base class, CDenseExactLogJob implementation of this and a unit-test.
- CIndependentComputationEngine base class, CSerialComputationEngine implementation, and a unit-test
- COperatorFunction<T> base class, CDenseMatrixExactLog implementation of this, a unit-test
- CTraceSampler base, CNormalSampler implementation, a unit-test
- CLogDetEstimator class and a unit-test
I have a small program ready (very much similar to the CLogDetEstimator unit-test) with this which I’ll add to the libshogun examples. This shows how the framework works –
- The sample method of CLogDetEstimator computes a number of log-det estimates of a linear operator function
- In sample, it first computes the stuff that are required for the rest of the computation job (call init on COperatorFunction and CTraceSampler). For example, CDenseMatrixExactLog implementation computes the log of the dense matrix and sets that as the linear operator in the operator function in its init. CLogRationalApproximation will compute complex shifts, weights etc using Jacobi elliptic functions in init, etc. CNormalSampler initializes the number of samples it should generate per log-det estimate (which is 1 in this case). CProbingSampler implementation will use graph coloring of the sparse linear operator and set the number of colors as the number of estimates, all in its init.
- Then for the required number of log-det estimates, sample uses the submit_job method of COperatorFunction to create a number of jobs per sample, and keeps the JobResultAggregators with itself. submit_job internally creates a number of jobs (based on the implementation) with one aggregator, submits the jobs to the computation engine, which may or may not start computing those job immediately (based on implementation) and passes the aggregator to sample.
- In serial implementation of the computation engine, it starts computing the jobs immediately as soon as they are submitted (call compute method on the job) and blocks until the computation is done. For CDenseExactLogJob, the computation task is to compute first (apply linear operator on vector, result is a vector), and then (vector-vector dot product). The vector is then safely discarded.
- sample then waits for all computation jobs to be completed with the engine’s wait_for_all method. Serial implementation returns immediately in this case since all jobs are already computed.
- sample calls finalize on all the job result aggregators which shapes the final aggregation of the job results into CScalarResults. It then computes an average of the estimates and returns.
Plan for next week
To implement the CLogRationalApproximation class for dense matrix that computes weights, shifts in init (requires the eigenvalues, currently thought of using Eigen3 for that), and then CLogRationalApproximationIndividual, which creates num_shifts jobs for each trace sample and moves the shifts into the operator. A CLogRationalApproximationIndividualJob will use Eigen3’s complex solver for direct solving. I planned to keep CLinearSolver as a base class for all the solvers, which I’ll try to implement next week, along with its DirectLinearSolver implementation. CStoreVectorAggregator has to be implemented as well.
Ah, its been too long for a report! I just realized it! :)
That’s all folks! See you next week!