chem/KPP/kpp/kpp-2.1/util/sutil.c

   1
   2 int KppDecomp( double *JVS )
   3 {
   4 double W[KPP_NVAR];
   5 double a;
   6 int k, kk, j, jj;
   7
   8   for( k = 0; k < KPP_NVAR; k++ ) {
   9     if( JVS[ LU_DIAG[k] ] == 0.0 ) return k+1;
  10     for( kk = LU_CROW[k]; kk < LU_CROW[k+1]; kk++ )
  11       W[ LU_ICOL[kk] ] = JVS[kk];
  12     for( kk = LU_CROW[k]; kk < LU_DIAG[k]; kk++ ) {
  13       j = LU_ICOL[kk];
  14       a = -W[j] / JVS[ LU_DIAG[j] ];
  15       W[j] = -a;
  16       for( jj = LU_DIAG[j]+1; jj < LU_CROW[j+1]; jj++ )
  17         W[ LU_ICOL[jj] ] += a*JVS[jj];
  18     }
  19     for( kk = LU_CROW[k]; kk < LU_CROW[k+1]; kk++ )
  20       JVS[kk] = W[ LU_ICOL[kk] ];
  21   }
  22   return 0;
  23 }