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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1811 - (show annotations)
Thu Sep 25 23:11:13 2008 UTC (11 years, 1 month ago) by ksteube
File MIME type: text/plain
File size: 8698 byte(s)
Copyright updated in all files

1
2 /*******************************************************
3 *
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
14
15 /**************************************************************/
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
37 double time0;
38 Paso_SystemMatrix*out=NULL;
39 dim_t n_norm,i;
40 Paso_SystemMatrixType pattern_format_out;
41
42 Paso_resetError();
43 time0=Paso_timer();
44 out=MEMALLOC(1,Paso_SystemMatrix);
45 if (! Paso_checkPtr(out)) {
46 out->type=type;
47 out->pattern=NULL;
48 out->row_distribution=NULL;
49 out->col_distribution=NULL;
50 out->mpi_info=Paso_MPIInfo_getReference(pattern->mpi_info);
51 out->row_coupler=NULL;
52 out->col_coupler=NULL;
53 out->mainBlock=NULL;
54 out->row_coupleBlock=NULL;
55 out->col_coupleBlock=NULL;
56 out->normalizer_is_valid=FALSE;
57 out->normalizer=NULL;
58 out->solver_package=PASO_PASO;
59 out->solver=NULL;
60 out->trilinos_data=NULL;
61 out->reference_counter=1;
62
63
64 pattern_format_out= (type & MATRIX_FORMAT_OFFSET1)? PATTERN_FORMAT_OFFSET1: PATTERN_FORMAT_DEFAULT;
65 /* ====== compressed sparse columns === */
66 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 if ((type & MATRIX_FORMAT_BLK1) || row_block_size!=col_block_size || col_block_size>3) {
71 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,pattern_format_out,col_block_size,row_block_size);
72 out->row_block_size=1;
73 out->col_block_size=1;
74 } else {
75 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,pattern_format_out,1,1);
76 out->row_block_size=row_block_size;
77 out->col_block_size=col_block_size;
78 }
79 }
80 out->row_distribution=Paso_Distribution_getReference(out->pattern->input_distribution);
81 out->col_distribution=Paso_Distribution_getReference(out->pattern->output_distribution);
82 } else {
83 /* ====== compressed sparse row === */
84 if (type & MATRIX_FORMAT_SYM) {
85 Paso_setError(TYPE_ERROR,"Generation of matrix pattern for symmetric CSR is not implemented yet.");
86 } else {
87 if ((type & MATRIX_FORMAT_BLK1) || row_block_size!=col_block_size || col_block_size>3) {
88 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,pattern_format_out,row_block_size,col_block_size);
89 out->row_block_size=1;
90 out->col_block_size=1;
91 } else {
92 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,pattern_format_out,1,1);
93 out->row_block_size=row_block_size;
94 out->col_block_size=col_block_size;
95 }
96 }
97 out->row_distribution=Paso_Distribution_getReference(out->pattern->output_distribution);
98 out->col_distribution=Paso_Distribution_getReference(out->pattern->input_distribution);
99 }
100 out->logical_row_block_size=row_block_size;
101 out->logical_col_block_size=col_block_size;
102 out->logical_block_size=out->logical_row_block_size*out->logical_block_size;
103 out->block_size=out->row_block_size*out->col_block_size;
104 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 /* 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 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 /* 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 }
130 /* all done: */
131 if (! Paso_noError()) {
132 Paso_SystemMatrix_free(out);
133 return NULL;
134 } else {
135 #ifdef Paso_TRACE
136 printf("timing: system matrix %.4e sec\n",Paso_timer()-time0);
137 printf("Paso_SystemMatrix_alloc: system matrix has been allocated.\n");
138 #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 return in;
148 }
149
150 /* deallocates a SystemMatrix: */
151
152 void Paso_SystemMatrix_free(Paso_SystemMatrix* in) {
153 if (in!=NULL) {
154 in->reference_counter--;
155 if (in->reference_counter<=0) {
156 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 Paso_Coupler_free(in->row_coupler);
161 Paso_Coupler_free(in->col_coupler);
162 Paso_SparseMatrix_free(in->mainBlock);
163 Paso_SparseMatrix_free(in->col_coupleBlock);
164 Paso_SparseMatrix_free(in->row_coupleBlock);
165 MEMFREE(in->normalizer);
166 Paso_solve_free(in);
167 #ifdef TRILINOS
168 Paso_TRILINOS_free(in->trilinos_data);
169 #endif
170 MEMFREE(in);
171 #ifdef Paso_TRACE
172 printf("Paso_SystemMatrix_free: system matrix as been deallocated.\n");
173 #endif
174 }
175 }
176 }
177 void Paso_SystemMatrix_startCollect(Paso_SystemMatrix* A,const double* in)
178 {
179 Paso_SystemMatrix_startColCollect(A,in);
180 }
181 double* Paso_SystemMatrix_finishCollect(Paso_SystemMatrix* A)
182 {
183 return Paso_SystemMatrix_finishColCollect(A);
184 }
185
186 void Paso_SystemMatrix_startColCollect(Paso_SystemMatrix* A,const double* in)
187 {
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 void Paso_SystemMatrix_startRowCollect(Paso_SystemMatrix* A,const double* in)
196 {
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 dim_t Paso_SystemMatrix_getTotalNumRows(const Paso_SystemMatrix* A){
206 return A->mainBlock->numRows * A->row_block_size;
207 }
208
209 dim_t Paso_SystemMatrix_getTotalNumCols(const Paso_SystemMatrix* A){
210 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 dim_t Paso_SystemMatrix_getNumOutput(Paso_SystemMatrix* A) {
228 return Paso_SystemMatrixPattern_getNumOutput(A->pattern);
229 }
230

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26