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

Diff of /branches/RW_WIN32/finley/src/CPPAdapter/SystemMatrixAdapter.cpp

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 108 by jgs, Thu Jan 27 06:21:59 2005 UTC revision 150 by jgs, Thu Sep 15 03:44:45 2005 UTC
# Line 13  Line 13 
13   $Id$   $Id$
14  */  */
15  extern "C" {  extern "C" {
16  #include "finley/finleyC/System.h"  #include "paso/SystemMatrix.h"
17    #include "paso/Options.h"
18  }  }
19  #include "escript/Data/Data.h"  #include "escript/Data/Data.h"
20    #include "escript/Data/UtilC.h"
21  #include "finley/CPPAdapter/SystemMatrixAdapter.h"  #include "finley/CPPAdapter/SystemMatrixAdapter.h"
22  #include "finley/CPPAdapter/FinleyAdapterException.h"  #include "finley/CPPAdapter/FinleyAdapterException.h"
23  #include "finley/CPPAdapter/FinleyError.h"  #include "finley/CPPAdapter/FinleyError.h"
# Line 38  SystemMatrixAdapter::SystemMatrixAdapter Line 40  SystemMatrixAdapter::SystemMatrixAdapter
40     throw FinleyAdapterException("Error - Illegal to generate default SystemMatrixAdapter.");     throw FinleyAdapterException("Error - Illegal to generate default SystemMatrixAdapter.");
41  }  }
42    
43  SystemMatrixAdapter::SystemMatrixAdapter(Finley_SystemMatrix* system_matrix,  SystemMatrixAdapter::SystemMatrixAdapter(Paso_SystemMatrix* system_matrix,
44                                           const int row_blocksize,                                           const int row_blocksize,
45                                           const escript::FunctionSpace& row_functionspace,                                           const escript::FunctionSpace& row_functionspace,
46                                           const int column_blocksize,                                           const int column_blocksize,
# Line 51  AbstractSystemMatrix(row_blocksize,row_f Line 53  AbstractSystemMatrix(row_blocksize,row_f
53  SystemMatrixAdapter::~SystemMatrixAdapter()  SystemMatrixAdapter::~SystemMatrixAdapter()
54  {  {
55      if (m_system_matrix.unique()) {      if (m_system_matrix.unique()) {
56          Finley_SystemMatrix* mat=m_system_matrix.get();          Paso_SystemMatrix* mat=m_system_matrix.get();
57          Finley_SystemMatrix_dealloc(mat);          Paso_SystemMatrix_dealloc(mat);
58      }      }
59  }  }
60    
61  Finley_SystemMatrix* SystemMatrixAdapter::getFinley_SystemMatrix() const  Paso_SystemMatrix* SystemMatrixAdapter::getPaso_SystemMatrix() const
62  {  {
63     return m_system_matrix.get();     return m_system_matrix.get();
64  }  }
65    
66  void SystemMatrixAdapter::ypAx(escript::Data& y, const escript::Data& x) const  void SystemMatrixAdapter::ypAx(escript::Data& y, const escript::Data& x) const
67  {  {
68     Finley_SystemMatrix* mat=getFinley_SystemMatrix();     Paso_SystemMatrix* mat=getPaso_SystemMatrix();
69     Finley_SystemMatrixVector(&(y.getDataC()),mat,&(x.getDataC()));  
70     checkFinleyError();    if (!x.isExpanded()) {
71       throw FinleyAdapterException("matrix vector product : input Data object has to be expanded");
72      } else if (!y.isExpanded()) {
73       throw FinleyAdapterException("matrix vector product : output Data object has to be expanded");
74      } else if ( x.getDataPointSize()  != getColumnBlockSize()) {
75       throw FinleyAdapterException("matrix vector product : column block size does not match the number of components in input.");
76      } else if (y.getDataPointSize() != getRowBlockSize()) {
77       throw FinleyAdapterException("matrix vector product : row block size does not match the number of components in output.");
78      } else if ( x.getFunctionSpace()  != getColumnFunctionSpace()) {
79       throw FinleyAdapterException("matrix vector product : column function space and function space of input don't match.");
80      } else if (y.getFunctionSpace() != getRowFunctionSpace()) {
81       throw FinleyAdapterException("matrix vector product : row function space and function space of output don't match.");
82      }
83      double* x_dp=((escript::Data)x).getSampleData(0);
84      double* y_dp=y.getSampleData(0);
85      Paso_SystemMatrix_MatrixVector(1., mat,x_dp, 1.,y_dp);
86      checkPasoError();
87    }
88    
89    int SystemMatrixAdapter::mapOptionToPaso(const int option)  {
90       switch (option) {
91           case  ESCRIPT_DEFAULT:
92              return PASO_DEFAULT;
93           case  ESCRIPT_DIRECT:
94              return PASO_DIRECT;
95           case  ESCRIPT_CHOLEVSKY:
96              return PASO_CHOLEVSKY;
97           case  ESCRIPT_PCG:
98              return PASO_PCG;
99           case  ESCRIPT_CR:
100              return PASO_CR;
101           case  ESCRIPT_CGS:
102              return PASO_CGS;
103           case  ESCRIPT_BICGSTAB:
104              return PASO_BICGSTAB;
105           case  ESCRIPT_SSOR:
106              return PASO_SSOR;
107           case  ESCRIPT_ILU0:
108              return PASO_ILU0;
109           case  ESCRIPT_ILUT:
110              return PASO_ILUT;
111           case  ESCRIPT_JACOBI:
112              return PASO_JACOBI;
113           case  ESCRIPT_GMRES:
114              return PASO_GMRES;
115           case  ESCRIPT_PRES20:
116              return PASO_PRES20;
117           case  ESCRIPT_NO_REORDERING:
118              return PASO_NO_REORDERING;
119           case  ESCRIPT_MINIMUM_FILL_IN:
120              return PASO_MINIMUM_FILL_IN;
121           case  ESCRIPT_NESTED_DISSECTION:
122              return PASO_NESTED_DISSECTION;
123           case  ESCRIPT_SCSL:
124              return PASO_SCSL;
125           case  ESCRIPT_MKL:
126              return PASO_MKL;
127           case  ESCRIPT_UMFPACK:
128              return PASO_UMFPACK;
129           case  ESCRIPT_ITERATIVE:
130              return PASO_ITERATIVE;
131           case  ESCRIPT_PASO:
132              return PASO_PASO;
133           case  ESCRIPT_LUMPING:
134              return PASO_LUMPING;
135           default:
136               stringstream temp;
137               temp << "Error - Cannot map option value "<< option << " onto Paso";
138               throw FinleyAdapterException(temp.str());
139        }
140  }  }
141    
142  void SystemMatrixAdapter::setToSolution(escript::Data& out, const escript::Data& in, const boost::python::dict& options) const  void SystemMatrixAdapter::setToSolution(escript::Data& out, const escript::Data& in, const boost::python::dict& options) const
143  {  {
144      Finley_SolverOptions finley_options;      Paso_SystemMatrix* mat=getPaso_SystemMatrix();
145      Finley_SystemMatrix_setDefaults(&finley_options);      Paso_Options paso_options;
146        Paso_Options_setDefaults(&paso_options);
147      // extract options      // extract options
148      #define EXTRACT(__key__,__val__,__type__) if ( options.has_key(__key__)) finley_options.__val__=boost::python::extract<__type__>(options.get(__key__))      #define EXTRACT(__key__,__val__,__type__) if ( options.has_key(__key__)) paso_options.__val__=boost::python::extract<__type__>(options.get(__key__))
149        #define EXTRACT_OPTION(__key__,__val__,__type__) if ( options.has_key(__key__)) paso_options.__val__=mapOptionToPaso(boost::python::extract<__type__>(options.get(__key__)))
150      EXTRACT("verbose",verbose,int);      EXTRACT("verbose",verbose,int);
151      EXTRACT("reordering",reordering,int);      EXTRACT_OPTION("reordering",reordering,int);
152      EXTRACT(ESCRIPT_TOLERANCE_KEY,tolerance,double);      EXTRACT(ESCRIPT_TOLERANCE_KEY,tolerance,double);
153      EXTRACT(ESCRIPT_METHOD_KEY,method,int);      EXTRACT_OPTION(ESCRIPT_METHOD_KEY,method,int);
154      EXTRACT(ESCRIPT_SYMMETRY_KEY,symmetric,int);      EXTRACT(ESCRIPT_SYMMETRY_KEY,symmetric,int);
155      EXTRACT("preconditioner",preconditioner,int);      EXTRACT_OPTION(ESCRIPT_PACKAGE_KEY,package,int);
156        EXTRACT_OPTION("preconditioner",preconditioner,int);
157      EXTRACT("iter_max",iter_max,int);      EXTRACT("iter_max",iter_max,int);
158      EXTRACT("drop_tolerance",drop_tolerance,double);      EXTRACT("drop_tolerance",drop_tolerance,double);
159      EXTRACT("drop_storage",drop_storage,double);      EXTRACT("drop_storage",drop_storage,double);
160      EXTRACT("truncation",truncation,double);      EXTRACT("truncation",truncation,double);
161      EXTRACT("restart",restart,double);      EXTRACT("restart",restart,double);
162      #undef EXTRACT      #undef EXTRACT
163      Finley_SystemMatrix_solve(getFinley_SystemMatrix(),&(out.getDataC()),&(in.getDataC()),&finley_options);      #undef EXTRACT_OPTION
164      checkFinleyError();      if (! out.isExpanded()) {
165         throw FinleyAdapterException("solve : result Data object has to be expanded");
166        } else if (!in.isExpanded()) {
167         throw FinleyAdapterException("solve : right hand side Data object has to be expanded");
168        } else if ( out.getDataPointSize()  != getColumnBlockSize()) {
169         throw FinleyAdapterException("solve : column block size does not match the number of components of solution.");
170        } else if ( in.getDataPointSize() != getRowBlockSize()) {
171         throw FinleyAdapterException("solve : row block size does not match the number of components of  right hand side.");
172        } else if ( out.getFunctionSpace()  != getColumnFunctionSpace()) {
173         throw FinleyAdapterException("solve : column function space and function space of solution don't match.");
174        } else if (in.getFunctionSpace() != getRowFunctionSpace()) {
175         throw FinleyAdapterException("solve : row function space and function space of right hand side don't match.");
176        }
177        double* out_dp=out.getSampleData(0);
178        double* in_dp=((escript::Data) in).getSampleData(0);
179        Paso_solve(mat,out_dp,in_dp,&paso_options);
180        checkPasoError();
181  }  }
182    
183  void SystemMatrixAdapter::nullifyRowsAndCols(const escript::Data& row_q, const escript::Data& col_q, const double mdv) const  void SystemMatrixAdapter::nullifyRowsAndCols(const escript::Data& row_q, const escript::Data& col_q, const double mdv) const
184  {  {
185      Finley_SystemMatrix* system_matrix_ptr = getFinley_SystemMatrix();      Paso_SystemMatrix* mat = getPaso_SystemMatrix();
186      escriptDataC row_qC = row_q.getDataC();      if (! col_q.isExpanded()) {
187      escriptDataC col_qC = col_q.getDataC();       throw FinleyAdapterException("nullifyRowsAndCols : column mask has to be expanded");
188        } else if (!row_q.isExpanded()) {
189      Finley_SystemMatrix_nullifyRowsAndCols(system_matrix_ptr, &row_qC, &col_qC, mdv);       throw FinleyAdapterException("nullifyRowsAndCols : row mask has to be expanded");
190        } else if ( col_q.getDataPointSize()  != getColumnBlockSize()) {
191         throw FinleyAdapterException("nullifyRowsAndCols : column block size does not match the number of components of column mask.");
192        } else if ( row_q.getDataPointSize() != getRowBlockSize()) {
193         throw FinleyAdapterException("nullifyRowsAndCols : row block size does not match the number of components of row mask.");
194        } else if ( col_q.getFunctionSpace()  != getColumnFunctionSpace()) {
195         throw FinleyAdapterException("nullifyRowsAndCols : column function space and function space of column mask don't match.");
196        } else if (row_q.getFunctionSpace() != getRowFunctionSpace()) {
197         throw FinleyAdapterException("nullifyRowsAndCols : row function space and function space of row mask don't match.");
198        }
199        double* row_q_dp=((escript::Data) row_q).getSampleData(0);
200        double* col_q_dp=((escript::Data) col_q).getSampleData(0);
201        Paso_SystemMatrix_nullifyRowsAndCols(mat,row_q_dp,col_q_dp, mdv);
202        checkPasoError();
203  }  }
204    
205  void SystemMatrixAdapter::saveMM(const std::string& fileName) const  void SystemMatrixAdapter::saveMM(const std::string& fileName) const
206  {  {
207      char fName[fileName.size()+1];      char fName[fileName.size()+1];
208      strcpy(fName,fileName.c_str());      strcpy(fName,fileName.c_str());
209      Finley_SystemMatrix* system_matrix_ptr = getFinley_SystemMatrix();      Paso_SystemMatrix* mat = getPaso_SystemMatrix();
210      Finley_SystemMatrix_saveMM(system_matrix_ptr,fName);      Paso_SystemMatrix_saveMM(mat,fName);
211      checkFinleyError();      checkPasoError();
212  }  }
213    
214  void SystemMatrixAdapter::setValue(const double value) const  void SystemMatrixAdapter::saveHB(const std::string& fileName) const
215  {  {
216     Finley_SystemMatrix* system_matrix_ptr = getFinley_SystemMatrix();      char fName[fileName.size()+1];
217     Finley_SystemMatrix_setValues(system_matrix_ptr,value);      strcpy(fName,fileName.c_str());
218     checkFinleyError();      Paso_SystemMatrix* mat = getPaso_SystemMatrix();
219        Paso_SystemMatrix_saveHB(mat,fName);
220        checkPasoError();
221  }  }
222    
223  void SystemMatrixAdapter::resetSolver() const  void SystemMatrixAdapter::resetValues() const
224  {  {
225     Finley_SystemMatrix* system_matrix_ptr = getFinley_SystemMatrix();     Paso_SystemMatrix* mat = getPaso_SystemMatrix();
226     Finley_SystemMatrix_solve_free(system_matrix_ptr);     Paso_SystemMatrix_setValues(mat,0.);
227     checkFinleyError();     Paso_solve_free(mat);
228       checkPasoError();
229  }  }
                                                                                                                                                       
230    
231  }  // end of namespace  }  // end of namespace

Legend:
Removed from v.108  
changed lines
  Added in v.150

  ViewVC Help
Powered by ViewVC 1.1.26