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

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

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

revision 1311 by ksteube, Wed Feb 14 04:40:49 2007 UTC revision 1312 by ksteube, Mon Sep 24 06:18:44 2007 UTC
# Line 1  Line 1 
1    
2  /* $Id$ */  /* $Id$ */
3    
4  /*  /*******************************************************
5  ********************************************************************************   *
6  *               Copyright   2006 by ACcESS MNRF                                *   *           Copyright 2003-2007 by ACceSS MNRF
7  *                                                                              *   *       Copyright 2007 by University of Queensland
8  *                 http://www.access.edu.au                                     *   *
9  *           Primary Business: Queensland, Australia                            *   *                http://esscc.uq.edu.au
10  *     Licensed under the Open Software License version 3.0             *   *        Primary Business: Queensland, Australia
11  *        http://www.opensource.org/licenses/osl-3.0.php                        *   *  Licensed under the Open Software License version 3.0
12  ********************************************************************************   *     http://www.opensource.org/licenses/osl-3.0.php
13  */   *
14     *******************************************************/
15    
16  /**************************************************************/  /**************************************************************/
17    
# Line 35  Line 37 
37    
38  void Paso_SystemMatrix_nullifyRowsAndCols(Paso_SystemMatrix* A, double* mask_row, double* mask_col, double main_diagonal_value) {  void Paso_SystemMatrix_nullifyRowsAndCols(Paso_SystemMatrix* A, double* mask_row, double* mask_col, double main_diagonal_value) {
39    
40    index_t index_offset=(A->type & MATRIX_FORMAT_OFFSET1 ? 1:0);    double *remote_values=NULL;
41    index_t ir,icol,iptr,icb,irb,irow,ic,l;    Paso_MPIInfo *mpi_info=A->mpi_info;
42      if (mpi_info->size>1) {
43    if (A ->col_block_size==1 && A ->row_block_size ==1) {       if (A ->col_block_size==1 && A ->row_block_size ==1) {
44      if (A->type & MATRIX_FORMAT_CSC) {         if (A->type & MATRIX_FORMAT_CSC) {
45        #pragma omp parallel for private(irow, iptr,icol) schedule(static)             Paso_setError(SYSTEM_ERROR,"Paso_SystemMatrix_nullifyRowsAndCols: CSC is not supported by MPI.");
46        for (icol=0;icol< A->pattern->n_ptr;icol++) {             return;
47      for (iptr=A->pattern->ptr[icol]-index_offset;iptr<A->pattern->ptr[icol+1]-index_offset; iptr++) {         } else if (A->type & MATRIX_FORMAT_TRILINOS_CRS) {
48        irow=A->pattern->index[iptr]-index_offset;             Paso_setError(SYSTEM_ERROR,"Paso_SystemMatrix_nullifyRowsAndCols: TRILINOS is not supported with MPI.");
49        if (mask_col[icol]>0. || mask_row[irow]>0. ) {             return;
50              if (irow==icol && mask_col[icol]>0. && mask_row[irow]>0.) {         } else {
51            A->val[iptr]=main_diagonal_value;           Paso_SystemMatrix_allocBuffer(A);
52              } else {           if (Paso_noError()) {
53            A->val[iptr]=0;              Paso_SystemMatrix_startCollect(A,mask_col) ;
54              }              Paso_SparseMatrix_nullifyRowsAndCols_CSR_BLK1(A->mainBlock,mask_row,mask_col,main_diagonal_value);
55        }              remote_values=Paso_SystemMatrix_finishCollect(A);
56      }              Paso_SparseMatrix_nullifyRowsAndCols_CSR_BLK1(A->coupleBlock,mask_row,remote_values,0.);
57        }           }
58      } else {           Paso_SystemMatrix_freeBuffer(A);
59        #pragma omp parallel for private(irow, iptr,icol) schedule(static)         }
60        for (irow=0;irow< A->pattern->n_ptr;irow++) {       } else {
61      /* TODO: test mask_row here amd not inside every loop */         if (A->type & MATRIX_FORMAT_CSC) {
62      for (iptr=A->pattern->ptr[irow]-index_offset;iptr<A->pattern->ptr[irow+1]-index_offset; iptr++) {             Paso_setError(SYSTEM_ERROR,"Paso_SystemMatrix_nullifyRowsAndCols: CSC is not supported by MPI.");
63        icol=A->pattern->index[iptr]-index_offset;             return;
64        if (mask_col[icol]>0. || mask_row[irow]>0. ) {         } else if (A->type & MATRIX_FORMAT_TRILINOS_CRS) {
65              if (irow==icol && mask_col[icol]>0. && mask_row[irow]>0.) {             Paso_setError(SYSTEM_ERROR,"Paso_SystemMatrix_nullifyRowsAndCols: TRILINOS is not supported with MPI.");
66            A->val[iptr]=main_diagonal_value;             return;
67              } else {         } else {
68            A->val[iptr]=0;           Paso_SystemMatrix_allocBuffer(A);
69              }           if (Paso_noError()) {
70        }              Paso_SystemMatrix_startCollect(A,mask_col) ;
71      }              Paso_SparseMatrix_nullifyRowsAndCols_CSR(A->mainBlock,mask_row,mask_col,main_diagonal_value);
72        }              remote_values=Paso_SystemMatrix_finishCollect(A);
73      }              Paso_SparseMatrix_nullifyRowsAndCols_CSR(A->coupleBlock,mask_row,remote_values,0.);
74    } else {           }
75      if (A->type & MATRIX_FORMAT_CSC) {           Paso_SystemMatrix_freeBuffer(A);
76        #pragma omp parallel for private(l,irow, iptr,icol,ic,irb,icb) schedule(static)         }
77        for (ic=0;ic< A->pattern->n_ptr;ic++) {       }
78      for (iptr=A->pattern->ptr[ic]-index_offset;iptr<A->pattern->ptr[ic+1]-index_offset; iptr++) {    } else {
79        for (irb=0;irb< A->row_block_size;irb++) {       if (A ->col_block_size==1 && A ->row_block_size ==1) {
80          irow=irb+A->row_block_size*(A->pattern->index[iptr]-index_offset);         if (A->type & MATRIX_FORMAT_CSC) {
81          for (icb=0;icb< A->col_block_size;icb++) {             Paso_SparseMatrix_nullifyRowsAndCols_CSC_BLK1(A->mainBlock,mask_row,mask_col,main_diagonal_value);
82            icol=icb+A->col_block_size*ic;         } else if (A->type & MATRIX_FORMAT_TRILINOS_CRS) {
83            if (mask_col[icol]>0. || mask_row[irow]>0. ) {             Paso_setError(SYSTEM_ERROR,"Paso_SystemMatrix_nullifyRowsAndCols: TRILINOS is not supported.");
84                  l=iptr*A->block_size+irb+A->row_block_size*icb;         } else {
85          if (irow==icol && mask_col[icol]>0. && mask_row[irow]>0.) {             Paso_SparseMatrix_nullifyRowsAndCols_CSR_BLK1(A->mainBlock,mask_row,mask_col,main_diagonal_value);
86            A->val[l]=main_diagonal_value;         }
87          } else {       } else {
88            A->val[l]=0;         if (A->type & MATRIX_FORMAT_CSC) {
89          }             Paso_SparseMatrix_nullifyRowsAndCols_CSC(A->mainBlock,mask_row,mask_col,main_diagonal_value);
90            }         } else if (A->type & MATRIX_FORMAT_TRILINOS_CRS) {
91          }             Paso_setError(SYSTEM_ERROR,"Paso_SystemMatrix_nullifyRowsAndCols: TRILINOS is not supported with MPI.");
92        }         } else {
93      }             Paso_SparseMatrix_nullifyRowsAndCols_CSR(A->mainBlock,mask_row,mask_col,main_diagonal_value);
94        }         }
95      } else {       }
       #pragma omp parallel for private(l,irow, iptr,icol,ir,irb,icb) schedule(static)  
       for (ir=0;ir< A->pattern->n_ptr;ir++) {  
     for (iptr=A->pattern->ptr[ir]-index_offset;iptr<A->pattern->ptr[ir+1]-index_offset; iptr++) {  
       for (irb=0;irb< A->row_block_size;irb++) {  
         irow=irb+A->row_block_size*ir;  
         for (icb=0;icb< A->col_block_size;icb++) {  
           icol=icb+A->col_block_size*(A->pattern->index[iptr]-index_offset);  
           if (mask_col[icol]>0. || mask_row[irow]>0. ) {  
                 l=iptr*A->block_size+irb+A->row_block_size*icb;  
         if (irow==icol && mask_col[icol]>0. && mask_row[irow]>0.) {  
           A->val[l]=main_diagonal_value;  
         } else {  
           A->val[l]=0;  
         }  
           }  
         }  
       }  
     }  
       }  
    }  
96    }    }
97    return;    return;
98  }  }
   
 /*  
  * $Log$  
  * Revision 1.2  2005/09/15 03:44:39  jgs  
  * Merge of development branch dev-02 back to main trunk on 2005-09-15  
  *  
  * Revision 1.1.2.1  2005/09/05 06:29:48  gross  
  * These files have been extracted from finley to define a stand alone libray for iterative  
  * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but  
  * has not been tested yet.  
  *  
  *  
  */  

Legend:
Removed from v.1311  
changed lines
  Added in v.1312

  ViewVC Help
Powered by ViewVC 1.1.26