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

Annotation of /trunk/paso/src/SystemMatrixPattern_getSubpattern.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 415 - (hide annotations)
Wed Jan 4 05:37:33 2006 UTC (14 years, 6 months ago) by gross
File MIME type: text/plain
File size: 2953 byte(s)
a better way of representing the matrix format type is introduced. this is needed for the Paradiso and UMFPACK interface
1 jgs 150 /* $Id$ */
2    
3     /**************************************************************/
4    
5     /* Paso: SystemMatrixPatternPattern */
6    
7     /**************************************************************/
8    
9     /* Copyrights by ACcESS Australia 2003, 2004, 2005 */
10     /* Author: gross@access.edu.au */
11    
12     /**************************************************************/
13    
14     #include "Paso.h"
15     #include "Util.h"
16     #include "SystemMatrixPattern.h"
17    
18     /**************************************************************/
19    
20     /* creates SystemMatrixPattern */
21    
22     Paso_SystemMatrixPattern* Paso_SystemMatrixPattern_getSubpattern(Paso_SystemMatrixPattern* pattern, \
23     int new_n_rows, index_t* row_list,index_t* new_col_index) {
24 gross 415 index_t index_offset=(pattern->type & PATTERN_FORMAT_OFFSET1 ? 1:0);
25 jgs 150 Paso_SystemMatrixPattern*out=NULL;
26     index_t *ptr=NULL,*index=NULL,k,j,subpattern_row,tmp;
27     dim_t i;
28     Paso_resetError();
29    
30     ptr=MEMALLOC(new_n_rows+1,index_t);
31     if (! Paso_checkPtr(ptr)) {
32     #pragma omp parallel
33     {
34     #pragma omp for private(i) schedule(static)
35     for (i=0;i<new_n_rows+1;++i) ptr[i]=0;
36    
37     /* find the number column entries in each row */
38     #pragma omp for private(i,k,j,subpattern_row) schedule(static)
39     for (i=0;i<new_n_rows;++i) {
40     j=0;
41     subpattern_row=row_list[i];
42 gross 415 for (k=pattern->ptr[subpattern_row]-index_offset;k<pattern->ptr[subpattern_row+1]-index_offset;++k)
43     if (new_col_index[pattern->index[k]-index_offset]>-1) j++;
44 jgs 150 ptr[i]=j;
45     }
46     }
47     /* accummulate ptr */
48     ptr[new_n_rows]=Paso_Util_cumsum(new_n_rows,ptr);
49     index=MEMALLOC(ptr[new_n_rows],index_t);
50     if (Paso_checkPtr(index)) {
51     MEMFREE(ptr);
52     } else {
53     /* find the number column entries in each row */
54     #pragma omp parallel for private(i,k,j,subpattern_row,tmp) schedule(static)
55     for (i=0;i<new_n_rows;++i) {
56     j=ptr[i];
57     subpattern_row=row_list[i];
58 gross 415 for (k=pattern->ptr[subpattern_row]-index_offset;k<pattern->ptr[subpattern_row+1]-index_offset;++k) {
59     tmp=new_col_index[pattern->index[k]-index_offset];
60 jgs 150 if (tmp>-1) {
61     index[j]=tmp;
62     ++j;
63     }
64     }
65     }
66     /* create return value */
67 gross 415 out=Paso_SystemMatrixPattern_alloc(pattern->type,new_n_rows,ptr,index);
68 jgs 150 if (! Paso_noError()) {
69     MEMFREE(index);
70     MEMFREE(ptr);
71     }
72     }
73     }
74     return out;
75     }
76     /*
77     * $Log$
78     * Revision 1.2 2005/09/15 03:44:39 jgs
79     * Merge of development branch dev-02 back to main trunk on 2005-09-15
80     *
81     * Revision 1.1.2.1 2005/09/05 06:29:47 gross
82     * These files have been extracted from finley to define a stand alone libray for iterative
83     * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
84     * has not been tested yet.
85     *
86     *
87     */

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26