gpt4 book ai didi

c++ - 如何有效地组装 FEM 稀疏矩阵

转载 作者:太空狗 更新时间:2023-10-29 21:34:16 25 4
gpt4 key购买 nike

全部处理,感谢您花时间阅读我的问题。我正在使用 Eigen3.3.4( http://eigen.tuxfamily.org/index.php?title=Main_Page ) 编写一些 FEM 代码。

我阅读了 Eigen3.3.4 文档,在这个网站(http://eigen.tuxfamily.org/dox/TopicFunctionTakingEigenTypes.html)中,它说我们应该使用 Ref<MatrixBase>避免额外的拷贝并获得高性能。

所以在我的 FEM 代码中,对于稀疏矩阵组装部分,假设函数是:

FormFE(const Ref<VectorXd> &U,const Ref<VectorXd> &V,
Ref<SparseMatrix<double> > AMATRIX,Ref<VectorXd> RHS)

其中 U 代表位移,V 代表速度项。 AMATRIX 是我的稀疏矩阵,RHS 是残差项。

然后我尝试在组装之前首先初始化我的 AMATRIX(我有一个包含所有非零元素及其值的 tripletList(我将值设置为零以进行初始化))所以我尝试了:

AMATRIX.setFromTriplets(ZeroTripList.begin(),ZeroTripList.end());

但是我有一个错误:

class Eigen::Ref<Eigen::SparseMatrix<double, 0, int> >’ has no member named ‘setFromTriplets

那么我该如何解决这个问题呢?

我的解决方案之一是使用:

FormFE(const Ref<VectorXd> &U,const Ref<VectorXd> &V,
SparseMatrix<double> &AMATRIX,Ref<VectorXd> RHS)

这工作得很好,但我不确定它是否有效。我不太擅长 cpp :P

其实我的问题是:

  1. 如何有效地使用 Eigen(尤其是 FEM 计算),我几乎在每个与 FEM 相关的函数中都使用了 Eigen 的 VectorXd 和 MatrixXd。
  2. 如何有效地组装一个稀疏矩阵?
  3. 是否可以为 FEM 组装做一些 OpenMP 并行化?
  4. 欢迎对基于 C++ 的 FEM 编码提出任何有用的建议(库推荐或任何有用的想法)!

谢谢。最好的问候。

最佳答案

是的,传递一个 SparseMatrix<double> &在这里做正确的事。 Ref<SparseMatrix>的目的|传递类似于 SparseMatrix 的组装对象,就像一个次稀疏矩阵,Map<SparseMatrix> ...

使用 setFromTriplets 也是确保获得良好性能的正确做法。使用 mat.insert(i,j) = val; 直接插入元素如果操作正确(即正确调用保留和正确的插入顺序),可能会快 2 倍。但如果你弄错了,它也会慢 100 倍......请参阅文档。与 SparseMatrix::insert也可以使用 OpenMP 填充矩阵,但这需要更加小心和严格,这是一个典型的模式:

int n_cols = ??, n_rows = ??;
std::vector<int> nnz_per_col(n_cols);
// set each nnz_per_col[j] to the exact number
// of non-zero entries in the j-th column (or more, but NOT less)
SparseMatrix<double> mat(n_rows, n_cols);
#pragma omp parallel for
for(int j=0; j<cols; ++j) {
for each non zero entry i in the j-th column {
// preferably with increasing i
double val_i_j = ...;
mat.insert(i,j) = val_i_j;
}
}

当然,如果对您来说更容易,您也可以按行工作。在这种情况下使用 SparseMatrix<double,RowMajor> .当然,您可以调整此模式以处理列/行 block 等。

如果对于装配,您需要使用一些密集矩阵/vector ,那么我相信它们非常小且大小固定。那么不要使用 MatrixXd/VectorXd,最好使用 Matrix<double,N,M>Matrix<double,N,1>静态分配的类型。这将防止大量内存分配/取消分配。

最后,最重要的建议:如果您关心性能,请不要忘记在优化代码时花时间和精力进行分析之前剖析您的代码。此外,始终在 bench/profile 中启用编译器优化。

关于c++ - 如何有效地组装 FEM 稀疏矩阵,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/47093836/

25 4 0
Copyright 2021 - 2024 cfsdn All Rights Reserved 蜀ICP备2022000587号
广告合作:1813099741@qq.com 6ren.com