我已经问了一些类似的问题,但这次我会说得更具体。
我需要在 1)Apache数学库 2)并行Colt库 3)JLapack库 在上述三种情况中的任何一种情况下,如果与MATLAB相比,所需的时间都非常长。 因此,我想知道在Java是否有任何用于Cholesky分解的高度优化的外部工具:例如,我一直在考虑CHOLMOD算法,它实际上在 我真的很感激你能就这件事给我一个全面的反馈。
下面是一些用于Java的BLAS库的一个很好的总结:Java-Matrix-Math-Libraries的Performance-of-Java-Matrix-Math-Libraries。您还可以在java-matrix-benchmark中看到许多这些库的基准测试。
然而,在我的经验中,这些库中的大多数似乎都不适合解决大型稀疏矩阵。在我的例子中,我所做的是通过JNI使用特征来实现求解。
Eigen对它的线性求解器有很好的讨论,包括Cholmod上的一个。
对于我的8860x8860稀疏矩阵,通过JNI使用Eigen's求解器比并行colt快20倍,比我自己的密集求解器快10倍。更重要的是,它看起来是以
实际上有一个带有Java的Eigen包装器,叫做JEigen,它使用JNI。然而,它没有实现稀疏矩阵的求解,所以它不能包装所有的东西。
我最初使用JNA,但对开销并不满意。关于如何使用JNI,维基百科有一个很好的例子。编写函数声明并用
例如用于
//Cholesky.java
package cfd.optimisation;
//ri, ci, v : matrix row indices, column indices, and values
//y = Ax where A is a nxn matrix with nnz non-zero values
public class Cholesky {
private static native void solve_eigenLDLTx(int[] ri, int[] ci, double[] v, double[] x, double[] y, int n, int nnz);
}
使用
JNIEXPORT void JNICALL Java_cfd_optimisation_Cholesky_solve_1eigenLDLTx
(JNIEnv *, jclass, jintArray, jintArray, jdoubleArray, jdoubleArray, jdoubleArray, jint, jint);
下面是我如何实现求解器的
JNIEXPORT void JNICALL Java_cfd_optimisation_Cholesky_solve_1eigenLDLTx(JNIEnv *env, jclass obj, jintArray arrri, jintArray arrci, jdoubleArray arrv, jdoubleArray arrx, jdoubleArray arry, jint jn, jint jnnz) {
int n = jn;
int *ri = (int*)env->GetPrimitiveArrayCritical(arrri, 0);
int *ci = (int*)env->GetPrimitiveArrayCritical(arrci, 0);
double *v = (double*)env->GetPrimitiveArrayCritical(arrv, 0);
int nnz = jnnz;
double *x = (double*)env->GetPrimitiveArrayCritical(arrx, 0);
double *y = (double*)env->GetPrimitiveArrayCritical(arry, 0);
Eigen::SparseMatrix<double> A = colt2eigen(ri, ci, v, nnz, n);
//Eigen::MappedSparseMatrix<double> A(n, n, nnz, ri, ci, v);
Eigen::VectorXd a(n), b(n);
for (int i = 0; i < n; i++) a(i) = x[i];
//a = Eigen::Map<Eigen::VectorXd>(x, n).cast<double>();
Eigen::SimplicialCholesky<Eigen::SparseMatrix<double> > solver;
solver.setMode(Eigen::SimplicialCholeskyLDLT);
b = solver.compute(A).solve(a);
for (int i = 0; i < n; i++) y[i] = b(i);
env->ReleasePrimitiveArrayCritical(arrri, ri, 0);
env->ReleasePrimitiveArrayCritical(arrci, ci, 0);
env->ReleasePrimitiveArrayCritical(arrv, v, 0);
env->ReleasePrimitiveArrayCritical(arrx, x, 0);
env->ReleasePrimitiveArrayCritical(arry, y, 0);
}
函数
Eigen::SparseMatrix<double> colt2eigen(int *ri, int *ci, double* v, int nnz, int n) {
std::vector<Eigen::Triplet<double>> tripletList;
for (int i = 0; i < nnz; i++) {
tripletList.push_back(Eigen::Triplet<double>(ri[i], ci[i], v[i]));
}
Eigen::SparseMatrix<double> m(n, n);
m.setFromTriplets(tripletList.begin(), tripletList.end());
return m;
}
其中一个棘手的部分是从Java和柯尔特那里得到这些数组。为了做这件事我做了这件事
//y = A x: x and y are double[] arrays and A is DoubleMatrix2D
int nnz = A.cardinality();
DoubleArrayList v = new DoubleArrayList(nnz);
IntArrayList ci = new IntArrayList(nnz);
IntArrayList ri = new IntArrayList(nnz);
A.forEachNonZero((row, column, value) -> {
v.add(value); ci.add(column); ri.add(row); return value;}
);
Cholesky.solve_eigenLDLTx(ri.elements(), ci.elements(), v.elements(), x, y, n, nnz);
我没有使用过这些工具,但我怀疑您正在受到这样一个事实的打击,即Java在某些版本中/在某些平台上没有使用原生处理器浮点平方根指令。
请参阅:在哪里可以找到Java's平方根函数的源代码?
如果绝对精确度不是一个要求,您可以尝试切换上述实现之一,使用平方根的近似值(参见Java中以牺牲精确度为代价的快速sqrt),这应该会更快一些。