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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 971 - (hide annotations)
Wed Feb 14 04:40:49 2007 UTC (12 years, 7 months ago) by ksteube
File MIME type: text/plain
File size: 7083 byte(s)
Had to undo commit to new MPI branch. The changes went into the original and
not the branch. The files committed here are exactly the same as revision 969.


1 jgs 150 /* $Id$ */
2    
3 dhawcroft 631
4     /*
5     ********************************************************************************
6 dhawcroft 633 * Copyright 2006 by ACcESS MNRF *
7 dhawcroft 631 * *
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 jgs 150 /**************************************************************/
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 "SystemMatrix.h"
28    
29     /**************************************************************/
30    
31     /* allocates a SystemMatrix of type type using the given matrix pattern
32     if type is UNKOWN CSR is used.
33     if CSC or CSC_BLK1 is used pattern has to give the CSC pattern.
34     if CSR or CSR_BLK1 is used pattern has to give the CSR pattern.
35     Values are initialized by zero. */
36    
37     Paso_SystemMatrix* Paso_SystemMatrix_alloc(Paso_SystemMatrixType type,Paso_SystemMatrixPattern *pattern, int row_block_size, int col_block_size) {
38 gross 432
39 jgs 150 double time0;
40     Paso_SystemMatrix*out=NULL;
41     Paso_resetError();
42     time0=Paso_timer();
43     dim_t n_norm,i;
44     out=MEMALLOC(1,Paso_SystemMatrix);
45     if (! Paso_checkPtr(out)) {
46     out->pattern=NULL;
47 gross 425 out->solver_package=PASO_PASO;
48     out->solver=NULL;
49 jgs 150 out->val=NULL;
50     out->reference_counter=1;
51 gross 415 out->type=type;
52     if (type & MATRIX_FORMAT_CSC) {
53     if (type & MATRIX_FORMAT_SYM) {
54     Paso_setError(TYPE_ERROR,"Generation of matrix pattern for symmetric CSC is not implemented yet.");
55     return NULL;
56     } else {
57 gross 425 if ((type & MATRIX_FORMAT_BLK1) || row_block_size!=col_block_size || col_block_size>3) {
58 gross 415 if (type & MATRIX_FORMAT_OFFSET1) {
59     out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,PATTERN_FORMAT_OFFSET1,col_block_size,row_block_size);
60     } else {
61     out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,PATTERN_FORMAT_DEFAULT,col_block_size,row_block_size);
62     }
63 jgs 150 out->row_block_size=1;
64     out->col_block_size=1;
65 gross 415 } else {
66     if ( (type & MATRIX_FORMAT_OFFSET1) ==(pattern->type & PATTERN_FORMAT_OFFSET1)) {
67     out->pattern=Paso_SystemMatrixPattern_reference(pattern);
68     } else {
69     out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,(type & MATRIX_FORMAT_OFFSET1)? PATTERN_FORMAT_OFFSET1: PATTERN_FORMAT_DEFAULT,1,1);
70     }
71 jgs 150 out->row_block_size=row_block_size;
72     out->col_block_size=col_block_size;
73 gross 415 }
74     }
75     out->num_rows=out->pattern->n_index;
76     out->num_cols=out->pattern->n_ptr;
77     n_norm = out->num_cols * out->col_block_size;
78     } else {
79     if (type & MATRIX_FORMAT_SYM) {
80     Paso_setError(TYPE_ERROR,"Generation of matrix pattern for symmetric CSR is not implemented yet.");
81     return NULL;
82     } else {
83 gross 425 if ((type & MATRIX_FORMAT_BLK1) || row_block_size!=col_block_size || col_block_size>3) {
84 gross 415 if (type & MATRIX_FORMAT_OFFSET1) {
85     out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,PATTERN_FORMAT_OFFSET1,row_block_size,col_block_size);
86     } else {
87     out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,PATTERN_FORMAT_DEFAULT,row_block_size,col_block_size);
88     }
89     out->row_block_size=1;
90     out->col_block_size=1;
91     } else {
92     if ((type & MATRIX_FORMAT_OFFSET1)==(pattern->type & PATTERN_FORMAT_OFFSET1)) {
93     out->pattern=Paso_SystemMatrixPattern_reference(pattern);
94     } else {
95     out->pattern=Paso_SystemMatrixPattern_unrollBlocks(pattern,(type & MATRIX_FORMAT_OFFSET1)? PATTERN_FORMAT_OFFSET1: PATTERN_FORMAT_DEFAULT,1,1);
96     }
97     out->row_block_size=row_block_size;
98     out->col_block_size=col_block_size;
99     }
100     }
101     out->num_rows=out->pattern->n_index;
102     out->num_cols=out->pattern->n_ptr;
103     n_norm = out->num_cols * out->col_block_size;
104     out->num_rows=out->pattern->n_ptr;
105     out->num_cols=out->pattern->n_index;
106     n_norm = out->num_rows * out->row_block_size;
107 jgs 150 }
108     out->logical_row_block_size=row_block_size;
109     out->logical_col_block_size=col_block_size;
110     out->logical_block_size=out->logical_row_block_size*out->logical_block_size;
111     out->block_size=out->row_block_size*out->col_block_size;
112     out->len=(size_t)(out->pattern->len)*(size_t)(out->block_size);
113     /* allocate memory for matrix entries */
114     out->val=MEMALLOC(out->len,double);
115     out->normalizer=MEMALLOC(n_norm,double);
116     out->normalizer_is_valid=FALSE;
117     if (! Paso_checkPtr(out->val)) {
118     Paso_SystemMatrix_setValues(out,DBLE(0));
119     }
120     if (! Paso_checkPtr(out->normalizer)) {
121     #pragma omp parallel for private(i) schedule(static)
122     for (i=0;i<n_norm;++i) out->normalizer[i]=0.;
123     }
124     }
125     /* all done: */
126     if (! Paso_noError()) {
127     Paso_SystemMatrix_dealloc(out);
128     return NULL;
129     } else {
130     #ifdef Paso_TRACE
131     printf("timing: system matrix %.4e sec\n",Paso_timer()-time0);
132     printf("Paso_SystemMatrix_alloc: %ld x %ld system matrix has been allocated.\n",(long)out->num_rows,(long)out->num_cols);
133     #endif
134     return out;
135     }
136     }
137    
138     /* returns a reference to Paso_SystemMatrix in */
139    
140     Paso_SystemMatrix* Paso_SystemMatrix_reference(Paso_SystemMatrix* in) {
141     if (in!=NULL) ++(in->reference_counter);
142     return NULL;
143     }
144    
145     /* deallocates a SystemMatrix: */
146    
147     void Paso_SystemMatrix_dealloc(Paso_SystemMatrix* in) {
148     if (in!=NULL) {
149     in->reference_counter--;
150     if (in->reference_counter<=0) {
151     MEMFREE(in->val);
152     MEMFREE(in->normalizer);
153     Paso_SystemMatrixPattern_dealloc(in->pattern);
154     Paso_solve_free(in);
155     MEMFREE(in);
156     #ifdef Paso_TRACE
157     printf("Paso_SystemMatrix_dealloc: system matrix as been deallocated.\n");
158     #endif
159     }
160     }
161     }
162     /*
163     * $Log$
164     * Revision 1.2 2005/09/15 03:44:38 jgs
165     * Merge of development branch dev-02 back to main trunk on 2005-09-15
166     *
167     * Revision 1.1.2.2 2005/09/07 00:59:08 gross
168     * some inconsistent renaming fixed to make the linking work.
169     *
170     * Revision 1.1.2.1 2005/09/05 06:29:47 gross
171     * These files have been extracted from finley to define a stand alone libray for iterative
172     * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
173     * has not been tested yet.
174     *
175     *
176     */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26