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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1639 - (show annotations)
Mon Jul 14 08:55:25 2008 UTC (11 years, 2 months ago) by gross
File MIME type: text/plain
File size: 8735 byte(s)


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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26