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

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

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

revision 414 by jgs, Wed Nov 9 02:02:19 2005 UTC revision 415 by gross, Wed Jan 4 05:37:33 2006 UTC
# Line 24  double* Paso_SystemMatrix_borrowNormaliz Line 24  double* Paso_SystemMatrix_borrowNormaliz
24     index_t iptr,icol_failed,irow_failed;     index_t iptr,icol_failed,irow_failed;
25     double fac;     double fac;
26     if (!A->normalizer_is_valid) {     if (!A->normalizer_is_valid) {
27        switch(A->type) {        if ((A->type & MATRIX_FORMAT_CSC) || (A->type & MATRIX_FORMAT_SYM) || (A->type & MATRIX_FORMAT_OFFSET1)) {
28          case CSR:          Paso_setError(TYPE_ERROR,"No normalization available for compressed sparse column, symmetric storage scheme or index offset 1.");
29          } else {
30            irow_failed=-1;            irow_failed=-1;
31            #pragma omp parallel for private(ir,irb,irow,fac,iptr,icb) schedule(static)            #pragma omp parallel for private(ir,irb,irow,fac,iptr,icb) schedule(static)
32            for (ir=0;ir< A->pattern->n_ptr;ir++) {            for (ir=0;ir< A->pattern->n_ptr;ir++) {
33          for (irb=0;irb< A->row_block_size;irb++) {          for (irb=0;irb< A->row_block_size;irb++) {
34            irow=irb+A->row_block_size*ir;            irow=irb+A->row_block_size*ir;
35                fac=0.;                fac=0.;
36            for (iptr=A->pattern->ptr[ir]-PTR_OFFSET;iptr<A->pattern->ptr[ir+1]-PTR_OFFSET; iptr++) {            for (iptr=A->pattern->ptr[ir];iptr<A->pattern->ptr[ir+1]; iptr++) {
37                for (icb=0;icb< A->col_block_size;icb++)                for (icb=0;icb< A->col_block_size;icb++)
38                      fac+=ABS(A->val[iptr*A->block_size+irb+A->row_block_size*icb]);                      fac+=ABS(A->val[iptr*A->block_size+irb+A->row_block_size*icb]);
39                }                }
# Line 44  double* Paso_SystemMatrix_borrowNormaliz Line 45  double* Paso_SystemMatrix_borrowNormaliz
45                }                }
46              }              }
47            }            }
48            if (irow_failed>=0) {            if (irow_failed>=0)
49               Paso_setError(ZERO_DIVISION_ERROR,"There is a row containing zero entries only.");               Paso_setError(ZERO_DIVISION_ERROR,"There is a row containing zero entries only.");
50            }            A->normalizer_is_valid=TRUE;
51            break;        }
         case CSC:  
           icol_failed=-1;  
           #pragma omp parallel for private(ic,icb,icol,fac,iptr,irb) schedule(static)  
           for (ic=0;ic< A->pattern->n_ptr;ic++) {  
              for (icb=0;icb< A->col_block_size;icb++) {  
            icol=icb+A->col_block_size*ic;  
                fac=0.;  
            for (iptr=A->pattern->ptr[ic]-PTR_OFFSET;iptr<A->pattern->ptr[ic+1]-PTR_OFFSET; iptr++) {  
               for (irb=0;irb< A->row_block_size;irb++)  
                      fac+=ABS(A->val[iptr*A->block_size+irb+A->row_block_size*icb]);  
                }  
                if (fac>0) {  
                   A->normalizer[icol]=1./fac;  
                } else {  
                   A->normalizer[icol]=1.;  
                   icol_failed=icol;  
                }  
              }  
            }  
            if (icol_failed>=0) {  
               Paso_setError(ZERO_DIVISION_ERROR,"There is a column containing zero entries only.");  
            }  
            break;  
       default:  
         Paso_setError(TYPE_ERROR,"No normalization available for this matrix type.");  
       } /* switch A->type */  
       A->normalizer_is_valid=TRUE;  
52     }     }
53     return A->normalizer;     return A->normalizer;
54  }  }

Legend:
Removed from v.414  
changed lines
  Added in v.415

  ViewVC Help
Powered by ViewVC 1.1.26