Fix QRDecomposition

Bug #1209830 reported by Richard Gomes
6
This bug affects 1 person
Affects Status Importance Assigned to Milestone
JQuantLib
New
Medium
Unassigned

Bug Description

QRDecomposition needs debugging.

QuantLib calls a QR decomposition implemented by MINPACK.
MINPACK was initially implemented in FORTRAN and then translated to Java.

Our code is based on MINPACK/J but adapted to look like other decompositions, which I adapted from JAMA.

Class MatrixTest has a test case for QRDecomposition.

Probably a good strategy to solve this issue would be calling MINPACK/J, mimicking what QuantLib does and compare against our test case. It would free us of having to run QuantLib inside CDT.

=============
Relationships
=============
related to http://bugs.launchpad.net/bugs/jquantlib-232

Tags: code-review
Revision history for this message
Richard Gomes (frgomes) wrote :

QRDecomposition is important to perform inversion of matrices.

I've adapted QRDecomposition from Apache Commons Math v2.0-RC4

MatrixTest has a number of test cases involving matrix decompositions.
QRDecomposition pass but QRsolve is still failing.

Revision history for this message
Richard Gomes (frgomes) wrote :

back to the pot

Revision history for this message
Richard Gomes (frgomes) wrote :

Kicked to v0.1.4

Let me know if there are dependencies in the current v0.1.3 which requires we solve this issue before the v0.1.3 release.

Thanks

Richard Gomes

Revision history for this message
Richard Gomes (frgomes) wrote :

In the next few days I will be commiting my changes in classes Cells, Array and Matrix. I've removed the support for FORTRAN indexing which was poluting the code and making it more complex, more difficult to understand and maintain.

Now it is possible to associate an AddressPolicy to classes Array and Matrix. In practice it means that we can determine a 'projection' over the underlying memory buffer which holds the data. For example: we can obtain a sub-Matrix from an existing Matrix just creating a new Matrix object and passing a AddressPolicy to it which tells how cell <0,0> can be obtained given the existing underlying Matrix. So, imagine that we have an existing 5x5 matrix and we are interested only on the last 2 rows of it, namely, rows 4 and 5. All we need to do it tell that row 0 on the new matrix is equivalent to row 4 on the existing matrix and row 1 on the new matrix is equivalent on row 1 on the existing matrix.

Back to the problem of QR decomposition.

MINPACK/J was always the solution I found the cleanest one. The only problem was it was using FORTRAN addressing, which means that row and column indexes start from 1 (like in FORTRAN) instead of zero (as C/C++/Java do). It was obliging us to create a new matrix with a void row zero and a void column 0 just to allow cell <1,1> match our previously existing cell <0,0>. It's certainly an ugly solution and imposes the creation of an intermediary matrix, consuming memory and precious CPU cycles.

With AddressPolicy, this problem is solved: all we need to do is create a new Matrix object (backed by the already existing data buffer) which tells that indexes need to be subtracted by 1 before being used to return data from the underlying data buffer.

This solution means that no temporary data storage is needed. The only price we will have to pay would be 2 additional arithmetic operations when calculating the address of a certain cell. Instead of calculating like this:

    return row * cols + col;

we will calculate as:

    return (row-1) * cols + (col -1);

This can be easily done by creating a specific AddressPolicy which implements the concept of FORTRAN addressing.

To post a comment you must log in.
This report contains Public information  
Everyone can see this information.

Other bug subscribers

Remote bug watches

Bug watches keep track of this bug in other bug trackers.