FEM 희소 행렬을 효율적으로 조립하는 방법

Question

모두 거래, 제 질문을 읽어 주셔서 감사합니다. 일부 FEM 코드를 작성하기 위해 Eigen3.3.4 (http://eigen.tuxfamily.org/index.php?title=Main_Page)를 사용하고 있습니다.

Eigen3.3.4 문서를 읽었으며이 웹 사이트 (http://eigen.tuxfamily.org/dox/TopicFunctionTakingEigenTypes.html)에 이라고 말하면 추가 사본을 피하고 고성능을 얻으려면 Ref<MatrixBase>을 사용해야합니다. 일부 조립 스파 스 매트릭스 내 FEM 코드 그래서

,의 함수가 가정 해 봅시다 : U는 변위를 나타냅니다

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

는, V는 속도 기간을 나타냅니다. AMATRIX는 내 희소 행렬이고 RHS는 잔여 항입니다.

그럼 내가 먼저 내 A- 매트릭스 전에 (내가 모든 비 - 제로 요소와 값 (I 초기화 0으로 값을 설정)을 포함하는 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. 효율적으로 (특히 FEM 계산을위한) 아이겐를 사용하는 방법, 나는 각각의 거의 모든 곳에서 아이겐의 VectorXd 및 MatrixXd 사용 내 FEM 관련 기능.
2. SparseMatrix를 효율적으로 조립하는 방법은 무엇입니까?
3. FEM 어셈블을 위해 OpenMP 병렬 처리를 수행 할 수 있습니까?
4. C++ 기반 FEM 코딩에 대한 유용한 제안 (라이브러리 권장 사항 또는 유용한 아이디어)을 환영합니다!

감사합니다. 감사합니다.

Answer 1

예, SparseMatrix<double> &을 전달하는 것이 올바른 방법입니다. Ref<SparseMatrix>의 목적은 Map<SparseMatrix> ...

사용 setFromTriplets도 OK 성능을해야 할 옳은 일이며, 하위 희소 행렬처럼 SparseMatrix에 ressemble 조립 객체를 전달합니다. mat.insert(i,j) = val;을 사용하여 직접 요소를 삽입하는 것은 올바르게 수행하면 (예 : 적절한 예약 및 적절한 삽입 순서) 최대 2 배까지 빨라질 수 있습니다. 하지만 잘못 이해하면 x100 배 더 느려질 수도 있습니다 ... 문서를 참조하십시오.물론

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::insert로는 OpenMP를를 사용하여 매트릭스를 작성하지만이 훨씬 더 신중하고 엄격해야 할 수도 있습니다, 여기에 대한 일반적인 패턴은 그것은 당신을 위해 더 쉽습니다. 이 경우 SparseMatrix<double,RowMajor>을 사용하십시오. 그리고 물론이 패턴을 조정하여 열/행 블록 등에서 작업 할 수 있습니다.

조립품의 경우 일부 고밀도 매트릭스/벡터를 사용해야하는 경우 고정 크기 . 그런 다음 MatrixXd/VectorXd를 사용하는 대신 정적으로 할당 된 Matrix<double,N,M> 및 Matrix<double,N,1> 유형을 사용하십시오. 이렇게하면 많은 메모리 할당/할당 해제를 방지 할 수 있습니다.

마지막으로, 가장 중요한 추천 : 당신이 성능에 대해 걱정하는 경우, 코드를 최적화하는 경우 시간과 노력을 조사하기 전에 코드 프로필 하는 것을 잊지 마십시오. 또한 항상 컴파일러 최적화를 사용하는 벤치/프로파일.