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

Contents of /trunk-mpi-branch/paso/src/SystemMatrix.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1096 - (show annotations)
Mon Apr 16 22:59:33 2007 UTC (14 years, 1 month ago) by ksteube
File MIME type: text/plain
File size: 8057 byte(s)
MPI implicit solver example run_simplesolve.py now compiling and
running successfully on one CPU of ess.

Adjusted SConscript, removed some debug print statements and removed
some partially implemented TRILINOS calls.


1 /* $Id$ */
2
3
4 /*
5 ********************************************************************************
6 * Copyright 2006 by ACcESS MNRF *
7 * *
8 * http://www.access.edu.au *
9 * Primary Business: Queensland, Australia *
10 * Licensed under the Open Software License version 3.0 *
11 * http://www.opensource.org/licenses/osl-3.0.php *
12 ********************************************************************************
13 */
14
15 /**************************************************************/
16
17 /* Paso: SystemMatrix */
18
19 /**************************************************************/
20
21 /* Copyrights by ACcESS Australia 2003, 2004,2005 */
22 /* Author: gross@access.edu.au */
23
24 /**************************************************************/
25
26 #include "Paso.h"
27 #include "Paso_MPI.h"
28 #include "SystemMatrix.h"
29 #include "TRILINOS.h"
30
31 /**************************************************************/
32
33 /* allocates a SystemMatrix of type type using the given matrix pattern
34 if type is UNKOWN CSR is used.
35 if CSC or CSC_BLK1 is used pattern has to give the CSC pattern.
36 if CSR or CSR_BLK1 is used pattern has to give the CSR pattern.
37 Values are initialized by zero. */
38
39 Paso_SystemMatrix* Paso_SystemMatrix_alloc(Paso_SystemMatrixType type,Paso_SystemMatrixPattern *pattern, int row_block_size, int col_block_size) {
40
41 double time0;
42 Paso_SystemMatrix*out=NULL;
43 Paso_resetError();
44 time0=Paso_timer();
45 dim_t n_norm,i;
46 out=MEMALLOC(1,Paso_SystemMatrix);
47 if (! Paso_checkPtr(out)) {
48 out->pattern=NULL;
49 out->solver_package=PASO_PASO;
50 out->solver=NULL;
51 out->val=NULL;
52 out->reference_counter=1;
53 out->type=type;
54 /* ====== compressed sparse columns === */
55 if (type & MATRIX_FORMAT_CSC) {
56 if (type & MATRIX_FORMAT_SYM) {
57 Paso_setError(TYPE_ERROR,"Generation of matrix pattern for symmetric CSC is not implemented yet.");
58 return NULL;
59 } else {
60 if ((type & MATRIX_FORMAT_BLK1) || row_block_size!=col_block_size || col_block_size>3) {
61 if (type & MATRIX_FORMAT_OFFSET1) {
62 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,PATTERN_FORMAT_OFFSET1,col_block_size,row_block_size);
63 } else {
64 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,PATTERN_FORMAT_DEFAULT,col_block_size,row_block_size);
65 }
66 out->row_block_size=1;
67 out->col_block_size=1;
68 } else {
69 if ( (type & MATRIX_FORMAT_OFFSET1) ==(pattern->type & PATTERN_FORMAT_OFFSET1)) {
70 out->pattern=Paso_SystemMatrixPattern_reference(pattern);
71 } else {
72 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,(type & MATRIX_FORMAT_OFFSET1)? PATTERN_FORMAT_OFFSET1: PATTERN_FORMAT_DEFAULT,1,1);
73 }
74 out->row_block_size=row_block_size;
75 out->col_block_size=col_block_size;
76 }
77 }
78 out->row_distribution=Paso_Distribution_getReference(out->pattern->input_distribution);
79 out->col_distribution=Paso_Distribution_getReference(out->pattern->output_distribution);
80 } else {
81 /* ====== compressed sparse row === */
82 if (type & MATRIX_FORMAT_SYM) {
83 Paso_setError(TYPE_ERROR,"Generation of matrix pattern for symmetric CSR is not implemented yet.");
84 return NULL;
85 } else {
86 if ((type & MATRIX_FORMAT_BLK1) || row_block_size!=col_block_size || col_block_size>3) {
87 if (type & MATRIX_FORMAT_OFFSET1) {
88 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,PATTERN_FORMAT_OFFSET1,row_block_size,col_block_size);
89 } else {
90 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,PATTERN_FORMAT_DEFAULT,row_block_size,col_block_size);
91 }
92 out->row_block_size=1;
93 out->col_block_size=1;
94 } else {
95 if ((type & MATRIX_FORMAT_OFFSET1)==(pattern->type & PATTERN_FORMAT_OFFSET1)) {
96 out->pattern=Paso_SystemMatrixPattern_reference(pattern);
97 } else {
98 out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,(type & MATRIX_FORMAT_OFFSET1)? PATTERN_FORMAT_OFFSET1: PATTERN_FORMAT_DEFAULT,1,1);
99 }
100 out->row_block_size=row_block_size;
101 out->col_block_size=col_block_size;
102 }
103 }
104 out->row_distribution=Paso_Distribution_getReference(out->pattern->output_distribution);
105 out->col_distribution=Paso_Distribution_getReference(out->pattern->input_distribution);
106 }
107 out->logical_row_block_size=row_block_size;
108 out->logical_col_block_size=col_block_size;
109 out->logical_block_size=out->logical_row_block_size*out->logical_block_size;
110 out->block_size=out->row_block_size*out->col_block_size;
111 out->myLen=(size_t)(out->pattern->myLen)*(size_t)(out->block_size);
112
113 out->numRows = out->row_distribution->numComponents;
114 out->myNumRows = out->row_distribution->myNumComponents;
115 out->myFirstRow = out->row_distribution->myFirstComponent;
116 out->maxNumRows = out->row_distribution->maxNumComponents;
117 out->numCols = out->col_distribution->numComponents;
118 out->myNumCols = out->col_distribution->myNumComponents;
119 out->myFirstCol = out->col_distribution->myFirstComponent;
120 out->maxNumCols = out->col_distribution->maxNumComponents;
121 out->mpi_info = Paso_MPIInfo_getReference(out->pattern->mpi_info);
122 /* allocate memory for matrix entries */
123 if (type & MATRIX_FORMAT_CSC) {
124 n_norm = out->maxNumCols * out->col_block_size;
125 } else {
126 n_norm = out->maxNumRows * out->row_block_size;
127 }
128 out->val=MEMALLOC(out->myLen,double);
129 out->normalizer=MEMALLOC(n_norm,double);
130 out->normalizer_is_valid=FALSE;
131 if (! Paso_checkPtr(out->val)) Paso_SystemMatrix_setValues(out,DBLE(0));
132 if (! Paso_checkPtr(out->normalizer)) {
133 #pragma omp parallel for private(i) schedule(static)
134 for (i=0;i<n_norm;++i) out->normalizer[i]=0.;
135 }
136 }
137 /* all done: */
138 if (! Paso_noError()) {
139 Paso_SystemMatrix_dealloc(out);
140 return NULL;
141 } else {
142 #ifdef Paso_TRACE
143 printf("timing: system matrix %.4e sec\n",Paso_timer()-time0);
144 printf("Paso_SystemMatrix_alloc: %ld x %ld system matrix has been allocated.\n",(long)out->numRows,(long)out->numCols);
145 #endif
146 return out;
147 }
148 }
149
150 /* returns a reference to Paso_SystemMatrix in */
151
152 Paso_SystemMatrix* Paso_SystemMatrix_reference(Paso_SystemMatrix* in) {
153 if (in!=NULL) ++(in->reference_counter);
154 return NULL;
155 }
156
157 /* deallocates a SystemMatrix: */
158
159 void Paso_SystemMatrix_dealloc(Paso_SystemMatrix* in) {
160 if (in!=NULL) {
161 in->reference_counter--;
162 if (in->reference_counter<=0) {
163 Paso_SystemMatrixPattern_dealloc(in->pattern);
164 Paso_Distribution_free(in->row_distribution);
165 Paso_Distribution_free(in->col_distribution);
166 Paso_MPIInfo_dealloc(in->mpi_info);
167 MEMFREE(in->val);
168 MEMFREE(in->normalizer);
169 Paso_solve_free(in);
170 MEMFREE(in);
171 #ifdef Paso_TRACE
172 printf("Paso_SystemMatrix_dealloc: system matrix as been deallocated.\n");
173 #endif
174 }
175 }
176 }
177 /*
178 * $Log$
179 * Revision 1.2 2005/09/15 03:44:38 jgs
180 * Merge of development branch dev-02 back to main trunk on 2005-09-15
181 *
182 * Revision 1.1.2.2 2005/09/07 00:59:08 gross
183 * some inconsistent renaming fixed to make the linking work.
184 *
185 * Revision 1.1.2.1 2005/09/05 06:29:47 gross
186 * These files have been extracted from finley to define a stand alone libray for iterative
187 * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
188 * has not been tested yet.
189 *
190 *
191 */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26