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

Diff of /trunk/finley/src/CPPAdapter/MeshAdapter.cpp

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

trunk/esys2/finley/src/CPPAdapter/MeshAdapter.cpp revision 97 by jgs, Tue Dec 14 05:39:33 2004 UTC trunk/finley/src/CPPAdapter/MeshAdapter.cpp revision 547 by gross, Tue Feb 21 06:10:54 2006 UTC
# Line 12  Line 12 
12   *                                                                            *   *                                                                            *
13   ******************************************************************************   ******************************************************************************
14  */  */
15  extern "C" {  
16  #include "finley/finleyC/Finley.h"  #include "MeshAdapter.h"
17  #include "finley/finleyC/Assemble.h"  
18  #include "finley/finleyC/Mesh.h"  #include "Data.h"
19  #include "finley/finleyC/Finley.h"  #include "DataFactory.h"
 #include "finley/finleyC/System.h"  
 }  
 #include "finley/CPPAdapter/SystemMatrixAdapter.h"  
 #include "finley/CPPAdapter/MeshAdapter.h"  
 #include "finley/CPPAdapter/FinleyError.h"  
 #include "finley/CPPAdapter/FinleyAdapterException.h"  
 #include "escript/Data/FunctionSpaceFactory.h"  
 #include "escript/Data/Data.h"  
 #include "escript/Data/DataArrayView.h"  
 #include "escript/Data/FunctionSpace.h"  
 #include "escript/Data/DataFactory.h"  
 #include <iostream>  
 #include <vector>  
 #include <sstream>  
20    
21  using namespace std;  using namespace std;
22  using namespace escript;  using namespace escript;
23    
24  namespace finley {  namespace finley {
25    
 struct null_deleter  
 {  
   void operator()(void const *ptr) const  
   {  
   }  
 };  
   
26  //  //
27  // define the statics  // define the static constants
28  MeshAdapter::FunctionSpaceNamesMapType MeshAdapter::m_functionSpaceTypeNames;  MeshAdapter::FunctionSpaceNamesMapType MeshAdapter::m_functionSpaceTypeNames;
29  const int MeshAdapter::DegreesOfFreedom=FINLEY_DEGREES_OF_FREEDOM;  const int MeshAdapter::DegreesOfFreedom=FINLEY_DEGREES_OF_FREEDOM;
30  const int MeshAdapter::ReducedDegreesOfFreedom=FINLEY_REDUCED_DEGREES_OF_FREEDOM;  const int MeshAdapter::ReducedDegreesOfFreedom=FINLEY_REDUCED_DEGREES_OF_FREEDOM;
# Line 64  MeshAdapter::MeshAdapter(Finley_Mesh* fi Line 43  MeshAdapter::MeshAdapter(Finley_Mesh* fi
43    // for us.    // for us.
44    m_finleyMesh.reset(finleyMesh,null_deleter());    m_finleyMesh.reset(finleyMesh,null_deleter());
45  }  }
46    
47  //  //
48  // The copy constructor should just increment the use count  // The copy constructor should just increment the use count
49  MeshAdapter::MeshAdapter(const MeshAdapter& in):  MeshAdapter::MeshAdapter(const MeshAdapter& in):
# Line 96  void MeshAdapter::write(const std::strin Line 76  void MeshAdapter::write(const std::strin
76    checkFinleyError();    checkFinleyError();
77  }  }
78    
 // void MeshAdapter::getTagList(int functionSpaceType,  
 //                  int* numTags) const  
 // {  
 //   Finley_Mesh_tagList(m_finleyMesh.get(),functionSpaceType,tagList,numTags);  
 //   return;  
 // }  
   
79  string MeshAdapter::getDescription() const  string MeshAdapter::getDescription() const
80  {  {
81    return "FinleyMesh";    return "FinleyMesh";
# Line 150  int MeshAdapter::getContinuousFunctionCo Line 123  int MeshAdapter::getContinuousFunctionCo
123  {  {
124    return Nodes;    return Nodes;
125  }  }
126    
127  int MeshAdapter::getFunctionCode() const  int MeshAdapter::getFunctionCode() const
128  {  {
129    return Elements;    return Elements;
130  }  }
131    
132  int MeshAdapter::getFunctionOnBoundaryCode() const  int MeshAdapter::getFunctionOnBoundaryCode() const
133  {  {
134    return FaceElements;    return FaceElements;
135  }  }
136    
137  int MeshAdapter::getFunctionOnContactZeroCode() const  int MeshAdapter::getFunctionOnContactZeroCode() const
138  {  {
139    return ContactElementsZero;    return ContactElementsZero;
140  }  }
141    
142  int MeshAdapter::getFunctionOnContactOneCode() const  int MeshAdapter::getFunctionOnContactOneCode() const
143  {  {
144    return ContactElementsOne;    return ContactElementsOne;
# Line 171  int MeshAdapter::getSolutionCode() const Line 148  int MeshAdapter::getSolutionCode() const
148  {  {
149    return DegreesOfFreedom;    return DegreesOfFreedom;
150  }  }
151    
152  int MeshAdapter::getReducedSolutionCode() const  int MeshAdapter::getReducedSolutionCode() const
153  {  {
154    return ReducedDegreesOfFreedom;    return ReducedDegreesOfFreedom;
155  }  }
156    
157  int MeshAdapter::getDiracDeltaFunctionCode() const  int MeshAdapter::getDiracDeltaFunctionCode() const
158  {  {
159    return Points;    return Points;
160  }  }
161    
 int MeshAdapter::getTagFromSampleNo(int functionSpaceType, int sampleNo) const  
 {  
   //  
   // It is assumed the sampleNo has been checked  
   // before calling this function.  
   int* tagList;  
   int numTags;  
   getTagList(functionSpaceType, &tagList, &numTags);  
   return tagList[sampleNo];  
 }  
   
 //  
 // returns a pointer to the tag list of samples of functionSpaceType  
 //  
 void MeshAdapter::getTagList(int functionSpaceType, int** tagList,  
                  int* numTags) const  
 {  
   *tagList=NULL;  
   *numTags=0;  
   Finley_Mesh* mesh=m_finleyMesh.get();  
   switch (functionSpaceType) {  
   case(Nodes):  
     if (mesh->Nodes!=NULL) {  
       *tagList=mesh->Nodes->Tag;  
       *numTags=mesh->Nodes->numNodes;  
     }  
     break;  
   case(Elements):  
     if (mesh->Elements!=NULL) {  
       *tagList=mesh->Elements->Tag;  
       *numTags=mesh->Elements->numElements;  
     }  
     break;  
   case(FaceElements):  
     if (mesh->FaceElements!=NULL) {  
       *tagList=mesh->FaceElements->Tag;  
       *numTags=mesh->FaceElements->numElements;  
     }  
     break;  
   case(Points):  
     if (mesh->Points!=NULL) {  
       *tagList=mesh->Points->Tag;  
       *numTags=mesh->Points->numElements;  
     }  
     break;  
   case(ContactElementsZero):  
     if (mesh->ContactElements!=NULL) {  
       *tagList=mesh->ContactElements->Tag;  
       *numTags=mesh->ContactElements->numElements;  
     }  
     break;  
   case(ContactElementsOne):  
     if (mesh->ContactElements!=NULL) {  
       *tagList=mesh->ContactElements->Tag;  
       *numTags=mesh->ContactElements->numElements;  
     }  
     break;  
   case(DegreesOfFreedom):  
     if (mesh->Nodes!=NULL) {  
       *tagList=NULL;  
       *numTags=0;  
     }  
     break;  
   case(ReducedDegreesOfFreedom):  
     if (mesh->Nodes!=NULL) {  
       *tagList=NULL;  
       *numTags=0;  
     }  
     break;  
   default:  
     stringstream temp;  
     temp << "Error - Invalid function space type: "  
      << functionSpaceType << " for domain: " << getDescription();  
     throw FinleyAdapterException(temp.str());  
     break;  
   }  
   return;  
 }  
   
162  //  //
163  // return the spatial dimension of the Mesh:  // return the spatial dimension of the Mesh:
164  //  //
# Line 268  int MeshAdapter::getDim() const Line 168  int MeshAdapter::getDim() const
168    checkFinleyError();    checkFinleyError();
169    return numDim;    return numDim;
170  }  }
171    
172  //  //
173  // return the number of data points per sample and the number of samples  // return the number of data points per sample and the number of samples
174  // needed to represent data on a parts of the mesh.  // needed to represent data on a parts of the mesh.
# Line 326  pair<int,int> MeshAdapter::getDataShape( Line 227  pair<int,int> MeshAdapter::getDataShape(
227             break;             break;
228        default:        default:
229          stringstream temp;          stringstream temp;
230          temp << "Error - Invalid function space type: "          temp << "Error - Invalid function space type: " << functionSpaceCode << " for domain: " << getDescription();
          << functionSpaceCode << " for domain: " << getDescription();  
231          throw FinleyAdapterException(temp.str());          throw FinleyAdapterException(temp.str());
232          break;          break;
233        }        }
234        return pair<int,int>(numDataPointsPerSample,numSamples);        return pair<int,int>(numDataPointsPerSample,numSamples);
235  }  }
236    
237  //  //
238  // adds linear PDE of second order into a given stiffness matrix and right hand side:  // adds linear PDE of second order into a given stiffness matrix and right hand side:
239  //  //
240  void MeshAdapter::addPDEToSystem(  void MeshAdapter::addPDEToSystem(
241                       SystemMatrixAdapter& mat, Data& rhs,                       SystemMatrixAdapter& mat, escript::Data& rhs,
242                       const Data& A, const Data& B, const Data& C,const  Data& D,const  Data& X,const  Data& Y,                       const escript::Data& A, const escript::Data& B, const escript::Data& C,const  escript::Data& D,const  escript::Data& X,const  escript::Data& Y,
243                       const Data& d, const Data& y,                       const escript::Data& d, const escript::Data& y,
244                       const Data& d_contact,const Data& y_contact) const                       const escript::Data& d_contact,const escript::Data& y_contact) const
245  {  {
246     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
247     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,mat.getFinley_SystemMatrix(),&(rhs.getDataC()),     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,mat.getPaso_SystemMatrix(),&(rhs.getDataC()),
248                         &(A.getDataC()),&(B.getDataC()),&(C.getDataC()),&(D.getDataC()),&(X.getDataC()),&(Y.getDataC()));                         &(A.getDataC()),&(B.getDataC()),&(C.getDataC()),&(D.getDataC()),&(X.getDataC()),&(Y.getDataC()));
249     checkFinleyError();     checkFinleyError();
250     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->FaceElements,     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->FaceElements,
251                    mat.getFinley_SystemMatrix(),                    mat.getPaso_SystemMatrix(),
252                    &(rhs.getDataC()),                    &(rhs.getDataC()),
253                                    &(d.getDataC()),&(y.getDataC()),                                    &(d.getDataC()),&(y.getDataC()),
254                                    Finley_Assemble_handelShapeMissMatch_Mean_out);                                    Finley_Assemble_handelShapeMissMatch_Mean_out);
255     checkFinleyError();     checkFinleyError();
256     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->FaceElements,     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->ContactElements,
257                    mat.getFinley_SystemMatrix(),                    mat.getPaso_SystemMatrix(),
258                    &(rhs.getDataC()),                    &(rhs.getDataC()),
259                                    &(d_contact.getDataC()),                                    &(d_contact.getDataC()),
260                    &(y_contact.getDataC()),                    &(y_contact.getDataC()),
261                                    Finley_Assemble_handelShapeMissMatch_Step_out);                                    Finley_Assemble_handelShapeMissMatch_Step_out);
262     checkFinleyError();     checkFinleyError();
263  }  }
264    
265  //  //
266  // adds linear PDE of second order into the right hand side only  // adds linear PDE of second order into the right hand side only
267  //  //
268  void MeshAdapter::addPDEToRHS( Data& rhs,  void MeshAdapter::addPDEToRHS( escript::Data& rhs,
269                       const  Data& X,const  Data& Y, const Data& y, const Data& y_contact) const                       const  escript::Data& X,const  escript::Data& Y, const escript::Data& y, const escript::Data& y_contact) const
270  {  {
271     Finley_Mesh* mesh=m_finleyMesh.get();     Finley_Mesh* mesh=m_finleyMesh.get();
272     Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,0,&(rhs.getDataC()),0,0,0,0,&(X.getDataC()),&(Y.getDataC()));  
273       // Finley_Assemble_PDE_RHS(mesh->Nodes,mesh->Elements,&(rhs.getDataC()),&(X.getDataC()),&(Y.getDataC()));
274       Finley_Assemble_PDE(mesh->Nodes,mesh->Elements,0,&(rhs.getDataC()),0,0,0,0,&(X.getDataC()),&(Y.getDataC()));
275     checkFinleyError();     checkFinleyError();
276     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->FaceElements,0,&(rhs.getDataC()),0,&(y.getDataC()),  
277                                    Finley_Assemble_handelShapeMissMatch_Mean_out);     // Finley_Assemble_RobinCondition_RHS(mesh->Nodes,mesh->FaceElements,&(rhs.getDataC()),&(y.getDataC()),Finley_Assemble_handelShapeMissMatch_Mean_out);
278  // cout << "Calling :addPDEToRHS." << endl;     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->FaceElements,0,&(rhs.getDataC()),0,&(y.getDataC()),Finley_Assemble_handelShapeMissMatch_Mean_out);
279    
280     checkFinleyError();     checkFinleyError();
281     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->FaceElements,0,&(rhs.getDataC()),0,&(y_contact.getDataC()),     Finley_Assemble_RobinCondition(mesh->Nodes,mesh->ContactElements,0,&(rhs.getDataC()),0,&(y_contact.getDataC()),Finley_Assemble_handelShapeMissMatch_Step_out);
282                                    Finley_Assemble_handelShapeMissMatch_Step_out);     // Finley_Assemble_RobinCondition_RHS(mesh->Nodes,mesh->ContactElements,&(rhs.getDataC()),&(y_contact.getDataC()),Finley_Assemble_handelShapeMissMatch_Step_out);
 // cout << "Calling :addPDEToRHS." << endl;  
283     checkFinleyError();     checkFinleyError();
284  }  }
285    
286  //  //
287  // interpolates data between different function spaces:  // interpolates data between different function spaces:
288  //  //
289  void MeshAdapter::interpolateOnDomain(Data& target,const Data& in) const  void MeshAdapter::interpolateOnDomain(escript::Data& target,const escript::Data& in) const
290  {  {
291    const MeshAdapter& inDomain=dynamic_cast<const MeshAdapter&>(in.getFunctionSpace().getDomain());    const MeshAdapter& inDomain=dynamic_cast<const MeshAdapter&>(in.getFunctionSpace().getDomain());
292    const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());    const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
# Line 415  void MeshAdapter::interpolateOnDomain(Da Line 320  void MeshAdapter::interpolateOnDomain(Da
320                 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));                 Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));
321                 break;                 break;
322             default:             default:
323                 Finley_ErrorCode=TYPE_ERROR;                 stringstream temp;
324                 sprintf(Finley_ErrorMsg,"Interpolation on Domain: Finley does not know anything about function space type %d",target.getFunctionSpace().getTypeCode());                 temp << "Error - Interpolation on Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
325                   throw FinleyAdapterException(temp.str());
326                 break;                 break;
327          }          }
328          break;          break;
# Line 424  void MeshAdapter::interpolateOnDomain(Da Line 330  void MeshAdapter::interpolateOnDomain(Da
330          if (target.getFunctionSpace().getTypeCode()==Elements) {          if (target.getFunctionSpace().getTypeCode()==Elements) {
331             Finley_Assemble_CopyElementData(mesh->Elements,&(target.getDataC()),&(in.getDataC()));             Finley_Assemble_CopyElementData(mesh->Elements,&(target.getDataC()),&(in.getDataC()));
332          } else {          } else {
333             Finley_ErrorCode=TYPE_ERROR;             throw FinleyAdapterException("Error - No interpolation with data on elements possible.");
            sprintf(Finley_ErrorMsg,"No interpolation with data on elements possible.");  
334          }          }
335          break;          break;
336       case(FaceElements):       case(FaceElements):
337          if (target.getFunctionSpace().getTypeCode()==FaceElements) {          if (target.getFunctionSpace().getTypeCode()==FaceElements) {
338             Finley_Assemble_CopyElementData(mesh->FaceElements,&(target.getDataC()),&(in.getDataC()));             Finley_Assemble_CopyElementData(mesh->FaceElements,&(target.getDataC()),&(in.getDataC()));
339          } else {          } else {
340             Finley_ErrorCode=TYPE_ERROR;             throw FinleyAdapterException("Error - No interpolation with data on face elements possible.");
            sprintf(Finley_ErrorMsg,"No interpolation with data on face elements possible.");  
341             break;             break;
342         }         }
343       case(Points):       case(Points):
344          if (target.getFunctionSpace().getTypeCode()==Points) {          if (target.getFunctionSpace().getTypeCode()==Points) {
345             Finley_Assemble_CopyElementData(mesh->Points,&(target.getDataC()),&(in.getDataC()));             Finley_Assemble_CopyElementData(mesh->Points,&(target.getDataC()),&(in.getDataC()));
346          } else {          } else {
347             Finley_ErrorCode=TYPE_ERROR;             throw FinleyAdapterException("Error - No interpolation with data on points possible.");
            sprintf(Finley_ErrorMsg,"No interpolation with data on points possible.");  
348             break;             break;
349          }          }
350          break;          break;
# Line 450  void MeshAdapter::interpolateOnDomain(Da Line 353  void MeshAdapter::interpolateOnDomain(Da
353          if (target.getFunctionSpace().getTypeCode()==ContactElementsZero || target.getFunctionSpace().getTypeCode()==ContactElementsOne) {          if (target.getFunctionSpace().getTypeCode()==ContactElementsZero || target.getFunctionSpace().getTypeCode()==ContactElementsOne) {
354             Finley_Assemble_CopyElementData(mesh->ContactElements,&(target.getDataC()),&(in.getDataC()));             Finley_Assemble_CopyElementData(mesh->ContactElements,&(target.getDataC()),&(in.getDataC()));
355          } else {          } else {
356             Finley_ErrorCode=TYPE_ERROR;             throw FinleyAdapterException("Error - No interpolation with data on contact elements possible.");
            sprintf(Finley_ErrorMsg,"No interpolation with data on contact elements possible.");  
357             break;             break;
358          }          }
359          break;          break;
# Line 476  void MeshAdapter::interpolateOnDomain(Da Line 378  void MeshAdapter::interpolateOnDomain(Da
378                Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));                Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));
379               break;               break;
380             default:             default:
381               Finley_ErrorCode=TYPE_ERROR;               stringstream temp;
382               sprintf(Finley_ErrorMsg,"Interpolation On Domain: Finley does not know anything about function space type %d",target.getFunctionSpace().getTypeCode());               temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
383                 throw FinleyAdapterException(temp.str());
384               break;               break;
385          }          }
386          break;          break;
387       case(ReducedDegreesOfFreedom):       case(ReducedDegreesOfFreedom):
388         switch(target.getFunctionSpace().getTypeCode()) {         switch(target.getFunctionSpace().getTypeCode()) {
389            case(ReducedDegreesOfFreedom):            case(ReducedDegreesOfFreedom):
           case(Nodes):  
390               Finley_Assemble_CopyNodalData(mesh->Nodes,&(target.getDataC()),&(in.getDataC()));               Finley_Assemble_CopyNodalData(mesh->Nodes,&(target.getDataC()),&(in.getDataC()));
391               break;               break;
392            case(Elements):            case(Elements):
# Line 500  void MeshAdapter::interpolateOnDomain(Da Line 402  void MeshAdapter::interpolateOnDomain(Da
402            case(ContactElementsOne):            case(ContactElementsOne):
403               Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));               Finley_Assemble_interpolate(mesh->Nodes,mesh->ContactElements,&(in.getDataC()),&(target.getDataC()));
404               break;               break;
405              case(Nodes):
406                 throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to mesh nodes.");
407                 break;
408            case(DegreesOfFreedom):            case(DegreesOfFreedom):
409               Finley_ErrorCode=TYPE_ERROR;               throw FinleyAdapterException("Error - Finley does not support interpolation from reduced degrees of freedom to degrees of freedom");
              sprintf(Finley_ErrorMsg,"Finley does not support interpolation from reduced degrees of freedom to degrees of freedom");  
410               break;               break;
411            default:            default:
412               Finley_ErrorCode=TYPE_ERROR;               stringstream temp;
413               sprintf(Finley_ErrorMsg,"Interpolation On Domain: Finley does not know anything about function space type %d",target.getFunctionSpace().getTypeCode());               temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << target.getFunctionSpace().getTypeCode();
414                 throw FinleyAdapterException(temp.str());
415               break;               break;
416         }         }
417         break;         break;
418       default:       default:
419          Finley_ErrorCode=TYPE_ERROR;          stringstream temp;
420          sprintf(Finley_ErrorMsg,"Interpolation On Domain: Finley does not know anything about function space type %d",in.getFunctionSpace().getTypeCode());          temp << "Error - Interpolation On Domain: Finley does not know anything about function space type %d" << in.getFunctionSpace().getTypeCode();
421            throw FinleyAdapterException(temp.str());
422          break;          break;
423    }    }
424    checkFinleyError();    checkFinleyError();
# Line 521  void MeshAdapter::interpolateOnDomain(Da Line 427  void MeshAdapter::interpolateOnDomain(Da
427  //  //
428  // copies the locations of sample points into x:  // copies the locations of sample points into x:
429  //  //
430  void MeshAdapter::setToX(Data& arg) const  void MeshAdapter::setToX(escript::Data& arg) const
431  {  {
432    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
433    if (argDomain!=*this)    if (argDomain!=*this)
# Line 531  void MeshAdapter::setToX(Data& arg) cons Line 437  void MeshAdapter::setToX(Data& arg) cons
437    if (arg.getFunctionSpace().getTypeCode()==Nodes) {    if (arg.getFunctionSpace().getTypeCode()==Nodes) {
438       Finley_Assemble_NodeCoordinates(mesh->Nodes,&(arg.getDataC()));       Finley_Assemble_NodeCoordinates(mesh->Nodes,&(arg.getDataC()));
439    } else {    } else {
440       Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);       escript::Data tmp_data=Vector(0.0,continuousFunction(asAbstractContinuousDomain()),true);
441       Finley_Assemble_NodeCoordinates(mesh->Nodes,&(tmp_data.getDataC()));       Finley_Assemble_NodeCoordinates(mesh->Nodes,&(tmp_data.getDataC()));
442       // this is then interpolated onto arg:       // this is then interpolated onto arg:
443       interpolateOnDomain(arg,tmp_data);       interpolateOnDomain(arg,tmp_data);
444    }    }
445    checkFinleyError();    checkFinleyError();
446  }  }
447    
448  //  //
449  // return the normal vectors at the location of data points as a Data object:  // return the normal vectors at the location of data points as a Data object:
450  //  //
451  void MeshAdapter::setToNormal(Data& normal) const  void MeshAdapter::setToNormal(escript::Data& normal) const
452  {  {
453    const MeshAdapter& normalDomain=dynamic_cast<const MeshAdapter&>(normal.getFunctionSpace().getDomain());    const MeshAdapter& normalDomain=dynamic_cast<const MeshAdapter&>(normal.getFunctionSpace().getDomain());
454    if (normalDomain!=*this)    if (normalDomain!=*this)
# Line 549  void MeshAdapter::setToNormal(Data& norm Line 456  void MeshAdapter::setToNormal(Data& norm
456    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
457    switch(normal.getFunctionSpace().getTypeCode()) {    switch(normal.getFunctionSpace().getTypeCode()) {
458      case(Nodes):      case(Nodes):
459        Finley_ErrorCode=TYPE_ERROR;        throw FinleyAdapterException("Error - Finley does not support surface normal vectors for nodes");
       sprintf(Finley_ErrorMsg,"Finley does not support surface normal vectors for nodes");  
460        break;        break;
461      case(Elements):      case(Elements):
462        Finley_ErrorCode=TYPE_ERROR;        throw FinleyAdapterException("Error - Finley does not support surface normal vectors for elements");
       sprintf(Finley_ErrorMsg,"Finley does not support surface normal vectors for elements");  
463        break;        break;
464      case (FaceElements):      case (FaceElements):
465        Finley_Assemble_setNormal(mesh->Nodes,mesh->FaceElements,&(normal.getDataC()));        Finley_Assemble_setNormal(mesh->Nodes,mesh->FaceElements,&(normal.getDataC()));
466        break;        break;
467      case(Points):      case(Points):
468        Finley_ErrorCode=TYPE_ERROR;        throw FinleyAdapterException("Error - Finley does not support surface normal vectors for point elements");
       sprintf(Finley_ErrorMsg,"Finley does not support surface normal vectors for point elements");  
469        break;        break;
470      case (ContactElementsOne):      case (ContactElementsOne):
471      case (ContactElementsZero):      case (ContactElementsZero):
472        Finley_Assemble_setNormal(mesh->Nodes,mesh->ContactElements,&(normal.getDataC()));        Finley_Assemble_setNormal(mesh->Nodes,mesh->ContactElements,&(normal.getDataC()));
473        break;        break;
474      case(DegreesOfFreedom):      case(DegreesOfFreedom):
475        Finley_ErrorCode=TYPE_ERROR;        throw FinleyAdapterException("Error - Finley does not support surface normal vectors for degrees of freedom.");
       sprintf(Finley_ErrorMsg,"Finley does not support surface normal vectors for degrees of freedom.");  
476        break;        break;
477      case(ReducedDegreesOfFreedom):      case(ReducedDegreesOfFreedom):
478        Finley_ErrorCode=TYPE_ERROR;        throw FinleyAdapterException("Error - Finley does not support surface normal vectors for reduced degrees of freedom.");
       sprintf(Finley_ErrorMsg,"Finley does not support surface normal vectors for reduced degrees of freedom.");  
479        break;        break;
480      default:      default:
481        Finley_ErrorCode=TYPE_ERROR;        stringstream temp;
482        sprintf(Finley_ErrorMsg,"Normal Vectors: Finley does not know anything about function space type %d",normal.getFunctionSpace().getTypeCode());        temp << "Error - Normal Vectors: Finley does not know anything about function space type " << normal.getFunctionSpace().getTypeCode();
483          throw FinleyAdapterException(temp.str());
484        break;        break;
485    }    }
486    checkFinleyError();    checkFinleyError();
487  }  }
488    
489  //  //
490  // interpolates data to other domain:  // interpolates data to other domain:
491  //  //
492  void MeshAdapter::interpolateACross(Data& target,const Data& source) const  void MeshAdapter::interpolateACross(escript::Data& target,const escript::Data& source) const
493  {  {
494    const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());    const MeshAdapter& targetDomain=dynamic_cast<const MeshAdapter&>(target.getFunctionSpace().getDomain());
495    if (targetDomain!=*this)    if (targetDomain!=*this)
496       throw FinleyAdapterException("Error - Illegal domain of interpolation target");       throw FinleyAdapterException("Error - Illegal domain of interpolation target");
497    
498    Finley_ErrorCode=SYSTEM_ERROR;    throw FinleyAdapterException("Error - Finley does not allow interpolation across domains yet.");
   sprintf(Finley_ErrorMsg,"Finley does not allow interpolation across domains yet.");  
   checkFinleyError();  
499  }  }
500    
501  //  //
502  // calculates the integral of a function defined of arg:  // calculates the integral of a function defined of arg:
503  //  //
504  void MeshAdapter::setToIntegrals(std::vector<double>& integrals,const Data& arg) const  void MeshAdapter::setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const
505  {  {
506    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
507    if (argDomain!=*this)    if (argDomain!=*this)
# Line 607  void MeshAdapter::setToIntegrals(std::ve Line 510  void MeshAdapter::setToIntegrals(std::ve
510    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
511    switch(arg.getFunctionSpace().getTypeCode()) {    switch(arg.getFunctionSpace().getTypeCode()) {
512       case(Nodes):       case(Nodes):
513          Finley_ErrorCode=TYPE_ERROR;          throw FinleyAdapterException("Error - Integral of data on nodes is not supported.");
         sprintf(Finley_ErrorMsg,"Integral of data on nodes is not supported.");  
514          break;          break;
515       case(Elements):       case(Elements):
516          Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&(arg.getDataC()),&integrals[0]);          Finley_Assemble_integrate(mesh->Nodes,mesh->Elements,&(arg.getDataC()),&integrals[0]);
# Line 617  void MeshAdapter::setToIntegrals(std::ve Line 519  void MeshAdapter::setToIntegrals(std::ve
519          Finley_Assemble_integrate(mesh->Nodes,mesh->FaceElements,&(arg.getDataC()),&integrals[0]);          Finley_Assemble_integrate(mesh->Nodes,mesh->FaceElements,&(arg.getDataC()),&integrals[0]);
520          break;          break;
521       case(Points):       case(Points):
522          Finley_ErrorCode=TYPE_ERROR;          throw FinleyAdapterException("Error - Integral of data on points is not supported.");
         sprintf(Finley_ErrorMsg,"Integral of data on points is not supported.");  
523          break;          break;
524       case(ContactElementsZero):       case(ContactElementsZero):
525          Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&(arg.getDataC()),&integrals[0]);          Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&(arg.getDataC()),&integrals[0]);
# Line 627  void MeshAdapter::setToIntegrals(std::ve Line 528  void MeshAdapter::setToIntegrals(std::ve
528          Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&(arg.getDataC()),&integrals[0]);          Finley_Assemble_integrate(mesh->Nodes,mesh->ContactElements,&(arg.getDataC()),&integrals[0]);
529          break;          break;
530       case(DegreesOfFreedom):       case(DegreesOfFreedom):
531          Finley_ErrorCode=TYPE_ERROR;          throw FinleyAdapterException("Error - Integral of data on degrees of freedom is not supported.");
         sprintf(Finley_ErrorMsg,"Integral of data on degrees of freedom is not supported.");  
532          break;          break;
533       case(ReducedDegreesOfFreedom):       case(ReducedDegreesOfFreedom):
534          Finley_ErrorCode=TYPE_ERROR;          throw FinleyAdapterException("Error - Integral of data on reduced degrees of freedom is not supported.");
         sprintf(Finley_ErrorMsg,"Integral of data on reduced degrees of freedom is not supported.");  
535          break;          break;
536       default:       default:
537          Finley_ErrorCode=TYPE_ERROR;          stringstream temp;
538          sprintf(Finley_ErrorMsg,"Integrals: Finley does not know anything about function space type %d",arg.getFunctionSpace().getTypeCode());          temp << "Error - Integrals: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
539            throw FinleyAdapterException(temp.str());
540          break;          break;
541    }    }
542    checkFinleyError();    checkFinleyError();
543  }  }
544    
545  //  //
546  // calculates the gradient of arg:  // calculates the gradient of arg:
547  //  //
548  void MeshAdapter::setToGradient(Data& grad,const Data& arg) const  void MeshAdapter::setToGradient(escript::Data& grad,const escript::Data& arg) const
549  {  {
550    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());    const MeshAdapter& argDomain=dynamic_cast<const MeshAdapter&>(arg.getFunctionSpace().getDomain());
551    if (argDomain!=*this)    if (argDomain!=*this)
# Line 656  void MeshAdapter::setToGradient(Data& gr Line 557  void MeshAdapter::setToGradient(Data& gr
557    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
558    switch(grad.getFunctionSpace().getTypeCode()) {    switch(grad.getFunctionSpace().getTypeCode()) {
559         case(Nodes):         case(Nodes):
560            Finley_ErrorCode=TYPE_ERROR;            throw FinleyAdapterException("Error - Gradient at nodes is not supported.");
           sprintf(Finley_ErrorMsg,"Gradient at nodes is not supported.");  
561            break;            break;
562         case(Elements):         case(Elements):
563            Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&(grad.getDataC()),&(arg.getDataC()));            Finley_Assemble_gradient(mesh->Nodes,mesh->Elements,&(grad.getDataC()),&(arg.getDataC()));
# Line 666  void MeshAdapter::setToGradient(Data& gr Line 566  void MeshAdapter::setToGradient(Data& gr
566            Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&(grad.getDataC()),&(arg.getDataC()));            Finley_Assemble_gradient(mesh->Nodes,mesh->FaceElements,&(grad.getDataC()),&(arg.getDataC()));
567            break;            break;
568         case(Points):         case(Points):
569            Finley_ErrorCode=TYPE_ERROR;            throw FinleyAdapterException("Error - Gradient at points is not supported.");
           sprintf(Finley_ErrorMsg,"Gradient at points is not supported.");  
570            break;            break;
571         case(ContactElementsZero):         case(ContactElementsZero):
572            Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&(arg.getDataC()));            Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&(arg.getDataC()));
# Line 676  void MeshAdapter::setToGradient(Data& gr Line 575  void MeshAdapter::setToGradient(Data& gr
575            Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&(arg.getDataC()));            Finley_Assemble_gradient(mesh->Nodes,mesh->ContactElements,&(grad.getDataC()),&(arg.getDataC()));
576            break;            break;
577         case(DegreesOfFreedom):         case(DegreesOfFreedom):
578            Finley_ErrorCode=TYPE_ERROR;            throw FinleyAdapterException("Error - Gradient at degrees of freedom is not supported.");
           sprintf(Finley_ErrorMsg,"Gradient at degrees of freedom is not supported.");  
579            break;            break;
580         case(ReducedDegreesOfFreedom):         case(ReducedDegreesOfFreedom):
581            Finley_ErrorCode=TYPE_ERROR;            throw FinleyAdapterException("Error - Gradient at reduced degrees of freedom is not supported.");
           sprintf(Finley_ErrorMsg,"Gradient at reduced degrees of freedom is not supported.");  
582            break;            break;
583         default:         default:
584            Finley_ErrorCode=TYPE_ERROR;            stringstream temp;
585            sprintf(Finley_ErrorMsg,"Gradient: Finley does not know anything about function space type %d",arg.getFunctionSpace().getTypeCode());            temp << "Error - Gradient: Finley does not know anything about function space type " << arg.getFunctionSpace().getTypeCode();
586              throw FinleyAdapterException(temp.str());
587            break;            break;
588    }    }
589    checkFinleyError();    checkFinleyError();
590  }  }
591    
592  //  //
593  // returns the size of elements:  // returns the size of elements:
594  //  //
595  void MeshAdapter::setToSize(Data& size) const  void MeshAdapter::setToSize(escript::Data& size) const
596  {  {
597    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
598    escriptDataC tmp=size.getDataC();    escriptDataC tmp=size.getDataC();
599    switch(size.getFunctionSpace().getTypeCode()) {    switch(size.getFunctionSpace().getTypeCode()) {
600         case(Nodes):         case(Nodes):
601            Finley_ErrorCode=TYPE_ERROR;            throw FinleyAdapterException("Error - Size of nodes is not supported.");
           sprintf(Finley_ErrorMsg,"Size of nodes is not supported.");  
602            break;            break;
603         case(Elements):         case(Elements):
604            Finley_Assemble_getSize(mesh->Nodes,mesh->Elements,&tmp);            Finley_Assemble_getSize(mesh->Nodes,mesh->Elements,&tmp);
# Line 709  void MeshAdapter::setToSize(Data& size) Line 607  void MeshAdapter::setToSize(Data& size)
607            Finley_Assemble_getSize(mesh->Nodes,mesh->FaceElements,&tmp);            Finley_Assemble_getSize(mesh->Nodes,mesh->FaceElements,&tmp);
608            break;            break;
609         case(Points):         case(Points):
610            Finley_ErrorCode=TYPE_ERROR;            throw FinleyAdapterException("Error - Size of point elements is not supported.");
           sprintf(Finley_ErrorMsg,"Size of point elements is not supported.");  
611            break;            break;
612         case(ContactElementsZero):         case(ContactElementsZero):
613         case(ContactElementsOne):         case(ContactElementsOne):
614            Finley_Assemble_getSize(mesh->Nodes,mesh->ContactElements,&tmp);            Finley_Assemble_getSize(mesh->Nodes,mesh->ContactElements,&tmp);
615            break;            break;
616         case(DegreesOfFreedom):         case(DegreesOfFreedom):
617            Finley_ErrorCode=TYPE_ERROR;            throw FinleyAdapterException("Error - Size of degrees of freedom is not supported.");
           sprintf(Finley_ErrorMsg,"Size of degrees of freedom is not supported.");  
618            break;            break;
619         case(ReducedDegreesOfFreedom):         case(ReducedDegreesOfFreedom):
620            Finley_ErrorCode=TYPE_ERROR;            throw FinleyAdapterException("Error - Size of reduced degrees of freedom is not supported.");
           sprintf(Finley_ErrorMsg,"Size of reduced degrees of freedom is not supported.");  
621            break;            break;
622         default:         default:
623            Finley_ErrorCode=TYPE_ERROR;            stringstream temp;
624            sprintf(Finley_ErrorMsg,"Element size: Finley does not know anything about function space type %d",size.getFunctionSpace().getTypeCode());            temp << "Error - Element size: Finley does not know anything about function space type " << size.getFunctionSpace().getTypeCode();
625              throw FinleyAdapterException(temp.str());
626            break;            break;
627    }    }
628    checkFinleyError();    checkFinleyError();
629  }  }
630    
631  // sets the location of nodes:  // sets the location of nodes:
632  void MeshAdapter::setNewX(const Data& new_x)  void MeshAdapter::setNewX(const escript::Data& new_x)
633  {  {
634    Finley_Mesh* mesh=m_finleyMesh.get();    Finley_Mesh* mesh=m_finleyMesh.get();
635    const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(new_x.getFunctionSpace().getDomain());    const MeshAdapter& newDomain=dynamic_cast<const MeshAdapter&>(new_x.getFunctionSpace().getDomain());
636    if (newDomain!=*this)    if (newDomain!=*this)
637       throw FinleyAdapterException("Error - Illegal domain of new point locations");       throw FinleyAdapterException("Error - Illegal domain of new point locations");
638    Finley_NodeFile_setCoordinates(mesh->Nodes,&(new_x.getDataC()));    Finley_Mesh_setCoordinates(mesh,&(new_x.getDataC()));
639    checkFinleyError();    checkFinleyError();
640  }  }
641    
642  // saves a data array in openDX format:  // saves a data array in openDX format:
643  void MeshAdapter::saveDX(const std::string& filename,const Data& arg) const  void MeshAdapter::saveDX(const std::string& filename,const boost::python::dict& arg) const
644  {  {
645    Finley_Mesh_saveDX(filename.c_str(),m_finleyMesh.get(),&(arg.getDataC()));      int MAX_namelength=256;
646    checkFinleyError();      const int num_data=boost::python::extract<int>(arg.attr("__len__")());
647        char names[num_data][MAX_namelength];
648        char* c_names[num_data];
649        escriptDataC data[num_data];
650        escriptDataC* ptr_data[num_data];
651    
652        boost::python::list keys=arg.keys();
653        for (int i=0;i<num_data;++i) {
654             escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
655             if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
656                 throw FinleyAdapterException("Error  in saveVTK: Data must be defined on same Domain");
657             data[i]=d.getDataC();
658             ptr_data[i]=&(data[i]);
659             std::string n=boost::python::extract<std::string>(keys[i]);
660             c_names[i]=&(names[i][0]);
661             if (n.length()>MAX_namelength-1) {
662                strncpy(c_names[i],n.c_str(),MAX_namelength-1);
663                c_names[i][MAX_namelength-1]='\0';
664             } else {
665                strcpy(c_names[i],n.c_str());
666             }
667        }
668        Finley_Mesh_saveDX(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
669        checkFinleyError();
670        return;
671    }
672    
673    // saves a data array in openVTK format:
674    void MeshAdapter::saveVTK(const std::string& filename,const boost::python::dict& arg) const
675    {
676        int MAX_namelength=256;
677        const int num_data=boost::python::extract<int>(arg.attr("__len__")());
678        char names[num_data][MAX_namelength];
679        char* c_names[num_data];
680        escriptDataC data[num_data];
681        escriptDataC* ptr_data[num_data];
682    
683        boost::python::list keys=arg.keys();
684        for (int i=0;i<num_data;++i) {
685             escript::Data& d=boost::python::extract<escript::Data&>(arg[keys[i]]);
686             if (dynamic_cast<const MeshAdapter&>(d.getFunctionSpace().getDomain()) !=*this)
687                 throw FinleyAdapterException("Error  in saveVTK: Data must be defined on same Domain");
688             data[i]=d.getDataC();
689             ptr_data[i]=&(data[i]);
690             std::string n=boost::python::extract<std::string>(keys[i]);
691             c_names[i]=&(names[i][0]);
692             if (n.length()>MAX_namelength-1) {
693                strncpy(c_names[i],n.c_str(),MAX_namelength-1);
694                c_names[i][MAX_namelength-1]='\0';
695             } else {
696                strcpy(c_names[i],n.c_str());
697             }
698        }
699        Finley_Mesh_saveVTK(filename.c_str(),m_finleyMesh.get(),num_data,c_names,ptr_data);
700        checkFinleyError();
701        return;
702  }  }
703                                                                                                                                                                            
704                                                                                                                                                                            
705  // creates a SystemMatrixAdapter stiffness matrix an initializes it with zeros:  // creates a SystemMatrixAdapter stiffness matrix an initializes it with zeros:
706  SystemMatrixAdapter MeshAdapter::newSystemMatrix(  SystemMatrixAdapter MeshAdapter::newSystemMatrix(
707                        const int row_blocksize,                        const int row_blocksize,
# Line 781  SystemMatrixAdapter MeshAdapter::newSyst Line 736  SystemMatrixAdapter MeshAdapter::newSyst
736      }      }
737      // generate matrix:      // generate matrix:
738            
739      Finley_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceRowOrder,reduceColOrder);      Paso_SystemMatrixPattern* fsystemMatrixPattern=Finley_getPattern(getFinley_Mesh(),reduceRowOrder,reduceColOrder);
     checkFinleyError();  
     Finley_SystemMatrix* fsystemMatrix=Finley_SystemMatrix_alloc(type,fsystemMatrixPattern,row_blocksize,column_blocksize);  
740      checkFinleyError();      checkFinleyError();
741        Paso_SystemMatrix* fsystemMatrix=Paso_SystemMatrix_alloc(type,fsystemMatrixPattern,row_blocksize,column_blocksize);
742        checkPasoError();
743        Paso_SystemMatrixPattern_dealloc(fsystemMatrixPattern);
744      return SystemMatrixAdapter(fsystemMatrix,row_blocksize,row_functionspace,column_blocksize,column_functionspace);      return SystemMatrixAdapter(fsystemMatrix,row_blocksize,row_functionspace,column_blocksize,column_functionspace);
745  }  }
746    
747  //  //
748  // vtkObject MeshAdapter::createVtkObject() const  // vtkObject MeshAdapter::createVtkObject() const
749  // TODO:  // TODO:
750  //  //
751    
752  //  //
753  // returns true if data at the atom_type is considered as being cell centered:  // returns true if data at the atom_type is considered as being cell centered:
754  bool MeshAdapter::isCellOriented(int functionSpaceCode) const  bool MeshAdapter::isCellOriented(int functionSpaceCode) const
# Line 809  bool MeshAdapter::isCellOriented(int fun Line 767  bool MeshAdapter::isCellOriented(int fun
767            return true;            return true;
768            break;            break;
769         default:         default:
770            Finley_ErrorCode=TYPE_ERROR;            stringstream temp;
771            sprintf(Finley_ErrorMsg,"Cell: Finley does not know anything about function space type %d",functionSpaceCode);            temp << "Error - Cell: Finley does not know anything about function space type " << functionSpaceCode;
772              throw FinleyAdapterException(temp.str());
773            break;            break;
774    }    }
775    checkFinleyError();    checkFinleyError();
776    return false;    return false;
777  }  }
778    
779  bool MeshAdapter::probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const  bool MeshAdapter::probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const
780  {  {
781    switch(functionSpaceType_source) {    switch(functionSpaceType_source) {
# Line 831  bool MeshAdapter::probeInterpolationOnDo Line 791  bool MeshAdapter::probeInterpolationOnDo
791             case(ContactElementsOne):             case(ContactElementsOne):
792                 return true;                 return true;
793             default:             default:
794                 Finley_ErrorCode=TYPE_ERROR;                 stringstream temp;
795                 sprintf(Finley_ErrorMsg,"Interpolation On Domain: Finley does not know anything about function space type %d",functionSpaceType_target);                 temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
796                   throw FinleyAdapterException(temp.str());
797          }          }
798          break;          break;
799       case(Elements):       case(Elements):
# Line 872  bool MeshAdapter::probeInterpolationOnDo Line 833  bool MeshAdapter::probeInterpolationOnDo
833             case(ContactElementsOne):             case(ContactElementsOne):
834                return true;                return true;
835             default:             default:
836               Finley_ErrorCode=TYPE_ERROR;               stringstream temp;
837               sprintf(Finley_ErrorMsg,"Interpolation On Domain: Finley does not know anything about function space type %d",functionSpaceType_target);               temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
838                 throw FinleyAdapterException(temp.str());
839          }          }
840          break;          break;
841       case(ReducedDegreesOfFreedom):       case(ReducedDegreesOfFreedom):
842         switch(functionSpaceType_target) {         switch(functionSpaceType_target) {
843            case(ReducedDegreesOfFreedom):            case(ReducedDegreesOfFreedom):
           case(Nodes):  
844            case(Elements):            case(Elements):
845            case(FaceElements):            case(FaceElements):
846            case(Points):            case(Points):
847            case(ContactElementsZero):            case(ContactElementsZero):
848            case(ContactElementsOne):            case(ContactElementsOne):
849                return true;                return true;
850              case(Nodes):
851            case(DegreesOfFreedom):            case(DegreesOfFreedom):
852               return false;               return false;
853            default:            default:
854               Finley_ErrorCode=TYPE_ERROR;               stringstream temp;
855               sprintf(Finley_ErrorMsg,"Interpolation On Domain: Finley does not know anything about function space type %d",functionSpaceType_target);               temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_target;
856                 throw FinleyAdapterException(temp.str());
857         }         }
858         break;         break;
859       default:       default:
860          Finley_ErrorCode=TYPE_ERROR;          stringstream temp;
861          sprintf(Finley_ErrorMsg,"Interpolation On Domain: Finley does not know anything about function space type %d",functionSpaceType_source);          temp << "Error - Interpolation On Domain: Finley does not know anything about function space type " << functionSpaceType_source;
862            throw FinleyAdapterException(temp.str());
863          break;          break;
864    }    }
865    checkFinleyError();    checkFinleyError();
866    return false;    return false;
867  }  }
868    
869  bool MeshAdapter::probeInterpolationACross(int functionSpaceType_source,const AbstractDomain& targetDomain, int functionSpaceType_target) const  bool MeshAdapter::probeInterpolationACross(int functionSpaceType_source,const AbstractDomain& targetDomain, int functionSpaceType_target) const
870  {  {
871      return false;      return false;
872  }  }
873  bool MeshAdapter::operator==(const MeshAdapter& other) const  
874    bool MeshAdapter::operator==(const AbstractDomain& other) const
875  {  {
876    return (m_finleyMesh==other.m_finleyMesh);    const MeshAdapter* temp=dynamic_cast<const MeshAdapter*>(&other);
877      if (temp!=0) {
878        return (m_finleyMesh==temp->m_finleyMesh);
879      } else {
880        return false;
881      }
882  }  }
883    
884  int MeshAdapter::getSystemMatrixTypeId(const int solver, const bool symmetry) const  bool MeshAdapter::operator!=(const AbstractDomain& other) const
885  {  {
886     int out=Finley_SystemMatrix_getSystemMatrixTypeId(solver,symmetry?1:0);    return !(operator==(other));
887     checkFinleyError();  }
888    
889    int MeshAdapter::getSystemMatrixTypeId(const int solver, const int package, const bool symmetry) const
890    {
891       int out=Paso_SystemMatrix_getSystemMatrixTypeId(SystemMatrixAdapter::mapOptionToPaso(solver),SystemMatrixAdapter::mapOptionToPaso(package),symmetry?1:0);
892       checkPasoError();
893     return out;     return out;
894  }  }
895  Data MeshAdapter::getX() const  
896    escript::Data MeshAdapter::getX() const
897  {  {
898    return continuousFunction(asAbstractContinuousDomain()).getX();    return continuousFunction(asAbstractContinuousDomain()).getX();
899  }  }
900  Data MeshAdapter::getNormal() const  
901    escript::Data MeshAdapter::getNormal() const
902  {  {
903    return functionOnBoundary(asAbstractContinuousDomain()).getNormal();    return functionOnBoundary(asAbstractContinuousDomain()).getNormal();
904  }  }
905  Data MeshAdapter::getSize() const  
906    escript::Data MeshAdapter::getSize() const
907  {  {
908    return function(asAbstractContinuousDomain()).getSize();    return function(asAbstractContinuousDomain()).getSize();
909  }  }
910    
911    int MeshAdapter::getTagFromSampleNo(int functionSpaceType, int sampleNo) const
 bool MeshAdapter::operator!=(const MeshAdapter& other) const  
912  {  {
913    return !operator==(other);    int out=0;
914      Finley_Mesh* mesh=m_finleyMesh.get();
915      switch (functionSpaceType) {
916      case(Nodes):
917        out=mesh->Nodes->Tag[sampleNo];
918        break;
919      case(Elements):
920        out=mesh->Elements->Tag[sampleNo];
921        break;
922      case(FaceElements):
923        out=mesh->FaceElements->Tag[sampleNo];
924        break;
925      case(Points):
926        out=mesh->Points->Tag[sampleNo];
927        break;
928      case(ContactElementsZero):
929        out=mesh->ContactElements->Tag[sampleNo];
930        break;
931      case(ContactElementsOne):
932        out=mesh->ContactElements->Tag[sampleNo];
933        break;
934      case(DegreesOfFreedom):
935        break;
936      case(ReducedDegreesOfFreedom):
937        break;
938      default:
939        stringstream temp;
940        temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
941        throw FinleyAdapterException(temp.str());
942        break;
943      }
944      return out;
945    }
946    int MeshAdapter::getReferenceNoFromSampleNo(int functionSpaceType, int sampleNo) const
947    {
948      int out=0,i;
949      Finley_Mesh* mesh=m_finleyMesh.get();
950      switch (functionSpaceType) {
951      case(Nodes):
952        if (mesh->Nodes!=NULL) {
953          out=mesh->Nodes->Id[sampleNo];
954          break;
955        }
956      case(Elements):
957        out=mesh->Elements->Id[sampleNo];
958        break;
959      case(FaceElements):
960        out=mesh->FaceElements->Id[sampleNo];
961        break;
962      case(Points):
963        out=mesh->Points->Id[sampleNo];
964        break;
965      case(ContactElementsZero):
966        out=mesh->ContactElements->Id[sampleNo];
967        break;
968      case(ContactElementsOne):
969        out=mesh->ContactElements->Id[sampleNo];
970        break;
971      case(DegreesOfFreedom):
972        for (i=0;i<mesh->Nodes->numNodes; ++i) {
973           if (mesh->Nodes->degreeOfFreedom[i]==sampleNo) {
974              out=mesh->Nodes->Id[i];
975              break;
976           }
977        }
978        break;
979      case(ReducedDegreesOfFreedom):
980        for (i=0;i<mesh->Nodes->numNodes; ++i) {
981           if (mesh->Nodes->reducedDegreeOfFreedom[i]==sampleNo) {
982              out=mesh->Nodes->Id[i];
983              break;
984           }
985        }
986        break;
987      default:
988        stringstream temp;
989        temp << "Error - Invalid function space type: " << functionSpaceType << " for domain: " << getDescription();
990        throw FinleyAdapterException(temp.str());
991        break;
992      }
993      return out;
994  }  }
995  // bool MeshAdapter::operator==(const AbstractDomain& other) const  
 // {  
   // try {  
     // const MeshAdapter& ref = dynamic_cast<const MeshAdapter&>(other);  
     // return (operator==(ref));  
   // }  
   // catch (bad_cast) {  
     // return false;  
   // }  
 // }  
996  }  // end of namespace  }  // end of namespace

Legend:
Removed from v.97  
changed lines
  Added in v.547

  ViewVC Help
Powered by ViewVC 1.1.26