/[escript]/trunk/paso/src/AMG.c
ViewVC logotype

Diff of /trunk/paso/src/AMG.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 3311 by gross, Mon Oct 25 04:33:31 2010 UTC revision 3312 by gross, Tue Oct 26 07:54:58 2010 UTC
# Line 109  Paso_Preconditioner_LocalAMG* Paso_Preco Line 109  Paso_Preconditioner_LocalAMG* Paso_Preco
109          */          */
110          n_F=Paso_Util_cumsum_maskedTrue(n,counter, split_marker);          n_F=Paso_Util_cumsum_maskedTrue(n,counter, split_marker);
111          n_C=n-n_F;          n_C=n-n_F;
112          if (verbose) printf("Paso AMG level %d: %d unknowns are flagged for elimination. %d left.\n",level,n_F,n-n_F);          if (verbose) printf("Paso: AMG level %d: %d unknowns are flagged for elimination. %d left.\n",level,n_F,n-n_F);
113            
114          if ( n_F == 0 ) {  /*  is a nasty case. a direct solver should be used, return NULL */          if ( n_F == 0 ) {  /*  is a nasty case. a direct solver should be used, return NULL */
115             out = NULL;             out = NULL;
# Line 165  Paso_Preconditioner_LocalAMG* Paso_Preco Line 165  Paso_Preconditioner_LocalAMG* Paso_Preco
165              if (SHOW_TIMING) printf("timing: level %d: Paso_SparseMatrix_getTranspose: %e\n",level,Esys_timer()-time0);              if (SHOW_TIMING) printf("timing: level %d: Paso_SparseMatrix_getTranspose: %e\n",level,Esys_timer()-time0);
166                                            
167              /*              /*
168              construct coarse level matrix (can we do this in one call?)              construct coarse level matrix:
169              */              */
170              time0=Esys_timer();              time0=Esys_timer();
171              Atemp=Paso_SparseMatrix_MatrixMatrix(A_p,out->P);              Atemp=Paso_SparseMatrix_MatrixMatrix(A_p,out->P);
# Line 195  Paso_Preconditioner_LocalAMG* Paso_Preco Line 195  Paso_Preconditioner_LocalAMG* Paso_Preco
195                    out->refinements = options->coarse_matrix_refinements;                    out->refinements = options->coarse_matrix_refinements;
196                    /* no coarse level matrix has been constructed. use direct solver */                    /* no coarse level matrix has been constructed. use direct solver */
197                    #ifdef MKL                    #ifdef MKL
198                      Atemp=Paso_SparseMatrix_unroll(A_C);                      out->A_C=Paso_SparseMatrix_unroll(MATRIX_FORMAT_BLK1 + MATRIX_FORMAT_OFFSET1, A_C);
199                      Paso_SparseMatrix_free(A_C);                      Paso_SparseMatrix_free(A_C);
                     out->A_C=Paso_SparseMatrix_alloc(MATRIX_FORMAT_BLK1 + MATRIX_FORMAT_OFFSET1, Atemp->pattern,1,1, FALSE);  
                     #pragma omp parallel for private(i) schedule(static)  
                     for (i=0;i<out->A->len;++i) {  
                        out->A_C->val[i]=Atemp->val[i];  
                     }  
                     Paso_SparseMatrix_free(Atemp);  
200                      out->A_C->solver_package = PASO_MKL;                      out->A_C->solver_package = PASO_MKL;
201                    #else                    #else
202                      #ifdef UMFPACK                      #ifdef UMFPACK
203                         out->A_C=Paso_SparseMatrix_unroll(A_C);                         out->A_C=Paso_SparseMatrix_unroll(MATRIX_FORMAT_BLK1 + MATRIX_FORMAT_CSC, A_C);
204                         Paso_SparseMatrix_free(A_C);                         Paso_SparseMatrix_free(A_C);
205                         out->A_C->solver_package = PASO_UMFPACK;                         out->A_C->solver_package = PASO_UMFPACK;
206                      #else                      #else
# Line 233  Paso_Preconditioner_LocalAMG* Paso_Preco Line 227  Paso_Preconditioner_LocalAMG* Paso_Preco
227    TMPMEMFREE(S);    TMPMEMFREE(S);
228    
229    if (Esys_noError()) {    if (Esys_noError()) {
       if (verbose) printf("AMG: level %d: %d unknowns eliminated.\n",level, n_F);  
230       return out;       return out;
231    } else  {    } else  {
232       Paso_Preconditioner_LocalAMG_free(out);       Paso_Preconditioner_LocalAMG_free(out);
# Line 258  void Paso_Preconditioner_LocalAMG_solve( Line 251  void Paso_Preconditioner_LocalAMG_solve(
251       if (amg->n_F < amg->n) { /* is there work on the coarse level? */       if (amg->n_F < amg->n) { /* is there work on the coarse level? */
252           time0=Esys_timer();           time0=Esys_timer();
253       Paso_Copy(n, amg->r, b);                            /*  r <- b */       Paso_Copy(n, amg->r, b);                            /*  r <- b */
254       Paso_SparseMatrix_MatrixVector_CSR_OFFSET0_DIAG(1.,amg->R,amg->r,0.,amg->b_C); /*r=r-Ax*/       Paso_SparseMatrix_MatrixVector_CSR_OFFSET0(-1.,A,x,1.,amg->r); /*r=r-Ax*/
255           Paso_SparseMatrix_MatrixVector_CSR_OFFSET0(1.,amg->R,amg->r,0.,amg->b_C);  /* b_c = R*r  */       Paso_SparseMatrix_MatrixVector_CSR_OFFSET0_DIAG(1.,amg->R,amg->r,0.,amg->b_C);  /* b_c = R*r  */
256           time0=Esys_timer()-time0;           time0=Esys_timer()-time0;
257       /* coarse level solve */       /* coarse level solve */
258       if ( amg->AMG_C == NULL) {       if ( amg->AMG_C == NULL) {
Legend:
Removed from v.3311  
changed lines
  Added in v.3312
  ViewVC Help
Powered by ViewVC 1.1.26