/[escript]/trunk/finley/src/CPPAdapter/SystemMatrixAdapter.cpp
ViewVC logotype

Annotation of /trunk/finley/src/CPPAdapter/SystemMatrixAdapter.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 817 - (hide annotations)
Sat Aug 26 03:08:52 2006 UTC (13 years, 11 months ago) by ksteube
File size: 9022 byte(s)
Can now compile and run with MPI on shake71


1 jgs 82 /*
2     ******************************************************************************
3     * *
4     * COPYRIGHT ACcESS 2004 - All Rights Reserved *
5     * *
6     * This software is the property of ACcESS. No part of this code *
7     * may be copied in any form or by any means without the expressed written *
8     * consent of ACcESS. Copying, use or modification of this software *
9     * by any unauthorised person is illegal unless that person has a software *
10     * license agreement with ACcESS. *
11     * *
12     ******************************************************************************
13 jgs 102 $Id$
14 jgs 82 */
15 jgs 472
16 ksteube 817 #ifdef PASO_MPI
17     #include <mpi.h>
18     #endif
19 jgs 203 #include "SystemMatrixAdapter.h"
20 jgs 82
21     using namespace std;
22    
23     namespace finley {
24    
25 jgs 102 struct null_deleter
26     {
27     void operator()(void const *ptr) const
28     {
29     }
30     };
31    
32    
33 jgs 82 SystemMatrixAdapter::SystemMatrixAdapter()
34     {
35     throw FinleyAdapterException("Error - Illegal to generate default SystemMatrixAdapter.");
36     }
37    
38 jgs 150 SystemMatrixAdapter::SystemMatrixAdapter(Paso_SystemMatrix* system_matrix,
39 jgs 82 const int row_blocksize,
40     const escript::FunctionSpace& row_functionspace,
41     const int column_blocksize,
42     const escript::FunctionSpace& column_functionspace):
43 jgs 102 AbstractSystemMatrix(row_blocksize,row_functionspace,column_blocksize,column_functionspace)
44 jgs 82 {
45 jgs 102 m_system_matrix.reset(system_matrix,null_deleter());
46 jgs 82 }
47    
48     SystemMatrixAdapter::~SystemMatrixAdapter()
49     {
50     if (m_system_matrix.unique()) {
51 jgs 150 Paso_SystemMatrix* mat=m_system_matrix.get();
52     Paso_SystemMatrix_dealloc(mat);
53 jgs 82 }
54     }
55    
56 jgs 150 Paso_SystemMatrix* SystemMatrixAdapter::getPaso_SystemMatrix() const
57 jgs 82 {
58     return m_system_matrix.get();
59     }
60    
61 jgs 153 void SystemMatrixAdapter::ypAx(escript::Data& y,escript::Data& x) const
62 jgs 82 {
63 jgs 150 Paso_SystemMatrix* mat=getPaso_SystemMatrix();
64    
65 jgs 153 if ( x.getDataPointSize() != getColumnBlockSize()) {
66 jgs 150 throw FinleyAdapterException("matrix vector product : column block size does not match the number of components in input.");
67     } else if (y.getDataPointSize() != getRowBlockSize()) {
68     throw FinleyAdapterException("matrix vector product : row block size does not match the number of components in output.");
69     } else if ( x.getFunctionSpace() != getColumnFunctionSpace()) {
70     throw FinleyAdapterException("matrix vector product : column function space and function space of input don't match.");
71     } else if (y.getFunctionSpace() != getRowFunctionSpace()) {
72     throw FinleyAdapterException("matrix vector product : row function space and function space of output don't match.");
73     }
74 jgs 153 x.expand();
75     y.expand();
76     double* x_dp=x.getSampleData(0);
77 jgs 150 double* y_dp=y.getSampleData(0);
78     Paso_SystemMatrix_MatrixVector(1., mat,x_dp, 1.,y_dp);
79     checkPasoError();
80 jgs 82 }
81    
82 jgs 150 int SystemMatrixAdapter::mapOptionToPaso(const int option) {
83     switch (option) {
84     case ESCRIPT_DEFAULT:
85     return PASO_DEFAULT;
86     case ESCRIPT_DIRECT:
87     return PASO_DIRECT;
88     case ESCRIPT_CHOLEVSKY:
89     return PASO_CHOLEVSKY;
90     case ESCRIPT_PCG:
91     return PASO_PCG;
92     case ESCRIPT_CR:
93     return PASO_CR;
94     case ESCRIPT_CGS:
95     return PASO_CGS;
96     case ESCRIPT_BICGSTAB:
97     return PASO_BICGSTAB;
98     case ESCRIPT_SSOR:
99     return PASO_SSOR;
100     case ESCRIPT_ILU0:
101     return PASO_ILU0;
102     case ESCRIPT_ILUT:
103     return PASO_ILUT;
104     case ESCRIPT_JACOBI:
105     return PASO_JACOBI;
106     case ESCRIPT_GMRES:
107     return PASO_GMRES;
108     case ESCRIPT_PRES20:
109     return PASO_PRES20;
110     case ESCRIPT_NO_REORDERING:
111     return PASO_NO_REORDERING;
112     case ESCRIPT_MINIMUM_FILL_IN:
113     return PASO_MINIMUM_FILL_IN;
114     case ESCRIPT_NESTED_DISSECTION:
115     return PASO_NESTED_DISSECTION;
116     case ESCRIPT_SCSL:
117     return PASO_SCSL;
118     case ESCRIPT_MKL:
119     return PASO_MKL;
120     case ESCRIPT_UMFPACK:
121     return PASO_UMFPACK;
122     case ESCRIPT_ITERATIVE:
123     return PASO_ITERATIVE;
124     case ESCRIPT_PASO:
125     return PASO_PASO;
126     case ESCRIPT_LUMPING:
127     return PASO_LUMPING;
128 gross 430 case ESCRIPT_AMG:
129     return PASO_AMG;
130     case ESCRIPT_RILU:
131     return PASO_RILU;
132 jgs 150 default:
133     stringstream temp;
134     temp << "Error - Cannot map option value "<< option << " onto Paso";
135     throw FinleyAdapterException(temp.str());
136     }
137     }
138    
139 jgs 153 void SystemMatrixAdapter::setToSolution(escript::Data& out,escript::Data& in, const boost::python::dict& options) const
140 jgs 82 {
141 jgs 150 Paso_SystemMatrix* mat=getPaso_SystemMatrix();
142     Paso_Options paso_options;
143     Paso_Options_setDefaults(&paso_options);
144 jgs 82 // extract options
145 jgs 150 #define EXTRACT(__key__,__val__,__type__) if ( options.has_key(__key__)) paso_options.__val__=boost::python::extract<__type__>(options.get(__key__))
146     #define EXTRACT_OPTION(__key__,__val__,__type__) if ( options.has_key(__key__)) paso_options.__val__=mapOptionToPaso(boost::python::extract<__type__>(options.get(__key__)))
147 jgs 82 EXTRACT("verbose",verbose,int);
148 jgs 150 EXTRACT_OPTION("reordering",reordering,int);
149 jgs 102 EXTRACT(ESCRIPT_TOLERANCE_KEY,tolerance,double);
150 jgs 150 EXTRACT_OPTION(ESCRIPT_METHOD_KEY,method,int);
151 jgs 102 EXTRACT(ESCRIPT_SYMMETRY_KEY,symmetric,int);
152 jgs 150 EXTRACT_OPTION(ESCRIPT_PACKAGE_KEY,package,int);
153     EXTRACT_OPTION("preconditioner",preconditioner,int);
154 jgs 82 EXTRACT("iter_max",iter_max,int);
155     EXTRACT("drop_tolerance",drop_tolerance,double);
156     EXTRACT("drop_storage",drop_storage,double);
157 jgs 102 EXTRACT("truncation",truncation,double);
158     EXTRACT("restart",restart,double);
159 jgs 82 #undef EXTRACT
160 jgs 150 #undef EXTRACT_OPTION
161 jgs 153 if ( out.getDataPointSize() != getColumnBlockSize()) {
162 jgs 150 throw FinleyAdapterException("solve : column block size does not match the number of components of solution.");
163     } else if ( in.getDataPointSize() != getRowBlockSize()) {
164     throw FinleyAdapterException("solve : row block size does not match the number of components of right hand side.");
165     } else if ( out.getFunctionSpace() != getColumnFunctionSpace()) {
166     throw FinleyAdapterException("solve : column function space and function space of solution don't match.");
167     } else if (in.getFunctionSpace() != getRowFunctionSpace()) {
168     throw FinleyAdapterException("solve : row function space and function space of right hand side don't match.");
169     }
170 jgs 153 out.expand();
171     in.expand();
172 jgs 150 double* out_dp=out.getSampleData(0);
173 jgs 153 double* in_dp=in.getSampleData(0);
174 jgs 150 Paso_solve(mat,out_dp,in_dp,&paso_options);
175     checkPasoError();
176 jgs 82 }
177    
178 jgs 153 void SystemMatrixAdapter::nullifyRowsAndCols(escript::Data& row_q,escript::Data& col_q, const double mdv) const
179 jgs 82 {
180 jgs 150 Paso_SystemMatrix* mat = getPaso_SystemMatrix();
181 jgs 153 if ( col_q.getDataPointSize() != getColumnBlockSize()) {
182 jgs 150 throw FinleyAdapterException("nullifyRowsAndCols : column block size does not match the number of components of column mask.");
183     } else if ( row_q.getDataPointSize() != getRowBlockSize()) {
184     throw FinleyAdapterException("nullifyRowsAndCols : row block size does not match the number of components of row mask.");
185     } else if ( col_q.getFunctionSpace() != getColumnFunctionSpace()) {
186     throw FinleyAdapterException("nullifyRowsAndCols : column function space and function space of column mask don't match.");
187     } else if (row_q.getFunctionSpace() != getRowFunctionSpace()) {
188     throw FinleyAdapterException("nullifyRowsAndCols : row function space and function space of row mask don't match.");
189     }
190 jgs 153 row_q.expand();
191     col_q.expand();
192     double* row_q_dp=row_q.getSampleData(0);
193     double* col_q_dp=col_q.getSampleData(0);
194 jgs 150 Paso_SystemMatrix_nullifyRowsAndCols(mat,row_q_dp,col_q_dp, mdv);
195     checkPasoError();
196 jgs 82 }
197    
198 jgs 102 void SystemMatrixAdapter::saveMM(const std::string& fileName) const
199     {
200 woo409 757 char *fName = (fileName.size()+1>0) ? TMPMEMALLOC(fileName.size()+1,char) : (char*)NULL;
201    
202     strcpy(fName,fileName.c_str());
203 jgs 150 Paso_SystemMatrix* mat = getPaso_SystemMatrix();
204     Paso_SystemMatrix_saveMM(mat,fName);
205     checkPasoError();
206 woo409 757 TMPMEMFREE(fName);
207    
208 jgs 102 }
209    
210 jgs 123 void SystemMatrixAdapter::saveHB(const std::string& fileName) const
211     {
212 woo409 757 char *fName = (fileName.size()+1>0) ? TMPMEMALLOC(fileName.size()+1,char) : (char*)NULL;
213    
214 jgs 123 strcpy(fName,fileName.c_str());
215 jgs 150 Paso_SystemMatrix* mat = getPaso_SystemMatrix();
216     Paso_SystemMatrix_saveHB(mat,fName);
217     checkPasoError();
218 woo409 757 TMPMEMFREE(fName);
219    
220 jgs 123 }
221    
222 jgs 149 void SystemMatrixAdapter::resetValues() const
223 jgs 102 {
224 jgs 150 Paso_SystemMatrix* mat = getPaso_SystemMatrix();
225     Paso_SystemMatrix_setValues(mat,0.);
226     Paso_solve_free(mat);
227     checkPasoError();
228 jgs 108 }
229 jgs 102
230 jgs 82 } // end of namespace

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26