/[escript]/trunk/finley/src/Assemble_addToSystemMatrix.c
ViewVC logotype

Annotation of /trunk/finley/src/Assemble_addToSystemMatrix.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1552 - (hide annotations)
Thu May 8 08:52:41 2008 UTC (11 years, 5 months ago) by gross
File MIME type: text/plain
File size: 9482 byte(s)
some changes to make the implementatiopn of a upwind MPI version easier
1 ksteube 1312
2     /* $Id$ */
3    
4     /*******************************************************
5     *
6     * Copyright 2003-2007 by ACceSS MNRF
7     * Copyright 2007 by University of Queensland
8     *
9     * http://esscc.uq.edu.au
10     * Primary Business: Queensland, Australia
11     * Licensed under the Open Software License version 3.0
12     * http://www.opensource.org/licenses/osl-3.0.php
13     *
14     *******************************************************/
15    
16 jgs 150 /**************************************************************/
17    
18     /* Finley: SystemMatrix and SystemVector */
19    
20     /* adds the matrix array[Equa,Sol,NN,NN] onto the matrix in. */
21     /* the rows/columns are given by */
22     /* i_Equa+Equa*Nodes_Equa[Nodes[j_Equa]] (i_Equa=0:Equa; j_Equa=0:NN_Equa). */
23     /* the routine has to be called from a parallel region */
24    
25     /* This routine assumes that in->Equa=in->Sol=1, i.e. */
26     /* array is fully packed. */
27     /* TODO: the case in->Equa!=1 */
28    
29     /**************************************************************/
30    
31     #include "Assemble.h"
32 ksteube 1312 #include "IndexList.h"
33 jgs 150
34     /**************************************************************/
35    
36     void Finley_Assemble_addToSystemMatrix(Paso_SystemMatrix* in,dim_t NN_Equa,index_t* Nodes_Equa, dim_t num_Equa,
37 ksteube 1312 dim_t NN_Sol,index_t* Nodes_Sol, dim_t num_Sol, double* array) {
38 gross 416 index_t index_offset=(in->type & MATRIX_FORMAT_OFFSET1 ? 1:0);
39 ksteube 1312 dim_t k_Equa,j_Equa,j_Sol,k_Sol,i_Equa,i_Sol,l_col,l_row,ic,ir,index,k,i_row, i_col;
40 gross 1552 index_t *mainBlock_ptr, *mainBlock_index, *col_coupleBlock_ptr, *col_coupleBlock_index, *row_coupleBlock_ptr, *row_coupleBlock_index;
41     double *mainBlock_val, *row_coupleBlock_val, *col_coupleBlock_val;
42 jgs 150 dim_t row_block_size=in->row_block_size;
43     dim_t col_block_size=in->col_block_size;
44     dim_t block_size=in->block_size;
45     dim_t num_subblocks_Equa=num_Equa/row_block_size;
46     dim_t num_subblocks_Sol=num_Sol/col_block_size;
47 ksteube 1312 dim_t numMyCols=in->pattern->mainPattern->numInput;
48     dim_t numMyRows=in->pattern->mainPattern->numOutput;
49 jgs 150
50 ksteube 1312
51 gross 416 if (in->type & MATRIX_FORMAT_CSC) {
52 ksteube 1312 /* MATRIX_FORMAT_CSC does not support MPI !!!!! */
53     mainBlock_ptr=in->mainBlock->pattern->ptr;
54     mainBlock_index=in->mainBlock->pattern->index;
55     mainBlock_val=in->mainBlock->val;
56 gross 416 for (k_Sol=0;k_Sol<NN_Sol;k_Sol++) {
57     j_Sol=Nodes_Sol[k_Sol];
58     for (l_col=0;l_col<num_subblocks_Sol;++l_col) {
59 ksteube 1312 i_col=j_Sol*num_subblocks_Sol+l_col;
60 gross 416 for (k_Equa=0;k_Equa<NN_Equa;k_Equa++) {
61     j_Equa=Nodes_Equa[k_Equa];
62     for (l_row=0;l_row<num_subblocks_Equa;++l_row) {
63 ksteube 1312 i_row=j_Equa*num_subblocks_Equa+index_offset+l_row;
64     for (k=mainBlock_ptr[i_col]-index_offset; k<mainBlock_ptr[i_col+1]-index_offset;++k) {
65     if (mainBlock_index[k] == i_row) {
66 jgs 150 for (ic=0;ic<col_block_size;++ic) {
67     i_Sol=ic+col_block_size*l_col;
68     for (ir=0;ir<row_block_size;++ir) {
69     i_Equa=ir+row_block_size*l_row;
70 ksteube 1312 mainBlock_val[k*block_size+ir+row_block_size*ic]+=
71 jgs 150 array[INDEX4(i_Equa,i_Sol,k_Equa,k_Sol,num_Equa,num_Sol,NN_Equa)];
72     }
73     }
74     break;
75     }
76     }
77     }
78     }
79 gross 416 }
80     }
81 ksteube 1312 } else if (in->type & MATRIX_FORMAT_TRILINOS_CRS) {
82     /* this needs to be modified */
83     #ifdef TRILINOS
84     for (k_Equa=0;k_Equa<NN_Equa;++k_Equa) { /* Down columns of array */
85     j_Equa=Nodes_Equa[k_Equa];
86     if (j_Equa < in->mainBlock->pattern->output_node_distribution->numLocal) {
87     for (k_Sol=0;k_Sol<NN_Sol;++k_Sol) { /* Across rows of array */
88     j_Sol=Nodes_Sol[k_Sol];
89     for (l_row=0;l_row<num_subblocks_Equa;++l_row) {
90     irow=j_Equa*row_block_size+l_row;
91     for (l_col=0;l_col<col_block_size;++l_col) {
92     icol=j_Sol*col_block_size+index_offset+l_col;
93 ksteube 1322 /* irow is local and icol is global */
94 ksteube 1312 Trilinos_SumIntoMyValues(in->trilinos_data, irow, icol, array[INDEX4(l_row,l_col,k_Equa,k_Sol,num_Equa,num_Sol,NN_Equa)]);
95     }
96     }
97     }
98     }
99     }
100     #endif
101 gross 416 } else {
102 ksteube 1312
103     mainBlock_ptr=in->mainBlock->pattern->ptr;
104     mainBlock_index=in->mainBlock->pattern->index;
105     mainBlock_val=in->mainBlock->val;
106 gross 1552 col_coupleBlock_ptr=in->col_coupleBlock->pattern->ptr;
107     col_coupleBlock_index=in->col_coupleBlock->pattern->index;
108     col_coupleBlock_val=in->col_coupleBlock->val;
109     row_coupleBlock_ptr=in->row_coupleBlock->pattern->ptr;
110     row_coupleBlock_index=in->row_coupleBlock->pattern->index;
111     row_coupleBlock_val=in->row_coupleBlock->val;
112 ksteube 1312
113     for (k_Equa=0;k_Equa<NN_Equa;++k_Equa) { /* Down columns of array */
114     j_Equa=Nodes_Equa[k_Equa];
115     for (l_row=0;l_row<num_subblocks_Equa;++l_row) {
116     i_row=j_Equa*num_subblocks_Equa+l_row;
117     /* only look at the matrix rows stored on this processor */
118     if (i_row < numMyRows) {
119     for (k_Sol=0;k_Sol<NN_Sol;++k_Sol) { /* Across rows of array */
120     j_Sol=Nodes_Sol[k_Sol];
121     for (l_col=0;l_col<num_subblocks_Sol;++l_col) {
122     /* only look at the matrix rows stored on this processor */
123     i_col=j_Sol*num_subblocks_Sol+index_offset+l_col;
124     if (i_col < numMyCols + index_offset ) {
125     for (k=mainBlock_ptr[i_row]-index_offset;k<mainBlock_ptr[i_row+1]-index_offset;++k) {
126     if (mainBlock_index[k]==i_col) {
127     /* Entry array(k_Sol, j_Equa) is a block (row_block_size x col_block_size) */
128     for (ic=0;ic<col_block_size;++ic) {
129     i_Sol=ic+col_block_size*l_col;
130     for (ir=0;ir<row_block_size;++ir) {
131     i_Equa=ir+row_block_size*l_row;
132     mainBlock_val[k*block_size+ir+row_block_size*ic]+=
133     array[INDEX4(i_Equa,i_Sol,k_Equa,k_Sol,num_Equa,num_Sol,NN_Equa)];
134     }
135     }
136     break;
137     }
138     }
139     } else {
140 gross 1552 for (k=col_coupleBlock_ptr[i_row]-index_offset;k<col_coupleBlock_ptr[i_row+1]-index_offset;++k) {
141     if (col_coupleBlock_index[k] == i_col-numMyCols) {
142 ksteube 1312 /* Entry array(k_Sol, j_Equa) is a block (row_block_size x col_block_size) */
143     for (ic=0;ic<col_block_size;++ic) {
144     i_Sol=ic+col_block_size*l_col;
145     for (ir=0;ir<row_block_size;++ir) {
146     i_Equa=ir+row_block_size*l_row;
147 gross 1552 col_coupleBlock_val[k*block_size+ir+row_block_size*ic]+=
148 ksteube 1312 array[INDEX4(i_Equa,i_Sol,k_Equa,k_Sol,num_Equa,num_Sol,NN_Equa)];
149     }
150     }
151     break;
152     }
153     }
154     }
155 jgs 150 }
156 ksteube 1312 }
157 gross 1552 } else {
158     for (k_Sol=0;k_Sol<NN_Sol;++k_Sol) { /* Across rows of array */
159     j_Sol=Nodes_Sol[k_Sol];
160     for (l_col=0;l_col<num_subblocks_Sol;++l_col) {
161     i_col=j_Sol*num_subblocks_Sol+index_offset+l_col;
162     if (i_col < numMyCols + index_offset ) {
163     for (k=row_coupleBlock_ptr[i_row-numMyRows]-index_offset;k<row_coupleBlock_ptr[i_row-numMyRows+1]-index_offset;++k) {
164     if (row_coupleBlock_index[k] == i_col) {
165     /* Entry array(k_Sol, j_Equa) is a block (row_block_size x col_block_size) */
166     for (ic=0;ic<col_block_size;++ic) {
167     i_Sol=ic+col_block_size*l_col;
168     for (ir=0;ir<row_block_size;++ir) {
169     i_Equa=ir+row_block_size*l_row;
170     row_coupleBlock_val[k*block_size+ir+row_block_size*ic]+=
171     array[INDEX4(i_Equa,i_Sol,k_Equa,k_Sol,num_Equa,num_Sol,NN_Equa)];
172     }
173     }
174     break;
175     }
176     }
177     }
178     }
179     }
180     }
181 ksteube 1312 }
182     }
183 jgs 150 }
184     }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26