/[escript]/trunk/paso/src/SparseMatrix_getSubmatrix.cpp
ViewVC logotype

Diff of /trunk/paso/src/SparseMatrix_getSubmatrix.cpp

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

temp_trunk_copy/paso/src/SparseMatrix_getSubmatrix.c revision 1384 by phornby, Fri Jan 11 02:29:38 2008 UTC trunk/paso/src/SparseMatrix_getSubmatrix.c revision 2881 by jfenwick, Thu Jan 28 02:03:15 2010 UTC
# Line 1  Line 1 
1    
 /* $Id: SparseMatrix_getSubmatrix.c 1306 2007-09-18 05:51:09Z ksteube $ */  
   
2  /*******************************************************  /*******************************************************
3   *  *
4   *           Copyright 2003-2007 by ACceSS MNRF  * Copyright (c) 2003-2010 by University of Queensland
5   *       Copyright 2007 by University of Queensland  * Earth Systems Science Computational Center (ESSCC)
6   *  * http://www.uq.edu.au/esscc
7   *                http://esscc.uq.edu.au  *
8   *        Primary Business: Queensland, Australia  * Primary Business: Queensland, Australia
9   *  Licensed under the Open Software License version 3.0  * Licensed under the Open Software License version 3.0
10   *     http://www.opensource.org/licenses/osl-3.0.php  * http://www.opensource.org/licenses/osl-3.0.php
11   *  *
12   *******************************************************/  *******************************************************/
13    
14    
15  /**************************************************************/  /**************************************************************/
16    
# Line 20  Line 19 
19  /**************************************************************/  /**************************************************************/
20    
21  /* Copyrights by ACcESS Australia 2003, 2004,2005 */  /* Copyrights by ACcESS Australia 2003, 2004,2005 */
22  /* Author: gross@access.edu.au */  /* Author: Lutz Gross, l.gross@uq.edu.au */
23    
24  /**************************************************************/  /**************************************************************/
25    
# Line 51  Paso_SparseMatrix* Paso_SparseMatrix_get Line 50  Paso_SparseMatrix* Paso_SparseMatrix_get
50           sub_pattern=Paso_Pattern_getSubpattern(A->pattern,n_row_sub,n_col_sub,row_list,new_col_index);           sub_pattern=Paso_Pattern_getSubpattern(A->pattern,n_row_sub,n_col_sub,row_list,new_col_index);
51           if (Paso_noError()) {           if (Paso_noError()) {
52              /* create the return object */              /* create the return object */
53              out=Paso_SparseMatrix_alloc(type,sub_pattern,A->row_block_size,A->col_block_size);              out=Paso_SparseMatrix_alloc(type,sub_pattern,A->row_block_size,A->col_block_size,TRUE);
54              if (Paso_noError()) {              if (Paso_noError()) {
55                   #pragma omp parallel for private(i,k,m,subpattern_row,tmp) schedule(static)                   #pragma omp parallel for private(i,k,m,subpattern_row,tmp) schedule(static)
56                   for (i=0;i<n_row_sub;++i) {                   for (i=0;i<n_row_sub;++i) {
# Line 62  Paso_SparseMatrix* Paso_SparseMatrix_get Line 61  Paso_SparseMatrix* Paso_SparseMatrix_get
61                             #pragma ivdep                             #pragma ivdep
62                             for (m=out->pattern->ptr[i]-index_offset;m<out->pattern->ptr[i+1]-index_offset;++m) {                             for (m=out->pattern->ptr[i]-index_offset;m<out->pattern->ptr[i+1]-index_offset;++m) {
63                                 if (out->pattern->index[m]==tmp+index_offset) {                                 if (out->pattern->index[m]==tmp+index_offset) {
64                                     Paso_copyDouble(A->block_size,&(A->val[k*A->block_size]),&(out->val[m*A->block_size]));                                     Paso_copyShortDouble(A->block_size,&(A->val[k*A->block_size]),&(out->val[m*A->block_size]));
65                                     break;                                     break;
66                                 }                                 }
67                             }                             }
# Line 75  Paso_SparseMatrix* Paso_SparseMatrix_get Line 74  Paso_SparseMatrix* Paso_SparseMatrix_get
74        }        }
75        return out;        return out;
76  }  }
77    
78    Paso_SparseMatrix* Paso_SparseMatrix_getBlock(Paso_SparseMatrix* A, int blockid){
79          
80          dim_t blocksize=A->row_block_size,i;
81          dim_t n=A->numRows;
82          index_t iptr;
83          
84          Paso_SparseMatrix* out=NULL;
85        
86          out=Paso_SparseMatrix_alloc(A->type,A->pattern, 1, 1, 0);
87          
88          if (blocksize==1) {
89                if (blockid==1) {
90                      #pragma omp parallel for private(i,iptr) schedule(static)
91                      for(i=0;i<n;++i) {
92                          for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {
93                                  out->val[iptr]=A->val[iptr];
94                          }
95                      }
96                }
97                else {
98                      Paso_setError(VALUE_ERROR,"Requested and actual block sizes do not match.");
99                }
100          }
101          else if (blocksize==2) {
102                if (blockid==1) {
103                      #pragma omp parallel for private(i,iptr) schedule(static)
104                      for(i=0;i<n;i++) {
105                            for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {
106                                  out->val[iptr]=A->val[4*iptr];
107                            }
108                      }
109                }
110                else if (blockid==2) {
111                      #pragma omp parallel for private(i,iptr) schedule(static)
112                      for(i=0;i<n;i++) {
113                            for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {
114                                  out->val[iptr]=A->val[4*iptr+3];
115                            }
116                      }
117                }
118                else {
119                      Paso_setError(VALUE_ERROR,"Requested and actual block sizes do not match.");
120                }
121          }
122          else if (blocksize==3) {
123                if (blockid==1) {
124                      #pragma omp parallel for private(i,iptr) schedule(static)
125                      for(i=0;i<n;i++) {
126                            for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {
127                                        out->val[iptr]=A->val[9*iptr];
128                            }
129                      }
130                }
131                else if (blockid==2) {
132                      #pragma omp parallel for private(i,iptr) schedule(static)
133                      for(i=0;i<n;i++) {
134                            for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {
135                                        out->val[iptr]=A->val[9*iptr+4];
136                            }
137                      }
138                }
139                else if (blockid==3) {
140                      #pragma omp parallel for private(i,iptr) schedule(static)
141                      for(i=0;i<n;i++) {
142                            for (iptr=A->pattern->ptr[i];iptr<A->pattern->ptr[i+1]; ++iptr) {
143                                        out->val[iptr]=A->val[9*iptr+8];
144                            }
145                      }
146                }
147                else {
148                      Paso_setError(VALUE_ERROR,"Requested and actual block sizes do not match.");
149                }
150          }
151                
152             return  out;
153    }

Legend:
Removed from v.1384  
changed lines
  Added in v.2881

  ViewVC Help
Powered by ViewVC 1.1.26