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

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

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

revision 1960 by artak, Tue Nov 4 05:50:41 2008 UTC revision 1975 by artak, Thu Nov 6 03:07:02 2008 UTC
# Line 81  Paso_Solver_AMG* Paso_Solver_getAMG(Paso Line 81  Paso_Solver_AMG* Paso_Solver_getAMG(Paso
81    dim_t i,k,j,j0;    dim_t i,k,j,j0;
82    Paso_SparseMatrix * schur=NULL;    Paso_SparseMatrix * schur=NULL;
83    Paso_SparseMatrix * schur_withFillIn=NULL;    Paso_SparseMatrix * schur_withFillIn=NULL;
84    double time0,time1,time2,S;    double time0=0,time1=0,time2=0,S;
85        
86    /* identify independend set of rows/columns */    /* identify independend set of rows/columns */
87    mis_marker=TMPMEMALLOC(n,index_t);    mis_marker=TMPMEMALLOC(n,index_t);
# Line 102  Paso_Solver_AMG* Paso_Solver_getAMG(Paso Line 102  Paso_Solver_AMG* Paso_Solver_getAMG(Paso
102    out->b_C=NULL;    out->b_C=NULL;
103    out->A=Paso_SparseMatrix_getReference(A_p);    out->A=Paso_SparseMatrix_getReference(A_p);
104    out->GS=Paso_Solver_getGS(A_p,verbose);    out->GS=Paso_Solver_getGS(A_p,verbose);
105    out->GS->sweeps=4;    out->GS->sweeps=2;
106    out->level=level;    out->level=level;
107        
108    if ( !(Paso_checkPtr(mis_marker) || Paso_checkPtr(out) || Paso_checkPtr(counter) ) ) {    if ( !(Paso_checkPtr(mis_marker) || Paso_checkPtr(out) || Paso_checkPtr(counter) ) ) {
# Line 144  Paso_Solver_AMG* Paso_Solver_getAMG(Paso Line 144  Paso_Solver_AMG* Paso_Solver_getAMG(Paso
144                for (i = 0; i < out->n_F; ++i) {                for (i = 0; i < out->n_F; ++i) {
145                  out->inv_A_FF[i]=0;                  out->inv_A_FF[i]=0;
146                  for (iPtr=A_p->pattern->ptr[out->rows_in_F[i]];iPtr<A_p->pattern->ptr[out->rows_in_F[i] + 1]; ++iPtr) {                  for (iPtr=A_p->pattern->ptr[out->rows_in_F[i]];iPtr<A_p->pattern->ptr[out->rows_in_F[i] + 1]; ++iPtr) {
147                   out->inv_A_FF[i]+=A_p->val[iPtr];                   j=A_p->pattern->index[iPtr];
148                     if (mis_marker[j])
149                         out->inv_A_FF[i]+=A_p->val[iPtr];
150                  }                  }
151                }                }
152                                
# Line 194  Paso_Solver_AMG* Paso_Solver_getAMG(Paso Line 196  Paso_Solver_AMG* Paso_Solver_getAMG(Paso
196             if( Paso_noError()) {             if( Paso_noError()) {
197                /* if there are no nodes in the coarse level there is no more work to do */                /* if there are no nodes in the coarse level there is no more work to do */
198                out->n_C=n-out->n_F;                out->n_C=n-out->n_F;
199                if (level<2) {                if (level<3) {
200                /*if (out->n_C>10) {*/                 /*if (out->n_C>10) {*/
201                     out->rows_in_C=MEMALLOC(out->n_C,index_t);                     out->rows_in_C=MEMALLOC(out->n_C,index_t);
202                     out->mask_C=MEMALLOC(n,index_t);                     out->mask_C=MEMALLOC(n,index_t);
203                     if (! (Paso_checkPtr(out->mask_C) || Paso_checkPtr(out->rows_in_C) ) ) {                     if (! (Paso_checkPtr(out->mask_C) || Paso_checkPtr(out->rows_in_C) ) ) {
# Line 296  Paso_Solver_AMG* Paso_Solver_getAMG(Paso Line 298  Paso_Solver_AMG* Paso_Solver_getAMG(Paso
298    if (Paso_noError()) {    if (Paso_noError()) {
299        if (verbose) {        if (verbose) {
300           printf("AMG: %d unknowns eliminated. %d left.\n",out->n_F,n-out->n_F);           printf("AMG: %d unknowns eliminated. %d left.\n",out->n_F,n-out->n_F);
301           if (level<2) {           if (level<3) {
302              printf("timing: AMG: MIS/reordering/elemination : %e/%e/%e\n",time2,time0,time1);              printf("timing: AMG: MIS/reordering/elemination : %e/%e/%e\n",time2,time0,time1);
303           } else {           } else {
304              printf("timing: AMG: MIS: %e\n",time2);              printf("timing: AMG: MIS: %e\n",time2);
# Line 341  void Paso_Solver_solveAMG(Paso_Solver_AM Line 343  void Paso_Solver_solveAMG(Paso_Solver_AM
343       double *bnew=MEMALLOC(amg->n,double);       double *bnew=MEMALLOC(amg->n,double);
344       double *x0=MEMALLOC(amg->n,double);       double *x0=MEMALLOC(amg->n,double);
345    
346       if (amg->level==2) {       if (amg->level==3) {
347       /*if (amg->n_C<=10) {*/       /*if (amg->n_C<=10) {*/
348          Paso_UMFPACK1(amg->A,x,b,0);          Paso_UMFPACK1(amg->A,x,b,0);
349       } else {       } else {
# Line 409  void Paso_Solver_solveAMG(Paso_Solver_AM Line 411  void Paso_Solver_solveAMG(Paso_Solver_AM
411       /*postsmoothing*/       /*postsmoothing*/
412       #pragma omp parallel for private(i) schedule(static)       #pragma omp parallel for private(i) schedule(static)
413       for (i=0;i<amg->n;++i) r[i]=b[i];       for (i=0;i<amg->n;++i) r[i]=b[i];
414        
      /*for (i=0;i<amg->n;++i) fprintf(stderr,"%lf %lf %d %d %d BBB\n",r[i],x[i],amg->n,amg->A->numCols, amg->level);*/  
       
      /*  
      for (i = 0; i < amg->A->numRows; ++i) {  
           for (iPtr=amg->A->pattern->ptr[i];iPtr<amg->A->pattern->ptr[i + 1]; ++iPtr) {  
               j=amg->A->pattern->index[iPtr];  
               fprintf(stderr,"A[%d,%d]=%lf ",i,j,amg->A->val[iPtr]);  
           }  
           fprintf(stderr,"\n");  
      }*/  
       
415       /*r=b-Ax*/       /*r=b-Ax*/
416       Paso_SparseMatrix_MatrixVector_CSR_OFFSET0(-1.,amg->A,x,1.,r);       Paso_SparseMatrix_MatrixVector_CSR_OFFSET0(-1.,amg->A,x,1.,r);
417       Paso_Solver_solveGS(amg->GS,x0,r);       Paso_Solver_solveGS(amg->GS,x0,r);
418            
419       #pragma omp parallel for private(i) schedule(static)       #pragma omp parallel for private(i) schedule(static)
420       for (i=0;i<amg->n;++i) {       for (i=0;i<amg->n;++i) {
       /*fprintf(stderr,"%lf %lf %lf\n",x[i],x0[i],x[i]+x0[i]);*/  
421        x[i]+=x0[i];        x[i]+=x0[i];
422       }       }
423            

Legend:
Removed from v.1960  
changed lines
  Added in v.1975

  ViewVC Help
Powered by ViewVC 1.1.26