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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2413 - (hide annotations)
Tue May 12 01:26:20 2009 UTC (10 years, 4 months ago) by artak
File MIME type: text/plain
File size: 8702 byte(s)
There was a bug in computing logical_block_size. It was always returning 0.
1 ksteube 1312
2     /*******************************************************
3 ksteube 1811 *
4     * Copyright (c) 2003-2008 by University of Queensland
5     * Earth Systems Science Computational Center (ESSCC)
6     * http://www.uq.edu.au/esscc
7     *
8     * Primary Business: Queensland, Australia
9     * Licensed under the Open Software License version 3.0
10     * http://www.opensource.org/licenses/osl-3.0.php
11     *
12     *******************************************************/
13 dhawcroft 631
14 ksteube 1811
15 jgs 150 /**************************************************************/
16    
17     /* Paso: SystemMatrix */
18    
19     /**************************************************************/
20    
21     /* Author: gross@access.edu.au */
22    
23     /**************************************************************/
24    
25     #include "SystemMatrix.h"
26    
27     /**************************************************************/
28    
29     /* allocates a SystemMatrix of type type using the given matrix pattern
30     if type is UNKOWN CSR is used.
31     if CSC or CSC_BLK1 is used pattern has to give the CSC pattern.
32     if CSR or CSR_BLK1 is used pattern has to give the CSR pattern.
33     Values are initialized by zero. */
34    
35     Paso_SystemMatrix* Paso_SystemMatrix_alloc(Paso_SystemMatrixType type,Paso_SystemMatrixPattern *pattern, int row_block_size, int col_block_size) {
36 gross 432
37 jgs 150 double time0;
38     Paso_SystemMatrix*out=NULL;
39 gross 1028 dim_t n_norm,i;
40 ksteube 1312 Paso_SystemMatrixType pattern_format_out;
41 gross 1028
42 jgs 150 Paso_resetError();
43     time0=Paso_timer();
44     out=MEMALLOC(1,Paso_SystemMatrix);
45     if (! Paso_checkPtr(out)) {
46 ksteube 1312 out->type=type;
47 jgs 150 out->pattern=NULL;
48 ksteube 1312 out->row_distribution=NULL;
49     out->col_distribution=NULL;
50     out->mpi_info=Paso_MPIInfo_getReference(pattern->mpi_info);
51 gross 1552 out->row_coupler=NULL;
52     out->col_coupler=NULL;
53 ksteube 1312 out->mainBlock=NULL;
54 gross 1552 out->row_coupleBlock=NULL;
55     out->col_coupleBlock=NULL;
56 ksteube 1312 out->normalizer_is_valid=FALSE;
57     out->normalizer=NULL;
58 gross 425 out->solver_package=PASO_PASO;
59     out->solver=NULL;
60 ksteube 1312 out->trilinos_data=NULL;
61 jgs 150 out->reference_counter=1;
62 ksteube 1312
63 gross 1552
64 ksteube 1312 pattern_format_out= (type & MATRIX_FORMAT_OFFSET1)? PATTERN_FORMAT_OFFSET1: PATTERN_FORMAT_DEFAULT;
65     /* ====== compressed sparse columns === */
66 gross 415 if (type & MATRIX_FORMAT_CSC) {
67     if (type & MATRIX_FORMAT_SYM) {
68     Paso_setError(TYPE_ERROR,"Generation of matrix pattern for symmetric CSC is not implemented yet.");
69     } else {
70 gross 425 if ((type & MATRIX_FORMAT_BLK1) || row_block_size!=col_block_size || col_block_size>3) {
71 ksteube 1312 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,pattern_format_out,col_block_size,row_block_size);
72 jgs 150 out->row_block_size=1;
73     out->col_block_size=1;
74 gross 415 } else {
75 ksteube 1312 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,pattern_format_out,1,1);
76 jgs 150 out->row_block_size=row_block_size;
77     out->col_block_size=col_block_size;
78 gross 415 }
79     }
80 ksteube 1312 out->row_distribution=Paso_Distribution_getReference(out->pattern->input_distribution);
81     out->col_distribution=Paso_Distribution_getReference(out->pattern->output_distribution);
82 gross 415 } else {
83 ksteube 1312 /* ====== compressed sparse row === */
84 gross 415 if (type & MATRIX_FORMAT_SYM) {
85     Paso_setError(TYPE_ERROR,"Generation of matrix pattern for symmetric CSR is not implemented yet.");
86     } else {
87 gross 425 if ((type & MATRIX_FORMAT_BLK1) || row_block_size!=col_block_size || col_block_size>3) {
88 ksteube 1312 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,pattern_format_out,row_block_size,col_block_size);
89 gross 415 out->row_block_size=1;
90     out->col_block_size=1;
91     } else {
92 ksteube 1312 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,pattern_format_out,1,1);
93 gross 415 out->row_block_size=row_block_size;
94     out->col_block_size=col_block_size;
95     }
96     }
97 ksteube 1312 out->row_distribution=Paso_Distribution_getReference(out->pattern->output_distribution);
98     out->col_distribution=Paso_Distribution_getReference(out->pattern->input_distribution);
99 jgs 150 }
100     out->logical_row_block_size=row_block_size;
101     out->logical_col_block_size=col_block_size;
102 artak 2413 out->logical_block_size=out->logical_row_block_size*out->logical_col_block_size;
103 jgs 150 out->block_size=out->row_block_size*out->col_block_size;
104 gross 1552 out->col_coupler=Paso_Coupler_alloc(pattern->col_connector,out->col_block_size);
105     out->row_coupler=Paso_Coupler_alloc(pattern->row_connector,out->row_block_size);
106 ksteube 1312 /* this should be bypassed if trilinos is used */
107     if (type & MATRIX_FORMAT_TRILINOS_CRS) {
108     #ifdef TRILINOS
109     out->trilinos_data=Paso_TRILINOS_alloc();
110     #endif
111     } else {
112     out->solver_package=PASO_PASO;
113     out->mainBlock=Paso_SparseMatrix_alloc(type,out->pattern->mainPattern,row_block_size,col_block_size);
114 gross 1552 out->col_coupleBlock=Paso_SparseMatrix_alloc(type,out->pattern->col_couplePattern,row_block_size,col_block_size);
115     out->row_coupleBlock=Paso_SparseMatrix_alloc(type,out->pattern->row_couplePattern,row_block_size,col_block_size);
116 ksteube 1312 /* allocate memory for matrix entries */
117     if (type & MATRIX_FORMAT_CSC) {
118     n_norm = out->mainBlock->numCols * out->col_block_size;
119     } else {
120     n_norm = out->mainBlock->numRows * out->row_block_size;
121     }
122     out->normalizer=MEMALLOC(n_norm,double);
123     out->normalizer_is_valid=FALSE;
124     if (! Paso_checkPtr(out->normalizer)) {
125     #pragma omp parallel for private(i) schedule(static)
126     for (i=0;i<n_norm;++i) out->normalizer[i]=0.;
127     }
128     }
129 jgs 150 }
130     /* all done: */
131     if (! Paso_noError()) {
132 ksteube 1312 Paso_SystemMatrix_free(out);
133 jgs 150 return NULL;
134     } else {
135     #ifdef Paso_TRACE
136     printf("timing: system matrix %.4e sec\n",Paso_timer()-time0);
137 ksteube 1312 printf("Paso_SystemMatrix_alloc: system matrix has been allocated.\n");
138 jgs 150 #endif
139     return out;
140     }
141     }
142    
143     /* returns a reference to Paso_SystemMatrix in */
144    
145     Paso_SystemMatrix* Paso_SystemMatrix_reference(Paso_SystemMatrix* in) {
146     if (in!=NULL) ++(in->reference_counter);
147 gross 1792 return in;
148 jgs 150 }
149    
150     /* deallocates a SystemMatrix: */
151    
152 ksteube 1312 void Paso_SystemMatrix_free(Paso_SystemMatrix* in) {
153 jgs 150 if (in!=NULL) {
154     in->reference_counter--;
155     if (in->reference_counter<=0) {
156 ksteube 1312 Paso_SystemMatrixPattern_free(in->pattern);
157     Paso_Distribution_free(in->row_distribution);
158     Paso_Distribution_free(in->col_distribution);
159     Paso_MPIInfo_free(in->mpi_info);
160 gross 1552 Paso_Coupler_free(in->row_coupler);
161     Paso_Coupler_free(in->col_coupler);
162 ksteube 1312 Paso_SparseMatrix_free(in->mainBlock);
163 gross 1552 Paso_SparseMatrix_free(in->col_coupleBlock);
164     Paso_SparseMatrix_free(in->row_coupleBlock);
165 jgs 150 MEMFREE(in->normalizer);
166     Paso_solve_free(in);
167 ksteube 1312 #ifdef TRILINOS
168     Paso_TRILINOS_free(in->trilinos_data);
169     #endif
170 jgs 150 MEMFREE(in);
171     #ifdef Paso_TRACE
172 ksteube 1312 printf("Paso_SystemMatrix_free: system matrix as been deallocated.\n");
173 jgs 150 #endif
174     }
175     }
176     }
177 gross 1639 void Paso_SystemMatrix_startCollect(Paso_SystemMatrix* A,const double* in)
178 ksteube 1312 {
179 gross 1552 Paso_SystemMatrix_startColCollect(A,in);
180 ksteube 1312 }
181     double* Paso_SystemMatrix_finishCollect(Paso_SystemMatrix* A)
182     {
183 gross 1552 return Paso_SystemMatrix_finishColCollect(A);
184 ksteube 1312 }
185    
186 gross 1639 void Paso_SystemMatrix_startColCollect(Paso_SystemMatrix* A,const double* in)
187 gross 1552 {
188     Paso_Coupler_startCollect(A->col_coupler, in);
189     }
190     double* Paso_SystemMatrix_finishColCollect(Paso_SystemMatrix* A)
191     {
192     Paso_Coupler_finishCollect(A->col_coupler);
193     return A->col_coupler->recv_buffer;
194     }
195 gross 1639 void Paso_SystemMatrix_startRowCollect(Paso_SystemMatrix* A,const double* in)
196 gross 1552 {
197     Paso_Coupler_startCollect(A->row_coupler, in);
198     }
199     double* Paso_SystemMatrix_finishRowCollect(Paso_SystemMatrix* A)
200     {
201     Paso_Coupler_finishCollect(A->row_coupler);
202     return A->row_coupler->recv_buffer;
203     }
204    
205 gross 1407 dim_t Paso_SystemMatrix_getTotalNumRows(const Paso_SystemMatrix* A){
206 ksteube 1312 return A->mainBlock->numRows * A->row_block_size;
207     }
208    
209 gross 1407 dim_t Paso_SystemMatrix_getTotalNumCols(const Paso_SystemMatrix* A){
210 ksteube 1312 return A->mainBlock->numCols * A->col_block_size;
211     }
212     dim_t Paso_SystemMatrix_getGlobalNumRows(Paso_SystemMatrix* A) {
213     if (A->type & MATRIX_FORMAT_CSC) {
214     return Paso_Distribution_getGlobalNumComponents(A->pattern->input_distribution);
215     } else {
216     return Paso_Distribution_getGlobalNumComponents(A->pattern->output_distribution);
217     }
218     }
219     dim_t Paso_SystemMatrix_getGlobalNumCols(Paso_SystemMatrix* A) {
220     if (A->type & MATRIX_FORMAT_CSC) {
221     return Paso_Distribution_getGlobalNumComponents(A->pattern->output_distribution);
222     } else {
223     return Paso_Distribution_getGlobalNumComponents(A->pattern->input_distribution);
224     }
225    
226     }
227 gross 1361 dim_t Paso_SystemMatrix_getNumOutput(Paso_SystemMatrix* A) {
228     return Paso_SystemMatrixPattern_getNumOutput(A->pattern);
229     }
230 ksteube 1312

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26