private void solveUsingTriangle(double real, int index, DMatrixRMaj r ) { for( int i = 0; i < index; i++ ) { implicit.A.add(i,i,-real); } SpecializedOps_DDRM.subvector(implicit.A,0,index,index,false,0,r); CommonOps_DDRM.changeSign(r); TriangularSolver_DDRM.solveU(implicit.A.data,r.data,implicit.A.numRows,0,index); for( int i = 0; i < index; i++ ) { implicit.A.add(i,i,real); } }
/** * Creates a random symmetric positive definite matrix. * * @param width The width of the square matrix it returns. * @param rand Random number generator used to make the matrix. * @return The random symmetric positive definite matrix. */ public static DMatrixRMaj symmetricPosDef(int width, Random rand) { // This is not formally proven to work. It just seems to work. DMatrixRMaj a = new DMatrixRMaj(width,1); DMatrixRMaj b = new DMatrixRMaj(width,width); for( int i = 0; i < width; i++ ) { a.set(i,0,rand.nextDouble()); } CommonOps_DDRM.multTransB(a,a,b); for( int i = 0; i < width; i++ ) { b.add(i,i,1); } return b; }
private void solveWithLU(double real, int index, DMatrixRMaj r ) { DMatrixRMaj A = new DMatrixRMaj(index,index); CommonOps_DDRM.extract(implicit.A,0,index,0,index,A,0,0); for( int i = 0; i < index; i++ ) { A.add(i,i,-real); } r.reshape(index,1, false); SpecializedOps_DDRM.subvector(implicit.A,0,index,index,false,0,r); CommonOps_DDRM.changeSign(r); // TODO this must be very inefficient if( !solver.setA(A)) throw new RuntimeException("Solve failed"); solver.solve(r,r); }