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

Diff of /trunk/paso/src/SparseMatrix.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 112  Paso_SparseMatrix* Paso_SparseMatrix_all Line 112  Paso_SparseMatrix* Paso_SparseMatrix_all
112    bool_t unroll=FALSE;    bool_t unroll=FALSE;
113    
114    if (patternIsUnrolled) {    if (patternIsUnrolled) {
115       if (! XNOR(type & MATRIX_FORMAT_OFFSET1, pattern->type & PATTERN_FORMAT_OFFSET1) ) {       if (! XNOR(type & MATRIX_FORMAT_OFFSET1, pattern->type & MATRIX_FORMAT_OFFSET1) ) {
116           Esys_setError(TYPE_ERROR,"Paso_SparseMatrix_alloc: requested offset and pattern offset does not match.");           Esys_setError(TYPE_ERROR,"Paso_SparseMatrix_alloc: requested offset and pattern offset does not match.");
117           return NULL;           return NULL;
118       }       }
# Line 126  Paso_SparseMatrix* Paso_SparseMatrix_all Line 126  Paso_SparseMatrix* Paso_SparseMatrix_all
126          /* or if lock size one requested and the block size is not 1 */          /* or if lock size one requested and the block size is not 1 */
127      ||  ((type & MATRIX_FORMAT_BLK1) &&  (col_block_size>1) )      ||  ((type & MATRIX_FORMAT_BLK1) &&  (col_block_size>1) )
128          /* offsets don't match */          /* offsets don't match */
129      || ( (type & MATRIX_FORMAT_OFFSET1) != ( pattern->type & PATTERN_FORMAT_OFFSET1) ) ;      || ( (type & MATRIX_FORMAT_OFFSET1) != ( pattern->type & MATRIX_FORMAT_OFFSET1) ) ;
130    
131    pattern_format_out= (type & MATRIX_FORMAT_OFFSET1)? PATTERN_FORMAT_OFFSET1:  PATTERN_FORMAT_DEFAULT;    pattern_format_out= (type & MATRIX_FORMAT_OFFSET1)? MATRIX_FORMAT_OFFSET1:  MATRIX_FORMAT_DEFAULT;
132    
133    Esys_resetError();    Esys_resetError();
134    out=MEMALLOC(1,Paso_SparseMatrix);    out=MEMALLOC(1,Paso_SparseMatrix);
# Line 324  Paso_SparseMatrix* Paso_SparseMatrix_loa Line 324  Paso_SparseMatrix* Paso_SparseMatrix_loa
324      }      }
325      row_ptr[M] = nz;      row_ptr[M] = nz;
326    
327          mainPattern=Paso_Pattern_alloc(PATTERN_FORMAT_DEFAULT,M,N,row_ptr,col_ind);          mainPattern=Paso_Pattern_alloc(MATRIX_FORMAT_DEFAULT,M,N,row_ptr,col_ind);
328      out  = Paso_SparseMatrix_alloc(MATRIX_FORMAT_DEFAULT, mainPattern, 1, 1, TRUE);      out  = Paso_SparseMatrix_alloc(MATRIX_FORMAT_DEFAULT, mainPattern, 1, 1, TRUE);
329      /* copy values and cleanup temps */      /* copy values and cleanup temps */
330      for( i=0; i<nz; i++ ) out->val[i] = val[i];      for( i=0; i<nz; i++ ) out->val[i] = val[i];
# Line 416  dim_t Paso_SparseMatrix_maxDeg(Paso_Spar Line 416  dim_t Paso_SparseMatrix_maxDeg(Paso_Spar
416  {  {
417     return Paso_Pattern_maxDeg(A_p->pattern);     return Paso_Pattern_maxDeg(A_p->pattern);
418  }  }
   
 Paso_SparseMatrix* Paso_SparseMatrix_unroll(Paso_SparseMatrix* A) {  
     Paso_SparseMatrix *out = NULL;  
     /*Paso_Pattern* mainPattern=NULL;*/  
     int i;  
     dim_t block_size=A->row_block_size;  
     index_t iptr, optr1=0, optr2=0, optr3=0;  
       
     /*mainPattern= Paso_Pattern_unrollBlocks(A->pattern, PATTERN_FORMAT_DEFAULT,  A->row_block_size,A->col_block_size);*/  
     out  = Paso_SparseMatrix_alloc(MATRIX_FORMAT_BLK1, A->pattern, A->row_block_size,A->col_block_size, FALSE);  
       
     if (block_size==1) {  
        #pragma omp parallel for private(i,iptr) schedule(static)  
        for (i=0;i<A->numRows;++i) {  
          for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {  
              out->val[iptr]=A->val[iptr];  
          }  
        }  
     } else if (block_size==2) {  
        #pragma omp parallel for private(i,iptr,optr1,optr2) schedule(static)  
        for (i=0;i<A->numRows;++i) {  
          optr1=out->pattern->ptr[i*block_size];  
          optr2=out->pattern->ptr[i*block_size+1];  
          for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {  
                 out->val[optr1]=A->val[iptr*block_size*block_size];  
             optr1++;  
             out->val[optr1]=A->val[iptr*block_size*block_size+2];  
             optr1++;  
             out->val[optr2]=A->val[iptr*block_size*block_size+1];  
             optr2++;  
             out->val[optr2]=A->val[iptr*block_size*block_size+3];  
             optr2++;  
          }  
        }  
           
     } else if (block_size==3) {  
        #pragma omp parallel for private(i,iptr,optr1,optr2,optr3) schedule(static)  
        for (i=0;i<A->numRows;++i) {  
             optr1=out->pattern->ptr[i*block_size];  
             optr2=out->pattern->ptr[i*block_size+1];  
             optr3=out->pattern->ptr[i*block_size+2];  
          for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {  
             out->val[optr1]=A->val[iptr*block_size*block_size];  
             optr1++;  
             out->val[optr1]=A->val[iptr*block_size*block_size+3];  
             optr1++;  
             out->val[optr1]=A->val[iptr*block_size*block_size+6];  
             optr1++;  
             out->val[optr2]=A->val[iptr*block_size*block_size+1];  
             optr2++;  
             out->val[optr2]=A->val[iptr*block_size*block_size+4];  
             optr2++;  
             out->val[optr2]=A->val[iptr*block_size*block_size+7];  
             optr2++;  
             out->val[optr3]=A->val[iptr*block_size*block_size+2];  
             optr3++;  
             out->val[optr3]=A->val[iptr*block_size*block_size+5];  
             optr3++;  
             out->val[optr3]=A->val[iptr*block_size*block_size+8];  
             optr3++;  
          }  
        }  
           
           
     }  
   
  return out;  
 }  
419  dim_t Paso_SparseMatrix_getTotalNumRows(const Paso_SparseMatrix* A){  dim_t Paso_SparseMatrix_getTotalNumRows(const Paso_SparseMatrix* A){
420     return A->numRows * A->row_block_size;     return A->numRows * A->row_block_size;
421  }  }
# Line 496  dim_t Paso_SparseMatrix_getNumRows(Paso_ Line 428  dim_t Paso_SparseMatrix_getNumRows(Paso_
428  }  }
429  dim_t Paso_SparseMatrix_getNumCols(Paso_SparseMatrix* A) {  dim_t Paso_SparseMatrix_getNumCols(Paso_SparseMatrix* A) {
430     return A->numCols;       return A->numCols;  
 }  
431    }

Legend:
Removed from v.3311  
changed lines
  Added in v.3312

  ViewVC Help
Powered by ViewVC 1.1.26