/[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 472 - (hide annotations)
Fri Jan 27 01:50:59 2006 UTC (14 years, 8 months ago) by jgs
File size: 8830 byte(s)
rationalise all #includes

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