public DoubleMatrix2D solve ( DoubleMatrix2D B ) { if ( B . rows ( ) != dim ) { log . error ( "wrong dimension of vector b: expected " + dim + ", actual " + B . rows ( ) ) ; throw new RuntimeException ( "wrong dimension of vector b: expected " + dim + ", actual " + B . rows ( ) ) ; } // with scaling, we must solve U.Q.U.z = U.b, after that we have x = U.z if ( this . rescaler != null ) { // B = ALG.mult(this.U, B); B = ColtUtils . diagonalMatrixMult ( this . U , B ) ; } int nOfColumns = B . columns ( ) ; // copy final double [ ] [ ] Y = B . copy ( ) . toArray ( ) ; // Solve LY = B (same as L.Yc = Bc for every column of Y and B) for ( int j = 0 ; j < dim ; j ++ ) { final double [ ] LTJ = LcolumnsValues [ j ] ; for ( int col = 0 ; col < nOfColumns ; col ++ ) { Y [ j ] [ col ] /= LTJ [ 0 ] ; // the diagonal of the matrix L final double YJCol = Y [ j ] [ col ] ; if ( Double . compare ( YJCol , 0. ) != 0 ) { for ( int i = j + 1 ; i < dim ; i ++ ) { Y [ i ] [ col ] -= YJCol * LTJ [ i - j ] ; } } } } // Solve L[T].X = Y (same as L[T].Xc = Yc for every column of X and Y) final DoubleMatrix2D X = F2 . make ( dim , nOfColumns ) ; for ( int i = dim - 1 ; i > - 1 ; i -- ) { final double [ ] LTI = LcolumnsValues [ i ] ; double [ ] sum = new double [ nOfColumns ] ; for ( int col = 0 ; col < nOfColumns ; col ++ ) { for ( int j = dim - 1 ; j > i ; j -- ) { sum [ col ] += LTI [ j - i ] * X . getQuick ( j , col ) ; } X . setQuick ( i , col , ( Y [ i ] [ col ] - sum [ col ] ) / LTI [ 0 ] ) ; } } if ( this . rescaler != null ) { // return ALG.mult(this.U, X); return ColtUtils . diagonalMatrixMult ( this . U , X ) ; } else { return X ; } }
Solves Q . X = B