1

일부 C++ 코드에서 희소 콜레스트 분해를 수행하기 위해 최근에 Cholmod을 설치했습니다.Cholmod와 Cholmod-Extra를 사용하여 희소 행렬의 역수를 계산합니다.

I를 계산하기 위해 원했던 : I는 (

d^T . (A^-1 + B^-1)^-1 . d 
d는 벡터 ^T 전조 나타낸다

, AB가 드문 드문 I는 다음과 같은 문제가있다)이 역행렬을 계산 부패를 사용하고자 B의 실제 역함수를 구한 다음 합계가 포함 된 솔루션을 선형으로 풀 수 있습니다. 내가 역을 계산하기 위해 의미되었다 https://cholmod-extra.readthedocs.org/en/latest/functions.html이 기능을 발견

cholmod_common Common, *cm ; 
cm = &Common ; 
cholmod_start (cm) ; 
cm->print=5; 
Common.supernodal = CHOLMOD_SUPERNODAL ; 
double *Tx, x; 
int *Ti, *Tj, *Rdeg, *Cdeg,j; 
cholmod_triplet *T ; 
size_t nrow;   
size_t ncol;   
size_t nnz ;    

nrow=Csize; 
ncol=Csize; 
nnz=numpulsars*numpulsars*numcoeff; 

T = cholmod_allocate_triplet(nrow,ncol,nnz,0,CHOLMOD_REAL, cm) ; 
Ti=(int*)T->i; 
Tj=(int*)T->j; 
Tx=(double*)T->x; 

for(int i=0;i<numpulsars;i++){  
      for(int j=0;j<numpulsars;j++){ 
        if(i==j){ 
          pcorr=1.0; 
     } 
        else{ 
          angle=pulsarCartesian[i][0]*pulsarCartesian[j][0] +pulsarCartesian[i][1]*pulsarCartesian[j][1]+pulsarCartesian[i][2]*pulsarCartesian[j][2]; 
          pcorr=(1.5*(0.5*(1-angle))*log(0.5*(1-angle)) - 0.25*0.5*(1-angle) + 0.5); 
        } 

        for(int c1=0; c1<numcoeff; c1++){ 
            val= pcorr*powercoeff[c1]; 
       Ti[T->nnz]=i*numcoeff+c1; 
       Tj[T->nnz]=j*numcoeff+c1; 
       Tx[T->nnz]=val; 
       (T->nnz)++; 


        } 
      } 
    } 


cholmod_sparse *PPFMSparse; 
cholmod_factor *L ; 
cholmod_dense *spinvK; 
PPFMSparse=cholmod_triplet_to_sparse(T,T->nnz,cm); 
// cholmod_print_sparse(PPFMSparse, "PPFMSparse", cm); 
L = cholmod_analyze (PPFMSparse, cm) ; 
cholmod_factorize (PPFMSparse, L, cm) ; 

cholmod_sparse *PPFMinv; 

PPFMinv=cholmod_spinv(L,cm); 
// cholmod_print_sparse(PPFMinv, "PPFMinv", cm); 
spinvK = cholmod_sparse_to_dense(PPFMinv, &Common) ; 
cholmod_print_dense(spinvK, "spinvK", cm); 
cholmod_free_sparse(&PPFMinv,cm); 
cholmod_free_factor (&L, cm) ; 
cholmod_free_sparse(&PPFMSparse,cm); 
cholmod_free_triplet (&T, cm) ; 
cholmod_free_dense (&spinvK, cm) ; 
cholmod_finish(cm); 

, 그러나 그것은 나에게 역의 제곱보다는 역에 관련된 무언가를 제공합니다 그것을 호출하는 코드는 다음과 같다. 난 누군가가 이것을 사용하여 어떤 경험을했는지 궁금 해서요, 또는 C + +에서 스파 스 매트릭스의 역행을 계산하는 equivilent.

건배 린

답변

1

나는 이전에 JAMA을 사용했다. 그것은 콜레 스키 분해입니다.

+0

JAMA는 스파 스 매트릭스를 지원합니까? 이 코드는 조밀 한 행렬에서 일반적인 콜레 스 decomp를 수행하는 것처럼 보입니다.이 경우에는 lapack이나 다른 것을 사용합니다. – LindleyLentati

+0

JAMA는 스파 스 매트릭스 인 http://math.nist.gov/tnt/overview.html#structures를 지원하는 TNT (Template Numerical Toolkit)를 사용합니다. – Smash